Package nifti :: Module niftiformat
[hide private]
[frames] | no frames]

Source Code for Module nifti.niftiformat

  1  #emacs: -*- mode: python-mode; py-indent-offset: 4; indent-tabs-mode: nil -*- 
  2  #ex: set sts=4 ts=4 sw=4 et: 
  3  ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ## 
  4  # 
  5  #   See COPYING file distributed along with the PyNIfTI package for the 
  6  #   copyright and license terms. 
  7  # 
  8  ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ## 
  9  """Python class representation of a NIfTI image header""" 
 10   
 11  __docformat__ = 'restructuredtext' 
 12   
 13   
 14  # the swig wrapper if the NIfTI C library 
 15  import nifticlib 
 16  from utils import nhdr2dict, updateNiftiHeaderFromDict, \ 
 17      Ndtype2niftidtype, nifti_xform_map 
 18  import numpy as N 
 19   
 20   
21 -class NiftiFormat(object):
22 """NIfTI header representation. 23 24 NIfTI header can be created by loading information from an existing NIfTI 25 file or by creating a matching NIfTI header for a ndarray. 26 27 In addition, a number of methods to manipulate the header information are 28 provided. However, this class is not able to write a NIfTI header back to 29 disk. Please refer to the NIfTIImage class for this functionality. 30 """
31 - def __init__(self, source, header=None):
32 """ 33 This method decides whether to load a nifti image header from file or 34 create one from ndarray data, depending on the datatype of `source`. 35 36 :Parameters: 37 source: str | ndarray 38 If source is a string, it is assumed to be a filename and an 39 attempt will be made to open the corresponding NIfTI file. 40 In case of an ndarray the array data will be used for the to be 41 created nifti image and a matching nifti header is generated. 42 If an object of a different type is supplied as 'source' a 43 ValueError exception will be thrown. 44 header: dict 45 Additional header data might be supplied. However, 46 dimensionality and datatype are determined from the ndarray and 47 not taken from a header dictionary. 48 """ 49 self.__nimg = None 50 51 if header == None: 52 header = {} 53 54 if type(source) == N.ndarray: 55 self.__newFromArray(source, header) 56 elif type(source) in (str, unicode): 57 self.__newFromFile(source) 58 else: 59 raise ValueError, \ 60 "Unsupported source type. Only NumPy arrays and filename " \ 61 + "string are supported."
62 63
64 - def __del__(self):
65 if self.__nimg: 66 nifticlib.nifti_image_free(self.__nimg)
67 68
69 - def __newFromArray(self, data, hdr = {}):
70 """Create a `nifti_image` struct from a ndarray. 71 72 :Parameters: 73 data: ndarray 74 Source ndarray. 75 hdr: dict 76 Optional dictionary with NIfTI header data. 77 """ 78 79 # check array 80 if len(data.shape) > 7: 81 raise ValueError, \ 82 "NIfTI does not support data with more than 7 dimensions." 83 84 # create template nifti header struct 85 niptr = nifticlib.nifti_simple_init_nim() 86 nhdr = nifticlib.nifti_convert_nim2nhdr(niptr) 87 88 # intermediate cleanup 89 nifticlib.nifti_image_free(niptr) 90 91 # convert virgin nifti header to dict to merge properties 92 # with supplied information and array properties 93 hdic = nhdr2dict(nhdr) 94 95 # copy data from supplied header dict 96 for k, v in hdr.iteritems(): 97 hdic[k] = v 98 99 # finally set header data that is determined by the data array 100 # convert NumPy to nifti datatype 101 hdic['datatype'] = Ndtype2niftidtype(data) 102 103 # make sure there are no zeros in the dim vector 104 # especially not in #4 as FSLView doesn't like that 105 hdic['dim'] = [ 1 for i in hdic['dim'] ] 106 107 # set number of dims 108 hdic['dim'][0] = len(data.shape) 109 110 # set size of each dim (and reverse the order to match nifti format 111 # requirements) 112 for i, s in enumerate(data.shape): 113 hdic['dim'][len(data.shape)-i] = s 114 115 # set magic field to mark as nifti file 116 hdic['magic'] = 'n+1' 117 118 # update nifti header with information from dict 119 updateNiftiHeaderFromDict(nhdr, hdic) 120 121 # convert nifti header to nifti image struct 122 self.__nimg = nifticlib.nifti_convert_nhdr2nim(nhdr, 'pynifti_none') 123 124 if not self.__nimg: 125 raise RuntimeError, "Could not create nifti image structure." 126 127 # kill filename for nifti images from arrays 128 self.__nimg.fname = '' 129 self.__nimg.iname = ''
130 131
132 - def __newFromFile(self, filename):
133 """Open a NIfTI file. 134 135 :Parameters: 136 filename: str 137 Filename of the to be opened image file. 138 """ 139 # do not load image data! 140 self.__nimg = nifticlib.nifti_image_read(filename, 0) 141 142 if not self.__nimg: 143 raise RuntimeError, "Error while opening nifti header."
144 145
146 - def getVoxDims(self):
147 """Returns a 3-tuple a voxel dimensions/size in (x,y,z). 148 149 The `voxdim` property is an alternative way to access this function. 150 """ 151 return (self.__nimg.dx, self.__nimg.dy, self.__nimg.dz)
152 153
154 - def setVoxDims(self, value):
155 """Set voxel dimensions/size. 156 157 The qform matrix and its inverse will be recalculated automatically. 158 159 :Parameter: 160 value: 3-tuple of floats 161 162 Besides reading it is also possible to set the voxel dimensions by 163 assigning to the `voxdim` property. 164 """ 165 if len(value) != 3: 166 raise ValueError, 'Requires 3-tuple.' 167 168 self.__nimg.dx = float(value[0]) 169 self.__nimg.dy = float(value[1]) 170 self.__nimg.dz = float(value[2]) 171 172 self.updateQFormFromQuaternion()
173 174
175 - def setPixDims(self, value):
176 """Set the pixel dimensions. 177 178 :Parameter: 179 value: sequence 180 Up to 7 values (max. number of dimensions supported by the 181 NIfTI format) are allowed in the sequence. 182 183 The supplied sequence can be shorter than seven elements. In 184 this case only present values are assigned starting with the 185 first dimension (spatial: x). 186 187 Calling `setPixDims()` with a length-3 sequence equals calling 188 `setVoxDims()`. 189 """ 190 if len(value) > 7: 191 raise ValueError, \ 192 'The Nifti format does not support more than 7 dimensions.' 193 194 pixdim = nifticlib.floatArray_frompointer( self.__nimg.pixdim ) 195 196 for i, val in enumerate(value): 197 pixdim[i+1] = float(val) 198 199 # The nifticlib uses dimension deltas (dx, dy, dz, dt...) to store 200 # the pixdim values (in addition to the pixdim array). When 201 # saving the image to a file, the deltas are used, not the pixdims. 202 # The nifti_update_dims_from_array sync's the deltas with the pixdims. 203 # (It also syncs the dim array with it's duplicate scalar variables.) 204 nifticlib.nifti_update_dims_from_array(self.__nimg)
205
206 - def getPixDims(self):
207 """Returns the pixel dimensions on all 7 dimensions. 208 209 The function is similar to `getVoxDims()`, but instead of the 3d 210 spatial dimensions of a voxel it returns the dimensions of an image 211 pixel on all 7 dimensions supported by the NIfTI dataformat. 212 """ 213 return \ 214 tuple([ nifticlib.floatArray_frompointer(self.__nimg.pixdim)[i] 215 for i in range(1,8) ] )
216 217
218 - def getExtent(self):
219 """Returns the shape of the dataimage. 220 221 :Returns: 222 Tuple with the size in voxel/timepoints. 223 224 The order of dimensions is (x,y,z,t,u,v,w). If the image has less 225 dimensions than 7 the return tuple will be shortened accordingly. 226 227 Please note that the order of dimensions is different from the tuple 228 returned by calling `NiftiImage.data.shape`! 229 230 See also `getVolumeExtent()` and `getTimepoints()`. 231 232 The `extent` property is an alternative way to access this function. 233 """ 234 # wrap dim array in nifti image struct 235 dims_array = nifticlib.intArray_frompointer(self.__nimg.dim) 236 dims = [ dims_array[i] for i in range(8) ] 237 238 return tuple( dims[1:dims[0]+1] )
239 240
241 - def getVolumeExtent(self):
242 """Returns the size/shape of the volume(s) in the image as a tuple. 243 244 :Returns: 245 Either a 3-tuple or 2-tuple or 1-tuple depending on the available 246 dimensions in the image. 247 248 The order of dimensions in the tuple is (x [, y [, z ] ] ). 249 250 The `volextent` property is an alternative way to access this function. 251 """ 252 253 # it is save to do this even if self.extent is shorter than 4 items 254 return self.extent[:3]
255 256
257 - def getTimepoints(self):
258 """Returns the number of timepoints in the image. 259 260 In case of a 3d (or less dimension) image this method returns 1. 261 262 The `timepoints` property is an alternative way to access this 263 function. 264 """ 265 266 if len(self.extent) < 4: 267 return 1 268 else: 269 return self.extent[3]
270 271
272 - def getRepetitionTime(self):
273 """Returns the temporal distance between the volumes in a timeseries. 274 275 The `rtime` property is an alternative way to access this function. 276 """ 277 return self.__nimg.dt
278 279
280 - def setRepetitionTime(self, value):
281 """Set the repetition time of a nifti image (dt). 282 """ 283 self.__nimg.dt = float(value)
284 285
286 - def asDict(self):
287 """Returns the header data of the `NiftiImage` in a dictionary. 288 289 Note, that modifications done to this dictionary do not cause any 290 modifications in the NIfTI image. Please use the `updateFromDict()` 291 method to apply changes to the image. 292 293 The `header` property is an alternative way to access this function. 294 But please note that the `header` property cannot be used like this:: 295 296 nimg.header['something'] = 'new value' 297 298 Instead one has to get the header dictionary, modify and later reassign 299 it:: 300 301 h = nimg.header 302 h['something'] = 'new value' 303 nimg.header = h 304 """ 305 # Convert nifti_image struct into nifti1 header struct. 306 # This get us all data that will actually make it into a 307 # NIfTI file. 308 nhdr = nifticlib.nifti_convert_nim2nhdr(self.__nimg) 309 310 return nhdr2dict(nhdr)
311 312
313 - def updateFromDict(self, hdrdict):
314 """Update NIfTI header information. 315 316 Updated header data is read from the supplied dictionary. One cannot 317 modify dimensionality and datatype of the image data. If such 318 information is present in the header dictionary it is removed before 319 the update. If resizing or datatype casting are required one has to 320 convert the image data into a separate array (`NiftiImage.assarray()`) 321 and perform resize and data manipulations on this array. When finished, 322 the array can be converted into a nifti file by calling the NiftiImage 323 constructor with the modified array as 'source' and the nifti header 324 of the original NiftiImage object as 'header'. 325 """ 326 # rebuild nifti header from current image struct 327 nhdr = nifticlib.nifti_convert_nim2nhdr(self.__nimg) 328 329 # remove settings from the hdrdict that are determined by 330 # the data set and must not be modified to preserve data integrity 331 if hdrdict.has_key('datatype'): 332 del hdrdict['datatype'] 333 if hdrdict.has_key('dim'): 334 del hdrdict['dim'] 335 336 # update the nifti header 337 updateNiftiHeaderFromDict(nhdr, hdrdict) 338 339 # if no filename was set already (e.g. image from array) set a temp 340 # name now, as otherwise nifti_convert_nhdr2nim will fail 341 have_temp_filename = False 342 if not self.filename: 343 self.__nimg.fname = 'pynifti_updateheader_temp_name' 344 self.__nimg.iname = 'pynifti_updateheader_temp_name' 345 have_temp_filename = True 346 347 # recreate nifti image struct 348 new_nimg = nifticlib.nifti_convert_nhdr2nim(nhdr, self.filename) 349 if not new_nimg: 350 raise RuntimeError, \ 351 "Could not recreate NIfTI image struct from updated header." 352 353 # replace old image struct by new one 354 # be careful with memory leak (still not checked whether successful) 355 356 # assign the new image struct 357 self.__nimg = new_nimg 358 359 # reset filename if temp name was set 360 if have_temp_filename: 361 self.__nimg.fname = '' 362 self.__nimg.iname = ''
363 364
365 - def setSlope(self, value):
366 """Set the slope attribute in the NIfTI header. 367 368 Besides reading it is also possible to set the slope by assigning 369 to the `slope` property. 370 """ 371 self.__nimg.scl_slope = float(value)
372 373
374 - def setIntercept(self, value):
375 """Set the intercept attribute in the NIfTI header. 376 377 Besides reading it is also possible to set the intercept by assigning 378 to the `intercept` property. 379 """ 380 self.__nimg.scl_inter = float(value)
381 382
383 - def setDescription(self, value):
384 """Set the description element in the NIfTI header. 385 386 :Parameter: 387 value: str 388 Description -- must not be longer than 79 characters. 389 390 Besides reading it is also possible to set the description by assigning 391 to the `description` property. 392 """ 393 if len(value) > 79: 394 raise ValueError, \ 395 "The NIfTI format only supports descriptions shorter than " \ 396 + "80 chars." 397 398 self.__nimg.descrip = value
399 400
401 - def setXFormCode(self, xform, code):
402 """Set the type of space described by the NIfTI transformations. 403 404 The NIfTI format defines five coordinate system types which are used 405 to describe the target space of a transformation (qform or sform). 406 Please note, that the last four transformation types are only available 407 in the NIfTI format and not when saving into ANALYZE. 408 409 'unkown', `NIFTI_XFORM_UNKNOWN`, 0: 410 Transformation is arbitrary. This is the ANALYZE compatibility 411 mode. In this case no *sform* matrix will be written, even when 412 stored in NIfTI and not in ANALYZE format. Additionally, only the 413 pixdim parts of the *qform* matrix will be saved (upper-left 3x3). 414 'scanner', `NIFTI_XFORM_SCANNER_ANAT`, 1: 415 Scanner-based anatomical coordinates. 416 'aligned', `NIFTI_XFORM_ALIGNED_ANAT`, 2: 417 Coordinates are aligned to another file's coordinate system. 418 'talairach', `NIFTI_XFORM_TALAIRACH`, 3: 419 Coordinate system is shifted to have its origin (0,0,0) at the 420 anterior commissure, as in the Talairach-Tournoux Atlas. 421 'mni152', `NIFTI_XFORM_MNI_152`, 4: 422 Coordinates are in MNI152 space. 423 424 :Parameters: 425 xform: str('qform' | 'sform') 426 Which of the two NIfTI transformations to set. 427 code: str | `NIFTI_XFORM_CODE` | int (0..4) 428 The Transformation code can be specified either by a string, the 429 `NIFTI_XFORM_CODE` defined in the nifti1.h header file (accessible 430 via the `nifticlib` module, or the corresponding integer value 431 (see above list for all possibilities). 432 """ 433 if isinstance(code, str): 434 code = nifti_xform_map[code] 435 436 if xform == 'qform': 437 self.raw_nimg.qform_code = code 438 elif xform == 'sform': 439 self.raw_nimg.sform_code = code 440 else: 441 raise ValueError, "Unkown transformation '%s'" % xform
442 443
444 - def getSForm(self):
445 """Returns the sform matrix. 446 447 Please note, that the returned SForm matrix is not bound to the 448 NiftiImage object. Therefore it cannot be successfully modified 449 in-place. Modifications to the SForm matrix can only be done by setting 450 a new SForm matrix either by calling `setSForm()` or by assigning it to 451 the `sform` attribute. 452 453 The `sform` property is an alternative way to access this function. 454 """ 455 return nifticlib.mat442array(self.__nimg.sto_xyz)
456 457
458 - def setSForm(self, m, code='mni152'):
459 """Sets the sform matrix. 460 461 The supplied value has to be a 4x4 matrix. The matrix elements will be 462 converted to floats. By definition the last row of the sform matrix has 463 to be (0,0,0,1). However, different values can be assigned, but will 464 not be stored when the niftifile is saved. 465 466 The inverse sform matrix will be automatically recalculated. 467 468 Besides reading it is also possible to set the sform matrix by 469 assigning to the `sform` property. 470 471 :Parameters: 472 m: 4x4 ndarray 473 The sform matrix. 474 code: str | `NIFTI_XFORM_CODE` | int (0..4) 475 The type of the coordinate system the sform matrix is describing. 476 By default this coordinate system is assumed to be the MNI152 space. 477 Please refer to the `setXFormCode()` method for a full list of 478 possible codes and their meaning. 479 """ 480 if m.shape != (4, 4): 481 raise ValueError, "SForm matrix has to be of size 4x4." 482 483 # make sure it is float 484 m = m.astype('float') 485 486 nifticlib.set_mat44( self.__nimg.sto_xyz, 487 m[0,0], m[0,1], m[0,2], m[0,3], 488 m[1,0], m[1,1], m[1,2], m[1,3], 489 m[2,0], m[2,1], m[2,2], m[2,3], 490 m[3,0], m[3,1], m[3,2], m[3,3] ) 491 492 # recalculate inverse 493 self.__nimg.sto_ijk = \ 494 nifticlib.nifti_mat44_inverse( self.__nimg.sto_xyz ) 495 496 # set sform code, which decides how the sform matrix is interpreted 497 self.setXFormCode('sform', code)
498 499
500 - def getInverseSForm(self):
501 """Returns the inverse sform matrix. 502 503 Please note, that the inverse SForm matrix cannot be modified in-place. 504 One needs to set a new SForm matrix instead. The corresponding inverse 505 matrix is then re-calculated automatically. 506 507 The `sform_inv` property is an alternative way to access this function. 508 """ 509 return nifticlib.mat442array(self.__nimg.sto_ijk)
510 511
512 - def getQForm(self):
513 """Returns the qform matrix. 514 515 Please note, that the returned QForm matrix is not bound to the 516 NiftiImage object. Therefore it cannot be successfully modified 517 in-place. Modifications to the QForm matrix can only be done by setting 518 a new QForm matrix either by calling `setQForm()` or by assigning it to 519 the `qform` property. 520 """ 521 return nifticlib.mat442array(self.__nimg.qto_xyz)
522 523
524 - def getInverseQForm(self):
525 """Returns the inverse qform matrix. 526 527 The `qform_inv` property is an alternative way to access this function. 528 529 Please note, that the inverse QForm matrix cannot be modified in-place. 530 One needs to set a new QForm matrix instead. The corresponding inverse 531 matrix is then re-calculated automatically. 532 """ 533 return nifticlib.mat442array(self.__nimg.qto_ijk)
534 535
536 - def setQForm(self, m, code='scanner'):
537 """Sets the qform matrix. 538 539 The supplied value has to be a 4x4 matrix. The matrix will be converted 540 to float. 541 542 The inverse qform matrix and the quaternion representation will be 543 automatically recalculated. 544 545 Besides reading it is also possible to set the qform matrix by 546 assigning to the `qform` property. 547 548 :Parameters: 549 m: 4x4 ndarray 550 The qform matrix. 551 code: str | `NIFTI_XFORM_CODE` | int (0..4) 552 The type of the coordinate system the qform matrix is describing. 553 By default this coordinate system is assumed to be the scanner 554 anatomical space. Please refer to the `setXFormCode()` method for 555 a full list of possible codes and their meaning. 556 """ 557 if m.shape != (4, 4): 558 raise ValueError, "QForm matrix has to be of size 4x4." 559 560 # make sure it is float 561 m = m.astype('float') 562 563 nifticlib.set_mat44( self.__nimg.qto_xyz, 564 m[0,0], m[0,1], m[0,2], m[0,3], 565 m[1,0], m[1,1], m[1,2], m[1,3], 566 m[2,0], m[2,1], m[2,2], m[2,3], 567 m[3,0], m[3,1], m[3,2], m[3,3] ) 568 569 # recalculate inverse 570 self.__nimg.qto_ijk = \ 571 nifticlib.nifti_mat44_inverse( self.__nimg.qto_xyz ) 572 573 # update quaternions 574 ( self.__nimg.quatern_b, self.__nimg.quatern_c, self.__nimg.quatern_d, 575 self.__nimg.qoffset_x, self.__nimg.qoffset_y, self.__nimg.qoffset_z, 576 self.__nimg.dx, self.__nimg.dy, self.__nimg.dz, 577 self.__nimg.qfac ) = \ 578 nifticlib.nifti_mat44_to_quatern( self.__nimg.qto_xyz ) 579 580 # set qform code, which decides how the qform matrix is interpreted 581 self.setXFormCode('qform', code)
582 583
585 """Recalculates the qform matrix (and the inverse) from the quaternion 586 representation. 587 """ 588 # recalculate qform 589 self.__nimg.qto_xyz = nifticlib.nifti_quatern_to_mat44 ( 590 self.__nimg.quatern_b, self.__nimg.quatern_c, self.__nimg.quatern_d, 591 self.__nimg.qoffset_x, self.__nimg.qoffset_y, self.__nimg.qoffset_z, 592 self.__nimg.dx, self.__nimg.dy, self.__nimg.dz, 593 self.__nimg.qfac ) 594 595 596 # recalculate inverse 597 self.__nimg.qto_ijk = \ 598 nifticlib.nifti_mat44_inverse( self.__nimg.qto_xyz )
599 600
601 - def setQuaternion(self, value, code='scanner'):
602 """Set Quaternion from 3-tuple (qb, qc, qd). 603 604 The qform matrix and it's inverse are re-computed automatically. 605 606 Besides reading it is also possible to set the quaternion by assigning 607 to the `quatern` property. 608 609 :Parameters: 610 value: length-3 sequence 611 qb, qc and qd quaternions. 612 code: str | `NIFTI_XFORM_CODE` | int (0..4) 613 The type of the coordinate system the corresponding qform matrix 614 is describing. By default this coordinate system is assumed to be 615 the scanner anatomical space. Please refer to the `setXFormCode()` 616 method for a full list of possible codes and their meaning. 617 """ 618 if len(value) != 3: 619 raise ValueError, 'Requires 3-tuple.' 620 621 self.__nimg.quatern_b = float(value[0]) 622 self.__nimg.quatern_c = float(value[1]) 623 self.__nimg.quatern_d = float(value[2]) 624 625 self.updateQFormFromQuaternion() 626 self.setXFormCode('qform', code)
627 628
629 - def getQuaternion(self):
630 """Returns a 3-tuple containing (qb, qc, qd). 631 632 The `quatern` property is an alternative way to access this function. 633 """ 634 return( ( self.__nimg.quatern_b, 635 self.__nimg.quatern_c, 636 self.__nimg.quatern_d ) )
637 638
639 - def setQOffset(self, value, code='scanner'):
640 """Set QOffset from 3-tuple (qx, qy, qz). 641 642 The qform matrix and its inverse are re-computed automatically. 643 644 Besides reading it is also possible to set the qoffset by assigning 645 to the `qoffset` property. 646 647 :Parameters: 648 value: length-3 sequence 649 qx, qy and qz offsets. 650 code: str | `NIFTI_XFORM_CODE` | int (0..4) 651 The type of the coordinate system the corresponding qform matrix 652 is describing. By default this coordinate system is assumed to be 653 the scanner anatomical space. Please refer to the `setXFormCode()` 654 method for a full list of possible codes and their meaning. 655 """ 656 if len(value) != 3: 657 raise ValueError, 'Requires 3-tuple.' 658 659 self.__nimg.qoffset_x = float(value[0]) 660 self.__nimg.qoffset_y = float(value[1]) 661 self.__nimg.qoffset_z = float(value[2]) 662 663 self.updateQFormFromQuaternion() 664 self.setXFormCode('qform', code)
665 666
667 - def getQOffset(self):
668 """Returns a 3-tuple containing (qx, qy, qz). 669 670 The `qoffset` property is an alternative way to access this function. 671 """ 672 return( ( self.__nimg.qoffset_x, 673 self.__nimg.qoffset_y, 674 self.__nimg.qoffset_z ) )
675 676
677 - def setQFac(self, value, code='scanner'):
678 """Set qfac (scaling factor of qform matrix). 679 680 The qform matrix and its inverse are re-computed automatically. 681 682 Besides reading it is also possible to set the qfac by assigning 683 to the `qfac` property. 684 685 :Parameters: 686 value: float 687 Scaling factor. 688 code: str | `NIFTI_XFORM_CODE` | int (0..4) 689 The type of the coordinate system the corresponding qform matrix 690 is describing. By default this coordinate system is assumed to be 691 the scanner anatomical space. Please refer to the `setXFormCode()` 692 method for a full list of possible codes and their meaning. 693 """ 694 self.__nimg.qfac = float(value) 695 self.updateQFormFromQuaternion() 696 self.setXFormCode('qform', code)
697 698
699 - def getQOrientation(self, as_string = False):
700 """Returns to orientation of the i, j and k axis as stored in the 701 qform matrix. 702 703 By default NIfTI orientation codes are returned, but if `as_string` is 704 set to true a string representation ala 'Left-to-right' is returned 705 instead. 706 """ 707 codes = nifticlib.nifti_mat44_to_orientation(self.__nimg.qto_xyz) 708 if as_string: 709 return [ nifticlib.nifti_orientation_string(i) for i in codes ] 710 else: 711 return codes
712 713
714 - def getSOrientation(self, as_string = False):
715 """Returns to orientation of the i, j and k axis as stored in the 716 sform matrix. 717 718 By default NIfTI orientation codes are returned, but if `as_string` is 719 set to true a string representation ala 'Left-to-right' is returned 720 instead. 721 """ 722 codes = nifticlib.nifti_mat44_to_orientation(self.__nimg.sto_xyz) 723 if as_string: 724 return [ nifticlib.nifti_orientation_string(i) for i in codes ] 725 else: 726 return codes
727 728
729 - def getFilename(self):
730 """Returns the filename. 731 732 To distinguish ANALYZE from 2-file NIfTI images the image filename is 733 returned for ANALYZE images while the header filename is returned for 734 NIfTI files. 735 736 The `filename` property is an alternative way to access this function. 737 """ 738 if self.__nimg.nifti_type == nifticlib.NIFTI_FTYPE_ANALYZE: 739 return self.__nimg.iname 740 else: 741 return self.__nimg.fname
742 743 # class properties 744 # read only 745 nvox = property(fget=lambda self: self.__nimg.nvox) 746 max = property(fget=lambda self: self.__nimg.cal_max) 747 min = property(fget=lambda self: self.__nimg.cal_min) 748 sform_inv = property(fget=getInverseSForm) 749 qform_inv = property(fget=getInverseQForm) 750 extent = property(fget=getExtent) 751 volextent = property(fget=getVolumeExtent) 752 timepoints = property(fget=getTimepoints) 753 raw_nimg = property(fget=lambda self: self.__nimg) 754 755 # read and write 756 filename = property(fget=getFilename) 757 slope = property(fget=lambda self: self.__nimg.scl_slope, 758 fset=setSlope) 759 intercept = property(fget=lambda self: self.__nimg.scl_inter, 760 fset=setIntercept) 761 voxdim = property(fget=getVoxDims, fset=setVoxDims) 762 pixdim = property(fget=getPixDims, fset=setPixDims) 763 description = property(fget=lambda self: self.__nimg.descrip, 764 fset=setDescription) 765 header = property(fget=asDict, fset=updateFromDict) 766 sform = property(fget=getSForm, fset=setSForm) 767 qform = property(fget=getQForm, fset=setQForm) 768 quatern = property(fget=getQuaternion, fset=setQuaternion) 769 qoffset = property(fget=getQOffset, fset=setQOffset) 770 qfac = property(fget=lambda self: self.__nimg.qfac, fset=setQFac) 771 rtime = property(fget=getRepetitionTime, fset=setRepetitionTime)
772