fix: add mean back for fsl 5.0.7 after highpass - test for shape in generating dof volumes

Author: satra

examples/fmri_ants_openfmri.py

@@ -471,23 +471,7 @@ def get_contrasts(contrast_file, task_id, conds):
                            name="modelspec")
     modelspec.inputs.input_units = 'secs'
 
-    '''
-    def check_behav_list(behav):
-        from nipype.external import six
-        out_behav = []
-        if isinstance(behav, six.string_types):
-            behav = [behav]
-        for val in behav:
-            if not isinstance(val, list):
-                out_behav.append([val])
-            else:
-                out_behav.append(val)
-        return out_behav
-
-    wf.connect(subjinfo, 'TR', modelspec, 'time_repetition')
-    wf.connect(datasource, ('behav', check_behav_list), modelspec, 'event_files')
-    '''
-    def check_behav_list(behav, conds):
+    def check_behav_list(behav, run_id, conds):
         from nipype.external import six
         import numpy as np
         num_conds = len(conds)
@@ -497,13 +481,14 @@ def check_behav_list(behav, conds):
         num_elements = behav_array.shape[0]
         return behav_array.reshape(num_elements/num_conds, num_conds).tolist()
 
-    reshape_behav = pe.Node(niu.Function(input_names=['behav', 'conds'],
+    reshape_behav = pe.Node(niu.Function(input_names=['behav', 'run_id', 'conds'],
                                        output_names=['behav'],
                                        function=check_behav_list),
                           name='reshape_behav')
 
     wf.connect(subjinfo, 'TR', modelspec, 'time_repetition')
     wf.connect(datasource, 'behav', reshape_behav, 'behav')
+    wf.connect(subjinfo, 'run_id', reshape_behav, 'run_id')
     wf.connect(subjinfo, 'conds', reshape_behav, 'conds')
     wf.connect(reshape_behav, 'behav', modelspec, 'event_files')
 

nipype/workflows/fmri/fsl/estimate.py

@@ -250,7 +250,10 @@ def get_dofvolumes(dof_files, cope_files):
         import nibabel as nb
         import numpy as np
         img = nb.load(cope_files[0])
-        out_data = np.zeros(img.get_shape())
+        if len(img.get_shape()) > 3:
+            out_data = np.zeros(img.get_shape())
+        else:
+            out_data = np.zeros(list(img.get_shape()) + [1])
         for i in range(out_data.shape[-1]):
             dof = np.loadtxt(dof_files[i])
             out_data[:, :, :, i] = dof

nipype/workflows/fmri/fsl/preprocess.py

@@ -8,6 +8,8 @@
 import nipype.interfaces.freesurfer as fs    # freesurfer
 import nipype.interfaces.spm as spm
 
+from nipype import LooseVersion
+
 from ...smri.freesurfer.utils import create_getmask_flow
 
 def getthreshop(thresh):
@@ -415,6 +417,11 @@ def create_featreg_preproc(name='featpreproc', highpass=True, whichvol='middle')
     >>> preproc.run() # doctest: +SKIP
     """
 
+    version = 0
+    if fsl.Info.version() and \
+        LooseVersion(fsl.Info.version()) > LooseVersion('5.0.6'):
+        version = 507
+
     featpreproc = pe.Workflow(name=name)
 
     """
@@ -658,17 +665,6 @@ def create_featreg_preproc(name='featpreproc', highpass=True, whichvol='middle')
 
     featpreproc.connect(medianval, ('out_stat', getmeanscale), meanscale, 'op_string')
 
-    """
-    Perform temporal highpass filtering on the data
-    """
-
-    if highpass:
-        highpass = pe.MapNode(interface=fsl.ImageMaths(suffix='_tempfilt'),
-                              iterfield=['in_file'],
-                              name='highpass')
-        featpreproc.connect(inputnode, ('highpass', highpass_operand), highpass, 'op_string')
-        featpreproc.connect(meanscale, 'out_file', highpass, 'in_file')
-        featpreproc.connect(highpass, 'out_file', outputnode, 'highpassed_files')
 
     """
     Generate a mean functional image from the first run
@@ -681,6 +677,37 @@ def create_featreg_preproc(name='featpreproc', highpass=True, whichvol='middle')
 
     featpreproc.connect(meanscale, ('out_file', pickfirst), meanfunc3, 'in_file')
     featpreproc.connect(meanfunc3, 'out_file', outputnode, 'mean')
+
+    """
+    Perform temporal highpass filtering on the data
+    """
+
+    if highpass:
+        highpass = pe.MapNode(interface=fsl.ImageMaths(suffix='_tempfilt'),
+                              iterfield=['in_file'],
+                              name='highpass')
+        featpreproc.connect(inputnode, ('highpass', highpass_operand), highpass, 'op_string')
+        featpreproc.connect(meanscale, 'out_file', highpass, 'in_file')
+
+        if version < 507:
+            featpreproc.connect(highpass, 'out_file', outputnode, 'highpassed_files')
+        else:
+            """
+            Add back the mean removed by the highpass filter operation as of FSL 5.0.7
+            """
+            meanfunc4 = pe.MapNode(interface=fsl.ImageMaths(op_string='-Tmean',
+                                                            suffix='_mean'),
+                                   iterfield=['in_file'],
+                                   name='meanfunc4')
+            
+            featpreproc.connect(meanscale, 'out_file', meanfunc4, 'in_file')
+            addmean = pe.MapNode(interface=fsl.BinaryMaths(operation='add'),
+                                 iterfield=['in_file', 'operand_file'],
+                                 name='addmean')
+            featpreproc.connect(highpass, 'out_file', addmean, 'in_file')
+            featpreproc.connect(meanfunc4, 'out_file', addmean, 'operand_file')
+            featpreproc.connect(addmean, 'out_file', outputnode, 'highpassed_files')
+
     return featpreproc