Welcome to Pirt’s documentation!

Pirt is the “Python image registration toolkit”. It is a library for (elastic, i.e. non-regid) image registration of 2D and 3D images with support for groupwise registration. It has support to constrain the deformations to be “diffeomorphic”, i.e. without folding or shearing, and thus invertable.

Pirt is written in pure Python and uses Numba for speed. It depends on Numpy, Scipy, Numba. It has an optional dependency on Visvis for visualization.

Pirt implements its own interpolation functions, which, incidentally, are faster than the corresponding functions in scipy and scikit-image (after Numba’s JIT warmup).

The functionality inside Pirt is implemented over a series of submodules, but (almost) all functions and classes are available in the main namespace.

The code lives on Bitbucket. Also check out the examples or the docs.

pirt.gaussfun - Gaussian diffusion

The gaussfun module implements functions for diffusion and Gaussian derivatives for data of any dimension.

Contents of this module:

  • gaussiankernel - Create a Gaussian kernel
  • gaussiankernel2 - Create a 2D Gaussian kernel
  • diffusionkernel - Create a discrete analog to the Gaussian kernel
  • gfilter - Filter data with Gaussian (derivative) kernels
  • diffuse - Filter data with true discrete diffusion kernels
  • gfilter2 - Filter data by specifying sigma in world coordinates
  • diffuse2 - Diffuse data by specifying sigma in world coordinates
pirt.gaussiankernel(sigma, order=0, N=None, returnt=False, warn=True)

Creates a 1D gaussian derivative kernel with the given sigma and the given order. (An order of 0 is a “regular” Gaussian.)

The returned kernel is a column vector, thus working in the first dimension (in images, this often is y).

The returned kernel is odd by default. Using N one can specify the full kernel size (if not int, the ceil operator is applied). By specifying a negative value for N, the tail length (number of elements on both sides of the center element) can be specified. The total kernel size than becomes ceil(-N)*2+1. Though the method to supply it may be a bit obscure, this measure can be handy, since the tail length if often related to the sigma. If not given, the optimal N is determined automatically, depending on sigma and order.

If the given scale is a small for the given order, a warning is produced (unless warn==True).

—– Used Literature:

Koenderink, J. J. The structure of images. Biological Cybernetics 50, 5 (1984), 363-370.

Lindeberg, T. Scale-space for discrete signals. IEEE Transactions on Pattern Analysis and Machine Intelligence 12, 3 (1990), 234-254.

Ter Haar Romeny, B. M., Niessen, W. J., Wilting, J., and Florack, L. M. J. Differential structure of images: Accuracy of representation. In First IEEE International Conference on Image Processing, (Austin, TX) (1994).

pirt.diffusionkernel(sigma, N=4, returnt=False)

A discrete analog to the continuous Gaussian kernel, as proposed by Toni Lindeberg.

N is the tail length factor (relative to sigma).

pirt.gfilter(L, sigma, order=0, mode='constant', warn=True)

Gaussian filterering and Gaussian derivative filters.

Parameters:
L : np.ndarray

The input data to filter

sigma : scalar or list-of-scalars

The smoothing parameter, can be given for each dimension

order : int or list-of-ints

The order of the derivative, can be given for each dimension

mode : {‘reflect’,’constant’,’nearest’,’mirror’, ‘wrap’}

Determines how edge effects are handled. (see scipy.ndimage.convolve1d)

warn : boolean

Whether to show a warning message if the sigma is too small to represent the required derivative.

Notes

Makes use of the seperability property of the Gaussian by convolving 1D kernels in each dimension.

pirt.gfilter2(L, scale, order=0, mode='reflect', warn=True)

Apply Gaussian filtering by specifying a scale in world coordinates rather than a sigma. This function determines the sigmas to apply, based on the sampling of the elements.

See gfilter for more information.

(If L is not an Aarray, this function yields the same result as gfilter.)

pirt.diffuse(L, sigma)

Diffusion using a discrete variant of the diffusion operator.

Parameters:
L : np.ndarray

The input data to filter

sigma : scalar or list-of-scalars

The smoothing parameter, can be given for each dimension
pirt.diffuse2(L, scale, mode='nearest')

Apply diffusion by specifying a scale in world coordinates rather than a sigma. This function determines the sigmas to apply, based on the sampling of the elements.

See diffuse for more information.

(If L is not an Aarray, this function yields the same result as diffuse.)

Technically implemented in pyramid.py, but very much related to Gaussian diffusion:

class pirt.ScaleSpacePyramid(data, min_scale=None, scale_offset=0, use_buffer=False, level_factor=2)

The scale space pyramid class provides a way to manage a scale space pyramid. Given an input image (of arbitrary dimension), it provides two simple methods to obtain the image at the a specified scale or level.

Parameters:
data : numpy array

An array of any dimension. Should preferably be of float type.

min_scale : scalar, optional

The minimum scale to sample from the pyramid. If not given, scale_offset is used. If larger than zero, the image is smoothed to this scale before creating the zeroth level. If the smoothness is sufficient, the data is also downsampled. This makes a registration algorithm much faster, because the image data for the final scales does not have a unnecessary high resolution.

scale_offset : scalar

The scale of the given data. Use this if the data is already smooth. Be careful not to set this value too high, as aliasing artifacts may be introduced. Default zero.

use_buffer : bool

Whether a result obtained with get_scale() is buffered for later use. Only one image is buffered. Default False.

level_factor : scalar

The scale distance between two levels. A larger number means saving a bit of memory in trade of speed. You’re probably fine with 2.0.

Notes

Note that this scale space representation handles anisotropic arrays and that scale is expressed in world units.

Note that images at higher levels do not always have a factor 2 sampling difference with the original! This is because the first and last pixel are kept the same, and the number of pixels is decreased with factors of two (or almost a factor of two if the number is uneven).

The images always have the same offset though.

We adopt the following concepts:
  • level: the level in the pyramid. Each level is a factor two smaller in size (in each dimension) than the previous.
  • scale: the scale in world coordinates
calculate(levels=None, min_shape=None)

Create the image pyramid now. Specify either the amount of levels, or the minimum shape component of the highest level. If neither levels nor min_shape is given, uses min_shape=8.

Returns (max_level, max_sigma) of the current pyramid.

get_level(self, level)

get_level(level):

Get the image at the specified (integer) level, zero being the lowest level (the original image).

Each level is approximately a factor two smaller in size that the previous level. All levels are buffered.

The returned image has two added properties: _pyramid_scale and _pyramid_level, wich specify the image scale and level in the pyramid.

get_scale(scale)

Get the image at the specified scale (expressed in world units). For higher scales, the image has a smaller shape than the original image. If min_scale and scale_offset are not used, a scale of 0 represents the original image.

To calculate the result, the image at the level corresponding to the nearest lower scale is obtained, and diffused an extra bit to obtain the requested scale.

The result is buffered (if the pyramid was instantiated with use_buffer=True), such that calling this function multiple times with the same scale is much faster. Only buffers the last used scale.

The returned image has two added properties: _pyramid_scale and _pyramid_level, wich specify the image scale and level in the pyramid.

pirt.interp - Interpolation functions

The interp module implements several functions for interpolation, implemented in Numba.

Interpolation by warping and projecting

pirt.warp(data, samples, order='linear', spline_type=0.0)

Interpolate (sample) data at the positions specified by samples (pixel coordinates).

Parameters:
data : array (float32 or float64)

Data to interpolate, can be 1D, 2D or 3D.

samples : tuple with numpy arrays

Each array specifies the sample position for one dimension (in x-y-z order). Can also be a stacked array as in skimage’s warp() (in z-y-x order).

order : integer or string

Order of interpolation. Can be 0:’nearest’, 1:’linear’, 2: ‘quadratic’, 3:’cubic’.

spline_type : float or string

Only for cubic interpolation. Specifies the type of spline. Can be ‘Basis’, ‘Hermite’, ‘Cardinal’, ‘Catmull-rom’, ‘Lagrange’, ‘Lanczos’, ‘quadratic’, or a float, specifying the tension parameter for the Cardinal spline. See the docs of get_cubic_spline_coefs() for more information.
Returns:
result : array

The result is of the same type as the data array, and of the same shape of the samples arrays, which can be of any shape. This flexibility makes this function suitable as a base function for higher level “sampling functions”.

Notes

The input data can have up to three dimensions. It can be of any dtype, but float32 or float64 is recommended in general.

An order of interpolation of 2 would naturally correspond to quadratic interpolation. However, due to its uneven coefficients it reques the same support (and speed) as a cubic interpolant. This implementation adds the two quadratic polynomials. Note that you can probably better use order=3 with a Catmull-Rom spline, which corresponds to the linear interpolation of the two quadratic polynomials.

It can be shown (see Thevenaz et al. 2000 “Interpolation Revisited”) that interpolation using a Cardinal spline is equivalent to interpolating-B-spline interpolation.

pirt.awarp(data, samples, order='linear', spline_type=0.0)

Interpolation in anisotropic array. Like warp(), but the samples are expressed in world coordimates.

pirt.project(data, samples)

Interpolate data to the positions specified by samples (pixel coordinates).

In contrast to warp(), the project() function applies forward deformations, moving the pixels/voxels to the given locations, rather than getting the pixel values from the given locations. Although this may feel closer to how one would like to think about deformations, this function is slower and has no options to determine the interpolation, because there is no interpolation, but splatting.

Parameters:
data : array (float32 or float64)

Data to interpolate, can be 1D, 2D or 3D.

samples : tuple with numpy arrays

Each array specifies the sample position for one dimension (in x-y-z order). In contrast to warp(), each array must have the same shape as data. Can also be a stacked array as in skimage’s warp() (in z-y-x order).
Returns:
result : array

The result is of the same type and shape as the data array.
pirt.aproject(data, samples)

Interpolation in anisotropic array. Like project(), but the samples are expressed in world coordimates.

Deforming

pirt.deform_backward(data, deltas, order=1, spline_type=0.0)

Interpolate data according to the deformations specified in deltas. Deltas should be a tuple of numpy arrays similar to ‘samples’ in the warp() function. They represent the relative sample positions expressed in world coordinates.

pirt.deform_forward(data, deltas)

Like deform_backward(), but applied to project (forward deformation).

High level functions

pirt.resize(data, new_shape, order=3, spline_type=0.0, prefilter=False, extra=False)

Resize the data to the specified new shape.

Parameters:
data : numpy array

The data to rezize.

new_shape : tuple

The new shape of the data (z-y-x order).

order : {0,1,3} or {‘nearest’, ‘linear’, ‘cubic’}

The interpolation order to use.

spline_type : float or string

Only for cubic interpolation. Specifies the type of spline. Can be ‘Basis’, ‘Hermite’, ‘Cardinal’, ‘Catmull-rom’, ‘Lagrange’, ‘Lanczos’, ‘quadratic’, or a float, specifying the tension parameter for the Cardinal spline. See the docs of get_cubic_spline_coefs() for more information.

prefilter : bool

Whether to apply (discrete Gaussian diffusion) anti-aliasing (when downampling). Default False.

extra : bool

Whether to extrapolate the data a bit. In this case, each datapoint is seen as spanning a space equal to the distance between the data points. This is the method used when you resize an image using e.g. paint.net or photoshop. If False, the first and last datapoint are exactly on top of the original first and last datapoint (like scipy.ndimage.zoom). Default False.
pirt.imresize(data, factor, order=3)

Convenience function to resize the image data (1D, 2D or 3D).

This function uses pirt.resize() with ‘prefilter’ and ‘extra’ set to True. This makes it more suitble for generic image resizing. Use pirt.resize() for more fine-grained control.

Parameters:
data : numpy array

The data to rezize.

new_shape : tuple

The new shape of the data (z-y-x order).

order : {0,1,3} or {‘nearest’, ‘linear’, ‘cubic’}

The interpolation order to use.
pirt.zoom(data, factor, order=3, spline_type=0.0, prefilter=False, extra=False)

Resize the data with the specified factor. The default behavior is the same as scipy.ndimage.zoom(), but three times faster.

Parameters:
data : numpy array

The data to rezize.

factor : scalar or tuple

The resize factor, optionally for each dimension (z-y-z order).

order : {0,1,3} or {‘nearest’, ‘linear’, ‘cubic’}

The interpolation order to use.

spline_type : float or string

Only for cubic interpolation. Specifies the type of spline. Can be ‘Basis’, ‘Hermite’, ‘Cardinal’, ‘Catmull-rom’, ‘Lagrange’, ‘Lanczos’, ‘quadratic’, or a float, specifying the tension parameter for the Cardinal spline. See the docs of get_cubic_spline_coefs() for more information.

prefilter : bool

Whether to apply (discrete Gaussian diffusion) anti-aliasing (when downampling). Default False.

extra : bool

Whether to extrapolate the data a bit. In this case, each datapoint is seen as spanning a space equal to the distance between the data points. This is the method used when you resize an image using e.g. paint.net or photoshop. If False, the first and last datapoint are exactly on top of the original first and last datapoint (like numpy.zoom). Default False.
pirt.imzoom(data, factor, order=3)

Convenience function to resize the image data (1D, 2D or 3D) with the specified factor.

This function uses pirt.interp.resize() with ‘prefilter’ and ‘extra’ set to True. This makes it more suitble for generic image resizing. Use pirt.resize() for more fine-grained control.

Parameters:
data : numpy array

The data to rezize.

factor : scalar or tuple

The resize factor, optionally for each dimension (z-y-x order).

order : {0,1,3} or {‘nearest’, ‘linear’, ‘cubic’}

The interpolation order to use.

Helper functions

pirt.get_cubic_spline_coefs(t, spline_type='Catmull-Rom')

Calculates the coefficients for a cubic spline and returns them as a tuple. t is the ratio between “left” point and “right” point on the lattice. If performance matters, consider using get_lut() instead.

spline_type can be (case insensitive):

<number between -1 and 1>: Gives a Cardinal spline with the specified number as its tension parameter. A Cardinal spline is a type of Hermite spline, where the tangents are calculated using points p0 and p3; the coefficients can be directly applied to p0 p1 p2 p3.

‘Catmull-Rom’ or ‘Cardinal0’: is a cardinal spline a tension of 0. An interesting note: if we would create two quadractic splines, that would fit a polynomial “f(t) = a*t*t + b*t + c” to the first three and last three knots respectively, and if we would then combine the two using linear interpolation, we would obtain a catmull-rom spline. I don’t know whether this is how the spline was designed, or if this is a “side-effect”.

‘B-pline’: Basis spline. Not an interpolating spline, but an approximating spline. Here too, the coeffs can be applied to p0-p3.

‘Hermite’: Gives the Hermite coefficients to be applied to p1 m1 p2 m2 (in this order), where p1 p2 are the closest knots and m1 m2 the tangents.

‘Lagrange’: The Lagrange spline or Lagrange polynomial is an interpolating spline. It is the same as Newton polynomials, but implemented in a different manner (wiki). The Newton implementation should be more efficient, but this implementation is much simpler, and is very similar to the B-spline implementation (only the coefficients are different!). Also, when for example interpolating an image, coefficients are reused often and can be precalculated to enhance speed.

‘Lanczos’: Lanczos interpolation (windowed sync funcion). Note that this is not really a spline, and that sum of the coefficients is not exactly 1. Often used in audio processing. The Lanczos spline is very similar to the Cardinal spline with a tension of -0.25.

‘quadratic’: Quadratic interpolation with a support of 4, essentially the addition of the two quadratic polynoms.

‘linear’: Linear interpolation. Effective support is 2. Added for completeness and testing.

‘nearest’: Nearest neighbour interpolation. Added for completeness and testing.

pirt.meshgrid(*args)

Meshgrid implementation for 1D, 2D, 3D and beyond.

meshgrid(nx, ny) will create meshgrids with the specified shapes.

meshgrid(ndarray) will create meshgrids corresponding to the shape of the given array (which must have 2 or 3 dimension).

meshgrid([2,3,4], [1,2,3], [4,1,2]) uses the supplied values to create the grids. These lists can also be numpy arrays.

Returns a tuple with the grids in x-y-z order, with arrays of type float32.

pirt.make_samples_absolute(samples)

Note: this function is intended for sampes that represent a deformation; the number of dimensions of each array should match the number of arrays.

Given a tuple of arrays that represent a relative deformation expressed in world coordinates (x,y,z order), returns a tuple of sample arrays that represents the absolute sample locations in pixel coordinates. It is assumed that the sampling of the data is the same as for the sample arrays. The origin property is not used.

This process can also be done with relative ease by adding a meshgrid and then using awarp() or aproject(). But by combining it in one step, this process becomes faster and saves memory. Note that the deform_*() functions use this function.

Sampling a 2D slice in a 3D volume

class pirt.interp.SliceInVolume(self, pos, normal=None, previous=None)

Defines a slice in a volume.

The two span vectors are in v and u respectively. In other words, vec1 is up, vec2 is right.

convert_local_to_global(self, p2d, p3d)

Convert local 2D points to global 3D points. UNTESTED

pirt.deform - Deformations classes

The deform module implements classes to represent deformations: The DeformationGrid represents a deformation in world coordinates using a spline grid, and the DeformationField represents a deformation in world coordinates using an array for each dimension; it describes the deformation for each pixel/voxel.

The aforementioned classes are actually base classes; one should use DeformationFieldBackward, DeformationFieldForward, DeformationGridBackward or DeformationGridForward.

The base deformation class

class pirt.Deformation

This class is an abstract base class for deformation grids and deformation fields.

A deformation maps one image (1D, 2D or 3D) to another. A deformation can be represented using a B-spline grid, or using a field (i.e. array). A deformation is either forward mapping or backward mapping.

add(other)

Combine two deformations by addition.

Returns a DeformationGrid instance if both deforms are grids. Otherwise returns deformation field. The mapping (forward/backward) is taken from the left deformation.

Notes

Note that the resulting deformation is not necesarily diffeomorphic even if the individual deformations are.

Although diffeomorphisms can in general not be averaged, the constraint of Choi used in this framework enables us to do so (add the individual deformations and scale with 1/n).

This function can also be invoked using the plus operator.

apply_deformation(data, interpolation=3)

Apply the deformation to the given data. Returns the deformed data.

Parameters:
data : numpy array

The data to deform

interpolation : {0,1,3}

The interpolation order (if backward mapping is used).
as_backward()

Returns the same deformation as a backward mapping. Returns the original if already in backward mapping.

as_backward_inverse()

Returns the inverse deformation as a forward mapping. Returns the inverse of the original if already in forward mapping. If in backward mapping, the data is the same, but wrapped in a DeformationFieldBackward instance.

Note: backward and forward mappings with the same data are each-others reverse.

as_deformation_field()

Obtain a deformation fields instance that represents the same deformation. If this is already a deformation field, returns self.

as_forward()

Returns the same deformation as a forward mapping. Returns the original if already in forward mapping.

as_forward_inverse()

Returns the inverse deformation as a forward mapping. Returns the inverse of the original if already in forward mapping. If in backward mapping, the data is the same, but wrapped in a Deformation{Field/Grid}Backward instance.

Note: backward and forward mappings with the same data are each-others reverse.

as_other(other)

Returns the deformation as a forward or backward mapping, so it matches the other deformations.

compose(def1, def2)

compose(other):

Combine two deformations by composition. The left is the “static” deformation, and the right is the “delta” deformation.

Always returns a DeformationField instance. The mapping (forward/backward) of the result is taken from the left deformation.

Notes

Let “h = f.compose(g)” and “o” the mathematical composition operator. Then mathematically “h(x) = g(f(x))” or “h = g o f”.

Practically, new deformation vectors are created by sampling in one deformation, at the locations specified by the vectors of the other.

For forward mapping we sample in g at the locations of f. For backward mapping we sample in f at the locations of g.

Since we cannot impose the field values for a B-spline grid without resorting to some multi-scale approach (which would use composition and is therefore infinitely recursive), the result is always a deformation field.

If the deformation to sample in (f for forward mapping, g for backward) is a B-spline grid, the composition does not introduce any errors; sampling in a field introduces interpolation errors. Since the static deformation f is often a DeformationField, forward mapping is preferred with regard to the accuracy of composition.

copy()

Copy this deformation instance (deep copy).

field_sampling

For each dim, the sampling (distance between pixels/voxels) of the field (all 1’s if isotropic).

field_shape

The shape of the deformation field.

forward_mapping

Returns True if this deformation is forward mapping.

get_deformation_locations()

Get a tuple of arrays (x,y,z order) that represent sample locations to apply the deformation. The locations are absolute and expressed in pixel coordinates. As such, they can be fed directly to interp() or project().

get_field(d)

Get the field corresponding to the given dimension.

get_field_in_points(pp, d, interpolation=1)

Obtain the field for dimension d in the specied points. The interpolation value is used only if this is a deformation field.

The points pp should be a point set (x-y-z order).

inverse()

Get the inverse deformation. This is only valid if the current deformation is diffeomorphic. The result is always a DeformationField instance.

is_identity

Returns True if this is an identity deform, or null deform; representing no deformation at all.

ndim

The number of dimensions of the deformation.

scale(factor)

Scale the deformation (in-place) with the given factor. Note that the result is diffeomorphic only if the original is diffeomorphic and the factor is between -1 and 1.

show(axes=None, axesAdjust=True)

Illustrates 2D deformations.

It does so by creating an image of a grid and deforming it. The image is displayed in the given (or current) axes. Returns the texture object of the grid image.

Requires visvis.

The base grid deformation class

class pirt.DeformationGrid(image, sampling=5)

A deformation grid represents a deformation using a spline grid.

The ‘grids’ property returns a tuple of SplineGrid instances (one for each dimension). These sub-grids can also obtained by indexing and are in z,y,x order.

Parameters:
image : shape-tuple, numpy-array, Aarray, or FieldDescription

A description of the field that this grid applies to. The image itself is not stored, only the field’s shape and sampling are of interest.

sampling : number

The spacing of the knots in the field. (For anisotropic fields, the spacing is expressed in world units.)
classmethod from_field(field, sampling, weights=None, injective=True, frozenedge=True, fd=None)

Create a DeformationGrid from the given deformation field (z-y-x order). Also see from_field_multiscale()

The optional weights array can be used to individually weight the field elements. Where the weight is zero, the values are not evaluated. The speed can therefore be significantly improved if there are relatively few nonzero elements.

Parameters:
field : list of numpy arrays

These arrays describe the deformation field (one per dimension).

sampling : scalar

The sampling of the returned grid

weights : numpy array

This array can be used to weigh the contributions of the individual elements.

injective : bool

Whether to prevent the grid from folding. This also penetalizes large relative deformations. An injective B-spline grid is diffeomorphic.

frozenedge : bool

Whether the edges should be frozen. This can help the registration process. Also, when used in conjunction with injective, a truly diffeomorphic deformation is obtained: every input pixel maps to a point within the image boundaries.

fd : field

Field description to describe the shape and sampling of the underlying field to be deformed.
classmethod from_field_multiscale(field, sampling, weights=None, fd=None)

Create a DeformationGrid from the given deformation field (z-y-x order). Applies from_field_multiscale() for each of its subgrids.

Parameters:
field : list of numpy arrays

These arrays describe the deformation field (one per dimension).

sampling : scalar

The sampling of the returned grid

weights : numpy array

This array can be used to weigh the contributions of the individual elements.

fd : field

Field description to describe the shape and sampling of the underlying field to be deformed.
classmethod from_points(image, sampling, pp1, pp2, injective=True, frozenedge=True)

Obtains the deformation field described by the two given sets of corresponding points. The deformation describes moving the points pp1 to points pp2. Note that backwards interpolation is used, so technically, the image is re-interpolated by sampling at the points in pp2 from pixels specified by the points in pp1.

Parameters:
image : numpy array or shape

The image (of any dimension) to which the deformation applies.

sampling : scalar

The sampling of the returned grid.

pp1 : PointSet, 2D ndarray

The base points.

pp2 : PointSet, 2D ndarray

The target points.

injective : bool

Whether to prevent the grid from folding. This also penetalizes large relative deformations. An injective B-spline grid is diffeomorphic.

frozenedge : bool

Whether the edges should be frozen. This can help the registration process. Also, when used in conjunction with injective, a truly diffeomorphic deformation is obtained: every input pixel maps to a point within the image boundaries.
classmethod from_points_multiscale(image, sampling, pp1, pp2)

Obtains the deformation field described by the two given sets of corresponding points. The deformation describes moving the points pp1 to points pp2. Applies from_points_multiscale() for each of its subgrids.

See DeformationField.from_points_multiscale() for a sound alternative.

Parameters:
image : numpy array or shape

The image (of any dimension) to which the deformation applies.

sampling : scalar

The sampling of the returned grid.

pp1 : PointSet, 2D ndarray

The base points.

pp2 : PointSet, 2D ndarray

The target points.
show(axes=None, axesAdjust=True, showGrid=True)

For 2D grids, illustrates the deformation and the knots of the grid. A grid image is made that is deformed and displayed in the given (or current) axes. By default the positions of the underlying knots are also shown using markers. Returns the texture object of the grid image.

Requires visvis.

The base field deformation class

class pirt.DeformationField(*fields)

A deformation field represents a deformation using an array for each dimension, thus specifying the deformation at each pixel/voxel. The fields should be in z-y-x order. The deformation is represented in world units (not pixels, unless pixel units are used).

Can be initialized with:
  • DeformationField(field_z, field_y, field_x)
  • DeformationField(image) # Null deformation
  • DeformationField(3) # Null deformation specifying only the ndims

This class has functionality to reshape the fields. This can be usefull during registration if using a scale space pyramid.

classmethod from_field_multiscale(field, sampling, weights=None, injective=True, frozenedge=True, fd=None)

Create a DeformationGrid from the given deformation field (z-y-x order).

Uses B-spline grids in a multi-scale approach to regularize the sparse known deformation. This produces a smooth field (minimal bending energy), similar to thin plate splines.

The optional weights array can be used to individually weight the field elements. Where the weight is zero, the values are not evaluated. The speed can therefore be significantly improved if there are relatively few nonzero elements.

Parameters:
field : list of numpy arrays

These arrays describe the deformation field (one per dimension).

sampling : scalar

The smallest sampling of the B-spline grid that is used to create the field.

weights : numpy array

This array can be used to weigh the contributions of the individual elements.

injective : bool or number

Whether to prevent the grid from folding. An injective B-spline grid is diffeomorphic. When a number between 0 and 1 is given, the unfold constraint can be tightened to obtain smoother deformations.

frozenedge : bool

Whether the edges should be frozen. This can help the registration process. Also, when used in conjunction with injective, a truly diffeomorphic deformation is obtained: every input pixel maps to a point within the image boundaries.

fd : field

Field description to describe the shape and sampling of the underlying field to be deformed.

Notes

The algorithmic is based on: Lee S, Wolberg G, Chwa K-yong, Shin SY. “Image Metamorphosis with Scattered Feature Constraints”. IEEE TRANSACTIONS ON VISUALIZATION AND COMPUTER GRAPHICS. 1996;2:337–354.

The injective constraint desctribed in this paper is not applied by this method, but by the DeformationGrid, since it is method specifically for deformations.

classmethod from_points_multiscale(image, sampling, pp1, pp2, injective=True, frozenedge=True)

Obtains the deformation field described by the two given sets of corresponding points. The deformation describes moving the points pp1 to points pp2. Note that backwards interpolation is used, so technically, the image is re-interpolated by sampling at the points in pp2 from pixels specified by the points in pp1.

Uses B-spline grids in a multi-scale approach to regularize the sparse known deformation. This produces a smooth field (minimal bending energy), similar to thin plate splines.

Parameters:
image : numpy array or shape

The image (of any dimension) to which the deformation applies.

sampling : scalar

The sampling of the smallest grid to describe the deform.

pp1 : PointSet, 2D ndarray

The base points.

pp2 : PointSet, 2D ndarray

The target points.

injective : bool or number

Whether to prevent the grid from folding. An injective B-spline grid is diffeomorphic. When a number between 0 and 1 is given, the unfold constraint can be tightened to obtain smoother deformations.

frozenedge : bool

Whether the edges should be frozen. This can help the registration process. Also, when used in conjunction with injective, a truly diffeomorphic deformation is obtained: every input pixel maps to a point within the image boundaries.

Notes

The algorithmic is based on: Lee S, Wolberg G, Chwa K-yong, Shin SY. “Image Metamorphosis with Scattered Feature Constraints”. IEEE TRANSACTIONS ON VISUALIZATION AND COMPUTER GRAPHICS. 1996;2:337–354.

The injective constraint desctribed in this paper is not applied by this method, but by the DeformationGrid, since it is method specifically for deformations.

resize_field(new_shape)

Create a new DeformationField instance, where the underlying field is resized.

The parameter new_shape can be anything that can be converted to a FieldDescription instance.

If the field is already of the correct size, returns self.

test_jacobian(show=True)

Test the determinand of the field’s Jacobian. It should be all positive for the field to be diffeomorphic.

Returns the number of pixels where the Jacobian <= 0. If show==True, will show a figure with these pixels indicated.

The concrete classes to use

class pirt.DeformationGridBackward(*args, **kwargs)

A deformation grid representing a backward mapping; the field represents where the pixels in the deformed image should be sampled to in the original image.

class pirt.DeformationGridForward(*args, **kwargs)

A deformation grid representing a forward mapping; to create the deformed image, the pixels are mapped to their new locations.

class pirt.DeformationFieldBackward(*fields)

A deformation field representing a backward mapping; the field represents where the pixels in the deformed image should be sampled to in the original image.

class pirt.DeformationFieldForward(*fields)

A deformation field representing a forward mapping; to create the deformed image, the pixels are mapped to their new locations.

class pirt.DeformationIdentity

Abstract identity deformation. It is not a grid nor a field, nor is it forward or backward mapping.

It is nothing more than a simple tool to initialize a deformation with.

pirt.splinegrid - Representing data using grids

The splinegrid module implements functionality for spline grids. A spline grid is used to representa field (i.e. data) using a grid. Essentially, spline grids allow the interpolation of sparse data in an optimally smooth way by adopting a multiscale approach. The only type of spline that makes sense for this purpose is the B-spline.

In Pirt, spline grids are used to represent deformations, but any kind of data can be represented, like e.g. image data.

The GridInterface class is the base class that enables basic grid properties and functionality. The SplineGrid class implements an actual grid that represents a scalar field. The GridContainer class can be used to wrap multiple SplineGrid instances in order to represent a vector/tensor field (such as color channels or deformations).

class pirt.GridInterface(field, sampling=5)

Abstract class to define the interface of a spline grid. Implemented by the Grid and GridContainer classes.

This class provides some generic methods and properties for grids in general. Most importantly, it handles initialization of the desctiption of the grid (dimensionality, shape and sampling of the underlying field, and the shape and sampling of the grid itself).

Parameters:
field : shape-tuple, numpy-array, Aarray, or FieldDescription

A description of the field that this grid applies to. The field itself is not stored, only the field’s shape and sampling are of interest.

sampling : number

The spacing of the knots in the field. (For anisotropic fields, the spacing is expressed in world units.)
add(other_grid)

Create a new grid by adding this grid and the given grid.

copy()

Return a deep copy of the grid.

field_sampling

For each dim, the sampling of the field, i.e. the distance (in world units) between pixels/voxels (all 1’s if isotropic).

field_shape

The shape of the underlying field. (i.e. the size in each dim.)

grid_sampling

A scalar indicating the spacing (in world units) between the knots.

grid_sampling_in_pixels

For each dim, the spacing (in sub-pixels) between the knots. A dimension that has a low field sampling will have a high grid sampling in pixels (since the pixels are small, more fit between two knots).

grid_shape

The shape of the grid. (i.e. the size in each dim.)

ndim

The number of dimensions of this grid.

refine()

Refine the grid, returning a new grid instance (of the same type) that represents the same field, but with half the grid_sampling.

resize_field(new_shape)

Create a new grid, where the underlying field is reshaped. The field is still the same; it only has a different shape and sampling.

The parameter new_shape can be anything that can be converted to a FieldDescription instance.

Note that the knots arrays are shallow copied (which makes this function very efficient).

class pirt.SplineGrid(field, sampling=5)

A SplineGrid is a representation of a scalar field in N dimensions. This field is represented in a sparse way using knots, which are distributed in a uniform grid.

The manner in which these knots describe the field depends on the underlying spline being used, which is a Cubic B-spline. This spline adopts a shape corresponding to minimum bending energy, which makes them the preferred choice for many interpolation tasks. (Earlier versions of Pirt allowed setting the spline types, but to make things easier, and because the B-spline is the only sensible choice, this option was removed.)

Parameters:
field : shape-tuple, numpy-array, Aarray, FieldDescription

A description of the field that this grid applies to. The field itself is not stored, only the field’s shape and sampling are of interest.

sampling : number

The spacing of the knots in the field. (For anisotropic fields, the spacing is expressed in world units.)
classmethod from_field(field, sampling, weights=None)

Create a SplineGrid from a given field. Note that the smoothness of the grid and the extent to which the grid follows the given values. Also see from_field_multiscale()

The optional weights array can be used to individually weight the field elements. Where the weight is zero, the values are not evaluated. The speed can therefore be significantly improved if there are relatively few nonzero elements.

Parameters:
field : numpy array or shape

The field to represent with this grid.

sampling : scalar

The sampling of the returned grid.

weights : (optional) numpy array

This array can be used to weigh the contributions of the individual elements.
classmethod from_field_multiscale(field, sampling, weights=None)

Create a SplineGrid from the given field. By performing a multi-scale approach the grid adopts a minimal bending to conform to the given field.

The optional weights array can be used to individually weight the field elements. Where the weight is zero, the values are not evaluated. The speed can therefore be significantly improved if there are relatively few nonzero elements.

Parameters:
field : numpy array or shape

The field to represent with this grid.

sampling : scalar

The sampling of the returned grid.

weights : (optional) numpy array

This array can be used to weigh the contributions of the individual elements.

Notes

The algorithmic is based on: Lee, Seungyong, George Wolberg, and Sung Yong Shin. 1997. “Scattered Data Interpolation with Multilevel B-splines”. IEEE TRANSACTIONS ON VISUALIZATION AND COMPUTER GRAPHICS 3 (3): 228-244.

classmethod from_points(field, sampling, pp, values)

Create a SplineGrid from the values specified at a set of points. Note that the smoothness of the grid and the extent to which the grid follows the given values. Also see from_points_multiscale()

Parameters:
field : numpy array or shape

The image (of any dimension) to which the grid applies.

sampling : scalar

The sampling of the returned grid.

pp : PointSet, 2D ndarray

The positions (in world coordinates) at which the values are given.

values : list or numpy array

The values specified at the given positions.
classmethod from_points_multiscale(field, sampling, pp, values)

Create a SplineGrid from the values specified at a set of points. By performing a multi-scale approach the grid adopts a minimal bending to conform to the given values.

Parameters:
field : numpy array or shape

The image (of any dimension) to which the grid applies.

sampling : scalar

The sampling of the returned grid.

pp : PointSet, 2D ndarray

The positions (in world coordinates) at which the values are given.

values : list or numpy array

The values specified at the given positions.

Notes

The algorithmic is based on: Lee, Seungyong, George Wolberg, and Sung Yong Shin. 1997. “Scattered Data Interpolation with Multilevel B-splines”. IEEE TRANSACTIONS ON VISUALIZATION AND COMPUTER GRAPHICS 3 (3): 228-244.

get_field()

Obtain the full field that this grid represents.

get_field_in_points(pp)

Obtain the field in the specied points (in world coordinates).

get_field_in_samples(pp)

Obtain the field in the specied samples (a tuple with pixel coordinates, in x-y-z order).

knots

A numpy array that represent the values of the knots.

show(axes=None, axesAdjust=True, showGrid=True)

For 2D grids, shows the field and the knots of the grid. The image is displayed in the given (or current) axes. By default the positions of the underlying knots are also shown using markers. Returns the texture object of the field image.

Requires visvis.

class pirt.GridContainer(field, sampling=5)

Abstract base class that represents multiple SplineGrid instances. Since each SplineGrid instance describes a field of scalar values, the GridContainer can be used to describe vectors/tensors. Examples are color and 2D/3D deformations.

The implementing class should:
  • instantiate SplineGrid instances and append them to ‘_grids’
  • implement methods to set the grid accordingly, probably using classmethods such as from_points, from_field, etc.
grids

A tuple of subgrids.

class pirt.FieldDescription(*args)

Describes the dimensions of a field (i.e. Aarray). It stores the following properties: shape, sampling, origin

This class can for example be used to instantiate a new grid without requiring the actual field.

This class can be instantiated with a shape and sampling tuple, or with any object that describes a field in a way that we know of (e.g. SplineGrid and DeformationField instances).

Examples

  • FieldDescription(shape, sampling=None, origin=None)
  • FieldDescription(grid)
  • FieldDescription(Aarray)
  • FieldDescription(np.ndarray) # assumes unit sampling and zero origin
defined_origin

Whether the origin was explicitly defined.

defined_sampling

Whether the sampling was explicitly defined.

ndim

The number of dimensions of the field.

origin

The origin (spatial offset) of the field.

sampling

The sampling between the pixels of the field.

shape

The shape of the field.

pirt.reg - The registration algorithms

The reg module implements the various registration algorithms.

Base registration classes

class pirt.AbstractRegistration(*images, makeFloat=True)

Base class for registration of 2 or more images. This class only provides a common interface for the user.

This base class can for example be inherited by registration classes that wrap an external registration algorithm, such as Elastix.

Also see pirt.BaseRegistration.

Implements:

  • progress, timer, and visualizer objects
  • properties to handle the mapping (forward or backward)
  • functionality to get and set parameters
  • basic functionality to get resulting deformations
  • functionality to show the result (2D only)
Parameters:
None
DeformationField

Depending on whether forward or backward mapping is used, returns the DeformationFieldForward or DeformationFieldBackward class.

DeformationGrid

Depending on whether forward or backward mapping is used, returns the DeformationGridForward or DeformationGridBackward class.

forward_mapping

Whether forward (True) or backward (False) mapping is to be used internally.

get_deform(i=0)

Get the DeformationField instance for image with index i. If groupwise registration was used, this deformation field maps image i to the mean shape.

get_final_deform(i=0, j=1, mapping=None)

Get the DeformationField instance that maps image with index i to the image with index j. If groupwise registration was used, the deform is a composition of deform ‘i’ with the inverse of deform ‘j’.

Parameters:
i : int

The source image

j : int

The target image

mapping : {‘forward’, ‘backward’, Deformation instance}

Whether the result should be a forward or backward deform. When specified here, the result can be calculated with less errors than for example using result.as_forward(). If a Deformation object is given, the mapping of that deform is used.
params

Get params structure (as a Parameters object).

progress

The progress object, can be used by the algorithm to indicate its progress.

register(verbose=1, fig=None)

Perform the registration process.

Parameters:
verbose : int

Verbosity level. 0 means silent, 1 means print some, 2 means print a lot.

fig : visvis figure or None

If given, will display the registration progress in the given figure.
classmethod register_and_get_object(*ims, **params)

Classmethod to register the given images with the given parameters, and return the resulting registration object (after the registration has been performed).

set_params(params=None, **kwargs)

Set any parameters. The parameters are updated with the given dict, Parameters object, and then with the parameters given via the keyword arguments.

Note that the parameter structure can also be accessed directly via the ‘params’ propery.

show_result(self, how=None, fig=None)

Convenience method to show the registration result. Only works for two dimensional data. Requires visvis.

timer

The timer object, can be used by the algorithm to measure the processing time of the different steps in the registration algorithm.

visualizer

The visualizer object, can be used by the algorithm to display the images as they are deformed during registration.

class pirt.BaseRegistration(*images)

Inherits from pirt.AbstractRegistration.

An abstract registration class that provides common functionality shared by almost all registration algorithms.

This class maintains an image pyramid for each image, implements methods to set the delta deform, and get the deformed image at a specified scale. Further, this class implements the high level aspect of the registration algorithm that iterates through scale space.

Parameters:
mapping : {‘forward’, ‘backward’}

Whether forward or backward mapping is used. Default forward.

combine_deforms : {‘compose’, ‘add’}

How deformations are combined. Default compose. While add is used in some (older) registration algorithms, it is a coarse approximation.

final_scale : scalar

The minimum scale used during the registration process. This is the scale at which the registration ends. Default 1.0. Because calculating differentials suffer from more errors as the scale decreases, the minimum value is limited at 0.5.

scale_levels : integer

The amount of scale levels to use during the registration. Each level represents a factor of two in scale. The default (4) works for most images, but for large images or large deformations a larger value can be used.

scale_sampling : scalar

The amount of iterations for each level (i.e. between each factor two in scale). What values are reasonable depends on the specific algorithm.

smooth_scale : bool

Whether a smooth scale space should be used (default) or the scale is reduced with a factor of two each scale_sampling iterations.
get_deformed_image(i, s=0)

Get the image i at scale s, deformed with its current deform. Mainly intended for the registration algorithms, but can be of interest during development/debugging.

class pirt.GDGRegistration(*images)

Inherits from pirt.BaseRegistration.

Generic Groupwise Diffeomorphic Registration. Abstract class that provides a generic way to perform diffeomorphic groupwise registration.

Parameters:
deform_wise : {‘groupwise’, ‘pairwise’}

Whether all images are deformed simultaneously, or only the first image is deformed. When registering more than 2 images, ‘groupwise’ registration should be used. Default ‘groupwise’.

injective : bool

Whether the injectivity constraint should be used. This value should only be set to False in specific (testing?) situation; the resulting deformations are only guaranteed to be diffeomorphic if injective=True.

frozenEdge : bool

Whether the deformation is set to zero at the edges. If True (default) the resulting deformation fields are fully diffeomorphic; no pixels are mapped from the inside to the outside nor vice versa.

final_grid_sampling : scalar

The grid sampling of the grid at the final level. During the registration process, the B-spine grid sampling scales along with the scale.

grid_sampling_factor : scalar between 0 and 1

To what extent the grid sampling scales with the scale. By making this value lower than 1, the grid is relatively fine at the the higher scales, allowing for more deformations. The default is 0.5. Note that setting this value to 1 when using ‘frozenedge’ can cause the image to be ‘stuck’ at higher scales.

deform_limit : scalar

If injective is True, the deformations at each iteration are constraint by a “magic” limit. By making this limit tighter (relative to the scale), the deformations stay in reasonable bounds. This feature helps a lot for convergence. Default value is 1.

Special classes

class pirt.NullRegistration(*images)

Inherits from pirt.AbstractRegistration.

A registration algorithm that does nothing. This can be usefull to test the result if no registration would be applied.

Parameters:
None

Demons registration

The demon’s algorithm is a simple yet often effective algorithm. We’ve improved upon the algorithm by making it diffeomorphic, greatly improving the realism of the produced deformations.

class pirt.OriginalDemonsRegistration(*images)

Inherits from pirt.BaseRegistration.

The original version of the Demons algorithm. Uses Gaussian diffusion to regularize the deformation field. This is the implementation as proposed by He Wang et al. in 2005 “Validation of an accelerated ‘demons’ algorithm for deformable image registration in radiation therapy”

See also pirt.DiffeomorphicDemonsRegistration.

The speed_factor and noise_factor parameters are specific to this algorithm. Other important parameters are also listed below.

Parameters:
speed_factor : scalar

The relative force of the transform. This one of the most important parameters to tune. Default 3.0.

noise_factor : scalar

The noise factor. Default 2.5.

final_scale : scalar

The minimum scale used during the registration process. This is the scale at which the registration ends. Default 1.0. Because calculating differentials suffer from more errors as the scale decreases, the minimum value is limited at 0.5.

scale_levels : integer

The amount of scale levels to use during the registration. Each level represents a factor of two in scale. The default (4) works for most images, but for large images or large deformations a larger value can be used.

scale_sampling : scalar

The amount of iterations for each level (i.e. between each factor two in scale). Values between 20 and 30 are reasonable in most situations. Default 25. Higher values yield better results in general. The speed of the algorithm scales linearly with this value.
class pirt.DiffeomorphicDemonsRegistration(*images)

Inherits from pirt.GDGRegistration.

A variant of the Demons algorithm that is diffeomorphic. Based on the generice diffeomorphic groupwise registration (GDGRegistration) method .

See also pirt.OriginalDemonsRegistration.

The speed_factor parameter is specific to this algorithm. The noise_factor works best set at 1.0, effectively disabling its use; it is made redundant by the B-spline based regularization. Other important parameters are also listed below.

Parameters:
speed_factor : scalar

The relative force of the transform. This one of the most important parameters to tune. Default 3.0.

mapping : {‘forward’, ‘backward’}

Whether forward or backward mapping is used. Default forward.

final_scale : scalar

The minimum scale used during the registration process. This is the scale at which the registration ends. Default 1.0. Because calculating differentials suffer from more errors as the scale decreases, the minimum value is limited at 0.5.

scale_levels : integer

The amount of scale levels to use during the registration. Each level represents a factor of two in scale. The default (4) works for most images, but for large images or large deformations a larger value can be used.

scale_sampling : scalar

The amount of iterations for each level (i.e. between each factor two in scale). Values between 20 and 30 are reasonable in most situations. Default 25. Higher values yield better results in general. The speed of the algorithm scales linearly with this value.

final_grid_sampling : scalar

The grid sampling of the grid at the final level. During the registration process, the B-spine grid sampling scales along with the scale. This parameter is usually best coupled to final_scale. (When increasing final scale, this value should often be increased accordingly.)

grid_sampling_factor : scalar between 0 and 1

To what extent the grid sampling scales with the scale. By making this value lower than 1, the grid is relatively fine at the the higher scales, allowing for more deformations. The default is 0.5. Note that setting this value to 1 when using ‘frozenedge’ can cause the image to be ‘stuck’ at higher scales.

Gravity registration

class pirt.GravityRegistration(*images)

Inherits from pirt.GDGRegistration

A registration algorithm based on attraction between masses in both images, which is robust for large differences between the images.

The most important parameters to tune the algorithm with are scale_sampling, speed_factor and final_grid_sampling.

The speed_factor and mass_transforms parameters are specific to this algorithm. Other important parameters are also listed below.

Parameters:
speed_factor : scalar

The relative force of the transform. This one of the most important parameters to tune. Typical values are between 1 and 5. Default 1.

mass_transforms : int or tuple of ints

How the image is transformed to obtain the mass image. The number refers to the order of differentiation; 1 and 2 are gradient magnitude and Laplacian respectively. 0 only performs normalization to subtract the background. Can be specified for all images or for each image individually. Default 1.

mapping : {‘forward’, ‘backward’}

Whether forward or backward mapping is used. Default forward.

final_scale : scalar

The minimum scale used during the registration process. This is the scale at which the registration ends. Default 1.0. Because calculating differentials suffer from more errors as the scale decreases, the minimum value is limited at 0.5.

scale_levels : integer

The amount of scale levels to use during the registration. Each level represents a factor of two in scale. The default (4) works for most images, but for large images or large deformations a larger value can be used.

scale_sampling : scalar

The amount of iterations for each level (i.e. between each factor two in scale). For the coarse implementation, this is the amount of iterations performed before moving to the next scale (which is always a factor of two smaller). Values between 10 and 20 are reasonable in most situations. Default 15. Higher values yield better results in general. The speed of the algorithm scales linearly with this value.

final_grid_sampling : scalar

The grid sampling of the grid at the final level. During the registration process, the B-spine grid sampling scales along with the scale. This parameter is usually best coupled to final_scale. (When increasing final scale, this value should often be increased accordingly.)

grid_sampling_factor : scalar between 0 and 1

To what extent the grid sampling scales with the scale. By making this value lower than 1, the grid is relatively fine at the the higher scales, allowing for more deformations. The default is 0.5. Note that setting this value to 1 when using ‘frozenedge’ can cause the image to be ‘stuck’ at higher scales.

Elastix registration

class pirt.ElastixRegistration(im1, im2)

Inherits from pirt.AbstractRegistration.

Registration class for registration using the Elastix toolkit. [http://elastix.isi.uu.nl/]

This class performs a bspline transformation by default. See also the convenience subclasses.

The params property this class returns a struct with a few common Elastix parameters. the params2 property contains another set of more advanced parameters. Note that any parameter known to elastix can be added to the parameter structures, which enables tuning the registration in great detail. Also note that two parameter structs can be combined by adding them.

Because the parameters directly represent the parameters for the Elastix toolkit, their names do not follow the style of most other registration objects in this package. Here we lists some of the common parameters, for more details we refer to the elastix manual.

Parameters:
FinalGridSpacingInPhysicalUnits : int

When using the BSplineTransform, the final spacing of the grid. This controls the smoothness of the final deformation.

NumberOfResolutions : int

Most registration algorithms adopt a multiresolution approach to direct the solution towards a global optimum and to speed up the process. This parameter specifies the number of scales to apply the registration at. (default 4)

MaximumNumberOfIterations : int

Maximum number of iterations in each resolution level. 200-2000 works usually fine for nonrigid registration. The more, the better, but the longer computation time. This is an important parameter! (default 500)
get_elastix_deformed_image()

Get the deformed input image found as Elastix created it. This should be the same (except for small interpolation errors) to the image obtained using reg.get_deformed_image(0).

class pirt.ElastixRegistration_rigid(*args)

Registration class for registration using the Elastix toolkit. [http://elastix.isi.uu.nl/]

This class performs a rigid transformation by default. See pirt.ElastixRegistration for more details.

Parameters:
AutomaticScalesEstimation : bool

When using a rigid or affine transform. Scales the affine matrix elements compared to the translations, to make sure they are in the same range. In general, it’s best to use automatic scales estimation.

AutomaticTransformInitialization : bool

When using a rigid or affine transform. Automatically guess an initial translation by aligning the geometric centers of the fixed and moving.

NumberOfResolutions : int

Most registration algorithms adopt a multiresolution approach to direct the solution towards a global optimum and to speed up the process. This parameter specifies the number of scales to apply the registration at. (default 4)

MaximumNumberOfIterations : int

Maximum number of iterations in each resolution level. 200-2000 works usually fine for nonrigid registration. The more, the better, but the longer computation time. This is an important parameter! (default 500)
class pirt.ElastixRegistration_affine(*args)

Registration class for registration using the Elastix toolkit. [http://elastix.isi.uu.nl/]

This class performs an affine transformation by default. See pirt.ElastixRegistration for more details.

Parameters:
AutomaticScalesEstimation : bool

When using a rigid or affine transform. Scales the affine matrix elements compared to the translations, to make sure they are in the same range. In general, it’s best to use automatic scales estimation.

AutomaticTransformInitialization : bool

When using a rigid or affine transform. Automatically guess an initial translation by aligning the geometric centers of the fixed and moving.

NumberOfResolutions : int

Most registration algorithms adopt a multiresolution approach to direct the solution towards a global optimum and to speed up the process. This parameter specifies the number of scales to apply the registration at. (default 4)

MaximumNumberOfIterations : int

Maximum number of iterations in each resolution level. 200-2000 works usually fine for nonrigid registration. The more, the better, but the longer computation time. This is an important parameter! (default 500)
class pirt.ElastixGroupwiseRegistration(*images)

Inherits from pirt.ElastixRegistration.

Registration class for registration using the Elastix toolkit. [http://elastix.isi.uu.nl/]

This variant uses the groupwise registration approach as proposed by Metz et al. “Nonrigid registration of dynamic medical imaging data using nD+t B-splines and a groupwise optimization approach”

The params property this class returns a struct with a few common Elastix parameters. the params2 property contains another set of more advanced parameters. Note that any parameter known to elastix can be added to the parameter structures, which enables tuning the registration in great detail. Also note that two parameter structs can be combined by adding them.

Because the parameters directly represent the parameters for the Elastix toolkit, their names do not follow the style of most other registration objects in this package.

Parameters:
FinalGridSpacingInPhysicalUnits : int

When using the BSplineTransform, the final spacing of the grid. This controls the smoothness of the final deformation.

NumberOfResolutions : int

Most registration algorithms adopt a multiresolution approach to direct the solution towards a global optimum and to speed up the process. This parameter specifies the number of scales to apply the registration at. (default 4)

MaximumNumberOfIterations : int

Maximum number of iterations in each resolution level. 200-2000 works usually fine for nonrigid registration. The more, the better, but the longer computation time. This is an important parameter! (default 500)

Anisotropic arrays and pointsets

Pirt uses two classes to represent anisotropic data and points. These allow taking anisotropy into account throughout the registration process (all the way to the interpolation). The visvis library is “aware” of the Aarray class (well, it checks arrays for a sampling and origin attribute, as duck-typing goes), and the PointSet class is really just a numpy nd array. Therefore (anisotropic) data represented with these classes can be correctly visualized with Visvis without applying any transformations.

class pirt.Aarray(shape_or_array, sampling=None, origin=None, fill=None, dtype='float32', **kwargs)

Anisotropic array; inherits from numpy.ndarray and adds a sampling and origin property which gives the sample distance and offset for each dimension.

Parameters:
shape_or_array : shape-tuple or numpy.ndarray

Specifies the shape of the produced array. If an array instance is given, the returned Aarray is a view of the same data (i.e. no data is copied).

sampling : tuple of ndim elements

Specifies the sample distance (i.e. spacing between elements) for each dimension. Default is all ones.

origin : tuple of ndim elements

Specifies the world coordinate at the first element for each dimension. Default is all zeros.

fill : scalar (optional)

If given, and the first argument is not an existing array, fills the array with this given value.

dtype : any valid numpy data type

The type of the data

All extra arguments are fed to the constructor of numpy.ndarray.
get_end()

Get the end of the array expressed in world coordinates.

get_size()

Get the size (as a vector) of the array expressed in world coordinates.

get_start()

Get the origin of the array expressed in world coordinates. Differs from the property ‘origin’ in that this method returns a point rather than indices z,y,x.

index_to_point(*index)

Given a multidimensional index, get the corresponding point in world coordinates.

origin

A tuple with the origin for each dimension.

point_to_index(point, non_on_index_error=False)

Given a point returns the sample index (z,y,x,..) closest to the given point. Returns a tuple with as many elements as there are dimensions.

If the point is outside the array an IndexError is raised by default, and None is returned when non_on_index_error == True.

sample(point, default=None)

Take a sample of the array, given the given point in world-coordinates, i.e. transformed using sampling. By default raises an IndexError if the point is not inside the array, and returns the value of “default” if it is given.

sampling

A tuple with the sample distance for each dimension.

class pirt.PointSet

The PointSet class can be used to represent sets of points or vectors, as well as singleton points. The dimensionality of the vectors in the pointset can be anything, and the dtype can be any of those supported by numpy.

This class inherits from np.ndarray, which makes it very flexible; you can threat it as a regular array, and also pass it to functions that require a numpy array. The shape of the array is NxD, with N the number of points, and D the dimensionality of each point.

This class has a __repr__ that displays a pointset-aware description. To see the underlying array, use print, or use pointset[…] to convert to a pure numpy array.

Parameters:
input : various

If input is in integer, it specifies the dimensionality of the array, and an empty pointset is created. If input is a list, it specifies a point, with which the pointset is initialized. If input is a numpy array, the pointset is a view on that array (ndim must be 2).

dtype : dtype descrsiption

The data type of the numpy array. If not given, the result will be float32.
angle(self, *p)

Calculate the angle (in radians) between two vectors. For 2D uses the arctan2 method so the angle has a sign. For 3D the angle is the smallest angles between the two vectors.

If no point is given, the angle is calculated relative to the positive x-axis.

angle2(self, *p)

Calculate the angle (in radians) of the vector between two points.

Say we have p1=(3,4) and p2=(2,1). p1.angle(p2) returns the difference of the angles of the two vectors: 0.142 = 0.927 - 0.785

p1.angle2(p2) returns the angle of the difference vector (1,3): p1.angle2(p2) == (p1-p2).angle()

append(self, *p)

Append a point to this pointset. One can give the elements of the points as separate arguments. Alternatively, a tuple or numpy array can be given.

can_resize

Whether points can be appended to/removed from this pointset. This can be False if the array does not own its own data or when it is not contiguous. In that case, one should make a copy first.

contains(self, *p)

Check whether the given point is already in this set.

cross(self, *p)

Calculate the cross product of two 3D vectors. Given two vectors, returns the vector that is orthogonal to both vectors. The right hand rule is applied; this vector is the middle finger, the argument the index finger, the returned vector points in the direction of the thumb.

distance(self, *p)

Calculate the Euclidian distance between two points or pointsets. Use norm() to calculate the length of a vector.

dot(self, *p)

Calculate the dot product of two pointsets. The dot product is the standard inner product of the orthonormal Euclidean space. The sizes of the point sets should match, or one point set should be singular.

extend(self, data)

Extend the point set with more points. The shape[1] of the given data must match with that of this array.

insert(self, index, *p)

Insert a point at the given index.

norm(self)

Calculate the norm (length) of the vector. This is the same as the distance to the origin, but implemented a bit faster.

normal(self)

Calculate the normalized normal of a vector. Use (p1-p2).normal() to calculate the normal of the line p1-p2. Only works on 2D points. For 3D points use cross().

normalize(self)

Return normalized vector (to unit length).

pop(self, index=-1)

Remove and returns a point from the pointset. Removes the last by default (which is more efficient than popping from anywhere else).

ravel([order])

Return a flattened array.

Refer to numpy.ravel for full documentation.

See also

numpy.ravel
equivalent function
ndarray.flat
a flat iterator on the array.
remove(self, *p, **kwargs)

Remove the given point from the point set. Produces an error if such a point is not present. If the keyword argument all is given and True, all occurances of that point are removed. Otherwise only the first occurance is removed.

pirt.apps - Simple demo apps

The apps module implements a few app classes for demo purposes. In some cases the apps can be used to collect user input.

To use, the apps should be imported explicitly.

Indices and tables