Skip to main content

cosmos_weighted

Function cosmos_weighted 

Source
pub fn cosmos_weighted<F: AsRef<[f64]>, W: AsRef<[f64]>>(
    local_fields: &[F],
    weights: &[W],
    bdirs: &[(f64, f64, f64)],
    mask: &[u8],
    grid: &Grid,
    params: &CosmosParams,
) -> Vec<f64>
Expand description

COSMOS dipole inversion with per-orientation spatial weighting.

This is the original formulation: each orientation’s residual is weighted by its signal magnitude, so voxels where that orientation had poor SNR (or was corrupted, or fell outside its own coverage) contribute less. Weighting is a spatial-domain operation and the dipole convolution is a k-space one, so the two no longer diagonalize together and the closed form of cosmos is unavailable. Instead the normal equations

    ( sum_t D_t W_t^2 D_t  +  lambda I ) chi  =  sum_t D_t W_t^2 f_t

are solved with conjugate gradient. The operator is symmetric positive semi-definite (D_t is a real, even k-space multiplier, hence self-adjoint), so CG applies directly. Each iteration costs two FFTs per orientation, making this substantially more expensive than cosmos — prefer the closed form unless the weighting earns its cost.

weights are multiplied by mask internally, so they need not be pre-masked. Typically they are the per-orientation magnitude images, normalized to a comparable scale across orientations (their relative size between orientations is what matters; a global factor only rescales lambda).

§Panics

As cosmos, and additionally if the weight count or any weight map size is wrong.