Source Code for Module mvpa.misc.fsl.flobs

  1  # emacs: -*- mode: python; py-indent-offset: 4; indent-tabs-mode: nil -*- 
  2  # vi: set ft=python sts=4 ts=4 sw=4 et: 
  3  ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ## 
  4  # 
  5  #   See COPYING file distributed along with the PyMVPA package for the 
  6  #   copyright and license terms. 
  7  # 
  8  ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ## 
  9  """Wrapper around FSLs halfcosbasis to generate HRF kernels""" 
 10   
 11  __docformat__ = 'restructuredtext' 
 12   
 13  import os 
 14  import tempfile 
 15  import shutil 
 16  import numpy as N 
 17  import math 
 18   
 19 -def makeFlobs(pre=0, rise=5, fall=5, undershoot=5, undershootamp=0.3,  
 20                nsamples=1, resolution=0.05, nsecs=-1, nbasisfns=2): 
 21      """Wrapper around the FSL tool halfcosbasis. 
 22   
 23      This function uses halfcosbasis to generate samples of HRF kernels. 
 24      Kernel parameters can be modified analogous to the Make_flobs GUI 
 25      which is part of FSL.  
 26   
 27      :: 
 28   
 29         ^         /-\\ 
 30         |        /   \\ 
 31         1       /     \\ 
 32         |      /       \\ 
 33         |     /         \\ 
 34         |    /           \\ 
 35        -----/             \\     /-----  | 
 36                            \\--/         |  undershootamp 
 37        |    |      |     |        | 
 38        |    |      |     |        | 
 39   
 40         pre   rise  fall  undershoot 
 41   
 42      Parameters 'pre', 'rise', 'fall', 'undershoot' and 'undershootamp' 
 43      can be specified as 2-tuples (min-max range for sampling) and single  
 44      value (setting exact values -- no sampling). 
 45   
 46      If 'nsec' is negative, the length of the samples is determined  
 47      automatically to include the whole kernel function (until it returns 
 48      to baseline). 'nsec' has to be an integer value and is set to the next  
 49      greater integer value if it is not. 
 50   
 51      All parameters except for 'nsamples' and 'nbasisfns' are in seconds. 
 52      """ 
 53      # create tempdir and temporary parameter file 
 54      pfile, pfilename = tempfile.mkstemp('pyflobs') 
 55      wdir = tempfile.mkdtemp('pyflobs') 
 56   
 57      # halfcosbasis can only handle >1 samples 
 58      # so we simply compute two and later ignore the other 
 59      if nsamples < 2: 
 60          rnsamples = 2 
 61      else: 
 62          rnsamples = nsamples 
 63   
 64      # make range tuples if not supplied 
 65      if not isinstance(pre, tuple): 
 66          pre = (pre, pre) 
 67      if not isinstance(rise, tuple): 
 68          rise = (rise, rise) 
 69      if not isinstance(fall, tuple): 
 70          fall = (fall, fall) 
 71      if not isinstance(undershoot, tuple): 
 72          undershoot = (undershoot, undershoot) 
 73      if not isinstance(undershootamp, tuple): 
 74          undershootamp = (undershootamp, undershootamp) 
 75   
 76      # calc minimum length of hrf if not specified 
 77      # looks like it has to be an integer 
 78      if nsecs < 0: 
 79          nsecs = int( math.ceil( pre[1] \ 
 80                                  + rise[1] \ 
 81                                  + fall[1] \ 
 82                                  + undershoot[1] \ 
 83                                  + resolution ) ) 
 84      else: 
 85          nsecs = math.ceil(nsecs) 
 86   
 87      # write parameter file 
 88      pfile = os.fdopen( pfile, 'w' ) 
 89   
 90      pfile.write(str(pre[0]) + ' ' + str(pre[1]) + '\n') 
 91      pfile.write(str(rise[0]) + ' ' + str(rise[1]) + '\n') 
 92      pfile.write(str(fall[0]) + ' ' + str(fall[1]) + '\n') 
 93      pfile.write(str(undershoot[0]) + ' ' + str(undershoot[1]) + '\n') 
 94      pfile.write('0 0\n0 0\n') 
 95      pfile.write(str(undershootamp[0]) + ' ' + str(undershootamp[1]) + '\n') 
 96      pfile.write('0 0\n') 
 97   
 98      pfile.close() 
 99   
100      # call halfcosbasis to generate the hrf samples 
101      tochild, fromchild, childerror = os.popen3('halfcosbasis'  
102                     + ' --hf=' + pfilename 
103                     + ' --nbfs=' + str(nbasisfns) 
104                     + ' --ns=' + str(nsecs) 
105                     + ' --logdir=' + os.path.join(wdir, 'out') 
106                     + ' --nhs=' + str(rnsamples) 
107                     + ' --res=' + str(resolution) ) 
108      err = childerror.readlines() 
109      if len(err) > 0: 
110          print err 
111          raise RuntimeError, "Problem while running halfcosbasis." 
112   
113      # read samples from file into an array 
114      hrfs = N.fromfile( os.path.join( wdir, 'out', 'hrfsamps.txt' ), 
115                         sep = ' ' ) 
116   
117      # reshape array to get one sample per row and 1d array only 
118      # for one sample hrf 
119      hrfs = \ 
120          hrfs.reshape( len(hrfs)/rnsamples, rnsamples).T[:nsamples].squeeze() 
121   
122      # cleanup working dir (ignore errors) 
123      shutil.rmtree( wdir, True ) 
124      # remove paramter file 
125      os.remove( pfilename ) 
126   
127      # and return an array 
128      return( hrfs ) 
129