[RTM] Gifti surfaces

CHANGES.rst

@@ -1,3 +1,8 @@
+Next release
+============
+
+* [ENH] Generate GrayWhite, Pial, MidThickness and inflated surfaces (#398)
+
 0.3.2 (7th of April 2017)
 =========================
 

docs/links.rst

@@ -6,3 +6,5 @@
 .. _workflows: workflows.html
 .. _FreeSurfer: https://surfer.nmr.mgh.harvard.edu/
 .. _`submillimeter reconstruction`: https://surfer.nmr.mgh.harvard.edu/fswiki/SubmillimeterRecon
+.. _GIFTI: https://www.nitrc.org/projects/gifti/
+.. _`Connectome Workbench`: https://www.humanconnectome.org/software/connectome-workbench.html

docs/workflows.rst

@@ -105,14 +105,22 @@ Reconstructed white and pial surfaces are included in the report.
 
  Surface reconstruction (FreeSurfer)
 
-If T1w voxel sizes are less 1mm in all dimensions (rounding to nearest
-.1mm), `submillimeter reconstruction`_ is used.
+If T1w voxel sizes are less than 1mm in all dimensions (rounding to nearest
+.1mm), `submillimeter reconstruction`_ is used, unless disabled with
+``--no-submm-recon``.
 
 In order to bypass reconstruction in ``fmriprep``, place existing reconstructed
 subjects in ``<output dir>/freesurfer`` prior to the run.
 ``fmriprep`` will perform any missing ``recon-all`` steps, but will not perform
 any steps whose outputs already exist.
 
+``lh.midthickness`` and ``rh.midthickness`` surfaces are created in the subject
+``surf/`` directory, corresponding to the surface half-way between the gray/white
+boundary and the pial surface.
+The ``smoothwm``, ``midthickness``, ``pial`` and ``inflated`` surfaces are also
+converted to GIFTI_ format and adjusted to be compatible with multiple software
+packages, including FreeSurfer and the `Connectome Workbench`_.
+
 
 BOLD preprocessing
 ------------------
@@ -301,6 +309,10 @@ Derivatives related to t1w files are in the ``anat`` subfolder:
 - ``*T1w_space-MNI152NLin2009cAsym_brainmask.nii.gz`` Same as above, but in MNI space.
 - ``*T1w_dtissue.nii.gz`` Tissue class map derived using FAST.
 - ``*T1w_preproc.nii.gz`` Bias field corrected t1w file, using ANTS' N4BiasFieldCorrection
+- ``*T1w_smoothwm.[LR].surf.gii`` Smoothed GrayWhite surfaces
+- ``*T1w_pial.[LR].surf.gii`` Pial surfaces
+- ``*T1w_midthickness.[LR].surf.gii`` MidThickness surfaces
+- ``*T1w_inflated.[LR].surf.gii`` FreeSurfer inflated surfaces for visualization
 - ``*T1w_space-MNI152NLin2009cAsym_preproc.nii.gz`` Same as above, but in MNI space
 - ``*T1w_space-MNI152NLin2009cAsym_class-CSF_probtissue.nii.gz``
 - ``*T1w_space-MNI152NLin2009cAsym_class-GM_probtissue.nii.gz``

fmriprep/workflows/anatomical.py

@@ -13,6 +13,7 @@
 from nipype.interfaces import ants
 from nipype.interfaces import freesurfer
 from nipype.interfaces import utility as niu
+from nipype.interfaces import io as nio
 from nipype.pipeline import engine as pe
 
 from niworkflows.interfaces.registration import RobustMNINormalizationRPT
@@ -32,7 +33,9 @@ def t1w_preprocessing(settings, name='t1w_preprocessing'):
 
  workflow = pe.Workflow(name=name)
 
- inputnode = pe.Node(niu.IdentityInterface(fields=['t1w', 't2w', 'subjects_dir']), name='inputnode')
+ inputnode = pe.Node(
+ niu.IdentityInterface(fields=['t1w', 't2w', 'subjects_dir']),
+ name='inputnode')
  outputnode = pe.Node(niu.IdentityInterface(
  fields=['t1_seg', 't1_tpms', 'bias_corrected_t1', 't1_brain', 't1_mask',
  't1_2_mni', 't1_2_mni_forward_transform',
@@ -80,7 +83,6 @@ def t1w_preprocessing(settings, name='t1w_preprocessing'):
  nthreads = settings['nthreads']
 
  def detect_inputs(t1w_list, t2w_list=[], hires_enabled=True):
- import os
  from nipype.interfaces.base import isdefined
  from nipype.utils.filemanip import filename_to_list
  from nipype.interfaces.traits_extension import Undefined
@@ -135,7 +137,7 @@ def bidsinfo(in_file):
  flags='-noskullstrip',
  openmp=nthreads,
  parallel=True),
- name='Reconstruction')
+ name='AutoRecon1')
  autorecon1.interface._can_resume = False
  autorecon1.interface.num_threads = nthreads
 
@@ -176,12 +178,12 @@ def inject_skullstripped(subjects_dir, subject_id, skullstripped):
  parallel=True,
  out_report='reconall.svg',
  generate_report=True),
- name='Reconstruction2')
+ name='ReconAll')
  reconall.interface.num_threads = nthreads
 
  fs_transform = pe.Node(
- freesurfer.utils.Tkregister2(fsl_out='freesurfer2subT1.mat',
-  reg_header=True),
+ freesurfer.Tkregister2(fsl_out='freesurfer2subT1.mat',
+ reg_header=True),
  name='FreeSurferTransform')
 
  recon_report = pe.Node(
@@ -190,6 +192,95 @@ def inject_skullstripped(subjects_dir, subject_id, skullstripped):
  name='ReconAll_Report'
  )
 
+ midthickness = pe.MapNode(
+ freesurfer.MRIsExpand(thickness=True, distance=0.5,
+ out_name='midthickness'),
+ iterfield='in_file',
+ name='MidThickness')
+
+ save_midthickness = pe.Node(nio.DataSink(parameterization=False),
+ name='SaveMidthickness')
+ surface_list = pe.Node(niu.Merge(4), name='SurfaceList')
+ gifticonv = pe.MapNode(freesurfer.MRIsConvert(out_datatype='gii'),
+ iterfield='in_file', name='GiftiSurfaces')
+
+ def get_gifti_name(in_file):
+ import os
+ import re
+ in_format = re.compile(r'(?P<LR>[lr])h.(?P<surf>.+)_converted.gii')
+ name = os.path.basename(in_file)
+ info = in_format.match(name).groupdict()
+ info['LR'] = info['LR'].upper()
+ return '{surf}.{LR}.surf'.format(**info)
+
+ name_surfs = pe.MapNode(
+ niu.Function(
+ function=get_gifti_name,
+ input_names=['in_file'],
+ output_names=['normalized']),
+ iterfield='in_file',
+ name='NameSurfs'
+ )
+
+ def normalize_surfs(in_file):
+ """ Re-center GIFTI coordinates to fit align to native T1 space
+
+ For midthickness surfaces, add MidThickness metadata
+
+ Coordinate update based on:
+ https://github.com/Washington-University/workbench/blob/1b79e56/src/Algorithms/AlgorithmSurfaceApplyAffine.cxx#L73-L91
+ and
+ https://github.com/Washington-University/Pipelines/blob/ae69b9a/PostFreeSurfer/scripts/FreeSurfer2CaretConvertAndRegisterNonlinear.sh#L147
+ """
+ import os
+ import numpy as np
+ import nibabel as nib
+ img = nib.load(in_file)
+ pointset = img.get_arrays_from_intent('NIFTI_INTENT_POINTSET')[0]
+ coords = pointset.data
+ c_ras_keys = ('VolGeomC_R', 'VolGeomC_A', 'VolGeomC_S')
+ ras = np.array([float(pointset.metadata[key])
+ for key in c_ras_keys])
+ # Apply C_RAS translation to coordinates
+ pointset.data = (coords + ras).astype(coords.dtype)
+
+ secondary = nib.gifti.GiftiNVPairs('AnatomicalStructureSecondary',
+ 'MidThickness')
+ geom_type = nib.gifti.GiftiNVPairs('GeometricType', 'Anatomical')
+ has_ass = has_geo = False
+ for nvpair in pointset.meta.data:
+ # Remove C_RAS translation from metadata to avoid double-dipping in FreeSurfer
+ if nvpair.name in c_ras_keys:
+ nvpair.value = '0.000000'
+ # Check for missing metadata
+ elif nvpair.name == secondary.name:
+ has_ass = True
+ elif nvpair.name == geom_type.name:
+ has_geo = True
+ fname = os.path.basename(in_file)
+ # Update metadata for MidThickness/graymid surfaces
+ if 'midthickness' in fname.lower() or 'graymid' in fname.lower():
+ if not has_ass:
+ pointset.meta.data.insert(1, secondary)
+ if not has_geo:
+ pointset.meta.data.insert(2, geom_type)
+ img.to_filename(fname)
+ return os.path.abspath(fname)
+
+ fix_surfs = pe.MapNode(
+ niu.Function(
+ function=normalize_surfs,
+ input_names=['in_file'],
+ output_names=['out_file']),
+ iterfield='in_file',
+ name='FixSurfs')
+
+ ds_surfs = pe.MapNode(
+ DerivativesDataSink(base_directory=settings['output_dir']),
+ iterfield=['in_file', 'suffix'],
+ name='DerivSurfs'
+ )
+
  # Resample the brain mask and the tissue probability maps into mni space
  bmask_mni = pe.Node(
  ants.ApplyTransforms(dimension=3, default_value=0, float=True,
@@ -283,6 +374,20 @@ def inject_skullstripped(subjects_dir, subject_id, skullstripped):
  (reconall, recon_report, [('out_report', 'in_file')]),
  (reconall, outputnode, [('subject_id', 'subject_id')]),
  (fs_transform, outputnode, [('fsl_file', 'fs_2_t1_transform')]),
+ (reconall, midthickness, [('smoothwm', 'in_file')]),
+ (reconall, save_midthickness, [('subjects_dir', 'base_directory'),
+ ('subject_id', 'container')]),
+ (midthickness, save_midthickness, [('out_file', 'surf.@graymid')]),
+ (reconall, surface_list, [('smoothwm', 'in1'),
+ ('pial', 'in2'),
+ ('inflated', 'in3')]),
+ (save_midthickness, surface_list, [('out_file', 'in4')]),
+ (surface_list, gifticonv, [('out', 'in_file')]),
+ (gifticonv, name_surfs, [('converted', 'in_file')]),
+ (gifticonv, fix_surfs, [('converted', 'in_file')]),
+ (inputnode, ds_surfs, [(('t1w', fix_multi_T1w_source_name), 'source_file')]),
+ (name_surfs, ds_surfs, [('normalized', 'suffix')]),
+ (fix_surfs, ds_surfs, [('out_file', 'in_file')]),
  ])
 
  # Write corrected file in the designated output dir

requirements.txt

@@ -1,2 +1,2 @@
-git+https://github.com/nipy/nipype.git@e904aec996a50b4c95f3d1c4ddc81c37aae59583#egg=nipype
+git+https://github.com/nipy/nipype.git@c15df990c1e0849aa789bc15aae286b0a6c16b6c#egg=nipype
 git+https://github.com/poldracklab/niworkflows.git@b0b3260d796826fd12d2b8c75e0c3dc78c83e701#egg=niworkflows