update dwi workflows

Author: josephmje

dmriprep/workflows/dwi/__init__.py

@@ -6,20 +6,23 @@
 ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
 
 .. automodule:: dmriprep.workflows.dwi.base
-.. automodule:: dmriprep.workflows.dwi.prep_dwi
+.. automodule:: dmriprep.workflows.dwi.artifacts
+.. automodule:: dmriprep.workflows.dwi.eddy
 .. automodule:: dmriprep.workflows.dwi.tensor
 .. automodule:: dmriprep.workflows.dwi.outputs
 
 """
 
 from .base import init_dwi_preproc_wf
-from .prep_dwi import init_prep_dwi_wf
-from .tensor import init_tensor_wf
-from .outputs import init_output_wf
+from .artifacts import init_dwi_artifacts_wf
+from .eddy import init_dwi_eddy_wf
+from .tensor import init_dwi_tensor_wf
+from .outputs import init_dwi_derivatives_wf
 
 __all__ = [
     "init_dwi_preproc_wf",
-    "init_prep_dwi_wf",
-    "init_tensor_wf",
-    "init_output_wf"
+    "init_dwi_artifacts_wf",
+    "init_dwi_eddy_wf",
+    "init_dwi_tensor_wf",
+    "init_dwi_derivatives_wf"
 ]

dmriprep/workflows/dwi/artifacts.py

@@ -0,0 +1,99 @@
+# -*- coding: utf-8 -*-
+
+"""
+Artefact removal
+^^^^^^^^^^^^^^^^
+
+.. autofunction:: init_dwi_artifacts_wf
+
+"""
+
+from nipype.pipeline import engine as pe
+from nipype.interfaces import mrtrix3, utility as niu
+
+
+def init_dwi_artifacts_wf(ignore):
+    """
+    This workflow performs denoising and unringing on the input dwi image.
+
+    Denoising is done using Mrtrix3's implementation of the MP-PCA
+    algorithm `[Veraart2016a_]` and `[Veraart2016b_]`.
+    Unringing is done using Mrtrix3  `[Kellner2016_]`.
+
+    .. workflow::
+        :graph2use: orig
+        :simple_form: yes
+
+        from dmriprep.workflows.dwi import init_dwi_artifacts_wf
+        wf = init_dwi_artifacts_wf(ignore=[])
+
+    **Parameters**
+
+        ignore : list
+            List of artefact removal steps to skip (default: [])
+
+    **Inputs**
+
+        dwi_file
+            dwi NIfTI file
+
+    **Outputs**
+
+        out_file
+            dwi NIfTI file after artefact removal
+
+    References
+    ^^^^^^^^^^
+
+    .. [Veraart2016a] Veraart, J.; Novikov, D.S.; Christiaens, D.; Ades-aron, B.;
+                      Sijbers, J. & Fieremans, E. Denoising of diffusion MRI
+                      using random matrix theory. NeuroImage, 2016, 142, 394-406,
+                      doi: `10.1016/j.neuroimage.2016.08.016 <https://doi.org/10.1016/j.neuroimage.2016.08.016>`_.
+
+    .. [Veraart2016b] Veraart, J.; Fieremans, E. & Novikov, D.S. Diffusion MRI
+                      noise mapping using random matrix theory. Magn. Res. Med.,
+                      2016, 76(5), 1582-1593,
+                      doi: `10.1002/mrm.26059 <https://doi.org/10.1002/mrm.26059>`_.
+
+    .. [Kellner2016] Kellner, E; Dhital, B; Kiselev, V.G & Reisert, M.
+                     Gibbs-ringing artifact removal based on local subvoxel-shifts.
+                     Magnetic Resonance in Medicine, 2016, 76, 1574–1581,
+                     doi: `10.1002/mrm.26054 <https://doi.org/10.1002/mrm.26054>`_.
+
+    """
+
+    wf = pe.Workflow(name="dwi_artifacts_wf")
+
+    inputnode = pe.Node(niu.IdentityInterface(fields=["dwi_file"]), name="inputnode")
+
+    outputnode = pe.Node(niu.IdentityInterface(fields=["out_file"]), name="outputnode")
+
+    denoise = pe.Node(mrtrix3.DWIDenoise(), name="denoise")
+
+    unring = pe.Node(mrtrix3.MRDeGibbs(), name="unring")
+
+    if ignore == ["denoise"]:
+        wf.connect([
+            (inputnode, unring, [("dwi_file", "in_file")]),
+            (unring, outputnode, [("out_file", "out_file")])
+        ])
+
+    elif ignore == ["unring"]:
+        wf.connect([
+            (inputnode, denoise, [("dwi_file", "in_file")]),
+            (denoise, outputnode, [("out_file", "out_file")])
+        ])
+
+    elif ["denoise", "unring"] in ignore:
+        wf.connect([
+            (inputnode, outputnode, "dwi_file", "out_file")
+        ])
+
+    else:
+        wf.connect([
+            (inputnode, denoise, [("dwi_file", "in_file")]),
+            (denoise, unring, [("out_file", "in_file")]),
+            (unring, outputnode, [("out_file", "out_file")])
+        ])
+
+    return wf

dmriprep/workflows/dwi/base.py

@@ -10,49 +10,61 @@
 
 from bids import BIDSLayout
 from nipype.pipeline import engine as pe
-from nipype.interfaces import ants, fsl, mrtrix3, utility as niu
+from nipype.interfaces import fsl, mrtrix3, utility as niu
 from numba import cuda
 
-from .prep_dwi import init_prep_dwi_wf
-from .tensor import init_tensor_wf
+from .artifacts import init_dwi_artifacts_wf
+from .tensor import init_dwi_tensor_wf
 
-# from ..fieldmap.base import init_sdc_prep_wf
+from ..fieldmap.base import init_sdc_prep_wf
 
 FMAP_PRIORITY = {"epi": 0, "fieldmap": 1, "phasediff": 2, "phase": 3, "syn": 4}
 
 
-def init_dwi_preproc_wf(subject_id, dwi_file, metadata, parameters):
+def init_dwi_preproc_wf(
+    subject_id,
+    dwi_file,
+    metadata,
+    layout,
+    ignore,
+    b0_thresh,
+    output_resolution,
+    bet_dwi,
+    bet_mag,
+    omp_nthreads,
+    synb0_dir
+):
+
+    fmaps = []
+    synb0 = ""
+
+    # If use_synb0 set, get synb0 from files
+    if synb0_dir:
+        synb0_layout = BIDSLayout(
+            synb0_dir, validate=False, derivatives=True
+        )
+        synb0 = synb0_layout.get(subject=subject_id, return_type="file")[0]
+    else:
+        fmaps = layout.get_fieldmap(dwi_file, return_list=True)
+        if not fmaps:
+            raise Exception(
+                "No fieldmaps found for participant {}. "
+                "All workflows require fieldmaps".format(subject_id)
+            )
 
-    # fmaps = []
-    # synb0 = ""
-    #
-    # # If use_synb0 set, get synb0 from files
-    # if parameters.synb0_dir:
-    #     synb0_layout = BIDSLayout(
-    #         parameters.synb0_dir, validate=False, derivatives=True
-    #     )
-    #     synb0 = synb0_layout.get(subject=subject_id, return_type="file")[0]
-    # else:
-    #     fmaps = parameters.layout.get_fieldmap(dwi_file, return_list=True)
-    #     if not fmaps:
-    #         raise Exception(
-    #             "No fieldmaps found for participant {}. "
-    #             "All workflows require fieldmaps".format(subject_id)
-    #         )
-    #
-    #     for fmap in fmaps:
-    #         fmap["metadata"] = parameters.layout.get_metadata(
-    #             fmap[fmap["suffix"]]
-    #         )
-    #
-    # sdc_wf = init_sdc_prep_wf(
-    #     subject_id,
-    #     fmaps,
-    #     metadata,
-    #     parameters.layout,
-    #     parameters.bet_mag,
-    #     synb0,
-    # )
+        for fmap in fmaps:
+            fmap["metadata"] = layout.get_metadata(
+                fmap[fmap["suffix"]]
+            )
+
+    sdc_wf = init_sdc_prep_wf(
+        subject_id,
+        fmaps,
+        metadata,
+        layout,
+        bet_mag,
+        synb0,
+    )
 
     multiple_runs = isinstance(dwi_file, list)
 
@@ -102,7 +114,7 @@ def init_dwi_preproc_wf(subject_id, dwi_file, metadata, parameters):
         name="outputnode",
     )
 
-    dwi_prep_wf = init_prep_dwi_wf(parameters.ignore, parameters.output_resolution)
+    dwi_artifacts_wf = init_dwi_artifacts_wf(ignore)
 
     def gen_index(in_file):
         """
@@ -118,6 +130,7 @@ def gen_index(in_file):
 
         out_file
             The output index file.
+
         """
 
         import os
@@ -143,6 +156,23 @@ def gen_index(in_file):
     )
 
     def gen_acqparams(in_file, metadata):
+        """
+        Create an acquisition parameters file for ``eddy``
+
+        **Inputs**
+
+        in_file
+            The dwi file
+        metadata
+            The BIDS metadata of the dwi file
+
+        **Outputs**
+
+        out_file
+            The output acquisition parameters file
+
+        """
+
         import os
         from nipype.utils.filemanip import fname_presuffix
 
@@ -189,7 +219,7 @@ def gen_acqparams(in_file, metadata):
 
     def b0_average(in_dwi, in_bval, b0_thresh, out_file=None):
         """
-        A function that averages the *b0* volumes from a DWI dataset.
+        Averages the *b0* volumes from a DWI dataset.
         As current dMRI data are being acquired with all b-values > 0.0,
         the *lowb* volumes are selected by specifying the parameter b0_thresh.
         .. warning:: *b0* should be already registered (head motion artifact
@@ -229,16 +259,16 @@ def b0_average(in_dwi, in_bval, b0_thresh, out_file=None):
         name="b0_avg_pre",
     )
 
-    avg_b0_0.inputs.b0_thresh = parameters.b0_thresh
+    avg_b0_0.inputs.b0_thresh = b0_thresh
 
     bet_dwi0 = pe.Node(
-        fsl.BET(frac=parameters.bet_dwi, mask=True, robust=True), name="bet_dwi_pre"
+        fsl.BET(frac=bet_dwi, mask=True, robust=True), name="bet_dwi_pre"
     )
 
     ecc = pe.Node(fsl.Eddy(repol=True, cnr_maps=True, residuals=True), name="fsl_eddy")
 
-    if parameters.omp_nthreads:
-        ecc.inputs.num_threads = parameters.omp_nthreads
+    if omp_nthreads:
+        ecc.inputs.num_threads = omp_nthreads
 
     try:
         if cuda.gpus:
@@ -253,15 +283,6 @@ def b0_average(in_dwi, in_bval, b0_thresh, out_file=None):
     get_path = lambda x: x.split(".nii.gz")[0].split("_fix")[0]
     get_qc_path = lambda x: x.split(".nii.gz")[0] + ".qc"
 
-    if set(parameters.ignore) == set(["denoise", "unring"]):
-        wf.connect(
-            [
-                (inputnode, avg_b0_0, [("dwi_file", "in_file")])(
-                    inputnode, ecc, [("dwi_file", "in_file")]
-                )
-            ]
-        )
-
     fslroi = pe.Node(fsl.ExtractROI(t_min=0, t_size=1), name="fslroi")
 
     bias_correct = pe.Node(mrtrix3.DWIBiasCorrect(use_ants=True), name="bias_correct")
@@ -342,7 +363,7 @@ def get_b0_mask_fn(b0_file):
     )
 
     # # If synb0 is meant to be used
-    # if parameters.synb0_dir:
+    # if synb0_dir:
     #     dwi_wf.connect(
     #         [
     #             (
@@ -360,7 +381,7 @@ def get_b0_mask_fn(b0_file):
     #             ),
     #         ]
     #     )
-    #     ecc.inputs.in_acqp = parameters.acqp_file
+    #     ecc.inputs.in_acqp = acqp_file
     # else:
     #     # Decide what ecc will take: topup or fmap
     #     fmaps.sort(key=lambda fmap: FMAP_PRIORITY[fmap["suffix"]])