Live notebook
You can run this notebook in a live session 
or view it on Github.
Extrapolate projection centers from a mean#
In this tutorial, we will extrapolate a plane of projection centers (PCs) from a mean PC, as determined from patterns spread out across the sample region of interest. This is an alternative to fitting a plane to PCs, as done in the tutorial Fit a plane to selected projection centers. As a verification of the extrapolated PCs, we will compare them to the PCs obtained from fitting a plane to the PCs as is done in that tutorial.
We’ll start by importing the necessary libraries
[1]:
%matplotlib inline
import matplotlib.pyplot as plt
import numpy as np
from diffsims.crystallography import ReciprocalLatticeVector
import hyperspy.api as hs
import kikuchipy as kp
from orix.crystal_map import PhaseList
plt.rcParams.update(
{
"figure.facecolor": "w",
"figure.dpi": 75,
"figure.figsize": (8, 8),
"font.size": 15,
}
)
Load and inspect data#
We will use nine calibration patterns from recrystallized nickel. The patterns have the full (480, 480) px\(^2\) resolution of the NORDIF UF-1100 detector they are acquired on. These patterns are always acquired from spread out sample positions across the region of interest (ROI) in order to calibrate the PCs prior to indexing the full dataset.
Read the calibration patterns
[2]:
s_cal = kp.data.ni_gain_calibration(1, allow_download=True)
s_cal
2023-04-07 12:00:40,602 - numexpr.utils - INFO - NumExpr defaulting to 4 threads.
[2]:
<EBSD, title: Calibration patterns, dimensions: (9|480, 480)>
Get information read with NORDIF calibration patterns
[3]:
omd = s_cal.original_metadata.as_dictionary()
Plot coordinates of calibration patterns on the secondary electron area overview image (part of the dataset), highlighting the region of interest (ROI)
[4]:
kp.draw.plot_pattern_positions_in_map(
rc=omd["calibration_patterns"]["indices_scaled"],
roi_shape=omd["roi"]["shape_scaled"],
roi_origin=omd["roi"]["origin_scaled"],
area_shape=omd["area"]["shape_scaled"],
area_image=omd["area_image"],
color="w",
)
Plot the nine calibration patterns
[5]:
fig = plt.figure(figsize=(10, 10))
_ = hs.plot.plot_images(
s_cal,
per_row=3,
axes_decor=None,
colorbar=False,
label=None,
fig=fig,
)
for i, ax in enumerate(fig.axes):
ax.axis("off")
ax.text(5, 10, str(i), va="top", ha="left", c="w")
fig.subplots_adjust(wspace=0.01, hspace=0.01)
We know that all patterns are of nickel. To get a description of nickel, we could create a Phase manually. However, we will later on use a dynamically simulated EBSD master pattern of nickel (created with EMsoft), which is loaded with a Phase. We will use this in the remaining analysis.
[6]:
mp = kp.data.ebsd_master_pattern(
"ni", allow_download=True, projection="lambert", energy=20
)
mp
[6]:
<EBSDMasterPattern, title: ni_mc_mp_20kv, dimensions: (|1001, 1001)>
Extract the phase, and change the lattice parameter unit from nm to Ångström
[7]:
phase = mp.phase
lat = phase.structure.lattice
lat.setLatPar(lat.a * 10, lat.b * 10, lat.c * 10)
<name: ni. space group: Fm-3m. point group: m-3m. proper point group: 432. color: tab:blue>
lattice=Lattice(a=3.5236, b=3.5236, c=3.5236, alpha=90, beta=90, gamma=90)
28 0.000000 0.000000 0.000000 1.0000
Estimate PCs with Hough indexing#
We will estimate PCs for the nine calibration patterns using PyEBSDIndex. See the Hough indexing tutorial for more details on the Hough indexing related commands below.
Note
PyEBSDIndex is an optional dependency of kikuchipy, and can be installed with both pip and conda (from conda-forge). To install PyEBSDIndex, see their installation instructions.
PyEBSDIndex supports indexing face centered and body centered cubic (FCC and BCC) materials.
Improve signal-to-noise ratio by removing the static and dynamic background
[9]:
s_cal.remove_static_background()
s_cal.remove_dynamic_background()
[########################################] | 100% Completed | 101.73 ms
[########################################] | 100% Completed | 101.63 ms
We need an EBSDIndexer to use PyEBSDIndex. We can obtain an indexer by passing a PhaseList to EBSDDetector.get_indexer(). Therefore, we need an initial EBSD detector
EBSDDetector (480, 480), px_size 1 um, binning 1, tilt 0.0, azimuthal 0.0, pc (0.5, 0.5, 0.5)
70.0
[11]:
phase_list = PhaseList(phase)
phase_list
[11]:
Id Name Space group Point group Proper point group Color
0 ni Fm-3m m-3m 432 tab:blue
[12]:
indexer = det_cal.get_indexer(phase_list)
print(indexer.phaselist)
['FCC']
We estimate the PC of each pattern with Hough indexing, and plot both the mean and standard deviation of the resulting PCs. Note that Bruker’s PC convention is used in kikuchipy. (We will “overwrite” the existing detector variable.)
[13]:
[0.42077662 0.2182365 0.49993497]
[0.0072175 0.00801777 0.00784516]
Plot the PCs
[14]:
det_cal.plot_pc("scatter", annotate=True)
Unfortunately, we do not recognize the spatial distribution from the overview image above, and the expected inverse relation between (PCz, PCy) is not present… We need better indexing of the patterns, which we get by refining the PCs using pattern matching with dynamical simulations. These simulations are created with EMsoft.
First, we need an initial guess for the orientations, using Hough indexing via EBSD.hough_indexing(). We will use the mean PC for all patterns
[15]:
indexer.PC = det_cal.pc_average
[16]:
xmap_hi = s_cal.hough_indexing(
phase_list=phase_list, indexer=indexer, verbose=2
)
Hough indexing with PyEBSDIndex information:
PyOpenCL: True
Projection center (Bruker): (0.4208, 0.2182, 0.4999)
Indexing 9 pattern(s) in 1 chunk(s)
Radon Time: 0.03620907699996678
Convolution Time: 0.0030053609998503816
Peak ID Time: 0.001584530999934941
Band Label Time: 0.04056105099994056
Total Band Find Time: 0.08138581000002887
Band Vote Time: 0.021481833000052575
Indexing speed: 29.15122 patterns/s
Refine PCs with pattern matching#
Refine the PCs (and orientations) using the Nelder-Mead implementation from NLopt
[17]:
Refinement information:
Method: LN_NELDERMEAD (local) from NLopt
Trust region (+/-): [5. 5. 5. 0.05 0.05 0.05]
Relative tolerance: 1e-07
Refining 9 orientation(s) and projection center(s):
[########################################] | 100% Completed | 84.17 ss
Refinement speed: 0.10692 patterns/s
Inspect some refinement statistics
[18]:
print("Score mean: ", xmap_ref.scores.mean())
print("Mean number of evaluations:", xmap_ref.num_evals.mean())
print("PC mean: ", det_ref.pc.mean(0))
print("PC std: ", det_ref.pc.std(0))
print("PC max. diff:", abs(det_cal.pc - det_ref.pc).max(0))
ori_diff = xmap_hi.orientations.angle_with(
xmap_ref.orientations, degrees=True
).max()
print("Ori. max. diff [deg]:", ori_diff)
Score mean: 0.5265933639473386
Mean number of evaluations: 536.7777777777778
PC mean: [0.41998311 0.21323872 0.50126773]
PC std: [0.00269245 0.00188403 0.00060997]
PC max. diff: [0.00819502 0.01585595 0.01511784]
Ori. max. diff [deg]: 1.401081894104408
[19]:
det_ref.plot_pc("scatter", annotate=True)
The diagonals 5-1-0-3-7 and 8-4-0-2-6 seen in the overview image above should be replicated in the (PCx, PCy) scatter plot. We see that 5-1-0-3-7 and 8-0-6 align as expected, but the PC values from the 2nd and 4th patterns do not lie in the expected range. Fitting a plane to all these values might not work to our satisfaction, so we will exclude the 2nd and 4th PC values when fitting a plane to the seven remaining PC values. The plane will have PC values for all points in the ROI in the overview image above.
Get PC plane by fitting#
The fitting is done by finding a transformation function which takes 2D sample coordinates and gives PC values for those coordinates. Both an affine and a projective transformation function is supported, following [Winkelmann et al., 2020]. By passing 2D indices of all ROI map points and of the points where the nine calibration patterns were obtained, EBSDDetector.fit_pc() returns a new detector with PC values for all map points. We will use the maximum difference between the above refined PC values and the corresponding fitted PC values as a measure of how good the fitted PC values are.
[21]:
pc_indices = omd["calibration_patterns"]["indices_scaled"]
pc_indices -= omd["roi"]["origin_scaled"]
pc_indices = pc_indices.T
map_indices = np.indices(omd["roi"]["shape_scaled"])
print("Full map shape (n rows, n columns):", omd["roi"]["shape_scaled"])
Full map shape (n rows, n columns): (149, 200)
Fit PC values using the affine transformation function
[22]:
[0.41979185 0.21355389 0.50146853]
72.41196011494577
Fit PC values using the projective transformation function
[23]:
[0.4197901 0.21355393 0.50146771]
73.03021604645001
Compare PC differences of all but the PC values from the 2nd and 4th patterns
[24]:
pc_indices2 = pc_indices.T[~is_outlier].T
# Refined PC values as a reference (ground truth)
pc_ref = det_ref.pc[~is_outlier]
# Difference in PC values from the affine transformation function
pc_diff_aff = det_fit_aff.pc[tuple(pc_indices2)] - pc_ref
pc_diff_aff_max = abs(pc_diff_aff).max(axis=0)
print(pc_diff_aff_max)
# Difference in PC values from the projective transformation function
pc_diff_proj = det_fit_proj.pc[tuple(pc_indices2)] - pc_ref
pc_diff_proj_max = abs(pc_diff_proj).max(axis=0)
print(pc_diff_proj_max)
# Is the difference from the affine transformation function greater than
# that from the projective one for (PCx, PCy, PCz)?
print(pc_diff_aff_max > pc_diff_proj_max)
[0.00063385 0.00036217 0.0002709 ]
[0.00067374 0.00033036 0.00029519]
[False True False]
Get PC plane by extrapolation#
Instead of fitting a plane to several PCs, we can extrapolate an average PC using EBSDDetector.extrapolate_pc(). To do this we need to know the detector pixel size and map step sizes, both given in the same unit.
[25]:
[0.41887152 0.21397264 0.50132107]
70.0
Difference in PC values between the refined PCs and the corresponding extrapolated PCs
[26]:
[0.00163081 0.00080355 0.0005684 ]
[ True True True]
The extrapolated PCs deviate more from the refined PCs than the PCs obtained from fitting, although the difference is small.
[27]:
det_ext.plot_pc()
Validate extrapolated PCs#
As a final check of the difference in PCs, we can plot the geometrical simulations on top of the patterns using the refined orientation but the three different PCs.
[28]:
rlv = ReciprocalLatticeVector.from_min_dspacing(phase.deepcopy())
rlv.sanitise_phase() # "Fill atoms in unit cell"
rlv.calculate_structure_factor()
structure_factor = abs(rlv.structure_factor)
rlv = rlv[structure_factor > 0.5 * structure_factor.max()]
rlv.print_table()
h k l d |F|_hkl |F|^2 |F|^2_rel Mult
1 1 1 2.034 11.8 140.0 100.0 8
2 0 0 1.762 10.4 108.2 77.3 6
2 2 0 1.246 7.4 55.0 39.3 12
3 1 1 1.062 6.2 38.6 27.6 24
[29]:
We will use a signal mask during refinement to exclude pattern intensities in the detector corners
[30]:
signal_mask = ~kp.filters.Window("circular", det_cal.shape).astype(bool)
Using PCs from the affine transformation function
[31]:
det_aff_cal = det_fit_aff.deepcopy()
det_aff_cal.pc = det_aff_cal.pc[tuple(pc_indices)]
[32]:
Refinement information:
Method: LN_NELDERMEAD (local) from NLopt
Trust region (+/-): [5 5 5]
Relative tolerance: 1e-07
Refining 9 orientation(s):
[########################################] | 100% Completed | 14.47 ss
Refinement speed: 0.62159 patterns/s
[33]:
sim_aff = simulator.on_detector(det_aff_cal, xmap_aff.rotations)
Finding bands that are in some pattern:
[########################################] | 100% Completed | 101.85 ms
Finding zone axes that are in some pattern:
[########################################] | 100% Completed | 102.18 ms
Calculating detector coordinates for bands and zone axes:
[########################################] | 100% Completed | 101.81 ms
Using PCs from the projective transformation function
[34]:
det_proj_cal = det_fit_proj.deepcopy()
det_proj_cal.pc = det_proj_cal.pc[tuple(pc_indices)]
[35]:
Refinement information:
Method: LN_NELDERMEAD (local) from NLopt
Trust region (+/-): [5 5 5]
Relative tolerance: 1e-07
Refining 9 orientation(s):
[########################################] | 100% Completed | 15.32 ss
Refinement speed: 0.58728 patterns/s
[36]:
sim_proj = simulator.on_detector(det_proj_cal, xmap_proj.rotations)
Finding bands that are in some pattern:
[########################################] | 100% Completed | 101.73 ms
Finding zone axes that are in some pattern:
[########################################] | 100% Completed | 104.01 ms
Calculating detector coordinates for bands and zone axes:
[########################################] | 100% Completed | 101.81 ms
Using the extrapolated PCs
[37]:
det_ext_cal = det_ext.deepcopy()
det_ext_cal.pc = det_ext_cal.pc[tuple(pc_indices)]
[38]:
Refinement information:
Method: LN_NELDERMEAD (local) from NLopt
Trust region (+/-): [5 5 5]
Relative tolerance: 1e-07
Refining 9 orientation(s):
[########################################] | 100% Completed | 14.67 ss
Refinement speed: 0.61334 patterns/s
[39]:
sim_ext = simulator.on_detector(det_ext_cal, xmap_ext.rotations)
Finding bands that are in some pattern:
[########################################] | 100% Completed | 102.02 ms
Finding zone axes that are in some pattern:
[########################################] | 100% Completed | 101.76 ms
Calculating detector coordinates for bands and zone axes:
[########################################] | 100% Completed | 102.76 ms
Compare normalized cross-correlation scores and number of evaluations (iterations)
[40]:
print("Scores\n------")
print(f"Affine: {xmap_aff.scores.mean():.7f}")
print(f"Projective: {xmap_proj.scores.mean():.7f}")
print(f"Extrapolated: {xmap_ext.scores.mean():.7f}\n")
print("Number of evaluations\n---------------------")
print(f"Affine: {xmap_aff.num_evals.mean():.1f}")
print(f"Projective: {xmap_proj.num_evals.mean():.1f}")
print(f"Extrapolated: {xmap_ext.num_evals.mean():.1f}")
Scores
------
Affine: 0.5692966
Projective: 0.5692996
Extrapolated: 0.5685377
Number of evaluations
---------------------
Affine: 168.4
Projective: 176.7
Extrapolated: 169.3
Plot Kikuchi bands on top of patterns for the solutions using the affine transformed PCs (red), projective transformed PCs (blue) and extrapolated PCs (white)
[41]:
fig, axes = plt.subplots(ncols=3, nrows=3, figsize=(12, 12))
for i, ax in enumerate(axes.ravel()):
ax.imshow(s_cal.inav[i].data, cmap="gray")
ax.axis("off")
# Affine
lines = sim_aff.as_collections(
i, lines_kwargs=dict(linewidth=4, alpha=0.4)
)[0]
ax.add_collection(lines)
# Projective
lines = sim_proj.as_collections(
i, lines_kwargs=dict(color="b", linewidth=3, alpha=0.4)
)[0]
ax.add_collection(lines)
# Extrapolated
lines = sim_ext.as_collections(
i, lines_kwargs=dict(color="w", linewidth=1, alpha=0.4)
)[0]
ax.add_collection(lines)
ax.text(5, 10, i, c="w", va="top", ha="left", fontsize=20)
fig.subplots_adjust(wspace=0.01, hspace=0.01)
As expected from the intermediate results above (similar PCs and mean NCC score), all PCs produce visually identical geometrical simulations. However, the orientations are still slightly different, as we see by calculating the misorientation angles between the orientations obtained using the extrapolated PCs and the fitted PCs
[42]:
xmap_proj.orientations.angle_with(xmap_ext.orientations, degrees=True)
[42]:
array([3.00016517, 2.99419537, 2.99722527, 3.00370345, 3.00205733,
2.99066348, 2.99660149, 3.00870014, 3.01203182])
The misorientation of 3\(^{\circ}\) most likely comes from the difference in applied sample tilt between the detector used for the projected and extrapolated PCs. Since these are experimental data, it’s difficult to say which sample tilt is more correct, although the nominal sample tilt from the microscope was 70\(^{\circ}\).