NiBabel

Access a cacophony of neuro-imaging file formats

Table Of Contents

Previous topic

parrec2nii_cmd

Next topic

py3k

Reggie -- the one

processing

Image processing functions for:

  • smoothing
  • resampling
  • converting sd to and from FWHM

Smoothing and resampling routines need scipy

adapt_affine(affine, n_dim) Adapt input / output dimensions of spatial affine for n_dims
fwhm2sigma(fwhm) Convert a FWHM value to sigma in a Gaussian kernel.
resample_from_to(from_img, to_vox_map[, ...]) Resample image from_img to mapped voxel space to_vox_map
resample_to_output(in_img[, voxel_sizes, ...]) Resample image in_img to output voxel axes (world space)
sigma2fwhm(sigma) Convert a sigma in a Gaussian kernel to a FWHM value
smooth_image(img, fwhm[, mode, cval, out_class]) Smooth image img along voxel axes by FWHM fwhm millimeters

adapt_affine

nibabel.processing.adapt_affine(affine, n_dim)

Adapt input / output dimensions of spatial affine for n_dims

Adapts a spatial (4, 4) affine that is being applied to an image with fewer than 3 spatial dimensions, or more than 3 dimensions. If there are more than three dimensions, assume an identity transformation for these dimensions.

Parameters:

affine : array-like

affine transform. Usually shape (4, 4). For what follows N, M = affine.shape

n_dims : int

Number of dimensions of underlying array, and therefore number of input dimensions for affine.

Returns:

adapted : shape (M, n_dims+1) array

Affine array adapted to number of input dimensions. Columns of the affine corresponding to missing input dimensions have been dropped, columns corresponding to extra input dimensions have an extra identity column added

fwhm2sigma

nibabel.processing.fwhm2sigma(fwhm)

Convert a FWHM value to sigma in a Gaussian kernel.

Parameters:

fwhm : array-like

FWHM value or values

Returns:

sigma : array or float

sigma values corresponding to fwhm values

Examples

>>> sigma = fwhm2sigma(6)
>>> sigmae = fwhm2sigma([6, 7, 8])
>>> sigma == sigmae[0]
True

resample_from_to

nibabel.processing.resample_from_to(from_img, to_vox_map, order=3, mode='constant', cval=0.0, out_class=<class 'nibabel.nifti1.Nifti1Image'>)

Resample image from_img to mapped voxel space to_vox_map

Resample using N-d spline interpolation.

Parameters:

from_img : object

Object having attributes dataobj, affine, header and shape. If out_class is not None, img.__class__ should be able to construct an image from data, affine and header.

to_vox_map : image object or length 2 sequence

If object, has attributes shape giving input voxel shape, and affine giving mapping of input voxels to output space. If length 2 sequence, elements are (shape, affine) with same meaning as above. The affine is a (4, 4) array-like.

order : int, optional

The order of the spline interpolation, default is 3. The order has to be in the range 0-5 (see scipy.ndimage.affine_transform)

mode : str, optional

Points outside the boundaries of the input are filled according to the given mode (‘constant’, ‘nearest’, ‘reflect’ or ‘wrap’). Default is ‘constant’ (see scipy.ndimage.affine_transform)

cval : scalar, optional

Value used for points outside the boundaries of the input if mode='constant'. Default is 0.0 (see scipy.ndimage.affine_transform)

out_class : None or SpatialImage class, optional

Class of output image. If None, use from_img.__class__.

Returns:

out_img : object

Image of instance specified by out_class, containing data output from resampling from_img into axes aligned to the output space of from_img.affine

resample_to_output

nibabel.processing.resample_to_output(in_img, voxel_sizes=None, order=3, mode='constant', cval=0.0, out_class=<class 'nibabel.nifti1.Nifti1Image'>)

Resample image in_img to output voxel axes (world space)

Parameters:

in_img : object

Object having attributes dataobj, affine, header. If out_class is not None, img.__class__ should be able to construct an image from data, affine and header.

voxel_sizes : None or sequence

Gives the diagonal entries of out_img.affine` (except the trailing 1 for the homogenous coordinates) (``out_img.affine == np.diag(voxel_sizes + [1])). If None, return identity out_img.affine. If scalar, interpret as vector [voxel_sizes] * len(in_img.shape).

order : int, optional

The order of the spline interpolation, default is 3. The order has to be in the range 0-5 (see scipy.ndimage.affine_transform).

mode : str, optional

Points outside the boundaries of the input are filled according to the given mode (‘constant’, ‘nearest’, ‘reflect’ or ‘wrap’). Default is ‘constant’ (see scipy.ndimage.affine_transform).

cval : scalar, optional

Value used for points outside the boundaries of the input if mode='constant'. Default is 0.0 (see scipy.ndimage.affine_transform).

out_class : None or SpatialImage class, optional

Class of output image. If None, use in_img.__class__.

Returns:

out_img : object

Image of instance specified by out_class, containing data output from resampling in_img into axes aligned to the output space of in_img.affine

sigma2fwhm

nibabel.processing.sigma2fwhm(sigma)

Convert a sigma in a Gaussian kernel to a FWHM value

Parameters:

sigma : array-like

sigma value or values

Returns:

fwhm : array or float

fwhm values corresponding to sigma values

Examples

>>> fwhm = sigma2fwhm(3)
>>> fwhms = sigma2fwhm([3, 4, 5])
>>> fwhm == fwhms[0]
True

smooth_image

nibabel.processing.smooth_image(img, fwhm, mode='nearest', cval=0.0, out_class=<class 'nibabel.nifti1.Nifti1Image'>)

Smooth image img along voxel axes by FWHM fwhm millimeters

Parameters:

img : object

Object having attributes dataobj, affine, header and shape. If out_class is not None, img.__class__ should be able to construct an image from data, affine and header.

fwhm : scalar or length 3 sequence

FWHM in mm over which to smooth. The smoothing applies to the voxel axes, not to the output axes, but is in millimeters. The function adjusts the FWHM to voxels using the voxel sizes calculated from the affine. A scalar implies the same smoothing across the spatial dimensions of the image, but 0 smoothing over any further dimensions such as time. A vector should be the same length as the number of image dimensions.

mode : str, optional

Points outside the boundaries of the input are filled according to the given mode (‘constant’, ‘nearest’, ‘reflect’ or ‘wrap’). Default is ‘nearest’. This is different from the default for scipy.ndimage.affine_transform, which is ‘constant’. ‘nearest’ might be a better choice when smoothing to the edge of an image where there is still strong brain signal, otherwise this signal will get blurred towards zero.

cval : scalar, optional

Value used for points outside the boundaries of the input if mode='constant'. Default is 0.0 (see scipy.ndimage.affine_transform).

out_class : None or SpatialImage class, optional

Class of output image. If None, use img.__class__.

Returns:

smoothed_img : object

Image of instance specified by out_class, containing data output from smoothing img data by given FWHM kernel.