Merge pull request #646 from satra/enh/restingwf

Example: Resting preprocessing workflow

Author: satra

examples/rsfmri_preprocessing.py

@@ -342,9 +342,9 @@ def __init__(self, **kwargs):
 class BinarizeInputSpec(FSTraitedSpec):
     in_file = File(exists=True, argstr='--i %s', mandatory=True,
                   copyfile=False, desc='input volume')
-    min = traits.Float(argstr='--min %f',
+    min = traits.Float(argstr='--min %f', xor=['wm_ven_csf'],
                        desc='min thresh')
-    max = traits.Float(argstr='--max %f',
+    max = traits.Float(argstr='--max %f', xor=['wm_ven_csf'],
                        desc='max thresh')
     rmin = traits.Float(argstr='--rmin %f',
                         desc='compute min based on rmin*globalmean')
@@ -356,7 +356,7 @@ class BinarizeInputSpec(FSTraitedSpec):
          desc='set match vals to 2 and 41 (aseg for cerebral WM)')
     ventricles = traits.Bool(argstr='--ventricles',
          desc='set match vals those for aseg ventricles+choroid (not 4th)')
-    wm_ven_csf = traits.Bool(argstr='--wm+vcsf',
+    wm_ven_csf = traits.Bool(argstr='--wm+vcsf', xor=['min', 'max'],
           desc='WM and ventricular CSF, including choroid (not 4th)')
     binary_file = File(argstr='--o %s', genfile=True,
                   desc='binary output volume')
@@ -430,7 +430,7 @@ def _list_outputs(self):
                 outfile = fname_presuffix(self.inputs.in_file,
                                           newpath=os.getcwd(),
                                           suffix='_thresh')
-        outputs['binary_file'] = outfile
+        outputs['binary_file'] = os.path.abspath(outfile)
         value = self.inputs.count_file
         if isdefined(value):
             if isinstance(value, bool):

nipype/interfaces/freesurfer/model.py

@@ -1804,7 +1804,6 @@ def _list_outputs(self):
         return outputs
 
     def _gen_filename(self, name):
-        if name in ('out_file'):
-            return fname_presuffix(self.inputs.in_file,
-                                   suffix='_glm.txt', use_ext=False)
+        if name in ['out_file']:
+            return self._gen_fname(self.inputs.in_file, suffix='_glm')
         return None

nipype/interfaces/fsl/model.py

@@ -214,7 +214,7 @@ class DataSink(IOBase):
     input_spec = DataSinkInputSpec
     output_spec = DataSinkOutputSpec
 
-    def __init__(self, infields=None, **kwargs):
+    def __init__(self, infields=None, force_run=True, **kwargs):
         """
         Parameters
         ----------
@@ -232,6 +232,8 @@ def __init__(self, infields=None, **kwargs):
                 self.inputs._outputs[key] = Undefined
                 undefined_traits[key] = Undefined
         self.inputs.trait_set(trait_change_notify=False, **undefined_traits)
+        if force_run:
+            self._always_run = True
 
     def _get_dst(self, src):
         ## If path is directory with trailing os.path.sep,

nipype/interfaces/io.py

@@ -1,3 +1,3 @@
 from .model import FitGLM, EstimateContrast
-from .preprocess import ComputeMask, FmriRealign4d
+from .preprocess import ComputeMask, FmriRealign4d, SpaceTimeRealigner
 from .utils import Similarity

nipype/interfaces/nipy/__init__.py

@@ -24,6 +24,8 @@
     from nipy.labs.mask import compute_mask
     from nipy.algorithms.registration import FmriRealign4d as FR4d
     from nipy import save_image, load_image
+    from nipy.algorithms.registration import SpaceTimeRealign
+    from nipy.algorithms.registration.groupwise_registration import SpaceRealign
 
 from ..base import (TraitedSpec, BaseInterface, traits,
                     BaseInterfaceInputSpec, isdefined, File,
@@ -82,7 +84,7 @@ class FmriRealign4dInputSpec(BaseInterfaceInputSpec):
                              desc="File to realign")
     tr = traits.Float(desc="TR in seconds",
                       mandatory=True)
-    slice_order = traits.List(traits.Int(), maxver=0.3,
+    slice_order = traits.List(traits.Int(),
             desc=('0 based slice order. This would be equivalent to entering'
                   'np.argsort(spm_slice_order) for this field. This effects'
                   'interleaved acquisition. This field will be deprecated in'
@@ -117,6 +119,7 @@ class FmriRealign4dOutputSpec(TraitedSpec):
                                desc="Motion parameter files")
 
 
+@np.deprecate_with_doc('Please use SpaceTimeRealign instead')
 class FmriRealign4d(BaseInterface):
     """Simultaneous motion and slice timing correction algorithm
 
@@ -196,6 +199,131 @@ def _list_outputs(self):
         return outputs
 
 
+class SpaceTimeRealignerInputSpec(BaseInterfaceInputSpec):
+
+    in_file = InputMultiPath(exists=True,
+                             mandatory=True,
+                             desc="File to realign")
+    tr = traits.Float(desc="TR in seconds", requires=['slice_times'])
+    slice_times = traits.Either(traits.List(traits.Float()),
+                                traits.Enum('asc_alt_2', 'asc_alt_2_1',
+                                            'asc_alt_half', 'asc_alt_siemens',
+                                            'ascending', 'desc_alt_2',
+                                            'desc_alt_half', 'descending'),
+                                desc=('Actual slice acquisition times.'))
+    slice_info = traits.Either(traits.Int,
+                               traits.List(min_len=2, max_len=2),
+                               desc=('Single integer or length 2 sequence '
+                                     'If int, the axis in `images` that is the '
+                                     'slice axis.  In a 4D image, this will '
+                                     'often be axis = 2.  If a 2 sequence, then'
+                                     ' elements are ``(slice_axis, '
+                                     'slice_direction)``, where ``slice_axis`` '
+                                     'is the slice axis in the image as above, '
+                                     'and ``slice_direction`` is 1 if the '
+                                     'slices were acquired slice 0 first, slice'
+                                     ' -1 last, or -1 if acquired slice -1 '
+                                     'first, slice 0 last.  If `slice_info` is '
+                                     'an int, assume '
+                                     '``slice_direction`` == 1.'),
+                               requires=['slice_times'],
+                               )
+
+
+class SpaceTimeRealignerOutputSpec(TraitedSpec):
+    out_file = OutputMultiPath(File(exists=True),
+                               desc="Realigned files")
+    par_file = OutputMultiPath(File(exists=True),
+                               desc=("Motion parameter files. Angles are not "
+                                     "euler angles"))
+
+
+class SpaceTimeRealigner(BaseInterface):
+    """Simultaneous motion and slice timing correction algorithm
+
+    If slice_times is not specified, this algorithm performs spatial motion
+    correction
+
+    This interface wraps nipy's SpaceTimeRealign algorithm [1]_ or simply the
+    SpatialRealign algorithm when timing info is not provided.
+
+    Examples
+    --------
+    >>> from nipype.interfaces.nipy import SpaceTimeRealigner
+    >>> #Run spatial realignment only
+    >>> realigner = SpaceTimeRealigner()
+    >>> realigner.inputs.in_file = ['functional.nii']
+    >>> res = realigner.run() # doctest: +SKIP
+
+    >>> realigner = SpaceTimeRealigner()
+    >>> realigner.inputs.in_file = ['functional.nii']
+    >>> realigner.inputs.tr = 2
+    >>> realigner.inputs.slice_times = range(0, 3, 67)
+    >>> realigner.inputs.slice_info = 2
+    >>> res = realigner.run() # doctest: +SKIP
+
+
+    References
+    ----------
+    .. [1] Roche A. A four-dimensional registration algorithm with \
+       application to joint correction of motion and slice timing \
+       in fMRI. IEEE Trans Med Imaging. 2011 Aug;30(8):1546-54. DOI_.
+
+    .. _DOI: http://dx.doi.org/10.1109/TMI.2011.2131152
+
+    """
+
+    input_spec = SpaceTimeRealignerInputSpec
+    output_spec = SpaceTimeRealignerOutputSpec
+    keywords = ['slice timing', 'motion correction']
+
+    def _run_interface(self, runtime):
+
+        all_ims = [load_image(fname) for fname in self.inputs.in_file]
+
+        if not isdefined(self.inputs.slice_times):
+            R = SpaceRealign(all_ims)
+        else:
+            R = SpaceTimeRealign(all_ims,
+                                 tr=self.inputs.tr,
+                                 slice_times=self.inputs.slice_times,
+                                 slice_info=self.inputs.slice_info,
+                                 )
+
+        R.estimate(refscan=None)
+
+        corr_run = R.resample()
+        self._out_file_path = []
+        self._par_file_path = []
+
+        for j, corr in enumerate(corr_run):
+            self._out_file_path.append(os.path.abspath('corr_%s.nii.gz' %
+                                      (split_filename(self.inputs.in_file[j])[1])))
+            save_image(corr, self._out_file_path[j])
+
+            self._par_file_path.append(os.path.abspath('%s.par' %
+                                      (os.path.split(self.inputs.in_file[j])[1])))
+            mfile = open(self._par_file_path[j], 'w')
+            motion = R._transforms[j]
+            # nipy does not encode euler angles. return in original form of
+            # translation followed by rotation vector see:
+            # http://en.wikipedia.org/wiki/Rodrigues'_rotation_formula
+            for i, mo in enumerate(motion):
+                params = ['%.10f' % item for item in np.hstack((mo.translation,
+                                                                mo.rotation))]
+                string = ' '.join(params) + '\n'
+                mfile.write(string)
+            mfile.close()
+
+        return runtime
+
+    def _list_outputs(self):
+        outputs = self._outputs().get()
+        outputs['out_file'] = self._out_file_path
+        outputs['par_file'] = self._par_file_path
+        return outputs
+
+
 class TrimInputSpec(BaseInterfaceInputSpec):
     in_file = File(
         exists=True, mandatory=True,

nipype/interfaces/nipy/preprocess.py

@@ -1211,6 +1211,7 @@ def hash_exists(self, updatehash=False):
         outdir = self.output_dir()
         hashfiles = glob(os.path.join(outdir, '_0x*.json'))
         if len(hashfiles) > 1:
+            logger.info(hashfiles)
             logger.info('Removing multiple hashfiles and forcing node to rerun')
             for hashfile in hashfiles:
                 os.unlink(hashfile)

nipype/pipeline/engine.py

@@ -75,8 +75,14 @@ def _get_result(self, taskid):
 
     def _submit_job(self, node, updatehash=False):
         self._taskid += 1
+        try:
+            if node.inputs.terminal_output == 'stream':
+                node.inputs.terminal_output = 'file'
+        except:
+            pass
         self._taskresult[self._taskid] = self.pool.apply_async(run_node,
-                                                               (node, updatehash,))
+                                                               (node,
+                                                                updatehash,))
         return self._taskid
 
     def _report_crash(self, node, result=None):