Skip to main content

Module crop

Module crop 

Source
Expand description

Scan geometry from the NIfTI affine: B0 direction, obliquity, and resampling to an axial grid.

The dipole kernel lives in the voxel grid, so an oblique acquisition must either supply the true B0 direction or be resampled to a cardinal-aligned grid. Wrapped phase has to be resampled in the complex domain — see geometry::resample_complex_to_axial. Cropping reconstruction to the region that carries signal, and putting the answer back.

FFT-based stages cost O(N log N) in the whole grid, not in the brain. See crop for why the box is also rounded up to FFT-friendly sizes, and for the wrap-around caveat. Cropping a volume to the region that actually carries signal, and putting it back.

Every FFT-based stage — background removal, dipole inversion, Laplacian unwrapping — costs O(N log N) in the whole grid, not in the brain. A brain typically fills a fifth of an acquired volume and a tenth of one that has been resampled to a cardinal grid (resampling wraps a box around a tilted slab, so half the result is empty corner). Reconstructing inside a box around the mask and padding the answer back afterwards is therefore close to free, and on a UK Biobank SWI it takes the resampled grid from 4.5× the cost of the acquired grid down to parity.

Two things make the difference, and they are worth separating because their risk differs:

  • Fewer voxels ([crop_box_for_mask]). Mask bounding box plus a margin, rather than the full field of view. This moves the FFT’s periodic boundary closer to the object.
  • Sizes an FFT likes ([fft_pad_box]). An extent with a large prime factor pushes rustfft onto Bluestein’s algorithm. An axially-resampled UK Biobank grid comes out 272×339×77 — that is 2⁴·17, 3·113 and 7·11, awkward on every axis, which is what resampling to a bounding box tends to produce. Padding to 280×343×80 costs 8% more voxels and takes the transform from 131 ms to 71 ms, a 1.85× speedup. Padding moves the boundary further from the object, so unlike cropping it carries no wrap-around risk.

Both are expressed as a [CropBox], which may sit inside the grid (cropping), extend beyond it (padding), or do both on different axes.

§This changes the numbers, not just the speed

This caveat applies to cropping, not to padding. FFT-based reconstruction is periodic, so moving the boundary closer to the object brings wrap-around with it, and the dipole kernel has infinite support. Measured on a UK Biobank acquisition, a crop that actually removed voxels changed χ by ~0.6% of its dynamic range at the median and ~4% at the 99th percentile. [margin_voxels] takes the margin in millimetres so anisotropic voxels get a geometrically equal margin on every side, but a caller should still validate a cropped reconstruction against an uncropped one rather than assume the two agree.

[fft_pad_box] has no such caveat: it discards nothing and only moves the boundary outward.

Structs§

CropBox
A box within a larger grid: where reconstruction actually happens.

Functions§

crop_box_for_mask
The box to reconstruct in: the mask’s bounding box, grown by margin_mm on every side, each axis then rounded up to an FFT-friendly size and clamped to the grid.
crop_volume
Copy the box out of a full-grid volume.
crop_volume_with
Copy the box out of a full-grid volume, filling anything outside the grid with fill.
fft_pad_box
A box covering the whole grid, each axis grown outward to an FFT-friendly size.
margin_voxels
Margin in voxels per axis for a margin given in millimetres, at least one voxel where the margin is positive. Anisotropic voxels get a geometrically equal margin rather than an equal voxel count — 8 voxels is 6.4 mm in-plane but 24 mm through-plane at 0.8 × 0.8 × 3 mm.
next_fft_friendly_size
Smallest size >= n whose prime factors are all at most 7 — the radices rustfft has dedicated butterflies for. Sizes with a large prime factor fall back to Bluestein’s algorithm and cost several times more despite holding the same data.
uncrop_volume
Put a cropped volume back into a full-grid volume, filling everything outside the box with fill (zero for field maps and χ, which are undefined outside the mask anyway).