This notebook is part of the kikuchipy
documentation https://kikuchipy.org.
Links to the documentation won't work from the notebook.
Note
This tutorial was given at the ESTEEM3 2022 workshop entitled Electron diffraction for solving engineering problems, held at NTNU in Trondheim, Norway in June 2022.
The tutorial has been updated to work with the current release of kikuchipy and to fit our documentation format. The dataset is available from Zenodo at Ånes et al. (2019), and is scan number 7 (17 dB) out of the series of 10 (22 dB) scans taken with increasing gain (0-22 dB).
In this tutorial we will inspect a small (100 MB) EBSD data of polycrystalline recrystallized nickel by (Hough and dictionary) indexing and inspect the results using geometrical EBSD simulations.
Steps:
Tools:
Import libraries (replace inline
with qt5
plotting backend for interactive
plots in separate windows)
%matplotlib inline
from pathlib import Path
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 CrystalMap, Phase, PhaseList
from orix.quaternion import Orientation, symmetry
from orix import io, plot, sampling
from orix.vector import Vector3d
from pyebsdindex import ebsd_index, pcopt
plt.rcParams.update(
{"figure.facecolor": "w", "font.size": 15, "figure.dpi": 75}
)
data_path = Path("../../../kikuchipy_test/esteem")
Load data
s = kp.data.ni_gain(7, allow_download=True) # External download
s
<EBSD, title: Pattern, dimensions: (200, 149|60, 60)>
Plot a pattern from a particular grain to show improvement in signal-to-noise ratio
s.axes_manager.indices = (156, 80)
Plot data
s.plot()
Mean intensity in each pattern
s_mean = s.mean(axis=(2, 3))
s_mean
<BaseSignal, title: Pattern, dimensions: (200, 149|)>
Save map
# plt.imsave(data_path / "maps_mean.png", s_mean.data, cmap="gray")
Set up VBSE image generator
vbse_imager = kp.imaging.VirtualBSEImager(s)
vbse_imager
VirtualBSEImager for <EBSD, title: Pattern, dimensions: (200, 149|60, 60)>
Plot VBSE grid
vbse_imager.plot_grid()
<EBSD, title: Pattern, dimensions: (|60, 60)>
Specify RGB colour channels
r = (2, 1)
b = (2, 2)
g = (2, 3)
Plot coloured grid tiles
vbse_imager.plot_grid(rgb_channels=[r, g, b])
<EBSD, title: Pattern, dimensions: (|60, 60)>
Get VBSE RGB image
vbse_rgb = vbse_imager.get_rgb_image(r, g, b)
s.plot(vbse_rgb)
Save RGB image
# vbse_rgb.save(data_path / "maps_vbse_rgb.png")
Remove static (constant) background
s.remove_static_background()
[########################################] | 100% Completed | 705.70 ms
Inspect statically corrected patterns
s.plot(vbse_rgb)
Remove dynamic (per pattern) background
s.remove_dynamic_background()
[########################################] | 100% Completed | 2.92 ss
Inspect dynamically corrected patterns
s.plot(vbse_rgb)
Average each patterns with the four nearest neighbours
s.average_neighbour_patterns()
[########################################] | 100% Completed | 1.11 sms
Inspect average patterns
s.plot(vbse_rgb)
Save corrected patterns
# s.save(data_path / "patterns_sda.h5")
Get image quality $Q$ map (not image quality IQ from Hough indexing!)
maps_iq = s.get_image_quality()
[########################################] | 100% Completed | 1.51 ss
Navigate patterns in $Q$ map
s_iq = hs.signals.Signal2D(maps_iq)
s.plot(s_iq)
Save $Q$ map to file
# plt.imsave(data_path / "maps_iq.png", maps_iq, cmap="gray")
Get average neighbour dot product map (ADP)
maps_adp = s.get_average_neighbour_dot_product_map()
[########################################] | 100% Completed | 6.65 ss
Navigate patterns in the ADP map
s_adp = hs.signals.Signal2D(maps_adp)
s.plot(s_adp)
Save ADP map to file
# plt.imsave(data_path / "maps_adp.png", maps_adp, cmap="gray")
Orientation of detector with respect to the sample:
Load calibration patterns to get a mean PC for the data
s_cal = kp.data.ni_gain_calibration(7)
s_cal
<EBSD, title: Calibration patterns, dimensions: (9|480, 480)>
s_cal.remove_static_background()
s_cal.remove_dynamic_background()
[########################################] | 100% Completed | 101.02 ms [########################################] | 100% Completed | 102.32 ms
s_cal.plot(navigator="none")
Generate an indexer instance used to optimize the PC and to perform Hough indexing
sig_shape_cal = s_cal.axes_manager.signal_shape[::-1]
indexer_cal = ebsd_index.EBSDIndexer(
phaselist=["FCC"],
vendor="KIKUCHIPY",
sampleTilt=70,
camElev=0,
patDim=sig_shape_cal,
)
Given an initial guess of the PC and optimize for the first calibration pattern,
looking for convergence by updating pc0
manually with the printed PC a couple
of times (the first run of optimize()
takes longer since PyEBSDIndex has to
compile some code before running; consecutive runs are much quicker)
pc0 = [0.42, 0.21, 0.50]
pc = pcopt.optimize(s_cal.inav[0].data, indexer=indexer_cal, PC0=pc0)
print(pc)
[0.40690636 0.22685186 0.50142137]
Optimize for all calibration patterns
pc_all = pcopt.optimize(s_cal.data, indexer_cal, pc, batch=True)
print(pc_all)
[[0.40690636 0.22685186 0.50142137] [0.4080457 0.23103236 0.50566783] [0.42190216 0.23017428 0.50047604] [0.41001409 0.2248648 0.51609248] [0.42017003 0.22805524 0.48807502] [0.42541351 0.21285912 0.50680637] [0.40812411 0.24206145 0.49607673] [0.41075829 0.22226264 0.4940148 ] [0.42612194 0.22455536 0.50103979]]
Calculate the mean PC
pc_mean = pc_all.mean(axis=0)
print(pc_mean)
print(pc_all.std(0))
[0.41527291 0.22696857 0.50107449] [0.00752652 0.00735793 0.00762872]
Index the calibration patterns
plt.figure()
hi_res, *_ = indexer_cal.index_pats(s_cal.data, PC=pc_mean, verbose=2)
Radon Time: 0.032236748000002535 Convolution Time: 0.00548116299998469 Peak ID Time: 0.0021082410000019536 Band Label Time: 0.04232218699999635 Total Band Find Time: 0.08219512799999507 Band Vote Time: 0.008918966000010187
<Figure size 480x360 with 0 Axes>
Plot the returned orientations in the inverse pole figure (IPF), showing the crystal direction $[uvw]$ parallel to the out-of-plane direction (Z)
G_hi = Orientation(hi_res[-1]["quat"], symmetry.Oh)
G_hi.scatter("ipf")
We can determine whether the indexed orientations are correct by projecting geometrical simulations (bands and zone axes) onto the patterns
Load a phase description from a CIF file
phase = Phase.from_cif(str(data_path / "ni.cif"))
phase
<name: ni. space group: Fm-3m. point group: m-3m. proper point group: 432. color: tab:blue>
phase.structure
[Ni 0.000000 0.000000 0.000000 1.0000, Ni 0.000000 0.500000 0.500000 1.0000, Ni 0.500000 0.000000 0.500000 1.0000, Ni 0.500000 0.500000 0.000000 1.0000]
phase.structure.lattice
Lattice(a=3.52387, b=3.52387, c=3.52387, alpha=90, beta=90, gamma=90)
Generate reflectors and filter the list on on minimum structure factor
ref = ReciprocalLatticeVector.from_min_dspacing(phase, 1)
ref.calculate_structure_factor()
F = abs(ref.structure_factor)
ref = ref[F > 0.5 * F.max()]
ref.print_table()
h k l d |F|_hkl |F|^2 |F|^2_rel Mult 1 1 1 2.035 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
Give each distinct set of $\{hkl\}$ a colour
hkl_sets = ref.get_hkl_sets()
hkl_sets
defaultdict(tuple, {(2.0, 0.0, 0.0): (array([ 6, 22, 24, 25, 27, 43]),), (2.0, 2.0, 0.0): (array([ 4, 5, 7, 8, 21, 23, 26, 28, 41, 42, 44, 45]),), (1.0, 1.0, 1.0): (array([12, 13, 16, 17, 32, 33, 36, 37]),), (3.0, 1.0, 1.0): (array([ 0, 1, 2, 3, 9, 10, 11, 14, 15, 18, 19, 20, 29, 30, 31, 34, 35, 38, 39, 40, 46, 47, 48, 49]),)})
hkl_rgb = np.zeros((ref.size, 3))
rgb = [[1, 0, 0], [0, 1, 0], [0, 0, 1], [0.75, 0, 0.75]]
for i, idx in enumerate(hkl_sets.values()):
hkl_rgb[idx] = rgb[i]
Plot each $(hkl)$ in the stereographic projection
ref.scatter(c=hkl_rgb, grid=True, s=100, axes_labels=["e1", "e2"])
Plot the plane trace of each $(hkl)$ in the stereographic projection
ref.draw_circle(color=hkl_rgb, axes_labels=["e1", "e2"])
Set up geometrical simulations
simulator = kp.simulations.KikuchiPatternSimulator(ref)
Calculate Bragg angles and plot bands in the stereographic projection
simulator.reflectors.calculate_theta(20e3)
fig = simulator.plot(hemisphere="upper", mode="bands", return_figure=True)
fig.axes[0].scatter(ref, c=hkl_rgb, s=100)
Specify the detector-sample geometry and perform geometrical simulations on the detector for each orientation found from Hough indexing
det_cal = kp.detectors.EBSDDetector(
shape=sig_shape_cal, pc=pc_mean, sample_tilt=70
)
det_cal
EBSDDetector (480, 480), px_size 1 um, binning 1, tilt 0, azimuthal 0, pc (0.415, 0.227, 0.501)
sim = simulator.on_detector(det_cal, G_hi)
Finding bands that are in some pattern: [########################################] | 100% Completed | 101.39 ms Finding zone axes that are in some pattern: [########################################] | 100% Completed | 103.17 ms Calculating detector coordinates for bands and zone axes: [########################################] | 100% Completed | 101.66 ms
Plot the simulations as markers
sim_markers = sim.as_markers(zone_axes_labels=True)
s_cal.add_marker(sim_markers, permanent=True, plot_signal=False)
s_cal.plot(navigator="none")
# To delete previously added permanent markers, do
del s_cal.metadata.Markers
Create a new indexer and Hough index all patterns
sig_shape = s.axes_manager.signal_shape[::-1]
indexer = ebsd_index.EBSDIndexer(
phaselist=["FCC"],
vendor="KIKUCHIPY",
sampleTilt=det_cal.sample_tilt,
camElev=det_cal.tilt,
patDim=sig_shape,
)
plt.figure()
hi_res_all, *_ = indexer.index_pats(
s.data.reshape((-1,) + sig_shape), PC=det_cal.pc, verbose=2
)
Radon Time: 2.3059811180000622 Convolution Time: 3.065812205000043 Peak ID Time: 2.4583082690002698 Band Label Time: 4.4460913610001 Total Band Find Time: 12.276845957000006 Band Vote Time: 18.648167346999998
<Figure size 480x360 with 0 Axes>
Inspect the first output object (the only one kept)
hi_res_all.dtype
dtype([('quat', '<f8', (4,)), ('iq', '<f4'), ('pq', '<f4'), ('cm', '<f4'), ('phase', '<i4'), ('fit', '<f4'), ('nmatch', '<i4'), ('matchattempts', '<i4', (4,)), ('totvotes', '<i4')])
Create map coordinate arrays
step_size = s.axes_manager["x"].scale
nav_shape = s.axes_manager.navigation_shape[::-1]
i, j = np.indices(nav_shape) * step_size
Contain indexing results in a crystal map for easy plotting and saving
xmap_hi = CrystalMap(
rotations=Orientation(hi_res_all[-1]["quat"]),
phase_list=PhaseList(phase),
x=j.ravel(),
y=i.ravel(),
prop={
"pq": hi_res_all[-1]["pq"],
"cm": hi_res_all[-1]["cm"],
"fit": hi_res_all[-1]["fit"],
},
scan_unit="um",
)
xmap_hi
Phase Orientations Name Space group Point group Proper point group Color 0 29800 (100.0%) ni Fm-3m m-3m 432 tab:blue Properties: pq, cm, fit Scan unit: um
Plot pattern quality (PC) and confidence metric (CM) maps
fig = xmap_hi.plot(
"pq", colorbar=True, colorbar_label="pq", return_figure=True
)
fig.savefig(data_path / "maps_hi_pq.png")
fig = xmap_hi.plot(
"cm", colorbar=True, colorbar_label="cm", return_figure=True
)
fig.savefig(data_path / "maps_hi_cm.png")
Plot (IPF-X) orientation map
ipfkey = plot.IPFColorKeyTSL(
xmap_hi.phases[0].point_group, direction=Vector3d.xvector()
)
fig = ipfkey.plot(return_figure=True)
# fig.savefig(data_path / "ipfkey.png")
G_hi = xmap_hi.orientations
rgb_hi = ipfkey.orientation2color(G_hi)
fig = xmap_hi.plot(rgb_hi, return_figure=True)
# fig.savefig(data_path / "maps_hi_ipfz.png")
fig = xmap_hi.plot(rgb_hi, overlay=maps_iq.ravel(), return_figure=True)
# fig.savefig(data_path / "maps_hi_ipfz_iq.png")
Save Hough indexing results to file (.ang file readable by MTEX, EDAX TSL OIM Analysis etc., HDF5 file can be read back in into Python)
# io.save(data_path / "xmap_hi.ang", xmap_hi)
# io.save(data_path / "xmap_hi.h5", xmap_hi)
Create a new detector with a shape corresponding to the experimental patterns
det = det_cal.deepcopy()
det.shape = sig_shape
Get geometrical simulations on this detector for all orientations from the map
G_hi_2d = G_hi.reshape(*xmap_hi.shape)
sim_hi = simulator.on_detector(det, G_hi_2d)
Finding bands that are in some pattern: [########################################] | 100% Completed | 101.53 ms Finding zone axes that are in some pattern: [########################################] | 100% Completed | 203.94 ms Calculating detector coordinates for bands and zone axes: [########################################] | 100% Completed | 304.15 ms
Plot the first simulated pattern
sim_hi.plot()
rgb_hi_2d = rgb_hi.reshape(xmap_hi.shape + (3,))
s_rgb_hi = kp.draw.get_rgb_navigator(rgb_hi_2d)
Add the geometrical simulations
s.add_marker(sim_hi.as_markers(), permanent=True, plot_signal=False)
s.plot(s_rgb_hi)
# To delete previously added permanent markers, do
del s.metadata.Markers
Improve indexing quality by using dynamical simulations. Start by inspecting the dynamically simulated master pattern in the familiar (?) stereographic projection
mp_sp = kp.data.nickel_ebsd_master_pattern_small(projection="stereographic")
mp_sp.plot()
Plot our geometrical simulation on top of the dynamical simulation
fig = simulator.plot(hemisphere="upper", mode="lines", return_figure=True, color="w")
fig.axes[0].scatter(ref, c=hkl_rgb, s=50)
fig.axes[0].imshow(mp_sp.data, cmap="gray", extent=(-1, 1, -1, 1));
Re-load the master pattern but in the Lambert projection, used to generate the dictionary of dynamically simulated patterns
mp = kp.data.nickel_ebsd_master_pattern_small(projection="lambert")
Quickly inspect the Lambert master pattern
mp.plot()
Discretely sample the complete orientation space of point group $m\bar{3}m$ (Oh) with an average misorientation of about 2.4$^{\circ}$ between points
Gr_sample = sampling.get_sample_fundamental(
resolution=2.4, point_group=mp.phase.point_group
)
G_sample = Orientation(Gr_sample, symmetry=mp.phase.point_group)
G_sample
Orientation (58453,) m-3m [[ 0.8614 -0.3311 -0.3311 -0.1968] [ 0.8614 -0.3349 -0.3349 -0.1836] [ 0.8614 -0.3384 -0.3384 -0.1703] ... [ 0.8614 0.3384 0.3384 0.1703] [ 0.8614 0.3349 0.3349 0.1836] [ 0.8614 0.3311 0.3311 0.1968]]
Plot 1000 sampled orientations in axis-angle space
G_sample.get_random_sample(1000).scatter()
Plot all sampled orientations
G_sample.scatter()
Plot all sampled orientations in the IPFs X, Y, and Z directions
directions = Vector3d(((1, 0, 0), (0, 1, 0), (0, 0, 1)))
G_sample.scatter("ipf", direction=directions, c="C0", s=5)
Set up generation of the dictionary of dynamically simulated patterns (with 2 000 patterns per chunk)
s_dict = mp.get_patterns(G_sample, det, energy=20, chunk_shape=2000)
Generate the five first patterns and plot them
s_dict.inav[:5].plot(navigator="none")
We will use a signal mask so that only pixels with the strongest signal are used during dictionary indexing and refinement
signal_mask = ~kp.filters.Window("circular", det.shape).astype(bool)
p = s.inav[0, 0].data
fig, ax = plt.subplots(ncols=2, figsize=(10, 5))
ax[0].imshow(p * signal_mask, cmap="gray")
ax[0].set_title("Not used in matching")
ax[1].imshow(p * ~signal_mask, cmap="gray")
ax[1].set_title("Used in matching");
Perform dictionary indexing by generating 2 000 simulated patterns at a time and compare them to all the experimental patterns
xmap_di = s.dictionary_indexing(s_dict, signal_mask=signal_mask)
Dictionary indexing information: Phase name: ni Matching 29800 experimental pattern(s) to 58453 dictionary pattern(s) NormalizedCrossCorrelationMetric: float32, greater is better, rechunk: False, navigation mask: False, signal mask: True
100%|███████████████████████████████████████████████████████████████████████████████████████████████████████| 30/30 [02:17<00:00, 4.60s/it]
Indexing speed: 216.13963 patterns/s, 12634009.67826 comparisons/s
xmap_di.scan_unit = "um"
xmap_di
Phase Orientations Name Space group Point group Proper point group Color 0 29800 (100.0%) ni Fm-3m m-3m 432 tab:blue Properties: scores, simulation_indices Scan unit: um
xmap_di.scores.shape
(29800, 20)
Plot similarity scores (normalized cross-correlation, NCC) between best matching experimental and simulated patterns
fig = xmap_di.plot(
xmap_di.scores[:, 0],
colorbar=True,
colorbar_label="NCC",
return_figure=True,
)
# fig.savefig(data_path / "maps_di_ncc.png")
Plot IPF-X orientation map
rgb_di = ipfkey.orientation2color(xmap_di.orientations)
fig = xmap_di.plot(rgb_di, return_figure=True)
# fig.savefig(data_path / "maps_di_ipfz.png")
Save dictionary indexing results to file
# io.save(data_path / "xmap_di.ang", xmap_di)
# io.save(data_path / "xmap_di.h5", xmap_di)
Refine dictionary indexing results by re-generating dynamically simulated patterns iteratively by letting the discretly sampled orientations vary to find the best match
xmap_ref = s.refine_orientation(
xmap=xmap_di,
detector=det,
master_pattern=mp,
energy=20,
signal_mask=signal_mask,
method="LN_NELDERMEAD",
rtol=1e-3,
)
Refinement information: Method: LN_NELDERMEAD (local) from NLopt Trust region (+/-): None Relative tolerance: 0.001 Refining 29800 orientation(s): [########################################] | 100% Completed | 225.72 s Refinement speed: 131.99694 patterns/s
xmap_ref
Phase Orientations Name Space group Point group Proper point group Color 0 29800 (100.0%) ni Fm-3m m-3m 432 tab:blue Properties: scores, num_evals Scan unit: um
Plot refined NCC scores
fig = xmap_ref.plot(
xmap_ref.scores, colorbar=True, colorbar_label="NCC", return_figure=True
)
# fig.savefig(data_path / "maps_di_ref_ncc.png")
Plot refined IPF-X orientation map
rgb_ref = ipfkey.orientation2color(xmap_ref.orientations)
fig = xmap_ref.plot(rgb_ref, return_figure=True)
# fig.savefig(data_path / "maps_di_ref_ipfz.png")
fig = xmap_ref.plot(rgb_ref, return_figure=True, overlay=xmap_ref.scores)
# fig.savefig(data_path / "maps_di_ref_ipfz_ncc.png")
Save refined results to file
# io.save(data_path / "xmap_di_ref.ang", xmap_ref)
# io.save(data_path / "xmap_di_ref.h5", xmap_ref)