WIP

Author: effigies

smriprep/workflows/anatomical.py

@@ -57,6 +57,8 @@
 import smriprep
 
 from ..interfaces import DerivativesDataSink
+
+# from ..interfaces.calc import T1T2Ratio
 from ..utils.misc import apply_lut as _apply_bids_lut
 from ..utils.misc import fs_isRunning as _fs_isRunning
 from .fit.registration import init_register_template_wf
@@ -645,6 +647,7 @@ def init_anat_fit_wf(
                 'sphere_reg_fsLR',
                 'sphere_reg_msm',
                 'anat_ribbon',
+                't1t2_ratio',
                 # Reverse transform; not computable from forward transform
                 'std2anat_xfm',
                 # Metadata
@@ -1308,6 +1311,16 @@ def init_anat_fit_wf(
         LOGGER.info('ANAT Stage 8a: Found pre-computed cortical ribbon mask')
         outputnode.inputs.anat_ribbon = precomputed['anat_ribbon']
 
+    if t2w and 't1t2ratio' not in precomputed:
+        LOGGER.info('ANAT Stage 8b: Creating T1w/T2w ratio map')
+        # Commented out to pacify linter.
+        # t1t2_ratio = pe.Node(T1T2Ratio(), name='t1t2_ratio')
+
+    elif not t2w:
+        LOGGER.info('ANAT No T2w images provided - skipping Stage 8b')
+    else:
+        LOGGER.info('ANAT Found precomputed T1w/T2w ratio map - skipping Stage 8b')
+
     # Stage 9: Baseline fsLR registration
     if len(precomputed.get('sphere_reg_fsLR', [])) < 2:
         LOGGER.info('ANAT Stage 9: Creating fsLR registration sphere')

smriprep/workflows/surfaces.py

@@ -30,6 +30,7 @@
 
 import typing as ty
 
+from nibabel.processing import fwhm2sigma
 from nipype.interfaces import freesurfer as fs
 from nipype.interfaces import io as nio
 from nipype.interfaces import utility as niu
@@ -54,6 +55,7 @@
 )
 
 import smriprep
+from smriprep.interfaces.cifti import GenerateDScalar
 from smriprep.interfaces.surf import MakeRibbon
 from smriprep.interfaces.workbench import SurfaceResample
 
@@ -1449,8 +1451,6 @@ def init_morph_grayords_wf(
     """
     from niworkflows.engine.workflows import LiterateWorkflow as Workflow
 
-    from smriprep.interfaces.cifti import GenerateDScalar
-
     workflow = Workflow(name=name)
     workflow.__desc__ = f"""\
 *Grayordinate* "dscalar" files containing {grayord_density} samples were
@@ -1495,7 +1495,14 @@ def init_morph_grayords_wf(
         name='outputnode',
     )
 
-    for metric in ('curv', 'sulc', 'thickness'):
+    metrics = ['curv', 'sulc', 'thickness']
+    select_surfaces = pe.Node(
+        KeySelect(fields=metrics, keys=['L', 'R']),
+        name='select_surfaces',
+        run_without_submitting=True,
+    )
+
+    for metric in metrics:
         resample_and_mask_wf = init_resample_and_mask_wf(
             grayord_density=grayord_density,
             omp_nthreads=omp_nthreads,
@@ -1510,14 +1517,16 @@ def init_morph_grayords_wf(
         )
 
         workflow.connect([
+            (inputnode, select_surfaces, [(metric, metric)]),
+            (hemisource, select_surfaces, [('hemi', 'key')]),
             (inputnode, resample_and_mask_wf, [
-                (metric, 'inputnode.in_file'),
                 ('midthickness', 'inputnode.midthickness'),
                 ('midthickness_fsLR', 'inputnode.midthickness_fsLR'),
                 ('sphere_reg_fsLR', 'inputnode.sphere_reg_fsLR'),
                 ('roi', 'inputnode.cortex_mask'),
             ]),
             (hemisource, resample_and_mask_wf, [('hemi', 'inputnode.hemi')]),
+            (select_surfaces, resample_and_mask_wf, [(metric, 'inputnode.in_file')]),
             (resample_and_mask_wf, cifti_metric, [('outputnode.out_file', 'scalar_surfs')]),
             (cifti_metric, outputnode, [
                 ('out_file', f'{metric}_fsLR'),
@@ -1528,6 +1537,156 @@ def init_morph_grayords_wf(
     return workflow
 
 
+def init_myelinmap_fsLR_wf(
+    grayord_density: ty.Literal['91k', '170k'],
+    omp_nthreads: int,
+    mem_gb: float,
+    name: str = 'myelinmap_fsLR_wf',
+):
+    """Resample myelinmap volume to fsLR surface.
+
+    Workflow Graph
+        .. workflow::
+            :graph2use: colored
+            :simple_form: yes
+
+            from smriprep.workflows.surfaces import init_myelinmap_fsLR_wf
+            wf = init_myelinmap_fsLR_wf(grayord_density='91k', omp_nthreads=1, mem_gb=1)
+
+    Parameters
+    ----------
+    grayord_density : :class:`str`
+        Either ``"91k"`` or ``"170k"``, representing the total *grayordinates*.
+    omp_nthreads : :class:`int`
+        Maximum number of threads an individual process may use
+    mem_gb : :class:`float`
+        Size of BOLD file in GB
+    name : :class:`str`
+        Name of workflow (default: ``"myelinmap_fsLR_wf"``)
+
+    Inputs
+    ------
+    in_file : :class:`str`
+        Path to the myelin map in subject volume space
+    thickness : :class:`list` of :class:`str`
+        Path to left and right hemisphere thickness GIFTI shape files
+    midthickness : :class:`list` of :class:`str`
+        Path to left and right hemisphere midthickness GIFTI surface files
+    midthickness_fsLR : :class:`list` of :class:`str`
+        Path to left and right hemisphere midthickness GIFTI surface files in fsLR space
+    sphere_reg_fsLR : :class:`list` of :class:`str`
+        Path to left and right hemisphere sphere.reg GIFTI surface files,
+        mapping from subject to fsLR
+    cortex_mask : :class:`list` of :class:`str`
+        Path to left and right hemisphere cortex mask GIFTI files
+
+    Outputs
+    -------
+    out_fsLR : :class:`str`
+        Path to the resampled myelin map in fsLR space
+
+    """
+    from niworkflows.engine.workflows import LiterateWorkflow as Workflow
+    from niworkflows.interfaces.utility import KeySelect
+    from niworkflows.interfaces.workbench import VolumeToSurfaceMapping
+
+    workflow = Workflow(name=name)
+
+    inputnode = pe.Node(
+        niu.IdentityInterface(
+            fields=[
+                'in_file',
+                'thickness',
+                'midthickness',
+                'midthickness_fsLR',
+                'sphere_reg_fsLR',
+                'cortex_mask',
+                'volume_roi',
+            ]
+        ),
+        name='inputnode',
+    )
+
+    outputnode = pe.Node(
+        niu.IdentityInterface(fields=['out_file', 'out_metadata']),
+        name='outputnode',
+    )
+
+    hemisource = pe.Node(
+        niu.IdentityInterface(fields=['hemi']),
+        name='hemisource',
+        iterables=[('hemi', ['L', 'R'])],
+    )
+
+    select_surfaces = pe.Node(
+        KeySelect(
+            fields=[
+                'thickness',
+                'midthickness',
+            ],
+            keys=['L', 'R'],
+        ),
+        name='select_surfaces',
+        run_without_submitting=True,
+    )
+
+    volume_to_surface = pe.Node(
+        VolumeToSurfaceMapping(method='myelin-style', sigma=fwhm2sigma(5)),
+        name='volume_to_surface',
+        mem_gb=mem_gb * 3,
+        n_procs=omp_nthreads,
+    )
+    # smooth = pe.Node(
+    #     MetricSmooth(sigma=fwhm2sigma(4), nearest=True),
+    #     name='metric_dilate',
+    #     mem_gb=1,
+    #     n_procs=omp_nthreads,
+    # )
+    resample_and_mask_wf = init_resample_and_mask_wf(
+        grayord_density=grayord_density,
+        omp_nthreads=omp_nthreads,
+        mem_gb=mem_gb,
+    )
+    cifti_myelinmap = pe.JoinNode(
+        GenerateDScalar(grayordinates=grayord_density, scalar_name='MyelinMap'),
+        name='cifti_myelinmap',
+        joinfield=['scalar_surfs'],
+        joinsource='hemisource',
+    )
+
+    workflow.connect([
+        (inputnode, select_surfaces, [
+            ('midthickness', 'midthickness'),
+            ('thickness', 'thickness'),
+        ]),
+        (hemisource, select_surfaces, [('hemi', 'key')]),
+        # Resample volume to native surface
+        (inputnode, volume_to_surface, [
+            ('in_file', 'volume_file'),
+            ('ribbon_file', 'ribbon_roi'),
+        ]),
+        (select_surfaces, volume_to_surface, [
+            ('midthickness', 'surface_file'),
+            ('thickness', 'thickness'),
+        ]),
+        (inputnode, resample_and_mask_wf, [
+            ('midthickness', 'inputnode.midthickness'),
+            ('midthickness_fsLR', 'inputnode.midthickness_fsLR'),
+            ('sphere_reg_fsLR', 'inputnode.sphere_reg_fsLR'),
+            ('cortex_mask', 'inputnode.cortex_mask'),
+        ]),
+        (hemisource, resample_and_mask_wf, [('hemi', 'inputnode.hemi')]),
+        (volume_to_surface, resample_and_mask_wf, [('out_file', 'inputnode.in_file')]),
+        (resample_and_mask_wf, cifti_myelinmap, [('outputnode.out_file', 'scalar_surfs')]),
+        (cifti_myelinmap, outputnode, [
+            ('out_file', 'out_file'),
+            ('out_metadata', 'out_metadata'),
+        ]),
+    ])  # fmt:skip
+
+    return workflow
+
+
 def init_resample_and_mask_wf(
     grayord_density: ty.Literal['91k', '170k'],
     omp_nthreads: int,
@@ -1561,17 +1720,23 @@ def init_resample_and_mask_wf(
 
     Inputs
     ------
+    in_file : :class:`str`
+        Path to metric file in subject space
+    hemi : :class:`str`
+        Hemisphere identifier (``"L"`` or ``"R"``)
     midthickness : :class:`list` of :class:`str`
         Path to left and right hemisphere midthickness GIFTI surfaces.
     midthickness_fsLR : :class:`list` of :class:`str`
         Path to left and right hemisphere midthickness GIFTI surfaces in fsLR space.
     sphere_reg_fsLR : :class:`list` of :class:`str`
         Path to left and right hemisphere sphere.reg GIFTI surfaces, mapping from subject to fsLR
+    cortex_mask : :class:`list` of :class:`str`
+        Path to left and right hemisphere cortex mask GIFTI files
 
     Outputs
     -------
-    metric_fsLR : :class:`list` of :class:`str`
-        Path to metrics resampled as GIFTI files in fsLR space
+    metric_fsLR : :class:`str`
+        Path to metric resampled as GIFTI file in fsLR space
 
     """
     import templateflow.api as tf
@@ -1604,7 +1769,6 @@ def init_resample_and_mask_wf(
     select_surfaces = pe.Node(
         KeySelect(
             fields=[
-                'in_file',
                 'midthickness',
                 'midthickness_fsLR',
                 'sphere_reg_fsLR',
@@ -1646,15 +1810,14 @@ def init_resample_and_mask_wf(
 
     workflow.connect([
         (inputnode, select_surfaces, [
-            ('in_file', 'in_file'),
             ('midthickness', 'midthickness'),
             ('midthickness_fsLR', 'midthickness_fsLR'),
             ('sphere_reg_fsLR', 'sphere_reg_fsLR'),
             ('cortex_mask', 'cortex_mask'),
             ('hemi', 'key'),
         ]),
+        (inputnode, resample_to_fsLR, [('in_file', 'in_file')]),
         (select_surfaces, resample_to_fsLR, [
-            ('in_file', 'in_file'),
             ('sphere_reg_fsLR', 'current_sphere'),
             ('template_sphere', 'new_sphere'),
             ('midthickness', 'current_area'),