Work on CIFTI conversion workflow (#10)

.zenodo.json

@@ -33,11 +33,6 @@
       "identifier": "https://smripost_linc.org",
       "relation": "documents",
       "scheme": "url"
-    },
-    {
-      "identifier": "10.1038/s41592-018-0235-4",
-      "relation": "isPartOf",
-      "scheme": "doi"
     }
   ],
   "upload_type": "software"

src/smripost_linc/interfaces/freesurfer.py

@@ -2,6 +2,7 @@
 
 import os
 import shutil
+from glob import glob
 
 from nipype.interfaces.base import (
     Directory,
@@ -120,3 +121,48 @@ def _run_interface(self, runtime):
         shutil.copyfile(self.inputs.in_file, out_file)
 
         return runtime
+
+
+class _CollectFSAverageSurfacesInputSpec(TraitedSpec):
+    freesurfer_dir = Directory(
+        exists=True,
+        mandatory=True,
+        desc='FreeSurfer directory',
+    )
+
+
+class _CollectFSAverageSurfacesOutputSpec(TraitedSpec):
+    lh_fsaverage_files = traits.List(
+        File(exists=True),
+        desc='Left-hemisphere fsaverage-space surfaces',
+    )
+    rh_fsaverage_files = traits.List(
+        File(exists=True),
+        desc='Right-hemisphere fsaverage-space surfaces',
+    )
+    names = traits.List(
+        traits.Str,
+        desc='Names of collected surfaces',
+    )
+
+
+class CollectFSAverageSurfaces(SimpleInterface):
+    input_spec = _CollectFSAverageSurfacesInputSpec
+    output_spec = _CollectFSAverageSurfacesOutputSpec
+
+    def _run_interface(self, runtime):
+        in_dir = os.path.join(
+            self.inputs.freesurfer_dir,
+            'surf',
+        )
+        lh_mgh_files = sorted(glob(os.path.join(in_dir, 'lh.*.fsaverage.mgh')))
+        self._results['lh_fsaverage_files'] = lh_mgh_files
+        self._results['names'] = []
+        self._results['rh_fsaverage_files'] = []
+        for lh_file in lh_mgh_files:
+            name = os.path.basename(lh_file).split('.')[1]
+            self._results['names'].append(name)
+            rh_file = os.path.join(in_dir, f'rh.{name}.fsaverage.mgh')
+            self._results['rh_fsaverage_files'].append(rh_file)
+
+        return runtime

src/smripost_linc/interfaces/misc.py

@@ -1,4 +1,6 @@
-"""Miscellaneous interfaces for fmriprep-aroma."""
+"""Miscellaneous interfaces for sMRIPost-LINC."""
+
+import os
 
 import numpy as np
 from nipype.interfaces.base import (
@@ -17,6 +19,8 @@
     _FixTraitApplyTransformsInputSpec,
 )
 
+from smripost_linc.utils.utils import split_filename
+
 
 class _ApplyTransformsInputSpec(_FixTraitApplyTransformsInputSpec):
     # Nipype's version doesn't have GenericLabel
@@ -156,6 +160,98 @@ class CiftiSeparateMetric(WBCommand):
     _cmd = 'wb_command  -cifti-separate'
 
 
+class _CiftiCreateDenseScalarInputSpec(_WBCommandInputSpec):
+    """Input specification for the CiftiSeparateVolumeAll command."""
+
+    out_file = File(
+        exists=False,
+        mandatory=False,
+        genfile=True,
+        argstr='%s',
+        position=0,
+        desc='The CIFTI output.',
+    )
+    left_metric = File(
+        exists=True,
+        mandatory=False,
+        argstr='-left-metric %s',
+        position=1,
+        desc='The input surface data from the left hemisphere.',
+    )
+    right_metric = File(
+        exists=True,
+        mandatory=False,
+        argstr='-right-metric %s',
+        position=2,
+        desc='The input surface data from the right hemisphere.',
+    )
+    volume_data = File(
+        exists=True,
+        mandatory=False,
+        argstr='-volume %s',
+        position=3,
+        desc='The input volumetric data.',
+    )
+    structure_label_volume = File(
+        exists=True,
+        mandatory=False,
+        argstr='%s',
+        position=4,
+        desc='A label file indicating the structure of each voxel in volume_data.',
+    )
+
+
+class _CiftiCreateDenseScalarOutputSpec(TraitedSpec):
+    """Output specification for the CiftiCreateDenseScalar command."""
+
+    out_file = File(exists=True, desc='output CIFTI file')
+
+
+class CiftiCreateDenseScalar(WBCommand):
+    """Extract volumetric data from CIFTI file (.dtseries).
+
+    Other structures can also be extracted.
+    The input cifti file must have a brain models mapping on the chosen
+    dimension, columns for .dtseries,
+
+    Examples
+    --------
+    >>> cifticreatedensescalar = CiftiCreateDenseScalar()
+    >>> cifticreatedensescalar.inputs.out_file = 'sub_01_task-rest.dscalar.nii'
+    >>> cifticreatedensescalar.inputs.left_metric = 'sub_01_task-rest_hemi-L.func.gii'
+    >>> cifticreatedensescalar.inputs.left_metric = 'sub_01_task-rest_hemi-R.func.gii'
+    >>> cifticreatedensescalar.inputs.volume_data = 'sub_01_task-rest_subcortical.nii.gz'
+    >>> cifticreatedensescalar.inputs.structure_label_volume = 'sub_01_task-rest_labels.nii.gz'
+    >>> cifticreatedensescalar.cmdline
+    wb_command -cifti-create-dense-scalar 'sub_01_task-rest.dscalar.nii' \
+        -left-metric 'sub_01_task-rest_hemi-L.func.gii' \
+        -right-metric 'sub_01_task-rest_hemi-R.func.gii' \
+        -volume-data 'sub_01_task-rest_subcortical.nii.gz' 'sub_01_task-rest_labels.nii.gz'
+    """
+
+    input_spec = _CiftiCreateDenseScalarInputSpec
+    output_spec = _CiftiCreateDenseScalarOutputSpec
+    _cmd = 'wb_command -cifti-create-dense-scalar'
+
+    def _gen_filename(self, name):
+        if name != 'out_file':
+            return None
+
+        if isdefined(self.inputs.out_file):
+            return self.inputs.out_file
+        elif isdefined(self.inputs.volume_data):
+            _, fname, _ = split_filename(self.inputs.volume_data)
+        else:
+            _, fname, _ = split_filename(self.inputs.left_metric)
+
+        return f'{fname}_converted.dscalar.nii'
+
+    def _list_outputs(self):
+        outputs = self.output_spec().get()
+        outputs['out_file'] = os.path.abspath(self._gen_filename('out_file'))
+        return outputs
+
+
 class _ParcellationStats2TSVInputSpec(DynamicTraitedSpec):
     in_file = File(exists=True, mandatory=True, desc='parcellated data')
     hemisphere = traits.Enum('lh', 'rh', usedefault=True, desc='hemisphere')

src/smripost_linc/utils/utils.py

@@ -241,3 +241,82 @@ def list_to_str(lst):
         return ' and '.join(lst_str)
     else:
         return f"{', '.join(lst_str[:-1])}, and {lst_str[-1]}"
+
+
+def split_filename(fname):
+    """Split a filename into parts: path, base filename and extension.
+
+    Parameters
+    ----------
+    fname : :obj:`str`
+        file or path name
+
+    Returns
+    -------
+    pth : :obj:`str`
+        base path from fname
+    fname : :obj:`str`
+        filename from fname, without extension
+    ext : :obj:`str`
+        file extension from fname
+
+    Examples
+    --------
+    >>> from nipype.utils.filemanip import split_filename
+    >>> pth, fname, ext = split_filename('/home/data/subject.nii.gz')
+    >>> pth
+    '/home/data'
+
+    >>> fname
+    'subject'
+
+    >>> ext
+    '.nii.gz'
+    """
+    # TM 07152022 - edited to add cifti and workbench extensions
+    special_extensions = [
+        '.nii.gz',
+        '.tar.gz',
+        '.niml.dset',
+        '.dconn.nii',
+        '.dlabel.nii',
+        '.dpconn.nii',
+        '.dscalar.nii',
+        '.dtseries.nii',
+        '.fiberTEMP.nii',
+        '.trajTEMP.wbsparse',
+        '.pconn.nii',
+        '.pdconn.nii',
+        '.plabel.nii',
+        '.pscalar.nii',
+        '.ptseries.nii',
+        '.sdseries.nii',
+        '.label.gii',
+        '.label.gii',
+        '.func.gii',
+        '.shape.gii',
+        '.rgba.gii',
+        '.surf.gii',
+        '.dpconn.nii',
+        '.dtraj.nii',
+        '.pconnseries.nii',
+        '.pconnscalar.nii',
+        '.dfan.nii',
+        '.dfibersamp.nii',
+        '.dfansamp.nii',
+    ]
+
+    pth = op.dirname(fname)
+    fname = op.basename(fname)
+
+    ext = None
+    for special_ext in special_extensions:
+        ext_len = len(special_ext)
+        if (len(fname) > ext_len) and (fname[-ext_len:].lower() == special_ext.lower()):
+            ext = fname[-ext_len:]
+            fname = fname[:-ext_len]
+            break
+    if not ext:
+        fname, ext = op.splitext(fname)
+
+    return pth, fname, ext

src/smripost_linc/workflows/base.py

@@ -405,6 +405,8 @@ def init_single_run_wf(anat_file, atlases):
         ]),
     ])  # fmt:skip
 
+    # Calculate myelin map if both T1w and T2w are available
+
     # Fill-in datasinks seen so far
     for node in workflow.list_node_names():
         node_name = node.split('.')[-1]