Function decompose 

pub fn decompose(
    chi_per_echo: &[f64],
    magnitude: &[f64],
    echo_times: &[f64],
    mask: &[u8],
    params: &DecomposeParams,
    progress: impl FnMut(usize, usize),
) -> Result<ChiSepMaps, InputError>

DECOMPOSE source separation from per-echo QSM and the multi-echo magnitude.

Arguments

• chi_per_echo — QSM of each echo in ppm, flattened as (n_voxels, n_echoes) (echo fastest per voxel, like magnitude): the authors’ input (QSM_multi, a STAR-QSM per echo in DECOMPOSE_analysis.m). It sets each echo’s phase, −(2/3)·χ(TE)·γ·B0·TE.
• magnitude — Multi-echo magnitude, flattened as (n_voxels, n_echoes) in row-major order (echo fastest per voxel) — the same layout as crate::separation::r2star_qsm_from_magnitude. Normalised internally by its global maximum.
• echo_times — Echo times in seconds (n_echoes).
• mask — Binary brain mask (n_voxels, 1 = inside).
• params — See DecomposeParams.
• progress — Progress callback (voxels_done, voxels_total), called periodically over the fitted (masked) voxels.

Returns

(chi_pos, chi_neg, chi_total) in ppm, restricted to mask — matching the chi_sep_ilsqr/chi_sep_medi convention: chi_pos ≥ 0 (paramagnetic), chi_neg ≤ 0 (diamagnetic, signed), and chi_total = chi_pos + chi_neg.