Source Code for Module mvpa.atlases.transformation

  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  """Coordinate transformations""" 
 10   
 11  import numpy as N 
 12   
 13  if __debug__: 
 14      from mvpa.base import debug 
 15   
 16 -class TypeProxy: 
 17      """ 
 18      Simple class to convert from and then back to original type 
 19      working with list, tuple, ndarray and having 
 20   
 21      XXX Obsolete functionality ?? 
 22      """ 
 23 -    def __init__(self, value, toType=N.array): 
 24          if   isinstance(value, list): self.__type = list 
 25          elif isinstance(value, tuple): self.__type = tuple 
 26          elif isinstance(value, N.ndarray): self.__type = N.array 
 27          else: 
 28              raise IndexError("Not understood format of coordinates '%s' for the transformation" % `coord`) 
 29   
 30 -    def __call__(self, value):    return self.__type(value) 
 31  #   def __getitem__(self, value): return self.__type(value) 
 32   
 33   
 34 -class TransformationBase: 
 35      """ 
 36      Basic class to describe a transformation. Pretty much an interface 
 37      """ 
 38   
 39 -    def __init__(self, previous=None): 
 40          self.previous = previous 
 41   
 42 -    def __getitem__(self, icoord): 
 43          """ 
 44          Obtain coordinates, apply the transformation and spit out in the same 
 45          format (list, tuple, numpy.array) 
 46          """ 
 47   
 48          # remember original type 
 49          #speed origType = TypeProxy(coord) 
 50   
 51          # just in case it is not an ndarray, and to provide a copy to manipulate with 
 52          coord = N.array(icoord) 
 53   
 54          # apply previous transformation if such defined 
 55          if self.previous: 
 56              # if __debug__: debug('ATL__', "Applying previous transformation on `%s`" % `coord`) 
 57              coord = self.previous[coord] 
 58   
 59          #speed if __debug__: debug('ATL__', "Applying main transformation on `%s`" % `coord`) 
 60          # apply transformation 
 61          coord_out = self.apply(coord) 
 62          #speed if __debug__: debug('ATL__', "Applied and got `%s`" % `coord_out`) 
 63   
 64          #speed return origType(coord_out) 
 65          return coord_out 
 66   
 67 -    def __call__(self, coord): 
 68          return self[coord] 
 69   
 70 -    def apply(self, coord): 
 71          return coord 
 72   
 73   
 74 -class SpaceTransformation(TransformationBase): 
 75      """ 
 76      To perform transformation from Voxel into Real Space. 
 77      Simple one -- would subtract the origin and multiply by voxelSize. 
 78      if toRealSpace is True then on call/getitem converts to RealSpace 
 79      """ 
 80 -    def __init__(self, voxelSize=None, origin=None, toRealSpace=True, 
 81                   *args, **kwargs): 
 82   
 83          TransformationBase.__init__(self, *args, **kwargs) 
 84   
 85          if not voxelSize is None: self.voxelSize = N.asarray(voxelSize) 
 86          else: self.voxelSize = 1 
 87   
 88          if not origin is None: self.origin = N.asarray(origin) 
 89          else: self.origin = 0 
 90   
 91          if toRealSpace: 
 92              self.apply = self.toRealSpace 
 93          else: 
 94              self.apply = self.toVoxelSpace 
 95   
 96 -    def toRealSpace(self, coord): 
 97          #speed if not self.origin is None: 
 98          coord -= self.origin 
 99          #speed if not self.voxelSize is None: 
100          coord *= self.voxelSize 
101          return coord 
102   
103 -    def toVoxelSpace(self, coord): 
104          #speed if not self.voxelSize is None: 
105          coord /= self.voxelSize 
106          #speed if not self.origin is None: 
107          coord += self.origin 
108          return map(lambda x:int(round(x)), coord) 
109   
110   
111 -class Linear(TransformationBase): 
112      """ 
113      Simple linear transformation defined by a matrix 
114      """ 
115 -    def __init__(self, transf=N.eye(4), **kwargs): 
116          transf = N.asarray(transf)  # assure that we have arrays not matrices 
117          prev = kwargs.get('previous', None) 
118          if prev is not None and isinstance(prev, Linear): 
119              if prev.N == transf.shape[0] -1: 
120                  if __debug__: debug('ATL__', "Colliding 2 linear transformations into 1") 
121                  transf = N.dot(transf, prev.M) 
122                  # reassign previous transformation to the current one 
123                  kwargs['previous'] = prev.previous 
124          TransformationBase.__init__(self, **kwargs) 
125          self.M = transf 
126          self.N = self.M.shape[0] - 1 
127   
128 -    def apply(self, coord): 
129          #speed if len(coord) != self.__N: 
130          #speed  raise ValueError("Transformation operates on %dD coordinates" \ 
131          #speed                   % self.__N ) 
132          #speed if __debug__: debug('ATL__', "Applying linear coord transformation + %s" % self.__M) 
133          # Might better come up with a linear transformation 
134          coord_ = N.r_[coord, [1.0]] 
135          result = N.dot(self.M, coord_) 
136          return result[0:-1] 
137   
138   
139 -class MNI2Tal_MatthewBrett(TransformationBase): 
140      """ 
141      Transformation to bring MNI coordinates into MNI space 
142   
143      Apparently it is due to Matthew Brett 
144      http://imaging.mrc-cbu.cam.ac.uk/imaging/MniTalairach 
145      """ 
146   
147 -    def __init__(self, *args, **kwargs): 
148          TransformationBase.__init__(self, *args, **kwargs) 
149          self.__upper = Linear( N.array([ [0.9900, 0, 0, 0 ], 
150                                            [0, 0.9688, 0.0460, 0 ], 
151                                            [0,-0.0485, 0.9189, 0 ], 
152                                            [0, 0, 0, 1.0000] ] ) ) 
153   
154          self.__lower = Linear(N.array( [ [0.9900, 0, 0, 0 ], 
155                                            [0, 0.9688, 0.0420, 0 ], 
156                                            [0,-0.0485, 0.8390, 0 ], 
157                                            [0, 0, 0, 1.0000] ] ) ) 
158   
159 -    def apply(self, coord): 
160          return {True: self.__upper, 
161                  False: self.__lower}[coord[2]>=0][coord] 
162   
163   
164 -def MNI2Tal_MeyerLindenberg98 (*args, **kwargs): 
165      """ 
166      Due to Andreas Meyer-Lindenberg 
167      Taken from 
168      http://imaging.mrc-cbu.cam.ac.uk/imaging/MniTalairach 
169      """ 
170   
171      return Linear( N.array([ 
172          [    0.88,   0,  0,  -0.8], 
173          [    0,   0.97,  0,  -3.32], 
174          [    0,   0.05,  0.88,   -0.44], 
175          [    0.00000,   0.00000,   0.00000,   1.00000] ]), *args, **kwargs ) 
176   
177   
178 -def MNI2Tal_YOHflirt (*args, **kwargs): 
179      """Transformations obtained using flirt from Talairach to Standard 
180   
181      Transformations were obtained by registration of 
182      grey/white matter image from talairach atlas to FSL's standard 
183      volume. Following sequence of commands was used: 
184   
185      fslroi /usr/share/rumba/atlases/data/talairach_atlas.nii.gz talairach_graywhite.nii.gz 3 1 
186      flirt -in talairach_graywhite.nii.gz \ 
187       -ref /usr/apps/fsl.4.1/data/standard/MNI152_T1_1mm_brain.nii.gz \ 
188       -out talairach2mni.nii.gz -omat talairach2mni.mat \ 
189       -searchrx -20 20 -searchry -20 20 -searchrz -20 20 -coarsesearch 10 -finesearch 6 -v 
190      flirt -datatype float -in talairach_graywhite.nii.gz -init talairach2mni.mat \ 
191       -ref /usr/apps/fsl.4.1/data/standard/MNI152_T1_1mm_brain.nii.gz \ 
192       -out talairach2mni_fine1.nii.gz -omat talairach2mni_fine1.mat \ 
193       -searchrx -10 10 -searchry -10 10 -searchrz -10 10 -coarsesearch 5 -finesearch 1 -v 
194      convert_xfm -inverse -omat mni2talairach.mat talairach2mni_fine1.mat 
195      """ 
196      return Linear( 
197          t=N.array([ 
198          [ 1.00448,  -0.00629625,  0.00741359,  0.70565,  ], 
199          [ 0.0130797,  0.978238,  0.0731315,  -3.8354,  ], 
200          [ 0.000248407,  -0.0374777,  0.838311,  18.6202,  ], 
201          [ 0,  0,  0,  1,  ], 
202          ]) 
203                     , *args, **kwargs ) 
204   
205   
206 -def Tal2MNI_YOHflirt (*args, **kwargs): 
207      """See MNI2Tal_YOHflirt doc 
208      """ 
209      return Linear( N.array([ 
210          [    1.00452,    0.00441281,  -0.011011,  -0.943886], 
211          [   -0.0141149,  1.00867,     -0.169177,  14.7016], 
212          [    0.00250222, 0.0920984,    1.18656,  -33.922], 
213          [    0.00000,   0.00000,   0.00000,   1.00000] ]), *args, **kwargs ) 
214   
215   
216   
217 -def MNI2Tal_Lancaster07FSL (*args, **kwargs): 
218      return Linear( N.array([ 
219          [  0.9464, 0.0034, -0.0026, -1.0680], 
220          [ -0.0083, 0.9479, -0.0580, -1.0239], 
221          [  0.0053, 0.0617,  0.9010,  3.1883], 
222          [  0.0000, 0.0000,  0.0000,  1.0000] ]), *args, **kwargs ) 
223   
224   
225 -def Tal2MNI_Lancaster07FSL (*args, **kwargs): 
226      return Linear( N.array([ 
227          [ 1.056585, -0.003972,  0.002793,  1.115461], 
228          [ 0.008834,  1.050528,  0.067651,  0.869379], 
229          [-0.00682 , -0.071916,  1.105229, -3.60472 ], 
230          [ 0.      ,  0.      ,  0.      ,  1.      ]]), *args, **kwargs ) 
231   
232   
233 -def MNI2Tal_Lancaster07pooled (*args, **kwargs): 
234      return Linear( N.array([ 
235          [    0.93570,   0.00290,  -0.00720,  -1.04230], 
236          [   -0.00650,   0.93960,  -0.07260,  -1.39400], 
237          [    0.01030,   0.07520,   0.89670,   3.64750], 
238          [    0.00000,   0.00000,   0.00000,   1.00000] ]), *args, **kwargs ) 
239   
240   
241 -def Tal2MNI_Lancaster07pooled (*args, **kwargs): 
242      return Linear( N.array([ 
243          [  1.06860,  -0.00396,   0.00826,   1.07816], 
244          [  0.00640,   1.05741,   0.08566,   1.16824], 
245          [ -0.01281,  -0.08863,   1.10792,  -4.17805], 
246          [  0.00000,   0.00000,   0.00000,   1.00000] ]), *args, **kwargs ) 
247   
248   
249  if __name__ == '__main__': 
250      #t = Tal2Mni 
251      tl = Tal2MNI_Lancaster07FSL() 
252      tli = MNI2Tal_Lancaster07FSL() 
253      tml = MNI2Tal_MeyerLindenberg98() 
254      #print t[1,3,2] 
255      print tl[(1,3,2)] 
256      print tli[[1,3,2]] 
257      print tml[[1,3,2]] 
258  #   print t[(1,3,2,2)] 
259