Show code cell source

Hide code cell source

import ase
import matplotlib.pyplot as plt
from ase.io import read

import abtem

SAED quickstart

This notebook demonstrates a basic simulation of selected area electron diffraction of SiO₂ zeolite.

Configuration

We start by (optionally) setting our configuration. See documentation for details.

abtem.config.set(
    {
        "device": "cpu",
        "fft": "numpy",
        "diagnostics.task_progress": True,
        "diagnostics.progress_bar": "tqdm",
    }
);

Atomic model

We import a unit cell from a .cif-file. See our walkthough or our tutorial on atomic models.

sio2 = read("data/SiO2.cif")

abtem.show_atoms(sio2, plane="xy", legend=True);

The structure is rotated such that the previous \(xy\)-plane is rotated into the \(xz\)-plane. The structure is repeated along \(z\) to obtain the required thickness. If frozen phonons are required the structure should also be repeated in the \(xy\) to describe a thermal ensemble.

sio2_rotated = abtem.atoms.rotate_atoms_to_plane(sio2, "xz")

sio2_repeated = sio2_rotated * (3, 2, 30)

fig, (ax1, ax2) = plt.subplots(1, 2, figsize=(12, 4))
abtem.show_atoms(sio2_repeated, ax=ax1, plane="xy", title="Beam view")
abtem.show_atoms(sio2_repeated, ax=ax2, plane="yz", title="Side view", linewidth=0);

Potential

We create an ensemble of potentials using the frozen phonon model. See our walkthrough on frozen phonons.

frozen_phonons = abtem.FrozenPhonons(sio2_repeated, 8, sigmas=0.1)

We create a potential from the frozen phonons model, see walkthrough on potentials.

potential = abtem.Potential(
    frozen_phonons,
    gpts=512,
    projection="infinite",
    slice_thickness=2,
    exit_planes=5,
)

len(potential)

300

Wave function

We create a plane wave function at an energy of \(200 \ \mathrm{keV}\). See our walkthrough on wave functions.

wave = abtem.PlaneWave(energy=200e3)

Multislice

We run the multislice algorithm to calculate the exit waves, see our walkthrough on multislice.

exit_waves = wave.multislice(potential, lazy=True)

We calculate the diffraction patterns up to a scattering angle of \(12 \ \mathrm{mrad}\).

measurement_ensemble = exit_waves.diffraction_patterns(max_angle=12)

measurement_ensemble.shape

(8, 61, 41, 41)

The result is an ensemble of images, one for each frozen phonon and exit plane, we average the ensemble across the phonon dimensions obtain the final thickness series.

measurement = measurement_ensemble.mean(0)

measurement.compute();

<abtem.measurements.DiffractionPatterns at 0x1c6a91e50>

Visualize results

We show the thickness series as an exploded plot. Before plotting the direct beam is blocked as it typically has a much higher intensity than the scattered beams.

visualization = (
    measurement[::15]
    .block_direct()
    .show(
        explode=True,
        figsize=(12, 5),
        cbar=True,
        common_color_scale=True,
    )
)

Index diffraction patterns

Instead of the pixelated representation above, the diffraction patterns can be represented as a collection indexed diffraction intensities.

# the conventional unit cell is given here
indexed_spots = measurement.index_diffraction_spots(cell=sio2.cell)

We show the last diffraction spots in the thickness series with an overlay showing the miller indices.

indexed_spots[-1].block_direct().show(
    scale=.01,
    figsize=(5, 5),
    annotation_kwargs={"threshold": 1e-2},
);

The diffraction spots can be converted to a pandas dataframe.

import pandas as pd

pd.set_option("display.max_rows", 8)

df = indexed_spots.block_direct().remove_low_intensity(1e-2).to_dataframe()

df

hkl	0 2 0	0 -2 0	1 1 0	1 5 0	1 -5 0	1 -1 0	2 0 0	2 2 0	2 8 0	2 -8 0	...	-3 -3 0	-2 0 0	-2 2 0	-2 8 0	-2 -8 0	-2 -2 0	-1 1 0	-1 5 0	-1 -5 0	-1 -1 0
z [Å]
0.000000	0.000000	0.000000	0.000000	0.000000	0.000000	0.000000	0.000000	0.000000	0.000000	0.000000	...	0.000000	0.000000	0.000000	0.000000	0.000000	0.000000	0.000000	0.000000	0.000000	0.000000
9.997273	0.000096	0.000098	0.000089	0.000280	0.000279	0.000094	0.000040	0.000022	0.000005	0.000005	...	0.000164	0.000039	0.000024	0.000006	0.000005	0.000022	0.000089	0.000282	0.000279	0.000094
19.994545	0.000383	0.000382	0.000360	0.001093	0.001107	0.000361	0.000161	0.000087	0.000023	0.000025	...	0.000653	0.000161	0.000091	0.000026	0.000024	0.000088	0.000359	0.001095	0.001108	0.000360
29.991818	0.000851	0.000872	0.000785	0.002477	0.002464	0.000830	0.000363	0.000200	0.000068	0.000063	...	0.001486	0.000360	0.000203	0.000072	0.000063	0.000199	0.000786	0.002482	0.002463	0.000826
...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...
569.844534	0.027980	0.028138	0.028099	0.001694	0.001758	0.028080	0.021415	0.011171	0.017014	0.016410	...	0.002691	0.021716	0.011047	0.016673	0.017149	0.010764	0.028736	0.001664	0.001691	0.028447
579.841807	0.026394	0.026363	0.027709	0.001028	0.001176	0.027259	0.021265	0.010764	0.017964	0.018010	...	0.001587	0.021627	0.010662	0.017654	0.018640	0.010368	0.028314	0.001000	0.001141	0.027687
589.839079	0.024416	0.024519	0.026684	0.000880	0.000822	0.026589	0.020916	0.010285	0.019563	0.019063	...	0.000771	0.021311	0.010238	0.019313	0.019520	0.009997	0.027276	0.000856	0.000823	0.026940
599.836352	0.022392	0.022462	0.025663	0.000970	0.000978	0.025312	0.020371	0.009851	0.020544	0.020768	...	0.000249	0.020885	0.009815	0.020315	0.021055	0.009521	0.026240	0.000963	0.001024	0.025660

61 rows × 36 columns

df[["4 0 0", "3 3 0", "1 5 0", "1 1 0"]].plot()

<Axes: xlabel='z [Å]'>