QSM.rs

A Rust library for Quantitative Susceptibility Mapping (QSM) of the brain.

QSM.rs provides a complete set of algorithms for reconstructing magnetic susceptibility maps from MRI phase data, including brain extraction, phase unwrapping, background field removal, and dipole inversion.

Website · API Documentation

QSM.rs is the shared reconstruction engine behind the wider QSM ecosystem: the QSMxT command-line pipeline and the QSMbly in-browser app both call directly into this library, and QSM-CI benchmarks these methods against phantom ground truth. See the ecosystem hub for an overview.

Installation

Add qsm-core to your Cargo.toml:

[dependencies]
qsm-core = { git = "https://github.com/astewartau/QSM.rs" }

Feature flags

Feature	Enables
parallel	Rayon multi-threading (and multithreaded ONNX kernels with onnx)
simd	SIMD-accelerated inner loops via wide
onnx	Deep-learning methods, run through the pure-Rust tract engine (native and WASM)
download	Native download and SHA-256-checked caching of model weights (see qsm_core::models)
introspection	Serialize on the *Params structs, for tooling that inspects parameter names

[dependencies]
qsm-core = { git = "https://github.com/astewartau/QSM.rs", features = ["parallel", "onnx"] }

Model weights are not bundled: deep-learning functions take the ONNX bytes from the caller.

Usage

There are two ways to use the crate. Both are demonstrated as runnable examples (examples/pipeline_highlevel.rs, examples/pipeline_lowlevel.rs):

cargo run --release --example pipeline_highlevel
cargo run --release --example pipeline_lowlevel

High-level pipeline (recommended)

Describe the scan once with a ScanMetadata, then run the stages. The run_* functions dispatch to the configured algorithm and handle unit conversions internally.

use qsm_core::pipeline::{
    ScanMetadata, FieldMappingConfig, BgRemovalConfig, InversionConfig, QsmReference,
    run_field_mapping, run_bg_removal, run_dipole_inversion, apply_reference,
};

# fn run(phase: Vec<f64>, magnitude: Vec<f64>, mask: Vec<u8>) -> Result<(), qsm_core::pipeline::PipelineError> {
let meta = ScanMetadata {
    dims: (128, 128, 64),
    voxel_size: (1.0, 1.0, 1.0),
    echo_times: vec![0.020],        // seconds
    field_strength: 3.0,            // Tesla
    b0_direction: (0.0, 0.0, 1.0),
};

let phases: Vec<&[f64]> = vec![&phase];
let mags: Vec<&[f64]> = vec![&magnitude];

let field = run_field_mapping(&phases, Some(&mags), &mask, &meta,
    &FieldMappingConfig::default(), &mut |_, _| {})?;
let bg = run_bg_removal(&field.b0_field_ppm, &mask, &meta,
    &BgRemovalConfig::default(), Some(&magnitude), &mut |_, _| {})?;
let chi = run_dipole_inversion(&bg.local_field_ppm, &bg.eroded_mask, &meta,
    &InversionConfig::default(), Some(&magnitude), &mut |_, _| {})?;
let chi = apply_reference(&chi, &bg.eroded_mask, QsmReference::Mean);
# let _ = chi; Ok(())
# }

Low-level building blocks

Call the individual algorithm functions directly when you want to wire stages together yourself. Each takes a Grid, a *Params struct (all implement Default), and — for iterative methods — a progress callback.

use qsm_core::{Grid, bet, unwrap, bgremove, inversion};
use qsm_core::bet::BetParams;
use qsm_core::bgremove::VsharpParams;
use qsm_core::inversion::TvParams;

# fn run(phase: &[f64], magnitude: &[f64]) {
let grid = Grid::new(128, 128, 64, 1.0, 1.0, 1.0);
let bdir = (0.0, 0.0, 1.0); // B0 direction

let mask = bet::run_bet(magnitude, &grid, &BetParams::default(), |_, _| {});
let unwrapped = unwrap::laplacian_unwrap(phase, &mask, &grid, unwrap::LaplacianSolver::Dct);
let (local, eroded) = bgremove::vsharp(&unwrapped, &mask, &grid, &VsharpParams::default(), |_, _| {});
let chi = inversion::tv_admm(&local, &eroded, &grid, bdir, &TvParams::default(), |_, _| {});
# let _ = chi;
# }

Load and save NIfTI volumes with qsm_core::io.

Algorithms

Methods marked † need the onnx feature and model weights. Parity: agreement with the original implementation, checked by QSM-CI.

Brain Extraction and Segmentation

Algorithm	Description	Reference	Parity
BET	Brain extraction by surface-mesh evolution (run_bet)	Smith, S.M. (2002). "Fast robust automated brain extraction." Human Brain Mapping, 17(3):143-155. DOI
HD-BET †	Deep-learning brain extraction from a magnitude image (hd_bet)	Isensee, F., et al. (2019). "Automated brain extraction of multisequence MRI using artificial neural networks." Human Brain Mapping, 40(17):4952-4964. DOI
RS2-Net †	Deep-learning rodent brain extraction (rs2_net)	Lin, Y., et al. (2024). "RS2-Net: An end-to-end deep learning framework for rodent skull stripping in multi-center brain MRI." NeuroImage, 298:120769. DOI
MindGrab †	Deep-learning brain extraction for any contrast (mindgrab)	Fani, A., et al. (2025). "MindGrab for BrainChop: Fast and Accurate Skull Stripping for Command Line and Browser." arXiv:2506.11860. DOI
SynthSeg †	Contrast-agnostic whole-brain segmentation into FreeSurfer labels (segment::synthseg)	Billot, B., et al. (2023). "SynthSeg: Segmentation of brain MRI scans of any contrast and resolution without retraining." , 86:102789.

Mask Operations

Algorithm	Description	Reference	Parity
Erosion / dilation	Iterated 6-connected erosion and dilation (erode_mask, dilate_mask)	—
Closing	Dilation then erosion with a spherical structuring element (morphological_close)	—
Hole filling	Fills enclosed background components up to a size limit (fill_holes)	—
Gaussian smoothing	Smooths the mask and re-thresholds it at 0.5 (MaskOp::GaussianSmooth)	—
Signal-gated erosion	Mask refinement that erodes only low-signal boundary voxels (signal_gated_erosion)	QSM-CI harmonization masking, QSMxT/QSM-CI
Connected components	6-connected component labelling, as scipy.ndimage.label (internal: signal-gated erosion, PRELUDE, SynthSeg)	—

Coil Combination and Field Mapping

Algorithm	Description	Reference	Parity
Real/imaginary input	Converts real/imaginary echoes to magnitude/phase, detecting and undoing alternate-slice π chopping with one decision for all echoes, and optionally negating the phase for GE's opposite sign convention (complex_to_mag_phase, ComplexToMagPhaseParams)	—
MCPC-3D-S	Coil combination and phase-offset removal (mcpc3ds_combine, phase_offset_removal)	Eckstein, K., et al. (2018). "Computationally Efficient Combination of Multi-channel Phase Data From Multi-echo Acquisitions (ASPIRE)." Magnetic Resonance in Medicine, 79:2996-3006. DOI
Bipolar correction	Removes the linear phase artefact of bipolar readout gradients; needs three or more echoes (bipolar_correction)	Eckstein, K. (2021). PhD thesis, TU Wien, Section 3.1.3. DOI
B0 weighted average	Weighted average of unwrapped phase / TE over echoes; weights by phase SNR, phase variance, magnitude, TE, uniform, or an assumed T2* decay (calculate_b0_weighted, B0WeightType)	—
Linear fit	Magnitude-weighted least-squares fit of phase against TE, with optional phase offset and reliability mask (multi_echo_linear_fit)	—
QSMART echo fit

Phase Unwrapping

Algorithm	Description	Reference	Parity
ROMEO	Quality-guided region growing (unwrap_romeo)	Dymerska, B., et al. (2021). "Phase unwrapping with a rapid opensource minimum spanning tree algorithm (ROMEO)." Magnetic Resonance in Medicine, 85(4):2294-2308. DOI
Laplacian	Poisson-equation unwrapping of the total field, DCT, zero-padded FFT or mirrored FFT (lapunwrap.m, as in QSMART) solver (laplacian_unwrap)	Ghiglia, D.C., Romero, L.A. (1994). JOSA A, 11(1):107-117. DOI; Schofield, M.A., Zhu, Y. (2003). "Fast phase unwrapping algorithm for interferometric applications." Optics Letters, 28(14):1194-1196. DOI
Best path	3D-SRNCP reliability-ordered edge merging (unwrap_bestpath)	Abdul-Rahman, H., et al. (2007). "Fast and robust three-dimensional best path phase unwrapping algorithm." Applied Optics, 46(26):6623-6635. DOI
PRELUDE	Region merging over phase-range bins, with regions labelled slice by slice (FSL prelude's default) or over the whole volume (prelude -f) (unwrap_prelude)	Jenkinson, M. (2003). "Fast, automated, N-dimensional phase-unwrapping algorithm." Magnetic Resonance in Medicine, 49(1):193-197. DOI

Background Field Removal

Algorithm	Description	Reference	Parity
V-SHARP	Variable-radius SHARP	Wu, B., et al. (2012). "Whole brain susceptibility mapping using compressed sensing." Magnetic Resonance in Medicine, 67(1):137-147. DOI
SHARP	Spherical-mean-value deconvolution with TSVD	Schweser, F., et al. (2011). "Quantitative imaging of intrinsic magnetic tissue properties using MRI signal phase." NeuroImage, 54(4):2789-2807. DOI
RESHARP	SHARP with Tikhonov regularization	Sun, H., Wilman, A.H. (2014). "Background field removal using spherical mean value filtering and Tikhonov regularization." Magnetic Resonance in Medicine, 71(3):1151-1157. DOI
SMV	Spherical mean value filtering	Schweser, F., et al. (2011). NeuroImage, 54(4):2789-2807. DOI
iSMV	Iterative spherical mean value	Wen, Y., et al. (2014). "An iterative spherical mean value method for background field removal in MRI." , 72(4):1065-1071.

Combined Phase Unwrapping and Background Removal

Algorithm	Description	Reference	Parity
HARPERELLA	Integrated Laplacian unwrapping and background removal, wrapped phase or total field (ppm) → local field (harperella, HarperellaInput)	Li, W., et al. (2014). "Integrated Laplacian-based phase unwrapping and background phase removal for quantitative susceptibility mapping." NMR in Biomedicine, 27(2):219-227. DOI
iHARPERELLA	Improved HARPERELLA (iharperella)	Li, W., Wu, B., Liu, C. (2015). "iHARPERELLA: an improved method for integrated 3D phase unwrapping and background phase removal." Proc. ISMRM 23, p.3313.
iQFM †	Deep-learning wrapped phase → local field (run_iqfm)	Gao, Y., et al. (2022). "Instant tissue field and magnetic susceptibility mapping from MRI raw phase using Laplacian enhanced deep neural networks." NeuroImage, 259:119410. DOI

Dipole Inversion

Algorithm	Description	Reference	Parity
TKD	Truncated k-space division	Shmueli, K., et al. (2009). "Magnetic susceptibility mapping of brain tissue in vivo using MRI phase data." Magnetic Resonance in Medicine, 62(6):1510-1522. DOI
TSVD	Truncated SVD of the dipole kernel	Shmueli, K., et al. (2009). Magnetic Resonance in Medicine, 62(6):1510-1522. DOI
Tikhonov	Closed-form L2 regularization (identity, gradient or Laplacian)	Bilgic, B., et al. (2014). "Fast image reconstruction with L2-regularization." Journal of Magnetic Resonance Imaging, 40(1):181-191. DOI
TV	Total variation via ADMM (tv_admm)	Bilgic, B., et al. (2014). "Fast quantitative susceptibility mapping with L1-regularization and automatic parameter selection." Magnetic Resonance in Medicine, 72(5):1444-1459. DOI
NLTV	Nonlinear-data-term TV via ADMM (QSM.jl nltv)	Milovic, C., Bilgic, B., Zhao, B., Acosta-Cabronero, J., Tejos, C. (2018). "Fast nonlinear susceptibility inversion with variational regularization." Magnetic Resonance in Medicine, 80(2):814-821. DOI

Single-Step QSM

Algorithm	Input	Description	Reference	Parity
TGV	Total field (ppm) or one wrapped echo, with TE and B0	Single-step total generalized variation (tgv_qsm, TgvInput)	Langkammer, C., et al. (2015). "Fast quantitative susceptibility mapping using 3D EPI and total generalized variation." NeuroImage, 111:622-630. DOI
TFI	Total field (ppm)	Preconditioned total field inversion	Liu, Z., et al. (2017). "Preconditioned total field inversion (TFI) method for quantitative susceptibility mapping." Magnetic Resonance in Medicine, 78(1):303-315. DOI
QSMART	Total field (ppm), optional reliability map R_0	SDF background removal, inversion, and vasculature separation	Yaghmaie, N., Syeda, W., et al. (2021). NeuroImage, 231:117701. DOI
AutoQSM †	Total field (ppm)	Patch V-Net from total field to χ	Wei, H., et al. (2019). "Learning-based single-step quantitative susceptibility mapping reconstruction without brain extraction." NeuroImage, 202:116064. DOI
NeXtQSM †	Total field (ppm)	Background-removal U-Net followed by a variational inversion network

Multi-Orientation

Algorithm	Description	Reference	Parity
COSMOS	Closed-form multi-orientation inversion (cosmos, cosmos_weighted)	Liu, T., et al. (2009). "Calculation of susceptibility through multiple orientation sampling (COSMOS)." Magnetic Resonance in Medicine, 61(1):196-204. DOI
STI	Susceptibility tensor imaging from six or more orientations (sti)	Liu, C. (2010). "Susceptibility tensor imaging." Magnetic Resonance in Medicine, 63(6):1471-1477. DOI

Relaxometry

Algorithm	Description	Reference	Parity
R2*/T2* (ARLO)	R2* mapping from multi-echo magnitude (arlo, a port of the MEDI toolbox's arlo.m)	Pei, M., et al. (2015). "Algorithm for fast monoexponential fitting based on Auto-Regression on Linear Operations (ARLO) of data." Magnetic Resonance in Medicine, 73(2):843-850. DOI
R2 (EPG)	R2 mapping from multi-echo spin echo by extended-phase-graph matching (r2_epg)	Weigel, M. (2015). "Extended phase graphs: dephasing, RF pulses, and echoes — pure and simple." Journal of Magnetic Resonance Imaging, 41(2):266-295. DOI
R2′	R2* − R2 from paired gradient- and spin-echo data (r2prime)	Yablonskiy, D.A., Haacke, E.M. (1994). "Theory of NMR signal behavior in magnetically inhomogeneous tissues: the static dephasing regime." Magnetic Resonance in Medicine, 32(6):749-763. DOI
R2PRIMEnet †	Deep-learning R2* → R2′ for GRE-only data (3 T)	Kim, M., et al. (2025). "χ-sepnet: Deep neural network for magnetic susceptibility source separation." Human Brain Mapping, 46(4):e70136. DOI

Susceptibility Source Separation

Algorithm	Description	Reference	Parity
χ-separation	Paramagnetic/diamagnetic separation from local field, R2′ and magnitude: the SNU χ-separation toolbox's MEDI variant (chi_sep_medi, which also takes N_std, a QSM and a CSF mask) and an iLSQR-style variant (chi_sep_ilsqr)	Shin, H., et al. (2021). "χ-separation: Magnetic susceptibility source separation toward iron and myelin mapping in the brain." NeuroImage, 240:118371. DOI
HC-ChiSep	Hollow-cylinder separation with fibre angle estimated from multi-echo magnitude	Wharton, S., Bowtell, R. (2012). "Fiber orientation-dependent white matter contrast in gradient echo MRI." PNAS, 109(45):18559-18564. DOI
R2*-QSM	Closed-form separation from QSM and R2*	Dimov, A.V., et al. (2022). "Magnetic susceptibility source separation solely from gradient echo data: histological validation." Tomography, 8(3):1544-1551. DOI
WaveSep	Wavelet-domain L1 separation from QSM and R2′	Fang, Z., et al. (2023). "Wavelet-based single-step susceptibility source separation." Proc. ISMRM.
DECOMPOSE	Per-voxel three-compartment fit of the complex multi-echo signal	Chen, J., et al. (2021). "Decompose quantitative susceptibility mapping (QSM) to sub-voxel diamagnetic and paramagnetic components based on gradient-echo MRI data." , 242:118477.

SWI Processing

Algorithm	Description	Reference	Parity
CLEAR-SWI	Multi-echo SWI: combined, bias-corrected magnitude × phase mask from the echo-combined, high-pass filtered phase; minimum intensity projection	Eckstein, K., et al. (2021). "Improved susceptibility weighted imaging at ultra-high field using bipolar multi-echo acquisition and optimized image processing: CLEAR-SWI." NeuroImage, 237:118175. doi:10.1016/j.neuroimage.2021.118175
SMWI	Magnitude weighted by a susceptibility-map mask	Gho, S.-M., et al. (2014). "Susceptibility map-weighted imaging (SMWI) for neuroimaging." Magnetic Resonance in Medicine, 72:337-346. DOI

Preprocessing and Artefact Reduction

Algorithm	Description	Reference	Parity
Bias correction	Receive-field homogeneity correction (makehomogeneous)	Eckstein, K., Trattnig, S., Robinson, S.D. (2019). "A Simple Homogeneity Correction for Neuroimaging at 7T." Proc. ISMRM 27th Annual Meeting.
N4 bias correction	N4ITK bias field correction as ANTs' N4BiasFieldCorrection (n4_bias_correction)	Tustison, N.J., et al. (2010). "N4ITK: Improved N3 Bias Correction." IEEE Transactions on Medical Imaging, 29(6):1310-1320. DOI
MP-PCA denoising	Marchenko–Pastur PCA denoising across echoes or volumes (mppca_denoise)	Veraart, J., et al. (2016). "Denoising of diffusion MRI using random matrix theory." NeuroImage, 142:394-406. DOI
Gibbs unringing	Local subvoxel-shift Gibbs-ringing removal (gibbs_unring): 3-D for 3-D acquisitions (default), 2-D slice-wise for 2-D multi-slice acquisitions (MRtrix3 mrdegibbs; DIPY as an option)	Kellner, E., et al. (2016). "Gibbs-ringing artifact removal based on local subvoxel-shifts." Magnetic Resonance in Medicine, 76(5):1574-1581. DOI; Bautista, T., et al. (2021). "Removal of Gibbs ringing artefacts for 3D acquisitions using subvoxel shifts." Proc. ISMRM 29:3535.
Tukey windowing	The SNU χ-separation toolbox's k-space Tukey window on each echo's complex image before processing (); taper ratio 0.4 for Siemens and GE, 0 for Philips

Utilities

Algorithm	Description	Reference	Parity
Frangi filter	Multi-scale 3D vesselness (frangi_filter_3d): Kroon's FrangiFilter3D in single precision, as QSMART runs it	Frangi, A.F., et al. (1998). "Multiscale vessel enhancement filtering." MICCAI'98, LNCS vol 1496, 130-137. DOI
Vasculature mask	Bottom-hat filtering then Frangi vesselness, for QSMART (generate_vasculature_mask), on QSMART's mean-echo, N4-corrected magnitude (qsmart_vessel_magnitude)	Yaghmaie, N., Syeda, W., et al. (2021). NeuroImage, 231:117701. DOI
Surface curvature	Discrete differential-geometry operators on triangle meshes (calculate_gaussian_curvature)	Meyer, M., et al. (2003). "Discrete Differential-Geometry Operators for Triangulated 2-Manifolds." Visualization and Mathematics III, 35-57. DOI
Rigid registration	6-DOF registration by normalised cross-correlation, with the B0 direction of each orientation (register_rigid)	—
Resampling / reorientation	Affine-derived B0 direction and obliquity, reorientation, and resampling (complex-domain for phase) to an axial or another grid (geometry)	—

Reference Implementations

This library was developed with reference to the following open-source implementations:

Repository	Algorithms	Language
mSMV	iSMV (iSMV.m, tdcsmv.m)	MATLAB
QSM.jl	SHARP, V-SHARP, SMV, iSMV, Laplacian unwrap, TKD, TSVD, Tikhonov, TV, RTS, NLTV, multi-echo linear fit	Julia
FANSI-toolbox	NDI, FANSI (nlTV/nlTGV), L1-QSM, WH-QSM	MATLAB
HD-QSM	HD-QSM	MATLAB
QSM_AMP_PE	AMP-PE	MATLAB
qsm_heidi	LSQR, HEIDI	MATLAB
QuantitativeSusceptibilityMappingTGV.jl	TGV	Julia
MriResearchTools.jl	ROMEO, MCPC-3D-S, B0 weighted average, R2*/T2*, bias correction, hole filling	Julia
MEDI toolbox (release 2024.11.26)	MEDI, MEDI+0, PDF, LBV, TFI, R2*/T2* (ARLO), complex fit	MATLAB
FSL-BET2	BET	C++
Rodent-Skull-Stripping

License

This project is licensed under the MIT License.