1
2
3
4
5
6
7
8
9 """Dataset that gets its samples from a NIfTI file"""
10
11 __docformat__ = 'restructuredtext'
12
13 from mvpa.base import externals
14
15 import sys
16 import numpy as N
17 from mvpa.support.copy import deepcopy
18
19 if __debug__:
20 from mvpa.base import debug
21
22 if externals.exists('nifti', raiseException=True):
23 if sys.version_info[:2] >= (2, 5):
24
25 NiftiImage = __import__('nifti', globals(), locals(), [], 0).NiftiImage
26 else:
27
28
29 oldname = __name__
30
31 __name__ = 'iaugf9zrkjsbdv89'
32 from nifti import NiftiImage
33
34 __name__ = oldname
35
36 from mvpa.datasets.base import Dataset
37 from mvpa.datasets.mapped import MappedDataset
38 from mvpa.datasets.event import EventDataset
39 from mvpa.mappers.base import CombinedMapper
40 from mvpa.mappers.metric import DescreteMetric, cartesianDistance
41 from mvpa.mappers.array import DenseArrayMapper
42 from mvpa.base import warning
43
44
46 """Load/access NIfTI data from files or instances.
47
48 :Parameters:
49 src: str | NiftiImage
50 Filename of a NIfTI image or a `NiftiImage` instance.
51 ensure : bool
52 If True, through ValueError exception if cannot be loaded.
53 enforce_dim : int or None
54 If not None, it is the dimensionality of the data to be enforced,
55 commonly 4D for the data, and 3D for the mask in case of fMRI.
56
57 :Returns:
58 NiftiImage | None
59 If the source is not supported None is returned.
60 """
61 nifti = None
62
63
64 if isinstance(src, str):
65
66 try:
67 nifti = NiftiImage(src)
68 except RuntimeError, e:
69 warning("ERROR: NiftiDatasets: Cannot open NIfTI file %s" \
70 % src)
71 raise e
72 elif isinstance(src, NiftiImage):
73
74 nifti = src
75 elif (isinstance(src, list) or isinstance(src, tuple)) \
76 and len(src)>0 \
77 and (isinstance(src[0], str) or isinstance(src[0], NiftiImage)):
78
79 if enforce_dim is not None: enforce_dim_ = enforce_dim - 1
80 else: enforce_dim_ = None
81 srcs = [getNiftiFromAnySource(s, ensure=ensure,
82 enforce_dim=enforce_dim_)
83 for s in src]
84 if __debug__:
85
86 shapes = [s.data.shape for s in srcs]
87 if not N.all([s == shapes[0] for s in shapes]):
88 raise ValueError, \
89 "Input volumes contain variable number of dimensions:" \
90 " %s" % (shapes,)
91
92 nifti = NiftiImage(N.array([s.asarray() for s in srcs]),
93 srcs[0].header)
94 elif ensure:
95 raise ValueError, "Cannot load NIfTI from %s" % (src,)
96
97 if nifti is not None and enforce_dim is not None:
98 shape, new_shape = nifti.data.shape, None
99 lshape = len(shape)
100
101
102 if lshape < enforce_dim:
103
104 new_shape = (1,)*(enforce_dim-lshape) + shape
105 elif lshape > enforce_dim:
106
107 bogus_dims = lshape - enforce_dim
108 if shape[:bogus_dims] != (1,)*bogus_dims:
109 raise ValueError, \
110 "Cannot enforce %dD on data with shape %s" \
111 % (enforce_dim, shape)
112 new_shape = shape[bogus_dims:]
113
114
115 if new_shape is not None:
116 if __debug__:
117 debug('DS_NIFTI', 'Enforcing shape %s for %s data from %s' %
118 (new_shape, shape, src))
119 nifti.data.shape = new_shape
120
121 return nifti
122
123
125 """Convenience function to extract the data array from a NiftiImage
126
127 This function will make use of advanced features of PyNIfTI to prevent
128 unnecessary copying if a sufficent version is available.
129 """
130 if externals.exists('nifti ge 0.20090205.1'):
131 return nim.data
132 else:
133 return nim.asarray()
134
135
137 """Dataset loading its samples from a NIfTI image or file.
138
139 Samples can be loaded from a NiftiImage instance or directly from a NIfTI
140 file. This class stores all relevant information from the NIfTI file header
141 and provides information about the metrics and neighborhood information of
142 all voxels.
143
144 Most importantly it allows to map data back into the original data space
145 and format via :meth:`~mvpa.datasets.nifti.NiftiDataset.map2Nifti`.
146
147 This class allows for convenient pre-selection of features by providing a
148 mask to the constructor. Only non-zero elements from this mask will be
149 considered as features.
150
151 NIfTI files are accessed via PyNIfTI. See
152 http://niftilib.sourceforge.net/pynifti/ for more information about
153 pynifti.
154 """
155
156
157 - def __init__(self, samples=None, mask=None, dsattr=None,
158 enforce_dim=4, **kwargs):
159 """
160 :Parameters:
161 samples: str | NiftiImage
162 Filename of a NIfTI image or a `NiftiImage` instance.
163 mask: str | NiftiImage | ndarray
164 Filename of a NIfTI image or a `NiftiImage` instance or an ndarray
165 of appropriate shape.
166 enforce_dim : int or None
167 If not None, it is the dimensionality of the data to be enforced,
168 commonly 4D for the data, and 3D for the mask in case of fMRI.
169 """
170
171 if not dsattr is None and dsattr.has_key('mapper'):
172 MappedDataset.__init__(self,
173 samples=samples,
174 dsattr=dsattr,
175 **kwargs)
176 return
177
178
179
180
181
182
183
184 niftisamples = getNiftiFromAnySource(samples, ensure=True,
185 enforce_dim=enforce_dim)
186 samples = niftisamples.data
187
188
189
190
191
192
193 dsattr = {'niftihdr': niftisamples.header}
194
195
196
197
198 niftimask = getNiftiFromAnySource(mask)
199 if niftimask is None:
200 pass
201 elif isinstance(niftimask, N.ndarray):
202 mask = niftimask
203 else:
204 mask = getNiftiData(niftimask)
205
206
207
208
209
210
211
212 elementsize = [i for i in reversed(niftisamples.voxdim)]
213 mapper = DenseArrayMapper(mask=mask, shape=samples.shape[1:],
214 metric=DescreteMetric(elementsize=elementsize,
215 distance_function=cartesianDistance))
216
217 MappedDataset.__init__(self,
218 samples=samples,
219 mapper=mapper,
220 dsattr=dsattr,
221 **kwargs)
222
223
225 """Maps a data vector into the dataspace and wraps it with a
226 NiftiImage. The header data of this object is used to initialize
227 the new NiftiImage.
228
229 :Parameters:
230 data : ndarray or Dataset
231 The data to be wrapped into NiftiImage. If None (default), it
232 would wrap samples of the current dataset. If it is a Dataset
233 instance -- takes its samples for mapping
234 """
235 if data is None:
236 data = self.samples
237 elif isinstance(data, Dataset):
238
239 data = data.samples
240 dsarray = self.mapper.reverse(data)
241 return NiftiImage(dsarray, self.niftihdr)
242
243
245 """Return the temporal distance of two samples/volumes.
246
247 This method tries to be clever and always returns `dt` in seconds, by
248 using unit information from the NIfTI header. If such information is
249 not present the assumed unit will also be `seconds`.
250 """
251
252 hdr = self.niftihdr
253 TR = hdr['pixdim'][4]
254
255
256 scale = 1.0
257
258
259 if hdr.has_key('time_unit'):
260 unit_code = hdr['time_unit'] / 8
261 elif hdr.has_key('xyzt_unit'):
262 unit_code = int(hdr['xyzt_unit']) / 8
263 else:
264 warning("No information on time units is available. Assuming "
265 "seconds")
266 unit_code = 0
267
268
269
270
271
272 if unit_code in [0, 1, 2, 3]:
273 if unit_code == 0:
274 warning("Time units were not specified in NiftiImage. "
275 "Assuming seconds.")
276 scale = [ 1.0, 1.0, 1e-3, 1e-6 ][unit_code]
277 else:
278 warning("Time units are incorrectly coded: value %d whenever "
279 "allowed are 8 (sec), 16 (millisec), 24 (microsec). "
280 "Assuming seconds." % (unit_code * 8,)
281 )
282 return TR * scale
283
284
285 niftihdr = property(fget=lambda self: self._dsattr['niftihdr'],
286 doc='Access to the NIfTI header dictionary.')
287
288 dt = property(fget=getDt,
289 doc='Time difference between two samples (in seconds). '
290 'AKA TR in fMRI world.')
291
292 samplingrate = property(fget=lambda self: 1.0 / self.dt,
293 doc='Sampling rate (based on .dt).')
294
295
297 """Dataset with event-defined samples from a NIfTI timeseries image.
298
299 This is a convenience dataset to facilitate the analysis of event-related
300 fMRI datasets. Boxcar-shaped samples are automatically extracted from the
301 full timeseries using :class:`~mvpa.misc.support.Event` definition lists.
302 For each event all volumes covering that particular event in time
303 (including partial coverage) are used to form the corresponding sample.
304
305 The class supports the conversion of events defined in 'realtime' into the
306 descrete temporal space defined by the NIfTI image. Moreover, potentially
307 varying offsets between true event onset and timepoint of the first selected
308 volume can be stored as an additional feature in the dataset.
309
310 Additionally, the dataset supports masking. This is done similar to the
311 masking capabilities of :class:`~mvpa.datasets.nifti.NiftiDataset`. However,
312 the mask can either be of the same shape as a single NIfTI volume, or
313 can be of the same shape as the generated boxcar samples, i.e.
314 a samples consisting of three volumes with 24 slices and 64x64 inplane
315 resolution needs a mask with shape (3, 24, 64, 64). In the former case the
316 mask volume is automatically expanded to be identical in a volumes of the
317 boxcar.
318 """
319 - def __init__(self, samples=None, events=None, mask=None, evconv=False,
320 storeoffset=False, tr=None, enforce_dim=4, **kwargs):
321 """
322 :Paramaters:
323 mask: str | NiftiImage | ndarray
324 Filename of a NIfTI image or a `NiftiImage` instance or an ndarray
325 of appropriate shape.
326 evconv: bool
327 Convert event definitions using `onset` and `duration` in some
328 temporal unit into #sample notation.
329 storeoffset: Bool
330 Whether to store temproal offset information when converting
331 Events into descrete time. Only considered when evconv == True.
332 tr: float
333 Temporal distance of two adjacent NIfTI volumes. This can be used
334 to override the corresponding value in the NIfTI header.
335 enforce_dim : int or None
336 If not None, it is the dimensionality of the data to be enforced,
337 commonly 4D for the data, and 3D for the mask in case of fMRI.
338 """
339
340 if events is None:
341 EventDataset.__init__(self, samples=samples, events=events,
342 mask=mask, **kwargs)
343 return
344
345 nifti = getNiftiFromAnySource(samples, ensure=True,
346 enforce_dim=enforce_dim)
347
348 samples = nifti.data
349
350
351
352
353
354
355 dsattr = {'niftihdr': nifti.header}
356
357
358 dt = nifti.rtime
359
360 if not tr is None:
361 dt = tr
362
363
364
365
366 elementsize = [dt] + [i for i in reversed(nifti.voxdim)]
367
368
369
370 metric = DescreteMetric(elementsize=elementsize,
371 distance_function=cartesianDistance)
372
373
374 if evconv:
375 if dt == 0:
376 raise ValueError, "'dt' cannot be zero when converting Events"
377
378 events = [ev.asDescreteTime(dt, storeoffset) for ev in events]
379 else:
380
381 events = deepcopy(events)
382
383
384
385 for ev in events:
386 oldonset = ev['onset']
387 oldduration = ev['duration']
388 ev['onset'] = int(ev['onset'])
389 ev['duration'] = int(ev['duration'])
390 if not oldonset == ev['onset'] \
391 or not oldduration == ev['duration']:
392 warning("Loosing information during automatic integer "
393 "conversion of EVs. Consider an explicit conversion"
394 " by setting `evconv` in ERNiftiDataset().")
395
396
397 if mask is None:
398 pass
399 elif isinstance(mask, N.ndarray):
400
401 pass
402 else:
403 mask_nim = getNiftiFromAnySource(mask)
404 if not mask_nim is None:
405 mask = getNiftiData(mask_nim)
406 else:
407 raise ValueError, "Cannot load mask from '%s'" % mask
408
409
410 EventDataset.__init__(self, samples=samples, events=events,
411 mask=mask, dametric=metric, dsattr=dsattr,
412 **kwargs)
413
414
416 """Maps a data vector into the dataspace and wraps it with a
417 NiftiImage. The header data of this object is used to initialize
418 the new NiftiImage.
419
420 .. note::
421 Only the features corresponding to voxels are mapped back -- not
422 any additional features passed via the Event definitions.
423
424 :Parameters:
425 data : ndarray or Dataset
426 The data to be wrapped into NiftiImage. If None (default), it
427 would wrap samples of the current dataset. If it is a Dataset
428 instance -- takes its samples for mapping
429 """
430 if data is None:
431 data = self.samples
432 elif isinstance(data, Dataset):
433
434 data = data.samples
435
436 mr = self.mapper.reverse(data)
437
438
439 if isinstance(self.mapper, CombinedMapper):
440
441 mr = mr[0]
442 else:
443 pass
444
445 return NiftiImage(mr, self.niftihdr)
446
447
448 niftihdr = property(fget=lambda self: self._dsattr['niftihdr'],
449 doc='Access to the NIfTI header dictionary.')
450