Merge pull request #2530 from effigies/fix/restore_sdc

FIX: Restore SyN-SDC

.circleci/config.yml

@@ -462,6 +462,12 @@ jobs:
          paths:
             - /tmp/ds005/work
 
+      - run:
+          name: Attempt run without PE metadata (should fail)
+          no_output_timeout: 2h
+          command: |
+            echo "TODO"
+
       - run:
           name: Run full fMRIPrep on ds005 (LegacyMultiProc plugin)
           no_output_timeout: 2h
@@ -470,6 +476,12 @@ jobs:
             if [ -f /tmp/.nofasttrack ]; then
                 FASTRACK_ARG=""
             fi
+
+            # Inject pretend metadata
+            json_sidecar=/tmp/data/${DATASET}/task-mixedgamblestask_bold.json
+            awk 'NR==1{print; print "    \"TotalReadoutTime\": 0.05,"} NR!=1' ${json_sidecar} > tmp && mv tmp ${json_sidecar}
+            awk 'NR==1{print; print "    \"PhaseEncodingDirection\": \"j\","} NR!=1' ${json_sidecar} > tmp && mv tmp ${json_sidecar}
+
             fmriprep-docker -i nipreps/fmriprep:latest \
                 -e FMRIPREP_DEV 1 --user $(id -u):$(id -g) \
                 --network none \
@@ -863,6 +875,12 @@ jobs:
          paths:
             - /tmp/ds210/work
 
+      - run:
+          name: Attempt run without PE metadata (should fail)
+          no_output_timeout: 2h
+          command: |
+            echo "TODO"
+
       - run:
           name: Run full fMRIPrep on ds000210
           no_output_timeout: 2h
@@ -871,6 +889,13 @@ jobs:
             if [ -f /tmp/.nofasttrack ]; then
                 FASTRACK_ARG=""
             fi
+
+            # Inject pretend metadata for SDCFlows not to crash
+            # TODO / open question - do all echos need the metadata?
+            chmod +w /tmp/data/${DATASET}
+            echo '{"PhaseEncodingDirection": "j", "TotalReadoutTime": 0.058}' >> /tmp/data/${DATASET}/task-cuedSGT_bold.json
+            chmod -R -w /tmp/data/${DATASET}
+
             fmriprep-docker -i nipreps/fmriprep:latest \
                 -e FMRIPREP_DEV 1 --user $(id -u):$(id -g) \
                 --config $PWD/nipype.cfg -w /tmp/${DATASET}/work \

.circleci/ds005_bids_fasttrack_outputs.txt

@@ -9,6 +9,12 @@ bids/logs/CITATION.html
 bids/logs/CITATION.md
 bids/logs/CITATION.tex
 bids/sub-01
+bids/sub-01/fmap
+bids/sub-01/fmap/sub-01_fmapid-auto00000_desc-coeff0_fieldmap.nii.gz
+bids/sub-01/fmap/sub-01_fmapid-auto00000_desc-coeff1_fieldmap.nii.gz
+bids/sub-01/fmap/sub-01_fmapid-auto00000_desc-magnitude_fieldmap.nii.gz
+bids/sub-01/fmap/sub-01_fmapid-auto00000_desc-preproc_fieldmap.json
+bids/sub-01/fmap/sub-01_fmapid-auto00000_desc-preproc_fieldmap.nii.gz
 bids/sub-01/func
 bids/sub-01/func/sub-01_task-mixedgamblestask_run-1_desc-confounds_timeseries.json
 bids/sub-01/func/sub-01_task-mixedgamblestask_run-1_desc-confounds_timeseries.tsv

.circleci/ds005_bids_outputs.txt

@@ -40,6 +40,12 @@ bids/sub-01/anat/sub-01_space-MNI152NLin2009cAsym_dseg.nii.gz
 bids/sub-01/anat/sub-01_space-MNI152NLin2009cAsym_label-CSF_probseg.nii.gz
 bids/sub-01/anat/sub-01_space-MNI152NLin2009cAsym_label-GM_probseg.nii.gz
 bids/sub-01/anat/sub-01_space-MNI152NLin2009cAsym_label-WM_probseg.nii.gz
+bids/sub-01/fmap
+bids/sub-01/fmap/sub-01_fmapid-auto00000_desc-coeff0_fieldmap.nii.gz
+bids/sub-01/fmap/sub-01_fmapid-auto00000_desc-coeff1_fieldmap.nii.gz
+bids/sub-01/fmap/sub-01_fmapid-auto00000_desc-magnitude_fieldmap.nii.gz
+bids/sub-01/fmap/sub-01_fmapid-auto00000_desc-preproc_fieldmap.json
+bids/sub-01/fmap/sub-01_fmapid-auto00000_desc-preproc_fieldmap.nii.gz
 bids/sub-01/func
 bids/sub-01/func/sub-01_task-mixedgamblestask_run-1_desc-confounds_timeseries.json
 bids/sub-01/func/sub-01_task-mixedgamblestask_run-1_desc-confounds_timeseries.tsv

.circleci/ds005_legacy_fasttrack_outputs.txt

@@ -9,6 +9,12 @@ fmriprep/logs/CITATION.html
 fmriprep/logs/CITATION.md
 fmriprep/logs/CITATION.tex
 fmriprep/sub-01
+fmriprep/sub-01/fmap
+fmriprep/sub-01/fmap/sub-01_fmapid-auto00000_desc-coeff0_fieldmap.nii.gz
+fmriprep/sub-01/fmap/sub-01_fmapid-auto00000_desc-coeff1_fieldmap.nii.gz
+fmriprep/sub-01/fmap/sub-01_fmapid-auto00000_desc-magnitude_fieldmap.nii.gz
+fmriprep/sub-01/fmap/sub-01_fmapid-auto00000_desc-preproc_fieldmap.json
+fmriprep/sub-01/fmap/sub-01_fmapid-auto00000_desc-preproc_fieldmap.nii.gz
 fmriprep/sub-01/func
 fmriprep/sub-01/func/sub-01_task-mixedgamblestask_run-1_AROMAnoiseICs.csv
 fmriprep/sub-01/func/sub-01_task-mixedgamblestask_run-1_desc-confounds_timeseries.json

.circleci/ds005_legacy_outputs.txt

@@ -50,6 +50,12 @@ fmriprep/sub-01/anat/sub-01_space-MNI152NLin6Asym_dseg.nii.gz
 fmriprep/sub-01/anat/sub-01_space-MNI152NLin6Asym_label-CSF_probseg.nii.gz
 fmriprep/sub-01/anat/sub-01_space-MNI152NLin6Asym_label-GM_probseg.nii.gz
 fmriprep/sub-01/anat/sub-01_space-MNI152NLin6Asym_label-WM_probseg.nii.gz
+fmriprep/sub-01/fmap
+fmriprep/sub-01/fmap/sub-01_fmapid-auto00000_desc-coeff0_fieldmap.nii.gz
+fmriprep/sub-01/fmap/sub-01_fmapid-auto00000_desc-coeff1_fieldmap.nii.gz
+fmriprep/sub-01/fmap/sub-01_fmapid-auto00000_desc-magnitude_fieldmap.nii.gz
+fmriprep/sub-01/fmap/sub-01_fmapid-auto00000_desc-preproc_fieldmap.json
+fmriprep/sub-01/fmap/sub-01_fmapid-auto00000_desc-preproc_fieldmap.nii.gz
 fmriprep/sub-01/func
 fmriprep/sub-01/func/sub-01_task-mixedgamblestask_run-1_AROMAnoiseICs.csv
 fmriprep/sub-01/func/sub-01_task-mixedgamblestask_run-1_desc-confounds_timeseries.json

.circleci/ds005_legacy_partial_fasttrack_outputs.txt

@@ -9,6 +9,12 @@ fmriprep/logs/CITATION.html
 fmriprep/logs/CITATION.md
 fmriprep/logs/CITATION.tex
 fmriprep/sub-01
+fmriprep/sub-01/fmap
+fmriprep/sub-01/fmap/sub-01_run-2_fmapid-auto00000_desc-coeff0_fieldmap.nii.gz
+fmriprep/sub-01/fmap/sub-01_run-2_fmapid-auto00000_desc-coeff1_fieldmap.nii.gz
+fmriprep/sub-01/fmap/sub-01_run-2_fmapid-auto00000_desc-magnitude_fieldmap.nii.gz
+fmriprep/sub-01/fmap/sub-01_run-2_fmapid-auto00000_desc-preproc_fieldmap.json
+fmriprep/sub-01/fmap/sub-01_run-2_fmapid-auto00000_desc-preproc_fieldmap.nii.gz
 fmriprep/sub-01/func
 fmriprep/sub-01/func/sub-01_task-mixedgamblestask_run-2_AROMAnoiseICs.csv
 fmriprep/sub-01/func/sub-01_task-mixedgamblestask_run-2_desc-confounds_timeseries.json

.circleci/ds005_legacy_partial_outputs.txt

@@ -50,6 +50,12 @@ fmriprep/sub-01/anat/sub-01_space-MNI152NLin6Asym_dseg.nii.gz
 fmriprep/sub-01/anat/sub-01_space-MNI152NLin6Asym_label-CSF_probseg.nii.gz
 fmriprep/sub-01/anat/sub-01_space-MNI152NLin6Asym_label-GM_probseg.nii.gz
 fmriprep/sub-01/anat/sub-01_space-MNI152NLin6Asym_label-WM_probseg.nii.gz
+fmriprep/sub-01/fmap
+fmriprep/sub-01/fmap/sub-01_run-2_fmapid-auto00000_desc-coeff0_fieldmap.nii.gz
+fmriprep/sub-01/fmap/sub-01_run-2_fmapid-auto00000_desc-coeff1_fieldmap.nii.gz
+fmriprep/sub-01/fmap/sub-01_run-2_fmapid-auto00000_desc-magnitude_fieldmap.nii.gz
+fmriprep/sub-01/fmap/sub-01_run-2_fmapid-auto00000_desc-preproc_fieldmap.json
+fmriprep/sub-01/fmap/sub-01_run-2_fmapid-auto00000_desc-preproc_fieldmap.nii.gz
 fmriprep/sub-01/func
 fmriprep/sub-01/func/sub-01_task-mixedgamblestask_run-02_AROMAnoiseICs.csv
 fmriprep/sub-01/func/sub-01_task-mixedgamblestask_run-02_desc-confounds_timeseries.json

.circleci/ds210_fasttrack_outputs.txt

@@ -7,6 +7,12 @@ fmriprep/logs/CITATION.html
 fmriprep/logs/CITATION.md
 fmriprep/logs/CITATION.tex
 fmriprep/sub-02
+fmriprep/sub-02/fmap
+fmriprep/sub-02/fmap/sub-02_run-1_fmapid-auto00000_desc-coeff0_fieldmap.nii.gz
+fmriprep/sub-02/fmap/sub-02_run-1_fmapid-auto00000_desc-coeff1_fieldmap.nii.gz
+fmriprep/sub-02/fmap/sub-02_run-1_fmapid-auto00000_desc-magnitude_fieldmap.nii.gz
+fmriprep/sub-02/fmap/sub-02_run-1_fmapid-auto00000_desc-preproc_fieldmap.json
+fmriprep/sub-02/fmap/sub-02_run-1_fmapid-auto00000_desc-preproc_fieldmap.nii.gz
 fmriprep/sub-02/func
 fmriprep/sub-02/func/sub-02_task-cuedSGT_run-1_desc-confounds_timeseries.json
 fmriprep/sub-02/func/sub-02_task-cuedSGT_run-1_desc-confounds_timeseries.tsv

.circleci/ds210_outputs.txt

@@ -26,6 +26,12 @@ fmriprep/sub-02/anat/sub-02_space-MNI152NLin2009cAsym_dseg.nii.gz
 fmriprep/sub-02/anat/sub-02_space-MNI152NLin2009cAsym_label-CSF_probseg.nii.gz
 fmriprep/sub-02/anat/sub-02_space-MNI152NLin2009cAsym_label-GM_probseg.nii.gz
 fmriprep/sub-02/anat/sub-02_space-MNI152NLin2009cAsym_label-WM_probseg.nii.gz
+fmriprep/sub-02/fmap
+fmriprep/sub-02/fmap/sub-02_run-1_fmapid-auto00000_desc-coeff0_fieldmap.nii.gz
+fmriprep/sub-02/fmap/sub-02_run-1_fmapid-auto00000_desc-coeff1_fieldmap.nii.gz
+fmriprep/sub-02/fmap/sub-02_run-1_fmapid-auto00000_desc-magnitude_fieldmap.nii.gz
+fmriprep/sub-02/fmap/sub-02_run-1_fmapid-auto00000_desc-preproc_fieldmap.json
+fmriprep/sub-02/fmap/sub-02_run-1_fmapid-auto00000_desc-preproc_fieldmap.nii.gz
 fmriprep/sub-02/func
 fmriprep/sub-02/func/sub-02_task-cuedSGT_run-1_desc-confounds_timeseries.json
 fmriprep/sub-02/func/sub-02_task-cuedSGT_run-1_desc-confounds_timeseries.tsv

fmriprep/workflows/base.py

@@ -311,27 +311,48 @@ def init_single_subject_wf(subject_id):
     if anat_only:
         return workflow
 
+    from sdcflows import fieldmaps as fm
     fmap_estimators = None
-    # TODO 21.0.0: Implement SyN
-    if any((config.workflow.use_syn_sdc, config.workflow.force_syn)):
-        config.loggers.workflow.critical("SyN processing is not yet implemented.")
 
-    if "fieldmaps" not in config.workflow.ignore:
+    if any(("fieldmaps" not in config.workflow.ignore,
+            config.workflow.use_syn_sdc,
+            config.workflow.force_syn)):
         from sdcflows.utils.wrangler import find_estimators
 
         # SDC Step 1: Run basic heuristics to identify available data for fieldmap estimation
         # For now, no fmapless
         fmap_estimators = find_estimators(
             layout=config.execution.layout,
             subject=subject_id,
-            fmapless=False,  # config.workflow.use_syn_sdc,
-            force_fmapless=False,  # config.workflow.force_syn,
+            fmapless=config.workflow.use_syn_sdc,
+            force_fmapless=config.workflow.force_syn,
         )
 
-        config.loggers.workflow.debug(
-            f"{len(fmap_estimators)} fieldmap estimators found: "
-            f"{[e.method for e in fmap_estimators]}"
-        )
+        if config.workflow.use_syn_sdc and not fmap_estimators:
+            message = ("Fieldmap-less (SyN) estimation was requested, but "
+                       "PhaseEncodingDirection information appears to be "
+                       "absent.")
+            config.loggers.workflow.error(message)
+            raise ValueError(message)
+
+        if (
+            "fieldmaps" in config.workflow.ignore
+            and [f for f in fmap_estimators
+                 if f.method != fm.EstimatorType.ANAT]
+        ):
+            config.loggers.workflow.info(
+                'Option "--ignore fieldmaps" was set, but either "--use-syn-sdc" '
+                'or "--force-syn" were given, so fieldmap-less estimation will be executed.'
+            )
+            fmap_estimators = [f for f in fmap_estimators
+                               if f.method == fm.EstimatorType.ANAT]
+
+        if fmap_estimators:
+            config.loggers.workflow.info(
+                "B0 field inhomogeneity map will be estimated with "
+                f" the following {len(fmap_estimators)} estimators: "
+                f"{[e.method for e in fmap_estimators]}."
+            )
 
     # Append the functional section to the existing anatomical exerpt
     # That way we do not need to stream down the number of bold datasets
@@ -373,7 +394,6 @@ def init_single_subject_wf(subject_id):
         return workflow
 
     from sdcflows.workflows.base import init_fmap_preproc_wf
-    from sdcflows import fieldmaps as fm
 
     fmap_wf = init_fmap_preproc_wf(
         debug="fieldmaps" in config.execution.debug,
@@ -412,6 +432,8 @@ def init_single_subject_wf(subject_id):
         config.loggers.workflow.info(f"""\
 Setting-up fieldmap "{estimator.bids_id}" ({estimator.method}) with \
 <{', '.join(s.path.name for s in estimator.sources)}>""")
+
+        # Mapped and phasediff can be connected internally by SDCFlows
         if estimator.method in (fm.EstimatorType.MAPPED, fm.EstimatorType.PHASEDIFF):
             continue
 
@@ -424,14 +446,56 @@ def init_single_subject_wf(subject_id):
             getattr(fmap_wf.inputs, f"in_{estimator.bids_id}").metadata = [
                 s.metadata for s in estimator.sources
             ]
-            continue
 
-        if estimator.method == fm.EstimatorType.PEPOLAR:
+        elif estimator.method == fm.EstimatorType.PEPOLAR:
             raise NotImplementedError(
                 "Sophisticated PEPOLAR schemes are unsupported."
             )
-        # TODO: SyN fieldmap processing
 
+        elif estimator.method == fm.EstimatorType.ANAT:
+            from niworkflows.interfaces.utility import KeySelect
+            from sdcflows.workflows.fit.syn import init_syn_preprocessing_wf
+
+            sources = [str(s.path) for s in estimator.sources if s.suffix == "bold"]
+            source_meta = [s.metadata for s in estimator.sources if s.suffix == "bold"]
+            syn_preprocessing_wf = init_syn_preprocessing_wf(
+                omp_nthreads=config.nipype.omp_nthreads,
+                debug=config.execution.sloppy,
+                auto_bold_nss=True,
+                t1w_inversion=False,
+                name=f"syn_preprocessing_{estimator.bids_id}",
+            )
+            syn_preprocessing_wf.inputs.inputnode.in_epis = sources
+            syn_preprocessing_wf.inputs.inputnode.in_meta = source_meta
+
+            # Select "MNI152NLin2009cAsym" from standard references.
+            fmap_select_std = pe.Node(
+                KeySelect(fields=["std2anat_xfm"], key="MNI152NLin2009cAsym"),
+                name="fmap_select_std",
+                run_without_submitting=True,
+            )
+
+            # fmt:off
+            workflow.connect([
+                (anat_preproc_wf, fmap_select_std, [
+                    ("outputnode.std2anat_xfm", "std2anat_xfm"),
+                    ("outputnode.template", "keys")]),
+                (anat_preproc_wf, syn_preprocessing_wf, [
+                    ("outputnode.t1w_preproc", "inputnode.in_anat"),
+                    ("outputnode.t1w_mask", "inputnode.mask_anat"),
+                ]),
+                (fmap_select_std, syn_preprocessing_wf, [
+                    ("std2anat_xfm", "inputnode.std2anat_xfm"),
+                ]),
+                (syn_preprocessing_wf, fmap_wf, [
+                    ("outputnode.epi_ref", f"in_{estimator.bids_id}.epi_ref"),
+                    ("outputnode.epi_mask", f"in_{estimator.bids_id}.epi_mask"),
+                    ("outputnode.anat_ref", f"in_{estimator.bids_id}.anat_ref"),
+                    ("outputnode.anat_mask", f"in_{estimator.bids_id}.anat_mask"),
+                    ("outputnode.sd_prior", f"in_{estimator.bids_id}.sd_prior"),
+                ]),
+            ])
+            # fmt:on
     return workflow
 
 

fmriprep/workflows/bold/base.py

@@ -252,18 +252,24 @@ def init_func_preproc_wf(bold_file, has_fieldmap=False):
             from sdcflows.fieldmaps import get_identifier
 
             # Fallback to IntendedFor
-            bold_rel = re.sub(
-                r"^sub-[a-zA-Z0-9]*/", "", str(Path(bold_file).relative_to(layout.root))
+            intended_rel = re.sub(
+                r"^sub-[a-zA-Z0-9]*/",
+                "",
+                str(Path(
+                    bold_file if not multiecho else bold_file[0]
+                ).relative_to(layout.root))
             )
-            estimator_key = get_identifier(bold_rel)
+            estimator_key = get_identifier(intended_rel)
 
         if not estimator_key:
             has_fieldmap = False
             config.loggers.workflow.critical(
-                f"None of the available B0 fieldmaps are associated to <{bold_rel}>"
+                f"None of the available B0 fieldmaps are associated to <{bold_file}>"
             )
         else:
-            config.loggers.workflow.info(f"Found usable B0 fieldmap <{estimator_key}>")
+            config.loggers.workflow.info(
+                f"Found usable B0-map (fieldmap) estimator(s) <{', '.join(estimator_key)}> "
+                f"to correct <{bold_file}> for susceptibility-derived distortions.")
 
     # Check whether STC must/can be run
     run_stc = (
@@ -1038,7 +1044,7 @@ def init_func_preproc_wf(bold_file, has_fieldmap=False):
         debug="fieldmaps" in config.execution.debug,
         omp_nthreads=config.nipype.omp_nthreads,
     )
-    unwarp_wf.inputs.inputnode.metadata = layout.get_metadata(str(bold_file))
+    unwarp_wf.inputs.inputnode.metadata = metadata
 
     output_select = pe.Node(
         KeySelect(fields=["fmap", "fmap_ref", "fmap_coeff", "fmap_mask", "sdc_method"]),
@@ -1139,7 +1145,13 @@ def init_func_preproc_wf(bold_file, has_fieldmap=False):
             ]
         ),
         joinsource=("meepi_echos" if run_stc is True else "boldbuffer"),
-        joinfield=["bold_files"],
+        joinfield=[
+            "fieldmap",
+            "fieldwarp",
+            "corrected",
+            "corrected_ref",
+            "corrected_mask",
+        ],
         name="join_sdc_echos",
     )
 
@@ -1163,7 +1175,7 @@ def _dpop(list_of_lists):
             ("corrected", "inputnode.bold_file"),
         ]),
         (join_sdc_echos, bold_t2s_wf, [
-            ("corrected_mask", "inputnode.bold_mask"),
+            (("corrected_mask", pop_file), "inputnode.bold_mask"),
         ]),
         (join_sdc_echos, bold_t1_trans_wf, [
             # TEMPORARY: For the moment we can't use frame-wise fieldmaps