Skip to main content

solve_poisson_dirichlet_roi

Function solve_poisson_dirichlet_roi 

Source
fn solve_poisson_dirichlet_roi(
    f: &[f64],
    mask: &[u8],
    grid: &Grid,
    tol: f64,
    max_iter: usize,
) -> Vec<f64>
Expand description

Laplacian phase unwrapping combined with background field removal.

Solves the Poisson equation with the Laplacian zeroed outside mask. That discards the field sources outside the ROI, and a field generated outside the ROI is harmonic inside it — so the background component is removed along with the wraps.

The result is not a total field. It is unwrapped and partially background-removed, by an amount that depends on the mask and the field geometry. Following this with a separate background-removal stage (V-SHARP, PDF, …) removes background twice.

Use laplacian_unwrap to unwrap without removing background.

Because ∇²(harmonic) = 0, the discarded component leaves no trace in the input to the Poisson solve and cannot be restored afterwards.

§Intended input

One wrapped phase volume. On the project’s test data it then matches [crate::bgremove::lbv] on the same field (r = 0.887 against 0.909, identical residual smooth content). It takes wrapped phase, so in a multi-echo pipeline the only way to apply it is to each echo before combining; that usage is not what the algorithm describes, has not been validated, and leaves visibly more background than either unwrapping then a field-map background removal or this function on a single volume. For multi-echo data use laplacian_unwrap or ROMEO, combine, then a background removal from crate::bgremove.

§Echo time

Accuracy falls off with the amount of phase to unwrap. On the 7 T test data, against the ground-truth local field: r = 0.84 at TE = 4 ms, 0.70 at 8 ms, 0.50 at 12 ms — where unwrapping then [crate::bgremove::lbv] gives 0.86, 0.84, 0.69 and [crate::bgremove::vsharp] holds near 0.82 throughout. Prefer the earliest echo.

§Arguments

  • phase - Wrapped phase (nx * ny * nz)
  • mask - Binary mask (nx * ny * nz), 1 = inside ROI
  • grid - Volume grid (dimensions and voxel sizes)

§Returns

Unwrapped, partially background-removed phase, zero outside mask.

§References

Schofield & Zhu (2003) for the unwrapping; Zhou et al. (2014) for the boundary-value formulation of the background removal. See the module docs. Solve ∇²u = f inside mask with u = 0 outside it (homogeneous Dirichlet on the ROI), by Gauss-Seidel with successive over-relaxation.

Masking the source term is only half of the ROI formulation: the solution has to be constrained at the ROI boundary too, or it picks up an arbitrary harmonic component. Solving the masked source over the whole volume with a periodic FFT does not constrain it, and that component is large — which is what this replaces.

Mirrors the solver in [crate::bgremove::lbv], which solves the homogeneous case (f = 0) with boundary values taken from the field.