Module sdf 

Spatially Dependent Filtering (SDF) for QSMART

SDF is the background field removal method used in QSMART. It uses variable-radius Gaussian filtering where the kernel size depends on the proximity to the brain boundary. This allows for aggressive filtering in the brain interior while preserving details near the surface.

The algorithm includes optional curvature-based weighting to further reduce artifacts at highly curved brain regions.

Reference pipeline

QSMART runs SDF twice on its echo-fit total field (QSMART_SDF.m): stage 1 on the BET mask × R_0 with no vessels and sdf_sigma1 (10) passed as both widths (l. 12), so σ = √(10² + 10²) ≈ 14.1 (sdf_curvature.m l. 8); stage 2 on the field × mask × R_0 with the vessel mask, σ1 = 8, σ2 = 2 (l. 20). Both use lowerLim 0.6, curvature constant 500 and an indent mask from closing the brain mask with strel('sphere', 8) voxels (Demo_QSMART.m l. 51–57). The paper instead lists σ_C = 10, σ_V = 2, σ_I = 10, σ = 10, C_s = 500 (Yaghmaie et al. 2021, §3.4); SdfParams follows the toolbox. The mask is not eroded first (paper §1, §5). Parameters were tuned on 7T human data at 0.6 mm; the Gaussian widths are in voxels (imgaussfilt3), so they are not rescaled for other voxel sizes, and the cortical proximity map uses [σ1, 2σ1, 2σ1] (sdf_curvature.m l. 12): half the width along the first voxel axis, whatever the orientation.

Reference: Yaghmaie, N., Syeda, W., et al. (2021). “QSMART: Quantitative Susceptibility Mapping Artifact Reduction Technique.” NeuroImage, 231:117701. https://doi.org/10.1016/j.neuroimage.2020.117701

Reference implementation: https://github.com/wtsyeda/QSMART

Agreement with the original

tests/qsmart_qsmci_parity.rs (measured 2026-11-02) runs QSMART’s stages on the inputs the verified QSM-CI qsmart-matlab v4 container (QSMxT/QSM-CI#296, #355, bit-identical to licensed MATLAB) gave its own, on a 36×36×28 vessel phantom at 3 T; the fixture is regenerated in Actions by tests/data/qsmart_qsmci_small.py. Given the toolbox’s total field and mask·R_0 (and stage 2 its vessel mask), sdf is 3.8e-5 % from QSMART_SDF stage 1 (max |Δ| 5e-9 ppm) and 1.6e-5 % from stage 2 (max |Δ| 2.7e-8 ppm): the fixture’s float32 storage. The 10.9 % and 2.1 % measured before were the curvature step’s triangulation and arithmetic (astewartau/QSM.rs#320, fixed by #370; see [crate::utils::curvature]). Stage 1’s widths are the toolbox’s (σ1 = σ2 = 10, above); with σ2 = 0, as recon.m’s s1.sdf_sigma2, it is 149 % off.

At full size (QSM-CI’s simulated head, 164×205×205, and KKI, 256×256×128; measured 2026-11-03, astewartau/QSM.rs#409) the local fields are 0.075 % / 0.20 % (sim, stage 1 / 2) and 0.13 % / 0.73 % (KKI) from the original’s. An instrumented copy of sdf_curvature.m and calculate_curvature.m under licensed MATLAB R2023b reproduces the container to 5e-15 %, and given MATLAB’s curvature map curvI, sdf’s smoothing widths are MATLAB’s in every voxel of both stages and both cases. All of the difference is the curvature weight of flat surface vertices (see [crate::utils::curvature]): 321 of 6,011 (sim) and 590 of 9,003 (KKI) get 1 where MATLAB gives 0 or the reverse, which moves 764 / 2,518 voxels to a neighbouring width in stage 1. The proximity’s rounding to two decimals and the vessel-proximity step add nothing.

Structs

SdfParams

Parameters for SDF background field removal

Functions

convolve_1d_direction 🔒

1D convolution along specified axis with replicate padding Matches MATLAB’s imgaussfilt3 default behavior

convolve_1d_direction_sigma 🔒

1D convolution with specified sigma

gaussian_smooth_3d_masked_f64 🔒

Gaussian smoothing with anisotropic sigma and mask

gaussian_smooth_3d_with_filter_size 🔒

3D Gaussian smoothing with specified filter size

sdf

SDF background field removal

sdf_curvature

SDF with curvature-based weighting (full QSMART pipeline)

sdf_simple

Simple SDF without curvature (faster, for testing)

smoothing_width 🔒

The Gaussian width SDF smooths a voxel with: sigma*round(prox^n, 2), as sdf_curvature.m computes it.