Reference¶
LiberTEM-blobfinder is structured into three parts:
A “base” package with numerics functions that work independent of LiberTEM.
A “common” package that uses other “common” aspects of LiberTEM for convenience, but can be used independent of LiberTEM core facilities.
A “udf” package with classes and functions to use this functionality with full LiberTEM integration.
Basic numerics functions¶
These functions work independent of any LiberTEM infrastructure.
-
libertem_blobfinder.base.correlation.
allocate_crop_bufs
(crop_size, n_peaks, dtype, limit=524288)[source]¶ Allocate buffer for stack of cropped peaks
The size is optimized to fit within
limit
. An aligned buffer for the FFT back-end is created if possible.- Parameters
crop_size (int) – The cropped parts will have size (2 * crop-size, 2 * crop_size)
n_peaks (int) – Number of peaks
dtype (numpy.dtype) – dtype of the buffer
limit (int, optional) – Upper limit, default 1/2 MB to be L3 cache friendly
- Returns
crop_bufs – Shape (n, 2*crop_size, 2*crop_size)
- Return type
np.ndarray
-
libertem_blobfinder.base.correlation.
do_correlations
(template, crop_parts)[source]¶ Calculate the correlation of the pre-calculated template with a stack of cropped peaks using fast correlation.
- Parameters
template (numpy.ndarray) – Real Fourier transform of the correlation pattern. The source pattern should have the same size as the cropped parts. Please note that the real Fourier transform (fft.rfft2) of the source pattern has a different shape!
crop_parts (numpy.ndarray) – Stack of peaks cropped from the frame.
- Returns
corrs – Correlation of the correlation pattern and the peaks.
- Return type
-
libertem_blobfinder.base.correlation.
get_buf_count
(crop_size, n_peaks, dtype, limit=524288)[source]¶ Calculate the optimal number of peaks in a stack to fit within the limit.
- Parameters
crop_size (int) – The cropped parts will have size (2 * crop-size, 2 * crop_size)
n_peaks (int) – Number of peaks
dtype (numpy.dtype) – dtype of the data for size calculation
limit (int, optional) – Upper limit, default 1/2 MB to be L3 cache friendly
- Returns
- Return type
-
libertem_blobfinder.base.correlation.
peak_elevation
(center, corrmap, height, r_min=1.5, r_max=inf)[source]¶ Return the slope of the tightest cone around
center
with heightheight
that touchescorrmap
betweenr_min
andr_max
.The correlation of two disks – mask and perfect diffraction spot – has the shape of a cone. The function’s return value correlates with the quality of a correlation. Higher slope means a strong peak and no side maxima, while weak signal or side maxima lead to a flatter slope.
- Parameters
center (numpy.ndarray) – (y, x) coordinates of the center within the
corrmap
corrmap (numpy.ndarray) – Correlation map
height (float) – The height is provided as a parameter since center can be float values from refinement and the height value is conveniently available from the calling function.
r_min (float, optional) – Masks out a small local plateau around the peak that would distort and dominate the calculation.
r_max (float, optional) – Mask out neighboring peaks if a large area with several legitimate peaks is correlated.
- Returns
elevation – Elevation of the tightest cone that fits the correlation map within the given parameter range.
- Return type
-
libertem_blobfinder.base.correlation.
process_frame_fast
(template, crop_size, frame, peaks, out_centers, out_refineds, out_heights, out_elevations, crop_bufs)[source]¶ Find the parameters of peaks in a diffraction pattern by correlation with a template
This function is designed to be used in an optimized pipeline with a pre-calculated Fourier transform of the match pattern and optional pre-allocated buffers. It is the engine of the
libertem_blobfinder.udf.correlation.FastCorrelationUDF
for stand-alone use independent of LiberTEM.libertem_blobfinder.common.correlation.process_frames_fast()
offers a more convenient interface for batch processing.- Parameters
template (numpy.ndarray) – Real Fourier transform of the correlation pattern. The source pattern should have size (2 * crop_size, 2 * crop_size). Please note that the real Fourier transform (fft.rfft2) of the source pattern has a different shape!
crop_size (int) – Half the size of the correlation pattern. Given as a parameter since real Fourier transform changes the size.
frame (np.ndarray) – Frame data. Currently, only Real values are supported.
peaks (np.ndarray) – List of peaks of shape (n_peaks, 2)
out_centers (np.ndarray) – Output buffer for center positions of shape (n_peaks, 2) and integer dtype.
out_refineds (np.ndarray) – Output buffer for refined center positions of shape (n_peaks, 2) and float dtype.
out_heights (np.ndarray) – Output buffer for peak height in log scaled frame. Shape (n_peaks, ) and float dtype.
out_elevations (np.ndarray) – Output buffer for peak elevation in log scaled frame. Shape (n_peaks, ) and float dtype.
crop_bufs (np.ndarray) – Aligned buffer for pyfftw. Shape (n, 2 * crop_size, 2 * crop_size) and float dtype. n doesn’t have to match the number of peaks. Instead, it should be chosen for good L3 cache efficiency.
allocate_crop_bufs()
can be used to allocate this buffer.
- Returns
The values are placed in the provided output buffers.
- Return type
Example
>>> from libertem_blobfinder.common.patterns import RadialGradient >>> from libertem_blobfinder.base.correlation import allocate_crop_bufs >>> >>> frames, indices, peaks = libertem.utils.generate.cbed_frame(radius=4) >>> pattern = RadialGradient(radius=4) >>> crop_size = pattern.get_crop_size() >>> template = pattern.get_template(sig_shape=(2 * crop_size, 2 * crop_size)) >>> >>> centers = np.zeros((len(frames), len(peaks), 2), dtype=np.uint16) >>> refineds = np.zeros((len(frames), len(peaks), 2), dtype=np.float32) >>> heights = np.zeros((len(frames), len(peaks)), dtype=np.float32) >>> elevations = np.zeros((len(frames), len(peaks)), dtype=np.float32) >>> >>> crop_bufs = allocate_crop_bufs(crop_size, len(peaks), frames.dtype) >>> >>> for i, f in enumerate(frames): ... process_frame_fast( ... template=template, crop_size=crop_size, ... frame=f, peaks=peaks.astype(np.int32), ... out_centers=centers[i], out_refineds=refineds[i], ... out_heights=heights[i], out_elevations=elevations[i], ... crop_bufs=crop_bufs ... ) >>> assert np.allclose(refineds[0], peaks, atol=0.1)
-
libertem_blobfinder.base.correlation.
process_frame_full
(template, crop_size, frame, peaks, out_centers=None, out_refineds=None, out_heights=None, out_elevations=None, frame_buf=None, buf_count=None)[source]¶ Find the parameters of peaks in a diffraction pattern by correlation with a template
This function is designed to be used in an optimized pipeline with a pre-calculated Fourier transform of the match pattern and optional pre-allocated buffers. It is the engine of the
libertem_blobfinder.udf.correlation.FullFrameCorrelationUDF
for stand-alone use independent of LiberTEM.libertem_blobfinder.common.correlation.process_frames_full()
offers a more convenient interface for batch processing.- Parameters
template (numpy.ndarray) – Real Fourier transform of the correlation pattern. The source pattern should have size (2 * crop_size, 2 * crop_size). Please note that the real Fourier transform (fft.rfft2) of the source pattern has a different shape!
crop_size (int) – Half the size of the correlation pattern. Given as a parameter since real Fourier transform changes the size.
frame (np.ndarray) – Frame data. Currently, only real values are supported.
peaks (np.ndarray) – List of peaks of shape (n_peaks, 2)
out_centers (np.ndarray, optional) – Output buffer for center positions of shape (n_peaks, 2) and integer dtype. Will be allocated if needed.
out_refineds (np.ndarray, optional) – Output buffer for refined center positions of shape (n_peaks, 2) and float dtype. Will be allocated if needed.
out_heights (np.ndarray, optional) – Output buffer for peak height in log scaled frame. Shape (n_peaks, ) and float dtype. Will be allocated if needed.
out_elevations (np.ndarray, optional) – Output buffer for peak elevation in log scaled frame. Shape (n_peaks, ) and float dtype. Will be allocated if needed.
frame_buf (np.ndarray) – Aligned buffer for FFT back-end, such as pyfftw. Shape of a frame and float dtype.
libertem_blobfinder.base.correlation.zero()
can be used.buf_count (int) – Number of peaks to process per outer loop iteration. This allows optimization of L3 cache efficiency.
- Returns
The values are placed in the provided output buffers.
- Return type
Example
>>> from libertem_blobfinder.common.patterns import RadialGradient >>> from libertem_blobfinder.base.correlation import get_buf_count, zeros >>> >>> frames, indices, peaks = libertem.utils.generate.cbed_frame() >>> pattern = RadialGradient(radius=4) >>> crop_size = pattern.get_crop_size() >>> template = pattern.get_template(sig_shape=frames[0].shape) >>> >>> centers = np.zeros((len(frames), len(peaks), 2), dtype=np.uint16) >>> refineds = np.zeros((len(frames), len(peaks), 2), dtype=np.float32) >>> heights = np.zeros((len(frames), len(peaks)), dtype=np.float32) >>> elevations = np.zeros((len(frames), len(peaks)), dtype=np.float32) >>> >>> frame_buf = zeros(frames[0].shape, dtype=np.float32) >>> buf_count = get_buf_count(crop_size, len(peaks), frame_buf.dtype) >>> >>> for i, f in enumerate(frames): ... process_frame_full( ... template=template, crop_size=crop_size, ... frame=f, peaks=peaks.astype(np.int32), ... out_centers=centers[i], out_refineds=refineds[i], ... out_heights=heights[i], out_elevations=elevations[i], ... frame_buf=frame_buf, buf_count=buf_count ... ) >>> assert np.allclose(refineds[0], peaks, atol=0.1)
Common classes and functions¶
These functions and classes depend on other LiberTEM “common” packages, but can be used without the LiberTEM core infrastructure.
-
class
libertem_blobfinder.common.patterns.
BackgroundSubtraction
(radius, search=None, radius_outer=None)[source]¶ Solid circular disk surrounded with a balancing negative area
This pattern rejects background and avoids false positives at positions between peaks
-
__init__
(radius, search=None, radius_outer=None)[source]¶ - Parameters
radius (float) – Radius of the circular pattern in px
search (float, optional) – Range from the center point in px to include in the correlation.
max(2*radius, radius_outer)
by default. Defining the size of the square correlation pattern.radius_outer (float, optional) – Radius of the negative region in px. 1.5x radius by default.
-
-
class
libertem_blobfinder.common.patterns.
Circular
(radius, search=None)[source]¶ Circular pattern with radius
radius
.This pattern is useful for constructing feature vectors using
feature_vector()
.New in version 0.3.0.
-
class
libertem_blobfinder.common.patterns.
MatchPattern
(search)[source]¶ Abstract base class for correlation patterns.
This class provides an API to provide a template for fast correlation-based peak finding.
-
class
libertem_blobfinder.common.patterns.
RadialGradient
(radius, search=None)[source]¶ Radial gradient from zero in the center to one at
radius
.This pattern rejects the influence of internal intensity variations of the CBED disk.
-
class
libertem_blobfinder.common.patterns.
RadialGradientBackgroundSubtraction
(radius, search=None, radius_outer=None, delta=1, radial_map=None)[source]¶ Combination of radial gradient with background subtraction
-
__init__
(radius, search=None, radius_outer=None, delta=1, radial_map=None)[source]¶ See
radial_gradient_background_subtraction()
for details.- Parameters
radius (float) – Radius of the circular pattern in px
search (float, optional) – Range from the center point in px to include in the correlation.
max(2*radius, radius_outer)
by default Defining the size of the square correlation pattern.radius_outer (float, optional) – Radius of the negative region in px. 1.5x radius by default.
delta (float, optional) – Width of the transition region between positive and negative in px
radial_map (numpy.ndarray, optional) – Radius value of each pixel in px. This can be used to distort the shape as needed or work in physical coordinates instead of pixels. A suitable map can be generated with
libertem.masks.polar_map()
.
Example
>>> import matplotlib.pyplot as plt
>>> (radius, phi) = libertem.masks.polar_map( ... centerX=64, centerY=64, ... imageSizeX=128, imageSizeY=128, ... stretchY=2., angle=np.pi/4 ... )
>>> template = RadialGradientBackgroundSubtraction( ... radius=30, radial_map=radius)
>>> # This shows an elliptical template that is stretched >>> # along the 45 ° bottom-left top-right diagonal >>> plt.imshow(template.get_mask(sig_shape=(128, 128))) <matplotlib.image.AxesImage object at ...> >>> plt.show()
-
-
class
libertem_blobfinder.common.patterns.
UserTemplate
(template, search=None)[source]¶ User-defined template
-
__init__
(template, search=None)[source]¶ - Parameters
template (numpy.ndarray) – Correlation template as 2D numpy.ndarray
search (float, optional) – Range from the center point in px to include in the correlation. Half diagonal of the template by default. Defining the size of the square correlation pattern.
-
-
libertem_blobfinder.common.patterns.
feature_vector
(imageSizeX, imageSizeY, peaks, match_pattern: libertem_blobfinder.common.patterns.MatchPattern)[source]¶ This function generates a sparse mask stack to extract a feature vector.
A match template based on the parameters in
parameters
is placed at each peak position in an individual mask layer. This mask stack can then be used inlibertem.udf.masks.ApplyMasksUDF
to generate a feature vector for each frame.Summing up the mask stack along the first axis generates a mask that can be used for virtual darkfield imaging of all peaks together.
- Parameters
imageSizeX (int) – Frame size in px
imageSizeY (int) – Frame size in px
peaks (numpy.ndarray) – Peak positions in px as numpy.ndarray of shape (n, 2) with integer type
match_pattern (MatchPattern) – Instance of
MatchPattern
-
libertem_blobfinder.common.correlation.
get_correlation
(sum_result, match_pattern: libertem_blobfinder.common.patterns.MatchPattern)[source]¶ Calculate the correlation between
sum_result
andmatch_pattern
.New in version 0.4.0.dev0.
- Parameters
sum_result (numpy.ndarray) – 2D result frame as correlation input
match_pattern (MatchPattern) – Instance of
MatchPattern
to correlatesum_result
with
-
libertem_blobfinder.common.correlation.
get_peaks
(sum_result, match_pattern: libertem_blobfinder.common.patterns.MatchPattern, num_peaks)[source]¶ Find peaks of the correlation between
sum_result
andmatch_pattern
.The result can then be used as input to
full_match()
to extract grid parameters,run_fastcorrelation()
to find the position in each frame or to construct a mask to extract feature vectors withfeature_vector()
.- Parameters
sum_result (numpy.ndarray) – 2D result frame as correlation input
match_pattern (MatchPattern) – Instance of
MatchPattern
to correlatesum_result
withnum_peaks (int) – Number of peaks to find
Example
>>> frame, _, _ = libertem.utils.generate.cbed_frame(radius=4) >>> pattern = libertem_blobfinder.common.patterns.RadialGradient(radius=4) >>> peaks = get_peaks(frame[0], pattern, 7) >>> print(peaks) [[64 64] [64 80] [80 80] [80 64] [48 80] [48 64] [64 96]]
-
libertem_blobfinder.common.correlation.
process_frames_fast
(pattern: libertem_blobfinder.common.patterns.MatchPattern, frames, peaks)[source]¶ Find the parameters of peaks in a diffraction pattern by correlation with a match pattern.
This method crops regions of interest around the peaks from the frames before correlation, which is usually fastest for a moderate amount of moderately sized peaks per frame.
Note
FastCorrelationUDF
is a parallelized, distributed version for large-scale data.- Parameters
pattern (MatchPattern) – Pattern to correlate with.
frames (np.ndarray) – Frame data. Currently, only Real values are supported.
peaks (np.ndarray) – List of peaks of shape (n_peaks, 2)
- Returns
centers (np.ndarray) – Center positions of shape (n_peaks, 2) and integer dtype.
refineds (np.ndarray) – Refined center positions of shape (n_peaks, 2) and float dtype.
heights (np.ndarray) – Peak height in log scaled frame. Shape (n_peaks, ) and float dtype.
elevations (np.ndarray) – Peak elevation in log scaled frame. Shape (n_peaks, ) and float dtype
Example
>>> frames, indices, peaks = libertem.utils.generate.cbed_frame() >>> pattern = libertem_blobfinder.common.patterns.RadialGradient(radius=4) >>> (centers, refineds, heights, elevations) = process_frames_fast( ... pattern=pattern, ... frames=frames, ... peaks=peaks.astype(np.int32), ... ) >>> assert np.allclose(refineds[0], peaks, atol=0.1)
-
libertem_blobfinder.common.correlation.
process_frames_full
(pattern: libertem_blobfinder.common.patterns.MatchPattern, frames, peaks)[source]¶ Find the parameters of peaks in a diffraction pattern by correlation with a match pattern.
This method crops regions of interest around the peaks after correlation, which can be faster for many peaks on smaller frames.
Note
FullFrameCorrelationUDF
is a parallelized, distributed version for large-scale data.- Parameters
pattern (MatchPattern) – Pattern to correlate with.
frame (np.ndarray) – Frame data. Currently, only real values are supported.
peaks (np.ndarray) – List of peaks of shape (n_peaks, 2)
- Returns
centers (np.ndarray) – Center positions of shape (n_peaks, 2) and integer dtype.
refineds (np.ndarray) – Refined center positions of shape (n_peaks, 2) and float dtype.
heights (np.ndarray) – Peak height in log scaled frame. Shape (n_peaks, ) and float dtype.
elevations (np.ndarray) – Peak elevation in log scaled frame. Shape (n_peaks, ) and float dtype
Example
>>> frames, indices, peaks = libertem.utils.generate.cbed_frame(radius=4) >>> pattern = libertem_blobfinder.common.patterns.RadialGradient(radius=4) >>> (centers, refineds, heights, elevations) = process_frames_full( ... pattern=pattern, ... frames=frames, ... peaks=peaks.astype(np.int32) ... ) >>> assert np.allclose(refineds[0], peaks, atol=0.1)
User-defined functions¶
These functions and classes depend on LiberTEM core infrastructure.
Correlation¶
UDFs and utility functions to find peaks and refine their positions by using correlation.
-
class
libertem_blobfinder.udf.correlation.
CorrelationUDF
(peaks, zero_shift=None, *args, **kwargs)[source]¶ Bases:
libertem.udf.base.UDF
Base class for peak correlation implementations
-
__init__
(peaks, zero_shift=None, *args, **kwargs)[source]¶ - Parameters
peaks (numpy.ndarray) – Numpy array of (y, x) coordinates with peak positions in px to correlate
zero_shift (Union[AUXBufferWrapper, numpy.ndarray, None], optional) – Zero shift, for example descan error. Can be
None
,numpy.array((y, x))
or AUX data with(y, x)
for each frame.
-
get_result_buffers
()[source]¶ The common buffers for all correlation methods.
centers
:(y, x) integer positions.
refineds
:(y, x) positions with subpixel refinement.
peak_values
:Peak height in the log scaled frame.
peak_elevations
:Peak quality (result of
peak_elevation()
).
See source code for details of the buffer declaration.
-
output_buffers
()[source]¶ This function allows abstraction of the result buffers from the default implementation in
get_result_buffers()
.Override this function if you wish to redirect the results to different buffers, for example ragged arrays or binned processing.
-
-
class
libertem_blobfinder.udf.correlation.
FastCorrelationUDF
(peaks, match_pattern, zero_shift=None, *args, **kwargs)[source]¶ Bases:
libertem_blobfinder.udf.correlation.CorrelationUDF
Fourier-based fast correlation-based refinement of peak positions within a search frame for each peak.
-
__init__
(peaks, match_pattern, zero_shift=None, *args, **kwargs)[source]¶ - Parameters
peaks (numpy.ndarray) – Numpy array of (y, x) coordinates with peak positions in px to correlate
match_pattern (MatchPattern) – Instance of
MatchPattern
zero_shift (Union[AUXBufferWrapper, numpy.ndarray, None], optional) – Zero shift, for example descan error. Can be
None
,numpy.array((y, x))
or AUX data with(y, x)
for each frame.
-
-
class
libertem_blobfinder.udf.correlation.
FullFrameCorrelationUDF
(peaks, match_pattern, zero_shift=None, *args, **kwargs)[source]¶ Bases:
libertem_blobfinder.udf.correlation.CorrelationUDF
Fourier-based correlation-based refinement of peak positions within a search frame for each peak using a single correlation step. This can be faster for correlating a large number of peaks in small frames in comparison to
FastCorrelationUDF
. However, it is more sensitive to interference from strong peaks next to the peak of interest.New in version 0.3.0.
-
__init__
(peaks, match_pattern, zero_shift=None, *args, **kwargs)[source]¶ - Parameters
peaks (numpy.ndarray) – Numpy array of (y, x) coordinates with peak positions in px to correlate
match_pattern (MatchPattern) – Instance of
MatchPattern
zero_shift (Union[AUXBufferWrapper, numpy.ndarray, None], optional) – Zero shift, for example descan error. Can be
None
,numpy.array((y, x))
or AUX data with(y, x)
for each frame.
-
-
class
libertem_blobfinder.udf.correlation.
SparseCorrelationUDF
(peaks, match_pattern, steps, *args, **kwargs)[source]¶ Bases:
libertem_blobfinder.udf.correlation.CorrelationUDF
Direct correlation using sparse matrices
This method allows to adjust the number of correlation steps independent of the template size.
-
__init__
(peaks, match_pattern, steps, *args, **kwargs)[source]¶ - Parameters
peaks (numpy.ndarray) – Numpy array of (y, x) coordinates with peak positions in px to correlate
match_pattern (MatchPattern) – Instance of
MatchPattern
steps (int) – The template is correlated with 2 * steps + 1 symmetrically around the peak position in x and y direction. This defines the maximum shift that can be detected. The number of calculations grows with the square of this value, that means keeping this as small as the data allows speeds up the calculation.
-
get_result_buffers
()[source]¶ This method adds the
corr
buffer to the result ofCorrelationUDF.get_result_buffers()
. See source code for the exact buffer declaration.
-
-
libertem_blobfinder.udf.correlation.
run_blobfinder
(ctx, dataset, match_pattern: libertem_blobfinder.common.patterns.MatchPattern, num_peaks, roi=None, progress=False)[source]¶ Wrapper function to find peaks in a dataset and refine their position using
FastCorrelationUDF
- Parameters
ctx (libertem.api.Context) –
dataset (libertem.io.dataset.base.DataSet) –
match_pattern (libertem_blobfinder.patterns.MatchPattern) –
num_peaks (int) – Number of peaks to look for
roi (numpy.ndarray, optional) – Boolean mask of the navigation dimension to select region of interest (ROI)
progress (bool, optional) – Show progress bar
- Returns
sum_result (numpy.ndarray) – Log-scaled sum frame of the dataset/ROI
centers, refineds, peak_values, peak_elevations (libertem.common.buffers.BufferWrapper) – See
CorrelationUDF.get_result_buffers()
for details.peaks (numpy.ndarray) – List of found peaks with (y, x) coordinates
-
libertem_blobfinder.udf.correlation.
run_fastcorrelation
(ctx, dataset, peaks, match_pattern: libertem_blobfinder.common.patterns.MatchPattern, zero_shift=None, roi=None, progress=False)[source]¶ Wrapper function to construct and run a
FastCorrelationUDF
- Parameters
ctx (libertem.api.Context) –
dataset (libertem.io.dataset.base.DataSet) –
peaks (numpy.ndarray) – List of peaks with (y, x) coordinates
match_pattern (libertem_blobfinder.patterns.MatchPattern) –
zero_shift (Union[AUXBufferWrapper, numpy.ndarray, None], optional) – Zero shift, for example descan error. Can be
None
,numpy.array((y, x))
or AUX data with(y, x)
for each frame.roi (numpy.ndarray, optional) – Boolean mask of the navigation dimension to select region of interest (ROI)
progress (bool, optional) – Show progress bar
- Returns
buffers – See
CorrelationUDF.get_result_buffers()
for details.- Return type
Refinement¶
UDFs and utility functions to refine grid parameters from peak positions.
-
class
libertem_blobfinder.udf.refinement.
AffineMixin
(*args, **kwargs)[source]¶ Bases:
libertem_blobfinder.udf.refinement.RefinementMixin
Refinement using
affinematch()
-
__init__
(*args, **kwargs)[source]¶ - Parameters
matcher (libertem.analysis.gridmatching.Matcher) – Instance of
Matcher
indices (numpy.ndarray) – List of indices [(h1, k1), (h2, k2), …] of all peaks. The indices can be non-integer and relative to any base vectors, including virtual ones like (1, 0); (0, 1). See documentation of
affinematch()
for details.
-
-
class
libertem_blobfinder.udf.refinement.
FastmatchMixin
(*args, **kwargs)[source]¶ Bases:
libertem_blobfinder.udf.refinement.RefinementMixin
Refinement using
fastmatch()
-
__init__
(*args, **kwargs)[source]¶ - Parameters
matcher (libertem.analysis.gridmatching.Matcher) – Instance of
Matcher
start_zero (numpy.ndarray) – Approximate value (y, x) in px for “zero” point (origin, zero order peak)
start_a (numpy.ndarray) – Approximate value (y, x) in px for “a” vector.
start_b (numpy.ndarray) – Approximate value (y, x) in px for “b” vector.
-
-
class
libertem_blobfinder.udf.refinement.
RefinementMixin
[source]¶ Bases:
object
To be combined with a
libertem_blobfinder.CorrelationUDF
using multiple inheritance.The mixin must come before the UDF in the inheritance list.
The subclasses implement a
postprocess
method that calculates a refinement of start_zero, start_a and start_b based on the correlation result and populates the appropriate result buffers with this refinement result.This allows combining arbitrary implementations of correlation-based matching with arbitrary implementations of the refinement by declaring an ad-hoc class that inherits from one subclass of RefinementMixin and one subclass of CorrelationUDF.
-
apply_match
(index, match)[source]¶ Override this method to change how a match is saved in the result buffers, for example to support binned processing or ragged result arrays.
-
get_result_buffers
()[source]¶ This adds
zero
,a
,b
,selector
,error
to the superclass result buffer declaration.zero
,a
,b
:Grid refinement parameters for each frame.
selector
:Boolean mask of the peaks that were used in the fit.
error
:Residual of the fit.
See source code for the exact buffer declaration.
-
-
libertem_blobfinder.udf.refinement.
run_refine
(ctx, dataset, zero, a, b, match_pattern: libertem_blobfinder.common.patterns.MatchPattern, matcher: libertem.analysis.gridmatching.Matcher, correlation='fast', match='fast', indices=None, steps=5, zero_shift=None, roi=None, progress=False)[source]¶ Wrapper function to refine the given lattice for each frame by calculating approximate peak positions and refining them for each frame using a combination of
libertem_blobfinder.CorrelationUDF
andlibertem_blobfinder.RefinementMixin
.Changed in version 0.3.0: Support for
FullFrameCorrelationUDF
through parametercorrelation = 'fullframe'
- Parameters
ctx (libertem.api.Context) – Instance of a LiberTEM
Context
dataset (libertem.io.dataset.base.DataSet) – Instance of a
DataSet
zero (numpy.ndarray) – Approximate value for “zero” point (y, x) in px (origin, zero order peak)
a (numpy.ndarray) – Approximate value for “a” vector (y, x) in px.
b (numpy.ndarray) – Approximate value for “b” vector (y, x) in px.
match_pattern (MatchPattern) – Instance of
MatchPattern
matcher (libertem.analysis.gridmatching.Matcher) – Instance of
Matcher
to perform the matchingcorrelation ({'fast', 'sparse', 'fullframe'}, optional) – ‘fast’, ‘sparse’ or ‘fullframe’ to select
FastCorrelationUDF
,SparseCorrelationUDF
orFullFrameCorrelationUDF
match ({'fast', 'affine'}, optional) – ‘fast’ or ‘affine’ to select
FastmatchMixin
orAffineMixin
indices (numpy.ndarray, optional) – Indices to refine. This is trimmed down to positions within the frame. As a convenience, for the indices parameter this function accepts both shape (n, 2) and (2, n, m) so that numpy.mgrid[h:k, i:j] works directly to specify indices. This saves boilerplate code when using this function. Default: numpy.mgrid[-10:10, -10:10].
steps (int, optional) – Only for correlation == ‘sparse’: Correlation steps. See
__init__()
for details.zero_shift (Union[AUXBufferWrapper, numpy.ndarray, None], optional) – Zero shift, for example descan error. Can be
None
,numpy.array((y, x))
or AUX data with(y, x)
for each frame. Only supported for correlation methodsfast
and fullframe.roi (numpy.ndarray, optional) – ROI for
run_udf()
progress (bool, optional) – Show progress bar. Default
False
- Returns
result (Dict[str, BufferWrapper]) – Result buffers of the UDF. See
libertem_blobfinder.correlation.CorrelationUDF.get_result_buffers()
andRefinementMixin.get_result_buffers()
for details on the available buffers.used_indices (numpy.ndarray) – The peak indices that were within the frame.
Examples
>>> dataset = ctx.load( ... filetype="memory", ... data=np.zeros(shape=(2, 2, 128, 128), dtype=np.float32) ... ) >>> (result, used_indices) = run_refine( ... ctx, dataset, ... zero=(64, 64), a=(1, 0), b=(0, 1), ... match_pattern=libertem_blobfinder.common.patterns.RadialGradient(radius=4), ... matcher=grm.Matcher() ... ) >>> result['centers'].data array(...)
Integration¶
UDFs to integrate peak intensity with positions specified per frame. If the peak
positions are sufficiently similar for all frames, you can use
libertem_blobfinder.common.patterns.feature_vector()
together with
libertem.udf.masks.ApplyMasksUDF
instead.
-
class
libertem_blobfinder.udf.integration.
IntegrationUDF
(centers, pattern)[source]¶ Bases:
libertem.udf.base.UDF
-
__init__
(centers, pattern)[source]¶ Integrate peak intensity at positions that are specified for each frame.
- Parameters
centers (AUXBufferWrapper) – Peak positions (y, x) as AUX buffer wrapper of kind “nav”, extra_shape (num_peaks, 2) and integer dtype.
pattern (libertem_blobfinder.common.patterns.MatchPattern) – Match pattern with the weight for each pixels.
libertem_blobfinder.common.patterns.BackgroundSubtraction
orlibertem_blobfinder.common.patterns.Circular
can be good choices.
Example
>>> from libertem_blobfinder.udf.integration import IntegrationUDF >>> from libertem_blobfinder.common.patterns import BackgroundSubtraction
>>> nav_shape = tuple(dataset.shape.nav) >>> sig_shape = tuple(dataset.shape.sig) >>> extra_shape = (3, 2) # three peaks with coordinates (y, x) >>> peaks_shape = nav_shape + extra_shape
>>> # Generate some random positions as an example >>> peaks = np.random.randint(low=0, high=np.min(sig_shape), size=peaks_shape)
>>> # Create an AuxBufferWrapper for the peaks >>> centers = IntegrationUDF.aux_data( ... data=peaks, ... kind='nav', ... dtype=np.int, ... extra_shape=extra_shape ... )
>>> udf = IntegrationUDF( ... centers=centers, ... pattern=BackgroundSubtraction(radius=5, radius_outer=6) ... )
>>> res = ctx.run_udf(udf=udf, dataset=dataset)
>>> nav_shape (16, 16) >>> # Integration result for each frame and peak >>> res['integration'].data.shape (16, 16, 3)
-