pyextal.larbed.empad ==================== .. py:module:: pyextal.larbed.empad .. autoapi-nested-parse:: Raw EMPAD I/O, virtual-aperture integration, and reciprocal-grid construction. These are the low-level routines behind :class:`pyextal.larbed.LARBEDAnalysis`. They operate on a 4D-STEM stack indexed ``matrix[scan_x, scan_y, det_x, det_y]``, where the first two axes are tilt/scan positions and the last two are detector pixels. Note on axis naming: throughout this module ``x`` is the *first* (slow) axis of the detector and ``y`` the second, matching the original PyLARBED convention. It is not the screen-space x/y of a plotted image. Attributes ---------- .. autoapisummary:: pyextal.larbed.empad.DEFAULT_DETECTOR_SHAPE Functions --------- .. autoapisummary:: pyextal.larbed.empad._in_bounds pyextal.larbed.empad.ReadRaw pyextal.larbed.empad.FEM_StackAverage pyextal.larbed.empad._aperture_pixels pyextal.larbed.empad.BackgroundFitWeights pyextal.larbed.empad.ApertureSum pyextal.larbed.empad.ApertureSumVariance pyextal.larbed.empad._correlation_matrix pyextal.larbed.empad.EstimatePixelNoise pyextal.larbed.empad.AlignByBeam pyextal.larbed.empad.getGridVectors pyextal.larbed.empad.getGridVectorsHOLZ pyextal.larbed.empad.find_nearest_peak pyextal.larbed.empad.getGridVectors_Adjust pyextal.larbed.empad.AnnularSum pyextal.larbed.empad.AnnularVar pyextal.larbed.empad.Precession pyextal.larbed.empad.Precession_Angled pyextal.larbed.empad.radial_profile pyextal.larbed.empad.get_gvectors Module Contents --------------- .. py:data:: DEFAULT_DETECTOR_SHAPE :value: (124, 124) .. py:function:: _in_bounds(x, y, shape) Checks that a detector coordinate is inside the usable pixel range. A one-pixel border is excluded on every side, so a (124, 124) detector admits coordinates 1..122 inclusive. :param x: Coordinate along the first detector axis. :type x: float :param y: Coordinate along the second detector axis. :type y: float :param shape: Detector shape as (n_x, n_y). :type shape: tuple[int, int] :returns: True if the coordinate is usable. :rtype: bool .. py:function:: ReadRaw(fname, type=0) Reads an EMPAD ``.raw`` file into a 4D array. :param fname: Path to the ``.raw`` file. :type fname: str | pathlib.Path :param type: Detector layout. ``0`` reads the standard 128x130 frame and crops the two trailing metadata rows to 124x124. ``1`` reads an already-cropped 124x124 frame. Defaults to 0. :type type: int, optional :returns: A (scan_x, scan_y, det_x, det_y) array with NaNs zeroed, negatives clipped to zero, and a 1e-10 floor applied so logarithms are safe. :rtype: np.ndarray .. py:function:: FEM_StackAverage(matrix, center_x, center_y, Avg_img=None, rot_average=False) Computes a fluctuation electron microscopy variance image over the stack. :param matrix: 4D diffraction stack. :type matrix: np.ndarray :param center_x: Direct-beam position along the first detector axis. :type center_x: int :param center_y: Direct-beam position along the second detector axis. :type center_y: int :param Avg_img: Stack to average for the normalisation image. Defaults to None. :type Avg_img: np.ndarray, optional :param rot_average: Add a rotationally averaged contribution. Defaults to False. :type rot_average: bool, optional :returns: The FEM image. :rtype: np.ndarray .. py:function:: _aperture_pixels(center_x, center_y, radius, shape) Lists the on-detector pixels of a square aperture. This is the one aperture definition shared by :func:`ApertureSum`, :func:`ApertureSumVariance` and :func:`BackgroundFitWeights`. :param center_x: Aperture centre along the first detector axis. :type center_x: int :param center_y: Aperture centre along the second detector axis. :type center_y: int :param radius: Aperture half-width in pixels. :type radius: int :param shape: Detector shape as (n_x, n_y). :type shape: tuple[int, int] :returns: Pixel coordinates, in summation order. :rtype: list[tuple[int, int]] .. py:function:: BackgroundFitWeights(center_x, center_y, radius, ri, ro, order=2, shape=DEFAULT_DETECTOR_SHAPE, exclude=()) Linear weights giving the fitted background summed under an aperture. A polynomial surface of the given order is fitted by unweighted least squares to the pixels of the annulus ``ri <= d <= ro`` that are neither in the aperture itself nor within ``radius`` (square half-width) of a node in `exclude`. Because the fit is linear, the fitted surface summed over the aperture pixels is ``data[fit pixels] @ w`` for a fixed `w`, so applying it to the stack fits every frame independently and exactly at once. .. note:: On a symmetric annulus the plane terms cancel at the centre, so ``order=1`` only differs from ``order=0`` when masking or the detector edge makes the fit region asymmetric. ``order=2`` corrects the curvature of a neighbour's tail, but only if the annulus spans a range of radii: on a 1-px ring the isotropic curvature is nearly indistinguishable from the constant and the fit extrapolates noisily. The outer radius should still stay clear of the neighbouring cores, where a quadratic no longer describes the tail. For a Lorentzian tail from a reflection 10-20 px away, ``ri=3`` with an annulus of ``(4, 6)`` cut the bias of the annulus mean by 10-30x, while ``(4, 8)`` overshot when the neighbour was only 10 px away. :param center_x: Aperture centre along the first detector axis. :type center_x: int :param center_y: Aperture centre along the second detector axis. :type center_y: int :param radius: Aperture half-width in pixels. :type radius: int :param ri: Inner radius of the fit annulus in pixels. :type ri: int :param ro: Outer radius of the fit annulus in pixels. :type ro: int :param order: Polynomial order, 0 (constant), 1 (plane) or 2 (quadratic). Defaults to 2. :type order: int, optional :param shape: Detector shape as (n_x, n_y). Defaults to ``DEFAULT_DETECTOR_SHAPE``. :type shape: tuple[int, int], optional :param exclude: Pixel positions of other reflections whose surroundings are left out of the fit. Defaults to (). :type exclude: Iterable[tuple[float, float]], optional :returns: ``(fi, fj, w, gain)`` -- the fit-pixel coordinates along each detector axis, their weights, and the noise gain ``||w|| / (count / sqrt(n_fit))``, which is 1 for a plain mean over the same pixels and grows as the fit extrapolates. :rtype: tuple[np.ndarray, np.ndarray, np.ndarray, float] :raises ValueError: If `order` is not 0, 1 or 2, or the fit region has too few pixels to determine the polynomial. .. py:function:: ApertureSum(matrix, center_x, center_y, radius=0, bkg=(0, 0), fit=None, clip=True) Sums counts within a square aperture, optionally subtracting a background. :param matrix: 4D diffraction stack. :type matrix: np.ndarray :param center_x: Aperture centre along the first detector axis. :type center_x: int :param center_y: Aperture centre along the second detector axis. :type center_y: int :param radius: Aperture half-width in pixels. ``0`` gives a single-pixel (point) virtual detector. Defaults to 0. :type radius: int, optional :param bkg: Inner and outer radii of an annulus whose mean intensity is subtracted. ``(0, 0)`` disables background subtraction. Ignored when `fit` is given. Defaults to (0, 0). :type bkg: tuple[int, int], optional :param fit: ``(fi, fj, w, ...)`` from :func:`BackgroundFitWeights`; subtracts the fitted background instead of the annulus mean. Defaults to None. :type fit: tuple, optional :param clip: Set negative results to zero. Defaults to True. :type clip: bool, optional :returns: The background-subtracted rocking curve and the subtracted background, both shaped (scan_x, scan_y). :rtype: tuple[np.ndarray, np.ndarray] .. py:function:: ApertureSumVariance(matrix, center_x, center_y, g=1, m=1, DQE=1, radius=0, bkg=(0, 0), varB=0.0, fit=None, rho=None) Propagates counting variance through the same aperture as :func:`ApertureSum`. The per-pixel variance follows the detector noise model used elsewhere in pyextal (see :func:`pyextal.LucyRichardson.sgnoise`): a signal-dependent term ``g*m*I/DQE`` plus a constant readout term ``g*m*varB``. Summing over ``count`` independent pixels scales the constant term by ``count``. Charge sharing makes neighbouring pixels' noise correlated. Given the correlation kernel `rho` (see :func:`EstimatePixelNoise`), the variance is instead the full quadratic form ``u^T C u``, with ``u = +1`` on the aperture, ``-w`` on the fit pixels and ``C_kl = rho(k - l) * sigma_k * sigma_l``. This includes the covariance between the aperture and the fit pixels next to it, which the independent-pixel sum leaves out. :param matrix: 4D diffraction stack. :type matrix: np.ndarray :param center_x: Aperture centre along the first detector axis. :type center_x: int :param center_y: Aperture centre along the second detector axis. :type center_y: int :param g: Detector gain. Defaults to 1. :type g: float, optional :param m: Mixing factor. Defaults to 1. :type m: float, optional :param DQE: Detective quantum efficiency. Defaults to 1. :type DQE: float, optional :param radius: Aperture half-width in pixels. Defaults to 0. :type radius: int, optional :param bkg: Inner and outer radii of the background annulus, whose variance is propagated into the result. Ignored when `fit` is given. Defaults to (0, 0). :type bkg: tuple[int, int], optional :param varB: Per-pixel background (readout) variance. Defaults to 0.0, which reproduces the original PyLARBED behaviour. :type varB: float, optional :param fit: ``(fi, fj, w, ...)`` from :func:`BackgroundFitWeights`. The same noise model is propagated through the fit weights, ``sum(w**2 * var_pixel)``; the fit pixels are disjoint from the aperture, so without `rho` the two terms add. Defaults to None. :type fit: tuple, optional :param rho: Inter-pixel noise correlation, a ``(2R+1, 2R+1)`` kernel indexed by pixel offset with ``rho[R, R] = 1``. Pairs further apart than ``R`` are taken as uncorrelated. Needs either `fit` or no background, because the annulus mean's variance is estimated from its spatial spread rather than propagated. Defaults to None (independent pixels). :type rho: np.ndarray, optional :returns: The variance of the rocking curve, shaped (scan_x, scan_y). :rtype: np.ndarray :raises ValueError: If `rho` is combined with the annulus-mean background. .. py:function:: _correlation_matrix(pixels, rho) Pairwise noise correlation between pixels, looked up in a kernel. :param pixels: Pixel coordinates, shaped (n, 2). :type pixels: np.ndarray :param rho: Correlation kernel, shaped (2R+1, 2R+1). :type rho: np.ndarray :returns: The (n, n) correlation matrix; zero beyond the kernel. :rtype: np.ndarray .. py:function:: EstimatePixelNoise(matrix, max_offset=3, counts=(5, 100), axis=1, stride=4, n_samples=400000, seed=0) Measures the per-pixel noise model and the inter-pixel noise correlation. Neighbouring scan positions see nearly the same pattern, so the second difference ``I[s-1] - 2 I[s] + I[s+1]`` along the scan removes the signal and leaves noise with six times the per-pixel variance. Its spread against intensity gives ``var = gm * I``, and its correlation between nearby detector pixels gives ``rho``. Only pixels whose mean is between ``counts[0]`` and ``counts[1]`` effective counts (``I / gm``) are used. Below that, the zero-clipping in :func:`ReadRaw` truncates the noise; above it, real changes of the rocking curve between scan positions leak into the second difference. .. note:: Difference along the fast scan axis (the default, ``axis=1``), whose neighbours are adjacent in time. On EMPAD scans the slow axis carries line-to-line jumps that are not counting noise. :param matrix: 4D diffraction stack (scan_x, scan_y, det_x, det_y). :type matrix: np.ndarray :param max_offset: Largest pixel offset ``R`` of the correlation kernel. Defaults to 3. :type max_offset: int, optional :param counts: Range of effective counts per pixel used. Defaults to (5, 100). :type counts: tuple[float, float], optional :param axis: Scan axis to difference along. Defaults to 1. :type axis: int, optional :param stride: Use every `stride`-th scan line across the other axis, to bound the cost. Defaults to 4. :type stride: int, optional :param n_samples: Random samples per chunk for the variance-versus-intensity fit. Defaults to 400000. :type n_samples: int, optional :param seed: Seed for that subsampling. Defaults to 0. :type seed: int, optional :returns: ``gm`` (per-pixel variance over intensity, i.e. ``g*m``), ``rho`` (the ``(2R+1, 2R+1)`` correlation kernel, ``rho[R, R] = 1``), ``gain_sum`` (``gm * rho.sum()``, the variance per unit intensity of a sum over a region much larger than the kernel), and ``n_pixels`` (the number of pixel samples behind ``rho``). :rtype: dict :raises ValueError: If no intensity bin falls in the `counts` range. .. py:function:: AlignByBeam(matrix, x, y, crop=4, reiterate=1) Aligns every pattern in the stack on the beam nearest a given position. :param matrix: 4D diffraction stack. :type matrix: np.ndarray :param x: Target position along the first detector axis. :type x: int :param y: Target position along the second detector axis. :type y: int :param crop: Half-width of the search window around (x, y). Defaults to 4. :type crop: int, optional :param reiterate: Number of alignment passes. Defaults to 1. :type reiterate: int, optional :returns: The aligned stack. :rtype: np.ndarray .. py:function:: getGridVectors(x0, y0, x2, y2, x1, y1, order=3, shape=DEFAULT_DETECTOR_SHAPE) Builds a reciprocal-lattice grid from two basis vectors. :param x0: Direct-beam position along the first detector axis. :type x0: float :param y0: Direct-beam position along the second detector axis. :type y0: float :param x2: First basis-vector tip, first axis. :type x2: float :param y2: First basis-vector tip, second axis. :type y2: float :param x1: Second basis-vector tip, first axis. :type x1: float :param y1: Second basis-vector tip, second axis. :type y1: float :param order: Grid extent, giving indices -order..order along each basis vector. Defaults to 3. :type order: int, optional :param shape: Detector shape used to reject off-detector nodes. Defaults to (124, 124). :type shape: tuple[int, int], optional :returns: ``((i, j), (x, y))`` pairs of lattice index and pixel position. :rtype: list .. py:function:: getGridVectorsHOLZ(x0, y0, x1, y1, x2, y2, x3, y3, x4, y4, order1=3, order2=3, order3=3, order4=3, shape=DEFAULT_DETECTOR_SHAPE) Builds a reciprocal-lattice grid from four basis vectors, including HOLZ. :param x0: Direct-beam position along the first detector axis. :type x0: float :param y0: Direct-beam position along the second detector axis. :type y0: float :param x1: Tips of the four basis vectors. :type x1: float :param y1: Tips of the four basis vectors. :type y1: float :param x2: Tips of the four basis vectors. :type x2: float :param y2: Tips of the four basis vectors. :type y2: float :param x3: Tips of the four basis vectors. :type x3: float :param y3: Tips of the four basis vectors. :type y3: float :param x4: Tips of the four basis vectors. :type x4: float :param y4: Tips of the four basis vectors. :type y4: float :param order1: Extent along the first basis vector. Defaults to 3. :type order1: int, optional :param order2: Extent along the second basis vector. Defaults to 3. :type order2: int, optional :param order3: Extent along the third basis vector. Defaults to 3. :type order3: int, optional :param order4: Extent along the fourth basis vector. Defaults to 3. :type order4: int, optional :param shape: Detector shape used to reject off-detector nodes. Defaults to (124, 124). :type shape: tuple[int, int], optional :returns: ``((i, j, k, l), (x, y))`` pairs of lattice index and pixel position. :rtype: list .. py:function:: find_nearest_peak(vector, matrix, crop=4) Finds the brightest pixel within a window around a predicted position. :param vector: Predicted (x, y) position. :type vector: tuple[float, float] :param matrix: 2D image to search, typically the stack average. :type matrix: np.ndarray :param crop: Half-width of the search window. Defaults to 4. :type crop: int, optional :returns: Index of the brightest pixel. :rtype: tuple[int, int] .. py:function:: getGridVectors_Adjust(x0, y0, x1, y1, x2, y2, average_after, order=4, crop=2, shape=DEFAULT_DETECTOR_SHAPE) Builds a two-vector grid and snaps each node onto the nearest observed peak. :param x0: Direct-beam position along the first detector axis. :type x0: float :param y0: Direct-beam position along the second detector axis. :type y0: float :param x1: First basis-vector tip, first axis. :type x1: float :param y1: First basis-vector tip, second axis. :type y1: float :param x2: Second basis-vector tip, first axis. :type x2: float :param y2: Second basis-vector tip, second axis. :type y2: float :param average_after: Averaged 2D pattern used to locate peaks. :type average_after: np.ndarray :param order: Grid extent along each basis vector. Defaults to 4. :type order: int, optional :param crop: Half-width of the peak search window. Defaults to 2. :type crop: int, optional :param shape: Detector shape used to reject off-detector nodes. Defaults to (124, 124). :type shape: tuple[int, int], optional :returns: ``((i, j), (x, y))`` pairs, with (x, y) snapped to observed peaks. :rtype: list .. py:function:: AnnularSum(matrix, center_x, center_y, ri, ro) Averages counts within an annular aperture. :param matrix: 4D diffraction stack. :type matrix: np.ndarray :param center_x: Annulus centre along the first detector axis. :type center_x: int :param center_y: Annulus centre along the second detector axis. :type center_y: int :param ri: Inner radius in pixels. :type ri: int :param ro: Outer radius in pixels. :type ro: int :returns: Mean intensity per scan position, shaped (scan_x, scan_y). :rtype: np.ndarray .. py:function:: AnnularVar(matrix, center_x, center_y, g, m, DQE, ri, ro) Computes the variance of the intensity inside an annular aperture. :param matrix: 4D diffraction stack. :type matrix: np.ndarray :param center_x: Annulus centre along the first detector axis. :type center_x: int :param center_y: Annulus centre along the second detector axis. :type center_y: int :param g: Detector gain. Unused; kept for signature compatibility. :type g: float :param m: Mixing factor. Unused; kept for signature compatibility. :type m: float :param DQE: Detective quantum efficiency. Unused; kept for signature compatibility. :type DQE: float :param ri: Inner radius in pixels. :type ri: int :param ro: Outer radius in pixels. :type ro: int :returns: Variance across the annulus per scan position. :rtype: np.ndarray .. py:function:: Precession(matrix, center_x, center_y, ri, ro) Averages diffraction patterns over an annulus of scan positions. :param matrix: 4D diffraction stack. :type matrix: np.ndarray :param center_x: Annulus centre along the first scan axis. :type center_x: int :param center_y: Annulus centre along the second scan axis. :type center_y: int :param ri: Inner radius in scan steps. :type ri: int :param ro: Outer radius in scan steps. :type ro: int :returns: The averaged 2D diffraction pattern. :rtype: np.ndarray .. py:function:: Precession_Angled(matrix, center_x, center_y, ri, ro, angle1, angle2) Averages diffraction patterns over an angular sector of an annulus. :param matrix: 4D diffraction stack. :type matrix: np.ndarray :param center_x: Annulus centre along the first scan axis. :type center_x: int :param center_y: Annulus centre along the second scan axis. :type center_y: int :param ri: Inner radius in scan steps. :type ri: int :param ro: Outer radius in scan steps. :type ro: int :param angle1: Start of the angular range in degrees. :type angle1: float :param angle2: End of the angular range in degrees. :type angle2: float :returns: The averaged 2D diffraction pattern. :rtype: np.ndarray .. py:function:: radial_profile(data, center) Computes the azimuthally averaged radial profile of an image. :param data: 2D image. :type data: np.ndarray :param center: Centre of the profile. :type center: tuple[int, int] :returns: Mean intensity as a function of integer radius. :rtype: np.ndarray .. py:function:: get_gvectors(g1, g2, grid, g3=None, g4=None) Converts lattice indices on a pixel grid into Miller indices. :param g1: Miller index of the first basis vector. :type g1: tuple[int, int, int] :param g2: Miller index of the second basis vector. :type g2: tuple[int, int, int] :param grid: Grid produced by :func:`getGridVectors`, :func:`getGridVectors_Adjust`, or :func:`getGridVectorsHOLZ`. :type grid: list :param g3: Miller index of the third basis vector. Supply together with ``g4`` for a four-vector (HOLZ) grid. Defaults to None. :type g3: tuple[int, int, int], optional :param g4: Miller index of the fourth basis vector. Defaults to None. :type g4: tuple[int, int, int], optional :returns: One Miller index per grid node. :rtype: list[tuple[int, int, int]]