ENH: Add reporting workflow for BOLD fit

Author: effigies

fmriprep/workflows/bold/fit.py

@@ -48,7 +48,12 @@
 
 # BOLD workflows
 from .hmc import init_bold_hmc_wf
-from .outputs import init_ds_boldref_wf, init_ds_hmc_wf, init_ds_registration_wf
+from .outputs import (
+    init_ds_boldref_wf,
+    init_ds_hmc_wf,
+    init_ds_registration_wf,
+    init_func_fit_reports_wf,
+)
 from .reference import init_raw_boldref_wf
 from .registration import init_bold_reg_wf
 
@@ -183,12 +188,32 @@ def init_bold_fit_wf(
         niu.IdentityInterface(fields=["boldref", "boldmask"]), name="regref_buffer"
     )
 
+    func_fit_reports_wf = init_func_fit_reports_wf(
+        # TODO: Enable sdc report even if we find coregref
+        sdc_correction=not (have_coregref or fieldmap_id is None),
+        freesurfer=config.workflow.run_reconall,
+        output_dir=config.execution.fmriprep_dir,
+    )
+
     # fmt:off
     workflow.connect([
         (hmcref_buffer, outputnode, [("boldref", "hmc_boldref")]),
         (regref_buffer, outputnode, [("boldref", "coreg_boldref"),
                                      ("boldmask", "bold_mask")]),
         (hmc_buffer, outputnode, [("hmc_xforms", "motion_xfm")]),
+        (inputnode, func_fit_reports_wf, [
+            ("bold_file", "inputnode.source_file"),
+            ("t1w_preproc", "inputnode.t1w_preproc"),
+            # May not need all of these
+            ("t1w_mask", "inputnode.t1w_mask"),
+            ("t1w_dseg", "inputnode.t1w_dseg"),
+            ("subjects_dir", "inputnode.subjects_dir"),
+            ("subject_id", "inputnode.subject_id"),
+        ]),
+        (outputnode, func_fit_reports_wf, [
+            ("coreg_boldref", "inputnode.coreg_boldref"),
+            ("boldref2anat_xfm", "inputnode.boldref2anat_xfm"),
+        ]),
     ])
     # fmt:on
 
@@ -275,6 +300,7 @@ def init_bold_fit_wf(
             (fmapref_buffer, enhance_boldref_wf, [("out", "inputnode.in_file")]),
             (fmapref_buffer, ds_coreg_boldref_wf, [("out", "inputnode.source_files")]),
             (ds_coreg_boldref_wf, regref_buffer, [("outputnode.boldref", "boldref")]),
+            (fmapref_buffer, func_fit_reports_wf, [("out", "inputnode.sdc_boldref")]),
         ])
         # fmt:on
 

fmriprep/workflows/bold/outputs.py

@@ -28,6 +28,7 @@
 import numpy as np
 from nipype.interfaces import utility as niu
 from nipype.pipeline import engine as pe
+from niworkflows.interfaces.fixes import FixHeaderApplyTransforms as ApplyTransforms
 from smriprep.workflows.outputs import _bids_relative
 
 from fmriprep import config
@@ -140,6 +141,169 @@ def prepare_timing_parameters(metadata: dict):
     return timing_parameters
 
 
+def init_func_fit_reports_wf(
+    *,
+    sdc_correction: bool,
+    freesurfer: bool,
+    output_dir: str,
+    name="func_fit_reports_wf",
+) -> pe.Workflow:
+    """
+    Set up a battery of datasinks to store reports in the right location.
+
+    Parameters
+    ----------
+    freesurfer : :obj:`bool`
+        FreeSurfer was enabled
+    output_dir : :obj:`str`
+        Directory in which to save derivatives
+    name : :obj:`str`
+        Workflow name (default: anat_reports_wf)
+
+    Inputs
+    ------
+    source_file
+        Input BOLD images
+
+    std_t1w
+        T1w image resampled to standard space
+    std_mask
+        Mask of skull-stripped template
+    subject_dir
+        FreeSurfer SUBJECTS_DIR
+    subject_id
+        FreeSurfer subject ID
+    t1w_conform_report
+        Conformation report
+    t1w_preproc
+        The T1w reference map, which is calculated as the average of bias-corrected
+        and preprocessed T1w images, defining the anatomical space.
+    t1w_dseg
+        Segmentation in T1w space
+    t1w_mask
+        Brain (binary) mask estimated by brain extraction.
+    template
+        Template space and specifications
+
+    """
+    from niworkflows.interfaces.reportlets.masks import ROIsPlot
+    from niworkflows.interfaces.reportlets.registration import (
+        SimpleBeforeAfterRPT as SimpleBeforeAfter,
+    )
+
+    workflow = pe.Workflow(name=name)
+
+    inputfields = [
+        "source_file",
+        "sdc_boldref",
+        "coreg_boldref",
+        "boldref2anat_xfm",
+        "t1w_preproc",
+        "t1w_mask",
+        "t1w_dseg",
+        # May be missing
+        "subject_id",
+        "subjects_dir",
+    ]
+    inputnode = pe.Node(niu.IdentityInterface(fields=inputfields), name="inputnode")
+
+    # Reportlets follow the structure of init_bold_fit_wf stages
+    # - SDC1:
+    #       Before: Pre-SDC boldref
+    #       After: Fieldmap reference resampled on boldref
+    #       Three-way: Fieldmap resampled on boldref
+    # - SDC2:
+    #       Before: Pre-SDC boldref with white matter mask
+    #       After: Post-SDC boldref with white matter mask
+    # - EPI-T1 registration:
+    #       Before: T1w brain with white matter mask
+    #       After: Resampled boldref with white matter mask
+
+    if sdc_correction:
+        # SDC2
+        sdc_report = pe.Node(
+            SimpleBeforeAfter(
+                before_label="Distorted",
+                after_label="Corrected",
+                dismiss_affine=True,
+            ),
+            name="sdc_report",
+            mem_gb=0.1,
+        )
+
+        ds_sdc_report = pe.Node(
+            DerivativesDataSink(
+                base_directory=output_dir,
+                desc="sdc",
+                suffix="bold",
+                datatype="figures",
+                dismiss_entities=("echo",),
+            ),
+            name="ds_sdc_report",
+        )
+
+        # fmt:off
+        workflow.connect([
+            (inputnode, sdc_report, [
+                ('sdc_boldref', 'before'),
+                ('coreg_boldref', 'after'),
+            ]),
+            (inputnode, ds_sdc_report, [('source_file', 'source_file')]),
+            (sdc_report, ds_sdc_report, [('out_report', 'in_file')]),
+        ])
+        # fmt:on
+
+    # EPI-T1 registration
+    # Resample T1w image onto EPI-space
+    t1w_boldref = pe.Node(
+        ApplyTransforms(
+            dimension=3,
+            default_value=0,
+            float=True,
+            invert_transform_flags=[True],
+            interpolation="LanczosWindowedSinc",
+        ),
+        name="t1w_boldref",
+        mem_gb=1,
+    )
+
+    epi_t1_report = pe.Node(
+        SimpleBeforeAfter(
+            before_label="T1w",
+            after_label="EPI",
+        ),
+        name="epi_t1_report",
+        mem_gb=0.1,
+    )
+
+    ds_epi_t1_report = pe.Node(
+        DerivativesDataSink(
+            base_directory=output_dir,
+            desc="coreg",
+            suffix="bold",
+            datatype="figures",
+            dismiss_entities=("echo",),
+        ),
+        name="ds_epi_t1_report",
+    )
+
+    # fmt:off
+    workflow.connect([
+        (inputnode, t1w_boldref, [
+            ('t1w_preproc', 'input_image'),
+            ('coreg_boldref', 'reference_image'),
+            ('boldref2anat_xfm', 'transforms'),
+        ]),
+        (inputnode, epi_t1_report, [('coreg_boldref', 'after')]),
+        (t1w_boldref, epi_t1_report, [('output_image', 'before')]),
+        (inputnode, ds_epi_t1_report, [('source_file', 'source_file')]),
+        (epi_t1_report, ds_epi_t1_report, [('out_report', 'in_file')]),
+    ])
+    # fmt:on
+
+    return workflow
+
+
 def init_ds_boldref_wf(
     *,
     bids_root,