Module slicewise 

Slice-wise (2D) phase unwrapping for 2D multi-slice acquisitions.

In a 2D multi-slice acquisition each slice is excited separately and carries its own constant receive phase offset; interleaved acquisition adds slice-to-slice jumps on top. Every unwrapper in crate::unwrap reads across slice boundaries, and all of them are corrupted by that:

• unwrap_romeo and unwrap_bestpath grow regions in 3D, so one bad slice boundary carries a wrong offset across everything downstream of it.
• laplacian_unwrap differentiates along z straight across the jump. It degrades more gracefully — it is a global FFT solve, so the error stays diffuse rather than propagating — but it is still wrong.

On clean data the region growers look like they get away with it, and it is worth being precise about why they do not. Their quality weighting puts the incoherent z edges last, so a slice tends to be entered once and to come out a whole number of wraps from the truth — the error a multi-echo fit is supposed to absorb. The catch is which whole number: the jump the grower rounds is offset_step + field_step * TE, so it depends on the echo. The per-slice wrap then differs between echoes, the offset stops being echo-independent, and the fit that was going to cancel it is corrupted instead. correct_multi_echo_wraps cannot repair it, because it corrects the whole volume at once and the error is per slice. Measured on this module’s test phantom, neither 3D path survives both offset patterns on its own — the template path misreads random offsets, the individual path misreads interleaved ones.

unwrap_slicewise runs the chosen unwrapper on each slice independently, so z is never differentiated across. Each slice is handed to the same code with a grid of (d0, d1, 1); the ±z neighbours then fall outside the volume and are skipped by the existing bounds checks, so the unwrappers themselves are untouched.

Reference pipeline

There is no single published reference for this module; the measurements below are QSM.rs’s own. Related options in the ROMEO authors’ tools: CLEAR-SWI’s laplacianslice unwrapping runs Laplacian unwrapping and its high-pass filter per 2D slice (CLEARSWI.jl 0949bc5, src/phase_processing.jl:106-116; ext/ClearswiApp/argparse.jl:55-57), and the ROMEO CLI has --fix-ge-phase for GE phase with slice jumps (ROMEO.jl 08776d4, ext/RomeoApp/argparse.jl:116-119). The QSM pipeline’s field_mapping.unwrapping_algorithm does not offer slice-wise unwrapping.

What this costs

Unwrapping a slice on its own fixes it only up to a multiple of 2π: nothing in the slice says which multiple. That is inherent to the problem, not to this implementation — the information that ties slices together is exactly the information a 2D multi-slice acquisition does not record.

Multi-echo data recovers it. The slice offset is echo-independent, so in phi(TE) = phi0 + gamma*dB*TE it lands entirely in the intercept, and a fit across echoes returns the field from the slope with phi0 discarded — slice offset and leftover 2π alike. QSMxT exposes that as --b0-estimation linear-fit.

Where inter-echo consistency belongs

That rescue only works if the leftover 2π is the same for every echo of a slice. Unwrap each (slice, echo) independently and it is not: each one picks its own multiple, the intercept stops being echo-independent, and the fit is corrupted rather than saved.

So unwrap_slicewise_multi_echo repairs it here, in the unwrapper, rather than leaving it to the B0 fit. The reason is not preference — it is that the repair needs the wrapped input phase, and the unwrapper is the only stage that holds both that and its own output. Given both, the echo-to-echo evolution is wrap(p_e - p_{e-1}) with no unwrapping and no knowledge of the offset, which pins the relative 2π exactly (see enforce_inter_echo_consistency). A B0 fit sees only unwrapped phase and would have to infer the jumps from the fit residual, which is both weaker and specific to one estimator — QSMxT’s other B0 estimators would silently keep the corruption.

The division of labour that leaves is:

• this module guarantees the leftover 2*pi*k is the same for every echo of a slice;
• the B0 fit cancels it, together with the physical slice offset, because both are now echo-independent.

Which half does which job

Measured on the qsm-forward 2D phantom at 3 mm slices with interleaved offsets, the two halves fix different things, and it is worth being exact about which, because the obvious reading — that slice-wise unwrapping is what rescues the field map — is wrong.

The consistency pass fixes the fitted field, from any starting point. Correlation of a linear-fit B0 estimate against the field the signal was generated from: 3D ROMEO alone reaches 0.20, and the same 3D output put through enforce_inter_echo_consistency reaches 0.983. So does the raw wrapped phase with no spatial unwrapping at all. A fitted slope is blind to a constant per voxel, so whatever wrap state a voxel starts in lands in the intercept and the fit discards it; all the slope needs is that consecutive echoes differ by their measured evolution, which is exactly what the pass enforces. If a linear-fit B0 map is all you want, you do not need this module — you need the pass, which is public and works on any unwrapper’s output.

Slice-wise unwrapping fixes the phase itself. The clean way to see this is on the first echo, because enforce_inter_echo_consistency loops for e in 1..n_echoes and so provably cannot touch it — any difference there belongs to the spatial unwrapping and nothing else, with no null needed to rule the pass out. Counting in-plane neighbour pairs that jump more than π in echo 1 of that session: 1432 in the wrapped input, 1432 after the pass (identical, as the structure requires), 0 after slice-wise unwrapping.

By the last echo both mechanisms are in play and the gap narrows but holds: 1433 wrapped, 2142 with the pass alone, 796 with slice-wise. The pass alone is worse there than not unwrapping, because it propagates echo 1’s wrap state faithfully into every later echo.

Anything reading the unwrapped phase rather than its TE-slope — a single-echo field map, a non-linear B0 estimator, phase fed straight to background removal — needs the spatial unwrapping, and on 2D multi-slice data that means per slice.

Thick slices need this even without per-slice offsets

This module was built for slices carrying independent receive phase offsets, and it turns out not to need them. On the matched session with uniform receive phase, differing from a 3D acquisition only in being 3 mm thick, 3D ROMEO reaches 0.936 against slice-wise’s 0.983.

3 mm slices by themselves produce through-slice phase steps large enough to mislead a region grower, so the trigger is slice thickness rather than 2D multi-slice acquisition as such. Found by an assertion failing: the no-offsets control used to require 3D to reach 0.95, on the assumption that there is nothing to fix when the receive phase is uniform.

Slice-wise UnwrapMethod::Laplacian reaches only 0.81 on the fit — its Poisson solution is not a whole number of wraps from the truth to begin with, so the pass cannot re-seat it cleanly. Prefer ROMEO or best path.

Single-echo data has no second measurement and so no way to close the gap. unwrap_slicewise returns each slice unwrapped and internally consistent, and the between-slice offsets survive into the field map. That is a property of the acquisition; background field removal will absorb a smoothly varying part of it, but not jumps.

Slice gaps

This module does not care whether the slices are contiguous. Unwrapping is per slice and never relates one slice to another spatially, so it is valid on gapped data — unlike the dipole kernel, which is not (see crate::bgremove).

Structs

SliceWiseParams

Parameters for slice-wise unwrapping.

Constants

TWO_PI 🔒

Functions

enforce_inter_echo_consistency

Re-seat every echo on the first, so the echoes of a voxel differ by their measured phase evolution and nothing else.

label_components 🔒

Label the connected components of a slice mask, 4-connected.

reanchor_to_wrapped 🔒

Re-anchor a Laplacian solution so it differs from the wrapped phase by whole 2π.

unwrap_slicewise

Unwrap each slice independently, in plane.

unwrap_slicewise_multi_echo

Unwrap each echo slice-wise, then make each slice’s leftover 2π the same across echoes.