process

calibration

Functions related to elliptical calibration, such as fitting elliptical distortions.

The user-facing representation of ellipses is in terms of the following 5 :param x0: :param y0 the center of the ellipse: :param a the semimajor axis length: :param b the semiminor axis length: :param theta the: to the x-axis, in radians :type theta the: positive, right handed

More details about the elliptical parameterization used can be found in the module docstring for process/utils/elliptical_coords.py.

py4DSTEM.process.calibration.ellipse.fit_ellipse_1D(ar, center=None, fitradii=None, mask=None)

For a 2d array ar, fits a 1d elliptical curve to the data inside an annulus centered at center with inner and outer radii at fitradii. The data to fit make optionally be additionally masked with the boolean array mask. See module docstring for more info.

Parameters:
  • ar (ndarray) – array containing the data to fit

  • center (2-tuple of floats) – the center (x0,y0) of the annular fitting region

  • fitradii (2-tuple of floats) – inner and outer radii (ri,ro) of the fit region

  • mask (ar-shaped ndarray of bools) – ignore data wherever mask==True

Returns:

A 5-tuple containing the ellipse parameters:
  • x0: the center x-position

  • y0: the center y-position

  • a: the semimajor axis length

  • b: the semiminor axis length

  • theta: the tilt of the ellipse semimajor axis with respect to the x-axis, in radians

Return type:

(5-tuple of floats)

py4DSTEM.process.calibration.ellipse.ellipse_err(p, x, y, val)

For a point (x,y) in a 2d cartesian space, and a function taking the value val at point (x,y), and some 1d ellipse in this space given by

A(x-x0)^2 + B(x-x0)(y-y0) + C(y-y0)^2 = 1

this function computes the error associated with the function’s value at (x,y) given by its deviation from the ellipse times val.

Note that this function is for internal use, and uses ellipse parameters p given in canonical form (x0,y0,A,B,C), which is different from the ellipse parameterization used in all the user-facing functions, for reasons of numerical stability.

py4DSTEM.process.calibration.ellipse.fit_ellipse_amorphous_ring(data, center, fitradii, p0=None, mask=None)

Fit the amorphous halo of a diffraction pattern, including any elliptical distortion.

The fit function is:

f(x,y; I0,I1,sigma0,sigma1,sigma2,c_bkgd,x0,y0,A,B,C) =
    Norm(r; I0,sigma0,0) +
    Norm(r; I1,sigma1,R)*Theta(r-R)
    Norm(r; I1,sigma2,R)*Theta(R-r) + c_bkgd

where

  • (x,y) are cartesian coordinates,

  • r is the radial coordinate,

  • (I0,I1,sigma0,sigma1,sigma2,c_bkgd,x0,y0,R,B,C) are parameters,

  • Norm(x;I,s,u) is a gaussian in the variable x with maximum amplitude I, standard deviation s, and mean u

  • Theta(x) is a Heavyside step function

  • R is the radial center of the double sided gaussian, derived from (A,B,C) and set to the mean of the semiaxis lengths

The function thus contains a pair of gaussian-shaped peaks along the radial direction of a polar-elliptical parametrization of a 2D plane. The first gaussian is centered at the origin. The second gaussian is centered about some finite R, and is ‘two-faced’: it’s comprised of two half-gaussians of different standard deviations, stitched together at their mean value of R. This Janus (two-faced ;p) gaussian thus comprises an elliptical ring with different inner and outer widths.

The parameters of the fit function are

  • I0: the intensity of the first gaussian function

  • I1: the intensity of the Janus gaussian

  • sigma0: std of first gaussian

  • sigma1: inner std of Janus gaussian

  • sigma2: outer std of Janus gaussian

  • c_bkgd: a constant offset

  • x0,y0: the origin

  • A,B,C: The ellipse parameters, in the form Ax^2 + Bxy + Cy^2 = 1

Parameters:
  • data (2d array) – the data

  • center (2-tuple of numbers) – the center (x0,y0)

  • fitradii (2-tuple of numbers) – the inner and outer radii of the fitting annulus

  • p0 (11-tuple) – initial guess parameters. If p0 is None, the function will compute a guess at all parameters. If p0 is a 11-tuple it must be populated by some mix of numbers and None; any parameters which are set to None will be guessed by the function. The parameters are the 11 parameters of the fit function described above, p0 = (I0,I1,sigma0,sigma1,sigma2,c_bkgd,x0,y0,A,B,C). Note that x0,y0 are redundant; their guess values are the x0,y0 values passed to the main function, but if they are passed as elements of p0 these will take precendence.

  • mask (2d array of bools) – only fit to datapoints where mask is True

Returns:

Returns a 2-tuple.

The first element is the ellipse parameters need to elliptically parametrize diffraction space, and is itself a 5-tuple:

  • x0: x center

  • y0: y center,

  • a: the semimajor axis length

  • b: the semiminor axis length

  • theta: tilt of a-axis w.r.t x-axis, in radians

The second element is the full set of fit parameters to the double sided gaussian function, described above, and is an 11-tuple

Return type:

(2-tuple comprised of a 5-tuple and an 11-tuple)

py4DSTEM.process.calibration.ellipse.double_sided_gaussian_fiterr(p, x, y, val)

Returns the fit error associated with a point (x,y) with value val, given parameters p.

py4DSTEM.process.calibration.ellipse.double_sided_gaussian(p, x, y)

Return the value of the double-sided gaussian function at point (x,y) given parameters p, described in detail in the fit_ellipse_amorphous_ring docstring.

py4DSTEM.process.calibration.ellipse.constrain_degenerate_ellipse(data, p_ellipse, r_inner, r_outer, phi_known, fitrad=6)

When fitting an ellipse to data containing 4 diffraction spots in a narrow annulus about the central beam, the answer is degenerate: an infinite number of ellipses correctly fit this data. Starting from one ellipse in the degenerate family of ellipses, this function selects the ellipse which will yield a final angle of phi_known between a pair of the diffraction peaks after performing elliptical distortion correction.

Note that there are two possible angles which phi_known might refer to, because the angle of interest is well defined up to a complementary angle. This function is written such that phi_known should be the smaller of these two angles.

Parameters:
  • data (ndarray) –

  • p_ellipse (5-tuple) – the ellipse parameters (x0,y0,a,b,theta)

  • r_inner (float) – the fitting annulus inner radius

  • r_outer (float) – the fitting annulus outer radius

  • phi_known (float) – the known angle between a pair of diffraction peaks, in radians

  • fitrad (float) – the region about the fixed data point used to refine its position

Returns:

A 2-tuple containing:

  • a_constrained: (float) the first semiaxis of the selected ellipse

  • b_constrained: (float) the second semiaxis of the selected ellipse

Return type:

(2-tuple)

py4DSTEM.process.calibration.origin.fit_origin(data, mask=None, fitfunction='plane', returnfitp=False, robust=False, robust_steps=3, robust_thresh=2)

Fits the position of the origin of diffraction space to a plane or parabola, given some 2D arrays (qx0_meas,qy0_meas) of measured center positions, optionally masked by the Boolean array mask. The 2D data arrays may be passed directly as a 2-tuple to the arg data, or, if data is either a DataCube or Calibration instance, they will be retreived automatically. If a DataCube or Calibration are passed, fitted origin and residuals are stored there directly.

Parameters:
  • data (2-tuple of 2d arrays) – the measured origin position (qx0,qy0)

  • mask (2b boolean array, optional) – ignore points where mask=False

  • fitfunction (str, optional) – must be ‘plane’ or ‘parabola’ or ‘bezier_two’ or ‘constant’

  • returnfitp (bool, optional) – if True, returns the fit parameters

  • robust (bool, optional) – If set to True, fit will be repeated with outliers removed.

  • robust_steps (int, optional) – Optional parameter. Number of robust iterations performed after initial fit.

  • robust_thresh (int, optional) – Threshold for including points, in units of root-mean-square (standard deviations) error of the predicted values after fitting.

Returns:

Return value depends on returnfitp. If returnfitp==False (default), returns a 4-tuple containing:

  • qx0_fit: (ndarray) the fit origin x-position

  • qy0_fit: (ndarray) the fit origin y-position

  • qx0_residuals: (ndarray) the x-position fit residuals

  • qy0_residuals: (ndarray) the y-position fit residuals

If returnfitp==True, returns a 2-tuple. The first element is the 4-tuple described above. The second element is a 4-tuple (popt_x,popt_y,pcov_x,pcov_y) giving fit parameters and covariance matrices with respect to the chosen fitting function.

Return type:

(variable)

py4DSTEM.process.calibration.origin.get_origin_single_dp(dp, r, rscale=1.2)

Find the origin for a single diffraction pattern, assuming (a) there is no beam stop, and (b) the center beam contains the highest intensity.

Parameters:
  • dp (ndarray) – the diffraction pattern

  • r (number) – the approximate disk radius

  • rscale (number) – factor by which r is scaled to generate a mask

Returns:

The origin

Return type:

(2-tuple)

py4DSTEM.process.calibration.origin.get_origin(datacube, r=None, rscale=1.2, dp_max=None, mask=None, fast_center=False)

Find the origin for all diffraction patterns in a datacube, assuming (a) there is no beam stop, and (b) the center beam contains the highest intensity. Stores the origin positions in the Calibration associated with datacube, and optionally also returns them.

Parameters:
  • datacube (DataCube) – the data

  • r (number or None) – the approximate radius of the center disk. If None (default), tries to compute r using the get_probe_size method. The data used for this is controlled by dp_max.

  • rscale (number) – expand ‘r’ by this amount to form a mask about the center disk when taking its center of mass

  • dp_max (ndarray or None) –

    the diffraction pattern or dp-shaped array used to compute the center disk radius, if r is left unspecified. Behavior depends on type:

    • if dp_max==None (default), computes and uses the maximal diffraction pattern. Note that for a large datacube, this may be a slow operation.

    • otherwise, this should be a (Q_Nx,Q_Ny) shaped array

  • mask (ndarray or None) – if not None, should be an (R_Nx,R_Ny) shaped boolean array. Origin is found only where mask==True, and masked arrays are returned for qx0,qy0

  • fast_center – (bool) Skip the center of mass refinement step.

Returns:

the origin, (x,y) at each scan position

Return type:

(2-tuple of (R_Nx,R_Ny)-shaped ndarrays)

py4DSTEM.process.calibration.origin.get_origin_single_dp_beamstop(DP: ndarray, mask: ndarray, **kwargs)

Find the origin for a single diffraction pattern, assuming there is a beam stop.

Parameters:
  • DP (np array) – diffraction pattern

  • mask (np array) – boolean mask which is False under the beamstop and True in the diffraction pattern. One approach to generating this mask is to apply a suitable threshold on the average diffraction pattern and use binary opening/closing to remove and holes

Returns:

qx0, qy0 (tuple) measured center position of diffraction pattern

py4DSTEM.process.calibration.origin.get_origin_beamstop(datacube: DataCube, mask: ndarray, **kwargs)

Find the origin for each diffraction pattern, assuming there is a beam stop.

Parameters:
  • datacube (DataCube) –

  • mask (np array) – boolean mask which is False under the beamstop and True in the diffraction pattern. One approach to generating this mask is to apply a suitable threshold on the average diffraction pattern and use binary opening/closing to remove any holes

Returns:

qx0, qy0 (tuple of np arrays) measured center position of each diffraction pattern

py4DSTEM.process.calibration.probe.get_probe_size(DP, thresh_lower=0.01, thresh_upper=0.99, N=100)

Gets the center and radius of the probe in the diffraction plane.

The algorithm is as follows: First, create a series of N binary masks, by thresholding the diffraction pattern DP with a linspace of N thresholds from thresh_lower to thresh_upper, measured relative to the maximum intensity in DP. Using the area of each binary mask, calculate the radius r of a circular probe. Because the central disk is typically very intense relative to the rest of the DP, r should change very little over a wide range of intermediate values of the threshold. The range in which r is trustworthy is found by taking the derivative of r(thresh) and finding identifying where it is small. The radius is taken to be the mean of these r values. Using the threshold corresponding to this r, a mask is created and the CoM of the DP times this mask it taken. This is taken to be the origin x0,y0.

Parameters:
  • DP (2D array) – the diffraction pattern in which to find the central disk. A position averaged, or shift-corrected and averaged, DP works best.

  • thresh_lower (float, 0 to 1) – the lower limit of threshold values

  • thresh_upper (float, 0 to 1) – the upper limit of threshold values

  • N (int) – the number of thresholds / masks to use

Returns:

A 3-tuple containing:

  • r: (float) the central disk radius, in pixels

  • x0: (float) the x position of the central disk center

  • y0: (float) the y position of the central disk center

Return type:

(3-tuple)

py4DSTEM.process.calibration.qpixelsize.get_Q_pixel_size(q_meas, q_known, units='A')

Computes the size of the Q-space pixels.

Parameters:
  • q_meas (number) – a measured distance in q-space in pixels

  • q_known (number) – the corresponding known real space distance

  • unit (str) – the units of the real space value of q_known

Returns:

the detector pixel size, the associated units

Return type:

(number,str)

py4DSTEM.process.calibration.qpixelsize.get_dq_from_indexed_peaks(qs, hkl, a)

Get dq, the size of the detector pixels in the diffraction plane, in inverse length units, using a set of measured peak distances from the optic axis, their Miller indices, and the known unit cell size.

Parameters:
  • qs (array) – the measured peak positions

  • hkl (list/tuple of length-3 tuples) – the Miller indices of the peak positions qs. The length of qs and hkl must be the same. To ignore any peaks, for this peak set (h,k,l)=(0,0,0).

  • a (number) – the unit cell size

Returns:

A 4-tuple containing:

  • dq: (number) the detector pixel size

  • qs_fit: (array) the fit positions of the peaks

  • hkl_fit: (list/tuple of length-3 tuples) the Miller indices of the fit peaks

  • mask: (array of bools) False wherever hkl[i]==(0,0,0)

Return type:

(4-tuple)

py4DSTEM.process.calibration.rotation.compare_QR_rotation(im_R, im_Q, QR_rotation, R_rotation=0, R_position=None, Q_position=None, R_pos_anchor='center', Q_pos_anchor='center', R_length=0.33, Q_length=0.33, R_width=0.001, Q_width=0.001, R_head_length_adjust=1, Q_head_length_adjust=1, R_head_width_adjust=1, Q_head_width_adjust=1, R_color='r', Q_color='r', figsize=(10, 5), returnfig=False)

Visualize a rotational offset between an image in real space, e.g. a STEM virtual image, and an image in diffraction space, e.g. a defocused CBED shadow image of the same region, by displaying an arrow overlaid over each of these two images with the specified QR rotation applied. The QR rotation is defined as the counter-clockwise rotation from real space to diffraction space, in degrees.

Parameters:
  • im_R (numpy array or other 2D image-like object (e.g. a VirtualImage)) – A real space image, e.g. a STEM virtual image

  • im_Q (numpy array or other 2D image-like object) – A diffraction space image, e.g. a defocused CBED image

  • QR_rotation (number) – The counterclockwise rotation from real space to diffraction space, in degrees

  • R_rotation (number) – The orientation of the arrow drawn in real space, in degrees

  • R_position (None or 2-tuple) – The position of the anchor point for the R-space arrow. If None, defaults to the center of the image

  • Q_position (None or 2-tuple) – The position of the anchor point for the Q-space arrow. If None, defaults to the center of the image

  • R_pos_anchor ('center' or 'tail' or 'head') – The anchor point for the R-space arrow, i.e. the point being specified by the R_position parameter

  • Q_pos_anchor ('center' or 'tail' or 'head') – The anchor point for the Q-space arrow, i.e. the point being specified by the Q_position parameter

  • R_length (number or None) – The length of the R-space arrow, as a fraction of the mean size of the image

  • Q_length (number or None) – The length of the Q-space arrow, as a fraction of the mean size of the image

  • R_width (number) – The width of the R-space arrow

  • Q_width (number) – The width of the R-space arrow

  • R_head_length_adjust (number) – Scaling factor for the R-space arrow head length

  • Q_head_length_adjust (number) – Scaling factor for the Q-space arrow head length

  • R_head_width_adjust (number) – Scaling factor for the R-space arrow head width

  • Q_head_width_adjust (number) – Scaling factor for the Q-space arrow head width

  • R_color (color) – Color of the R-space arrow

  • Q_color (color) – Color of the Q-space arrow

  • figsize (2-tuple) – The figure size

  • returnfig (bool) – Toggles returning the figure and axes

py4DSTEM.process.calibration.rotation.get_Qvector_from_Rvector(vx, vy, QR_rotation)

For some vector (vx,vy) in real space, and some rotation QR between real and reciprocal space, determine the corresponding orientation in diffraction space. Returns both R and Q vectors, normalized.

Parameters:
  • vx (numbers) – the (x,y) components of a real space vector

  • vy (numbers) – the (x,y) components of a real space vector

  • QR_rotation (number) – the offset angle between real and reciprocal space.

  • Specifically

  • to (the counterclockwise rotation of real space with respect) –

  • degrees. (diffraction space. In) –

Returns:

4-tuple consisting of:

  • vx_R: the x component of the normalized real space vector

  • vy_R: the y component of the normalized real space vector

  • vx_Q: the x component of the normalized reciprocal space vector

  • vy_Q: the y component of the normalized reciprocal space vector

Return type:

(4-tuple)

py4DSTEM.process.calibration.rotation.get_Rvector_from_Qvector(vx, vy, QR_rotation)

For some vector (vx,vy) in diffraction space, and some rotation QR between real and reciprocal space, determine the corresponding orientation in diffraction space. Returns both R and Q vectors, normalized.

Parameters:
  • vx (numbers) – the (x,y) components of a reciprocal space vector

  • vy (numbers) – the (x,y) components of a reciprocal space vector

  • QR_rotation (number) – the offset angle between real and reciprocal space. Specifically, the counterclockwise rotation of real space with respect to diffraction space. In degrees.

Returns:

4-tuple consisting of:

  • vx_R: the x component of the normalized real space vector

  • vy_R: the y component of the normalized real space vector

  • vx_Q: the x component of the normalized reciprocal space vector

  • vy_Q: the y component of the normalized reciprocal space vector

Return type:

(4-tuple)

classification

class py4DSTEM.process.classification.braggvectorclassification.BraggVectorClassification(braggpeaks, Qx, Qy, X_is_boolean=True, max_dist=None)

A class for classifying 4D-STEM data based on which Bragg peaks are found at each diffraction pattern.

A BraggVectorClassification instance enables classification using several methods; a brief overview is provided here, with more details in each individual method’s documentation.

Initialization methods:

__init__:

Determine the initial classes. The approach here involves first segmenting diffraction space, using maxima of a Bragg vector map.

get_initial_classes_by_cooccurrence:

Class refinement methods: Each of these methods creates a new set of candidate classes, but does not yet overwrite the old classes. This enables the new classes to be viewed and compared to the old classes before deciding whether to accept or reject them. Thus running two of these methods in succession, without accepting changes in between, simply discards the first set of candidate classes.

nmf:

Nonnegative matrix factorization (X = WH) to refine the classes. Briefly, after constructing a matrix X which describes which Bragg peaks were observed in each diffraction pattern, we factor X into two smaller matrices, W and H. Physically, W and H describe a small set of classes, each of which corresponds to some subset of (or, more strictly, weights for) the Bragg peaks and the scan positions. We additionally impose the contraint that, on physical grounds, all the elements of X, W, and H must be nonnegative.

split:

If any classes contain multiple non-contiguous segments in real space, divide these into distinct classes.

merge:

If any classes contain sufficient overlap in both scan positions and BPs, merge them into a single class.

Accepting/rejecting changes:

accept:

Updates classes (the W and H matrices) with the current candidate classes.

reject:

Discard the current candidate classes.

Class examination methods:

get_class:

get a single class, returning both its BP weights and scan position weights

get_class_BPs:

get the BP weights for a single class

get_class_image:

get the image, i.e. scan position weights, associated with a single class

get_candidate_class:

as above, for the current candidate class

get_candidate_class_BPs:

as above, for the current candidate class

get_candidate_class_image:

as above, for the current candidate class

Parameters:
  • braggpeaks (PointListArray) – Bragg peaks; must have coords ‘qx’ and ‘qy’

  • Qx (ndarray of floats) – x-coords of the voronoi points

  • Qy (ndarray of floats) – y-coords of the voronoi points

  • X_is_boolean (bool) – if True, populate X with bools (BP is or is not present). if False, populate X with floats (BP c.c. intensities)

  • max_dist (None or number) – maximum distance from a given voronoi point a peak can be and still be associated with this label

__init__(braggpeaks, Qx, Qy, X_is_boolean=True, max_dist=None)

Initializes a BraggVectorClassification instance.

This method: 1. Gets integer labels for all of the detected Bragg peaks, according to which

(Qx,Qy) is closest, then generating a corresponding set of integers for each scan position. See get_braggpeak_labels_by_scan_position() docstring for more info.

  1. Generates the data matrix X. See the nmf() method docstring for more info.

This method should be followed by one of the methods which populates the initial classes - currently, either get_initial_classes_by_cooccurrence() or get_initial_classes_from_images. These methods generate the W and H matrices – i.e. the decompositions of the X matrix in terms of scan positions and Bragg peaks – which are necessary for any subsequent processing.

Parameters:
  • braggpeaks (PointListArray) – Bragg peaks; must have coords ‘qx’ and ‘qy’

  • Qx (ndarray of floats) – x-coords of the voronoi points

  • Qy (ndarray of floats) – y-coords of the voronoi points

  • X_is_boolean (bool) – if True, populate X with bools (BP is or is not present). if False, populate X with floats (BP c.c. intensities)

  • max_dist (None or number) – maximum distance from a given voronoi point a peak can be and still be associated with this label

R_Nx

shape of real space (x)

R_Ny

shape of real space (y)

Qx

x-coordinates of the voronoi points

Qy

y-coordinates of the voronoi points

braggpeak_labels

the sets of Bragg peaks present at each scan position

N_feat

first dimension of the data matrix; the number of bragg peaks

N_meas

second dimension of the data matrix; the number of scan positions

X

the data matrix

get_initial_classes_by_cooccurrence(thresh=0.3, BP_fraction_thresh=0.1, max_iterations=200, X_is_boolean=True, n_corr_init=2)

Populate the initial classes by finding sets of Bragg peaks that tend to co-occur in the same diffraction patterns.

Beginning from the sets of Bragg peaks labels for each scan position (determined in __init__), this method gets initial classes by determining which labels are most likely to co-occur with each other – see get_initial_classes() docstring for more info. Then the matrices W and H are generated – see nmf() doscstring for discussion.

Parameters:
  • thresh (float in [0,1]) – threshold for adding new BPs to a class

  • BP_fraction_thresh (float in [0,1]) – algorithm terminates if fewer than this fraction of the BPs have not been assigned to a class

  • max_iterations (int) – algorithm terminates after this many iterations

  • n_corr_init (int) – seed new classes by finding maxima of the n-point joint probability function. Must be 2 or 3.

get_initial_classes_from_images(class_images)

Populate the initial classes using a set of user-defined class images.

Parameters:

class_images (ndarray) – must have shape (R_Nx,R_Ny,N_c), where N_c is the number of classes, and class_images[:,:,i] is the image of class i.

nmf(max_iterations=1)

Nonnegative matrix factorization to refine the classes.

The data matrix X is factored into two smaller matrices, W and H:

X = WH

Here,

  • ``X``is the data matrix. It has shape (N_feat,N_meas), where N_feat is the number of Bragg peak integer labels (i.e. len(Qx)) and N_meas is the number of diffraction patterns (i.e. R_Nx*R_Ny). Element X[i,j] represents the value of the i’th BP in the j’th DP. The values depend on the flag datamatrix_is_boolean: if True, X[i,j] is 1 if this BP was present in this DP, or 0 if not; if False, X[i,j] is the cross correlation intensity of this BP in this DP.

  • W is the class matrix. It has shape (N_feat,N_c), where N_c is the number of classes. The i’th column vector, w_i = W[:,i], describes the weight of each Bragg peak in the i’th class. w_i has length N_feat, and w_i[j] describes how strongly the j’th BP is associated with the i’th class.

  • H is the coefficient matrix. It has shape (N_c,N_meas). The i’th column vector H[:,i] describes the contribution of each class to scan position i.

Alternatively, we can completely equivalently think of H as a class matrix, and W as a coeffient matrix. In this picture, the i’th row vector of H, h_i = H[i,:], describes the weight of each scan position in the i’th class. h_i has length N_meas, and h_i[j] describes how strongly the j’th scan position is associated with the i’th class. The row vector W[i,:] is then a coefficient vector, which gives the contributions each of the (H) classes to the measured values of the i’th BP. These pictures are related by a transpose: X = WH is equivalent to X.T = (H.T)(W.T).

In nonnegative matrix factorization we impose the constrain, here on physical grounds, that all elements of X, W, and H should be nonnegative.

The computation itself is performed using the sklearn nmf class. When this method is called, the three relevant matrices should already be defined. This method refines W and H, with up to max_iterations NMF steps.

Parameters:

max_iterations (int) – the maximum number of NMF steps to take

split(sigma=2, threshold_split=0.25, expand_mask=1, minimum_pixels=1)

If any classes contain multiple non-contiguous segments in real space, divide these regions into distinct classes.

Algorithm is as follows: First, an image of each class is obtained from its scan position weights. Then, the image is convolved with a gaussian of std sigma. This is then turned into a binary mask, by thresholding with threshold_split. Stray pixels are eliminated by performing a one pixel binary closing, then binary opening. The mask is then expanded by expand_mask pixels. Finally, the contiguous regions of the resulting mask are found. These become the new class components by scan position.

The splitting itself involves creating two classes - i.e. adding a column to W and a row to H. The new BP classes (W columns) have exactly the same values as the old BP class. The two new scan position classes (H rows) divide up the non-zero entries of the old scan position class into two or more non-intersecting subsets, each of which becomes its own new class.

Parameters:
  • sigma (float) – std of gaussian kernel used to smooth the class images before thresholding and splitting.

  • threshold_split (float) – used to threshold the class image to create a binary mask.

  • expand_mask (int) – number of pixels by which to expand the mask before separating into contiguous regions.

  • minimum_pixels (int) – if, after splitting, a potential new class contains fewer than this number of pixels, ignore it

merge(threshBPs=0.1, threshScanPosition=0.1, return_params=True)

If any classes contain sufficient overlap in both scan positions and BPs, merge them into a single class.

The algorithm is as follows: First, the Pearson correlation coefficient matrix is calculated for the classes according to both their diffraction space, Bragg peak representations (i.e. the correlations of the columns of W) and according to their real space, scan position representations (i.e. the correlations of the rows of H). Class pairs whose BP correlation coefficient exceeds threshBPs and whose scan position correlation coefficient exceed threshScanPosition are deemed ‘sufficiently overlapped’, and are marked as merge candidates. To account for intransitivity issues (e.g. class pairs 1/2 and 2/3 are merge candidates, but class pair 1/3 is not), merging is then performed beginning with candidate pairs with the greatest product of the two correlation coefficients, skipping later merge candidate pairs if one of the two classes has already been merged.

The algorithm can be looped until no more merge candidates satisfying the specified thresholds remain with the merge_iterative method.

The merging itself involves turning two classes into one by combining a pair of W columns (i.e. the BP representations of the classes) and the corresponding pair of H rows (i.e. the scan position representation of the class) into a single W column / H row. In terms of scan positions, the new row of H is generated by simply adding the two old H rows. In terms of Bragg peaks, the new column of W is generated by adding the two old columns of W, while weighting each by its total intensity in real space (i.e. the sum of its H row).

Parameters:
  • threshBPs (float) – the threshold for the bragg peaks correlation coefficient, above which the two classes are considered candidates for merging

  • threshScanPosition (float) – the threshold for the scan position correlation coefficient, above which two classes are considered candidates for merging

  • return_params (bool) – if True, returns W_corr, H_corr, and merge_candidates. Otherwise, returns nothing. Incompatible with iterative=True.

merge_by_class_index(i, j)

Merge classes i and j into a single class.

Columns i and j of W pair of W (i.e. the BP representations of the classes) and the corresponding pair of H rows (i.e. the scan position representation of the class) are mergedinto a single W column / H row. In terms of scan positions, the new row of H is generated by simply adding the two old H rows. In terms of Bragg peaks, the new column of W is generated by adding the two old columns of W, while weighting each by its total intensity in real space (i.e. the sum of its H row).

Parameters:
  • i (int) – index of the first class to merge

  • j (int) – index of the second class to merge

split_by_class_index(i, sigma=2, threshold_split=0.25, expand_mask=1, minimum_pixels=1)

If class i contains multiple non-contiguous segments in real space, divide these regions into distinct classes.

Algorithm is as described in the docstring for self.split.

Parameters:
  • i (int) – index of the class to split

  • sigma (float) – std of gaussian kernel used to smooth the class images before thresholding and splitting.

  • threshold_split (float) – used to threshold the class image to create a binary mask.

  • expand_mask (int) – number of pixels by which to expand the mask before separating into contiguous regions.

  • minimum_pixels (int) – if, after splitting, a potential new class contains fewer than this number of pixels, ignore it

remove_class(i)

Remove class i.

Parameters:

i (int) – index of the class to remove

merge_iterative(threshBPs=0.1, threshScanPosition=0.1)

If any classes contain sufficient overlap in both scan positions and BPs, merge them into a single class.

Identical to the merge method, with the addition of iterating until no new merge pairs are found.

Parameters:
  • threshBPs (float) – the threshold for the bragg peaks correlation coefficient, above which the two classes are considered candidates for merging

  • threshScanPosition (float) – the threshold for the scan position correlation coefficient, above which two classes are considered candidates for merging

accept()

Updates classes (the W and H matrices) with the current candidate classes.

reject()

Discard the current candidate classes.

get_class(i)

Get a single class, returning both its BP weights and scan position weights.

Parameters:

i (int) – the class index

Returns:

A 2-tuple containing:

  • class_BPs: (length N_feat array of floats) the weights of the N_feat Bragg peaks for this class

  • class_image: (shape (R_Nx,R_Ny) array of floats) the weights of each scan position in this class

Return type:

(2-tuple)

get_class_BPs(i)

Get a single class, returning its BP weights.

Parameters:

i (int) – the class index

Returns:

the weights of the N_feat Bragg peaks for this class

Return type:

(length N_feat array of floats)

get_class_image(i)

Get a single class, returning its scan position weights.

Parameters:

i (int) – the class index

Returns:

the weights of each scan position in this class

Return type:

(shape (R_Nx,R_Ny) array of floats)

get_candidate_class(i)

Get a single candidate class, returning both its BP weights and scan position weights.

Parameters:

i (int) –

Returns:

A 2-tuple containing:

  • class_BPs: (length N_feat array of floats) the weights of the N_feat Bragg peaks for this class

  • class_image: (shape (R_Nx,R_Ny) array of floats) the weights of each scan position in this class

Return type:

(2-tuple)

get_candidate_class_BPs(i)

Get a single candidate class, returning its BP weights.

Accepts:

i (int) the class index

Returns:

class_BPs (length N_feat array of floats) the weights of the N_feat Bragg peaks for

this class

get_candidate_class_image(i)

Get a single candidate class, returning its scan position weights.

Parameters:

i (int) – the class index

Returns:

the weights of each scan position in this class

Return type:

(shape (R_Nx,R_Ny) array of floats)

py4DSTEM.process.classification.braggvectorclassification.get_braggpeak_labels_by_scan_position(braggpeaks, Qx, Qy, max_dist=None)

For each scan position, gets a set of integers, specifying the bragg peaks at this scan position.

From a set of positions in diffraction space (Qx,Qy), assign each detected bragg peak in the PointListArray braggpeaks a label corresponding to the index of the closest position; thus for a bragg peak at (qx,qy), if the closest position in (Qx,Qy) is (Qx[i],Qy[i]), assign this peak the label i. This is equivalent to assigning each bragg peak (qx,qy) a label according to the Voronoi region it lives in, given a voronoi tesselation seeded from the points (Qx,Qy).

For each scan position, get the set of all indices i for all bragg peaks found at this scan position.

Parameters:
  • braggpeaks (PointListArray) – Bragg peaks; must have coords ‘qx’ and ‘qy’

  • Qx (ndarray of floats) – x-coords of the voronoi points

  • Qy (ndarray of floats) – y-coords of the voronoi points

  • max_dist (None or number) – maximum distance from a given voronoi point a peak can be and still be associated with this label

Returns:

(list of lists of sets) the labels found at each scan position. Scan position (Rx,Ry) is accessed via braggpeak_labels[Rx][Ry]

py4DSTEM.process.classification.braggvectorclassification.get_initial_classes(braggpeak_labels, N, thresh=0.3, BP_fraction_thresh=0.1, max_iterations=200, n_corr_init=2)

From the sets of Bragg peaks present at each scan position, get an initial guess classes at which Bragg peaks should be grouped together into classes.

The algorithm is as follows: 1. Calculate an n-point correlation function, i.e. the joint probability of any given n BPs coexisting in a diffraction pattern. n is controlled by n_corr_init, and must be 2 or 3. peaks i, j, and k are all in the same DP. 2. Find the BP triplet maximizing the 3-point function; include these three BPs in a class. 3. Get all DPs containing the class BPs. From these, find the next most likely BP to also be present. If its probability of coexisting with the known class BPs is greater than thresh, add it to the class and repeat this step. Otherwise, proceed to the next step. 4. Check: if the new class is the same as a class that has already been found, OR if the fraction of BPs which have not yet been placed in a class is less than BP_fraction_thresh, or more than max_iterations have been attempted, finish, returning all classes. Otherwise, set all slices of the 3-point function containing the BPs in the new class to zero, and begin a new iteration, starting at step 2 using the new, altered 3-point function.

Parameters:
  • N (int) – the total number of indexed Bragg peaks in the 4D-STEM dataset

  • braggpeak_labels (list of lists of sets) – the Bragg peak labels found at each scan position; see get_braggpeak_labels_by_scan_position().

  • thresh (float in [0,1]) – threshold for adding new BPs to a class

  • BP_fraction_thresh (float in [0,1]) – algorithm terminates if fewer than this fraction of the BPs have not been assigned to a class

  • max_iterations (int) – algorithm terminates after this many iterations

  • n_corr_init (int) – seed new classes by finding maxima of the n-point joint probability function. Must be 2 or 3.

Returns:

the sets of Bragg peaks constituting the classes

Return type:

(list of sets)

py4DSTEM.process.classification.classutils.get_class_DP(datacube, class_image, thresh=0.01, xshifts=None, yshifts=None, darkref=None, intshifts=True)

Get the average diffraction pattern for the class described in real space by class_image.

Parameters:
  • datacube (DataCube) – a datacube

  • class_image (2D array) – the weight of the class at each position in real space

  • thresh (float) – only include diffraction patterns for scan positions with a value greater than or equal to thresh in class_image

  • xshifts (2D array, or None) – the x diffraction shifts at each real space pixel. If None, no shifting is performed.

  • yshifts (2D array, or None) – the y diffraction shifts at each real space pixel. If None, no shifting is performed.

  • darkref (2D array, or None) – background to remove from each diffraction pattern

  • intshifts (bool) – if True, round shifts to the nearest integer to speed up computation

Returns:

the average diffraction pattern for the class

Return type:

(2D array)

py4DSTEM.process.classification.classutils.get_class_DP_without_Bragg_scattering(datacube, class_image, braggpeaks, radius, x0, y0, thresh=0.01, xshifts=None, yshifts=None, darkref=None, intshifts=True)

Get the average diffraction pattern, removing any Bragg scattering, for the class described in real space by class_image.

Bragg scattering is eliminated by masking circles of size radius about each of the detected peaks in braggpeaks in each diffraction pattern before adding to the average image. Importantly, braggpeaks refers to the peak positions in the raw data - i.e. BEFORE any shift correction is applied. Passing shifted Bragg peaks will yield incorrect results. For speed, the Bragg peaks are removed with a binary mask, rather than a continuous sigmoid, so selecting a radius that is slightly (~1 pix) larger than the disk size is recommended.

Parameters:
  • datacube (DataCube) – a datacube

  • class_image (2D array) – the weight of the class at each position in real space

  • braggpeaks (PointListArray) – the detected Bragg peak positions, with respect to the raw data (i.e. not diffraction shift or ellipse corrected)

  • radius (number) – the radius to mask about each detected Bragg peak - should be slightly larger than the disk radius

  • x0 (number) – x-position of the optic axis

  • y0 (number) – y-position of the optic axis

  • thresh (float) – only include diffraction patterns for scan positions with a value greater than or equal to thresh in class_image

  • xshifts (2D array, or None) – the x diffraction shifts at each real space pixel. If None, no shifting is performed.

  • yshifts (2D array, or None) – the y diffraction shifts at each real space pixel. If None, no shifting is performed.

  • darkref (2D array, or None) – background to remove from each diffraction pattern

  • intshifts (bool) – if True, round shifts to the nearest integer to speed up computation

Returns:

class_DP (2D array) the average diffraction pattern for the class

class py4DSTEM.process.classification.featurization.Featurization(features, R_Nx, R_Ny, name)

A class for feature selection, modification, and classification of 4D-STEM data based on a user defined array of input features for each pattern. Features are stored under Featurization. Features and can be used for a variety of unsupervised classification tasks.

Initialization methods:
__init__:

Creates instance of featurization

concatenate_features:

Creates instance of featurization from a list of featurization instances

from_braggvectors:

Creates instance of featurization from a BraggVectors instance

Feature Dictionary Modification Methods
add_feature:

Adds features to the features array

remove_feature:

Removes features to the features array

Feature Preprocessing Methods
MinMaxScaler:

Performs sklearn MinMaxScaler operation on features stored at a key

RobustScaler:

Performs sklearn RobustScaler operation on features stored at a key

mean_feature:

Takes the rowwise average of a matrix stored at a key, such that only one column is left, reducing a set of n features down to 1 feature per pattern.

median_feature:

Takes the rowwise median of a matrix stored at a key, such that only one column is left, reducing a set of n features down to 1 feature per pattern.

max_feature:

Takes the rowwise max of a matrix stored at a key, such that only one column is left, reducing a set of n features down to 1 feature per pattern.

Classification Methods
PCA:

Principal Component Analysis to refine features.

ICA:

Independent Component Analysis to refine features.

NMF:

Performs either traditional or iterative Nonnegative Matrix Factorization (NMF) to refine features.

GMM:

Gaussian mixture model to predict class labels. Fits a gaussian based on covariance of features.

Class Examination Methods
get_class_DPs:

Gets weighted class diffraction patterns (DPs) for an NMF or GMM operation

get_class_ims:

Gets weighted class images (ims) for an NMF or GMM operation

__init__(features, R_Nx, R_Ny, name)

Initializes classification instance.

This method: 1. Generates key:value pair to access input features 2. Initializes the empty dictionaries for feature modification and classification

Parameters:
  • features (list) – A list of ndarrays which will each be associated with value stored at the key in the same index within the list

  • R_Nx (int) – The real space x dimension of the dataset

  • R_Ny (int) – The real space y dimension of the dataset

  • name (str) – The name of the featurization object

Returns:

New Featurization instance

Return type:

new_instance

from_braggvectors(bins_x, bins_y, intensity_scale, name, mask=None)

Initialize a featurization instance from a BraggVectors instance

Parameters:
  • braggvectors (BraggVectors) – BraggVectors instance containing calibrations

  • bins_x (int) – Number of pixels per bin in x direction

  • bins_y (int) – Number of pixels per bin in y direction

  • intensity_scale (float) – Value to scale intensity of detected disks by

  • name (str) – Name of featurization instance

  • mask (bool) – Mask to remove disks in unwanted positions in diffraction space

Returns:

Featurization instance

Return type:

new_instance

Details:

Transforms the calibrated pointlistarray in BraggVectors instance into a numpy array that can be clustered using the methods in featurization.

concatenate_features(name)

Concatenates featurization instances (features) and outputs a new Featurization instance containing the concatenated features from each featurization instance. R_Nx, R_Ny will be inherited from the featurization instances and must be consistent across objects.

Parameters:
  • features (list) – A list of keys to be concatenated into one array

  • name (str) – The name of the featurization instance

Returns:

Featurization instance

Return type:

new_instance

add_features(feature)

Add a feature to the end of the features array

Parameters:
  • key (int, float, str) – A key in which a feature can be accessed from

  • feature (ndarray) – The feature associated with the key

delete_features(index)

Deletes feature columns from the feature array

Parameters:

index (int, list) – A key which will be removed

mean_feature(index)

Takes columnwise mean and replaces features in ‘index’.

Parameters:

index (list of int) – Indices of features to take the mean of. New feature array is placed in self.features.

median_feature(index)

Takes columnwise median and replaces features in ‘index’. New feature array is placed in self.features.

Parameters:

index (list of int) – Indices of features to take the median of.

max_feature(index)

Takes columnwise max and replaces features in ‘index’. New feature array is placed in self.features.

Parameters:

index (list of int) – Indices of features to take the max of.

MinMaxScaler(return_scaled=True)

Uses sklearn MinMaxScaler to scale a subset of the input features. Replaces a feature with the positive shifted array.

Parameters:

return_scaled (bool) – returns the scaled array

RobustScaler(return_scaled=True)

Uses sklearn RobustScaler to scale a subset of the input features. Replaces a feature with the positive shifted array.

Parameters:

return_scaled (bool) – returns the scaled array

shift_positive(return_scaled=True)

Replaces a feature with the positive shifted array.

Parameters:

return_scaled (bool) – returns the scaled array

PCA(components, return_results=False)

Performs PCA on features

Parameters:

components (list) – A list of ints for each key. This will be the output number of features

ICA(components, return_results=True)

Performs ICA on features

Parameters:

components (list) – A list of ints for each key. This will be the output number of features

NMF(max_components, num_models, merge_thresh=1, max_iterations=1, random_seed=None, save_all_models=True, return_results=False)

Performs either traditional Nonnegative Matrix Factoriation (NMF) or iteratively on input features. For Traditional NMF:

set either merge_threshold = 1, max_iterations = 1, or both. Default is to set

Parameters:
  • max_components (int) – Number of initial components to start the first NMF iteration

  • merge_thresh (float) – Correlation threshold to merge features

  • num_models (int) – Number of independent models to run (number of learners that will be combined in consensus).

  • max_iterations (int) – Number of iterations. Default 1, which runs traditional NMF

  • random_seed (int) – Random seed.

  • save_all_models (bool) – Whether or not to return all of the models - default is to return all outputs for consensus clustering. if False, will only return the model with the lowest NMF reconstruction error.

  • return_results (bool) – Whether or not to return the final class weights

Details:

This method may require trial and error for proper selection of parameters. To perform traditional NMF, the defaults should be used:

merge_thresh = 1 max_iterations = 1

Note that the max_components in this case will be equivalent to the number of classes the NMF model identifies.

Iterative NMF calculates the correlation between all of the output columns from an NMF iteration, merges the features correlated above the merge_thresh, and performs NMF until either max_iterations is reached or until no more columns are correlated above merge_thresh.

GMM(cv, components, num_models, random_seed=None, return_results=False)

Performs gaussian mixture model on input features

Parameters:
  • cv (str) – Covariance type - must be ‘spherical’, ‘tied’, ‘diag’, or ‘full’

  • components (int) – Number of components

  • num_models (int) – Number of models to run

  • random_seed (int) – Random seed

get_class_DPs(datacube, method, thresh)

Returns weighted class patterns based on classification instance datacube must be vectorized in real space (shape = (R_Nx * R_Ny, 1, Q_Nx, Q_Ny)

Parameters:
  • classification_method (str) – Either ‘nmf’ or ‘gmm’ - finds location of clusters

  • datacube (py4DSTEM datacube) – Vectorized in real space, with shape (R_Nx * R_Ny, Q_Nx, Q_Ny)

get_class_ims(classification_method)

Returns weighted class maps based on classification instance

Parameters:

classification_method (str) – Location to retrieve class images from - NMF, GMM, PCA, or ICA

spatial_separation(size, threshold=0, method=None, clean=True)

Identify spatially distinct regions from class images and separate based on a threshold and size.

Parameters:
  • size (int) – Number of pixels which is the minimum to keep a class - all spatially distinct regions with less than ‘size’ pixels will be removed

  • threshold (float) – Intensity weight of a component to keep

  • method (str) – (Optional) Filter method, default None. Accepts options ‘yen’ and ‘otsu’.

  • clean (bool) – Whether or not to ‘clean’ cluster sets based on overlap, i.e. remove clusters that do not have any unique components

consensus(threshold=0, location='spatially_separated_ims', split=0, method='mean', drop_bins=0)

Consensus Clustering takes the outcome of a prepared set of 2D images from each cluster and averages the outcomes.

Parameters:
  • threshold (float) – Threshold weights, default 0

  • location (str) – Where to get the consensus from - after spatial separation = ‘spatially_separated_ims’

  • split_value (float) – Threshold in which to separate classes during label correspondence (Default 0). This should be proportional to the expected class weights- the sum of the weights in the current class image that match nonzero values in each bin are calculated and then checked for splitting.

  • method (str) – Method in which to combine the consensus clusters - either mean or median.

  • drop_bins (int) – Number of clusters needed in each class to keep cluster set in the consensus. Default 0, meaning

Details:

This method involves 2 steps: Label correspondence and consensus clustering.

Label correspondence sorts the classes found by the independent models into bins based on class overlap in real space. Arguments related to label correspondence are the threshold and split_value. The threshold is related to the weights of the independent classes. If the weight of the observation in the class is less than the threshold, it will be set to 0. The split_value indicates the extent of similarity the independent classes must have before intializing a new bin. The default is 0 - this means if the class of interest has 0 overlap with the identified bins, a new bin will be created. The value is based on the sum of the weights in the current class image that match the nonzero values in the current bins.

Consensus clustering combines these sorted bin into 1 class based on the selected method (either ‘mean’ which takes the average of the bin, or ‘median’ which takes the median of the bin). Bins with less than the drop_bins value will not be included in the final results.

diffraction

py4DSTEM.process.diffraction.WK_scattering_factors.compute_WK_factor(g: ndarray, Z: int, accelerating_voltage: float, thermal_sigma: float | None = None, include_core: bool = True, include_phonon: bool = True, verbose=False) complex128

Compute the Weickenmeier-Kohl atomic scattering factors, using the parameterization of the elastic part and computation of the inelastic part found in EMsoftLib/others.f90. Return value should be in Å.

This implementation always returns the absorptive, relativistically corrected factors.

Currently this is mostly a direct translation of the Fortran code, along with the accompanying comments from the original in quotation marks. Colin Ophus vectorized it around v0.13.17. Currently it is only vectorized over g (i.e. Z and all other args must be a single value.)

This method uses an 8-parameter fit to the elastic form factors, and then computes the absorptive form factors using an analytic solution based on that fitting function.

Args: (note that these values cannot be arrays: the code is not vectorized)
g (float/ndarray): Scattering vector magnitude in the crystallographic/py4DSTEM

convention, 1/d_hkl in units of 1/Å

Z (int): Atomic number. Data are available for H thru Cf (1 thru 98) accelerating_voltage (float): Accelerating voltage in eV. thermal_sigma (float): RMS atomic displacement for TDS, in Å

(This is often written as 〈u〉in papers)

include_core (bool): If True, include the core loss contribution to the absorptive

form factors.

include_phonon (bool): If True, include the phonon/TDS contribution to the

absorptive form factors.

Returns:

The computed atomic form factor

Return type:

Fscatt (np.complex128)

py4DSTEM.process.diffraction.WK_scattering_factors.RIH2(X)

WERTET X*EXP(-X)*EI(X) AUS FUER GROSSE X DURCH INTERPOLATION DER TABELLE … AUS ABRAMOWITZ

class py4DSTEM.process.diffraction.crystal.Crystal(positions, numbers, cell, occupancy=None)

A class storing a single crystal structure, and associated diffraction data.

orientation_plan(zone_axis_range: ndarray = array([[0, 1, 1], [1, 1, 1]]), angle_step_zone_axis: float = 2.0, angle_coarse_zone_axis: float | None = None, angle_refine_range: float | None = None, angle_step_in_plane: float = 2.0, accel_voltage: float = 300000.0, corr_kernel_size: float = 0.08, radial_power: float = 1.0, intensity_power: float = 0.25, calculate_correlation_array=True, tol_peak_delete=None, tol_distance: float = 0.01, fiber_axis=None, fiber_angles=None, figsize: list | tuple | ndarray = (6, 6), CUDA: bool = False, progress_bar: bool = True)

Calculate the rotation basis arrays for an SO(3) rotation correlogram.

Parameters:
  • zone_axis_range (float) –

    Row vectors give the range for zone axis orientations. If user specifies 2 vectors (2x3 array), we start at [0,0,1]

    to make z-x-z rotation work.

    If user specifies 3 vectors (3x3 array), plan will span these vectors. Setting to ‘full’ as a string will use a hemispherical range. Setting to ‘half’ as a string will use a quarter sphere range. Setting to ‘fiber’ as a string will make a spherical cap around a given vector. Setting to ‘auto’ will use pymatgen to determine the point group symmetry

    of the structure and choose an appropriate zone_axis_range

  • angle_step_zone_axis (float) – Approximate angular step size for zone axis search [degrees]

  • angle_coarse_zone_axis (float) – Coarse step size for zone axis search [degrees]. Setting to None uses the same value as angle_step_zone_axis.

  • angle_refine_range (float) – Range of angles to use for zone axis refinement. Setting to None uses same value as angle_coarse_zone_axis.

  • angle_step_in_plane (float) – Approximate angular step size for in-plane rotation [degrees]

  • accel_voltage (float) – Accelerating voltage for electrons [Volts]

  • corr_kernel_size (float) – Correlation kernel size length. The size of the overlap kernel between the measured Bragg peaks and diffraction library Bragg peaks. [1/Angstroms]

  • radial_power (float) – Power for scaling the correlation intensity as a function of the peak radius

  • intensity_power (float) – Power for scaling the correlation intensity as a function of the peak intensity

  • calculate_correlation_array (bool) – Set to false to skip calculating the correlation array. This is useful when we only want the angular range / rotation matrices.

  • tol_peak_delete (float) – Distance to delete peaks for multiple matches. Default is kernel_size * 0.5

  • tol_distance (float) – Distance tolerance for radial shell assignment [1/Angstroms]

  • fiber_axis (float) – (3,) vector specifying the fiber axis

  • fiber_angles (float) – (2,) vector specifying angle range from fiber axis, and in-plane angular range [degrees]

  • cartesian_directions (bool) – When set to true, all zone axes and projection directions are specified in Cartesian directions.

  • figsize (float) – (2,) vector giving the figure size

  • CUDA (bool) – Use CUDA for the Fourier operations.

  • progress_bar (bool) – If false no progress bar is displayed

match_orientations(bragg_peaks_array: PointListArray, num_matches_return: int = 1, min_angle_between_matches_deg=None, min_number_peaks: int = 3, inversion_symmetry: bool = True, multiple_corr_reset: bool = True, return_orientation: bool = True, progress_bar: bool = True)
Parameters:
  • bragg_peaks_array (PointListArray) – PointListArray containing the Bragg peaks and intensities, with calibrations applied

  • num_matches_return (int) – return these many matches as 3th dim of orient (matrix)

  • min_angle_between_matches_deg (int) – Minimum angle between zone axis of multiple matches, in degrees. Note that I haven’t thought how to handle in-plane rotations, since multiple matches are possible.

  • min_number_peaks (int) – Minimum number of peaks required to perform ACOM matching

  • inversion_symmetry (bool) – check for inversion symmetry in the matches

  • multiple_corr_reset (bool) – keep original correlation score for multiple matches

  • return_orientation (bool) – Return orientation map from function for inspection. The map is always stored in the Crystal object.

  • progress_bar (bool) – Show or hide the progress bar

match_single_pattern(bragg_peaks: PointList, num_matches_return: int = 1, min_angle_between_matches_deg=None, min_number_peaks=3, inversion_symmetry=True, multiple_corr_reset=True, plot_polar: bool = False, plot_corr: bool = False, returnfig: bool = False, figsize: list | tuple | ndarray = (12, 4), verbose: bool = False)

Solve for the best fit orientation of a single diffraction pattern.

Parameters:
  • bragg_peaks (PointList) – numpy array containing the Bragg positions and intensities (‘qx’, ‘qy’, ‘intensity’)

  • num_matches_return (int) – return these many matches as 3th dim of orient (matrix)

  • min_angle_between_matches_deg (int) – Minimum angle between zone axis of multiple matches, in degrees. Note that I haven’t thought how to handle in-plane rotations, since multiple matches are possible.

  • min_number_peaks (int) – Minimum number of peaks required to perform ACOM matching

  • bool (multiple_corr_reset) – check for inversion symmetry in the matches

  • bool – keep original correlation score for multiple matches

  • subpixel_tilt (bool) – set to false for faster matching, returning the nearest corr point

  • plot_polar (bool) – set to true to plot the polar transform of the diffraction pattern

  • plot_corr (bool) – set to true to plot the resulting correlogram

  • returnfig (bool) – return figure handles

  • figsize (list) – size of figure

  • verbose (bool) – Print the fitted zone axes, correlation scores

  • CUDA (bool) – Enable CUDA for the FFT steps

Returns:

  • orientation (Orientation) – Orientation class containing all outputs

  • fig, ax (handles) – Figure handles for the plotting output

cluster_grains(threshold_add=1.0, threshold_grow=0.1, angle_tolerance_deg=5.0, progress_bar=True)

Cluster grains using rotation criterion, and correlation values.

Parameters:
  • threshold_add (float) – Minimum signal required for a probe position to initialize a cluster.

  • threshold_grow (float) – Minimum signal required for a probe position to be added to a cluster.

  • angle_tolerance_deg (float) – Rotation rolerance for clustering grains.

  • progress_bar (bool) – Turns on the progress bar for the polar transformation

cluster_orientation_map(stripe_width=(2, 2), area_min=2)

Produce a new orientation map from the clustered grains. Use a stripe pattern for the overlapping grains.

Parameters:
  • stripe_width ((int,int)) – Width of stripes for plotting maps with overlapping grains

  • area_min ((int)) – Minimum size of grains to include

Returns:

The clustered orientation map

Return type:

orientation_map

calculate_strain(bragg_peaks_array: PointListArray, orientation_map: OrientationMap, corr_kernel_size=None, sigma_excitation_error=0.02, tol_excitation_error_mult: float = 3, tol_intensity: float = 0.0001, k_max: float | None = None, min_num_peaks=5, rotation_range=None, mask_from_corr=True, corr_range=(0, 2), corr_normalize=True, progress_bar=True)

This function takes in both a PointListArray containing Bragg peaks, and a corresponding OrientationMap, and uses least squares to compute the deformation tensor which transforms the simulated diffraction pattern into the experimental pattern, for all probe positons.

TODO: add robust fitting?

Parameters:
  • bragg_peaks_array (PointListArray) – All Bragg peaks

  • orientation_map (OrientationMap) – Orientation map generated from ACOM

  • corr_kernel_size (float) – Correlation kernel size - if user does not specify, uses self.corr_kernel_size.

  • sigma_excitation_error (float) – sigma value for envelope applied to s_g (excitation errors) in units of inverse Angstroms

  • tol_excitation_error_mult (float) – tolerance in units of sigma for s_g inclusion

  • tol_intensity (np float) – tolerance in intensity units for inclusion of diffraction spots

  • k_max (float) – Maximum scattering vector

  • min_num_peaks (int) – Minimum number of peaks required.

  • rotation_range (float) – Maximum rotation range in radians (for symmetry reduction).

  • progress_bar (bool) – Show progress bar

  • mask_from_corr (bool) – Use ACOM correlation signal for mask

  • corr_range (np.ndarray) – Range of correlation signals for mask

  • corr_normalize (bool) – Normalize correlation signal before masking

Returns:

strain tensor

Return type:

strain_map (RealSlice)

symmetry_reduce_directions(orientation, match_ind=0, plot_output=False, figsize=(15, 6), el_shift=0.0, az_shift=-30.0)

This function calculates the symmetry-reduced cartesian directions from and orientation matrix stored in orientation.matrix, and outputs them into orientation.family. It optionally plots the 3D output.

save_ang_file(file_name, orientation_map, ind_orientation=0, pixel_size=1.0, pixel_units='px', transpose_xy=True, flip_x=False)

This function outputs an ascii text file in the .ang format, containing the Euler angles of an orientation map.

Parameters:
  • file_name (str) – Path to save .ang file.

  • orientation_map (OrientationMap) – Class containing orientation matrices, correlation values, etc.

  • ind_orientation (int) – Which orientation match to plot if num_matches > 1

  • pixel_size (float) – Pixel size, if known.

  • pixel_units (str) – Units of the pixel size

  • transpose_xy (bool) – Transpose x and y pixel coordinates.

  • flip_x (bool) – Swap x direction pixels (after transpose).

Returns:

nothing

plot_structure(orientation_matrix: ndarray | None = None, zone_axis_lattice: ndarray | None = None, proj_x_lattice: ndarray | None = None, zone_axis_cartesian: ndarray | None = None, proj_x_cartesian: ndarray | None = None, size_marker: float = 400, tol_distance: float = 0.001, plot_limit: ndarray | None = None, camera_dist: float | None = None, show_axes: bool = False, perspective_axes: bool = True, figsize: tuple | list | ndarray = (8, 8), returnfig: bool = False)

Quick 3D plot of the untit cell /atomic structure.

Parameters:
  • orientation_matrix (array) – (3,3) orientation matrix, where columns represent projection directions.

  • zone_axis_lattice (array) – (3,) projection direction in lattice indices

  • proj_x_lattice (array) – (3,) x-axis direction in lattice indices

  • zone_axis_cartesian (array) – (3,) cartesian projection direction

  • proj_x_cartesian (array) – (3,) cartesian projection direction

  • scale_markers (float) – Size scaling for markers

  • tol_distance (float) – Tolerance for repeating atoms on edges on cell boundaries.

  • plot_limit (float) – (2,3) numpy array containing x y z plot min and max in columns. Default is 1.1* unit cell dimensions.

  • camera_dist (float) – Move camera closer to the plot (relative to matplotlib default of 10)

  • show_axes (bool) – Whether to plot axes or not.

  • perspective_axes (bool) – Select either perspective (true) or orthogonal (false) axes

  • figsize (2 element float) – Size scaling of figure axes.

  • returnfig (bool) – Return figure and axes handles.

Returns:

fig, ax (optional) figure and axes handles

plot_structure_factors(orientation_matrix: ndarray | None = None, zone_axis_lattice: ndarray | None = None, proj_x_lattice: ndarray | None = None, zone_axis_cartesian: ndarray | None = None, proj_x_cartesian: ndarray | None = None, scale_markers: float = 1000.0, plot_limit: list | tuple | ndarray | None = None, camera_dist: float | None = None, show_axes: bool = True, perspective_axes: bool = True, figsize: list | tuple | ndarray = (8, 8), returnfig: bool = False)

3D scatter plot of the structure factors using magnitude^2, i.e. intensity.

Parameters:
  • orientation_matrix (array) – (3,3) orientation matrix, where columns represent projection directions.

  • zone_axis_lattice (array) – (3,) projection direction in lattice indices

  • proj_x_lattice (array) – (3,) x-axis direction in lattice indices

  • zone_axis_cartesian (array) – (3,) cartesian projection direction

  • proj_x_cartesian (array) – (3,) cartesian projection direction

  • scale_markers (float) – size scaling for markers

  • plot_limit (float) – x y z plot limits, default is [-1 1]*self.k_max

  • camera_dist (float) – Move camera closer to the plot (relative to matplotlib default of 10)

  • show_axes (bool) – Whether to plot axes or not.

  • perspective_axes (bool) – Select either perspective (true) or orthogonal (false) axes

  • figsize (2 element float) – size scaling of figure axes

  • returnfig (bool) – set to True to return figure and axes handles

Returns:

fig, ax (optional) figure and axes handles

plot_scattering_intensity(k_min=0.0, k_max=None, k_step=0.001, k_broadening=0.0, k_power_scale=0.0, int_power_scale=0.5, int_scale=1.0, remove_origin=True, bragg_peaks=None, bragg_k_power=0.0, bragg_intensity_power=1.0, bragg_k_broadening=0.005, figsize: list | tuple | ndarray = (10, 4), returnfig: bool = False)

1D plot of the structure factors

Parameters:
  • k_min (float) – min k value for profile range.

  • k_max (float) – max k value for profile range.

  • k_step (float) – Step size of k in profile range.

  • k_broadening (float) – Broadening of simulated pattern.

  • k_power_scale (float) – Scale SF intensities by k**k_power_scale.

  • int_power_scale (float) – Scale SF intensities**int_power_scale.

  • int_scale (float) – Scale output profile by this value.

  • remove_origin (bool) – Remove origin from plot.

  • bragg_peaks (BraggVectors) – Passed in bragg_peaks for comparison with simulated pattern.

  • bragg_k_power (float) – bragg_peaks scaled by k**bragg_k_power.

  • bragg_intensity_power (float) – bragg_peaks scaled by intensities**bragg_intensity_power.

  • bragg_k_broadening (float) – Broadening applied to bragg_peaks.

  • figsize (list, tuple, np.ndarray) – Figure size for plot.

  • (bool) (returnfig) – Return figure and axes handles if this is True.

Returns:

figure and axes handles

Return type:

fig, ax (optional)

plot_orientation_zones(azim_elev: list | tuple | ndarray | None = None, proj_dir_lattice: list | tuple | ndarray | None = None, proj_dir_cartesian: list | tuple | ndarray | None = None, tol_den=10, marker_size: float = 20, plot_limit: list | tuple | ndarray = array([-1.1, 1.1]), figsize: list | tuple | ndarray = (8, 8), returnfig: bool = False)

3D scatter plot of the structure factors using magnitude^2, i.e. intensity.

Parameters:
  • azim_elev (array) – az and el angles for plot

  • proj_dir_lattice (array) – (3,) projection direction in lattice

  • proj_dir_cartesian – (array): (3,) projection direction in cartesian

  • tol_den (int) – tolerance for rational index denominator

  • dir_proj (float) – projection direction, either [elev azim] or normal vector Default is mean vector of self.orientation_zone_axis_range rows

  • marker_size (float) – size of markers

  • plot_limit (float) – x y z plot limits, default is [0, 1.05]

  • figsize (2 element float) – size scaling of figure axes

  • returnfig (bool) – set to True to return figure and axes handles

Returns:

fig, ax (optional) figure and axes handles

plot_orientation_plan(index_plot: int = 0, zone_axis_lattice: ndarray | None = None, zone_axis_cartesian: ndarray | None = None, figsize: list | tuple | ndarray = (14, 6), returnfig: bool = False)

3D scatter plot of the structure factors using magnitude^2, i.e. intensity.

Parameters:
  • index_plot (int) – which index slice to plot

  • zone_axis_plot (3 element float) – which zone axis slice to plot

  • figsize (2 element float) – size scaling of figure axes

  • returnfig (bool) – set to True to return figure and axes handles

Returns:

fig, ax (optional) figure and axes handles

plot_orientation_maps(orientation_map=None, orientation_ind: int = 0, dir_in_plane_degrees: float = 0.0, corr_range: ndarray = array([0, 5]), corr_normalize: bool = True, scale_legend: bool | None = None, figsize: list | tuple | ndarray = (16, 5), figbound: list | tuple | ndarray = (0.01, 0.005), show_axes: bool = True, camera_dist=None, plot_limit=None, plot_layout=0, swap_axes_xy_limits=False, returnfig: bool = False, progress_bar=False)

Plot the orientation maps.

Parameters:
  • orientation_map (OrientationMap) – Class containing orientation matrices, correlation values, etc. Optional - can reference internally stored OrientationMap.

  • orientation_ind (int) – Which orientation match to plot if num_matches > 1

  • dir_in_plane_degrees (float) – In-plane angle to plot in degrees. Default is 0 / x-axis / vertical down.

  • corr_range (np.ndarray) – Correlation intensity range for the plot

  • corr_normalize (bool) – If true, set mean correlation to 1.

  • scale_legend (float) – 2 elements, x and y scaling of legend panel

  • figsize (array) – 2 elements defining figure size

  • figbound (array) – 2 elements defining figure boundary

  • show_axes (bool) – Flag setting whether orienation map axes are visible.

  • camera_dist (float) – distance of camera from legend

  • plot_limit (array) – 2x3 array defining plot boundaries of legend

  • plot_layout (int) – subplot layout: 0 - 1 row, 3 col 1 - 3 row, 1 col

  • swap_axes_xy_limits (bool) – swap x and y boundaries for legend (not sure why we need this in some cases)

  • returnfig (bool) – set to True to return figure and axes handles

  • progress_bar (bool) – Enable progressbar when calculating orientation images.

Returns:

RGB images fig, axs (handles): Figure and axes handes for the

Return type:

images_orientation (int)

Note

Currently, no symmetry reduction. Therefore the x and y orientations are going to be correct only for [001][011][111] orientation triangle.

plot_fiber_orientation_maps(orientation_map, orientation_ind: int = 0, symmetry_order: int | None = None, symmetry_mirror: bool = False, dir_in_plane_degrees: float = 0.0, corr_range: ndarray = array([0, 2]), corr_normalize: bool = True, show_axes: bool = True, medfilt_size: int | None = None, cmap_out_of_plane: str = 'plasma', leg_size: int = 200, figsize: list | tuple | ndarray = (12, 8), figbound: list | tuple | ndarray = (0.005, 0.04), returnfig: bool = False)

Generate and plot the orientation maps from fiber texture plots.

Parameters:
  • orientation_map (OrientationMap) – Class containing orientation matrices, correlation values, etc.

  • orientation_ind (int) – Which orientation match to plot if num_matches > 1

  • dir_in_plane_degrees (float) – Reference in-plane angle (degrees). Default is 0 / x-axis / vertical down.

  • corr_range (np.ndarray) – Correlation intensity range for the plot

  • corr_normalize (bool) – If true, set mean correlation to 1.

  • show_axes (bool) – Flag setting whether orienation map axes are visible.

  • figsize (array) – 2 elements defining figure size

  • figbound (array) – 2 elements defining figure boundary

  • returnfig (bool) – set to True to return figure and axes handles

Returns:

RGB images fig, axs (handles): Figure and axes handes for the

Return type:

images_orientation (int)

Note

Currently, no symmetry reduction. Therefore the x and y orientations are going to be correct only for [001][011][111] orientation triangle.

plot_clusters(area_min=2, outline_grains=True, outline_thickness=1, fill_grains=0.25, smooth_grains=1.0, cmap='viridis', figsize=(8, 8), returnfig=False)

Plot the clusters as an image.

Parameters:
  • area_min (int (optional)) – Min cluster size to include, in units of probe positions.

  • outline_grains (bool (optional)) – Set to True to draw grains with outlines

  • outline_thickness (int (optional)) – Thickenss of the grain outline

  • fill_grains (float (optional)) – Outlined grains are filled with this value in pixels.

  • smooth_grains (float (optional)) – Grain boundaries are smoothed by this value in pixels.

  • figsize (tuple) – Size of the figure panel

  • returnfig (bool) – Setting this to true returns the figure and axis handles

Returns:

Figure and axes handles

Return type:

fig, ax (optional)

plot_cluster_size(area_min=None, area_max=None, area_step=1, weight_intensity=False, pixel_area=1.0, pixel_area_units='px^2', figsize=(8, 6), returnfig=False)

Plot the cluster sizes

Parameters:
  • area_min (int (optional)) – Min area to include in pixels^2

  • area_max (int (optional)) – Max area bin in pixels^2

  • area_step (int (optional)) – Step size of the histogram bin in pixels^2

  • weight_intensity (bool) – Weight histogram by the peak intensity.

  • pixel_area (float) – Size of pixel area unit square

  • pixel_area_units (string) – Units of the pixel area

  • figsize (tuple) – Size of the figure panel

  • returnfig (bool) – Setting this to true returns the figure and axis handles

Returns:

Figure and axes handles

Return type:

fig, ax (optional)

calibrate_pixel_size(bragg_peaks, scale_pixel_size=1.0, bragg_k_power=1.0, bragg_intensity_power=1.0, k_min=0.0, k_max=None, k_step=0.002, k_broadening=0.002, fit_all_intensities=True, set_calibration_in_place=False, verbose=True, plot_result=False, figsize: list | tuple | ndarray = (12, 6), returnfig=False)

Use the calculated structure factor scattering lengths to compute 1D diffraction patterns, and solve the best-fit relative scaling between them. Returns the fit pixel size in Å^-1.

Parameters:
  • bragg_peaks (BraggVectors) – Input Bragg vectors.

  • scale_pixel_size (float) – Initial guess for scaling of the existing pixel size If the pixel size is currently uncalibrated, this is a guess of the pixel size in Å^-1. If the pixel size is already (approximately) calibrated, this is the scaling factor to correct that existing calibration.

  • bragg_k_power (float) – Input Bragg peak intensities are multiplied by k**bragg_k_power to change the weighting of longer scattering vectors

  • bragg_intensity_power (float) – Input Bragg peak intensities are raised power **bragg_intensity_power.

  • k_min (float) – min k value for fitting range (Å^-1)

  • k_max (float) – max k value for fitting range (Å^-1)

  • k_step (float) step size of k in fitting range (Å^-1) –

  • k_broadening (float) – Initial guess for Gaussian broadening of simulated pattern (Å^-1)

  • fit_all_intensities (bool) – Set to true to allow all peak intensities to change independently False forces a single intensity scaling.

  • set_calibration (bool) – if True, set the fit pixel size to the calibration metadata, and calibrate bragg_peaks

  • verbose (bool) – Output the calibrated pixel size.

  • plot_result (bool) – Plot the resulting fit.

  • figsize (list, tuple, np.ndarray) – Figure size of the plot.

  • returnfig (bool) – Return handles figure and axis

Returns:

fig, ax – Figure and axis handles, if returnfig=True.

Return type:

handles, optional

calibrate_unit_cell(bragg_peaks, coef_index=None, coef_update=None, bragg_k_power=1.0, bragg_intensity_power=1.0, k_min=0.0, k_max=None, k_step=0.005, k_broadening=0.02, fit_all_intensities=True, verbose=True, plot_result=False, figsize: list | tuple | ndarray = (12, 6), returnfig=False)

Solve for the best fit scaling between the computed structure factors and bragg_peaks.

Parameters:
  • bragg_peaks (BraggVectors) – Input Bragg vectors.

  • coef_index (list of ints) – List of ints that act as pointers to unit cell parameters and angles to update.

  • coef_update (list of bool) – List of booleans to indicate whether or not to update the cell at that position

  • bragg_k_power (float) – Input Bragg peak intensities are multiplied by k**bragg_k_power to change the weighting of longer scattering vectors

  • bragg_intensity_power (float) – Input Bragg peak intensities are raised power **bragg_intensity_power.

  • k_min (float) – min k value for fitting range (Å^-1)

  • k_max (float) – max k value for fitting range (Å^-1)

  • k_step (float) – step size of k in fitting range (Å^-1)

  • k_broadening (float) – Initial guess for Gaussian broadening of simulated pattern (Å^-1)

  • fit_all_intensities (bool) – Set to true to allow all peak intensities to change independently False forces a single intensity scaling.

  • verbose (bool) – Output the calibrated pixel size.

  • plot_result (bool) – Plot the resulting fit.

  • figsize (list, tuple, np.ndarray) –

  • returnfig (bool) – Return handles figure and axis

Returns:

Optional figure and axis handles, if returnfig=True.

Return type:

fig, ax (handles)

Details: User has the option to define what is allowed to update in the unit cell using the arguments coef_index and coef_update. Each has 6 entries, corresponding to the a, b, c, alpha, beta, gamma parameters of the unit cell, in this order. The coef_update argument is a list of bools specifying whether or not the unit cell value will be allowed to change (True) or must maintain the original value (False) upon fitting. The coef_index argument provides a pointer to the index in which the code will update to.

For example, to update a, b, c, alpha, beta, gamma all independently of eachother, the following arguments should be used:

coef_index = [0, 1, 2, 3, 4, 5] coef_update = [True, True, True, True, True, True,]

The default is set to automatically define what can update in a unit cell based on the point group constraints. When either ‘coef_index’ or ‘coef_update’ are None, these constraints will be automatically pulled from the pointgroup.

For example, the default for cubic unit cells is:

coef_index = [0, 0, 0, 3, 3, 3] coef_update = [True, True, True, False, False, False]

Which allows a, b, and c to update (True in first 3 indices of coef_update) but b and c update based on the value of a (0 in the 1 and 2 list entries in coef_index) such that a = b = c. While coef_update is False for alpha, beta, and gamma (entries 3, 4, 5), no updates will be made to the angles.

The user has the option to predefine coef_index or coef_update to override defaults. In the coef_update list, there must be 6 entries and each are boolean. In the coef_index list, there must be 6 entries, with the first 3 entries being between 0 - 2 and the last 3 entries between 3 - 5. These act as pointers to pull the updated parameter from.

generate_dynamical_diffraction_pattern(beams: PointList, thickness: float | list | tuple | ndarray, zone_axis_lattice: ndarray | None = None, zone_axis_cartesian: ndarray | None = None, foil_normal_lattice: ndarray | None = None, foil_normal_cartesian: ndarray | None = None, verbose: bool = False, always_return_list: bool = False, dynamical_matrix_cache: DynamicalMatrixCache | None = None, return_complex: bool = False, return_eigenvectors: bool = False, return_Smatrix: bool = False) PointList | List[PointList]

Generate a dynamical diffraction pattern (or thickness series of patterns) using the Bloch wave method.

The beams to be included in the Bloch calculation must be pre-calculated and passed as a PointList containing at least (qx, qy, h, k, l) fields.

If thickness is a single value, one new PointList will be returned. If thickness is a sequence of values, a list of PointLists will be returned,

corresponding to each thickness value in the input.

Frequent reference will be made to “Introduction to conventional transmission electron microscopy”

by DeGraef, whose overall approach we follow here.

Parameters:
  • beams (PointList) – PointList from the kinematical diffraction generator which will define the beams included in the Bloch calculation

  • thickness (float or list/array) – The main Bloch calculation can be reused for multiple thicknesses without much overhead.

  • direction. (zone_axis & foil_normal Incident beam orientation and foil normal) – Each can be specified in the Cartesian or crystallographic basis, using e.g. zone_axis_lattice or zone_axis_cartesian. These are internally parsed by Crystal.parse_orientation

Less commonly used args:
always_return_list (bool): When True, the return is always a list of PointLists,

even for a single thickness

dynamical_matrix_cache: (DyanmicalMatrixCache) Dataclass used for caching of the

dynamical matrix. If the cached matrix does not exist, it is computed and stored. Subsequent calls will use the cached matrix for the off-diagonal components of the A matrix and overwrite the diagonal elements. This is used for CBED calculations.

return_complex (bool): When True, returns both the complex amplitude and intensity. Defaults to (False)

Returns:

Bragg peaks with fields [qx, qy, intensity, h, k, l]

or

[bragg_peaks,…] (PointList): If thickness is a list/array, or always_return_list is True,

a list of PointLists is returned.

if return_complex = True:
bragg_peaks (PointList): Bragg peaks with fields [qx, qy, intensity, amplitude, h, k, l]

or

[bragg_peaks,…] (PointList): If thickness is a list/array, or always_return_list is True,

a list of PointLists is returned.

if return_Smatrix = True:
[S_matrix, …], psi_0: Returns a list of S-matrices for each thickness (this is always a list),

and the vector representing the incident plane wave. The beams of the S-matrix have the same order as in the input beams.

Return type:

bragg_peaks (PointList)

generate_CBED(beams: ~emdfile.classes.pointlist.PointList, thickness: float | list | tuple | ~numpy.ndarray, alpha_mrad: float, pixel_size_inv_A: float, DP_size_inv_A: float | None = None, zone_axis_lattice: ~numpy.ndarray | None = None, zone_axis_cartesian: ~numpy.ndarray | None = None, foil_normal_lattice: ~numpy.ndarray | None = None, foil_normal_cartesian: ~numpy.ndarray | None = None, LACBED: bool = False, dtype: ~numpy.dtype = <class 'numpy.float32'>, verbose: bool = False, progress_bar: bool = True, return_mask: bool = False, two_beam_zone_axis_lattice: ~numpy.ndarray | None = None, return_probe: bool = False) ndarray | List[ndarray] | Dict[Tuple[int], ndarray]

Generate a dynamical CBED pattern using the Bloch wave method.

Parameters:
  • beams (PointList) – PointList from the kinematical diffraction generator which will define the beams included in the Bloch calculation

  • thickness (float or list/array) – The main Bloch calculation can be reused for multiple thicknesses without much overhead.

  • alpha_mrad (float) – Convergence angle for CBED pattern. Note that if disks in the calculation overlap, they will be added incoherently, and the resulting CBED will thus represent the average over the unit cell (i.e. a PACBED pattern, as described in LeBeau et al., Ultramicroscopy 110(2): 2010.)

  • pixel_size_inv_A (float) – CBED pixel size in 1/Å.

  • DP_size_inv_A (optional float) – If specified, defines the extents of the diffraction pattern. If left unspecified, the DP will be automatically scaled to fit all of the beams present in the input plus some small buffer.

  • zone_axis (np float vector) – 3 element projection direction for sim pattern Can also be a 3x3 orientation matrix (zone axis 3rd column)

  • foil_normal – 3 element foil normal - set to None to use zone_axis

  • LACBED (bool) – keyed by tuples of (h,k,l).

  • proj_x_axis (np float vector) – 3 element vector defining image x axis (vertical)

  • PointList (two_beam_zone_axis_lattice When only two beams are present in the "beams") – the computation of the projected crystallographic directions becomes ambiguous. In this case, you must specify the indices of the zone axis used to generate the beams.

:paramthe computation of the projected crystallographic directions

becomes ambiguous. In this case, you must specify the indices of the zone axis used to generate the beams.

Parameters:

return_probe (bool) – If True, the probe (np.ndarray) will be returned in additon to the CBED

Returns:

CBED pattern as np.ndarray If thickness is a sequence: CBED patterns for each thickness value as a list of np.ndarrays If LACBED is True and thickness is scalar: Dictionary with tuples of ints (h,k,l) as keys, mapping to np.ndarray. If LACBED is True and thickness is a sequence: List of dictionaries, structured as above. If return_probe is True: will return a tuple (<CBED/LACBED object>, Probe)

Return type:

If thickness is a scalar

calculate_dynamical_structure_factors(accelerating_voltage: float, method: str = 'WK-CP', k_max: float = 2.0, thermal_sigma: float | dict | None = None, tol_structure_factor: float = 0.0, recompute_kinematic_structure_factors=True, g_vec_precision=None)

Calculate and store the relativistic corrected structure factors used for Bloch computations in a dictionary for faster lookup.

Parameters:
  • accelerating_voltage (float) – accelerating voltage in eV

  • method (str) –

    Choose which parameterization of the structure factors to use: “Lobato”: Uses the kinematic structure factors from crystal.py, using the parameterization from

    Lobato & Van Dyck, Acta Cryst A 70:6 (2014)

    ”Lobato-absorptive”: Lobato factors plus an imaginary part

    equal to 0.1•f, as a simple but inaccurate way to include absorption, per Hashimoto, Howie, & Whelan, Proc R Soc Lond A 269:80-103 (1962)

    ”WK”: Uses the Weickenmeier-Kohl parameterization for

    the elastic form factors, including Debye-Waller factor, with no absorption, as described in Weickenmeier & Kohl, Acta Cryst A 47:5 (1991)

    ”WK-C”: WK form factors plus the “core” contribution to absorption

    following H. Rose, Optik 45:2 (1976)

    ”WK-P”: WK form factors plus the phonon/TDS absorptive contribution “WK-CP”: WK form factors plus core and phonon absorption (default)

  • k_max (float) – max scattering length to compute structure factors to. Setting this to 2x the k_max used in generating the beamsn included in a simulation will retain all possible couplings

  • thermal_sigma (float or dict{int->float}) – RMS atomic diplacement for attenuating form factors to account for thermal broadening of the potential, only used when a “WK” method is selected. Required when WK-P or WK-CP are selected. Units are Å. (This is often written as 〈u〉in papers) To specify different 〈u〉 for each element, pass a dictionary with Z as the key, mapping to the appropriate float value

  • tol_structure_factor (float) – tolerance for removing low-valued structure factors. Reflections with structure factor below the tolerance will have zero coupling in the dynamical calculations (i.e. they are the ignored weak beams)

  • recompute_kinematic_structure_factors (bool) – When True, recomputes the kinematic structure factors using the same tol_structure_factor, and with k_max set to half the k_max for the dynamical factors. The factor of half ensures that every beam in a simulation can couple to every other beam (no high-angle couplings in the Bloch matrix are set to zero.)

  • g_vec_precision (optional int) – If specified, rounds |g| to this many decimal places so that automatic caching of the atomic form factors is not slowed down due to floating point errors. Setting this to 3 can give substantial speedup at the cost of some reduced accuracy

  • factors. (See WK_scattering_factors.py for details on the Weickenmeier-Kohl form) –

__init__(positions, numbers, cell, occupancy=None)
Parameters:
  • positions (np.array) – fractional coordinates of each atom in the cell

  • numbers (np.array) – Z number for each atom in the cell, if one number passed it is used for all atom positions

  • cell (np.array) – specify the unit cell, using a variable number of parameters 1 number: the lattice parameter for a cubic cell 3 numbers: the three lattice parameters for an orthorhombic cell 6 numbers: the a,b,c lattice parameters and ɑ,β,ɣ angles for any cell 3x3 array: row vectors containing the (u,v,w) lattice vectors.

  • occupancy (np.array) – Partial occupancy values for each atomic site. Must match the length of positions

positions

fractional atomic coordinates

get_strained_crystal(exx=0.0, eyy=0.0, ezz=0.0, exy=0.0, exz=0.0, eyz=0.0, deformation_matrix=None, return_deformation_matrix=False)

This method returns new Crystal class with strain applied. The directions of (x,y,z) are with respect to the default Crystal orientation, which can be checked with print(Crystal.lat_real) applied to the original Crystal.

Strains are given in fractional values, so exx = 0.01 is 1% strain along the x direction. Deformation matrix should be of the form:

deformation_matrix = np.array([

[1.0+exx, 1.0*exy, 1.0*exz], [1.0*exy, 1.0+eyy, 1.0*eyz], [1.0*exz, 1.0*eyz, 1.0+ezz],

])

Parameters:
  • (float) (eyz) – fractional strain along the xx direction

  • (float) – fractional strain along the yy direction

  • (float) – fractional strain along the zz direction

  • (float) – fractional strain along the xy direction

  • (float) – fractional strain along the xz direction

  • (float) – fractional strain along the yz direction

  • (np.ndarray) (deformation_matrix) – 3x3 array describing deformation matrix

  • (bool) (return_deformation_matrix) – boolean switch to return deformation matrix

Returns:

  • return_deformation_matrix == False – strained_crystal (py4DSTEM.Crystal)

  • return_deformation_matrix == True – (strained_crystal, deformation_matrix)

static from_ase(atoms)

Create a py4DSTEM Crystal object from an ASE atoms object

Parameters:

atoms (ase.Atoms) – an ASE atoms object

static from_prismatic(filepath)

Create a py4DSTEM Crystal object from an prismatic style xyz co-ordinate file

Parameters:

filepath (str|Pathlib.Path) – path to the prismatic format xyz file

static from_CIF(CIF, primitive: bool = True, conventional_standard_structure: bool = True)

Create a Crystal object from a CIF file, using pymatgen to import the CIF

Note that pymatgen typically prefers to return primitive unit cells, which can be overridden by setting conventional_standard_structure=True.

Parameters:
  • CIF – (str or Path) path to the CIF File

  • conventional_standard_structure – (bool) if True, conventional standard unit cell will be returned instead of the primitive unit cell pymatgen typically returns

static from_pymatgen_structure(structure=None, formula=None, space_grp=None, MP_key=None, conventional_standard_structure=True)

Create a Crystal object from a pymatgen Structure object. If a Materials Project API key is installed, you may pass the Materials Project ID of a structure, which will be fetched through the MP API. For setup information see: https://pymatgen.org/usage.html#setting-the-pmg-mapi-key-in-the-config-file. Alternatively, Materials Porject API key can be pass as an argument through the function (MP_key). To get your API key, please visit Materials Project website and login/sign up using your email id. Once logged in, go to the dashboard to generate your own API key (https://materialsproject.org/dashboard).

Note that pymatgen typically prefers to return primitive unit cells, which can be overridden by setting conventional_standard_structure=True.

Parameters:
  • structure – (pymatgen Structure or str), if specified as a string, it will be considered as a Materials Project ID of a structure, otherwise it will accept only pymatgen Structure object. if None, MP database will be queried using the specified formula and/or space groups for the available structure

  • formula – (str), pretty formula to search in the MP database, (note that the forumlas in MP database are not always formatted in the conventional order. Please visit Materials Project website for information (https://materialsproject.org/) if None, structure argument must not be None

  • space_grp – (int) space group number of the forumula provided to query MP database. If None, MP will search for all the available space groups for the formula provided and will consider the one with lowest unit cell volume, only specify when using formula to search MP database

  • MP_key – (str) Materials Project API key

  • conventional_standard_structure – (bool) if True, conventional standard unit cell will be returned instead of the primitive unit cell pymatgen returns

static from_unitcell_parameters(latt_params, elements, positions, space_group=None, lattice_type='cubic', from_cartesian=False, conventional_standard_structure=True)

Create a Crystal using pymatgen to generate unit cell manually from user inputs

Parameters:
  • latt_params – (list of floats) list of lattice parameters. For example, for cubic: latt_params = [a], for hexagonal: latt_params = [a, c], for monoclinic: latt_params = [a,b,c,beta], and in general: latt_params = [a,b,c,alpha,beta,gamma]

  • elements – (list of strings) list of elements, for example for SnS: elements = [“Sn”, “S”]

  • positions – (list) list of (x,y,z) positions for each element present in the elements, default: fractional coord

  • space_group – (optional) (string or int) space group of the crystal system, if specified, unit cell will be created using pymatgen Structure.from_spacegroup function

  • lattice_type – (string) type of crystal family: cubic, hexagonal, triclinic etc; default: ‘cubic’

  • from_cartesian – (bool) if True, positions will be considered as cartesian, default: False

  • conventional_standard_structure – (bool) if True, conventional standard unit cell will be returned instead of the primitive unit cell pymatgen returns

Returns:

Crystal object

setup_diffraction(accelerating_voltage: float)

Set up attributes used for diffraction calculations without going through the full ACOM pipeline.

calculate_structure_factors(k_max: float = 2.0, tol_structure_factor: float = 0.0001, return_intensities: bool = False)

Calculate structure factors for all hkl indices up to max scattering vector k_max

Parameters:
  • k_max (float) – max scattering vector to include (1/Angstroms)

  • tol_structure_factor (float) – tolerance for removing low-valued structure factors

  • return_intensities (bool) – return the intensities and positions of all structure factor peaks.

Returns:

Tuple of the q vectors and intensities of each structure factor.

Return type:

(q_SF, I_SF)

generate_diffraction_pattern(orientation: Orientation | None = None, ind_orientation: int | None = 0, orientation_matrix: ndarray | None = None, zone_axis_lattice: ndarray | None = None, proj_x_lattice: ndarray | None = None, foil_normal_lattice: list | tuple | ndarray | None = None, zone_axis_cartesian: ndarray | None = None, proj_x_cartesian: ndarray | None = None, foil_normal_cartesian: list | tuple | ndarray | None = None, sigma_excitation_error: float = 0.02, tol_excitation_error_mult: float = 3, tol_intensity: float = 0.0001, k_max: float | None = None, keep_qz=False, return_orientation_matrix=False)

Generate a single diffraction pattern, return all peaks as a pointlist.

Parameters:
  • orientation (Orientation) – an Orientation class object

  • orientations (ind_orientation If input is an Orientation class object with multiple) – this input can be used to select a specific orientation.

:param : this input can be used to select a specific orientation. :param orientation_matrix: (3,3) orientation matrix, where columns represent projection directions. :type orientation_matrix: array :param zone_axis_lattice: (3,) projection direction in lattice indices :type zone_axis_lattice: array :param proj_x_lattice: (3,) x-axis direction in lattice indices :type proj_x_lattice: array :param zone_axis_cartesian: (3,) cartesian projection direction :type zone_axis_cartesian: array :param proj_x_cartesian: (3,) cartesian projection direction :type proj_x_cartesian: array :param foil_normal: 3 element foil normal - set to None to use zone_axis :param proj_x_axis: 3 element vector defining image x axis (vertical) :type proj_x_axis: np float vector :param accel_voltage: Accelerating voltage in Volts. If not specified,

we check to see if crystal already has voltage specified.

Parameters:
  • sigma_excitation_error (float) – sigma value for envelope applied to s_g (excitation errors) in units of inverse Angstroms

  • tol_excitation_error_mult (float) – tolerance in units of sigma for s_g inclusion

  • tol_intensity (np float) – tolerance in intensity units for inclusion of diffraction spots

  • k_max (float) – Maximum scattering vector

  • keep_qz (bool) – Flag to return out-of-plane diffraction vectors

  • return_orientation_matrix (bool) – Return the orientation matrix

Returns:

list of all Bragg peaks with fields [qx, qy, intensity, h, k, l] orientation_matrix (array): 3x3 orientation matrix (optional)

Return type:

bragg_peaks (PointList)

generate_ring_pattern(k_max=2.0, use_bloch=False, thickness=None, bloch_params=None, orientation_plan_params=None, sigma_excitation_error=0.02, tol_intensity=0.001, plot_rings=True, plot_params={}, return_calc=True)

Calculate polycrystalline diffraction pattern from structure

Parameters:
  • k_max (float) – Maximum scattering vector

  • use_bloch (bool) – if true, use dynamic instead of kinematic approach

  • thickness (float) – thickness in Ångström to evaluate diffraction patterns, only needed for dynamical calculations

  • bloch_params (dict) – optional, parameters to calculate dynamical structure factor, see calculate_dynamical_structure_factors doc strings

  • orientation_plan_params (dict) – optional, parameters to calculate orientation plan, see orientation_plan doc strings

  • sigma_excitation_error (float) – sigma value for envelope applied to s_g (excitation errors) in units of inverse Angstroms

  • tol_intensity (np float) – tolerance in intensity units for inclusion of diffraction spots

  • plot_rings (bool) – if true, plot diffraction rings with plot_ring_pattern

  • return_calc (bool) – return radii and intensities

Returns:

radii of ring pattern in units of scattering vector k intensity_unique (np array): intensity of rings weighted by frequency of diffraciton spots

Return type:

radii_unique (np array)

excitation_errors(g, foil_normal=None)

Calculate the excitation errors, assuming k0 = [0, 0, -1/lambda]. If foil normal is not specified, we assume it is [0,0,-1].

calculate_bragg_peak_histogram(bragg_peaks, bragg_k_power=1.0, bragg_intensity_power=1.0, k_min=0.0, k_max=None, k_step=0.005)

Prepare experimental bragg peaks for lattice parameter or unit cell fitting.

Parameters:
  • bragg_peaks (BraggVectors) – Input Bragg vectors.

  • bragg_k_power (float) – Input Bragg peak intensities are multiplied by k**bragg_k_power to change the weighting of longer scattering vectors

  • bragg_intensity_power (float) – Input Bragg peak intensities are raised power **bragg_intensity_power.

  • k_min (float) – min k value for fitting range (Å^-1)

  • k_max (float) – max k value for fitting range (Å^-1)

  • k_step (float) – step size of k in fitting range (Å^-1)

Returns:

Bragg vectors after calibration fig, ax (handles): Optional figure and axis handles, if returnfig=True.

Return type:

bragg_peaks_cali (BraggVectors)

py4DSTEM.process.diffraction.crystal.generate_moire_diffraction_pattern(bragg_peaks_0, bragg_peaks_1, thresh_0=0.0002, thresh_1=0.0002, exx_1=0.0, eyy_1=0.0, exy_1=0.0, phi_1=0.0, power=2.0)

Calculate a Moire lattice from 2 parent diffraction patterns. The second lattice can be rotated and strained with respect to the original lattice. Note that this strain is applied in real space, and so the inverse of the calculated infinitestimal strain tensor is applied.

Parameters:
  • bragg_peaks_0 (BraggVector) – Bragg vectors for parent lattice 0.

  • bragg_peaks_1 (BraggVector) – Bragg vectors for parent lattice 1.

  • thresh_0 (float) – Intensity threshold for structure factors from lattice 0.

  • thresh_1 (float) – Intensity threshold for structure factors from lattice 1.

  • exx_1 (float) – Strain of lattice 1 in x direction (vertical) in real space.

  • eyy_1 (float) – Strain of lattice 1 in y direction (horizontal) in real space.

  • exy_1 (float) – Shear strain of lattice 1 in (x,y) direction (diagonal) in real space.

  • phi_1 (float) – Rotation of lattice 1 in real space.

  • power (float) – Plotting power law (default is amplitude**2.0, i.e. intensity).

Returns:

parent_peaks_0, parent_peaks_1, moire_peaks – Bragg vectors for the rotated & strained parent lattices and the moire lattice

Return type:

BraggVectors

py4DSTEM.process.diffraction.crystal.plot_moire_diffraction_pattern(bragg_parent_0, bragg_parent_1, bragg_moire, int_range=(0, 0.005), k_max=1.0, plot_subpixel=True, labels=None, marker_size_parent=16, marker_size_moire=4, text_size_parent=10, text_size_moire=6, add_labels_parent=False, add_labels_moire=False, dist_labels=0.03, dist_check=0.06, sep_labels=0.03, figsize=(8, 6), returnfig=False)

Plot Moire lattice and parent lattices.

Parameters:
  • bragg_peaks_0 (BraggVector) – Bragg vectors for parent lattice 0.

  • bragg_peaks_1 (BraggVector) – Bragg vectors for parent lattice 1.

  • bragg_moire (BraggVector) – Bragg vectors for moire lattice.

  • int_range ((float, float)) – Plotting intensity range for the Moire peaks.

  • k_max (float) – Max k value of the plotted Moire lattice.

  • plot_subpixel (bool) – Apply subpixel corrections to the Bragg spot positions. Matplotlib default scatter plot rounds to the nearest pixel.

  • labels (list) – List of text labels for parent lattices

  • marker_size_parent (float) – Size of plot markers for the two parent lattices.

  • marker_size_moire (float) – Size of plot markers for the Moire lattice.

  • text_size_parent (float) – Label text size for parent lattice.

  • text_size_moire (float) – Label text size for Moire lattice.

  • add_labels_parent (bool) – Plot the parent lattice index labels.

  • add_labels_moire (bool) – Plot the parent lattice index labels for the Moire spots.

  • dist_labels (float) – Distance to move the labels off the spots.

  • dist_check (float) – Set to some distance to “push” the labels away from each other if they are within this distance.

  • sep_labels (float) – Separation distance for labels which are “pushed” apart.

  • figsize ((float,float)) – Size of output figure.

  • returnfig (bool) – Return the (fix,ax) handles of the plot.

Returns:

fig, ax – Figure and axes handles for the moire plot.

Return type:

matplotlib handles (optional)

py4DSTEM.process.diffraction.crystal_ACOM.orientation_plan(self, zone_axis_range: ndarray = array([[0, 1, 1], [1, 1, 1]]), angle_step_zone_axis: float = 2.0, angle_coarse_zone_axis: float | None = None, angle_refine_range: float | None = None, angle_step_in_plane: float = 2.0, accel_voltage: float = 300000.0, corr_kernel_size: float = 0.08, radial_power: float = 1.0, intensity_power: float = 0.25, calculate_correlation_array=True, tol_peak_delete=None, tol_distance: float = 0.01, fiber_axis=None, fiber_angles=None, figsize: list | tuple | ndarray = (6, 6), CUDA: bool = False, progress_bar: bool = True)

Calculate the rotation basis arrays for an SO(3) rotation correlogram.

Parameters:
  • zone_axis_range (float) –

    Row vectors give the range for zone axis orientations. If user specifies 2 vectors (2x3 array), we start at [0,0,1]

    to make z-x-z rotation work.

    If user specifies 3 vectors (3x3 array), plan will span these vectors. Setting to ‘full’ as a string will use a hemispherical range. Setting to ‘half’ as a string will use a quarter sphere range. Setting to ‘fiber’ as a string will make a spherical cap around a given vector. Setting to ‘auto’ will use pymatgen to determine the point group symmetry

    of the structure and choose an appropriate zone_axis_range

  • angle_step_zone_axis (float) – Approximate angular step size for zone axis search [degrees]

  • angle_coarse_zone_axis (float) – Coarse step size for zone axis search [degrees]. Setting to None uses the same value as angle_step_zone_axis.

  • angle_refine_range (float) – Range of angles to use for zone axis refinement. Setting to None uses same value as angle_coarse_zone_axis.

  • angle_step_in_plane (float) – Approximate angular step size for in-plane rotation [degrees]

  • accel_voltage (float) – Accelerating voltage for electrons [Volts]

  • corr_kernel_size (float) – Correlation kernel size length. The size of the overlap kernel between the measured Bragg peaks and diffraction library Bragg peaks. [1/Angstroms]

  • radial_power (float) – Power for scaling the correlation intensity as a function of the peak radius

  • intensity_power (float) – Power for scaling the correlation intensity as a function of the peak intensity

  • calculate_correlation_array (bool) – Set to false to skip calculating the correlation array. This is useful when we only want the angular range / rotation matrices.

  • tol_peak_delete (float) – Distance to delete peaks for multiple matches. Default is kernel_size * 0.5

  • tol_distance (float) – Distance tolerance for radial shell assignment [1/Angstroms]

  • fiber_axis (float) – (3,) vector specifying the fiber axis

  • fiber_angles (float) – (2,) vector specifying angle range from fiber axis, and in-plane angular range [degrees]

  • cartesian_directions (bool) – When set to true, all zone axes and projection directions are specified in Cartesian directions.

  • figsize (float) – (2,) vector giving the figure size

  • CUDA (bool) – Use CUDA for the Fourier operations.

  • progress_bar (bool) – If false no progress bar is displayed

py4DSTEM.process.diffraction.crystal_ACOM.match_orientations(self, bragg_peaks_array: PointListArray, num_matches_return: int = 1, min_angle_between_matches_deg=None, min_number_peaks: int = 3, inversion_symmetry: bool = True, multiple_corr_reset: bool = True, return_orientation: bool = True, progress_bar: bool = True)
Parameters:
  • bragg_peaks_array (PointListArray) – PointListArray containing the Bragg peaks and intensities, with calibrations applied

  • num_matches_return (int) – return these many matches as 3th dim of orient (matrix)

  • min_angle_between_matches_deg (int) – Minimum angle between zone axis of multiple matches, in degrees. Note that I haven’t thought how to handle in-plane rotations, since multiple matches are possible.

  • min_number_peaks (int) – Minimum number of peaks required to perform ACOM matching

  • inversion_symmetry (bool) – check for inversion symmetry in the matches

  • multiple_corr_reset (bool) – keep original correlation score for multiple matches

  • return_orientation (bool) – Return orientation map from function for inspection. The map is always stored in the Crystal object.

  • progress_bar (bool) – Show or hide the progress bar

py4DSTEM.process.diffraction.crystal_ACOM.match_single_pattern(self, bragg_peaks: PointList, num_matches_return: int = 1, min_angle_between_matches_deg=None, min_number_peaks=3, inversion_symmetry=True, multiple_corr_reset=True, plot_polar: bool = False, plot_corr: bool = False, returnfig: bool = False, figsize: list | tuple | ndarray = (12, 4), verbose: bool = False)

Solve for the best fit orientation of a single diffraction pattern.

Parameters:
  • bragg_peaks (PointList) – numpy array containing the Bragg positions and intensities (‘qx’, ‘qy’, ‘intensity’)

  • num_matches_return (int) – return these many matches as 3th dim of orient (matrix)

  • min_angle_between_matches_deg (int) – Minimum angle between zone axis of multiple matches, in degrees. Note that I haven’t thought how to handle in-plane rotations, since multiple matches are possible.

  • min_number_peaks (int) – Minimum number of peaks required to perform ACOM matching

  • bool (multiple_corr_reset) – check for inversion symmetry in the matches

  • bool – keep original correlation score for multiple matches

  • subpixel_tilt (bool) – set to false for faster matching, returning the nearest corr point

  • plot_polar (bool) – set to true to plot the polar transform of the diffraction pattern

  • plot_corr (bool) – set to true to plot the resulting correlogram

  • returnfig (bool) – return figure handles

  • figsize (list) – size of figure

  • verbose (bool) – Print the fitted zone axes, correlation scores

  • CUDA (bool) – Enable CUDA for the FFT steps

Returns:

  • orientation (Orientation) – Orientation class containing all outputs

  • fig, ax (handles) – Figure handles for the plotting output

py4DSTEM.process.diffraction.crystal_ACOM.cluster_grains(self, threshold_add=1.0, threshold_grow=0.1, angle_tolerance_deg=5.0, progress_bar=True)

Cluster grains using rotation criterion, and correlation values.

Parameters:
  • threshold_add (float) – Minimum signal required for a probe position to initialize a cluster.

  • threshold_grow (float) – Minimum signal required for a probe position to be added to a cluster.

  • angle_tolerance_deg (float) – Rotation rolerance for clustering grains.

  • progress_bar (bool) – Turns on the progress bar for the polar transformation

py4DSTEM.process.diffraction.crystal_ACOM.cluster_orientation_map(self, stripe_width=(2, 2), area_min=2)

Produce a new orientation map from the clustered grains. Use a stripe pattern for the overlapping grains.

Parameters:
  • stripe_width ((int,int)) – Width of stripes for plotting maps with overlapping grains

  • area_min ((int)) – Minimum size of grains to include

Returns:

The clustered orientation map

Return type:

orientation_map

py4DSTEM.process.diffraction.crystal_ACOM.calculate_strain(self, bragg_peaks_array: PointListArray, orientation_map: OrientationMap, corr_kernel_size=None, sigma_excitation_error=0.02, tol_excitation_error_mult: float = 3, tol_intensity: float = 0.0001, k_max: float | None = None, min_num_peaks=5, rotation_range=None, mask_from_corr=True, corr_range=(0, 2), corr_normalize=True, progress_bar=True)

This function takes in both a PointListArray containing Bragg peaks, and a corresponding OrientationMap, and uses least squares to compute the deformation tensor which transforms the simulated diffraction pattern into the experimental pattern, for all probe positons.

TODO: add robust fitting?

Parameters:
  • bragg_peaks_array (PointListArray) – All Bragg peaks

  • orientation_map (OrientationMap) – Orientation map generated from ACOM

  • corr_kernel_size (float) – Correlation kernel size - if user does not specify, uses self.corr_kernel_size.

  • sigma_excitation_error (float) – sigma value for envelope applied to s_g (excitation errors) in units of inverse Angstroms

  • tol_excitation_error_mult (float) – tolerance in units of sigma for s_g inclusion

  • tol_intensity (np float) – tolerance in intensity units for inclusion of diffraction spots

  • k_max (float) – Maximum scattering vector

  • min_num_peaks (int) – Minimum number of peaks required.

  • rotation_range (float) – Maximum rotation range in radians (for symmetry reduction).

  • progress_bar (bool) – Show progress bar

  • mask_from_corr (bool) – Use ACOM correlation signal for mask

  • corr_range (np.ndarray) – Range of correlation signals for mask

  • corr_normalize (bool) – Normalize correlation signal before masking

Returns:

strain tensor

Return type:

strain_map (RealSlice)

py4DSTEM.process.diffraction.crystal_ACOM.save_ang_file(self, file_name, orientation_map, ind_orientation=0, pixel_size=1.0, pixel_units='px', transpose_xy=True, flip_x=False)

This function outputs an ascii text file in the .ang format, containing the Euler angles of an orientation map.

Parameters:
  • file_name (str) – Path to save .ang file.

  • orientation_map (OrientationMap) – Class containing orientation matrices, correlation values, etc.

  • ind_orientation (int) – Which orientation match to plot if num_matches > 1

  • pixel_size (float) – Pixel size, if known.

  • pixel_units (str) – Units of the pixel size

  • transpose_xy (bool) – Transpose x and y pixel coordinates.

  • flip_x (bool) – Swap x direction pixels (after transpose).

Returns:

nothing

py4DSTEM.process.diffraction.crystal_ACOM.symmetry_reduce_directions(self, orientation, match_ind=0, plot_output=False, figsize=(15, 6), el_shift=0.0, az_shift=-30.0)

This function calculates the symmetry-reduced cartesian directions from and orientation matrix stored in orientation.matrix, and outputs them into orientation.family. It optionally plots the 3D output.

class py4DSTEM.process.diffraction.crystal_bloch.DynamicalMatrixCache(has_valid_cache: bool = False, cached_U_gmh: <built-in function array> = None)
__init__(has_valid_cache: bool = False, cached_U_gmh: array | None = None) None
py4DSTEM.process.diffraction.crystal_bloch.calculate_dynamical_structure_factors(self, accelerating_voltage: float, method: str = 'WK-CP', k_max: float = 2.0, thermal_sigma: float | dict | None = None, tol_structure_factor: float = 0.0, recompute_kinematic_structure_factors=True, g_vec_precision=None)

Calculate and store the relativistic corrected structure factors used for Bloch computations in a dictionary for faster lookup.

Parameters:
  • accelerating_voltage (float) – accelerating voltage in eV

  • method (str) –

    Choose which parameterization of the structure factors to use: “Lobato”: Uses the kinematic structure factors from crystal.py, using the parameterization from

    Lobato & Van Dyck, Acta Cryst A 70:6 (2014)

    ”Lobato-absorptive”: Lobato factors plus an imaginary part

    equal to 0.1•f, as a simple but inaccurate way to include absorption, per Hashimoto, Howie, & Whelan, Proc R Soc Lond A 269:80-103 (1962)

    ”WK”: Uses the Weickenmeier-Kohl parameterization for

    the elastic form factors, including Debye-Waller factor, with no absorption, as described in Weickenmeier & Kohl, Acta Cryst A 47:5 (1991)

    ”WK-C”: WK form factors plus the “core” contribution to absorption

    following H. Rose, Optik 45:2 (1976)

    ”WK-P”: WK form factors plus the phonon/TDS absorptive contribution “WK-CP”: WK form factors plus core and phonon absorption (default)

  • k_max (float) – max scattering length to compute structure factors to. Setting this to 2x the k_max used in generating the beamsn included in a simulation will retain all possible couplings

  • thermal_sigma (float or dict{int->float}) – RMS atomic diplacement for attenuating form factors to account for thermal broadening of the potential, only used when a “WK” method is selected. Required when WK-P or WK-CP are selected. Units are Å. (This is often written as 〈u〉in papers) To specify different 〈u〉 for each element, pass a dictionary with Z as the key, mapping to the appropriate float value

  • tol_structure_factor (float) – tolerance for removing low-valued structure factors. Reflections with structure factor below the tolerance will have zero coupling in the dynamical calculations (i.e. they are the ignored weak beams)

  • recompute_kinematic_structure_factors (bool) – When True, recomputes the kinematic structure factors using the same tol_structure_factor, and with k_max set to half the k_max for the dynamical factors. The factor of half ensures that every beam in a simulation can couple to every other beam (no high-angle couplings in the Bloch matrix are set to zero.)

  • g_vec_precision (optional int) – If specified, rounds |g| to this many decimal places so that automatic caching of the atomic form factors is not slowed down due to floating point errors. Setting this to 3 can give substantial speedup at the cost of some reduced accuracy

  • factors. (See WK_scattering_factors.py for details on the Weickenmeier-Kohl form) –

py4DSTEM.process.diffraction.crystal_bloch.generate_dynamical_diffraction_pattern(self, beams: PointList, thickness: float | list | tuple | ndarray, zone_axis_lattice: ndarray | None = None, zone_axis_cartesian: ndarray | None = None, foil_normal_lattice: ndarray | None = None, foil_normal_cartesian: ndarray | None = None, verbose: bool = False, always_return_list: bool = False, dynamical_matrix_cache: DynamicalMatrixCache | None = None, return_complex: bool = False, return_eigenvectors: bool = False, return_Smatrix: bool = False) PointList | List[PointList]

Generate a dynamical diffraction pattern (or thickness series of patterns) using the Bloch wave method.

The beams to be included in the Bloch calculation must be pre-calculated and passed as a PointList containing at least (qx, qy, h, k, l) fields.

If thickness is a single value, one new PointList will be returned. If thickness is a sequence of values, a list of PointLists will be returned,

corresponding to each thickness value in the input.

Frequent reference will be made to “Introduction to conventional transmission electron microscopy”

by DeGraef, whose overall approach we follow here.

Parameters:
  • beams (PointList) – PointList from the kinematical diffraction generator which will define the beams included in the Bloch calculation

  • thickness (float or list/array) – The main Bloch calculation can be reused for multiple thicknesses without much overhead.

  • direction. (zone_axis & foil_normal Incident beam orientation and foil normal) – Each can be specified in the Cartesian or crystallographic basis, using e.g. zone_axis_lattice or zone_axis_cartesian. These are internally parsed by Crystal.parse_orientation

Less commonly used args:
always_return_list (bool): When True, the return is always a list of PointLists,

even for a single thickness

dynamical_matrix_cache: (DyanmicalMatrixCache) Dataclass used for caching of the

dynamical matrix. If the cached matrix does not exist, it is computed and stored. Subsequent calls will use the cached matrix for the off-diagonal components of the A matrix and overwrite the diagonal elements. This is used for CBED calculations.

return_complex (bool): When True, returns both the complex amplitude and intensity. Defaults to (False)

Returns:

Bragg peaks with fields [qx, qy, intensity, h, k, l]

or

[bragg_peaks,…] (PointList): If thickness is a list/array, or always_return_list is True,

a list of PointLists is returned.

if return_complex = True:
bragg_peaks (PointList): Bragg peaks with fields [qx, qy, intensity, amplitude, h, k, l]

or

[bragg_peaks,…] (PointList): If thickness is a list/array, or always_return_list is True,

a list of PointLists is returned.

if return_Smatrix = True:
[S_matrix, …], psi_0: Returns a list of S-matrices for each thickness (this is always a list),

and the vector representing the incident plane wave. The beams of the S-matrix have the same order as in the input beams.

Return type:

bragg_peaks (PointList)

py4DSTEM.process.diffraction.crystal_bloch.generate_CBED(self, beams: ~emdfile.classes.pointlist.PointList, thickness: float | list | tuple | ~numpy.ndarray, alpha_mrad: float, pixel_size_inv_A: float, DP_size_inv_A: float | None = None, zone_axis_lattice: ~numpy.ndarray | None = None, zone_axis_cartesian: ~numpy.ndarray | None = None, foil_normal_lattice: ~numpy.ndarray | None = None, foil_normal_cartesian: ~numpy.ndarray | None = None, LACBED: bool = False, dtype: ~numpy.dtype = <class 'numpy.float32'>, verbose: bool = False, progress_bar: bool = True, return_mask: bool = False, two_beam_zone_axis_lattice: ~numpy.ndarray | None = None, return_probe: bool = False) ndarray | List[ndarray] | Dict[Tuple[int], ndarray]

Generate a dynamical CBED pattern using the Bloch wave method.

Parameters:
  • beams (PointList) – PointList from the kinematical diffraction generator which will define the beams included in the Bloch calculation

  • thickness (float or list/array) – The main Bloch calculation can be reused for multiple thicknesses without much overhead.

  • alpha_mrad (float) – Convergence angle for CBED pattern. Note that if disks in the calculation overlap, they will be added incoherently, and the resulting CBED will thus represent the average over the unit cell (i.e. a PACBED pattern, as described in LeBeau et al., Ultramicroscopy 110(2): 2010.)

  • pixel_size_inv_A (float) – CBED pixel size in 1/Å.

  • DP_size_inv_A (optional float) – If specified, defines the extents of the diffraction pattern. If left unspecified, the DP will be automatically scaled to fit all of the beams present in the input plus some small buffer.

  • zone_axis (np float vector) – 3 element projection direction for sim pattern Can also be a 3x3 orientation matrix (zone axis 3rd column)

  • foil_normal – 3 element foil normal - set to None to use zone_axis

  • LACBED (bool) – keyed by tuples of (h,k,l).

  • proj_x_axis (np float vector) – 3 element vector defining image x axis (vertical)

  • PointList (two_beam_zone_axis_lattice When only two beams are present in the "beams") – the computation of the projected crystallographic directions becomes ambiguous. In this case, you must specify the indices of the zone axis used to generate the beams.

:paramthe computation of the projected crystallographic directions

becomes ambiguous. In this case, you must specify the indices of the zone axis used to generate the beams.

Parameters:

return_probe (bool) – If True, the probe (np.ndarray) will be returned in additon to the CBED

Returns:

CBED pattern as np.ndarray If thickness is a sequence: CBED patterns for each thickness value as a list of np.ndarrays If LACBED is True and thickness is scalar: Dictionary with tuples of ints (h,k,l) as keys, mapping to np.ndarray. If LACBED is True and thickness is a sequence: List of dictionaries, structured as above. If return_probe is True: will return a tuple (<CBED/LACBED object>, Probe)

Return type:

If thickness is a scalar

py4DSTEM.process.diffraction.crystal_calibrate.calibrate_pixel_size(self, bragg_peaks, scale_pixel_size=1.0, bragg_k_power=1.0, bragg_intensity_power=1.0, k_min=0.0, k_max=None, k_step=0.002, k_broadening=0.002, fit_all_intensities=True, set_calibration_in_place=False, verbose=True, plot_result=False, figsize: list | tuple | ndarray = (12, 6), returnfig=False)

Use the calculated structure factor scattering lengths to compute 1D diffraction patterns, and solve the best-fit relative scaling between them. Returns the fit pixel size in Å^-1.

Parameters:
  • bragg_peaks (BraggVectors) – Input Bragg vectors.

  • scale_pixel_size (float) – Initial guess for scaling of the existing pixel size If the pixel size is currently uncalibrated, this is a guess of the pixel size in Å^-1. If the pixel size is already (approximately) calibrated, this is the scaling factor to correct that existing calibration.

  • bragg_k_power (float) – Input Bragg peak intensities are multiplied by k**bragg_k_power to change the weighting of longer scattering vectors

  • bragg_intensity_power (float) – Input Bragg peak intensities are raised power **bragg_intensity_power.

  • k_min (float) – min k value for fitting range (Å^-1)

  • k_max (float) – max k value for fitting range (Å^-1)

  • k_step (float) step size of k in fitting range (Å^-1) –

  • k_broadening (float) – Initial guess for Gaussian broadening of simulated pattern (Å^-1)

  • fit_all_intensities (bool) – Set to true to allow all peak intensities to change independently False forces a single intensity scaling.

  • set_calibration (bool) – if True, set the fit pixel size to the calibration metadata, and calibrate bragg_peaks

  • verbose (bool) – Output the calibrated pixel size.

  • plot_result (bool) – Plot the resulting fit.

  • figsize (list, tuple, np.ndarray) – Figure size of the plot.

  • returnfig (bool) – Return handles figure and axis

Returns:

fig, ax – Figure and axis handles, if returnfig=True.

Return type:

handles, optional

py4DSTEM.process.diffraction.crystal_calibrate.calibrate_unit_cell(self, bragg_peaks, coef_index=None, coef_update=None, bragg_k_power=1.0, bragg_intensity_power=1.0, k_min=0.0, k_max=None, k_step=0.005, k_broadening=0.02, fit_all_intensities=True, verbose=True, plot_result=False, figsize: list | tuple | ndarray = (12, 6), returnfig=False)

Solve for the best fit scaling between the computed structure factors and bragg_peaks.

Parameters:
  • bragg_peaks (BraggVectors) – Input Bragg vectors.

  • coef_index (list of ints) – List of ints that act as pointers to unit cell parameters and angles to update.

  • coef_update (list of bool) – List of booleans to indicate whether or not to update the cell at that position

  • bragg_k_power (float) – Input Bragg peak intensities are multiplied by k**bragg_k_power to change the weighting of longer scattering vectors

  • bragg_intensity_power (float) – Input Bragg peak intensities are raised power **bragg_intensity_power.

  • k_min (float) – min k value for fitting range (Å^-1)

  • k_max (float) – max k value for fitting range (Å^-1)

  • k_step (float) – step size of k in fitting range (Å^-1)

  • k_broadening (float) – Initial guess for Gaussian broadening of simulated pattern (Å^-1)

  • fit_all_intensities (bool) – Set to true to allow all peak intensities to change independently False forces a single intensity scaling.

  • verbose (bool) – Output the calibrated pixel size.

  • plot_result (bool) – Plot the resulting fit.

  • figsize (list, tuple, np.ndarray) –

  • returnfig (bool) – Return handles figure and axis

Returns:

Optional figure and axis handles, if returnfig=True.

Return type:

fig, ax (handles)

Details: User has the option to define what is allowed to update in the unit cell using the arguments coef_index and coef_update. Each has 6 entries, corresponding to the a, b, c, alpha, beta, gamma parameters of the unit cell, in this order. The coef_update argument is a list of bools specifying whether or not the unit cell value will be allowed to change (True) or must maintain the original value (False) upon fitting. The coef_index argument provides a pointer to the index in which the code will update to.

For example, to update a, b, c, alpha, beta, gamma all independently of eachother, the following arguments should be used:

coef_index = [0, 1, 2, 3, 4, 5] coef_update = [True, True, True, True, True, True,]

The default is set to automatically define what can update in a unit cell based on the point group constraints. When either ‘coef_index’ or ‘coef_update’ are None, these constraints will be automatically pulled from the pointgroup.

For example, the default for cubic unit cells is:

coef_index = [0, 0, 0, 3, 3, 3] coef_update = [True, True, True, False, False, False]

Which allows a, b, and c to update (True in first 3 indices of coef_update) but b and c update based on the value of a (0 in the 1 and 2 list entries in coef_index) such that a = b = c. While coef_update is False for alpha, beta, and gamma (entries 3, 4, 5), no updates will be made to the angles.

The user has the option to predefine coef_index or coef_update to override defaults. In the coef_update list, there must be 6 entries and each are boolean. In the coef_index list, there must be 6 entries, with the first 3 entries being between 0 - 2 and the last 3 entries between 3 - 5. These act as pointers to pull the updated parameter from.

class py4DSTEM.process.diffraction.crystal_phase.Crystal_Phase(crystals, orientation_maps, name)

A class storing multiple crystal structures, and associated diffraction data. Must be initialized after matching orientations to a pointlistarray???

__init__(crystals, orientation_maps, name)
Parameters:
  • crystals (list) – List of crystal instances

  • orientation_maps (list) – List of orientation maps

  • name (str) – Name of Crystal_Phase instance

plot_all_phase_maps(map_scale_values=None, index=0)

Visualize phase maps of dataset.

Parameters:

map_scale_values (float) – Value to scale correlations by

quantify_phase(pointlistarray, tolerance_distance=0.08, method='nnls', intensity_power=0, mask_peaks=None)

Quantification of the phase of a crystal based on the crystal instances and the pointlistarray.

Parameters:
  • pointlisarray (pointlistarray) – Pointlistarray to quantify phase of

  • tolerance_distance (float) – Distance allowed between a peak and match

  • method (str) – Numerical method used to quantify phase

  • intensity_power (float) –

  • mask_peaks (list, optional) – A pointer of which positions to mask peaks from

Details:

quantify_phase_pointlist(pointlistarray, position, method='nnls', tolerance_distance=0.08, intensity_power=0, mask_peaks=None)
Parameters:
  • pointlisarray (pointlistarray) – Pointlistarray to quantify phase of

  • position (tuple/list) – Position of pointlist in pointlistarray

  • tolerance_distance (float) – Distance allowed between a peak and match

  • method (str) – Numerical method used to quantify phase

  • intensity_power (float) –

  • mask_peaks (list, optional) – A pointer of which positions to mask peaks from

Returns:

Peak matches in the rows of array and the crystals in the columns phase_weights (np.ndarray): Weights of each phase phase_residuals (np.ndarray): Residuals crystal_identity (list): List of lists, where the each entry represents the position in the

crystal and orientation match that is associated with the phase weights. for example, if the output was [[0,0], [0,1], [1,0], [0,1]], the first entry [0,0] in phase weights is associated with the first crystal the first match within that crystal. [0,1] is the first crystal and the second match within that crystal.

Return type:

pointlist_peak_intensity_matches (np.ndarray)

py4DSTEM.process.diffraction.crystal_viz.plot_structure(self, orientation_matrix: ndarray | None = None, zone_axis_lattice: ndarray | None = None, proj_x_lattice: ndarray | None = None, zone_axis_cartesian: ndarray | None = None, proj_x_cartesian: ndarray | None = None, size_marker: float = 400, tol_distance: float = 0.001, plot_limit: ndarray | None = None, camera_dist: float | None = None, show_axes: bool = False, perspective_axes: bool = True, figsize: tuple | list | ndarray = (8, 8), returnfig: bool = False)

Quick 3D plot of the untit cell /atomic structure.

Parameters:
  • orientation_matrix (array) – (3,3) orientation matrix, where columns represent projection directions.

  • zone_axis_lattice (array) – (3,) projection direction in lattice indices

  • proj_x_lattice (array) – (3,) x-axis direction in lattice indices

  • zone_axis_cartesian (array) – (3,) cartesian projection direction

  • proj_x_cartesian (array) – (3,) cartesian projection direction

  • scale_markers (float) – Size scaling for markers

  • tol_distance (float) – Tolerance for repeating atoms on edges on cell boundaries.

  • plot_limit (float) – (2,3) numpy array containing x y z plot min and max in columns. Default is 1.1* unit cell dimensions.

  • camera_dist (float) – Move camera closer to the plot (relative to matplotlib default of 10)

  • show_axes (bool) – Whether to plot axes or not.

  • perspective_axes (bool) – Select either perspective (true) or orthogonal (false) axes

  • figsize (2 element float) – Size scaling of figure axes.

  • returnfig (bool) – Return figure and axes handles.

Returns:

fig, ax (optional) figure and axes handles

py4DSTEM.process.diffraction.crystal_viz.plot_structure_factors(self, orientation_matrix: ndarray | None = None, zone_axis_lattice: ndarray | None = None, proj_x_lattice: ndarray | None = None, zone_axis_cartesian: ndarray | None = None, proj_x_cartesian: ndarray | None = None, scale_markers: float = 1000.0, plot_limit: list | tuple | ndarray | None = None, camera_dist: float | None = None, show_axes: bool = True, perspective_axes: bool = True, figsize: list | tuple | ndarray = (8, 8), returnfig: bool = False)

3D scatter plot of the structure factors using magnitude^2, i.e. intensity.

Parameters:
  • orientation_matrix (array) – (3,3) orientation matrix, where columns represent projection directions.

  • zone_axis_lattice (array) – (3,) projection direction in lattice indices

  • proj_x_lattice (array) – (3,) x-axis direction in lattice indices

  • zone_axis_cartesian (array) – (3,) cartesian projection direction

  • proj_x_cartesian (array) – (3,) cartesian projection direction

  • scale_markers (float) – size scaling for markers

  • plot_limit (float) – x y z plot limits, default is [-1 1]*self.k_max

  • camera_dist (float) – Move camera closer to the plot (relative to matplotlib default of 10)

  • show_axes (bool) – Whether to plot axes or not.

  • perspective_axes (bool) – Select either perspective (true) or orthogonal (false) axes

  • figsize (2 element float) – size scaling of figure axes

  • returnfig (bool) – set to True to return figure and axes handles

Returns:

fig, ax (optional) figure and axes handles

py4DSTEM.process.diffraction.crystal_viz.plot_scattering_intensity(self, k_min=0.0, k_max=None, k_step=0.001, k_broadening=0.0, k_power_scale=0.0, int_power_scale=0.5, int_scale=1.0, remove_origin=True, bragg_peaks=None, bragg_k_power=0.0, bragg_intensity_power=1.0, bragg_k_broadening=0.005, figsize: list | tuple | ndarray = (10, 4), returnfig: bool = False)

1D plot of the structure factors

Parameters:
  • k_min (float) – min k value for profile range.

  • k_max (float) – max k value for profile range.

  • k_step (float) – Step size of k in profile range.

  • k_broadening (float) – Broadening of simulated pattern.

  • k_power_scale (float) – Scale SF intensities by k**k_power_scale.

  • int_power_scale (float) – Scale SF intensities**int_power_scale.

  • int_scale (float) – Scale output profile by this value.

  • remove_origin (bool) – Remove origin from plot.

  • bragg_peaks (BraggVectors) – Passed in bragg_peaks for comparison with simulated pattern.

  • bragg_k_power (float) – bragg_peaks scaled by k**bragg_k_power.

  • bragg_intensity_power (float) – bragg_peaks scaled by intensities**bragg_intensity_power.

  • bragg_k_broadening (float) – Broadening applied to bragg_peaks.

  • figsize (list, tuple, np.ndarray) – Figure size for plot.

  • (bool) (returnfig) – Return figure and axes handles if this is True.

Returns:

figure and axes handles

Return type:

fig, ax (optional)

py4DSTEM.process.diffraction.crystal_viz.plot_orientation_zones(self, azim_elev: list | tuple | ndarray | None = None, proj_dir_lattice: list | tuple | ndarray | None = None, proj_dir_cartesian: list | tuple | ndarray | None = None, tol_den=10, marker_size: float = 20, plot_limit: list | tuple | ndarray = array([-1.1, 1.1]), figsize: list | tuple | ndarray = (8, 8), returnfig: bool = False)

3D scatter plot of the structure factors using magnitude^2, i.e. intensity.

Parameters:
  • azim_elev (array) – az and el angles for plot

  • proj_dir_lattice (array) – (3,) projection direction in lattice

  • proj_dir_cartesian – (array): (3,) projection direction in cartesian

  • tol_den (int) – tolerance for rational index denominator

  • dir_proj (float) – projection direction, either [elev azim] or normal vector Default is mean vector of self.orientation_zone_axis_range rows

  • marker_size (float) – size of markers

  • plot_limit (float) – x y z plot limits, default is [0, 1.05]

  • figsize (2 element float) – size scaling of figure axes

  • returnfig (bool) – set to True to return figure and axes handles

Returns:

fig, ax (optional) figure and axes handles

py4DSTEM.process.diffraction.crystal_viz.plot_orientation_plan(self, index_plot: int = 0, zone_axis_lattice: ndarray | None = None, zone_axis_cartesian: ndarray | None = None, figsize: list | tuple | ndarray = (14, 6), returnfig: bool = False)

3D scatter plot of the structure factors using magnitude^2, i.e. intensity.

Parameters:
  • index_plot (int) – which index slice to plot

  • zone_axis_plot (3 element float) – which zone axis slice to plot

  • figsize (2 element float) – size scaling of figure axes

  • returnfig (bool) – set to True to return figure and axes handles

Returns:

fig, ax (optional) figure and axes handles

py4DSTEM.process.diffraction.crystal_viz.plot_diffraction_pattern(bragg_peaks: PointList, bragg_peaks_compare: PointList | None = None, scale_markers: float = 500, scale_markers_compare: float | None = None, power_markers: float = 1, plot_range_kx_ky: list | tuple | ndarray | None = None, add_labels: bool = True, shift_labels: float = 0.08, shift_marker: float = 0.005, min_marker_size: float = 1e-06, max_marker_size: float = 1000, figsize: list | tuple | ndarray = (12, 6), returnfig: bool = False, input_fig_handle=None)

2D scatter plot of the Bragg peaks

Parameters:
  • bragg_peaks (PointList) – numpy array containing (‘qx’, ‘qy’, ‘intensity’, ‘h’, ‘k’, ‘l’)

  • bragg_peaks_compare (PointList) – numpy array containing (‘qx’, ‘qy’, ‘intensity’)

  • scale_markers (float) – size scaling for markers

  • scale_markers_compare (float) – size scaling for markers of comparison

  • power_markers (float) – power law scaling for marks (default is 1, i.e. amplitude)

  • plot_range_kx_ky (float) – 2 element numpy vector giving the plot range

  • add_labels (bool) – flag to add hkl labels to peaks

  • min_marker_size (float) – minimum marker size for the comparison peaks

  • max_marker_size (float) – maximum marker size for the comparison peaks

  • figsize (2 element float) – size scaling of figure axes

  • returnfig (bool) – set to True to return figure and axes handles

  • input_fig_handle (fig,ax) –

py4DSTEM.process.diffraction.crystal_viz.plot_orientation_maps(self, orientation_map=None, orientation_ind: int = 0, dir_in_plane_degrees: float = 0.0, corr_range: ndarray = array([0, 5]), corr_normalize: bool = True, scale_legend: bool | None = None, figsize: list | tuple | ndarray = (16, 5), figbound: list | tuple | ndarray = (0.01, 0.005), show_axes: bool = True, camera_dist=None, plot_limit=None, plot_layout=0, swap_axes_xy_limits=False, returnfig: bool = False, progress_bar=False)

Plot the orientation maps.

Parameters:
  • orientation_map (OrientationMap) – Class containing orientation matrices, correlation values, etc. Optional - can reference internally stored OrientationMap.

  • orientation_ind (int) – Which orientation match to plot if num_matches > 1

  • dir_in_plane_degrees (float) – In-plane angle to plot in degrees. Default is 0 / x-axis / vertical down.

  • corr_range (np.ndarray) – Correlation intensity range for the plot

  • corr_normalize (bool) – If true, set mean correlation to 1.

  • scale_legend (float) – 2 elements, x and y scaling of legend panel

  • figsize (array) – 2 elements defining figure size

  • figbound (array) – 2 elements defining figure boundary

  • show_axes (bool) – Flag setting whether orienation map axes are visible.

  • camera_dist (float) – distance of camera from legend

  • plot_limit (array) – 2x3 array defining plot boundaries of legend

  • plot_layout (int) – subplot layout: 0 - 1 row, 3 col 1 - 3 row, 1 col

  • swap_axes_xy_limits (bool) – swap x and y boundaries for legend (not sure why we need this in some cases)

  • returnfig (bool) – set to True to return figure and axes handles

  • progress_bar (bool) – Enable progressbar when calculating orientation images.

Returns:

RGB images fig, axs (handles): Figure and axes handes for the

Return type:

images_orientation (int)

Note

Currently, no symmetry reduction. Therefore the x and y orientations are going to be correct only for [001][011][111] orientation triangle.

py4DSTEM.process.diffraction.crystal_viz.plot_fiber_orientation_maps(self, orientation_map, orientation_ind: int = 0, symmetry_order: int | None = None, symmetry_mirror: bool = False, dir_in_plane_degrees: float = 0.0, corr_range: ndarray = array([0, 2]), corr_normalize: bool = True, show_axes: bool = True, medfilt_size: int | None = None, cmap_out_of_plane: str = 'plasma', leg_size: int = 200, figsize: list | tuple | ndarray = (12, 8), figbound: list | tuple | ndarray = (0.005, 0.04), returnfig: bool = False)

Generate and plot the orientation maps from fiber texture plots.

Parameters:
  • orientation_map (OrientationMap) – Class containing orientation matrices, correlation values, etc.

  • orientation_ind (int) – Which orientation match to plot if num_matches > 1

  • dir_in_plane_degrees (float) – Reference in-plane angle (degrees). Default is 0 / x-axis / vertical down.

  • corr_range (np.ndarray) – Correlation intensity range for the plot

  • corr_normalize (bool) – If true, set mean correlation to 1.

  • show_axes (bool) – Flag setting whether orienation map axes are visible.

  • figsize (array) – 2 elements defining figure size

  • figbound (array) – 2 elements defining figure boundary

  • returnfig (bool) – set to True to return figure and axes handles

Returns:

RGB images fig, axs (handles): Figure and axes handes for the

Return type:

images_orientation (int)

Note

Currently, no symmetry reduction. Therefore the x and y orientations are going to be correct only for [001][011][111] orientation triangle.

py4DSTEM.process.diffraction.crystal_viz.plot_clusters(self, area_min=2, outline_grains=True, outline_thickness=1, fill_grains=0.25, smooth_grains=1.0, cmap='viridis', figsize=(8, 8), returnfig=False)

Plot the clusters as an image.

Parameters:
  • area_min (int (optional)) – Min cluster size to include, in units of probe positions.

  • outline_grains (bool (optional)) – Set to True to draw grains with outlines

  • outline_thickness (int (optional)) – Thickenss of the grain outline

  • fill_grains (float (optional)) – Outlined grains are filled with this value in pixels.

  • smooth_grains (float (optional)) – Grain boundaries are smoothed by this value in pixels.

  • figsize (tuple) – Size of the figure panel

  • returnfig (bool) – Setting this to true returns the figure and axis handles

Returns:

Figure and axes handles

Return type:

fig, ax (optional)

py4DSTEM.process.diffraction.crystal_viz.plot_cluster_size(self, area_min=None, area_max=None, area_step=1, weight_intensity=False, pixel_area=1.0, pixel_area_units='px^2', figsize=(8, 6), returnfig=False)

Plot the cluster sizes

Parameters:
  • area_min (int (optional)) – Min area to include in pixels^2

  • area_max (int (optional)) – Max area bin in pixels^2

  • area_step (int (optional)) – Step size of the histogram bin in pixels^2

  • weight_intensity (bool) – Weight histogram by the peak intensity.

  • pixel_area (float) – Size of pixel area unit square

  • pixel_area_units (string) – Units of the pixel area

  • figsize (tuple) – Size of the figure panel

  • returnfig (bool) – Setting this to true returns the figure and axis handles

Returns:

Figure and axes handles

Return type:

fig, ax (optional)

py4DSTEM.process.diffraction.crystal_viz.atomic_colors(Z, scheme='jmol')

Return atomic colors for Z.

Modes are “colin” and “jmol”. “colin” uses the handmade but incomplete scheme of Colin Ophus “jmol” uses the JMOL scheme, from http://jmol.sourceforge.net/jscolors

which includes all elements up to 109

py4DSTEM.process.diffraction.crystal_viz.plot_ring_pattern(radii, intensity, theta=[-3.141592653589793, 3.141592653589793, 200], intensity_scale=1, intensity_constant=False, color='k', figsize=(10, 10), returnfig=False, input_fig_handle=None, **kwargs)

2D plot of diffraction rings

Parameters:
  • radii (PointList) – 1D numpy array containing radii for diffraction rings

  • intensity (PointList) – 1D numpy array containing intensities for diffraciton rings

  • theta (3-tuple) – first two values specify angle range, and the last specifies the number of points used for plotting

  • intensity_scale (float) – size scaling for ring thickness

  • intensity_constant (bool) – if true, all rings are plotted with same line width

  • color (matplotlib color) – color of ring, any format recognized by matplotlib

  • figsize (2 element float) – size scaling of figure axes

  • returnfig (bool) – set to True to return figure and axes handles

  • input_fig_handle (fig,ax) –

py4DSTEM.process.diffraction.flowlines.make_orientation_histogram(bragg_peaks: PointListArray | None = None, radial_ranges: ndarray | None = None, orientation_map=None, orientation_ind: int = 0, orientation_growth_angles: array = 0.0, orientation_separate_bins: bool = False, orientation_flip_sign: bool = False, upsample_factor=4.0, theta_step_deg=1.0, sigma_x=1.0, sigma_y=1.0, sigma_theta=3.0, normalize_intensity_image: bool = False, normalize_intensity_stack: bool = True, progress_bar: bool = True)

Create an 3D or 4D orientation histogram from a braggpeaks PointListArray from user-specified radial ranges, or from the Euler angles from a fiber texture OrientationMap generated by the ACOM module of py4DSTEM.

Parameters:
  • bragg_peaks (PointListArray) – 2D of pointlists containing centered peak locations.

  • radial_ranges (np array) – Size (N x 2) array for N radial bins, or (2,) for a single bin.

  • orientation_map (OrientationMap) – Class containing the Euler angles to generate a flowline map.

  • orientation_ind (int) – Index of the orientation map (default 0)

  • orientation_growth_angles (array) – Angles to place into histogram, relative to orientation.

  • orientation_separate_bins (bool) – whether to place multiple angles into multiple radial bins.

  • upsample_factor (float) – Upsample factor

  • theta_step_deg (float) – Step size along annular direction in degrees

  • sigma_x (float) – Smoothing in x direction before upsample

  • sigma_y (float) – Smoothing in x direction before upsample

  • sigma_theta (float) – Smoothing in annular direction (units of bins, periodic)

  • normalize_intensity_image (bool) – Normalize to max peak intensity = 1, per image

  • normalize_intensity_stack (bool) – Normalize to max peak intensity = 1, all images

  • progress_bar (bool) – Enable progress bar

Returns:

4D array containing Bragg peak intensity histogram

[radial_bin x_probe y_probe theta]

Return type:

orient_hist (array)

py4DSTEM.process.diffraction.flowlines.make_flowline_map(orient_hist, thresh_seed=0.2, thresh_grow=0.05, thresh_collision=0.001, sep_seeds=None, sep_xy=6.0, sep_theta=5.0, sort_seeds='intensity', linewidth=2.0, step_size=0.5, min_steps=4, max_steps=1000, sigma_x=1.0, sigma_y=1.0, sigma_theta=2.0, progress_bar: bool = True)

Create an 3D or 4D orientation flowline map - essentially a pixelated “stream map” which represents diffraction data.

Parameters:
  • orient_hist (array) – Histogram of all orientations with coordinates [radial_bin x_probe y_probe theta] We assume theta bin ranges from 0 to 180 degrees and is periodic.

  • thresh_seed (float) – Threshold for seed generation in histogram.

  • thresh_grow (float) – Threshold for flowline growth in histogram.

  • thresh_collision (float) – Threshold for termination of flowline growth in histogram.

  • sep_seeds (float) – Initial seed separation in bins - set to None to use default value, which is equal to 0.5*sep_xy.

  • sep_xy (float) – Search radius for flowline direction in x and y.

  • = (sep_theta) – Search radius for flowline direction in theta.

  • sort_seeds (str) – How to sort the initial seeds for growth: None - no sorting ‘intensity’ - sort by histogram intensity ‘random’ - random order

  • linewidth (float) – Thickness of the flowlines in pixels.

  • step_size (float) – Step size for flowline growth in pixels.

  • min_steps (int) – Minimum number of steps for a flowline to be drawn.

  • max_steps (int) – Maximum number of steps for a flowline to be drawn.

  • sigma_x (float) – Weighted sigma in x direction for direction update.

  • sigma_y (float) – Weighted sigma in y direction for direction update.

  • sigma_theta (float) – Weighted sigma in theta for direction update.

  • progress_bar (bool) – Enable progress bar

Returns:

4D array containing flowlines

[radial_bin x_probe y_probe theta]

Return type:

orient_flowlines (array)

py4DSTEM.process.diffraction.flowlines.make_flowline_rainbow_image(orient_flowlines, int_range=[0, 0.2], sym_rotation_order=2, theta_offset=0.0, greyscale=False, greyscale_max=True, white_background=False, power_scaling=1.0, sum_radial_bins=False, plot_images=True, figsize=None)

Generate RGB output images from the flowline arrays.

Parameters:
  • orient_flowline (array) – Histogram of all orientations with coordinates [x y radial_bin theta] We assume theta bin ranges from 0 to 180 degrees and is periodic.

  • int_range (float) –

  • sym_rotation_order (int) – rotational symmety for colouring

  • theta_offset (float) – Offset the anglular coloring by this value in radians.

  • greyscale (bool) – Set to False for color output, True for greyscale output.

  • greyscale_max (bool) – If output is greyscale, use max instead of mean for overlapping flowlines.

  • white_background (bool) – For either color or greyscale output, switch to white background (from black).

  • power_scaling (float) – Power law scaling for flowline intensity output.

  • sum_radial_bins (bool) – Sum all radial bins (alternative is to output separate images).

  • plot_images (bool) – Plot the outputs for quick visualization.

  • figsize (2-tuple) – Size of output figure.

Returns:

3D or 4D array containing flowline images

Return type:

im_flowline (array)

py4DSTEM.process.diffraction.flowlines.make_flowline_rainbow_legend(im_size=array([256, 256]), sym_rotation_order=2, theta_offset=0.0, white_background=False, return_image=False, radial_range=array([0.45, 0.9]), plot_legend=True, figsize=(4, 4))

This function generates a legend for a the rainbow colored flowline maps, and returns it as an RGB image.

Parameters:
  • im_size (np.array) – Size of legend image in pixels.

  • sym_rotation_order (int) – rotational symmety for colouring

  • theta_offset (float) – Offset the anglular coloring by this value in radians.

  • white_background (bool) – For either color or greyscale output, switch to white background (from black).

  • return_image (bool) – Return the image array.

  • radial_range (np.array) – Inner and outer radius for the legend ring.

  • plot_legend (bool) – Plot the generated legend.

  • figsize (tuple or list) – Size of the plotted legend.

Returns:

Image array for the legend.

Return type:

im_legend (array)

py4DSTEM.process.diffraction.flowlines.make_flowline_combined_image(orient_flowlines, int_range=[0, 0.2], cvals=array([[0., 0.7, 0.], [1., 0., 0.], [0., 0.7, 1.]]), white_background=False, power_scaling=1.0, sum_radial_bins=True, plot_images=True, figsize=None)

Generate RGB output images from the flowline arrays.

Parameters:
  • orient_flowline (array) – Histogram of all orientations with coordinates [x y radial_bin theta] We assume theta bin ranges from 0 to 180 degrees and is periodic.

  • int_range (float) –

  • cvals (array) – Nx3 size array containing RGB colors for different radial ibns.

  • white_background (bool) – For either color or greyscale output, switch to white background (from black).

  • power_scaling (float) – Power law scaling for flowline intensities.

  • sum_radial_bins (bool) – Sum outputs over radial bins.

  • plot_images (bool) – Plot the output images for quick visualization.

  • figsize (2-tuple) – Size of output figure.

Returns:

flowline images

Return type:

im_flowline (array)

py4DSTEM.process.diffraction.flowlines.orientation_correlation(orient_hist, radius_max=None, progress_bar=True)

Take in the 4D orientation histogram, and compute the distance-angle (auto)correlations

Parameters:
  • orient_hist (array) – 3D or 4D histogram of all orientations with coordinates [x y radial_bin theta]

  • radius_max (float) – Maximum radial distance for correlogram calculation. If set to None, the maximum radius will be set to min(orient_hist.shape[0],orient_hist.shape[1])/2.

Returns:

3D or 4D array containing correlation images as function of (dr,dtheta)

Return type:

orient_corr (array)

py4DSTEM.process.diffraction.flowlines.plot_orientation_correlation(orient_corr, prob_range=[0.1, 10.0], calculate_coefs=False, fraction_coefs=0.5, length_fit_slope=10, plot_overlaid_coefs=True, inds_plot=None, pixel_size=None, pixel_units=None, fontsize=10, figsize=(8, 6), returnfig=False)

Plot the distance-angle (auto)correlations in orient_corr.

Parameters:
  • (array) (figsize) – 3D or 4D array containing correlation images as function of (dr,dtheta) 1st index represents each pair of rings.

  • (array) – Plotting range in units of “multiples of random distribution”.

  • (bool) (returnfig) – If this value is True, the 0.5 and 0.1 distribution fraction of the radial and annular correlations will be calculated and printed.

  • (float) (fontsize) – What fraction to calculate the correlation distribution coefficients for.

  • (int) (length_fit_slope) – Number of pixels to fit the slope of angular vs radial intercept.

  • (bool) – If this value is True, the 0.5 and 0.1 distribution fraction of the radial and annular correlations will be overlaid onto the plots.

  • (float) – Which indices to plot for orient_corr. Set to “None” to plot all pairs.

  • (float) – Pixel size for x axis.

  • (str) (pixel_units) – units of pixels.

  • (float) – Font size. Title will be slightly larger, axis slightly smaller.

  • (array) – Size of the figure panels.

  • (bool) – Set to True to return figure axes.

Returns:

Figure and axes handles (optional).

Return type:

fig, ax (handles)

This module provides access to some objects used or maintained by the interpreter and to functions that interact strongly with the interpreter.

Dynamic objects:

argv – command line arguments; argv[0] is the script pathname if known path – module search path; path[0] is the script directory, else ‘’ modules – dictionary of loaded modules

displayhook – called to show results in an interactive session excepthook – called to handle any uncaught exception other than SystemExit

To customize printing in an interactive session or to install a custom top-level exception handler, assign other functions to replace these.

stdin – standard input file object; used by input() stdout – standard output file object; used by print() stderr – standard error object; used for error messages

By assigning other file objects (or objects that behave like files) to these, it is possible to redirect all of the interpreter’s I/O.

last_type – type of last uncaught exception last_value – value of last uncaught exception last_traceback – traceback of last uncaught exception

These three are only available in an interactive session after a traceback has been printed.

Static objects:

builtin_module_names – tuple of module names built into this interpreter copyright – copyright notice pertaining to this interpreter exec_prefix – prefix used to find the machine-specific Python library executable – absolute path of the executable binary of the Python interpreter float_info – a named tuple with information about the float implementation. float_repr_style – string indicating the style of repr() output for floats hash_info – a named tuple with information about the hash algorithm. hexversion – version information encoded as a single integer implementation – Python implementation information. int_info – a named tuple with information about the int implementation. maxsize – the largest supported length of containers. maxunicode – the value of the largest Unicode code point platform – platform identifier prefix – prefix used to find the Python library thread_info – a named tuple with information about the thread implementation. version – the version of this interpreter as a string version_info – version information as a named tuple __stdin__ – the original stdin; don’t touch! __stdout__ – the original stdout; don’t touch! __stderr__ – the original stderr; don’t touch! __displayhook__ – the original displayhook; don’t touch! __excepthook__ – the original excepthook; don’t touch!

Functions:

displayhook() – print an object to the screen, and save it in builtins._ excepthook() – print an exception and its traceback to sys.stderr exc_info() – return thread-safe information about the current exception exit() – exit the interpreter by raising SystemExit getdlopenflags() – returns flags to be used for dlopen() calls getprofile() – get the global profiling function getrefcount() – return the reference count for an object (plus one :-) getrecursionlimit() – return the max recursion depth for the interpreter getsizeof() – return the size of an object in bytes gettrace() – get the global debug tracing function setdlopenflags() – set the flags to be used for dlopen() calls setprofile() – set the global profiling function setrecursionlimit() – set the max recursion depth for the interpreter settrace() – set the global debug tracing function

class py4DSTEM.process.diffraction.utils.Orientation(num_matches: int)

A class for storing output orientations, generated by fitting a Crystal class orientation plan or Bloch wave pattern matching to a PointList.

__init__(num_matches: int) None
class py4DSTEM.process.diffraction.utils.OrientationMap(num_x: int, num_y: int, num_matches: int)

A class for storing output orientations, generated by fitting a Crystal class orientation plan or Bloch wave pattern matching to a PointListArray.

__init__(num_x: int, num_y: int, num_matches: int) None
py4DSTEM.process.diffraction.utils.sort_orientation_maps(orientation_map, sort='intensity', cluster_thresh=0.1)

Sort the orientation maps along the ind_match direction, either by intensity or by clustering similar angles (greedily, in order of intensity).

Parameters:
  • OrientationMap (orientation_map Initial) –

  • sort (string) – “intensity” or “cluster” for sorting method.

  • cluster_thresh (float) – similarity threshold for clustering method

Returns:

orientation_sort Sorted OrientationMap

py4DSTEM.process.diffraction.utils.calc_1D_profile(k, g_coords, g_int, remove_origin=True, k_broadening=0.0, int_scale=None, normalize_intensity=True)

Utility function to calculate a 1D histogram from the diffraction vector lengths stored in a Crystal class.

Parameters:
  • k (np.array) – k coordinates.

  • g_coords (np.array) – Scattering vector lengths g.

  • bragg_intensity_power (np.array) – Scattering vector intensities.

  • remove_origin (bool) – Remove the origin peak from the profile.

  • k_broadening (float) – Broadening applied to full profile.

  • int_scale (np.array) – Either a scalar value mulitiplied into all peak intensities, or a vector with 1 value per peak to scale peaks individually.

  • normalize_intensity – Normalize maximum output value to 1.

diskdetection

fit

py4DSTEM.process.fit.fit.fit_1D_gaussian(xdata, ydata, xmin, xmax)

Fits a 1D gaussian to the subset of the 1D curve f(xdata)=ydata within the window (xmin,xmax). Returns A,mu,sigma. Retrieve the full curve with

>>> fit_gaussian = py4DSTEM.process.fit.gaussian(xdata,A,mu,sigma)
py4DSTEM.process.fit.fit.fit_2D(function, data, data_mask=None, popt=None, robust=False, robust_steps=3, robust_thresh=2)

Performs a 2D fit.

TODO: make returning the mask optional

Parameters:
  • function (callable) – Some function( xy, **p) where xy is a length 2 vector (1D np array) specifying the pixel position (x,y), and p is the function parameters

  • data (ndarray) – Some 2D array of any shape (n,m)

  • data_mask (None or boolean array of shape (n,m), optional) – If specified, fits only the pixels in data where this array is True

  • popt (dict) – Initial guess at the parameters p of function. Note that positions in pixels (i.e. the xy positions) are linearly scaled to the space [0,1]

  • robust (bool) – Toggles robust fitting, which iteratively rejects outlier data points which have a root-mean-square error beyond robust_thresh

  • robust_steps (int) – The number of robust fitting iterations to perform

  • robust_thresh (int) – The robust fitting cutoff

  • Returns

  • (popt (4-tuple) – The optimal fit parameters, the fitting covariance matrix, the the fit array with the returned popt params, and the mask

  • pcov (4-tuple) – The optimal fit parameters, the fitting covariance matrix, the the fit array with the returned popt params, and the mask

  • fit_at (4-tuple) – The optimal fit parameters, the fitting covariance matrix, the the fit array with the returned popt params, and the mask

  • mask) (4-tuple) – The optimal fit parameters, the fitting covariance matrix, the the fit array with the returned popt params, and the mask

py4DSTEM.process.fit.fit.fit_2D_polar_gaussian(data, mask=None, p0=None, robust=False, robust_steps=3, robust_thresh=2, constant_background=False)

NOTE - this cannot work without using pixel coordinates - something is wrong in the workflow.

Fits a 2D gaussian to the pixels in data which are set to True in mask.

The gaussian is anisotropic and oriented along (t,q), centered at (mu_t,mu_q), has standard deviations (sigma_t,sigma_q), maximum of I0, and an optional constant offset of C, and is periodic in t.

f(x,y) = I0 * exp( - (x-mu_x)^2/(2sig_x^2) + (y-mu_y)^2/(2sig_y^2) ) or f(x,y) = I0 * exp( - (x-mu_x)^2/(2sig_x^2) + (y-mu_y)^2/(2sig_y^2) ) + C

Parameters:
  • data (2d array) – the data to fit

  • p0 (6-tuple) – initial guess at fit parameters, (I0,mu_x,mu_y,sigma_x_sigma_y,C)

  • mask (2d boolean array) – ignore pixels where mask is False

  • robust (bool) – toggle robust fitting

  • robust_steps (int) – number of robust fit iterations

  • robust_thresh (number) – the robust fitting threshold

  • constant_background (bool) – whether or not to include constant background

Returns:

(popt,pcov,fit_ar) – the optimal fit parameters, the covariance matrix, and the fit array

Return type:

3-tuple

latticevectors

phase

py4DSTEM.process.phase.utils.polar_symbols = ('C10', 'C12', 'phi12', 'C21', 'phi21', 'C23', 'phi23', 'C30', 'C32', 'phi32', 'C34', 'phi34', 'C41', 'phi41', 'C43', 'phi43', 'C45', 'phi45', 'C50', 'C52', 'phi52', 'C54', 'phi54', 'C56', 'phi56')

Symbols for the polar representation of all optical aberrations up to the fifth order.

py4DSTEM.process.phase.utils.polar_aliases = {'C5': 'C50', 'Cs': 'C30', 'astigmatism': 'C12', 'astigmatism_angle': 'phi12', 'coma': 'C21', 'coma_angle': 'phi21', 'defocus': 'C10'}

Aliases for the most commonly used optical aberrations.

class py4DSTEM.process.phase.utils.ComplexProbe(energy: float, gpts: Tuple[int, int], sampling: Tuple[float, float], semiangle_cutoff: float = inf, rolloff: float = 2.0, vacuum_probe_intensity: ndarray | None = None, device: str = 'cpu', focal_spread: float = 0.0, angular_spread: float = 0.0, gaussian_spread: float = 0.0, phase_shift: float = 0.0, parameters: Mapping[str, float] | None = None, **kwargs)

Complex Probe Class.

Simplified version of CTF and Probe from abTEM: https://github.com/abTEM/abTEM/blob/master/abtem/transfer.py https://github.com/abTEM/abTEM/blob/master/abtem/waves.py

Parameters:
  • energy (float) – The electron energy of the wave functions this contrast transfer function will be applied to [eV].

  • semiangle_cutoff (float) – The semiangle cutoff describes the sharp Fourier space cutoff due to the objective aperture [mrad].

  • gpts (Tuple[int,int]) – Number of grid points describing the wave functions.

  • sampling (Tuple[float,float]) – Lateral sampling of wave functions in Å

  • device (str, optional) – Device to perform calculations on. Must be either ‘cpu’ or ‘gpu’

  • rolloff (float, optional) – Tapers the cutoff edge over the given angular range [mrad].

  • vacuum_probe_intensity (np.ndarray, optional) – Squared of corner-centered aperture amplitude to use, instead of semiangle_cutoff + rolloff

  • focal_spread (float, optional) – The 1/e width of the focal spread due to chromatic aberration and lens current instability [Å].

  • angular_spread (float, optional) – The 1/e width of the angular deviations due to source size [mrad].

  • gaussian_spread (float, optional) – The 1/e width image deflections due to vibrations and thermal magnetic noise [Å].

  • phase_shift (float, optional) – A constant phase shift [radians].

  • parameters (dict, optional) – Mapping from aberration symbols to their corresponding values. All aberration magnitudes should be given in Å and angles should be given in radians.

  • kwargs – Provide the aberration coefficients as keyword arguments.

__init__(energy: float, gpts: Tuple[int, int], sampling: Tuple[float, float], semiangle_cutoff: float = inf, rolloff: float = 2.0, vacuum_probe_intensity: ndarray | None = None, device: str = 'cpu', focal_spread: float = 0.0, angular_spread: float = 0.0, gaussian_spread: float = 0.0, phase_shift: float = 0.0, parameters: Mapping[str, float] | None = None, **kwargs)
set_parameters(parameters: dict)

Set the phase of the phase aberration. :param parameters: Mapping from aberration symbols to their corresponding values. :type parameters: dict

polar_coordinates(x, y)

Calculate a polar grid for a given Cartesian grid.

build()

Builds corner-centered complex probe in the center of the region of interest.

visualize(**kwargs)

Plots the probe intensity.

py4DSTEM.process.phase.utils.spatial_frequencies(gpts: ~typing.Tuple[int, int], sampling: ~typing.Tuple[float, float], xp=<module 'numpy' from '/home/docs/checkouts/readthedocs.org/user_builds/py4dstem/envs/latest/lib/python3.10/site-packages/numpy/__init__.py'>)

Calculate spatial frequencies of a grid.

Parameters:
  • gpts (tuple of int) – Number of grid points.

  • sampling (tuple of float) – Sampling of the potential [1 / Å].

Return type:

tuple of arrays

py4DSTEM.process.phase.utils.fourier_translation_operator(positions: ~numpy.ndarray, shape: tuple, xp=<module 'numpy' from '/home/docs/checkouts/readthedocs.org/user_builds/py4dstem/envs/latest/lib/python3.10/site-packages/numpy/__init__.py'>) ndarray

Create an array representing one or more phase ramp(s) for shifting another array.

Parameters:
  • positions (array of xy-positions) – Positions to calculate fourier translation operators for

  • shape (two int) – Array dimensions to be fourier-shifted

  • xp (Callable) – Array computing module

Return type:

Fourier translation operators

py4DSTEM.process.phase.utils.fft_shift(array, positions, xp=<module 'numpy' from '/home/docs/checkouts/readthedocs.org/user_builds/py4dstem/envs/latest/lib/python3.10/site-packages/numpy/__init__.py'>)

Fourier-shift array using positions.

Parameters:
  • array (np.ndarray) – Array to be shifted

  • positions (array of xy-positions) – Positions to fourier-shift array with

  • xp (Callable) – Array computing module

Return type:

Fourier-shifted array

py4DSTEM.process.phase.utils.subdivide_into_batches(num_items: int, num_batches: int | None = None, max_batch: int | None = None)

Split an n integer into m (almost) equal integers, such that the sum of smaller integers equals n.

Parameters:
  • n (int) – The integer to split.

  • m (int) – The number integers n will be split into.

Return type:

list of int

class py4DSTEM.process.phase.utils.AffineTransform(scale0: float = 1.0, scale1: float = 1.0, shear1: float = 0.0, angle: float = 0.0, t0: float = 0.0, t1: float = 0.0, dilation: float = 1.0)

Affine Transform Class.

Simplified version of AffineTransform from tike: https://github.com/AdvancedPhotonSource/tike/blob/f9004a32fda5e49fa63b987e9ffe3c8447d59950/src/tike/ptycho/position.py

AffineTransform() -> Identity

Parameters:
  • scale0 (float) – x-scaling

  • scale1 (float) – y-scaling

  • shear1 (float) – gamma shear

  • angle (float) – theta rotation angle

  • t0 (float) – x-translation

  • t1 (float) – y-translation

  • dilation (float) – Isotropic expansion (multiplies scale0 and scale1)

__init__(scale0: float = 1.0, scale1: float = 1.0, shear1: float = 0.0, angle: float = 0.0, t0: float = 0.0, t1: float = 0.0, dilation: float = 1.0)
classmethod fromarray(T: ndarray)

Return an Affine Transfrom from a 2x2 matrix. Use decomposition method from Graphics Gems 2 Section 7.1

asarray()

Return an 2x2 matrix of scale, shear, rotation. This matrix is scale @ shear @ rotate from left to right.

asarray3()

Return an 3x2 matrix of scale, shear, rotation, translation. This matrix is scale @ shear @ rotate from left to right. Expects a homogenous (z) coordinate of 1.

astuple()

Return the constructor parameters in a tuple.

py4DSTEM.process.phase.utils.estimate_global_transformation(positions0: ~numpy.ndarray, positions1: ~numpy.ndarray, origin: ~typing.Tuple[int, int] = (0, 0), translation_allowed: bool = True, xp=<module 'numpy' from '/home/docs/checkouts/readthedocs.org/user_builds/py4dstem/envs/latest/lib/python3.10/site-packages/numpy/__init__.py'>)

Use least squares to estimate the global affine transformation.

py4DSTEM.process.phase.utils.estimate_global_transformation_ransac(positions0: ~numpy.ndarray, positions1: ~numpy.ndarray, origin: ~typing.Tuple[int, int] = (0, 0), translation_allowed: bool = True, min_sample: int = 64, max_error: float = 16, min_consensus: float = 0.75, max_iter: int = 20, xp=<module 'numpy' from '/home/docs/checkouts/readthedocs.org/user_builds/py4dstem/envs/latest/lib/python3.10/site-packages/numpy/__init__.py'>)

Use RANSAC to estimate the global affine transformation.

py4DSTEM.process.phase.utils.fourier_ring_correlation(image_1, image_2, pixel_size=None, bin_size=None, sigma=None, align_images=False, upsample_factor=8, device='cpu', plot_frc=True, frc_color='red', half_bit_color='blue')

Computes fourier ring correlation (FRC) of 2 arrays. Arrays must bet the same size.

Parameters
image1: ndarray

first image for FRC

image2: ndarray

second image for FRC

pixel_size: tuple

size of pixels in A (x,y)

bin_size: float, optional

size of bins for ring profile

sigma: float, optional

standard deviation for Gaussian kernel

align_images: bool

if True, aligns images using DFT upsampling of cross correlation.

upsample factor: int

if align_images, upsampling for correlation. Must be greater than 2.

device: str, optional

calculation device will be perfomed on. Must be ‘cpu’ or ‘gpu’

plot_frc: bool, optional

if True, plots frc

frc_color: str, optional

color of FRC line in plot

half_bit_color: str, optional

color of half-bit line

Returns:

  • q_frc (ndarray) – spatial frequencies of FRC

  • frc (ndarray) – fourier ring correlation

  • half_bit (ndarray) – half-bit criteria

py4DSTEM.process.phase.utils.return_1D_profile(intensity, pixel_size=None, bin_size=None, sigma=None, device='cpu')

Return 1D radial profile from corner centered array

Parameters
intensity: ndarray

Array for computing 1D profile

pixel_size: tuple

Size of pixels in A (x,y)

bin_size: float, optional

Size of bins for ring profile

sigma: float, optional

standard deviation for Gaussian kernel

device: str, optional

calculation device will be perfomed on. Must be ‘cpu’ or ‘gpu’

Returns:

  • q_bins (ndarray) – spatial frequencies of bins

  • I_bins (ndarray) – Intensity of bins

  • n (ndarray) – Number of pixels in each bin

py4DSTEM.process.phase.utils.fourier_rotate_real_volume(array, angle, axes=(0, 1), xp=<module 'numpy' from '/home/docs/checkouts/readthedocs.org/user_builds/py4dstem/envs/latest/lib/python3.10/site-packages/numpy/__init__.py'>)

Rotates a 3D array using three Fourier-based shear operators.

Parameters
array: ndarray

3D array to rotate

angle: float

Angle in deg to rotate array by

axes: tuple, Optional

Axes defining plane in which to rotate about

xp: Callable, optional

Array computing module

Returns:

output_arr – Fourier-rotated array

Return type:

ndarray

py4DSTEM.process.phase.utils.array_slice(axis, ndim, start, end, step=1)

Returns array slice along dynamic axis

py4DSTEM.process.phase.utils.periodic_centered_difference(array, spacing, axis, xp=<module 'numpy' from '/home/docs/checkouts/readthedocs.org/user_builds/py4dstem/envs/latest/lib/python3.10/site-packages/numpy/__init__.py'>)

Computes second-order centered difference with periodic BCs

py4DSTEM.process.phase.utils.compute_divergence_periodic(vector_field, spacings, xp=<module 'numpy' from '/home/docs/checkouts/readthedocs.org/user_builds/py4dstem/envs/latest/lib/python3.10/site-packages/numpy/__init__.py'>)

Computes divergence of vector_field

py4DSTEM.process.phase.utils.compute_gradient_periodic(scalar_field, spacings, xp=<module 'numpy' from '/home/docs/checkouts/readthedocs.org/user_builds/py4dstem/envs/latest/lib/python3.10/site-packages/numpy/__init__.py'>)

Computes gradient of scalar_field

py4DSTEM.process.phase.utils.preconditioned_laplacian_periodic_3D(shape, xp=<module 'numpy' from '/home/docs/checkouts/readthedocs.org/user_builds/py4dstem/envs/latest/lib/python3.10/site-packages/numpy/__init__.py'>)

FFT eigenvalues

py4DSTEM.process.phase.utils.preconditioned_poisson_solver_periodic_3D(rhs, gauge=None, xp=<module 'numpy' from '/home/docs/checkouts/readthedocs.org/user_builds/py4dstem/envs/latest/lib/python3.10/site-packages/numpy/__init__.py'>)

FFT based poisson solver

py4DSTEM.process.phase.utils.project_vector_field_divergence_periodic_3D(vector_field, xp=<module 'numpy' from '/home/docs/checkouts/readthedocs.org/user_builds/py4dstem/envs/latest/lib/python3.10/site-packages/numpy/__init__.py'>)

Returns solenoidal part of vector field using projection:

f - grad{p} s.t. laplacian{p} = div{f}

py4DSTEM.process.phase.utils.cartesian_to_polar_transform_2Ddata(im_cart, xy_center, num_theta_bins=90, radius_max=None, corner_centered=False, xp=<module 'numpy' from '/home/docs/checkouts/readthedocs.org/user_builds/py4dstem/envs/latest/lib/python3.10/site-packages/numpy/__init__.py'>)

Quick cartesian to polar conversion.

py4DSTEM.process.phase.utils.polar_to_cartesian_transform_2Ddata(im_polar, xy_size, xy_center, corner_centered=False, xp=<module 'numpy' from '/home/docs/checkouts/readthedocs.org/user_builds/py4dstem/envs/latest/lib/python3.10/site-packages/numpy/__init__.py'>)

Quick polar to cartesian conversion.

py4DSTEM.process.phase.utils.regularize_probe_amplitude(probe_init, width_max_pixels=2.0, nearest_angular_neighbor_averaging=5, enforce_constant_intensity=True, corner_centered=False)

Fits sigmoid for each angular direction.

Parameters:
  • probe_init (np.array) – 2D complex image of the probe in Fourier space.

  • width_max_pixels (float) – Maximum edge width of the probe in pixels.

  • nearest_angular_neighbor_averaging (int) – Number of nearest angular neighbor pixels to average to make aperture less jagged.

  • enforce_constant_intensity (bool) – Set to true to make intensity inside the aperture constant.

  • corner_centered (bool) – If True, the probe is assumed to be corner-centered

Returns:

  • probe_corr (np.ndarray) – 2D complex image of the corrected probe in Fourier space.

  • coefs_all (np.ndarray) – coefficients for the sigmoid fits

py4DSTEM.process.phase.utils.interleave_ndarray_symmetrically(array_nd, axis, xp=<module 'numpy' from '/home/docs/checkouts/readthedocs.org/user_builds/py4dstem/envs/latest/lib/python3.10/site-packages/numpy/__init__.py'>)

[a,b,c,d,e,f] -> [a,c,e,f,d,b]

py4DSTEM.process.phase.utils.dct_II_using_FFT_base(array_nd, xp=<module 'numpy' from '/home/docs/checkouts/readthedocs.org/user_builds/py4dstem/envs/latest/lib/python3.10/site-packages/numpy/__init__.py'>)

FFT-based DCT-II

py4DSTEM.process.phase.utils.interleave_ndarray_symmetrically_inverse(array_nd, axis, xp=<module 'numpy' from '/home/docs/checkouts/readthedocs.org/user_builds/py4dstem/envs/latest/lib/python3.10/site-packages/numpy/__init__.py'>)

[a,c,e,f,d,b] -> [a,b,c,d,e,f]

py4DSTEM.process.phase.utils.idct_II_using_FFT_base(array_nd, xp=<module 'numpy' from '/home/docs/checkouts/readthedocs.org/user_builds/py4dstem/envs/latest/lib/python3.10/site-packages/numpy/__init__.py'>)

FFT-based IDCT-II

py4DSTEM.process.phase.utils.idct_II_using_FFT(array_nd, xp=<module 'numpy' from '/home/docs/checkouts/readthedocs.org/user_builds/py4dstem/envs/latest/lib/python3.10/site-packages/numpy/__init__.py'>)

FFT-based IDCT-II

py4DSTEM.process.phase.utils.preconditioned_laplacian_neumann_2D(shape, xp=<module 'numpy' from '/home/docs/checkouts/readthedocs.org/user_builds/py4dstem/envs/latest/lib/python3.10/site-packages/numpy/__init__.py'>)

DCT eigenvalues

py4DSTEM.process.phase.utils.preconditioned_poisson_solver_neumann_2D(rhs, gauge=None, xp=<module 'numpy' from '/home/docs/checkouts/readthedocs.org/user_builds/py4dstem/envs/latest/lib/python3.10/site-packages/numpy/__init__.py'>)

DCT based poisson solver

py4DSTEM.process.phase.utils.unwrap_phase_2d(array, weights=None, gauge=None, corner_centered=True, xp=<module 'numpy' from '/home/docs/checkouts/readthedocs.org/user_builds/py4dstem/envs/latest/lib/python3.10/site-packages/numpy/__init__.py'>)

Weigted phase unwrapping using DCT-based poisson solver

py4DSTEM.process.phase.utils.rotate_point(origin, point, angle)

Rotate a point (x1, y1) counterclockwise by a given angle around a given origin (x0, y0).

Parameters:
  • origin (2-tuple of floats) – (x0, y0)

  • point (2-tuple of floats) – (x1, y1)

  • angle (float (radians)) –

Return type:

rotated points (2-tuple)

py4DSTEM.process.phase.utils.bilinearly_interpolate_array(image, xa, ya, xp=<module 'numpy' from '/home/docs/checkouts/readthedocs.org/user_builds/py4dstem/envs/latest/lib/python3.10/site-packages/numpy/__init__.py'>)

Bilinear sampling of intensities from an image array and pixel positions.

Parameters:
  • image (np.ndarray) – Image array to sample from

  • xa (np.ndarray) – Vertical interpolation sampling positions of image array in pixels

  • ya (np.ndarray) – Horizontal interpolation sampling positions of image array in pixels

Returns:

intensities – Bilinearly-sampled intensities of array at (xa,ya) positions

Return type:

np.ndarray

py4DSTEM.process.phase.utils.lanczos_interpolate_array(image, xa, ya, alpha, xp=<module 'numpy' from '/home/docs/checkouts/readthedocs.org/user_builds/py4dstem/envs/latest/lib/python3.10/site-packages/numpy/__init__.py'>)

Lanczos sampling of intensities from an image array and pixel positions.

Parameters:
  • image (np.ndarray) – Image array to sample from

  • xa (np.ndarray) – Vertical Interpolation sampling positions of image array in pixels

  • ya (np.ndarray) – Horizontal interpolation sampling positions of image array in pixels

  • alpha (int) – Lanczos kernel order

Returns:

intensities – Lanczos-sampled intensities of array at (xa,ya) positions

Return type:

np.ndarray

py4DSTEM.process.phase.utils.pixel_rolling_kernel_density_estimate(stack, shifts, upsampling_factor, kde_sigma, lowpass_filter=False, xp=<module 'numpy' from '/home/docs/checkouts/readthedocs.org/user_builds/py4dstem/envs/latest/lib/python3.10/site-packages/numpy/__init__.py'>, gaussian_filter=<function gaussian_filter>)

kernel density estimate from a set coordinates (xa,ya) and intensity weights.

Parameters:
  • stack (np.ndarray) – Unshifted image stack, shape (N,P,S)

  • shifts (np.ndarray) – Shifts for each image in stack, shape: (N,2)

  • upsampling_factor (int) – Upsampling factor

  • kde_sigma (float) – KDE gaussian kernel bandwidth in upsampled pixels

  • lowpass_filter (bool, optional) – If True, the resulting KDE upsampled image is lowpass-filtered using a sinc-function

Returns:

pix_output – Upsampled intensity image

Return type:

np.ndarray

py4DSTEM.process.phase.utils.bilinear_kernel_density_estimate(xa, ya, intensities, output_shape, kde_sigma, lowpass_filter=False, xp=<module 'numpy' from '/home/docs/checkouts/readthedocs.org/user_builds/py4dstem/envs/latest/lib/python3.10/site-packages/numpy/__init__.py'>, gaussian_filter=<function gaussian_filter>)

kernel density estimate from a set coordinates (xa,ya) and intensity weights.

Parameters:
  • xa (np.ndarray) – Vertical positions of intensity array in pixels

  • ya (np.ndarray) – Horizontal positions of intensity array in pixels

  • intensities (np.ndarray) – Intensity array weights

  • output_shape ((int,int)) – Upsampled intensities shape

  • kde_sigma (float) – KDE gaussian kernel bandwidth in upsampled pixels

  • lowpass_filter (bool, optional) – If True, the resulting KDE upsampled image is lowpass-filtered using a sinc-function

Returns:

pix_output – Upsampled intensity image

Return type:

np.ndarray

py4DSTEM.process.phase.utils.lanczos_kernel_density_estimate(xa, ya, intensities, output_shape, kde_sigma, alpha, lowpass_filter=False, xp=<module 'numpy' from '/home/docs/checkouts/readthedocs.org/user_builds/py4dstem/envs/latest/lib/python3.10/site-packages/numpy/__init__.py'>, gaussian_filter=<function gaussian_filter>)

kernel density estimate from a set coordinates (xa,ya) and intensity weights.

Parameters:
  • xa (np.ndarray) – Vertical positions of intensity array in pixels

  • ya (np.ndarray) – Horizontal positions of intensity array in pixels

  • intensities (np.ndarray) – Intensity array weights

  • output_shape ((int,int)) – Upsampled intensities shape

  • kde_sigma (float) – KDE gaussian kernel bandwidth in upsampled pixels

  • alpha (int) – Lanczos kernel order

  • lowpass_filter (bool, optional) – If True, the resulting KDE upsampled image is lowpass-filtered using a sinc-function

Returns:

pix_output – Upsampled intensity image

Return type:

np.ndarray

py4DSTEM.process.phase.utils.vectorized_bilinear_resample(array, scale=None, output_size=None, mode='grid-wrap', grid_mode=True, xp=<module 'numpy' from '/home/docs/checkouts/readthedocs.org/user_builds/py4dstem/envs/latest/lib/python3.10/site-packages/numpy/__init__.py'>)

Resize an array along its final two axes. Note, this is vectorized and thus very memory-intensive.

The scaling of the array can be specified by passing either scale, which sets the scaling factor along both axes to be scaled; or by passing output_size, which specifies the final dimensions of the scaled axes.

Parameters:
  • array (np.ndarray) – Input array to be resampled

  • scale (float) – Scalar value giving the scaling factor for all dimensions

  • output_size ((int,int)) – Tuple of two values giving the output size for the final two axes

  • xp (Callable) – Array computing module

Returns:

resampled_array – Resampled array

Return type:

np.ndarray

py4DSTEM.process.phase.utils.vectorized_fourier_resample(array, scale=None, output_size=None, xp=<module 'numpy' from '/home/docs/checkouts/readthedocs.org/user_builds/py4dstem/envs/latest/lib/python3.10/site-packages/numpy/__init__.py'>)

Resize a 2D array along any dimension, using Fourier interpolation. For 4D input arrays, only the final two axes can be resized. Note, this is vectorized and thus very memory-intensive.

The scaling of the array can be specified by passing either scale, which sets the scaling factor along both axes to be scaled; or by passing output_size, which specifies the final dimensions of the scaled axes (and allows for different scaling along the x,y or kx,ky axes.)

Parameters:
  • array (np.ndarray) – Input 2D/4D array to be resampled

  • scale (float) – Scalar value giving the scaling factor for all dimensions

  • output_size ((int,int)) – Tuple of two values giving eith the (x,y) or (kx,ky) output size for 2D and 4D respectively.

  • xp (Callable) – Array computing module

Returns:

resampled_array – Resampled 2D/4D array

Return type:

np.ndarray

py4DSTEM.process.phase.utils.partition_list(lst, size)

Partitions lst into chunks of size. Returns a generator.

py4DSTEM.process.phase.utils.copy_to_device(array, device='cpu')

Copies array to device. Default allows one to use this as asnumpy()

probe

rdf

py4DSTEM.process.rdf.amorph.fit_stack(datacube, init_coefs, mask=None)

This will fit an ellipse using the polar elliptical transform code to all the diffraction patterns. It will take in a datacube and return a coefficient array which can then be used to map strain, fit the centers, etc.

Parameters:
  • datacute – a datacube of diffraction data

  • init_coefs – an initial starting guess for the fit

  • mask – a mask, either 2D or 4D, for either one mask for the whole stack, or one per pattern.

Returns:

an array of coefficients of the fit

py4DSTEM.process.rdf.amorph.calculate_coef_strain(coef_cube, r_ref)

This function will calculate the strains from a 3D matrix output by fit_stack

Coefs order:
  • I0 the intensity of the first gaussian function

  • I1 the intensity of the Janus gaussian

  • sigma0 std of first gaussian

  • sigma1 inner std of Janus gaussian

  • sigma2 outer std of Janus gaussian

  • c_bkgd a constant offset

  • R center of the Janus gaussian

  • x0,y0 the origin

  • B,C 1x^2 + Bxy + Cy^2 = 1

Parameters:
  • coef_cube – output from fit_stack

  • r_ref – a reference 0 strain radius - needed because we fit r as well as B and C

Returns:

  • exx: strain in the x axis direction in image coordinates

  • eyy: strain in the y axis direction in image coordinates

  • exy: shear

Return type:

(3-tuple) A 3-tuple containing

py4DSTEM.process.rdf.amorph.plot_strains(strains, cmap='RdBu_r', vmin=None, vmax=None, mask=None)

This function will plot strains with a unified color scale.

Parameters:
  • strains (3-tuple of arrays) – (exx, eyy, exy)

  • cmap – imshow parameters

  • vmin – imshow parameters

  • vmax – imshow parameters

  • mask – real space mask of values not to show (black)

py4DSTEM.process.rdf.amorph.convert_stack_polar(datacube, coef_cube)

This function will take the coef_cube from fit_stack and apply it to the image stack, to return polar transformed images.

Parameters:
  • datacube – data in datacube format

  • coef_cube – coefs from fit_stack

Returns:

polar transformed datacube

py4DSTEM.process.rdf.amorph.compute_polar_stack_symmetries(datacube_polar)

This function will take in a datacube of polar-transformed diffraction patterns, and do the autocorrelation, before taking the fourier transform along the theta direction, such that symmetries can be measured. They will be plotted by a different function

Parameters:

datacube_polar – diffraction pattern cube that has been polar transformed

Returns:

the normalized fft along the theta direction of the autocorrelated patterns in datacube_polar

py4DSTEM.process.rdf.amorph.plot_symmetries(datacube_symmetries, sym_order)

This function will take in a datacube from compute_polar_stack_symmetries and plot a specific symmetry order.

Parameters:
  • datacube_symmetries – result of compute_polar_stack_symmetries, the stack of fft’d autocorrelated diffraction patterns

  • sym_order – symmetry order desired to plot

Returns:

None

py4DSTEM.process.rdf.rdf.get_radial_intensity(polar_img, polar_mask)

Takes in a radial transformed image and the radial mask (if any) applied to that image. Designed to be compatible with polar-elliptical transforms from utils

py4DSTEM.process.rdf.rdf.fit_scattering_factor(scale, elements, composition, q_arr, units)

Scale is linear factor Elements is an 1D array of atomic numbers. Composition is a 1D array, same length as elements, describing the average atomic composition of the sample. If the Q_coords is a 1D array of Fourier coordinates, given in inverse Angstroms. Units is a string of ‘VA’ or ‘A’, which returns the scattering factor in volt angtroms or in angstroms.

py4DSTEM.process.rdf.rdf.get_phi(radialIntensity, scatter, q_arr)

ymean scale*scatter.fe**2

py4DSTEM.process.rdf.rdf.get_mask(left, right, midpoint, slopes, q_arr)

start is float stop is float midpoint is float slopes is [float,float]

py4DSTEM.process.rdf.rdf.get_rdf(phi, q_arr)

phi can be masked or not masked

utils

py4DSTEM.process.utils.cross_correlate.get_cross_correlation(ar, template, corrPower=1, _returnval='real')

Get the cross/phase/hybrid correlation of ar with template, where the latter is in real space.

If _returnval is ‘real’, returns the real-valued cross-correlation. Otherwise, returns the complex valued result.

py4DSTEM.process.utils.cross_correlate.get_cross_correlation_FT(ar, template_FT, corrPower=1, _returnval='real')

Get the cross/phase/hybrid correlation of ar with template_FT, where the latter is already in Fourier space (i.e. template_FT is np.conj(np.fft.fft2(template)).

If _returnval is ‘real’, returns the real-valued cross-correlation. Otherwise, returns the complex valued result.

py4DSTEM.process.utils.cross_correlate.get_shift(ar1, ar2, corrPower=1)

Determine the relative shift between a pair of arrays giving the best overlap.

Shift determination uses the brightest pixel in the cross correlation, and is

thus limited to pixel resolution. corrPower specifies the cross correlation power, with 1 corresponding to a cross correlation and 0 a phase correlation.

Args:

ar1,ar2 (2D ndarrays):

corrPower (float between 0 and 1, inclusive): 1=cross correlation, 0=phase

correlation

Returns:

(shiftx,shifty) - the relative image shift, in pixels

Return type:

(2-tuple)

py4DSTEM.process.utils.cross_correlate.align_images_fourier(G1, G2, upsample_factor, device='cpu')

Alignment of two images using DFT upsampling of cross correlation.

Parameters:
  • G1 (ndarray) – fourier transform of image 1

  • G2 (ndarray) – fourier transform of image 2

  • upsample_factor (float) – upsampling for correlation. Must be greater than 2.

  • device (str, optional) – calculation device will be perfomed on. Must be ‘cpu’ or ‘gpu’

  • Returns – xy_shift [pixels]

py4DSTEM.process.utils.cross_correlate.align_and_shift_images(image_1, image_2, upsample_factor, device='cpu')

Alignment of two images using DFT upsampling of cross correlation.

Parameters:
  • image_1 (ndarray) – image 1

  • image_2 (ndarray) – image 2

  • upsample_factor (float) – upsampling for correlation. Must be greater than 2.

  • device (str, optional) – calculation device will be perfomed on. Must be ‘cpu’ or ‘gpu’.

  • Returns – shifted image [pixels]

Contains functions relating to polar-elliptical calculations.

This includes
  • transforming data from cartesian to polar-elliptical coordinates

  • converting between ellipse representations

  • radial and polar-elliptical radial integration

Functions for measuring/fitting elliptical distortions are found in process/calibration/ellipse.py. Functions for computing radial and polar-elliptical radial backgrounds are found in process/preprocess/ellipse.py.

py4DSTEM uses 2 ellipse representations - one user-facing representation, and one internal representation. The user-facing represenation is in terms of the following 5 parameters:

x0,y0 the center of the ellipse a the semimajor axis length b the semiminor axis length theta the (positive, right handed) tilt of the a-axis

to the x-axis, in radians

Internally, fits are performed using the canonical ellipse parameterization, in terms of the parameters (x0,y0,A,B,C):

A(x-x0)^2 + B(x-x0)(y-y0) C(y-y0)^2 = 1

It is possible to convert between (a,b,theta) <–> (A,B,C) using the convert_ellipse_params() and convert_ellipse_params_r() methods.

Transformation from cartesian to polar-elliptical space is done using

x = x0 + a*r*cos(phi)*cos(theta) + b*r*sin(phi)*sin(theta) y = y0 + a*r*cos(phi)*sin(theta) - b*r*sin(phi)*cos(theta)

where (r,phi) are the polar-elliptical coordinates. All angular quantities are in radians.

py4DSTEM.process.utils.elliptical_coords.convert_ellipse_params(A, B, C)

Converts ellipse parameters from canonical form (A,B,C) into semi-axis lengths and tilt (a,b,theta). See module docstring for more info.

Parameters:
  • A (floats) – parameters of an ellipse in the form: Ax^2 + Bxy + Cy^2 = 1

  • B (floats) – parameters of an ellipse in the form: Ax^2 + Bxy + Cy^2 = 1

  • C (floats) – parameters of an ellipse in the form: Ax^2 + Bxy + Cy^2 = 1

Returns:

A 3-tuple consisting of:

  • a: (float) the semimajor axis length

  • b: (float) the semiminor axis length

  • theta: (float) the tilt of the ellipse semimajor axis with respect to the x-axis, in radians

Return type:

(3-tuple)

py4DSTEM.process.utils.elliptical_coords.convert_ellipse_params_r(a, b, theta)

Converts from ellipse parameters (a,b,theta) to (A,B,C). See module docstring for more info.

Parameters:
  • a (floats) – parameters of an ellipse, where a/b are the semimajor/semiminor axis lengths, and theta is the tilt of the semimajor axis with respect to the x-axis, in radians.

  • b (floats) – parameters of an ellipse, where a/b are the semimajor/semiminor axis lengths, and theta is the tilt of the semimajor axis with respect to the x-axis, in radians.

  • theta (floats) – parameters of an ellipse, where a/b are the semimajor/semiminor axis lengths, and theta is the tilt of the semimajor axis with respect to the x-axis, in radians.

Returns:

A 3-tuple consisting of (A,B,C), the ellipse parameters in

canonical form.

Return type:

(3-tuple)

py4DSTEM.process.utils.elliptical_coords.cartesian_to_polarelliptical_transform(cartesianData, p_ellipse, dr=1, dphi=0.03490658503988659, r_range=None, mask=None, maskThresh=0.99)

Transforms an array of data in cartesian coordinates into a data array in polar-elliptical coordinates.

Discussion of the elliptical parametrization used can be found in the docstring for the process.utils.elliptical_coords module.

Parameters:
  • cartesianData (2D float array) – the data in cartesian coordinates

  • p_ellipse (5-tuple) – specifies (qx0,qy0,a,b,theta), the parameters for the transformation. These are the same 5 parameters which are outputs of the elliptical fitting functions in the process.calibration module, e.g. fit_ellipse_amorphous_ring and fit_ellipse_1D. For more details, see the process.utils.elliptical_coords module docstring

  • dr (float) – sampling of the (r,phi) coords: the width of the bins in r

  • dphi (float) – sampling of the (r,phi) coords: the width of the bins in phi, in radians

  • r_range (number or length 2 list/tuple or None) –

    specifies the sampling of the (r,theta) coords. Precise behavior which depends on the parameter type:

    • if None, autoselects max r value

    • if r_range is a number, specifies the maximum r value

    • if r_range is a length 2 list/tuple, specifies the min/max r values

  • mask (2d array of bools) – shape must match cartesianData; where mask==False, ignore these datapoints in making the polarElliptical data array

  • maskThresh (float) – the final data mask is calculated by converting mask (above) from cartesian to polar elliptical coords. Due to interpolation, this results in some non-boolean values - this is converted back to a boolean array by taking polarEllipticalMask = polarTrans(mask) < maskThresh. Cells where polarTrans is less than 1 (i.e. has at least one masked NN) should generally be masked, hence the default value of 0.99.

Returns:

A 3-tuple, containing:

  • polarEllipticalData: (2D masked array) a masked array containing the data and the data mask, in polarElliptical coordinates

  • rr: (2D array) meshgrid of the r coordinates

  • pp: (2D array) meshgrid of the phi coordinates

Return type:

(3-tuple)

py4DSTEM.process.utils.elliptical_coords.elliptical_resample_datacube(datacube, p_ellipse, mask=None, maskThresh=0.99)

Perform elliptic resamplig on each diffraction pattern in a DataCube Detailed description of the args is found in elliptical_resample.

NOTE: Only use this function if you need to resample the raw data. If you only need for Bragg disk positions to be corrected, use the BraggVector calibration routines, as it is much faster to perform this on the peak positions than the entire datacube.

py4DSTEM.process.utils.elliptical_coords.elliptical_resample(data, p_ellipse, mask=None, maskThresh=0.99)

Resamples data with elliptic distortion to correct distortion of the input pattern.

Discussion of the elliptical parametrization used can be found in the docstring for the process.utils.elliptical_coords module.

Parameters:
  • data (2D float array) – the data in cartesian coordinates

  • p_ellipse (5-tuple) – specifies (qx0,qy0,a,b,theta), the parameters for the transformation. These are the same 5 parameters which are outputs of the elliptical fitting functions in the process.calibration module, e.g. fit_ellipse_amorphous_ring and fit_ellipse_1D. For more details, see the process.utils.elliptical_coords module docstring

  • dr (float) – sampling of the (r,phi) coords: the width of the bins in r

  • dphi (float) – sampling of the (r,phi) coords: the width of the bins in phi, in radians

  • r_range (number or length 2 list/tuple or None) –

    specifies the sampling of the (r,theta) coords. Precise behavior which depends on the parameter type:

    • if None, autoselects max r value

    • if r_range is a number, specifies the maximum r value

    • if r_range is a length 2 list/tuple, specifies the min/max r values

  • mask (2d array of bools) – shape must match cartesianData; where mask==False, ignore these datapoints in making the polarElliptical data array

  • maskThresh (float) – the final data mask is calculated by converting mask (above) from cartesian to polar elliptical coords. Due to interpolation, this results in some non-boolean values - this is converted back to a boolean array by taking polarEllipticalMask = polarTrans(mask) < maskThresh. Cells where polarTrans is less than 1 (i.e. has at least one masked NN) should generally be masked, hence the default value of 0.99.

Returns:

A 3-tuple, containing:

  • resampled_data: (2D masked array) a masked array containing the data and the data mask, in polarElliptical coordinates

Return type:

(3-tuple)

py4DSTEM.process.utils.elliptical_coords.radial_elliptical_integral(ar, dr, p_ellipse, rmax=None)

Computes the radial integral of array ar from center (x0,y0) with a step size in r of dr.

Parameters:
  • ar (2d array) – the data

  • dr (number) – the r sampling

  • p_ellipse (5-tuple) – the parameters (x0,y0,a,b,theta) for the ellipse

  • r_max (float) – maximum radial value

Returns:

A 2-tuple containing:

  • rbin_centers: (1d array) the bins centers of the radial integral

  • radial_integral: (1d array) the radial integral

radial_integral (1d array) the radial integral

Return type:

(2-tuple)

py4DSTEM.process.utils.elliptical_coords.radial_integral(ar, x0=None, y0=None, dr=0.1, rmax=None)

Computes the radial integral of array ar from center (x0,y0) with a step size in r of dr.

Parameters:
  • ar (2d array) – the data

  • x0 (floats) – the origin

  • y0 (floats) – the origin

  • dr (number) – radial step size

  • rmax (float) – maximum radial dimension

Returns:

A 2-tuple containing:

  • rbin_centers: (1d array) the bins centers of the radial integral

  • radial_integral: (1d array) the radial integral

Return type:

(2-tuple)

py4DSTEM.process.utils.masks.get_beamstop_mask(dp, qx0, qy0, theta, dtheta=1, w=10, r=10)

Generates a beamstop shaped mask.

Parameters:
  • dp (2d array) – a diffraction pattern

  • qx0 (numbers) – the center position of the beamstop

  • qy0 (numbers) – the center position of the beamstop

  • theta (number) – the orientation of the beamstop, in degrees

  • dtheta (number) – angular span of the wedge representing the beamstop, in degrees

  • w (integer) – half the width of the beamstop arm, in pixels

  • r (number) – the radius of a circle at the end of the beamstop, in pixels

Returns:

the mask

Return type:

(2d boolean array)

py4DSTEM.process.utils.masks.make_circular_mask(shape, qxy0, radius)

Create a hard circular mask, for use in DPC integration or or to use as a filter in diffraction or real space.

Parameters:
  • shape (2-tuple of ints) –

  • qxy0 (2-tuple of floats) center coordinates, in pixels. Must be in (row, column) –

  • radius (float) –

Returns:

mask (2D boolean array) the mask

loosely based on multicorr.py found at: https://github.com/ercius/openNCEM/blob/master/ncempy/algo/multicorr.py

modified by SEZ, May 2019 to integrate with py4DSTEM utility functions
  • rewrote upsampleFFT (previously did not work correctly)

  • modified upsampled_correlation to accept xyShift, the point around which to

upsample the DFT * eliminated the factor-2 FFT upsample step in favor of using parabolic for first-pass subpixel (since parabolic is so fast) * rewrote the matrix multiply DFT to be more pythonic

py4DSTEM.process.utils.multicorr.upsampled_correlation(imageCorr, upsampleFactor, xyShift, device='cpu')

Refine the correlation peak of imageCorr around xyShift by DFT upsampling.

There are two approaches to Fourier upsampling for subpixel refinement: (a) one can pad an (appropriately shifted) FFT with zeros and take the inverse transform, or (b) one can compute the DFT by matrix multiplication using modified transformation matrices. The former approach is straightforward but requires performing the FFT algorithm (which is fast) on very large data. The latter method trades one speedup for a slowdown elsewhere: the matrix multiply steps are expensive but we operate on smaller matrices. Since we are only interested in a very small region of the FT around a peak of interest, we use the latter method to get a substantial speedup and enormous decrease in memory requirement. This “DFT upsampling” approach computes the transformation matrices for the matrix- multiply DFT around a small 1.5px wide region in the original imageCorr.

Following the matrix multiply DFT we use parabolic subpixel fitting to get even more precision! (below 1/upsampleFactor pixels)

NOTE: previous versions of multiCorr operated in two steps: using the zero- padding upsample method for a first-pass factor-2 upsampling, followed by the DFT upsampling (at whatever user-specified factor). I have implemented it differently, to better support iterating over multiple peaks. The DFT is always upsampled around xyShift, which MUST be specified to HALF-PIXEL precision (no more, no less) to replicate the behavior of the factor-2 step. (It is possible to refactor this so that peak detection is done on a Fourier upsampled image rather than using the parabolic subpixel and rounding as now… I like keeping it this way because all of the parameters and logic will be identical to the other subpixel methods.)

Parameters:
  • imageCorr (complex valued ndarray) – Complex product of the FFTs of the two images to be registered i.e. m = np.fft.fft2(DP) * probe_kernel_FT; imageCorr = np.abs(m)**(corrPower) * np.exp(1j*np.angle(m))

  • upsampleFactor (int) – Upsampling factor. Must be greater than 2. (To do upsampling with factor 2, use upsampleFFT, which is faster.)

  • xyShift – Location in original image coordinates around which to upsample the FT. This should be given to exactly half-pixel precision to replicate the initial FFT step that this implementation skips

Returns:

Refined location of the peak in image coordinates.

Return type:

(2-element np array)

py4DSTEM.process.utils.multicorr.upsampleFFT(cc, device='cpu')

Zero-padding FFT upsampling. Returns the real IFFT of the input with 2x upsampling. This may have an error for matrices with an odd size. Takes a complex np array as input.

py4DSTEM.process.utils.multicorr.dftUpsample(imageCorr, upsampleFactor, xyShift, device='cpu')

This performs a matrix multiply DFT around a small neighboring region of the inital correlation peak. By using the matrix multiply DFT to do the Fourier upsampling, the efficiency is greatly improved. This is adapted from the subfuction dftups found in the dftregistration function on the Matlab File Exchange.

https://www.mathworks.com/matlabcentral/fileexchange/18401-efficient-subpixel-image-registration-by-cross-correlation

The matrix multiplication DFT is from:

Manuel Guizar-Sicairos, Samuel T. Thurman, and James R. Fienup, “Efficient subpixel image registration algorithms,” Opt. Lett. 33, 156-158 (2008). http://www.sciencedirect.com/science/article/pii/S0045790612000778

Parameters:
  • imageCorr (complex valued ndarray) – Correlation image between two images in Fourier space.

  • upsampleFactor (int) – Scalar integer of how much to upsample.

  • xyShift (list of 2 floats) – Coordinates in the UPSAMPLED GRID around which to upsample. These must be single-pixel IN THE UPSAMPLED GRID

Returns:

Upsampled image from region around correlation peak.

Return type:

(ndarray)

py4DSTEM.process.utils.utils.radial_reduction(ar, x0, y0, binsize=1, fn=<function mean>, coords=None)

Evaluate a reduction function on pixels within annular rings centered on (x0,y0), with a ring width of binsize.

By default, returns the mean value of pixels within each annulus. Some other useful reductions include: np.sum, np.std, np.count, np.median, …

When running in a loop, pre-compute the pixel coordinates and pass them in for improved performance, like so:

coords = np.mgrid[0:ar.shape[0],0:ar.shape[1]] radial_sums = radial_reduction(ar, x0,y0, coords=coords)

py4DSTEM.process.utils.utils.sector_mask(shape, centre, radius, angle_range=(0, 360))

Return a boolean mask for a circular sector. The start/stop angles in angle_range should be given in clockwise order.

Parameters:
  • shape – 2D shape of the mask

  • centre – 2D center of the circular sector

  • radius – radius of the circular mask

  • angle_range – angular range of the circular mask

py4DSTEM.process.utils.utils.get_qx_qy_1d(M, dx=[1, 1], fft_shifted=False)

Generates 1D Fourier coordinates for a (Nx,Ny)-shaped 2D array. Specifying the dx argument sets a unit size.

Parameters:
  • M – (2,) shape of the returned array

  • dx – (2,) tuple, pixel size

  • fft_shifted – True if result should be fft_shifted to have the origin in the center of the array

py4DSTEM.process.utils.utils.make_Fourier_coords2D(Nx, Ny, pixelSize=1)
Generates Fourier coordinates for a (Nx,Ny)-shaped 2D array.

Specifying the pixelSize argument sets a unit size.

py4DSTEM.process.utils.utils.get_CoM(ar, device='cpu', corner_centered=False)

Finds and returns the center of mass of array ar. If corner_centered is True, uses fftfreq for indices.

py4DSTEM.process.utils.utils.get_maxima_1D(ar, sigma=0, minSpacing=0, minRelativeIntensity=0, relativeToPeak=0)

Finds the indices where 1D array ar is a local maximum. Optional parameters allow blurring the array and filtering the output; setting each to 0 (default) turns off these functions.

Parameters:
  • ar (1D array) –

  • sigma (number) – gaussian blur std to apply to ar before finding maxima

  • minSpacing (number) – if two maxima are found within minSpacing, the dimmer one is removed

  • minRelativeIntensity (number) – maxima dimmer than minRelativeIntensity compared to the relativeToPeak’th brightest maximum are removed

  • relativeToPeak (int) – 0=brightest maximum. 1=next brightest, etc.

Returns:

An array of indices where ar is a local maximum, sorted by intensity.

Return type:

(array of ints)

py4DSTEM.process.utils.utils.linear_interpolation_1D(ar, x)

Calculates the 1D linear interpolation of array ar at position x using the two nearest elements.

py4DSTEM.process.utils.utils.add_to_2D_array_from_floats(ar, x, y, I)

Adds the values I to array ar, distributing the value between the four pixels nearest (x,y) using linear interpolation. Inputs (x,y,I) may be floats or arrays of floats.

Note that if the same [x,y] coordinate appears more than once in the input array, only the final value of I at that coordinate will get added.

py4DSTEM.process.utils.utils.get_voronoi_vertices(voronoi, nx, ny, dist=10)

From a scipy.spatial.Voronoi instance, return a list of ndarrays, where each array is shape (N,2) and contains the (x,y) positions of the vertices of a voronoi region.

The problem this function solves is that in a Voronoi instance, some vertices outside the field of view of the tesselated region are left unspecified; only the existence of a point beyond the field is referenced (which may or may not be ‘at infinity’). This function specifies all points, such that the vertices and edges of the tesselation may be directly laid over data.

Parameters:
  • voronoi (scipy.spatial.Voronoi) – the voronoi tesselation

  • nx (int) – the x field-of-view of the tesselated region

  • ny (int) – the y field-of-view of the tesselated region

  • dist (float, optional) – place new vertices by extending new voronoi edges outside the frame by a distance of this factor times the distance of its known vertex from the frame edge

Returns:

the (x,y) coords of the vertices of each voronoi region

Return type:

(list of ndarrays of shape (N,2))

py4DSTEM.process.utils.utils.get_ewpc_filter_function(Q_Nx, Q_Ny)

Returns a function for computing the exit wave power cepstrum of a diffraction pattern using a Hanning window. This can be passed as the filter_function in the Bragg disk detection functions (with the probe an array of ones) to find the lattice vectors by the EWPC method (but be careful as the lengths are now in realspace units!) See https://arxiv.org/abs/1911.00984

py4DSTEM.process.utils.utils.fourier_resample(array, scale=None, output_size=None, force_nonnegative=False, bandlimit_nyquist=None, bandlimit_power=2, dtype=<class 'numpy.float32'>)

Resize a 2D array along any dimension, using Fourier interpolation / extrapolation. For 4D input arrays, only the final two axes can be resized.

The scaling of the array can be specified by passing either scale, which sets the scaling factor along both axes to be scaled; or by passing output_size, which specifies the final dimensions of the scaled axes (and allows for different scaling along the x,y or kx,ky axes.)

Parameters:
  • array (2D/4D numpy array) – Input array, or 4D stack of arrays, to be resized.

  • scale (float) – scalar value giving the scaling factor for all dimensions

  • output_size (2-tuple of ints) – two values giving either the (x,y) output size for 2D, or (kx,ky) for 4D

  • force_nonnegative (bool) – Force all outputs to be nonnegative, after filtering

  • bandlimit_nyquist (float) – Gaussian filter information limit in Nyquist units (0.5 max in both directions)

  • bandlimit_power (float) – Gaussian filter power law scaling (higher is sharper)

  • dtype (numpy dtype) – datatype for binned array. default is single precision float

Returns:

the resized array (2D/4D numpy array)

virtualdiffraction

virtualimage

wholepatternfit

class py4DSTEM.process.wholepatternfit.wp_models.WPFModelType(value)

Flags to signify capabilities and other semantics of a Model

class py4DSTEM.process.wholepatternfit.wp_models.WPFModel(name: str, params: dict, model_type=WPFModelType.DUMMY)

Prototype class for a compent of a whole-pattern model. Holds the following:

name: human-readable name of the model params: a dict of names and initial (or returned) values of the model parameters func: a function that takes as arguments:

  • the diffraction pattern being built up, which the function should modify in place

  • positional arguments in the same order as the params dictionary

  • keyword arguments. this is to provide some pre-computed information for convenience
    kwargs will include:
    • xArray, yArray meshgrid of the x and y coordinates

    • global_x0 global x-coordinate of the pattern center

    • global_y0 global y-coordinate of the pattern center

jacobian: a function that takes as arguments:
  • the diffraction pattern being built up, which the function should modify in place

  • positional arguments in the same order as the params dictionary

  • offset: the first index (j) that values should be written into

    (the function should ONLY write into 0,1, and offset:offset+nParams) 0 and 1 are the entries for global_x0 and global_y0, respectively REMEMBER TO ADD TO 0 and 1 SINCE ALL MODELS CAN CONTRIBUTE TO THIS PARTIAL DERIVATIVE

  • keyword arguments. this is to provide some pre-computed information for convenience

__init__(name: str, params: dict, model_type=WPFModelType.DUMMY)
class py4DSTEM.process.wholepatternfit.wp_models.DCBackground(background_value=0.0, name='DC Background')

Model representing constant background intensity.

Parameters:

background_value

Background intensity value. Specified as initial_value, (initial_value, deviation), or

(initial_value, lower_bound, upper_bound). See Parameter documentation for details.

__init__(background_value=0.0, name='DC Background')
class py4DSTEM.process.wholepatternfit.wp_models.GaussianBackground(WPF, sigma, intensity, global_center=True, x0=0.0, y0=0.0, name='Gaussian Background')

Model representing a 2D Gaussian intensity distribution

Parameters:
  • WPF (WholePatternFit) – Parent WPF object

  • sigma

    parameter specifying width of the Gaussian Specified as initial_value, (initial_value, deviation), or

    (initial_value, lower_bound, upper_bound). See Parameter documentation for details.

  • intensity

    parameter specifying intensity of the Gaussian Specified as initial_value, (initial_value, deviation), or

    (initial_value, lower_bound, upper_bound). See Parameter documentation for details.

  • global_center (bool) – If True, uses same center coordinate as the global model If False, uses an independent center

  • x0

    Center coordinates of model for local origin Specified as initial_value, (initial_value, deviation), or

    (initial_value, lower_bound, upper_bound). See Parameter documentation for details.

  • y0

    Center coordinates of model for local origin Specified as initial_value, (initial_value, deviation), or

    (initial_value, lower_bound, upper_bound). See Parameter documentation for details.

__init__(WPF, sigma, intensity, global_center=True, x0=0.0, y0=0.0, name='Gaussian Background')
class py4DSTEM.process.wholepatternfit.wp_models.GaussianRing(WPF, radius, sigma, intensity, global_center=True, x0=0.0, y0=0.0, name='Gaussian Ring')

Model representing a halo with Gaussian falloff

Parameters:
  • WPF (WholePatternFit) – parent fitting object

  • radius

    radius of halo Specified as initial_value, (initial_value, deviation), or

    (initial_value, lower_bound, upper_bound). See Parameter documentation for details.

  • sigma

    width of Gaussian falloff Specified as initial_value, (initial_value, deviation), or

    (initial_value, lower_bound, upper_bound). See Parameter documentation for details.

  • intensity

    Intensity of the halo Specified as initial_value, (initial_value, deviation), or

    (initial_value, lower_bound, upper_bound). See Parameter documentation for details.

  • global_center (bool) – If True, uses same center coordinate as the global model If False, uses an independent center

  • x0

    Center coordinates of model for local origin Specified as initial_value, (initial_value, deviation), or

    (initial_value, lower_bound, upper_bound). See Parameter documentation for details.

  • y0

    Center coordinates of model for local origin Specified as initial_value, (initial_value, deviation), or

    (initial_value, lower_bound, upper_bound). See Parameter documentation for details.

__init__(WPF, radius, sigma, intensity, global_center=True, x0=0.0, y0=0.0, name='Gaussian Ring')
class py4DSTEM.process.wholepatternfit.wp_models.SyntheticDiskLattice(WPF, ux: float, uy: float, vx: float, vy: float, disk_radius: float, disk_width: float, u_max: int, v_max: int, intensity_0: float, refine_radius: bool = False, refine_width: bool = False, global_center: bool = True, x0: float = 0.0, y0: float = 0.0, exclude_indices: list = [], include_indices: list | None = None, name='Synthetic Disk Lattice', verbose=False)

Model representing a lattice of diffraction disks with a soft edge

Parameters:
  • WPF (WholePatternFit) – parent fitting object

  • ux

    x and y components of the lattice vectors u and v. Specified as initial_value, (initial_value, deviation), or

    (initial_value, lower_bound, upper_bound). See Parameter documentation for details.

  • uy

    x and y components of the lattice vectors u and v. Specified as initial_value, (initial_value, deviation), or

    (initial_value, lower_bound, upper_bound). See Parameter documentation for details.

  • vx

    x and y components of the lattice vectors u and v. Specified as initial_value, (initial_value, deviation), or

    (initial_value, lower_bound, upper_bound). See Parameter documentation for details.

  • vy

    x and y components of the lattice vectors u and v. Specified as initial_value, (initial_value, deviation), or

    (initial_value, lower_bound, upper_bound). See Parameter documentation for details.

  • disk_radius

    Radius of each diffraction disk. Specified as initial_value, (initial_value, deviation), or

    (initial_value, lower_bound, upper_bound). See Parameter documentation for details.

  • disk_width

    Width of the smooth falloff at the edge of the disk Specified as initial_value, (initial_value, deviation), or

    (initial_value, lower_bound, upper_bound). See Parameter documentation for details.

  • u_max – Maximum lattice indices to include in the pattern. Disks outside the pattern are automatically clipped.

  • v_max – Maximum lattice indices to include in the pattern. Disks outside the pattern are automatically clipped.

  • intensity_0

    Initial intensity for each diffraction disk. Each disk intensity is an independent fit variable in the final model Specified as initial_value, (initial_value, deviation), or

    (initial_value, lower_bound, upper_bound). See Parameter documentation for details.

  • refine_radius (bool) – Flag whether disk radius is made a fitting parameter

  • refine_width (bool) – Flag whether disk edge width is made a fitting parameter

  • global_center (bool) – If True, uses same center coordinate as the global model If False, uses an independent center

  • x0

    Center coordinates of model for local origin Specified as initial_value, (initial_value, deviation), or

    (initial_value, lower_bound, upper_bound). See Parameter documentation for details.

  • y0

    Center coordinates of model for local origin Specified as initial_value, (initial_value, deviation), or

    (initial_value, lower_bound, upper_bound). See Parameter documentation for details.

  • exclude_indices (list) – Indices to exclude from the pattern

  • include_indices (list) – If specified, only the indices in the list are added to the pattern

__init__(WPF, ux: float, uy: float, vx: float, vy: float, disk_radius: float, disk_width: float, u_max: int, v_max: int, intensity_0: float, refine_radius: bool = False, refine_width: bool = False, global_center: bool = True, x0: float = 0.0, y0: float = 0.0, exclude_indices: list = [], include_indices: list | None = None, name='Synthetic Disk Lattice', verbose=False)
class py4DSTEM.process.wholepatternfit.wp_models.SyntheticDiskMoire(WPF, lattice_a: SyntheticDiskLattice, lattice_b: SyntheticDiskLattice, intensity_0: float, decorated_peaks: list | None = None, link_moire_disk_intensities: bool = False, link_disk_parameters: bool = True, refine_width: bool = True, edge_width: list | None = None, refine_radius: bool = True, disk_radius: list | None = None, name: str = 'Moire Lattice')

Model of diffraction disks arising from interference between two lattices.

The Moire unit cell is determined automatically using the two input lattices.

Parameters:
  • WPF (WholePatternFit) – parent fitting object

  • lattice_a (SyntheticDiskLattice) – parent lattices for the Moire

  • lattice_b (SyntheticDiskLattice) – parent lattices for the Moire

  • intensity_0

    Initial guess of Moire disk intensity Specified as initial_value, (initial_value, deviation), or

    (initial_value, lower_bound, upper_bound). See Parameter documentation for details.

  • decorated_peaks (list) – When specified, only the reflections in the list are decorated with Moire spots If not specified, all peaks are decorated

  • link_moire_disk_intensities (bool) – When False, each Moire disk has an independently fit intensity When True, Moire disks arising from the same order of parent reflection share the same intensity

  • link_disk_parameters (bool) – When True, edge_width and disk_radius are inherited from lattice_a

  • refine_width (bool) – Flag whether disk edge width is a fit variable

  • edge_width

    Width of the soft edge of the diffraction disk. Specified as initial_value, (initial_value, deviation), or

    (initial_value, lower_bound, upper_bound). See Parameter documentation for details.

  • refine_radius (bool) – Flag whether disk radius is a fit variable

  • radius (disk) –

    Radius of the diffraction disks Specified as initial_value, (initial_value, deviation), or

    (initial_value, lower_bound, upper_bound). See Parameter documentation for details.

__init__(WPF, lattice_a: SyntheticDiskLattice, lattice_b: SyntheticDiskLattice, intensity_0: float, decorated_peaks: list | None = None, link_moire_disk_intensities: bool = False, link_disk_parameters: bool = True, refine_width: bool = True, edge_width: list | None = None, refine_radius: bool = True, disk_radius: list | None = None, name: str = 'Moire Lattice')
class py4DSTEM.process.wholepatternfit.wp_models.ComplexOverlapKernelDiskLattice(WPF, probe_kernel: ndarray, ux: float, uy: float, vx: float, vy: float, u_max: int, v_max: int, intensity_0: float, exclude_indices: list = [], global_center: bool = True, x0=0.0, y0=0.0, name='Complex Overlapped Disk Lattice', verbose=False)
__init__(WPF, probe_kernel: ndarray, ux: float, uy: float, vx: float, vy: float, u_max: int, v_max: int, intensity_0: float, exclude_indices: list = [], global_center: bool = True, x0=0.0, y0=0.0, name='Complex Overlapped Disk Lattice', verbose=False)
class py4DSTEM.process.wholepatternfit.wp_models.KernelDiskLattice(WPF, probe_kernel: ndarray, ux: float, uy: float, vx: float, vy: float, u_max: int, v_max: int, intensity_0: float, exclude_indices: list = [], global_center: bool = True, x0=0.0, y0=0.0, name='Custom Kernel Disk Lattice', verbose=False)
__init__(WPF, probe_kernel: ndarray, ux: float, uy: float, vx: float, vy: float, u_max: int, v_max: int, intensity_0: float, exclude_indices: list = [], global_center: bool = True, x0=0.0, y0=0.0, name='Custom Kernel Disk Lattice', verbose=False)
py4DSTEM.process.wholepatternfit.wpf_viz.show_lattice_points(self, im=None, vmin=None, vmax=None, power=None, show_vectors=True, crop_to_pattern=False, returnfig=False, moire_origin_idx=[0, 0, 0, 0], *args, **kwargs)

Plotting utility to show the initial lattice points.

Parameters:
  • im (np.ndarray) – Optional: Image to show, defaults to mean CBED

  • vmin (float) – Intensity ranges for plotting im

  • vmax (float) – Intensity ranges for plotting im

  • power (float) – Gamma level for showing im

  • show_vectors (bool) – Flag to plot the lattice vectors

  • crop_to_pattern (bool) – Flag to limit the field of view to the pattern area. If False, spots outside the pattern are shown

  • returnfig (bool) – If True, (fig,ax) are returned and plt.show() is not called

  • moire_origin_idx (list of length 4) – Indices of peak on which to draw Moire vectors, written as [a_u, a_v, b_u, b_v]

  • args – Passed to plt.subplots

  • kwargs – Passed to plt.subplots

Returns:

fig,ax

Return type:

If returnfig=True