Function is_nyquist_bin 

fn is_nyquist_bin(
    i: usize,
    j: usize,
    k: usize,
    dims: (usize, usize, usize),
) -> bool

Whether a bin sits on a Nyquist plane, which only exists on even-length axes.

These bins are excluded from the reconstruction, for the same reason DC is: the model does not determine them.

At a Nyquist bin +f and -f alias onto the same sample, so fftfreq must pick one (it reports the negative) and the stored frequency no longer flips sign under the Hermitian pairing k -> -k. That matters because the STI row is not even in k component-wise — (H.k)/|k|^2 is linear in k, so flipping one axis while leaving the others alone changes the row. Away from Nyquist this never bites, since the Hermitian partner of a bin negates every axis and the row is even under a full negation. At a Nyquist bin the partner flips only the other axes, evenness fails, and the forward projection stops being Hermitian — after which the imaginary part discarded by the inverse FFT is real signal.

Restoring evenness by averaging the row over both signs does make the projection Hermitian, but the terms linear in that axis then cancel, leaving the six tensor components no longer separable at those bins — the system becomes rank-deficient rather than merely awkward. The information genuinely is not there, so these bins are zeroed rather than guessed at. They are the outermost sampled frequencies, where a tissue field map carries almost no energy; on a 164x205x205 volume exactly one axis is even, so this is one plane in 164.