ENH: Merge T2w images and coregister to T1w template

smriprep/workflows/anatomical.py

@@ -36,6 +36,7 @@
 from nipype.interfaces.ants import N4BiasFieldCorrection
 
 from niworkflows.engine.workflows import LiterateWorkflow as Workflow
+from niworkflows.interfaces.fixes import FixHeaderApplyTransforms as ApplyTransforms
 from niworkflows.interfaces.freesurfer import (
     StructuralReference,
     PatchedLTAConvert as LTAConvert,
@@ -62,6 +63,7 @@ def init_anat_preproc_wf(
     hires,
     longitudinal,
     t1w,
+    t2w,
     omp_nthreads,
     output_dir,
     skull_strip_mode,
@@ -98,6 +100,7 @@ def init_anat_preproc_wf(
                 hires=True,
                 longitudinal=False,
                 t1w=['t1w.nii.gz'],
+                t2w=[],
                 omp_nthreads=1,
                 output_dir='.',
                 skull_strip_mode='force',
@@ -228,7 +231,8 @@ def init_anat_preproc_wf(
 
     outputnode = pe.Node(
         niu.IdentityInterface(
-            fields=["template", "subjects_dir", "subject_id"] + get_outputnode_spec()
+            fields=["template", "subjects_dir", "subject_id", "t2w_preproc"]
+            + get_outputnode_spec()
         ),
         name="outputnode",
     )
@@ -454,6 +458,7 @@ def _check_img(img):
         bids_root=bids_root,
         freesurfer=freesurfer,
         num_t1w=num_t1w,
+        t2w=t2w,
         output_dir=output_dir,
         spaces=spaces,
     )
@@ -474,6 +479,7 @@ def _check_img(img):
             ('t1w_mask', 'inputnode.t1w_mask'),
             ('t1w_dseg', 'inputnode.t1w_dseg'),
             ('t1w_tpms', 'inputnode.t1w_tpms'),
+            ('t2w_preproc', 'inputnode.t2w_preproc'),
         ]),
     ])
     # fmt:on
@@ -509,6 +515,7 @@ def _check_img(img):
                 (('outputnode.out_file', _pop), 't1w_brain'),
                 ('outputnode.out_mask', 't1w_mask')]),
         ])
+        # fmt:on
         return workflow
 
     # check for older IsRunning files and remove accordingly
@@ -522,6 +529,55 @@ def _check_img(img):
         name="surface_recon_wf", omp_nthreads=omp_nthreads, hires=hires
     )
     applyrefined = pe.Node(fsl.ApplyMask(), name="applyrefined")
+
+    if t2w:
+        t2w_template_wf = init_anat_template_wf(
+            longitudinal=longitudinal,
+            omp_nthreads=omp_nthreads,
+            num_t1w=len(t2w),
+            name="t2w_template_wf",
+        )
+        bbreg = pe.Node(
+            fs.BBRegister(
+                contrast_type="t2",
+                init="coreg",
+                dof=6,
+                out_lta_file=True,
+                args="--gm-proj-abs 2 --wm-proj-abs 1",
+            ),
+            name="bbreg",
+        )
+        coreg_xfms = pe.Node(niu.Merge(2), name="merge_xfms", run_without_submitting=True)
+        t2tot1_xfm = pe.Node(ConcatenateXFMs(), name="t2tot1_xfm", run_without_submitting=True)
+        t2w_resample = pe.Node(
+            ApplyTransforms(
+                dimension=3,
+                default_value=0,
+                float=True,
+                interpolation="LanczosWindowedSinc",
+            ),
+            name="t2w_resample",
+        )
+        # fmt:off
+        workflow.connect([
+            (inputnode, t2w_template_wf, [('t2w', 'inputnode.t1w')]),
+            (t2w_template_wf, bbreg, [('outputnode.t1w_ref', 'source_file')]),
+            (surface_recon_wf, bbreg, [
+                ('outputnode.subject_id', 'subject_id'),
+                ('outputnode.subjects_dir', 'subjects_dir'),
+            ]),
+            (bbreg, coreg_xfms, [('out_lta_file', 'in1')]),
+            (surface_recon_wf, coreg_xfms, [('outputnode.fsnative2t1w_xfm', 'in2')]),
+            (coreg_xfms, t2tot1_xfm, [('out', 'in_xfms')]),
+            (t2w_template_wf, t2w_resample, [('outputnode.t1w_ref', 'input_image')]),
+            (brain_extraction_wf, t2w_resample, [
+                (('outputnode.bias_corrected', _pop), 'reference_image'),
+            ]),
+            (t2tot1_xfm, t2w_resample, [('out_xfm', 'transforms')]),
+            (t2w_resample, outputnode, [('output_image', 't2w_preproc')]),
+        ])
+        # fmt:on
+
     # fmt:off
     workflow.connect([
         (inputnode, fs_isrunning, [

smriprep/workflows/base.py

@@ -308,6 +308,7 @@ def init_single_subject_wf(
         # for documentation purposes
         subject_data = {
             "t1w": ["/completely/made/up/path/sub-01_T1w.nii.gz"],
+            "t2w": [],
         }
     else:
         subject_data = collect_data(layout, subject_id, bids_filters=bids_filters)[0]
@@ -406,6 +407,7 @@ def init_single_subject_wf(
         longitudinal=longitudinal,
         name="anat_preproc_wf",
         t1w=subject_data["t1w"],
+        t2w=subject_data["t2w"],
         omp_nthreads=omp_nthreads,
         output_dir=output_dir,
         skull_strip_fixed_seed=skull_strip_fixed_seed,

smriprep/workflows/outputs.py

@@ -201,6 +201,7 @@ def init_anat_derivatives_wf(
     bids_root,
     freesurfer,
     num_t1w,
+    t2w,
     output_dir,
     spaces,
     name="anat_derivatives_wf",
@@ -241,6 +242,8 @@ def init_anat_derivatives_wf(
         Segmentation in T1w space
     t1w_tpms
         Tissue probability maps in T1w space
+    t2w_preproc
+        The preprocessed T2w image, bias-corrected and resampled into anatomical space.
     anat2std_xfm
         Nonlinear spatial transform to resample imaging data given in anatomical space
         into standard space.
@@ -284,6 +287,7 @@ def init_anat_derivatives_wf(
                 "t1w_mask",
                 "t1w_dseg",
                 "t1w_tpms",
+                "t2w_preproc",
                 "anat2std_xfm",
                 "std2anat_xfm",
                 "t1w2fsnative_xfm",
@@ -558,6 +562,27 @@ def init_anat_derivatives_wf(
     if not freesurfer:
         return workflow
 
+    # T2w coregistration requires FreeSurfer surfaces, so only try to save if freesurfer
+    if t2w:
+        ds_t2w_preproc = pe.Node(
+            DerivativesDataSink(
+                base_directory=output_dir,
+                desc="preproc",
+                suffix="T2w",
+                compress=True,
+            ),
+            name="ds_t2w_preproc",
+            run_without_submitting=True,
+        )
+        ds_t2w_preproc.inputs.SkullStripped = False
+        ds_t2w_preproc.inputs.source_file = t2w
+
+        # fmt:off
+        workflow.connect([
+            (inputnode, ds_t2w_preproc, [('t2w_preproc', 'in_file')]),
+        ])
+        # fmt:on
+
     from niworkflows.interfaces.nitransforms import ConcatenateXFMs
     from niworkflows.interfaces.surf import Path2BIDS