Package mvpa :: Package atlases :: Module fsl

[hide private]
[frames] | no frames]

# Source Code for Module mvpa.atlases.fsl

      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  """FSL atlases interfaces
     10
     11  """
     12
     13  from mvpa.base import warning, externals
     14
     15  if externals.exists('nifti', raiseException=True):
     16      from nifti import NiftiImage
     17
     18  import os, re
     19  import numpy as N
     20
     21  from mvpa.misc.support import reuseAbsolutePath
     22  from mvpa.base.dochelpers import enhancedDocString
     23
     24  from mvpa.atlases.base import XMLBasedAtlas, LabelsLevel
     25
     26  if __debug__:
     27          from mvpa.base import debug
     28
     29  #
     30  # Atlases from FSL
     31  #
     32
     33 -class FSLAtlas(XMLBasedAtlas):
     34      """Base class for FSL atlases
     35
     36      """
     37      source = 'FSL'
     38
     39
     40 -    def __init__(self, *args, **kwargs):
     41          """
     42          """
     43          XMLBasedAtlas.__init__(self, *args, **kwargs)
     44          self.space = 'MNI'
     45
     46
     47      __doc__ = enhancedDocString('FSLAtlas', locals(), XMLBasedAtlas)
     48
     49
     50 -    def _loadImages(self):
     51          resolution = self._resolution
     52          header = self.header
     53          images = header.images
     54          # Load present images
     55          # XXX might be refactored to avoid duplication of
     56          #     effort with PyMVPAAtlas
     57          ni_image = None
     58          resolutions = []
     59          if self._force_image_file is None:
     60              imagefile_candidates = [
     61                  reuseAbsolutePath(self._filename, i.imagefile.text, force=True)
     62                  for i in images]
     63          else:
     64              imagefile_candidates = [self._force_image_file]
     65
     66          for imagefilename in imagefile_candidates:
     67              try:
     68                  ni_image_  = NiftiImage(imagefilename, load=False)
     69              except RuntimeError, e:
     70                  raise RuntimeError, " Cannot open file " + imagefilename
     71
     72              resolution_ = ni_image_.pixdim[0]
     73              if resolution is None:
     74                  # select this one if the best
     75                  if ni_image is None or \
     76                         resolution_ < ni_image.pixdim[0]:
     77                      ni_image = ni_image_
     78                      self._image_file = imagefilename
     79              else:
     80                  if resolution_ == resolution:
     81                      ni_image = ni_image_
     82                      self._image_file = imagefilename
     83                      break
     84                  else:
     85                      resolutions += [resolution_]
     86              # TODO: also make use of summaryimagefile may be?
     87
     88          if ni_image is None:
     89              msg = "Could not find an appropriate atlas among %d atlases."
     90              if resolution is not None:
     91                  msg += " Atlases had resolutions %s" % \
     92                        (resolutions,)
     93              raise RuntimeError, msg
     94          if __debug__:
     95              debug('ATL__', "Loading atlas data from %s" % self._image_file)
     96          self._image = ni_image
     97          self._resolution = ni_image.pixdim[0]
     98          self._origin = N.abs(ni_image.header['qoffset']) * 1.0  # XXX
     99          self._data   = self._image.data
    100
    101
    102 -    def _loadData(self):
    103          """   """
    104          # Load levels
    105          self._levels_dict = {}
    106          # preprocess labels for different levels
    107          self.Nlevels = 1
    108          #level = Level.generateFromXML(self.data, levelType='label')
    109          level = LabelsLevel.generateFromXML(self.data)#, levelType='label')
    110          level.description = self.header.name.text
    111          self._levels_dict = {0: level}
    112          #for index, child in enumerate(self.data.getchildren()):
    113          #   if child.tag == 'level':
    114          #       level = Level.generateFromXML(child)
    115          #       self._levels_dict[level.description] = level
    116          #       try:
    117          #           self._levels_dict[level.index] = level
    118          #       except:
    119          #           pass
    120          #   else:
    121          #       raise XMLAtlasException("Unknown child '%s' within data" % child.tag)
    122          #   self.Nlevels += 1
    123          #pass
    124
    125
    126      @staticmethod
    127 -    def _checkVersion(version):
    128          return re.search('^[0-9]+\.[0-9]', version) is not None
    129
    130
    131 -class FSLProbabilisticAtlas(FSLAtlas):
    132      """Probabilistic FSL atlases
    133      """
    134
    135 -    def __init__(self, thr=0.0, strategy='all', sort=True, *args, **kwargs):
    136          """
    137
    138          :Parameters:
    139            thr : float
    140              Value to threshold at
    141            strategy : basestring
    142              Possible values
    143                all - all entries above thr
    144                max - entry with maximal value
    145            sort : bool
    146              Either to sort entries for 'all' strategy according to
    147              probability
    148          """
    149
    150          FSLAtlas.__init__(self, *args, **kwargs)
    151          self.thr = thr
    152          self.strategy = strategy
    153          self.sort = sort
    154
    155      __doc__ = enhancedDocString('FSLProbabilisticAtlas', locals(), FSLAtlas)
    156
    157 -    def labelVoxel(self, c, levels=None):
    158          """Return labels for the voxel
    159
    160          :Parameters:
    161            - c : tuple of coordinates (xyz)
    162            - levels : just for API consistency (heh heh). Must be 0 for FSL atlases
    163          """
    164
    165          if levels is not None and not (levels in [0, [0], (0,)]):
    166              raise ValueError, \
    167                    "I guess we don't support levels other than 0 in FSL atlas"
    168
    169          # check range
    170          c = self._checkRange(c)
    171
    172          # XXX think -- may be we should better assign each map to a
    173          # different level
    174          level = 0
    175          resultLabels = []
    176          for index, area in enumerate(self._levels_dict[level]):
    177              prob =  int(self._data[index, c[2], c[1], c[0]])
    178              if prob > self.thr:
    179                  resultLabels += [dict(index=index,
    180                                        #id=
    181                                        label=area.text,
    182                                        prob=prob)]
    183
    184          if self.sort or self.strategy == 'max':
    185              resultLabels.sort(cmp=lambda x,y: cmp(x['prob'], y['prob']),
    186                                reverse=True)
    187
    188          if self.strategy == 'max':
    189              resultLabels = resultLabels[:1]
    190          elif self.strategy == 'all':
    191              pass
    192          else:
    193              raise ValueError, 'Unknown strategy %s' % self.strategy
    194
    195          result = {'voxel_queried' : c,
    196                    # in the list since we have only single level but
    197                    # with multiple entries
    198                    'labels': [resultLabels]}
    199
    200          return result
    201
    202 -    def find(self, *args, **kwargs):
    203          """Just a shortcut to the only level.
    204
    205          See :class:`~mvpa.atlases.base.Level.find` for more info
    206          """
    207          return self.levels_dict[0].find(*args, **kwargs)
    208
    209 -    def getMap(self, target, strategy='unique'):
    210          """Return a probability map
    211
    212          :Parameters:
    213            target : int or str or re._pattern_type
    214              If int, map for given index is returned. Otherwise, .find is called
    215              with unique=True to find matching area
    216            strategy : str in ('unique', 'max')
    217              If 'unique', then if multiple areas match, exception would be raised.
    218              In case of 'max', each voxel would get maximal value of probabilities
    219              from all matching areas
    220          """
    221          if isinstance(target, int):
    222              return self._data[target]
    223          else:
    224              lev = self.levels_dict[0]       # we have just 1 here
    225              if strategy == 'unique':
    226                  return self.getMap(lev.find(target, unique=True).index)
    227              else:
    228                  maps = N.array(self.getMaps(target))
    229                  return N.max(maps, axis=0)
    230
    231 -    def getMaps(self, target):
    232          """Return a list of probability maps for the target
    233
    234          :Parameters:
    235            target : str or re._pattern_type
    236              .find is called with a target and unique=False to find all matches
    237          """
    238          lev = self.levels_dict[0]       # we have just 1 here
    239          return [self.getMap(l.index) for l in lev.find(target, unique=False)]
    240
    241
    242 -class FSLLabelsAtlas(XMLBasedAtlas):
    243      """Not sure what this one was for"""
    244 -    def __init__(self, *args, **kwargs):
    245          """not implemented"""
    246          FSLAtlas.__init__(self, *args, **kwargs)
    247          raise NotImplementedError
    248