FIX: correctness fixes for MEDIC pipeline#1
Conversation
`warpkit.api.medic` exposes three outputs: `fieldmap` (undistorted-space, the consumer-facing output), `fieldmap_native` (a distorted-space debug output), and `displacement_map`. The Stage-3 wiring used `fieldmap_native`, which makes fmriprep's pull-resampler look up the field at the wrong physical location at each target voxel (it's off by the SDC displacement itself). `resample_vol` documents fmap_hz as "sampled to the target space, in Hz" (interfaces/resampling.py:272-273) — i.e., the undistorted/output grid. That's `fieldmap`, not `fieldmap_native`. Verified empirically by running `resample_series_async` against both fields on two datasets: * warpkit's 3.5mm fixture: RMS Δ ≈ 358, max|Δ| ≈ 17.5k vs mean signal ~7050 * ds007637 sub-04 at 2mm: RMS Δ ≈ 211, max|Δ| ≈ 21k vs mean signal ~4175 Relative error grows with resolution, as expected — a fixed physical displacement spans more voxels in a finer grid.
`func_fit_reports_wf.inputnode.fieldmap` flows into FieldmapReportlet (workflows/bold/outputs.py:373-401), which expects a 3D static field. The Stage-3 warpkit branch was passing the 4D HMC-aligned per-frame field directly, so the reportlet either crashes on the 4D NIfTI or silently renders frame 0 only. Pipe the temporal mean (already computed by `warpkit_mean` on line 526) instead. The 4D series remains the outputnode.fieldmap that the bold native workflow consumes for per-frame correction.
`mem_gb['resampled']` is the memory footprint of one resampled BOLD volume (`filesize * 4` in `estimate_bold_mem_usage`). MEDIC actually holds all N echoes of phase + magnitude in memory as float64, plus warpkit's internal state (unwrapped phases, masks, per-frame fields). For typical multi-echo data this underestimates real usage by an order of magnitude. ds007637 sub-04 is 5 echoes × ~320 MB compressed — close to 3 GB resident before warpkit's working set. The nipype scheduler will happily oversubscribe and OOM. Use `2 * n_echoes * filesize` instead, which roughly doubles the input footprint to cover outputs and warpkit's internal arrays.
The Stage-3 warpkit branch cleared `fmapref_buffer.sbref_files` and gated `raw_sbref_wf` with `not warpkit_enabled`. SBRef is acquired separately from the BOLD run and is orthogonal to MEDIC — dropping it silently is a regression for sites that acquire it. If there's a concrete incompatibility between SBRef and the MEDIC field-application path, it should be a stated limitation with a clear error, not an invisible behavior change. Removing the gate so SBRef flows through the same way it does on every other distortion-correction path.
`bold_series` is sorted by EchoTime in `workflows/base.py:250` and `phase_files` is independently sorted by EchoTime in `collect_bold_part_files`. They *should* line up, but a missing or zero-default EchoTime sidecar on one side and not the other silently transposes the pair set fed to warpkit. Add an explicit cross-check so that misalignment fails fast with a readable message instead of producing garbage fieldmaps.
vanandrew
changed the title
vanandrew
changed the title
|
@DVSneuro — these are five correctness fixes I found while testing your MEDIC add-on for fMRIPrep (nipreps#3645). Please give them a spin on your end before folding them into your PR; you have more real-world MEDIC runs queued up than I do and I'd rather you catch any regression before this lands upstream. Heads-up on scope: I originally had this branch also re-routing MEDIC through sdcflows' new |
|
@vanandrew -- Thanks for fixing these issues! I rebuilt and tested your branch on the multi-echo datasets I’ve been using for validation, and I did not spot any regressions. In this round, MEDIC ran successfully on ds002278/sub-Blossom, ds005085/sub-10017, ds006131/sub-20188, ds006926/sub-a01, and ds007637/sub-01. For ds005085/sub-10017 specifically, the single-echo runs followed the standard PHASEDIFF path, while the eligible 4-echo runs used MEDIC/warpkit. Thanks also for thinking carefully about the sdcflows scope. It sounds like routing MEDIC through sdcflows is likely the long-term architecture, but I agree it makes sense to keep these correctness fixes separate while nipreps/sdcflows#541 is still in flux. I’m happy to defer to you and the NiPreps maintainers on when/how to revisit that integration. |
2 participants
This PR is opened against your
warpkit-medicbranch (which is fMRIPrep PR nipreps#3645). It collects five independent correctness fixes to the existing MEDIC pipeline. Each commit is bisectable / revertable on its own.Bug fixes
FIX: feed undistorted MEDIC fieldmap (fmap) to the resamplerThe Stage-3 wiring fed
WarpkitMEDIC.fieldmap_nativeintoResampleSeries. That's warpkit's debug output — the field lives in distorted/native space, not the target/undistorted grid that fmriprep's pull-resampler expects (interfaces/resampling.py:272-273documentsfmap_hzas "sampled to the target space, in Hz").FIX: feed 3D temporal-mean field to FieldmapReportletfunc_fit_reports_wf.inputnode.fieldmapflows intoFieldmapReportlet, which expects a 3D field. Passing the 4D HMC-aligned series would either crash or silently render volume 0. Usewarpkit_mean.out_file(already computed) instead.FIX: scale warpkit_medic mem_gb by echo countmem_gb['resampled']is the size of one resampled BOLD volume. MEDIC holds N echoes × (phase + magnitude) in memory as float64 plus warpkit's internal state — ds007637 sub-04 is ~3 GB before warpkit's working set. Bump to2 * n_echoes * filesize.FIX: stop silently disabling SBRef when MEDIC is enabledfmapref_buffer.sbref_files = []and thenot warpkit_enabledgate onraw_sbref_wfsilently drop SBRef on MEDIC runs. SBRef is independent of MEDIC; if there's a real incompatibility, please make it a stated limitation rather than an invisible behavior change.FIX: assert echo ordering matches between magnitude and phasebold_seriesandphase_filesare sorted independently. A missing/zero EchoTime sidecar on one side silently transposes the pair set fed to warpkit. Add an explicit cross-check so misalignment fails fast.All five are independent of any sdcflows changes and ready to land on top of PR nipreps#3645 as-is.