recovar.reconstruction¶
3-D reconstruction algorithms, noise estimation, and regularization.
homogeneous¶
Homogeneous (mean) 3-D reconstruction via direct Fourier inversion.
recovar.reconstruction.homogeneous
¶
Homogeneous (mean) 3-D reconstruction via direct Fourier inversion.
MeanEstimate(combined, corrected0, corrected1, corrected0reg, corrected1reg, lhs, prior, cubic_coeffs=None)
dataclass
¶
Result of half-set mean reconstruction.
Attributes¶
combined : ndarray — average of the two half-set reconstructions. corrected0, corrected1 : ndarray — unregularized half-set means. corrected0reg, corrected1reg : ndarray — regularized half-set means. lhs : ndarray — average Hermitian (CTF²) filter. prior : ndarray — FSC-derived regularization prior.
get_mean_conformation_relion(dataset, batch_size, noise_variance=None, use_regularization=False, upsampling_factor=2)
¶
Compute the mean conformation using RELION-style half-set reconstruction.
Runs one accumulation pass per half-set, computes the FSC-based prior, then optionally re-post-processes with regularization.
Parameters¶
dataset : CryoEMDataset (with halfset_indices set)
batch_size : int
noise_variance : array, optional
use_regularization : bool
If True, means["combined"] is the average of the per-half
regularized reconstructions; otherwise it is the average of the
unregularized reconstructions.
upsampling_factor : int
Volume oversampling factor (default 2).
Returns¶
means : MeanEstimate mean_prior : ndarray — FSC-derived regularization prior fsc : ndarray
noise¶
Noise variance estimation from cryo-EM image residuals.
recovar.reconstruction.noise
¶
Noise variance estimation from cryo-EM image residuals.
RadialNoiseModel(noise_variance_radial, image_shape=None)
¶
Per-dataset radial noise profile with lazy full/half expansions.
The previous implementation eagerly materialized both full-spectrum and half-spectrum arrays during construction. Most call sites consume only one layout, so we now build each cached representation on first use to reduce peak GPU memory without changing outputs.
get_half(*args, **kwargs)
¶
Return noise in half-spectrum (rfft-packed) layout (1, H*(W//2+1)).
VariableRadialNoiseModel(noise_variance_radials, dose_indices, image_shape=None)
¶
Radial noise profile indexed by dose/tilt level with lazy expansions.
get_half(indices, *args, **kwargs)
¶
Return noise in half-spectrum (rfft-packed) layout (B, H*(W//2+1)).
ConstantNoiseModel(noise_array)
¶
Wraps a pre-expanded noise array that is the same for every image.
Satisfies the same get / get_half interface as the radial models so
callers can treat all noise sources uniformly.
as_noise_model(cov_noise, image_shape)
¶
to_batched_pixel_noise(noise_variances, image_shape, batch_size=None)
¶
Normalize noise variance into shape broadcastable with (B, D*D).
Returns (1, DD) for uniform noise (callers rely on XLA broadcasting), or (B, DD) when the input is already per-image. Accepts common forms used in the codebase: - scalar (white noise) - (D, D) - (1, D, D) or (B, D, D) - (DD,) - (1, DD) or (B, D*D)
to_batched_half_pixel_noise(noise_variances, image_shape, batch_size=None)
¶
Normalize noise variance into half-spectrum layout broadcastable with (B, H*(W//2+1)).
Returns (1, N_half) for uniform noise or (B, N_half) for per-image noise. Callers should rely on XLA broadcasting rather than expecting exact (B, N_half).
predict_noise_variance(noise_variance, CTF_params, voxel_size, ctf, image_masks, image_shape, radial=True, premultiplied_ctf=False, upsample_factor=1)
¶
Predict noise variance in images, optionally handling upsampling.
Parameters:
| Name | Type | Description | Default |
|---|---|---|---|
noise_variance
|
Base noise variance (radial or scalar) |
required | |
CTF_params
|
CTF parameters |
required | |
voxel_size
|
Voxel size |
required | |
ctf
|
Function to compute CTF |
required | |
image_masks
|
Image masks |
required | |
image_shape
|
Image shape |
required | |
radial
|
Whether noise is radial |
True
|
|
premultiplied_ctf
|
Whether CTF is premultiplied |
False
|
|
upsample_factor
|
Factor to upsample by (default 1 for no upsampling) |
1
|
fit_noise_model_to_images(experiment_dataset, volume_mask, mean_estimate, image_subset, batch_size, invert_mask, disc_type='linear_interp', use_batch_solver=True, tilt_dose_inner=False, image_n_iter=10000.0)
¶
Fit noise model to images, handling tilt series data specially.
Parameters:
| Name | Type | Description | Default |
|---|---|---|---|
experiment_dataset
|
Dataset containing images and metadata |
required | |
volume_mask
|
Mask for the volume |
required | |
mean_estimate
|
Estimate of mean volume |
required | |
image_subset
|
Subset of images to use, or None for all |
required | |
batch_size
|
Batch size for processing |
required | |
invert_mask
|
Whether to invert the mask |
required | |
disc_type
|
Type of discretization to use |
'linear_interp'
|
|
use_batch_solver
|
Whether to use batch solver vs full dataset |
True
|
|
tilt_dose_inner
|
Whether this is an inner call for tilt series |
False
|
Returns:
| Name | Type | Description |
|---|---|---|
|
For tilt series: Array of noise variances per tilt |
||
Otherwise |
Single noise variance and initial estimate |
estimate_noise_variance(experiment_dataset, batch_size, max_images=10000)
¶
Estimate per-image noise variance from corner pixels.
Computes the noise power spectrum from image regions outside the particle mask, subsampling to at most max_images for efficiency.
Parameters:
| Name | Type | Description | Default |
|---|---|---|---|
experiment_dataset
|
A |
required | |
batch_size
|
Number of images to process per GPU batch. |
required | |
max_images
|
Maximum number of images to use for estimation. |
10000
|
Returns:
| Type | Description |
|---|---|
|
Tuple |
|
|
is a scalar noise variance and radial_noise_profile is the |
|
|
averaged radial power spectrum of the noise. |
make_radial_noise_half(average_image_PS, image_shape)
¶
Expand radial noise profile directly to half-image (rfft-packed) layout.
Output shape: (H * (W//2+1),) instead of (H * W,).
batch_basis_times_coords2(basis, coords)
¶
Compute basis @ coords.T reshaped for batched images, memory-efficiently.
average_residual_square(config, images, model, image_mask, contrasts, basis_coordinates, rotation_matrices, translations, ctf_params)
¶
Compute average residual squared — Equinox API.
Replaces the 19-param get_average_residual_square_inner_v2.
regularization¶
Fourier-shell regularization priors and FSC computation.
recovar.reconstruction.regularization
¶
Fourier-shell regularization priors and FSC computation.
compute_relion_prior(halfset_datasets, cov_noise, image0, image1, batch_size, estimate_merged_SNR=False, noise_level=None)
¶
Compute a RELION-style spectral prior from two half-set reconstructions.
Parameters:
| Name | Type | Description | Default |
|---|---|---|---|
halfset_datasets
|
Pair of half-set datasets. |
required | |
cov_noise
|
Scalar noise variance. |
required | |
image0
|
First half-map (Fourier coefficients). |
required | |
image1
|
Second half-map (Fourier coefficients). |
required | |
batch_size
|
GPU batch size for noise estimation. |
required | |
estimate_merged_SNR
|
Estimate SNR from merged map. |
False
|
|
noise_level
|
Pre-computed noise level (skips estimation if given). |
None
|
Returns:
| Type | Description |
|---|---|
|
Tuple |
|
|
curve, and averaged prior. |
get_fsc(vol1, vol2, volume_shape, substract_shell_mean=False, frequency_shift=0)
¶
Compute the Fourier Shell Correlation between two volumes.
Parameters:
| Name | Type | Description | Default |
|---|---|---|---|
vol1
|
First volume (flattened Fourier coefficients). |
required | |
vol2
|
Second volume (flattened Fourier coefficients). |
required | |
volume_shape
|
Tuple |
required | |
substract_shell_mean
|
Subtract per-shell mean before correlating. |
False
|
|
frequency_shift
|
Shift applied to frequency indices. |
0
|
Returns:
| Type | Description |
|---|---|
|
1-D array of FSC values, one per radial shell. |
relion_functions¶
RELION-compatible reconstruction and Wiener filtering routines.
recovar.reconstruction.relion_functions
¶
RELION-compatible reconstruction and Wiener filtering routines.
griddingCorrect(vol_in, ori_size, padding_factor, order=0)
¶
Radial sinc gridding correction.
griddingCorrect_square(vol_in, ori_size, padding_factor, order=0)
¶
Per-axis sinc product gridding correction (Fourier transform of trilinear interpolator).
relion_style_triangular_kernel(experiment_dataset, cov_noise, batch_size, disc_type='linear_interp', index_subset=None, upsampling_factor=None, by_image=True)
¶
RELION-style triangular kernel reconstruction.
Accumulates weighted back-projection (Ft_y) and CTF weight sum (Ft_ctf) in half-volume layout, then expands to full volume before returning.
Parameters¶
experiment_dataset : CryoEMDataset
cov_noise : array or None
Radial shell variances or pre-expanded noise. When None, noise is
drawn per-batch from experiment_dataset.noise.
batch_size : int
disc_type : str
index_subset : array-like, optional
If given, only iterate over this subset in the domain selected by
by_image.
upsampling_factor : int, optional
Defaults to 1.
by_image : bool, default True
When True, interpret index_subset as image indices.
When False, interpret it as group / particle indices and iterate
with the grouped dataset iterator. Needed for tilt-series particle
subsets such as junk-detection halfmap reconstruction.
relion_kernel_batch(config, images, ctf_params, rotation_matrices, translations, noise_variance, Ft_y=None, Ft_ctf=None)
¶
RELION-style triangular kernel batch for raw real-space images.
Applies pad + rfft2 internally, then backprojects with half_image and half_volume layouts for maximum memory efficiency.
relion_kernel_batch_from_fft(config, images, ctf_params, rotation_matrices, translations, noise_variance, Ft_y=None, Ft_ctf=None)
¶
RELION-style triangular kernel batch for pre-FFTed complex images.
Extracts the half-spectrum from full-spectrum images, then backprojects using half_image and half_volume layouts.
residual_relion_kernel_trilinear(config, mean_estimate, images, ctf_params, rotation_matrices, translations, cov_noise, Ft_y=None, Ft_ctf=None)
¶
Residual RELION-style kernel (trilinear), accumulating into half-volume layout.
Parameters¶
Ft_y, Ft_ctf : optional accumulator volumes (half-volume layout) to add into.
residual_relion_style_triangular_kernel(experiment_dataset, mean_estimate, cov_noise, batch_size, index_subset=None)
¶
Residual RELION-style triangular kernel reconstruction.
adjust_regularization_relion_style(filter, volume_shape, tau=None, padding_factor=1, max_res_shell=None, half_volume=False)
¶
Adjust the RELION-style regularization filter.
Adds 1/tau to the filter (Wiener denominator) and floors small values at 1/1000 of the spherically-averaged filter to avoid division by zero. See RELION backprojector.cpp for the original algorithm.
post_process_from_filter(cryo, Ft_ctf, F_ty, tau=None, disc_type='nearest', use_spherical_mask=True, grid_correct=True, gridding_correct='square', kernel_width=1)
¶
Post-process RELION-style reconstruction from filter weights.
Thin wrapper around post_process_from_filter_v2 that extracts the
necessary geometry from a dataset object (cryo).
post_process_from_filter_v2(Ft_ctf, F_ty, og_volume_shape, volume_upsampling_factor, tau=None, kernel='triangular', use_spherical_mask=True, grid_correct=True, gridding_correct='square', kernel_width=1, volume_mask=None, return_real_space=False, return_half_volume=False, input_half_volume=None)
¶
Post-process RELION-style reconstruction from filter weights.
Steps: regularize -> iDFT -> crop -> spherical mask -> grid correct -> DFT.
Supports both full Fourier inputs (N0*N1*N2,) and packed half-volume
inputs (N0*N1*(N2//2+1),).
relion_reconstruct(cryo, noise_variance, batch_size=100, disc_type='linear_interp', use_spherical_mask=True, upsampling_factor=2, grid_correct=True, gridding_correct='square', tau=None)
¶
Full mean reconstruction pipeline: accumulate → post-process.
ewald¶
Ewald sphere curvature correction for high-resolution reconstruction.
recovar.reconstruction.ewald
¶
Ewald Sphere Code: models, multipliers, etc. This code represents generally the most updated models for forward and backwards PART 1: Coordinate/Ewald Sphere Functions PART 2: CTF Explicit Functions PART 3: FORWARD/BACKWARD IMPLICIT MODELS PART 4: ADAPTATIONS FOR CG INPUT (i.e/mask/unmask)
logger = logging.getLogger(__name__)
module-attribute
¶
PART 1: Coordinate/Ewald Sphere Functions
get_good_idx_mask(volume_shape)
¶
Build real/imag frequency masks. Real part includes DC; imaginary excludes it.
volt_to_wavelength(volt)
¶
Convert accelerating voltage (kV) to electron wavelength (Angstroms).