This page was generated from doc/tutorials/esteem2022_diffraction_workshop.ipynb. Interactive online version: Binder badge.

ESTEEM3 workshop#

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:

  1. Data overview (virtual backscatter electron imaging)

  2. Enhance Kikuchi pattern (background correction)

  3. Data overview (“feature maps”)

  4. Hough indexing

  5. Verification of results (geometrical simulations)

  6. Dictionary indexing

  7. Orientation refinement

  8. Summary

Tools:

Import libraries (replace inline with qt5 plotting backend for interactive plots in separate windows)

[1]:
%matplotlib inline

import os

from diffsims.crystallography import ReciprocalLatticeVector
import hyperspy.api as hs
import kikuchipy as kp
import matplotlib.pyplot as plt
import numpy as np
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})

Load data

[2]:
dir_data = "../../../kikuchipy_test/esteem"
s = kp.load(os.path.join(dir_data, "ebsd_data/Pattern.dat"))
[3]:
s
[3]:
<EBSD, title: Pattern, dimensions: (200, 149|60, 60)>

Plot a pattern from a particular grain to show improvement in signal-to-noise ratio

[4]:
s.axes_manager.indices = (156, 80)

Plot data

[5]:
s.plot()
../_images/tutorials_esteem2022_diffraction_workshop_10_0.png
../_images/tutorials_esteem2022_diffraction_workshop_10_1.png

Data overview (virtual backscatter electron [VBSE] imaging)#

Mean intensity in each pattern

[6]:
s_mean = s.mean(axis=(2, 3))
s_mean
[6]:
<BaseSignal, title: Pattern, dimensions: (200, 149|)>

Save map

[7]:
# plt.imsave(os.path.join(dir_data, "maps_mean.png"), s_mean.data, cmap="gray")

Set up VBSE generator

[8]:
VirtualBSEGenerator for <EBSD, title: Pattern, dimensions: (200, 149|60, 60)>

Plot VBSE grid

[9]:
vbse_gen.plot_grid()
[9]:
<EBSD, title: Pattern, dimensions: (|60, 60)>
../_images/tutorials_esteem2022_diffraction_workshop_19_1.png

Specify RGB colour channels

[10]:
r = (2, 1)
b = (2, 2)
g = (2, 3)

Plot coloured grid tiles

[11]:
vbse_gen.plot_grid(rgb_channels=[r, g, b])
[11]:
<EBSD, title: Pattern, dimensions: (|60, 60)>
../_images/tutorials_esteem2022_diffraction_workshop_23_1.png

Get VBSE RGB image

[12]:
vbse_rgb = vbse_gen.get_rgb_image(r, g, b)
[13]:
s.plot(vbse_rgb)
../_images/tutorials_esteem2022_diffraction_workshop_26_0.png
../_images/tutorials_esteem2022_diffraction_workshop_26_1.png

Save RGB image

[14]:
# vbse_rgb.save(os.path.join(dir_data, "maps_vbse_rgb.png"))

Enhance Kikuchi pattern (background correction)#

Remove static (constant) background

[15]:
s.remove_static_background()
[########################################] | 100% Completed | 809.69 ms

Inspect statically corrected patterns

[16]:
s.plot(vbse_rgb)
../_images/tutorials_esteem2022_diffraction_workshop_33_0.png
../_images/tutorials_esteem2022_diffraction_workshop_33_1.png

Remove dynamic (per pattern) background

[17]:
s.remove_dynamic_background()
[########################################] | 100% Completed | 3.02 ss

Inspect dynamically corrected patterns

[18]:
s.plot(vbse_rgb)
../_images/tutorials_esteem2022_diffraction_workshop_37_0.png
../_images/tutorials_esteem2022_diffraction_workshop_37_1.png

Average each patterns with the four nearest neighbours

[19]:
s.average_neighbour_patterns()
[########################################] | 100% Completed | 1.12 sms

Inspect average patterns

[20]:
s.plot(vbse_rgb)
../_images/tutorials_esteem2022_diffraction_workshop_41_0.png
../_images/tutorials_esteem2022_diffraction_workshop_41_1.png

Save corrected patterns

[21]:
# s.save(os.path.join(dir_data, "patterns_sda.h5"))

Data overview (“feature maps”)#

Get image quality \(Q\) map (not image quality IQ from Hough indexing!)

[22]:
iq = s.get_image_quality()
[########################################] | 100% Completed | 1.61 ss

Navigate patterns in \(Q\) map

[24]:
s.plot(s_iq)
../_images/tutorials_esteem2022_diffraction_workshop_49_0.png
../_images/tutorials_esteem2022_diffraction_workshop_49_1.png

Save \(Q\) map to file

[25]:
# plt.imsave(os.path.join(dir_data, "maps_iq.png"), s_iq.data, cmap="gray")

Get average neighbour dot product map (ADP)

[26]:
adp = s.get_average_neighbour_dot_product_map()
[########################################] | 100% Completed | 7.37 ss

Navigate patterns in the ADP map

[27]:
[28]:
s.plot(s_adp)
../_images/tutorials_esteem2022_diffraction_workshop_56_0.png
../_images/tutorials_esteem2022_diffraction_workshop_56_1.png

Save ADP map to file

[29]:
# plt.imsave(os.path.join(dir_data, "maps_adp.png"), s_adp.data, cmap="gray")

Hough indexing (of calibration patterns)#

Orientation of detector with respect to the sample: * Known: * Sample tilt (about microscope X) * Camera tilt (about microscope X) * Unknown: * Projection/pattern centre (PCx, PCy, PCz): Shortest distance from source point to detector

Load calibration patterns to get a mean PC for the data

[30]:
s_cal = kp.load(os.path.join(dir_data, "ebsd_data/Setting.txt"))
[31]:
s_cal
[31]:
<EBSD, title: Calibration patterns, dimensions: (9|480, 480)>
[32]:
s_cal.remove_static_background()
s_cal.remove_dynamic_background()
[########################################] | 100% Completed | 102.44 ms
[########################################] | 100% Completed | 102.48 ms
[33]:
s_cal.plot(navigator="none")
../_images/tutorials_esteem2022_diffraction_workshop_63_0.png

Generate an indexer instance used to optimize the PC and to perform Hough indexing

[34]:
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)

[35]:
pc0 = [0.42, 0.21, 0.50]
pc = pcopt.optimize(s_cal.inav[0].data, indexer=indexer_cal, PC0=pc0)
print(pc)
[0.41857685 0.2271845  0.50193261]

Optimize for all calibration patterns

[36]:
pc_all = np.zeros((9, 3))
for i in range(len(s_cal)):
    pc_all[i] = pcopt.optimize(s_cal.inav[i].data, indexer=indexer_cal, PC0=pc)

# Replace the above with the following after PyEBSDIndex > 0.1.0 is released
# pc_all = pcopt.optimize(s_cal.data, indexer_cal, pc, batch=True)

print(pc_all)
[[0.42033693 0.22717537 0.50521349]
 [0.42671022 0.21948888 0.49898992]
 [0.42756565 0.23393217 0.50839931]
 [0.41205772 0.21906086 0.50599244]
 [0.42011539 0.21376857 0.49933781]
 [0.44451362 0.23659495 0.49572229]
 [0.41587398 0.22465204 0.49432407]
 [0.41274765 0.23215346 0.50806298]
 [0.42889034 0.23029534 0.49482328]]

Calculate the mean PC

[37]:
pc_mean = pc_all.mean(axis=0)
pc_mean
[37]:
array([0.42320128, 0.22634685, 0.50120729])

Index the calibration patterns

[38]:
hi_res = indexer_cal.index_pats(s_cal.data, PC=pc_mean, verbose=2)[0]
Radon Time: 0.18021900999883655
Convolution Time: 0.0047501619992544875
Peak ID Time: 0.004734218993689865
Band Label Time: 0.00044612400233745575
Total Band Find Time: 0.19017342998995446
../_images/tutorials_esteem2022_diffraction_workshop_73_1.png
Band Vote Time:  0.023265301002538763

Plot the returned orientations in the inverse pole figure (IPF), showing the crystal direction \([uvw]\) parallel to the out-of-plane direction (Z)

[39]:
g = Orientation(hi_res[-1]["quat"], symmetry.Oh)
g.scatter("ipf")
../_images/tutorials_esteem2022_diffraction_workshop_75_0.png

We can determine whether the indexed orientations are correct by projecting geometrical simulations (bands and zone axes) onto the patterns

Geometrical simulations#

Load a phase description from a CIF file

[40]:
phase = Phase.from_cif(os.path.join(dir_data, "ni.cif"))
[41]:
[41]:
<name: ni. space group: Fm-3m. point group: m-3m. proper point group: 432. color: tab:blue>
[42]:
phase.structure
[42]:
[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]
[43]:
phase.structure.lattice
[43]:
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

[44]:
ref = ReciprocalLatticeVector.from_min_dspacing(phase, 1)

ref.calculate_structure_factor()

f_real = abs(ref.structure_factor)
ref = ref[f_real > 0.5 * f_real.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

[45]:
hkl_sets = ref.get_hkl_sets()
[46]:
hkl_sets
[46]:
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]),)})
[47]:
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

[48]:
ref.scatter(c=hkl_rgb, grid=True, s=100, axes_labels=["e1", "e2"])
../_images/tutorials_esteem2022_diffraction_workshop_89_0.png

Plot the plane trace of each \((hkl)\) in the stereographic projection

[49]:
ref.draw_circle(color=hkl_rgb, axes_labels=["e1", "e2"])
../_images/tutorials_esteem2022_diffraction_workshop_91_0.png

Set up geometrical simulations

Calculate Bragg angles and plot bands in the stereographic projection

[51]:
simulator.reflectors.calculate_theta(20e3)
[52]:
fig = simulator.plot(hemisphere="upper", mode="bands", return_figure=True)
fig.axes[0].scatter(ref, c=hkl_rgb, s=100)
../_images/tutorials_esteem2022_diffraction_workshop_96_0.png

Specify the detector-sample geometry and perform geometrical simulations on the detector for each orientation found from Hough indexing

[53]:
detector_cal = kp.detectors.EBSDDetector(
    shape=sig_shape_cal, pc=pc_mean, sample_tilt=70
)
detector_cal
[53]:
EBSDDetector (480, 480), px_size 1 um, binning 1, tilt 0, azimuthal 0, pc (0.423, 0.226, 0.501)
[54]:
sim = simulator.on_detector(detector_cal, g)
Finding bands that are in some pattern:
[########################################] | 100% Completed | 101.29 ms
Finding zone axes that are in some pattern:
[########################################] | 100% Completed | 103.05 ms
Calculating detector coordinates for bands and zone axes:
[########################################] | 100% Completed | 102.26 ms

Plot the simulations as markers

[55]:
sim_markers = sim.as_markers(zone_axes_labels=True)
[56]:
s_cal.add_marker(sim_markers, permanent=True, plot_signal=False)
[57]:
s_cal.plot(navigator="none")
../_images/tutorials_esteem2022_diffraction_workshop_103_0.png
[58]:
# To delete previously added permanent markers, do
del s_cal.metadata.Markers

Hough indexing of all patterns#

Create a new indexer and Hough index all patterns

[59]:
sig_shape = s.axes_manager.signal_shape[::-1]
indexer = ebsd_index.EBSDIndexer(
    phaselist=["FCC"],
    vendor="KIKUCHIPY",
    sampleTilt=70,
    camElev=0,
    patDim=sig_shape,
)
[60]:
hi_res_all = indexer.index_pats(
    s.data.reshape((-1,) + sig_shape), PC=detector_cal.pc, verbose=2
)[0]
Radon Time: 75.71765887404035
Convolution Time: 23.454073408967815
Peak ID Time: 17.955805922014406
Band Label Time: 11.42529403300432
Total Band Find Time: 128.55397094899672
../_images/tutorials_esteem2022_diffraction_workshop_107_1.png
Band Vote Time:  17.4025655919977

Inspect the first output object (the only one kept)

[61]:
hi_res_all.dtype
[61]:
dtype([('quat', '<f8', (4,)), ('iq', '<f4'), ('pq', '<f4'), ('cm', '<f4'), ('phase', '<i4'), ('fit', '<f4'), ('nmatch', '<i4'), ('matchattempts', '<i4', (2,)), ('totvotes', '<i4')])

Create map coordinate arrays

[62]:
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

[63]:
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",
)
[64]:
[64]:
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

[65]:
fig = xmap_hi.plot("pq", colorbar=True, colorbar_label="pq", return_figure=True)
fig.savefig(os.path.join(dir_data, "maps_hi_pq.png"))
../_images/tutorials_esteem2022_diffraction_workshop_116_0.png
[66]:
fig = xmap_hi.plot("cm", colorbar=True, colorbar_label="cm", return_figure=True)
fig.savefig(os.path.join(dir_data, "maps_hi_cm.png"))
../_images/tutorials_esteem2022_diffraction_workshop_117_0.png

Plot (IPF-Z) orientation map

[67]:
ipfkey = plot.IPFColorKeyTSL(xmap_hi.phases[0].point_group)
[68]:
fig = ipfkey.plot(return_figure=True)
../_images/tutorials_esteem2022_diffraction_workshop_120_0.png
[69]:
# fig.savefig(os.path.join(dir_data, "ipfkey.png"))
[70]:
ori_hi = xmap_hi.orientations
[71]:
ori_rgb = ipfkey.orientation2color(ori_hi)
[72]:
fig = xmap_hi.plot(ori_rgb, return_figure=True)
# fig.savefig(os.path.join(dir_data, "maps_hi_ipfz.png"))
../_images/tutorials_esteem2022_diffraction_workshop_124_0.png
[73]:
fig = xmap_hi.plot(ori_rgb, overlay=iq.ravel(), return_figure=True)
# fig.savefig(os.path.join(dir_data, "maps_hi_ipfz_iq.png"))
../_images/tutorials_esteem2022_diffraction_workshop_125_0.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)

[74]:
# io.save(os.path.join(dir_data, "xmap_hi.ang"), xmap_hi)
[75]:
# io.save(os.path.join(dir_data, "xmap_hi.h5"), xmap_hi)

Create a new detector with a shape corresponding to the experimental patterns

[76]:
detector = detector_cal.deepcopy()
detector.shape = sig_shape

Get geometrical simulations on this detector for all orientations from the map

[77]:
ori_hi_2d = ori_hi.reshape(*xmap_hi.shape)
sim_hi = simulator.on_detector(detector, ori_hi_2d)
Finding bands that are in some pattern:
[########################################] | 100% Completed | 101.56 ms
Finding zone axes that are in some pattern:
[########################################] | 100% Completed | 202.87 ms
Calculating detector coordinates for bands and zone axes:
[########################################] | 100% Completed | 309.11 ms

Plot the first simulated pattern

[78]:
sim_hi.plot()
../_images/tutorials_esteem2022_diffraction_workshop_134_0.png
[79]:
ori_rgb_2d = ori_rgb.reshape(xmap_hi.shape + (3,))
s_rgb_hi = kp.draw.get_rgb_navigator(ori_rgb_2d)

Add the geometrical simulations

[80]:
s.add_marker(sim_hi.as_markers(), permanent=True, plot_signal=False)
[81]:
s.plot(s_rgb_hi)
../_images/tutorials_esteem2022_diffraction_workshop_138_0.png
../_images/tutorials_esteem2022_diffraction_workshop_138_1.png
[82]:
# To delete previously added permanent markers, do
del s.metadata.Markers

Dictionary indexing#

Improve indexing quality by using dynamical simulations. Start by inspecting the dynamically simulated master pattern in the familiar (?) stereographic projection

[83]:
mp_sp = kp.data.nickel_ebsd_master_pattern_small(projection="stereographic")
[84]:
../_images/tutorials_esteem2022_diffraction_workshop_142_0.png

Plot our geometrical simulation on top of the dynamical simulation

[85]:
fig = simulator.plot(hemisphere="upper", mode="lines", return_figure=True)
fig.axes[0].scatter(ref, c=hkl_rgb, s=50)

fig.axes[0].imshow(mp_sp.data, cmap="gray", extent=(-1, 1, -1, 1));
../_images/tutorials_esteem2022_diffraction_workshop_144_0.png

Re-load the master pattern but in the Lambert projection, used to generate the dictionary of dynamically simulated patterns

[86]:

Quickly inspect the Lambert master pattern

[87]:
../_images/tutorials_esteem2022_diffraction_workshop_148_0.png

Discretely sample the complete orientation space of point group \(m\bar{3}m\) (Oh) with an average misorientation of 2.5\(^{\circ}\) between points

[88]:
g_sample = sampling.get_sample_fundamental(
    resolution=2.5, point_group=mp.phase.point_group
)
g_sample = Orientation(g_sample, symmetry=mp.phase.point_group)
[89]:
[89]:
Orientation (52607,) m-3m
[[ 0.8642 -0.3165 -0.3165 -0.2297]
 [ 0.8642 -0.3211 -0.3211 -0.2167]
 [ 0.8642 -0.3254 -0.3254 -0.2034]
 ...
 [ 0.8642  0.3254  0.3254  0.2034]
 [ 0.8642  0.3211  0.3211  0.2167]
 [ 0.8642  0.3165  0.3165  0.2297]]

Plot 1000 sampled orientations in axis-angle space

[90]:
../_images/tutorials_esteem2022_diffraction_workshop_153_0.png

Plot all sampled orientations

[91]:
g_sample.scatter()
../_images/tutorials_esteem2022_diffraction_workshop_155_0.png

Plot all sampled orientations in the IPFs X, Y and Z directions

[92]:
directions = Vector3d(((1, 0, 0), (0, 1, 0), (0, 0, 1)))
g_sample.scatter("ipf", direction=directions, c="C0", s=5)
../_images/tutorials_esteem2022_diffraction_workshop_157_0.png

Set up generation of the dictionary of dynamically simulated patterns

[93]:
s_dict = mp.get_patterns(g_sample, detector, energy=20)

Generate the five first patterns and plot them

[94]:
s_dict.inav[:5].plot(navigator="none")
../_images/tutorials_esteem2022_diffraction_workshop_161_0.png

Perform dictionary indexing by generating 2000 simulated patterns at a time and compare them to all the experimental patterns

[95]:
xmap_di = s.dictionary_indexing(s_dict, n_per_iteration=2000)
Dictionary indexing information:
        Phase name: ni/ni
        Matching 29800 experimental pattern(s) to 52607 dictionary pattern(s)
        NormalizedCrossCorrelationMetric: float32, greater is better, rechunk: False, signal mask: False
100%|████████████████████████████████████████████████████████████████████| 27/27 [02:59<00:00,  6.66s/it]
        Indexing speed: 165 patterns/s, 8717798 comparisons/s
[96]:
xmap_di.scan_unit = "um"
[97]:
xmap_di
[97]:
Phase    Orientations   Name  Space group  Point group  Proper point group     Color
    0  29800 (100.0%)  ni/ni        Fm-3m         m-3m                 432  tab:blue
Properties: scores, simulation_indices
Scan unit: um
[98]:
xmap_di.scores.shape
[98]:
(29800, 20)

Plot similarity scores (normalized cross-correlation, NCC) between best matching experimental and simulated patterns

[99]:
fig = xmap_di.plot(
    xmap_di.scores[:, 0], colorbar=True, colorbar_label="NCC", return_figure=True
)
# fig.savefig(os.path.join(dir_data, "maps_di_ncc.png"))
../_images/tutorials_esteem2022_diffraction_workshop_168_0.png

Plot IPF-Z orientation map

[100]:
ipfkey = plot.IPFColorKeyTSL(xmap_di.phases[0].point_group)
g_rgb_di = ipfkey.orientation2color(xmap_di.orientations)
[101]:
fig = xmap_di.plot(g_rgb_di, return_figure=True)
# fig.savefig(os.path.join(dir_data, "maps_di_ipfz.png"))
../_images/tutorials_esteem2022_diffraction_workshop_171_0.png

Save dictionary indexing results to file

[102]:
# io.save(os.path.join(dir_data, "xmap_di.ang"), xmap_di)
# io.save(os.path.join(dir_data, "xmap_di.h5"), xmap_di)

Orientation refinement#

Refine dictionary indexing results by re-generating dynamically simulated patterns iteratively by letting the discretly sampled orientations vary within +/- 1.5\(^{\circ}\) and finding the best match

[103]:
xmap_ref = s.refine_orientation(
    xmap=xmap_di,
    detector=detector,
    master_pattern=mp,
    energy=20,
    trust_region=[1.5, 1.5, 1.5],
    method_kwargs=dict(options=dict(fatol=0.001)),
)
Refinement information:
        Local optimization method: Nelder-Mead (minimize)
        Keyword arguments passed to method: {'options': {'fatol': 0.001}, 'method': 'Nelder-Mead'}
        Trust region: [1.5, 1.5, 1.5]
Refining 29800 orientation(s):
[########################################] | 100% Completed | 468.64 s
Refinement speed: 63 patterns/s
[104]:
xmap_ref
[104]:
Phase    Orientations   Name  Space group  Point group  Proper point group     Color
    0  29800 (100.0%)  ni/ni        Fm-3m         m-3m                 432  tab:blue
Properties: scores
Scan unit: um

Plot refined NCC scores

[105]:
fig = xmap_ref.plot(
    xmap_ref.scores, colorbar=True, colorbar_label="NCC", return_figure=True
)
# fig.savefig(os.path.join(dir_data, "maps_di_ref_ncc.png"))
../_images/tutorials_esteem2022_diffraction_workshop_178_0.png

Plot refined IPF-Z orientation map

[106]:
g_rgb_ref = ipfkey.orientation2color(xmap_ref.orientations)
[107]:
fig = xmap_ref.plot(g_rgb_ref, return_figure=True)
# fig.savefig(os.path.join(dir_data, "maps_di_ref_ipfz.png"))
../_images/tutorials_esteem2022_diffraction_workshop_181_0.png
[108]:
fig = xmap_ref.plot(g_rgb_ref, return_figure=True, overlay=xmap_ref.scores)
# fig.savefig(os.path.join(dir_data, "maps_di_ref_ipfz_ncc.png"))
../_images/tutorials_esteem2022_diffraction_workshop_182_0.png

Save refined results to file

[109]:
# io.save(os.path.join(dir_data, "xmap_di_ref.ang"), xmap_ref)
# io.save(os.path.join(dir_data, "xmap_di_ref.h5"), xmap_ref)