gradunwarp wip

Author: bpinsard

fmriprep/cli/parser.py

@@ -279,6 +279,13 @@ def _slice_time_ref(value, parser):
         help="Ignore selected aspects of the input dataset to disable corresponding "
         "parts of the workflow (a space delimited list)",
     )
+    g_conf.add_argument(
+        "--gradunwarp-file",
+        metavar="PATH",
+        type=Path,
+        help="Path to vendor file for gradunwarp gradient distortion "
+        "correction.",
+    )
     g_conf.add_argument(
         "--output-spaces",
         nargs="*",

fmriprep/config.py

@@ -536,7 +536,9 @@ class workflow(_Config):
     ignore = None
     """Ignore particular steps for *fMRIPrep*."""
     longitudinal = False
-    """Run FreeSurfer ``recon-all`` with the ``-logitudinal`` flag."""
+    """Run FreeSurfer ``recon-all`` with the ``-longitudinal`` flag."""
+    gradunwarp_file = None
+    """Run GradUnwarp using the provided gradient or coefficient file"""
     medial_surface_nan = None
     """Fill medial surface with :abbr:`NaNs (not-a-number)` when sampling."""
     project_goodvoxels = False

fmriprep/workflows/base.py

@@ -319,6 +319,7 @@ def init_single_subject_wf(subject_id: str):
         spaces=spaces,
         t1w=subject_data['t1w'],
         t2w=subject_data['t2w'],
+        gradunwarp_file=config.workflow.gradunwarp_file,
         cifti_output=config.workflow.cifti_output,
     )
     # fmt:off

fmriprep/workflows/bold/base.py

@@ -111,6 +111,8 @@ def init_func_preproc_wf(bold_file, has_fieldmap=False):
         LTA-style affine matrix translating from T1w to FreeSurfer-conformed subject space
     fsnative2t1w_xfm
         LTA-style affine matrix translating from FreeSurfer-conformed subject space to T1w
+    gradunwarp_file
+        Vendor provided .coeff or .grad file for gradient distortion correction.
 
     Outputs
     -------
@@ -187,6 +189,7 @@ def init_func_preproc_wf(bold_file, has_fieldmap=False):
 
     """
     from niworkflows.engine.workflows import LiterateWorkflow as Workflow
+    from niworkflows.workflows.gradunwarp import init_gradunwarp_wf
     from niworkflows.func.util import init_bold_reference_wf
     from niworkflows.interfaces.nibabel import ApplyMask
     from niworkflows.interfaces.reportlets.registration import (
@@ -478,6 +481,14 @@ def init_func_preproc_wf(bold_file, has_fieldmap=False):
     # Top-level BOLD splitter
     bold_split = pe.Node(FSLSplit(dimension="t"), name="bold_split", mem_gb=mem_gb["filesize"] * 3)
 
+    # Gradient unwarping
+    if config.workflow.gradunwarp_file:
+        gradunwarp_wf = init_gradunwarp_wf()
+        gradunwarp_wf.inputs.inputnode.grad_file = gradunwarp_file
+        workflow.connect(
+            [(initial_boldref_wf, gradunwarp_wf, [('raw_ref_image', 'in_file')])]
+        )
+
     # HMC on the BOLD
     bold_hmc_wf = init_bold_hmc_wf(
         name="bold_hmc_wf", mem_gb=mem_gb["filesize"], omp_nthreads=omp_nthreads
@@ -505,6 +516,7 @@ def init_func_preproc_wf(bold_file, has_fieldmap=False):
         use_compression=False,
     )
     bold_t1_trans_wf.inputs.inputnode.fieldwarp = "identity"
+    bold_t1_trans_wf.inputs.inputnode.gradient_warp = "identity"
 
     # get confounds
     bold_confounds_wf = init_bold_confs_wf(
@@ -779,6 +791,7 @@ def init_func_preproc_wf(bold_file, has_fieldmap=False):
             use_compression=not config.execution.low_mem,
         )
         bold_std_trans_wf.inputs.inputnode.fieldwarp = "identity"
+        bold_std_trans_wf.inputs.inputnode.gradient_warp = "identity"
 
         # fmt:off
         workflow.connect([
@@ -1010,6 +1023,7 @@ def _last(inlist):
             name="bold_bold_trans_wf",
         )
         bold_bold_trans_wf.inputs.inputnode.fieldwarp = "identity"
+        bold_bold_trans_wf.inputs.inputnode.gradient_warp = "identity"
 
         # fmt:off
         workflow.connect([

fmriprep/workflows/bold/resampling.py

@@ -912,6 +912,7 @@ def init_bold_std_trans_wf(
                 "bold_split",
                 "t2star",
                 "fieldwarp",
+                "gradient_warp",
                 "hmc_xforms",
                 "itk_bold_to_t1",
                 "name_source",
@@ -964,7 +965,7 @@ def init_bold_std_trans_wf(
     )
 
     merge_xforms = pe.Node(
-        niu.Merge(4),
+        niu.Merge(5),
         name="merge_xforms",
         run_without_submitting=True,
         mem_gb=DEFAULT_MEMORY_MIN_GB,
@@ -994,7 +995,8 @@ def init_bold_std_trans_wf(
                                  ("templates", "keys")]),
         (inputnode, mask_std_tfm, [("bold_mask", "input_image")]),
         (inputnode, gen_ref, [(("bold_split", _first), "moving_image")]),
-        (inputnode, merge_xforms, [("hmc_xforms", "in4"),
+        (inputnode, merge_xforms, [("gradient_warp", "in5"),
+                                   ("hmc_xforms", "in4"),
                                    ("fieldwarp", "in3"),
                                    (("itk_bold_to_t1", _aslist), "in2")]),
         (inputnode, merge, [("name_source", "header_source")]),
@@ -1170,7 +1172,9 @@ def init_bold_preproc_trans_wf(
     )
 
     inputnode = pe.Node(
-        niu.IdentityInterface(fields=["name_source", "bold_file", "hmc_xforms", "fieldwarp"]),
+        niu.IdentityInterface(
+            fields=["name_source", "bold_file", "hmc_xforms", "fieldwarp", "gradient_warp"]
+        ),
         name="inputnode",
     )
 
@@ -1203,7 +1207,8 @@ def init_bold_preproc_trans_wf(
     # fmt:off
     workflow.connect([
         (inputnode, merge_xforms, [("fieldwarp", "in1"),
-                                   ("hmc_xforms", "in2")]),
+                                   ("hmc_xforms", "in2"),
+                                   ("gradient_warp", "in3")]),
         (inputnode, bold_transform, [("bold_file", "input_image"),
                                      (("bold_file", _first), "reference_image")]),
         (inputnode, merge, [("name_source", "header_source")]),