NiBabel

Access a cacophony of neuro-imaging file formats

Previous topic

nibabel.nifti1.Nifti1Image

Next topic

nibabel.nifti1.Nifti1PairHeader

Reggie -- the one

nibabel.nifti1.Nifti1Pair

digraph inheritance5815df6a85 { rankdir=LR; ratio=compress; fontsize=14; size="6.0, 8.0"; "AnalyzeImage" [shape=ellipse,URL="nibabel.analyze.AnalyzeImage.html#nibabel.analyze.AnalyzeImage",fontname=Vera Sans, DejaVu Sans, Liberation Sans, Arial, Helvetica, sans,fontsize=14,color=dodgerblue1,style=filled,height=0.75]; "SpatialImage" -> "AnalyzeImage" [arrowsize=0.5,style="setlinewidth(0.5)"]; "SpatialImage" [shape=ellipse,URL="nibabel.spatialimages.SpatialImage.html#nibabel.spatialimages.SpatialImage",fontname=Vera Sans, DejaVu Sans, Liberation Sans, Arial, Helvetica, sans,fontsize=14,color=dodgerblue1,style=filled,height=0.75]; "Nifti1Pair" [shape=ellipse,URL="#nibabel.nifti1.Nifti1Pair",fontname=Vera Sans, DejaVu Sans, Liberation Sans, Arial, Helvetica, sans,fontsize=14,color=dodgerblue1,style=filled,height=0.75]; "AnalyzeImage" -> "Nifti1Pair" [arrowsize=0.5,style="setlinewidth(0.5)"]; }

class nibabel.nifti1.Nifti1Pair(data, affine, header=None, extra=None, file_map=None)

Initialize image

The image is a combination of (array, affine matrix, header), with optional metadata in extra, and filename / file-like objects contained in the file_map mapping.

Parameters :

data : object

image data. It should be some object that retuns an array from np.asanyarray. It should have a shape attribute or property

affine : None or (4,4) array-like

homogenous affine giving relationship between voxel coordinates and world coordinates. Affine can also be None. In this case, obj.get_affine() also returns None, and the affine as written to disk will depend on the file format.

header : None or mapping or header instance, optional

metadata for this image format

extra : None or mapping, optional

metadata to associate with image that cannot be stored in the metadata of this image type

file_map : mapping, optional

mapping giving file information for this image format

get_qform(coded=False)

Return 4x4 affine matrix from qform parameters in header

Parameters :

coded : bool, optional

If True, return {affine or None}, and qform code. If False, just return affine. {affine or None} means, return None if qform code == 0, and affine otherwise.

Returns :

affine : None or (4,4) ndarray

If coded is False, always return affine reconstructed from qform quaternion. If coded is True, return None if qform code is 0, else return the affine.

code : int

Qform code. Only returned if coded is True.

get_sform(coded=False)

Return 4x4 affine matrix from sform parameters in header

Parameters :

coded : bool, optional

If True, return {affine or None}, and sform code. If False, just return affine. {affine or None} means, return None if sform code == 0, and affine otherwise.

Returns :

affine : None or (4,4) ndarray

If coded is False, always return affine from sform fields. If coded is True, return None if sform code is 0, else return the affine.

code : int

Sform code. Only returned if coded is True.

header_class

alias of Nifti1PairHeader

set_qform(affine, code=None, strip_shears=True, **kwargs)

Set qform header values from 4x4 affine

Parameters :

hdr : nifti1 header

affine : None or 4x4 array

affine transform to write into sform. If None, only set code.

code : None, string or integer

String or integer giving meaning of transform in affine. The default is None. If code is None, then: * If affine is None, code-> 0 * If affine not None and sform code in header == 0, code-> 2

(aligned)

  • If affine not None and sform code in header != 0, code-> sform code in header

strip_shears : bool, optional

Whether to strip shears in affine. If True, shears will be silently stripped. If False, the presence of shears will raise a HeaderDataError

update_affine : bool, optional

Whether to update the image affine from the header best affine after setting the qform. Must be keyword argumemt (because of different position in set_qform). Default is True

Examples

>>> data = np.arange(24).reshape((2,3,4))
>>> aff = np.diag([2, 3, 4, 1])
>>> img = Nifti1Pair(data, aff)
>>> img.get_qform()
array([[ 2.,  0.,  0.,  0.],
       [ 0.,  3.,  0.,  0.],
       [ 0.,  0.,  4.,  0.],
       [ 0.,  0.,  0.,  1.]])
>>> img.get_qform(coded=True)
(None, 0)
>>> aff2 = np.diag([3, 4, 5, 1])
>>> img.set_qform(aff2, 'talairach')
>>> qaff, code = img.get_qform(coded=True)
>>> np.all(qaff == aff2)
True
>>> int(code)
3
set_sform(affine, code=None, **kwargs)

Set sform transform from 4x4 affine

Parameters :

hdr : nifti1 header

affine : None or 4x4 array

affine transform to write into sform. If None, only set code

code : None, string or integer

String or integer giving meaning of transform in affine. The default is None. If code is None, then: * If affine is None, code-> 0 * If affine not None and sform code in header == 0, code-> 2

(aligned)

  • If affine not None and sform code in header != 0, code-> sform code in header

update_affine : bool, optional

Whether to update the image affine from the header best affine after setting the qform. Must be keyword argumemt (because of different position in set_qform). Default is True

Examples

>>> data = np.arange(24).reshape((2,3,4))
>>> aff = np.diag([2, 3, 4, 1])
>>> img = Nifti1Pair(data, aff)
>>> img.get_sform()
array([[ 2.,  0.,  0.,  0.],
       [ 0.,  3.,  0.,  0.],
       [ 0.,  0.,  4.,  0.],
       [ 0.,  0.,  0.,  1.]])
>>> saff, code = img.get_sform(coded=True)
>>> saff
array([[ 2.,  0.,  0.,  0.],
       [ 0.,  3.,  0.,  0.],
       [ 0.,  0.,  4.,  0.],
       [ 0.,  0.,  0.,  1.]])
>>> int(code)
2
>>> aff2 = np.diag([3, 4, 5, 1])
>>> img.set_sform(aff2, 'talairach')
>>> saff, code = img.get_sform(coded=True)
>>> np.all(saff == aff2)
True
>>> int(code)
3
update_header()

Harmonize header with image data and affine

See AnalyzeImage.update_header for more examples

Examples

>>> data = np.zeros((2,3,4))
>>> affine = np.diag([1.0,2.0,3.0,1.0])
>>> img = Nifti1Image(data, affine)
>>> hdr = img.get_header()
>>> np.all(hdr.get_qform() == affine)
True
>>> np.all(hdr.get_sform() == affine)
True