align
#


VerbosityLevels This enum defines the four levels of verbosity we use in the align module. 
Module: align._public
#
Registration API: simplified API for registration of MRI data and of streamlines.

Register a 2D/3D source image (moving) to a 2D/3D target image (static). 

Register DWI data to a template through the B0 volumes. 

Write out a syn registration mapping to a nifti file. 

Read a syn registration mapping from a nifti file. 

Resample an image (moving) from one space to another (static). 

Find the affine transformation between two 3D images. 

Implements a center of mass transform. 

Implements a translation transform. 

Implements a rigid transform. 

Implements a rigid isoscaling transform. 

Implements a rigid scaling transform. 

Implements an affine transform. 
dict() > new empty dictionary dict(mapping) > new dictionary initialized from a mapping object's (key, value) pairs dict(iterable) > new dictionary initialized as if via: d = {} for k, v in iterable: d[k] = v dict(**kwargs) > new dictionary initialized with the name=value pairs in the keyword argument list. For example: dict(one=1, two=2). 


Register a series to a reference image. 

Register a DWI series to the mean of the B0 images in that series. 

Apply a motion correction to a DWI dataset (BetweenVolumes Motion correction) 

Register two collections of streamlines ('bundles') to each other. 
Module: align.bundlemin
#

Determine the effective number of threads to be used for OpenMP calls 

Minimum direct flipped distance matrix between two streamline sets 
Module: align.cpd
#
Note#
This file is copied (possibly with major modifications) from the sources of the pycpd project  siavashk/pycpd. It remains licensed as the rest of PyCPD (MIT license as of October 2010).
# ## ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ## # # See COPYING file distributed along with the PyCPD package for the # copyright and license terms. # # ## ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ##

Deformable point cloud registration. 



Calculate num_eig eigenvectors and eigenvalues of gaussian matrix G. 

Initialize the variance (sigma2). 

Calculate eigenvectors and eigenvalues of gaussian matrix G. 
Module: align.crosscorr
#
Utility functions used by the Cross Correlation (CC) metric

Gradient of the CC Metric w.r.t. 

Gradient of the CC Metric w.r.t. 

Gradient of the CC Metric w.r.t. 

Gradient of the CC Metric w.r.t. 

Precomputations to quickly compute the gradient of the CC Metric 

Precomputations to quickly compute the gradient of the CC Metric 

Precomputations to quickly compute the gradient of the CC Metric 

Precomputations to quickly compute the gradient of the CC Metric 
Module: align.expectmax
#

Demons step for EM metric in 2D Computes the demons step [Vercauteren09] for SSDdriven registration ( eq. 

Demons step for EM metric in 3D Computes the demons step [Vercauteren09] for SSDdriven registration ( eq. 

Computes the mean and std. 

Computes the mean and std. 

Quantizes a 2D image to num_levels quantization levels 

Quantizes a 3D volume to num_levels quantization levels 
Module: align.imaffine
#
Affine image registration module consisting of the following classes:
 AffineMap: encapsulates the necessary information to perform affine
transforms between two domains, defined by a static and a moving image. The domain of the transform is the set of points in the static image’s grid, and the codomain is the set of points in the moving image. When we call the transform method, AffineMap maps each point x of the domain (static grid) to the codomain (moving grid) and interpolates the moving image at that point to obtain the intensity value to be placed at x in the resulting grid. The transform_inverse method performs the opposite operation mapping points in the codomain to points in the domain.
 ParzenJointHistogram: computes the marginal and joint distributions of
intensities of a pair of images, using Parzen windows [Parzen62] with a cubic spline kernel, as proposed by Mattes et al. [Mattes03]. It also computes the gradient of the joint histogram w.r.t. the parameters of a given transform.
 MutualInformationMetric: computes the value and gradient of the mutual
information metric the way Optimizer needs them. That is, given a set of transform parameters, it will use ParzenJointHistogram to compute the value and gradient of the joint intensity histogram evaluated at the given parameters, and evaluate the value and gradient of the histogram’s mutual information.
 AffineRegistration: it runs the multiresolution registration, putting
all the pieces together. It needs to create the scale space of the images and run the multiresolution registration by using the Metric and the Optimizer at each level of the Gaussian pyramid. At each level, it will setup the metric to compute value and gradient of the metric with the input images with different levels of smoothing.
References#
 [Parzen62] E. Parzen. On the estimation of a probability density
function and the mode. Annals of Mathematical Statistics, 33(3), 10651076, 1962.
 [Mattes03] Mattes, D., Haynor, D. R., Vesselle, H., Lewellen, T. K.,
& Eubank, W. PETCT image registration in the chest using freeform deformations. IEEE Transactions on Medical Imaging, 22(1), 1208, 2003.






dict() > new empty dictionary dict(mapping) > new dictionary initialized from a mapping object's (key, value) pairs dict(iterable) > new dictionary initialized as if via: d = {} for k, v in iterable: d[k] = v dict(**kwargs) > new dictionary initialized with the name=value pairs in the keyword argument list. For example: dict(one=1, two=2). 


Transformation to align the center of mass of the input images. 

Transformation to align the geometric center of the input images. 

Transformation to align the origins of the input images. 
Module: align.imwarp
#
Classes and functions for Symmetric Diffeomorphic Registration






Registration Stages 

Instances of the Logger class represent a single logging channel. 


Returns the matrix product A.dot(B) considering None as the identity 

Extracts the rotational and spacing components from a matrix 
Module: align.metrics
#
Metrics for Symmetric Diffeomorphic Registration









Multiresolution GaussSeidel solver using Vtype cycles 

Multiresolution GaussSeidel solver using Vtype cycles 
Module: align.parzenhist
#



Computes the mutual information and its gradient (if requested) 

Evaluates the cubic spline at a set of values 
Evaluates the cubic spline derivative at a set of values 


Take floor(total_voxels/k) samples from a (2D or 3D) grid 
Module: align.reslice
#

Reslice data with new voxel resolution defined by 
Module: align.scalespace
#




Instances of the Logger class represent a single logging channel. 
Module: align.streamlinear
#



Bundlebased Minimum Distance aka BMD 

Bundlebased Minimum Distance aka BMD 

Asymmetric Bundlebased Minimum distance. 

Bundlebased Sum Distance aka BMD 

Bundlebased Minimum Distance for joint optimization. 






Instances of the Logger class represent a single logging channel. 


MDF distance optimization function (SUM). 

MDFbased pairwise distance optimization function (MIN). 

MDFbased pairwise distance optimization function (MIN). 
MDFbased pairwise distance optimization function (MIN). 




Progressive SLR. 

Utility function for registering large tractograms. 

Function to perform unbiased groupwise bundle registration. 

Make unique pairs from n_bundle bundles. 

Compose a 4x4 transformation matrix. 

Given a 4x4 homogeneous matrix return the parameter vector. 
Module: align.streamwarp
#

Find average Euclidean length of the bundle in mm. 

Find unmatched streamline indices in moving bundle. 

Register two bundles using nonlinear method. 

Calculate vector fields. 

Calculate bundle shape difference profile. 
Module: align.sumsqdiff
#
Utility functions used by the Sum of Squared Differences (SSD) metric

Sum of squared differences between two 2D images 

Sum of squared differences between two 3D volumes 
The residual displacement field to be fit on the next iteration 

The residual displacement field to be fit on the next iteration 


Demons step for 2D SSDdriven registration Computes the demons step for SSDdriven registration ( eq. 

Demons step for 3D SSDdriven registration Computes the demons step for SSDdriven registration ( eq. 
One iteration of a large linear system solver for 2D SSD registration 

One iteration of a large linear system solver for 3D SSD registration 


Solves a 2variable symmetric positivedefinite linear system Solves the symmetric positivedefinite linear system \(Mx = y\) given by:: M = [[A[0], A[1]], [A[1], A[2]]] Parameters  A : array, shape (3,) the array containing the entries of the symmetric 2x2 matrix y : array, shape (2,) righthand side of the system to be solved Returns  out : array, shape (2,) the array the output will be stored in 

Solves a 3variable symmetric positivedefinite linear system Solves the symmetric semipositivedefinite linear system \(Mx = y\) given by \(M = (g g^{T} + \tau I)\). 
Module: align.transforms
#
Base class (contract) for all transforms for affine image registration Each transform must define the following (fast, nogil) methods: 

Module: align.vector_fields
#

Computes the composition of two 2D displacement fields 

Computes the composition of two 3D displacement fields 

Create a binary 2D image where pixel values are 1 iff their distance to the center of the image is less than or equal to radius. 

Creates an invertible 2D displacement field 

Creates an invertible 3D displacement field 

Creates a random 2D displacement 'exactly' mapping points of two grids 

Creates a random 3D displacement 'exactly' mapping points of two grids 

Create a binary 3D image where voxel values are 1 iff their distance to the center of the image is less than or equal to radius. 
Downsamples the 2D input vector field by a factor of 2 Downsamples the input vector field by a factor of 2. 

Downsamples the input 3D vector field by a factor of 2 


Downsamples the input 2D image by a factor of 2 

Downsamples the input volume by a factor of 2 

Gradient of an image in physical space 

Computes the inverse of a 2D displacement fields 

Computes the inverse of a 3D displacement fields 



Linearly transforms all vectors of a 2D displacement field 

Linearly transforms all vectors of a 3D displacement field 

Resamples a 2D vector field to a custom target shape 

Resamples a 3D vector field to a custom target shape 

Simplifies a nonlinear warping function combined with an affine transform 

Simplifies a nonlinear warping function combined with an affine transform 

Gradient of an image in physical space 

Transforms a 2D image by an affine transform with bilinear interp. 

Transforms a 2D image by an affine transform with NN interpolation 

Transforms a 3D volume by an affine transform with trilinear interp. 

Transforms a 3D volume by an affine transform with NN interpolation 

Warps a 2D image using bilinear interpolation 

Warps a 2D image using nearest neighbor interpolation 

Warps a 3D volume using trilinear interpolation 

Warps a 3D volume using using nearestneighbor interpolation 

Parameters points : array, shape (n, 2) d1 : array, shape (S, R, C, 2) in2world : array, shape (3, 3) world2out : array, shape (3, 3) field_world2grid : array, shape (3, 3) 

Parameters points : array, shape (n, 3) d1 : array, shape (S, R, C, 3) in2world : array, shape (4, 4) world2out : array, shape (4, 4) field_world2grid : array, shape (4, 4) 
Bunch
#
VerbosityLevels#
 dipy.align.VerbosityLevels()#
VerbosityLevels This enum defines the four levels of verbosity we use in the align module. NONE : do not print anything STATUS : print information about the current status of the algorithm DIAGNOSE : print high level information of the components involved in the registration that can be used to detect a failing component. DEBUG : print as much information as possible to isolate the cause of a bug.
syn_registration#
 dipy.align._public.syn_registration(moving, static, moving_affine=None, static_affine=None, step_length=0.25, metric='CC', dim=3, level_iters=None, prealign=None, **metric_kwargs)#
Register a 2D/3D source image (moving) to a 2D/3D target image (static).
Parameters#
 moving, staticarray or nib.Nifti1Image or str.
Either as a 2D/3D array or as a nifti image object, or as a string containing the full path to a nifti file.
 moving_affine, static_affine4x4 array, optional.
Must be provided for data provided as an array. If provided together with Nifti1Image or str data, this input will override the affine that is stored in the data input. Default: use the affine stored in data.
 metricstring, optional
The metric to be optimized. One of CC, EM, SSD, Default: ‘CC’ => CCMetric.
 dim: int (either 2 or 3), optional
The dimensions of the image domain. Default: 3
 level_iterslist of int, optional
the number of iterations at each level of the Gaussian Pyramid (the length of the list defines the number of pyramid levels to be used). Default: [10, 10, 5].
 metric_kwargsdict, optional
Parameters for initialization of the metric object. If not provided, uses the default settings of each metric.
Returns#
 warped_movingndarray
The data in moving, warped towards the static data.
 forwardndarray (…, 3)
The vector field describing the forward warping from the source to the target.
 backwardndarray (…, 3)
The vector field describing the backward warping from the target to the source.
register_dwi_to_template#
 dipy.align._public.register_dwi_to_template(dwi, gtab, dwi_affine=None, template=None, template_affine=None, reg_method='syn', **reg_kwargs)#
Register DWI data to a template through the B0 volumes.
Parameters#
 dwi4D array, nifti image or str
Containing the DWI data, or full path to a nifti file with DWI.
 gtabGradientTable or sequence of strings
The gradients associated with the DWI data, or a sequence with (fbval, fbvec), full paths to bvals and bvecs files.
 dwi_affine4x4 array, optional
An affine transformation associated with the DWI. Required if data is provided as an array. If provided together with nifti/path, will override the affine that is in the nifti.
 template3D array, nifti image or str
Containing the data for the template, or full path to a nifti file with the template data.
 template_affine4x4 array, optional
An affine transformation associated with the template. Required if data is provided as an array. If provided together with nifti/path, will override the affine that is in the nifti.
 reg_methodstr,
One of “syn” or “aff”, which designates which registration method is used. Either syn, which uses the
syn_registration()
function oraffine_registration()
function. Default: “syn”. reg_kwargskeyword arguments for
syn_registration()
or
Returns#
warped_b0, mapping: The fist is an array with the b0 volume warped to the template. If reg_method is “syn”, the second is a DiffeomorphicMap class instance that can be used to transform between the two spaces. Otherwise, if reg_method is “aff”, this is a 4x4 matrix encoding the affine transform.
Notes#
This function assumes that the DWI data is already internally registered. See
register_dwi_series()
.
write_mapping#
 dipy.align._public.write_mapping(mapping, fname)#
Write out a syn registration mapping to a nifti file.
Parameters#
mapping : a DiffeomorphicMap object derived from
syn_registration()
fname : strFull path to the nifti file storing the mapping
Notes#
The data in the file is organized with shape (X, Y, Z, 3, 2), such that the forward mapping in each voxel is in data[i, j, k, :, 0] and the backward mapping in each voxel is in data[i, j, k, :, 1].
read_mapping#
 dipy.align._public.read_mapping(disp, domain_img, codomain_img, prealign=None)#
Read a syn registration mapping from a nifti file.
Parameters#
 dispstr or Nifti1Image
A file of image containing the mapping displacement field in each voxel Shape (x, y, z, 3, 2)
domain_img : str or Nifti1Image
codomain_img : str or Nifti1Image
Returns#
A
DiffeomorphicMap
object.Notes#
See
write_mapping()
for the data format expected.
resample#
 dipy.align._public.resample(moving, static, moving_affine=None, static_affine=None, between_affine=None)#
Resample an image (moving) from one space to another (static).
Parameters#
 movingarray, nifti image or str
Containing the data for the moving object, or full path to a nifti file with the moving data.
 moving_affine4x4 array, optional
An affine transformation associated with the moving object. Required if data is provided as an array. If provided together with nifti/path, will override the affine that is in the nifti.
 staticarray, nifti image or str
Containing the data for the static object, or full path to a nifti file with the moving data.
 static_affine4x4 array, optional
An affine transformation associated with the static object. Required if data is provided as an array. If provided together with nifti/path, will override the affine that is in the nifti.
 between_affine: 4x4 array, optional
If an additional affine is needed between the two spaces. Default: identity (no additional registration).
Returns#
A Nifti1Image class instance with the data from the moving object resampled into the space of the static object.
affine_registration#
 dipy.align._public.affine_registration(moving, static, moving_affine=None, static_affine=None, pipeline=None, starting_affine=None, metric='MI', level_iters=None, sigmas=None, factors=None, ret_metric=False, moving_mask=None, static_mask=None, **metric_kwargs)#
Find the affine transformation between two 3D images. Alternatively, find the combination of several linear transformations.
Parameters#
 movingarray, nifti image or str
Containing the data for the moving object, or full path to a nifti file with the moving data.
 staticarray, nifti image or str
Containing the data for the static object, or full path to a nifti file with the moving data.
 moving_affine4x4 array, optional
An affine transformation associated with the moving object. Required if data is provided as an array. If provided together with nifti/path, will override the affine that is in the nifti.
 static_affine4x4 array, optional
An affine transformation associated with the static object. Required if data is provided as an array. If provided together with nifti/path, will override the affine that is in the nifti.
 pipelinelist of str, optional
Sequence of transforms to use in the gradual fitting. Default: gradual fit of the full affine (executed from left to right):
["center_of_mass", "translation", "rigid", "affine"]
Alternatively, any other combination of the following registration methods might be used: center_of_mass, translation, rigid, rigid_isoscaling, rigid_scaling and affine. starting_affine: 4x4 array, optional
Initial guess for the transformation between the spaces. Default: identity.
 metricstr, optional.
Currently only supports ‘MI’ for MutualInformationMetric.
 level_iterssequence, optional
AffineRegistration keyword argument: the number of iterations at each scale of the scale space. level_iters[0] corresponds to the coarsest scale, level_iters[1] the finest, where n is the length of the sequence. By default, a 3level scale space with iterations sequence equal to [10000, 1000, 100] will be used.
 sigmassequence of floats, optional
AffineRegistration keyword argument: custom smoothing parameter to build the scale space (one parameter for each scale). By default, the sequence of sigmas will be [3, 1, 0].
 factorssequence of floats, optional
AffineRegistration keyword argument: custom scale factors to build the scale space (one factor for each scale). By default, the sequence of factors will be [4, 2, 1].
 ret_metricboolean, optional
Set it to True to return the value of the optimized coefficients and the optimization quality metric.
 moving_maskarray, shape (S’, R’, C’) or (R’, C’), optional
moving image mask that defines which pixels in the moving image are used to calculate the mutual information.
 static_maskarray, shape (S, R, C) or (R, C), optional
static image mask that defines which pixels in the static image are used to calculate the mutual information.
 nbinsint, optional
MutualInformationMetric keyword argument: the number of bins to be used for computing the intensity histograms. The default is 32.
 sampling_proportionNone or float in interval (0, 1], optional
MutualInformationMetric keyword argument: There are two types of sampling: dense and sparse. Dense sampling uses all voxels for estimating the (joint and marginal) intensity histograms, while sparse sampling uses a subset of them. If sampling_proportion is None, then dense sampling is used. If sampling_proportion is a floating point value in (0,1] then sparse sampling is used, where sampling_proportion specifies the proportion of voxels to be used. The default is None (dense sampling).
Returns#
transformed : array with moving data resampled to the static space after computing the affine transformation affine : the affine 4x4 associated with the transformation. xopt : the value of the optimized coefficients. fopt : the value of the optimization quality metric.
Notes#
Performs a gradual registration between the two inputs, using a pipeline that gradually approximates the final registration. If the final default step (affine) is omitted, the resulting affine may not have all 12 degrees of freedom adjusted.
center_of_mass#
 dipy.align._public.center_of_mass(moving, static, moving_affine=None, static_affine=None, *, pipeline=['center_of_mass'], starting_affine=None, metric='MI', level_iters=None, sigmas=None, factors=None, ret_metric=False, moving_mask=None, static_mask=None, **metric_kwargs)#
Implements a center of mass transform. Based on affine_registration().
translation#
 dipy.align._public.translation(moving, static, moving_affine=None, static_affine=None, *, pipeline=['translation'], starting_affine=None, metric='MI', level_iters=None, sigmas=None, factors=None, ret_metric=False, moving_mask=None, static_mask=None, **metric_kwargs)#
Implements a translation transform. Based on affine_registration().
rigid#
 dipy.align._public.rigid(moving, static, moving_affine=None, static_affine=None, *, pipeline=['rigid'], starting_affine=None, metric='MI', level_iters=None, sigmas=None, factors=None, ret_metric=False, moving_mask=None, static_mask=None, **metric_kwargs)#
Implements a rigid transform. Based on affine_registration().
rigid_isoscaling#
 dipy.align._public.rigid_isoscaling(moving, static, moving_affine=None, static_affine=None, *, pipeline=['rigid_isoscaling'], starting_affine=None, metric='MI', level_iters=None, sigmas=None, factors=None, ret_metric=False, moving_mask=None, static_mask=None, **metric_kwargs)#
Implements a rigid isoscaling transform. Based on affine_registration().
rigid_scaling#
 dipy.align._public.rigid_scaling(moving, static, moving_affine=None, static_affine=None, *, pipeline=['rigid_scaling'], starting_affine=None, metric='MI', level_iters=None, sigmas=None, factors=None, ret_metric=False, moving_mask=None, static_mask=None, **metric_kwargs)#
Implements a rigid scaling transform. Based on affine_registration().
affine#
 dipy.align._public.affine(moving, static, moving_affine=None, static_affine=None, *, pipeline=['affine'], starting_affine=None, metric='MI', level_iters=None, sigmas=None, factors=None, ret_metric=False, moving_mask=None, static_mask=None, **metric_kwargs)#
Implements an affine transform. Based on affine_registration().
_METHOD_DICT#
 dipy.align._public._METHOD_DICT()#
dict() > new empty dictionary dict(mapping) > new dictionary initialized from a mapping object’s
(key, value) pairs
 dict(iterable) > new dictionary initialized as if via:
d = {} for k, v in iterable:
d[k] = v
 dict(**kwargs) > new dictionary initialized with the name=value pairs
in the keyword argument list. For example: dict(one=1, two=2)
register_series#
 dipy.align._public.register_series(series, ref, pipeline=None, series_affine=None, ref_affine=None, static_mask=None)#
Register a series to a reference image.
Parameters#
 series4D array or nib.Nifti1Image class instance or str
The data is 4D with the last dimension separating different 3D volumes
 refint or 3D array or nib.Nifti1Image class instance or str
If this is an int, this is the index of the reference image within the series. Otherwise it is an array of data to register with (associated with a ref_affine required) or a nifti img or full path to a file containing one.
 pipelinesequence, optional
Sequence of transforms to do for each volume in the series. Default: (executed from left to right): [center_of_mass, translation, rigid, affine]
 series_affine, ref_affine4x4 arrays, optional.
The affine. If provided, this input will override the affine provided together with the nifti img or file.
 static_maskarray, shape (S, R, C) or (R, C), optional
static image mask that defines which pixels in the static image are used to calculate the mutual information.
Returns#
xformed, affines : 4D array with transformed data and a (4,4,n) array with 4x4 matrices associated with each of the volumes of the input moving data that was used to transform it into register with the static data.
register_dwi_series#
 dipy.align._public.register_dwi_series(data, gtab, affine=None, b0_ref=0, pipeline=None, static_mask=None)#
Register a DWI series to the mean of the B0 images in that series.
all first registered to the first B0 volume
Parameters#
 data4D array or nibabel Nifti1Image class instance or str
Diffusion data. Either as a 4D array or as a nifti image object, or as a string containing the full path to a nifti file.
 gtaba GradientTable class instance or tuple of strings
If provided as a tuple of strings, these are assumed to be full paths to the bvals and bvecs files (in that order).
 affine4x4 array, optional.
Must be provided for data provided as an array. If provided together with Nifti1Image or str data, this input will override the affine that is stored in the data input. Default: use the affine stored in data.
 b0_refint, optional.
Which b0 volume to use as reference. Default: 0
 pipelinelist of callables, optional.
The transformations to perform in sequence (from left to right): Default:
[center_of_mass, translation, rigid, affine]
 static_maskarray, shape (S, R, C) or (R, C), optional
static image mask that defines which pixels in the static image are used to calculate the mutual information.
Returns#
xform_img, affine_array: a Nifti1Image containing the registered data and using the affine of the original data and a list containing the affine transforms associated with each of the
motion_correction#
 dipy.align._public.motion_correction(data, gtab, affine=None, b0_ref=0, *, pipeline=['center_of_mass', 'translation', 'rigid', 'affine'], static_mask=None)#
Apply a motion correction to a DWI dataset (BetweenVolumes Motion correction)
Parameters#
 data4D array or nibabel Nifti1Image class instance or str
Diffusion data. Either as a 4D array or as a nifti image object, or as a string containing the full path to a nifti file.
 gtaba GradientTable class instance or tuple of strings
If provided as a tuple of strings, these are assumed to be full paths to the bvals and bvecs files (in that order).
 affine4x4 array, optional.
Must be provided for data provided as an array. If provided together with Nifti1Image or str data, this input will override the affine that is stored in the data input. Default: use the affine stored in data.
 b0_refint, optional.
Which b0 volume to use as reference. Default: 0
 pipelinelist of callables, optional.
The transformations to perform in sequence (from left to right): Default:
[center_of_mass, translation, rigid, affine]
 static_maskarray, shape (S, R, C) or (R, C), optional
static image mask that defines which pixels in the static image are used to calculate the mutual information.
Returns#
xform_img, affine_array: a Nifti1Image containing the registered data and using the affine of the original data and a list containing the affine transforms associated with each of the
streamline_registration#
 dipy.align._public.streamline_registration(moving, static, n_points=100, native_resampled=False)#
Register two collections of streamlines (‘bundles’) to each other.
Parameters#
 moving, staticlists of 3 by n, or str
The two bundles to be registered. Given either as lists of arrays with 3D coordinates, or strings containing full paths to these files.
 n_pointsint, optional
How many points to resample to. Default: 100.
 native_resampledbool, optional
Whether to return the moving bundle in the original space, but resampled in the static space to n_points.
Returns#
 alignedlist
Streamlines from the moving group, moved to be closely matched to the static group.
 matrixarray (4, 4)
The affine transformation that takes us from ‘moving’ to ‘static’
determine_num_threads#
 dipy.align.bundlemin.determine_num_threads(num_threads)#
Determine the effective number of threads to be used for OpenMP calls
For
num_threads = None
,  if theOMP_NUM_THREADS
environment variable is set, return that value  otherwise, return the maximum number of cores retrieved byopenmp.opm_get_num_procs()
.For
num_threads > 0
, return this value.For
num_threads < 0
, return the maximal number of threads minusnum_threads + 1
. In particularnum_threads = 1
will use as many threads as there are available cores on the machine.For
num_threads = 0
a ValueError is raised.
Parameters#
 num_threadsint or None
Desired number of threads to be used.
distance_matrix_mdf#
 dipy.align.bundlemin.distance_matrix_mdf(streamlines_a, streamlines_b)#
Minimum direct flipped distance matrix between two streamline sets
All streamlines need to have the same number of points
Parameters#
 streamlines_asequence
of streamlines as arrays, [(N, 3) .. (N, 3)]
 streamlines_bsequence
of streamlines as arrays, [(N, 3) .. (N, 3)]
Returns#
 DMarray, shape (len(streamlines_a), len(streamlines_b))
distance matrix
DeformableRegistration
#
 class dipy.align.cpd.DeformableRegistration(X, Y, sigma2=None, alpha=None, beta=None, low_rank=False, num_eig=100, max_iterations=None, tolerance=None, w=None, *args, **kwargs)#
Bases:
object
Deformable point cloud registration.
Attributes#
 X: numpy array
NxD array of target points.
 Y: numpy array
MxD array of source points.
 TY: numpy array
MxD array of transformed source points.
 sigma2: float (positive)
Initial variance of the Gaussian mixture model.
 N: int
Number of target points.
 M: int
Number of source points.
 D: int
Dimensionality of source and target points
 iteration: int
The current iteration throughout registration.
 max_iterations: int
Registration will terminate once the algorithm has taken this many iterations.
 tolerance: float (positive)
Registration will terminate once the difference between consecutive objective function values falls within this tolerance.
 w: float (between 0 and 1)
Contribution of the uniform distribution to account for outliers. Valid values span 0 (inclusive) and 1 (exclusive).
 q: float
The objective function value that represents the misalignment between source and target point clouds.
 diff: float (positive)
The absolute difference between the current and previous objective function values.
 P: numpy array
MxN array of probabilities. P[m, n] represents the probability that the mth source point corresponds to the nth target point.
 Pt1: numpy array
Nx1 column array. Multiplication result between the transpose of P and a column vector of all 1s.
 P1: numpy array
Mx1 column array. Multiplication result between P and a column vector of all 1s.
 Np: float (positive)
The sum of all elements in P.
 alpha: float (positive)
Represents the tradeoff between the goodness of maximum likelihoo fit and regularization.
 beta: float(positive)
Width of the Gaussian kernel.
 low_rank: bool
Whether to use low rank approximation.
 num_eig: int
Number of eigenvectors to use in lowrank calculation.
 __init__(X, Y, sigma2=None, alpha=None, beta=None, low_rank=False, num_eig=100, max_iterations=None, tolerance=None, w=None, *args, **kwargs)#
 expectation()#
Compute the expectation step of the EM algorithm.
 get_registration_parameters()#
Return the current estimate of the deformable transformation parameters.
Returns#
 self.G: numpy array
Gaussian kernel matrix.
 self.W: numpy array
Deformable transformation matrix.
 iterate()#
Perform one iteration of the EM algorithm.
 maximization()#
Compute the maximization step of the EM algorithm.
 register(callback=<function DeformableRegistration.<lambda>>)#
Perform the EM registration.
Parameters#
 callback: function
A function that will be called after each iteration. Can be used to visualize the registration process.
Returns#
 self.TY: numpy array
MxD array of transformed source points.
 registration_parameters:
Returned params dependent on registration method used.
 transform_point_cloud(Y=None)#
Update a point cloud using the new estimate of the deformable transformation.
Parameters#
 Y: numpy array, optional
Array of points to transform  use to predict on new set of points. Best for predicting on new points not used to run initial registration. If None, self.Y used.
Returns#
If Y is None, returns None. Otherwise, returns the transformed Y.
 update_transform()#
Calculate a new estimate of the deformable transformation. See Eq. 22 of https://arxiv.org/pdf/0905.2635.pdf.
 update_variance()#
Update the variance of the mixture model.
This is using the new estimate of the deformable transformation. See the update rule for sigma2 in Eq. 23 of of https://arxiv.org/pdf/0905.2635.pdf.
gaussian_kernel#
 dipy.align.cpd.gaussian_kernel(X, beta, Y=None)#
low_rank_eigen#
 dipy.align.cpd.low_rank_eigen(G, num_eig)#
Calculate num_eig eigenvectors and eigenvalues of gaussian matrix G.
Enables lower dimensional solving.
initialize_sigma2#
lowrankQS#
 dipy.align.cpd.lowrankQS(G, beta, num_eig, eig_fgt=False)#
Calculate eigenvectors and eigenvalues of gaussian matrix G.
!!! This function is a placeholder for implementing the fast gauss transform. It is not yet implemented. !!!
Parameters#
 G: numpy array
Gaussian kernel matrix.
 beta: float
Width of the Gaussian kernel.
 num_eig: int
Number of eigenvectors to use in lowrank calculation of G
 eig_fgt: bool
If True, use fast gauss transform method to speed up.
compute_cc_backward_step_2d#
 dipy.align.crosscorr.compute_cc_backward_step_2d(grad_moving, factors, radius)#
Gradient of the CC Metric w.r.t. the backward transformation
Computes the gradient of the Cross Correlation metric for symmetric registration (SyN) [Avants2008] w.r.t. the displacement associated to the static image (‘forward’ step) as in [Avants2011]
Parameters#
 grad_movingarray, shape (R, C, 2)
the gradient of the moving image
 factorsarray, shape (R, C, 5)
the precomputed cross correlation terms obtained via precompute_cc_factors_2d
Returns#
 outarray, shape (R, C, 2)
the gradient of the cross correlation metric with respect to the displacement associated to the static image
energy : the cross correlation energy (data term) at this iteration
References#
compute_cc_backward_step_3d#
 dipy.align.crosscorr.compute_cc_backward_step_3d(grad_moving, factors, radius)#
Gradient of the CC Metric w.r.t. the backward transformation
Computes the gradient of the Cross Correlation metric for symmetric registration (SyN) [Avants08] w.r.t. the displacement associated to the static volume (‘backward’ step) as in [Avants11]
Parameters#
 grad_movingarray, shape (S, R, C, 3)
the gradient of the moving volume
 factorsarray, shape (S, R, C, 5)
the precomputed cross correlation terms obtained via precompute_cc_factors_3d
 radiusint
the radius of the neighborhood used for the CC metric when computing the factors. The returned vector field will be zero along a boundary of width radius voxels.
Returns#
 outarray, shape (S, R, C, 3)
the gradient of the cross correlation metric with respect to the displacement associated to the static volume
energy : the cross correlation energy (data term) at this iteration
References#
 [Avants08] Avants, B. B., Epstein, C. L., Grossman, M., & Gee, J. C. (2008)
Symmetric Diffeomorphic Image Registration with CrossCorrelation: Evaluating Automated Labeling of Elderly and Neurodegenerative Brain, Med Image Anal. 12(1), 2641.
 [Avants11] Avants, B. B., Tustison, N., & Song, G. (2011).
Advanced Normalization Tools (ANTS), 135.
compute_cc_forward_step_2d#
 dipy.align.crosscorr.compute_cc_forward_step_2d(grad_static, factors, radius)#
Gradient of the CC Metric w.r.t. the forward transformation
Computes the gradient of the Cross Correlation metric for symmetric registration (SyN) [Avants2008] w.r.t. the displacement associated to the moving image (‘backward’ step) as in [Avants2011]
Parameters#
 grad_staticarray, shape (R, C, 2)
the gradient of the static image
 factorsarray, shape (R, C, 5)
the precomputed cross correlation terms obtained via precompute_cc_factors_2d
Returns#
 outarray, shape (R, C, 2)
the gradient of the cross correlation metric with respect to the displacement associated to the moving image
energy : the cross correlation energy (data term) at this iteration
Notes#
Currently, the gradient of the static image is not being used, but some authors suggest that symmetrizing the gradient by including both, the moving and static gradients may improve the registration quality. We are leaving this parameter as a placeholder for future investigation
References#
compute_cc_forward_step_3d#
 dipy.align.crosscorr.compute_cc_forward_step_3d(grad_static, factors, radius)#
Gradient of the CC Metric w.r.t. the forward transformation
Computes the gradient of the Cross Correlation metric for symmetric registration (SyN) [Avants2008] w.r.t. the displacement associated to the moving volume (‘forward’ step) as in [Avants2011]
Parameters#
 grad_staticarray, shape (S, R, C, 3)
the gradient of the static volume
 factorsarray, shape (S, R, C, 5)
the precomputed cross correlation terms obtained via precompute_cc_factors_3d
 radiusint
the radius of the neighborhood used for the CC metric when computing the factors. The returned vector field will be zero along a boundary of width radius voxels.
Returns#
 outarray, shape (S, R, C, 3)
the gradient of the cross correlation metric with respect to the displacement associated to the moving volume
energy : the cross correlation energy (data term) at this iteration
References#
precompute_cc_factors_2d#
 dipy.align.crosscorr.precompute_cc_factors_2d(static, moving, radius)#
Precomputations to quickly compute the gradient of the CC Metric
Precomputes the separate terms of the cross correlation metric [Avants2008] and image norms at each voxel considering a neighborhood of the given radius to efficiently [Avants2011] compute the gradient of the metric with respect to the deformation field.
Parameters#
 staticarray, shape (R, C)
the static volume, which also defines the reference registration domain
 movingarray, shape (R, C)
the moving volume (notice that both images must already be in a common reference domain, i.e. the same R, C)
radius : the radius of the neighborhood(square of (2*radius + 1)^2 voxels)
Returns#
 factorsarray, shape (R, C, 5)
the precomputed cross correlation terms: factors[:,:,0] : static minus its mean value along the neighborhood factors[:,:,1] : moving minus its mean value along the neighborhood factors[:,:,2] : sum of the pointwise products of static and moving
along the neighborhood
factors[:,:,3] : sum of sq. values of static along the neighborhood factors[:,:,4] : sum of sq. values of moving along the neighborhood
References#
precompute_cc_factors_2d_test#
 dipy.align.crosscorr.precompute_cc_factors_2d_test(static, moving, radius)#
Precomputations to quickly compute the gradient of the CC Metric
This version of precompute_cc_factors_2d is for testing purposes, it directly computes the local crosscorrelation without any optimization.
precompute_cc_factors_3d#
 dipy.align.crosscorr.precompute_cc_factors_3d(static, moving, radius, num_threads=None)#
Precomputations to quickly compute the gradient of the CC Metric
Precomputes the separate terms of the cross correlation metric and image norms at each voxel considering a neighborhood of the given radius to efficiently compute the gradient of the metric with respect to the deformation field [Ocegueda2016] [Avants2008] [Avants2011].
Parameters#
 staticarray, shape (S, R, C)
the static volume, which also defines the reference registration domain
 movingarray, shape (S, R, C)
the moving volume (notice that both images must already be in a common reference domain, i.e. the same S, R, C)
radius : the radius of the neighborhood (cube of (2 * radius + 1)^3 voxels)
Returns#
 factorsarray, shape (S, R, C, 5)
the precomputed cross correlation terms: factors[:,:,:,0] : static minus its mean value along the neighborhood factors[:,:,:,1] : moving minus its mean value along the neighborhood factors[:,:,:,2] : sum of the pointwise products of static and moving
along the neighborhood
factors[:,:,:,3] : sum of sq. values of static along the neighborhood factors[:,:,:,4] : sum of sq. values of moving along the neighborhood
References#
precompute_cc_factors_3d_test#
 dipy.align.crosscorr.precompute_cc_factors_3d_test(static, moving, radius)#
Precomputations to quickly compute the gradient of the CC Metric
This version of precompute_cc_factors_3d is for testing purposes, it directly computes the local crosscorrelation factors without any optimization, so it is less errorprone than the accelerated version.
compute_em_demons_step_2d#
 dipy.align.expectmax.compute_em_demons_step_2d(delta_field, sigma_sq_field, gradient_moving, sigma_sq_x, out)#
Demons step for EM metric in 2D Computes the demons step [Vercauteren09] for SSDdriven registration ( eq. 4 in [Vercauteren09] ) using the EM algorithm [Arce14] to handle multimodality images. In this case, \(\sigma_i\) in eq. 4 of [Vercauteren] is estimated using the EM algorithm, while in the original version of diffeomorphic demons it is estimated by the difference between the image values at each pixel. Parameters ——— delta_field : array, shape (R, C) contains, at each pixel, the difference between the moving image (warped under the current deformation s(. , .) ) J and the static image I: delta_field[i,j] = J(s(i,j))  I(i,j). The order is important, changing to delta_field[i,j] = I(i,j)  J(s(i,j)) yields the backward demons step warping the static image towards the moving, which may not be the intended behavior unless the ‘gradient_moving’ passed corresponds to the gradient of the static image sigma_sq_field : array, shape (R, C) contains, at each pixel (i, j), the estimated variance (not std) of the hidden variable associated to the intensity at static[i,j] (which must have been previously quantized) gradient_moving : array, shape (R, C, 2) the gradient of the moving image sigma_sq_x : float parameter controlling the amount of regularization. It corresponds to \(\sigma_x^2\) in algorithm 1 of Vercauteren et al.[2] out : array, shape (R, C, 2) the resulting demons step will be written to this array Returns —— demons_step : array, shape (R, C, 2) the demons step to be applied for updating the current displacement field energy : float the current em energy (before applying the returned demons_step) References ——— [Arce14] Arcesantana, E., Camposdelgado, D. U., & Viguerasg, F. (2014). Nonrigid Multimodal Image Registration Based on the ExpectationMaximization Algorithm, (168140), 3647. [Vercauteren09] Vercauteren, T., Pennec, X., Perchant, A., & Ayache, N. (2009). Diffeomorphic demons: efficient nonparametric image registration. NeuroImage, 45(1 Suppl), S6172. doi:10.1016/j.neuroimage.2008.10.040
compute_em_demons_step_3d#
 dipy.align.expectmax.compute_em_demons_step_3d(delta_field, sigma_sq_field, gradient_moving, sigma_sq_x, out)#
Demons step for EM metric in 3D Computes the demons step [Vercauteren09] for SSDdriven registration ( eq. 4 in [Vercauteren09] ) using the EM algorithm [Arce14] to handle multimodality images. In this case, \(\sigma_i\) in eq. 4 of [Vercauteren09] is estimated using the EM algorithm, while in the original version of diffeomorphic demons it is estimated by the difference between the image values at each pixel. Parameters ——— delta_field : array, shape (S, R, C) contains, at each pixel, the difference between the moving image (warped under the current deformation s ) J and the static image I: delta_field[k,i,j] = J(s(k,i,j))  I(k,i,j). The order is important, changing to delta_field[k,i,j] = I(k,i,j)  J(s(k,i,j)) yields the backward demons step warping the static image towards the moving, which may not be the intended behavior unless the ‘gradient_moving’ passed corresponds to the gradient of the static image sigma_sq_field : array, shape (S, R, C) contains, at each pixel (k, i, j), the estimated variance (not std) of the hidden variable associated to the intensity at static[k,i,j] (which must have been previously quantized) gradient_moving : array, shape (S, R, C, 2) the gradient of the moving image sigma_sq_x : float parameter controlling the amount of regularization. It corresponds to \(\sigma_x^2\) in algorithm 1 of Vercauteren et al.[2]. out : array, shape (S, R, C, 2) the resulting demons step will be written to this array Returns —— demons_step : array, shape (S, R, C, 3) the demons step to be applied for updating the current displacement field energy : float the current em energy (before applying the returned demons_step) References ——— [Arce14] Arcesantana, E., Camposdelgado, D. U., & Viguerasg, F. (2014). Nonrigid Multimodal Image Registration Based on the ExpectationMaximization Algorithm, (168140), 3647. [Vercauteren09] Vercauteren, T., Pennec, X., Perchant, A., & Ayache, N. (2009). Diffeomorphic demons: efficient nonparametric image registration. NeuroImage, 45(1 Suppl), S6172. doi:10.1016/j.neuroimage.2008.10.040
compute_masked_class_stats_2d#
 dipy.align.expectmax.compute_masked_class_stats_2d(mask, v, num_labels, labels)#
Computes the mean and std. for each quantization level.
Computes the mean and standard deviation of the intensities in ‘v’ for each corresponding label in ‘labels’. In other words, for each label L, it computes the mean and standard deviation of the intensities in ‘v’ at pixels whose label in ‘labels’ is L. This is used by the EM metric to compute statistics for each hidden variable represented by the labels.
Parameters#
 maskarray, shape (R, C)
the mask of pixels that will be taken into account for computing the statistics. All zero pixels in mask will be ignored
 varray, shape (R, C)
the image which the statistics will be computed from
 num_labelsint
the number of different labels in ‘labels’ (equal to the number of hidden variables in the EM metric)
 labelsarray, shape (R, C)
the label assigned to each pixel
Returns#
 meansarray, shape (num_labels,)
means[i], 0<=i<num_labels will be the mean intensity in v of all voxels labeled i, or 0 if no voxels are labeled i
 variancesarray, shape (num_labels,)
variances[i], 0<=i<num_labels will be the standard deviation of the intensities in v of all voxels labeled i, or infinite if less than 2 voxels are labeled i.
compute_masked_class_stats_3d#
 dipy.align.expectmax.compute_masked_class_stats_3d(mask, v, num_labels, labels)#
Computes the mean and std. for each quantization level.
Computes the mean and standard deviation of the intensities in ‘v’ for each corresponding label in ‘labels’. In other words, for each label L, it computes the mean and standard deviation of the intensities in ‘v’ at voxels whose label in ‘labels’ is L. This is used by the EM metric to compute statistics for each hidden variable represented by the labels.
Parameters#
 maskarray, shape (S, R, C)
the mask of voxels that will be taken into account for computing the statistics. All zero voxels in mask will be ignored
 varray, shape (S, R, C)
the volume which the statistics will be computed from
 num_labelsint
the number of different labels in ‘labels’ (equal to the number of hidden variables in the EM metric)
 labelsarray, shape (S, R, C)
the label assigned to each pixel
Returns#
 meansarray, shape (num_labels,)
means[i], 0<=i<num_labels will be the mean intensity in v of all voxels labeled i, or 0 if no voxels are labeled i
 variancesarray, shape (num_labels,)
variances[i], 0<=i<num_labels will be the standard deviation of the intensities in v of all voxels labeled i, or infinite if less than 2 voxels are labeled i.
quantize_positive_2d#
 dipy.align.expectmax.quantize_positive_2d(v, num_levels)#
Quantizes a 2D image to num_levels quantization levels
Quantizes the input image at num_levels intensity levels considering <=0 as a special value. Those input pixels <=0, and only those, will be assigned a quantization level of 0. The positive values are divided into the remaining num_levels1 uniform quantization levels.
The following are undefined, and raise a ValueError: * Quantizing at zero levels because at least one level must be assigned * Quantizing at one level because positive values should be assigned a
level different from the secial level 0 (at least 2 levels are needed)
Parameters#
 varray, shape (R, C)
the image to be quantized
 num_levelsint
the number of levels
Returns#
 outarray, shape (R, C), same shape as v
the quantized image
 levels: array, shape (num_levels,)
the quantization values: levels[0]=0, and levels[i] is the midpoint of the interval of intensities that are assigned to quantization level i, i=1, …, num_levels1.
 hist: array, shape (num_levels,)
histogram: the number of pixels that were assigned to each quantization level
quantize_positive_3d#
 dipy.align.expectmax.quantize_positive_3d(v, num_levels)#
Quantizes a 3D volume to num_levels quantization levels
Quantizes the input volume at num_levels intensity levels considering <=0 as a special value. Those input voxels <=0, and only those, will be assigned a quantization level of 0. The positive values are divided into the remaining num_levels1 uniform quantization levels.
The following are undefined, and raise a ValueError: * Quantizing at zero levels because at least one level must be assigned * Quantizing at one level because positive values should be assigned a
level different from the secial level 0 (at least 2 levels are needed)
Parameters#
 varray, shape (S, R, C)
the volume to be quantized
 num_levelsint
the number of levels
Returns#
 outarray, shape (S, R, C), same shape as v
the quantized volume
 levels: array, shape (num_levels,)
the quantization values: levels[0]=0, and levels[i] is the midpoint of the interval of intensities that are assigned to quantization level i, i=1, …, num_levels1.
 hist: array, shape (num_levels,)
histogram: the number of voxels that were assigned to each quantization level
AffineInversionError
#
AffineInvalidValuesError
#
AffineMap
#
 class dipy.align.imaffine.AffineMap(affine, domain_grid_shape=None, domain_grid2world=None, codomain_grid_shape=None, codomain_grid2world=None)#
Bases:
object
 __init__(affine, domain_grid_shape=None, domain_grid2world=None, codomain_grid_shape=None, codomain_grid2world=None)#
AffineMap.
Implements an affine transformation whose domain is given by domain_grid and domain_grid2world, and whose codomain is given by codomain_grid and codomain_grid2world.
The actual transform is represented by the affine matrix, which operate in world coordinates. Therefore, to transform a moving image towards a static image, we first map each voxel (i,j,k) of the static image to world coordinates (x,y,z) by applying domain_grid2world. Then we apply the affine transform to (x,y,z) obtaining (x’, y’, z’) in moving image’s world coordinates. Finally, (x’, y’, z’) is mapped to voxel coordinates (i’, j’, k’) in the moving image by multiplying (x’, y’, z’) by the inverse of codomain_grid2world. The codomain_grid_shape is used analogously to transform the static image towards the moving image when calling transform_inverse.
If the domain/codomain information is not provided (None) then the sampling information needs to be specified each time the transform or transform_inverse is called to transform images. Note that such sampling information is not necessary to transform points defined in physical space, such as stream lines.
Parameters#
 affinearray, shape (dim + 1, dim + 1)
the matrix defining the affine transform, where dim is the dimension of the space this map operates in (2 for 2D images, 3 for 3D images). If None, then self represents the identity transformation.
 domain_grid_shapesequence, shape (dim,), optional
the shape of the default domain sampling grid. When transform is called to transform an image, the resulting image will have this shape, unless a different sampling information is provided. If None, then the sampling grid shape must be specified each time the transform method is called.
 domain_grid2worldarray, shape (dim + 1, dim + 1), optional
the gridtoworld transform associated with the domain grid. If None (the default), then the gridtoworld transform is assumed to be the identity.
 codomain_grid_shapesequence of integers, shape (dim,)
the shape of the default codomain sampling grid. When transform_inverse is called to transform an image, the resulting image will have this shape, unless a different sampling information is provided. If None (the default), then the sampling grid shape must be specified each time the transform_inverse method is called.
 codomain_grid2worldarray, shape (dim + 1, dim + 1)
the gridtoworld transform associated with the codomain grid. If None (the default), then the gridtoworld transform is assumed to be the identity.
 get_affine()#
Return the value of the transformation, not a reference.
Returns#
 affinendarray
Copy of the transform, not a reference.
 set_affine(affine)#
Set the affine transform (operating in physical space).
Also sets self.affine_inv  the inverse of affine, or None if there is no inverse.
Parameters#
 affinearray, shape (dim + 1, dim + 1)
the matrix representing the affine transform operating in physical space. The domain and codomain information remains unchanged. If None, then self represents the identity transformation.
 transform(image, interpolation='linear', image_grid2world=None, sampling_grid_shape=None, sampling_grid2world=None, resample_only=False)#
Transform the input image from codomain to domain space.
By default, the transformed image is sampled at a grid defined by self.domain_shape and self.domain_grid2world. If such information was not provided then sampling_grid_shape is mandatory.
Parameters#
 image2D or 3D array
the image to be transformed
 interpolationstring, either ‘linear’ or ‘nearest’
the type of interpolation to be used, either ‘linear’ (for klinear interpolation) or ‘nearest’ for nearest neighbor
 image_grid2worldarray, shape (dim + 1, dim + 1), optional
the gridtoworld transform associated with image. If None (the default), then the gridtoworld transform is assumed to be the identity.
 sampling_grid_shapesequence, shape (dim,), optional
the shape of the grid where the transformed image must be sampled. If None (the default), then self.codomain_shape is used instead (which must have been set at initialization, otherwise an exception will be raised).
 sampling_grid2worldarray, shape (dim + 1, dim + 1), optional
the gridtoworld transform associated with the sampling grid (specified by sampling_grid_shape, or by default self.codomain_shape). If None (the default), then the gridtoworld transform is assumed to be the identity.
 resample_onlyBoolean, optional
If False (the default) the affine transform is applied normally. If True, then the affine transform is not applied, and the input image is just resampled on the domain grid of this transform.
Returns#
 transformedarray, shape sampling_grid_shape or
self.codomain_shape
the transformed image, sampled at the requested grid
 transform_inverse(image, interpolation='linear', image_grid2world=None, sampling_grid_shape=None, sampling_grid2world=None, resample_only=False)#
Transform the input image from domain to codomain space.
By default, the transformed image is sampled at a grid defined by self.codomain_shape and self.codomain_grid2world. If such information was not provided then sampling_grid_shape is mandatory.
Parameters#
 image2D or 3D array
the image to be transformed
 interpolationstring, either ‘linear’ or ‘nearest’
the type of interpolation to be used, either ‘linear’ (for klinear interpolation) or ‘nearest’ for nearest neighbor
 image_grid2worldarray, shape (dim + 1, dim + 1), optional
the gridtoworld transform associated with image. If None (the default), then the gridtoworld transform is assumed to be the identity.
 sampling_grid_shapesequence, shape (dim,), optional
the shape of the grid where the transformed image must be sampled. If None (the default), then self.codomain_shape is used instead (which must have been set at initialization, otherwise an exception will be raised).
 sampling_grid2worldarray, shape (dim + 1, dim + 1), optional
the gridtoworld transform associated with the sampling grid (specified by sampling_grid_shape, or by default self.codomain_shape). If None (the default), then the gridtoworld transform is assumed to be the identity.
 resample_onlyBoolean, optional
If False (the default) the affine transform is applied normally. If True, then the affine transform is not applied, and the input image is just resampled on the domain grid of this transform.
Returns#
 transformedarray, shape sampling_grid_shape or
self.codomain_shape
the transformed image, sampled at the requested grid
MutualInformationMetric
#
 class dipy.align.imaffine.MutualInformationMetric(nbins=32, sampling_proportion=None)#
Bases:
object
 __init__(nbins=32, sampling_proportion=None)#
Initialize an instance of the Mutual Information metric.
This class implements the methods required by Optimizer to drive the registration process.
Parameters#
 nbinsint, optional
the number of bins to be used for computing the intensity histograms. The default is 32.
 sampling_proportionNone or float in interval (0, 1], optional
There are two types of sampling: dense and sparse. Dense sampling uses all voxels for estimating the (joint and marginal) intensity histograms, while sparse sampling uses a subset of them. If sampling_proportion is None, then dense sampling is used. If sampling_proportion is a floating point value in (0,1] then sparse sampling is used, where sampling_proportion specifies the proportion of voxels to be used. The default is None.
Notes#
Since we use linear interpolation, images are not, in general, differentiable at exact voxel coordinates, but they are differentiable between voxel coordinates. When using sparse sampling, selected voxels are slightly moved by adding a small random displacement within one voxel to prevent sampling points from being located exactly at voxel coordinates. When using dense sampling, this random displacement is not applied.
 distance(params)#
Numeric value of the negative Mutual Information.
We need to change the sign so we can use standard minimization algorithms.
Parameters#
 paramsarray, shape (n,)
the parameter vector of the transform currently used by the metric (the transform name is provided when self.setup is called), n is the number of parameters of the transform
Returns#
 neg_mifloat
the negative mutual information of the input images after transforming the moving image by the currently set transform with params parameters
 distance_and_gradient(params)#
Numeric value of the metric and its gradient at given parameters.
Parameters#
 paramsarray, shape (n,)
the parameter vector of the transform currently used by the metric (the transform name is provided when self.setup is called), n is the number of parameters of the transform
Returns#
 neg_mifloat
the negative mutual information of the input images after transforming the moving image by the currently set transform with params parameters
 neg_mi_gradarray, shape (n,)
the gradient of the negative Mutual Information
 gradient(params)#
Numeric value of the metric’s gradient at the given parameters.
Parameters#
 paramsarray, shape (n,)
the parameter vector of the transform currently used by the metric (the transform name is provided when self.setup is called), n is the number of parameters of the transform
Returns#
 gradarray, shape (n,)
the gradient of the negative Mutual Information
 setup(transform, static, moving, static_grid2world=None, moving_grid2world=None, starting_affine=None, static_mask=None, moving_mask=None)#
Prepare the metric to compute intensity densities and gradients.
The histograms will be setup to compute probability densities of intensities within the minimum and maximum values of static and moving
Parameters#
 transform: instance of Transform
the transformation with respect to whose parameters the gradient must be computed
 staticarray, shape (S, R, C) or (R, C)
static image
 movingarray, shape (S’, R’, C’) or (R’, C’)
moving image. The dimensions of the static (S, R, C) and moving (S’, R’, C’) images do not need to be the same.
 static_grid2worldarray (dim+1, dim+1), optional
the gridtospace transform of the static image. The default is None, implying the transform is the identity.
 moving_grid2worldarray (dim+1, dim+1)
the gridtospace transform of the moving image. The default is None, implying the spacing along all axes is 1.
 starting_affinearray, shape (dim+1, dim+1), optional
the prealigning matrix (an affine transform) that roughly aligns the moving image towards the static image. If None, no prealignment is performed. If a prealignment matrix is available, it is recommended to provide this matrix as starting_affine instead of manually transforming the moving image to reduce interpolation artifacts. The default is None, implying no prealignment is performed.
 static_maskarray, shape (S, R, C) or (R, C), optional
static image mask that defines which pixels in the static image are used to calculate the mutual information.
 moving_maskarray, shape (S’, R’, C’) or (R’, C’), optional
moving image mask that defines which pixels in the moving image are used to calculate the mutual information.
AffineRegistration
#
 class dipy.align.imaffine.AffineRegistration(metric=None, level_iters=None, sigmas=None, factors=None, method='LBFGSB', ss_sigma_factor=None, options=None, verbosity=1)#
Bases:
object
 __init__(metric=None, level_iters=None, sigmas=None, factors=None, method='LBFGSB', ss_sigma_factor=None, options=None, verbosity=1)#
Initialize an instance of the AffineRegistration class. Parameters ——— metric : None or object, optional an instance of a metric. The default is None, implying the Mutual Information metric with default settings. level_iters : sequence, optional the number of iterations at each scale of the scale space. level_iters[0] corresponds to the coarsest scale, level_iters[1] the finest, where n is the length of the sequence. By default, a 3level scale space with iterations sequence equal to [10000, 1000, 100] will be used. sigmas : sequence of floats, optional custom smoothing parameter to build the scale space (one parameter for each scale). By default, the sequence of sigmas will be [3, 1, 0]. factors : sequence of floats, optional custom scale factors to build the scale space (one factor for each scale). By default, the sequence of factors will be [4, 2, 1]. method : string, optional optimization method to be used. If Scipy version < 0.12, then only LBFGSB is available. Otherwise, method can be any gradientbased method available in dipy.core.Optimize: CG, BFGS, NewtonCG, dogleg or trustncg. The default is ‘LBFGSB’. ss_sigma_factor : float, optional If None, this parameter is not used and an isotropic scale space with the given factors and sigmas will be built. If not None, an anisotropic scale space will be used by automatically selecting the smoothing sigmas along each axis according to the voxel dimensions of the given image. The ss_sigma_factor is used to scale the automatically computed sigmas. For example, in the isotropic case, the sigma of the kernel will be \(factor * (2 ^ i)\) where \(i = 1, 2, ..., n_scales  1\) is the scale (the finest resolution image \(i=0\) is never smoothed). The default is None. options : dict, optional extra optimization options. The default is None, implying no extra options are passed to the optimizer. verbosity: int (one of {0, 1, 2, 3}), optional Set the verbosity level of the algorithm: 0 : do not print anything 1 : print information about the current status of the algorithm 2 : print high level information of the components involved in the registration that can be used to detect a failing component. 3 : print as much information as possible to isolate the cause of a bug. Default: 1
 docstring_addendum = 'verbosity: int (one of {0, 1, 2, 3}), optional\n Set the verbosity level of the algorithm:\n 0 : do not print anything\n 1 : print information about the current status of the algorithm\n 2 : print high level information of the components involved in\n the registration that can be used to detect a failing\n component.\n 3 : print as much information as possible to isolate the cause\n of a bug.\n Default: 1\n '#
 optimize(static, moving, transform, params0, static_grid2world=None, moving_grid2world=None, starting_affine=None, ret_metric=False, static_mask=None, moving_mask=None)#
Start the optimization process.
Parameters#
 static2D or 3D array
the image to be used as reference during optimization.
 moving2D or 3D array
the image to be used as “moving” during optimization. It is necessary to prealign the moving image to ensure its domain lies inside the domain of the deformation fields. This is assumed to be accomplished by “prealigning” the moving image towards the static using an affine transformation given by the ‘starting_affine’ matrix
 transforminstance of Transform
the transformation with respect to whose parameters the gradient must be computed
 params0array, shape (n,)
parameters from which to start the optimization. If None, the optimization will start at the identity transform. n is the number of parameters of the specified transformation.
 static_grid2worldarray, shape (dim+1, dim+1), optional
the voxeltospace transformation associated with the static image. The default is None, implying the transform is the identity.
 moving_grid2worldarray, shape (dim+1, dim+1), optional
the voxeltospace transformation associated with the moving image. The default is None, implying the transform is the identity.
 starting_affinestring, or matrix, or None, optional
 If string:
‘mass’: align centers of gravity ‘voxelorigin’: align physical coordinates of voxel (0,0,0) ‘centers’: align physical coordinates of central voxels
 If matrix:
array, shape (dim+1, dim+1).
 If None:
Start from identity.
The default is None.
 ret_metricboolean, optional
if True, it returns the parameters for measuring the similarity between the images (default ‘False’). The metric containing optimal parameters and the distance between the images.
 static_maskarray, shape (S, R, C) or (R, C), optional
static image mask that defines which pixels in the static image are used to calculate the mutual information.
 moving_maskarray, shape (S’, R’, C’) or (R’, C’), optional
moving image mask that defines which pixels in the moving image are used to calculate the mutual information.
Returns#
 affine_mapinstance of AffineMap
the affine resulting affine transformation
 xoptoptimal parameters
the optimal parameters (translation, rotation shear etc.)
 foptSimilarity metric
the value of the function at the optimal parameters.
_transform_method#
 dipy.align.imaffine._transform_method()#
dict() > new empty dictionary dict(mapping) > new dictionary initialized from a mapping object’s
(key, value) pairs
 dict(iterable) > new dictionary initialized as if via:
d = {} for k, v in iterable:
d[k] = v
 dict(**kwargs) > new dictionary initialized with the name=value pairs
in the keyword argument list. For example: dict(one=1, two=2)
transform_centers_of_mass#
 dipy.align.imaffine.transform_centers_of_mass(static, static_grid2world, moving, moving_grid2world)#
Transformation to align the center of mass of the input images.
Parameters#
 staticarray, shape (S, R, C)
static image
 static_grid2worldarray, shape (dim+1, dim+1)
the voxeltospace transformation of the static image
 movingarray, shape (S, R, C)
moving image
 moving_grid2worldarray, shape (dim+1, dim+1)
the voxeltospace transformation of the moving image
Returns#
 affine_mapinstance of AffineMap
the affine transformation (translation only, in this case) aligning the center of mass of the moving image towards the one of the static image
transform_geometric_centers#
 dipy.align.imaffine.transform_geometric_centers(static, static_grid2world, moving, moving_grid2world)#
Transformation to align the geometric center of the input images.
With “geometric center” of a volume we mean the physical coordinates of its central voxel
Parameters#
 staticarray, shape (S, R, C)
static image
 static_grid2worldarray, shape (dim+1, dim+1)
the voxeltospace transformation of the static image
 movingarray, shape (S, R, C)
moving image
 moving_grid2worldarray, shape (dim+1, dim+1)
the voxeltospace transformation of the moving image
Returns#
 affine_mapinstance of AffineMap
the affine transformation (translation only, in this case) aligning the geometric center of the moving image towards the one of the static image
transform_origins#
 dipy.align.imaffine.transform_origins(static, static_grid2world, moving, moving_grid2world)#
Transformation to align the origins of the input images.
With “origin” of a volume we mean the physical coordinates of voxel (0,0,0)
Parameters#
 staticarray, shape (S, R, C)
static image
 static_grid2worldarray, shape (dim+1, dim+1)
the voxeltospace transformation of the static image
 movingarray, shape (S, R, C)
moving image
 moving_grid2worldarray, shape (dim+1, dim+1)
the voxeltospace transformation of the moving image
Returns#
 affine_mapinstance of AffineMap
the affine transformation (translation only, in this case) aligning the origin of the moving image towards the one of the static image
DiffeomorphicMap
#
 class dipy.align.imwarp.DiffeomorphicMap(dim, disp_shape, disp_grid2world=None, domain_shape=None, domain_grid2world=None, codomain_shape=None, codomain_grid2world=None, prealign=None)#
Bases:
object
 __init__(dim, disp_shape, disp_grid2world=None, domain_shape=None, domain_grid2world=None, codomain_shape=None, codomain_grid2world=None, prealign=None)#
DiffeomorphicMap
Implements a diffeomorphic transformation on the physical space. The deformation fields encoding the direct and inverse transformations share the same domain discretization (both the discretization grid shape and voxeltospace matrix). The input coordinates (physical coordinates) are first aligned using prealign, and then displaced using the corresponding vector field interpolated at the aligned coordinates.
Parameters#
 dimint, 2 or 3
the transformation’s dimension
 disp_shapearray, shape (dim,)
the number of slices (if 3D), rows and columns of the deformation field’s discretization
 disp_grid2worldthe voxeltospace transform between the def. fields
grid and space
 domain_shapearray, shape (dim,)
the number of slices (if 3D), rows and columns of the default discretization of this map’s domain
 domain_grid2worldarray, shape (dim+1, dim+1)
the default voxeltospace transformation between this map’s discretization and physical space
 codomain_shapearray, shape (dim,)
the number of slices (if 3D), rows and columns of the images that are ‘normally’ warped using this transformation in the forward direction (this will provide default transformation parameters to warp images under this transformation). By default, we assume that the inverse transformation is ‘normally’ used to warp images with the same discretization and voxeltospace transformation as the deformation field grid.
 codomain_grid2worldarray, shape (dim+1, dim+1)
the voxeltospace transformation of images that are ‘normally’ warped using this transformation (in the forward direction).
 prealignarray, shape (dim+1, dim+1)
the linear transformation to be applied to align input images to the reference space before warping under the deformation field.
 allocate()#
Creates a zero displacement field
Creates a zero displacement field (the identity transformation).
 compute_inversion_error()#
Inversion error of the displacement fields
Estimates the inversion error of the displacement fields by computing statistics of the residual vectors obtained after composing the forward and backward displacement fields.
Returns#
 residualarray, shape (R, C) or (S, R, C)
the displacement field resulting from composing the forward and backward displacement fields of this transformation (the residual should be zero for a perfect diffeomorphism)
 statsarray, shape (3,)
statistics from the norms of the vectors of the residual displacement field: maximum, mean and standard deviation
Notes#
Since the forward and backward displacement fields have the same discretization, the final composition is given by
comp[i] = forward[ i + Dinv * backward[i]]
where Dinv is the spacetogrid transformation of the displacement fields
 expand_fields(expand_factors, new_shape)#
Expands the displacement fields from current shape to new_shape
Upsamples the discretization of the displacement fields to be of new_shape shape.
Parameters#
 expand_factorsarray, shape (dim,)
the factors scaling current spacings (voxel sizes) to spacings in the expanded discretization.
 new_shapearray, shape (dim,)
the shape of the arrays holding the upsampled discretization
 get_backward_field()#
Deformation field to transform an image in the backward direction
Returns the deformation field that must be used to warp an image under this transformation in the backward direction (note the ‘is_inverse’ flag).
 get_forward_field()#
Deformation field to transform an image in the forward direction
Returns the deformation field that must be used to warp an image under this transformation in the forward direction (note the ‘is_inverse’ flag).
 get_simplified_transform()#
Constructs a simplified version of this Diffeomorhic Map
The simplified version incorporates the prealign transform, as well as the domain and codomain affine transforms into the displacement field. The resulting transformation may be regarded as operating on the image spaces given by the domain and codomain discretization. As a result, self.prealign, self.disp_grid2world, self.domain_grid2world and self.codomain affine will be None (denoting Identity) in the resulting diffeomorphic map.
 interpret_matrix(obj)#
Try to interpret obj as a matrix
Some operations are performed faster if we know in advance if a matrix is the identity (so we can skip the actual matrixvector multiplication). This function returns None if the given object is None or the ‘identity’ string. It returns the same object if it is a numpy array. It raises an exception otherwise.
Parameters#
 objobject
any object
Returns#
 objobject
the same object given as argument if obj is None or a numpy array. None if obj is the ‘identity’ string.
 inverse()#
Inverse of this DiffeomorphicMap instance
Returns a diffeomorphic map object representing the inverse of this transformation. The internal arrays are not copied but just referenced.
Returns#
 invDiffeomorphicMap object
the inverse of this diffeomorphic map.
 shallow_copy()#
Shallow copy of this DiffeomorphicMap instance
Creates a shallow copy of this diffeomorphic map (the arrays are not copied but just referenced)
Returns#
 new_mapDiffeomorphicMap object
the shallow copy of this diffeomorphic map
 transform(image, interpolation='linear', image_world2grid=None, out_shape=None, out_grid2world=None)#
Warps an image in the forward direction
Transforms the input image under this transformation in the forward direction. It uses the “is_inverse” flag to switch between “forward” and “backward” (if is_inverse is False, then transform(…) warps the image forwards, else it warps the image backwards).
Parameters#
 imagearray, shape (s, r, c) if dim = 3 or (r, c) if dim = 2
the image to be warped under this transformation in the forward direction
 interpolationstring, either ‘linear’ or ‘nearest’
the type of interpolation to be used for warping, either ‘linear’ (for klinear interpolation) or ‘nearest’ for nearest neighbor
 image_world2gridarray, shape (dim+1, dim+1)
the transformation bringing world (space) coordinates to voxel coordinates of the image given as input
 out_shapearray, shape (dim,)
the number of slices, rows and columns of the desired warped image
 out_grid2worldthe transformation bringing voxel coordinates of the
warped image to physical space
Returns#
 warpedarray, shape = out_shape or self.codomain_shape if None
the warped image under this transformation in the forward direction
Notes#
See _warp_forward and _warp_backward documentation for further information.
 transform_inverse(image, interpolation='linear', image_world2grid=None, out_shape=None, out_grid2world=None)#
Warps an image in the backward direction
Transforms the input image under this transformation in the backward direction. It uses the “is_inverse” flag to switch between “forward” and “backward” (if is_inverse is False, then transform_inverse(…) warps the image backwards, else it warps the image forwards)
Parameters#
 imagearray, shape (s, r, c) if dim = 3 or (r, c) if dim = 2
the image to be warped under this transformation in the forward direction
 interpolationstring, either ‘linear’ or ‘nearest’
the type of interpolation to be used for warping, either ‘linear’ (for klinear interpolation) or ‘nearest’ for nearest neighbor
 image_world2gridarray, shape (dim+1, dim+1)
the transformation bringing world (space) coordinates to voxel coordinates of the image given as input
 out_shapearray, shape (dim,)
the number of slices, rows, and columns of the desired warped image
 out_grid2worldthe transformation bringing voxel coordinates of the
warped image to physical space
Returns#
 warpedarray, shape = out_shape or self.codomain_shape if None
warped image under this transformation in the backward direction
Notes#
See _warp_forward and _warp_backward documentation for further information.
 transform_points(points, coord2world=None, world2coord=None)#
Warp the list of points in the forward direction.
Applies this diffeomorphic map to the list of points (or streamlines) given by points. We assume that the points’ coordinates are mapped to world coordinates by applying the coord2world affine transform. The warped coordinates are given in world coordinates unless world2coord affine transform is specified, in which case the warped points will be transformed to the corresponding coordinate system.
Parameters#
points : array, shape (N, dim) or Streamlines object
 coord2worldarray, shape (dim+1, dim+1), optional
affine matrix mapping points to world coordinates
 world2coordarray, shape (dim+1, dim+1), optional
affine matrix mapping world coordinates to points
 transform_points_inverse(points, coord2world=None, world2coord=None)#
Warp the list of points in the backward direction.
Applies this diffeomorphic map to the list of points (or streamlines) given by points. We assume that the points’ coordinates are mapped to world coordinates by applying the coord2world affine transform. The warped coordinates are given in world coordinates unless world2coord affine transform is specified, in which case the warped points will be transformed to the corresponding coordinate system.
Parameters#
points : array, shape (N, dim) or Streamlines object
 coord2worldarray, shape (dim+1, dim+1), optional
affine matrix mapping points to world coordinates
 world2coordarray, shape (dim+1, dim+1), optional
affine matrix mapping world coordinates to points
 warp_endomorphism(phi)#
Composition of this DiffeomorphicMap with a given endomorphism
Creates a new DiffeomorphicMap C with the same properties as self and composes its displacement fields with phi’s corresponding fields. The resulting diffeomorphism is of the form C(x) = phi(self(x)) with inverse C^{1}(y) = self^{1}(phi^{1}(y)). We assume that phi is an endomorphism with the same discretization and domain affine as self to ensure that the composition inherits self’s properties (we also assume that the prealigning matrix of phi is None or identity).
Parameters#
 phiDiffeomorphicMap object
the endomorphism to be warped by this diffeomorphic map
Returns#
 compositionthe composition of this diffeomorphic map with the
endomorphism given as input
Notes#
The problem with our current representation of a DiffeomorphicMap is that the set of Diffeomorphism that can be represented this way (a prealigning matrix followed by a nonlinear endomorphism given as a displacement field) is not closed under the composition operation.
Supporting a general DiffeomorphicMap class, closed under composition, may be extremely costly computationally, and the kind of transformations we actually need for Avants’ midpoint algorithm (SyN) are much simpler.
DiffeomorphicRegistration
#
 class dipy.align.imwarp.DiffeomorphicRegistration(metric=None)#
Bases:
object
 __init__(metric=None)#
Diffeomorphic Registration
This abstract class defines the interface to be implemented by any optimization algorithm for diffeomorphic registration.
Parameters#
 metricSimilarityMetric object
the object measuring the similarity of the two images. The registration algorithm will minimize (or maximize) the provided similarity.
 abstract get_map()#
Returns the resulting diffeomorphic map after optimization
 abstract optimize()#
Starts the metric optimization
This is the main function each specialized class derived from this must implement. Upon completion, the deformation field must be available from the forward transformation model.
 set_level_iters(level_iters)#
Sets the number of iterations at each pyramid level
Establishes the maximum number of iterations to be performed at each level of the Gaussian pyramid, similar to ANTS.
Parameters#
 level_iterslist
the number of iterations at each level of the Gaussian pyramid. level_iters[0] corresponds to the finest level, level_iters[n1] the coarsest, where n is the length of the list
SymmetricDiffeomorphicRegistration
#
 class dipy.align.imwarp.SymmetricDiffeomorphicRegistration(metric, level_iters=None, step_length=0.25, ss_sigma_factor=0.2, opt_tol=1e05, inv_iter=20, inv_tol=0.001, callback=None)#
Bases:
DiffeomorphicRegistration
 __init__(metric, level_iters=None, step_length=0.25, ss_sigma_factor=0.2, opt_tol=1e05, inv_iter=20, inv_tol=0.001, callback=None)#
Symmetric Diffeomorphic Registration (SyN) Algorithm
Performs the multiresolution optimization algorithm for nonlinear registration using a given similarity metric.
Parameters#
 metricSimilarityMetric object
the metric to be optimized
 level_iterslist of int
the number of iterations at each level of the Gaussian Pyramid (the length of the list defines the number of pyramid levels to be used)
 opt_tolfloat
the optimization will stop when the estimated derivative of the energy profile w.r.t. time falls below this threshold
 inv_iterint
the number of iterations to be performed by the displacement field inversion algorithm
 step_lengthfloat
the length of the maximum displacement vector of the update displacement field at each iteration
 ss_sigma_factorfloat
parameter of the scalespace smoothing kernel. For example, the std. dev. of the kernel will be factor*(2^i) in the isotropic case where i = 0, 1, …, n_scales is the scale
 inv_tolfloat
the displacement field inversion algorithm will stop iterating when the inversion error falls below this threshold
 callbackfunction(SymmetricDiffeomorphicRegistration)
a function receiving a SymmetricDiffeomorphicRegistration object to be called after each iteration (this optimizer will call this function passing self as parameter)
 get_map()#
Return the resulting diffeomorphic map.
Returns the DiffeomorphicMap registering the moving image towards the static image.
 optimize(static, moving, static_grid2world=None, moving_grid2world=None, prealign=None)#
Starts the optimization
Parameters#
 staticarray, shape (S, R, C) or (R, C)
the image to be used as reference during optimization. The displacement fields will have the same discretization as the static image.
 movingarray, shape (S, R, C) or (R, C)
the image to be used as “moving” during optimization. Since the deformation fields’ discretization is the same as the static image, it is necessary to prealign the moving image to ensure its domain lies inside the domain of the deformation fields. This is assumed to be accomplished by “prealigning” the moving image towards the static using an affine transformation given by the ‘prealign’ matrix
 static_grid2worldarray, shape (dim+1, dim+1)
the voxeltospace transformation associated to the static image
 moving_grid2worldarray, shape (dim+1, dim+1)
the voxeltospace transformation associated to the moving image
 prealignarray, shape (dim+1, dim+1)
the affine transformation (operating on the physical space) prealigning the moving image towards the static
Returns#
 static_to_refDiffeomorphicMap object
the diffeomorphic map that brings the moving image towards the static one in the forward direction (i.e. by calling static_to_ref.transform) and the static image towards the moving one in the backward direction (i.e. by calling static_to_ref.transform_inverse).
 update(current_displacement, new_displacement, disp_world2grid, time_scaling)#
Composition of the current displacement field with the given field
Interpolates new displacement at the locations defined by current_displacement. Equivalently, computes the composition C of the given displacement fields as C(x) = B(A(x)), where A is current_displacement and B is new_displacement. This function is intended to be used with deformation fields of the same sampling (e.g. to be called by a registration algorithm).
Parameters#
 current_displacementarray, shape (R’, C’, 2) or (S’, R’, C’, 3)
the displacement field defining where to interpolate new_displacement
 new_displacementarray, shape (R, C, 2) or (S, R, C, 3)
the displacement field to be warped by current_displacement
 disp_world2gridarray, shape (dim+1, dim+1)
the spacetogrid transform associated with the displacements’ grid (we assume that both displacements are discretized over the same grid)
 time_scalingfloat
scaling factor applied to d2. The effect may be interpreted as moving d1 displacements along a factor (time_scaling) of d2.
Returns#
 updatedarray, shape (the same as new_displacement)
the warped displacement field
mean_norm : the mean norm of all vectors in current_displacement
RegistrationStages#
 dipy.align.imwarp.RegistrationStages()#
logger#
 dipy.align.imwarp.logger()#
Instances of the Logger class represent a single logging channel. A “logging channel” indicates an area of an application. Exactly how an “area” is defined is up to the application developer. Since an application can have any number of areas, logging channels are identified by a unique string. Application areas can be nested (e.g. an area of “input processing” might include subareas “read CSV files”, “read XLS files” and “read Gnumeric files”). To cater for this natural nesting, channel names are organized into a namespace hierarchy where levels are separated by periods, much like the Java or Python package namespace. So in the instance given above, channel names might be “input” for the upper level, and “input.csv”, “input.xls” and “input.gnu” for the sublevels. There is no arbitrary limit to the depth of nesting.
mult_aff#
 dipy.align.imwarp.mult_aff(A, B)#
Returns the matrix product A.dot(B) considering None as the identity
Parameters#
A : array, shape (n,k) B : array, shape (k,m)
Returns#
The matrix product A.dot(B). If any of the input matrices is None, it is treated as the identity matrix. If both matrices are None, None is returned
get_direction_and_spacings#
 dipy.align.imwarp.get_direction_and_spacings(affine, dim)#
Extracts the rotational and spacing components from a matrix
Extracts the rotational and spacing (voxel dimensions) components from a matrix. An image gradient represents the local variation of the image’s gray values per voxel. Since we are iterating on the physical space, we need to compute the gradients as variation per millimeter, so we need to divide each gradient’s component by the voxel size along the corresponding axis, that’s what the spacings are used for. Since the image’s gradients are oriented along the grid axes, we also need to reorient the gradients to be given in physical space coordinates.
Parameters#
 affinearray, shape (k, k), k = 3, 4
the matrix transforming grid coordinates to physical space.
Returns#
 directionarray, shape (k1, k1)
the rotational component of the input matrix
 spacingsarray, shape (k1,)
the scaling component (voxel size) of the matrix
SimilarityMetric
#
 class dipy.align.metrics.SimilarityMetric(dim)#
Bases:
object
 __init__(dim)#
Similarity Metric abstract class
A similarity metric is in charge of keeping track of the numerical value of the similarity (or distance) between the two given images. It also computes the update field for the forward and inverse displacement fields to be used in a gradientbased optimization algorithm. Note that this metric does not depend on any transformation (affine or nonlinear) so it assumes the static and moving images are already warped
Parameters#
 dimint (either 2 or 3)
the dimension of the image domain
 abstract compute_backward()#
Computes one step bringing the static image towards the moving.
Computes the backward update field to register the static image towards the moving image in a gradientbased optimization algorithm
 abstract compute_forward()#
Computes one step bringing the reference image towards the static.
Computes the forward update field to register the moving image towards the static image in a gradientbased optimization algorithm
 abstract free_iteration()#
Releases the resources no longer needed by the metric
This method is called by the RegistrationOptimizer after the required iterations have been computed (forward and / or backward) so that the SimilarityMetric can safely delete any data it computed as part of the initialization
 abstract get_energy()#
Numerical value assigned by this metric to the current image pair
Must return the numeric value of the similarity between the given static and moving images
 abstract initialize_iteration()#
Prepares the metric to compute one displacement field iteration.
This method will be called before any compute_forward or compute_backward call, this allows the Metric to precompute any useful information for speeding up the update computations. This initialization was needed in ANTS because the updates are called once per voxel. In Python this is unpractical, though.
 set_levels_above(levels)#
Informs the metric how many pyramid levels are above the current one
Informs this metric the number of pyramid levels above the current one. The metric may change its behavior (e.g. number of inner iterations) accordingly
Parameters#
 levelsint
the number of levels above the current Gaussian Pyramid level
 set_levels_below(levels)#
Informs the metric how many pyramid levels are below the current one
Informs this metric the number of pyramid levels below the current one. The metric may change its behavior (e.g. number of inner iterations) accordingly
Parameters#
 levelsint
the number of levels below the current Gaussian Pyramid level
 set_moving_image(moving_image, moving_affine, moving_spacing, moving_direction)#
Sets the moving image being compared against the static one.
Sets the moving image. The default behavior (of this abstract class) is simply to assign the reference to an attribute, but generalizations of the metric may need to perform other operations
Parameters#
 moving_imagearray, shape (R, C) or (S, R, C)
the moving image
 set_static_image(static_image, static_affine, static_spacing, static_direction)#
Sets the static image being compared against the moving one.
Sets the static image. The default behavior (of this abstract class) is simply to assign the reference to an attribute, but generalizations of the metric may need to perform other operations
Parameters#
 static_imagearray, shape (R, C) or (S, R, C)
the static image
 use_moving_image_dynamics(original_moving_image, transformation)#
This is called by the optimizer just after setting the moving image
This method allows the metric to compute any useful information from knowing how the current static image was generated (as the transformation of an original static image). This method is called by the optimizer just after it sets the static image. Transformation will be an instance of DiffeomorficMap or None if the original_moving_image equals self.moving_image.
Parameters#
 original_moving_imagearray, shape (R, C) or (S, R, C)
original image from which the current moving image was generated
 transformationDiffeomorphicMap object
the transformation that was applied to the original image to generate the current moving image
 use_static_image_dynamics(original_static_image, transformation)#
This is called by the optimizer just after setting the static image.
This method allows the metric to compute any useful information from knowing how the current static image was generated (as the transformation of an original static image). This method is called by the optimizer just after it sets the static image. Transformation will be an instance of DiffeomorficMap or None if the original_static_image equals self.moving_image.
Parameters#
 original_static_imagearray, shape (R, C) or (S, R, C)
original image from which the current static image was generated
 transformationDiffeomorphicMap object
the transformation that was applied to original image to generate the current static image
CCMetric
#
 class dipy.align.metrics.CCMetric(dim, sigma_diff=2.0, radius=4)#
Bases:
SimilarityMetric
 __init__(dim, sigma_diff=2.0, radius=4)#
Normalized CrossCorrelation Similarity metric.
Parameters#
 dimint (either 2 or 3)
the dimension of the image domain
 sigma_diffthe standard deviation of the Gaussian smoothing kernel to
be applied to the update field at each iteration
 radiusint
the radius of the squared (cubic) neighborhood at each voxel to be considered to compute the cross correlation
 compute_backward()#
Computes one step bringing the static image towards the moving.
Computes the update displacement field to be used for registration of the static image towards the moving image
 compute_forward()#
Computes one step bringing the moving image towards the static.
Computes the update displacement field to be used for registration of the moving image towards the static image
 free_iteration()#
Frees the resources allocated during initialization
 get_energy()#
Numerical value assigned by this metric to the current image pair
Returns the Cross Correlation (data term) energy computed at the largest iteration
 initialize_iteration()#
Prepares the metric to compute one displacement field iteration.
Precomputes the crosscorrelation factors for efficient computation of the gradient of the Cross Correlation w.r.t. the displacement field. It also precomputes the image gradients in the physical space by reorienting the gradients in the voxel space using the corresponding affine transformations.
EMMetric
#
 class dipy.align.metrics.EMMetric(dim, smooth=1.0, inner_iter=5, q_levels=256, double_gradient=True, step_type='gauss_newton')#
Bases:
SimilarityMetric
 __init__(dim, smooth=1.0, inner_iter=5, q_levels=256, double_gradient=True, step_type='gauss_newton')#
ExpectationMaximization Metric
Similarity metric based on the ExpectationMaximization algorithm to handle multimodal images. The transfer function is modeled as a set of hidden random variables that are estimated at each iteration of the algorithm.
Parameters#
 dimint (either 2 or 3)
the dimension of the image domain
 smoothfloat
smoothness parameter, the larger the value the smoother the deformation field
 inner_iterint
number of iterations to be performed at each level of the multi resolution GaussSeidel optimization algorithm (this is not the number of steps per Gaussian Pyramid level, that parameter must be set for the optimizer, not the metric)
 q_levelsnumber of quantization levels (equal to the number of hidden
variables in the EM algorithm)
 double_gradientboolean
if True, the gradient of the expected static image under the moving modality will be added to the gradient of the moving image, similarly, the gradient of the expected moving image under the static modality will be added to the gradient of the static image.
 step_typestring (‘gauss_newton’, ‘demons’)
the optimization schedule to be used in the multiresolution GaussSeidel optimization algorithm (not used if Demons Step is selected)
 compute_backward()#
Computes one step bringing the static image towards the moving.
Computes the update displacement field to be used for registration of the static image towards the moving image
 compute_demons_step(forward_step=True)#
Demons step for EM metric
Parameters#
 forward_stepboolean
if True, computes the Demons step in the forward direction (warping the moving towards the static image). If False, computes the backward step (warping the static image to the moving image)
Returns#
 displacementarray, shape (R, C, 2) or (S, R, C, 3)
the Demons step
 compute_forward()#
Computes one step bringing the reference image towards the static.
Computes the forward update field to register the moving image towards the static image in a gradientbased optimization algorithm
 compute_gauss_newton_step(forward_step=True)#
Computes the GaussNewton energy minimization step
Computes the Newton step to minimize this energy, i.e., minimizes the linearized energy function with respect to the regularized displacement field (this step does not require postsmoothing, as opposed to the demons step, which does not include regularization). To accelerate convergence we use the multigrid GaussSeidel algorithm proposed by Bruhn and Weickert et al [Bruhn05]
Parameters#
 forward_stepboolean
if True, computes the Newton step in the forward direction (warping the moving towards the static image). If False, computes the backward step (warping the static image to the moving image)
Returns#
 displacementarray, shape (R, C, 2) or (S, R, C, 3)
the Newton step
References#
 [Bruhn05] Andres Bruhn and Joachim Weickert, “Towards ultimate motion
estimation: combining highest accuracy with realtime performance”, 10th IEEE International Conference on Computer Vision, 2005. ICCV 2005.
 free_iteration()#
Frees the resources allocated during initialization
 get_energy()#
The numerical value assigned by this metric to the current image pair
Returns the EM (data term) energy computed at the largest iteration
 initialize_iteration()#
Prepares the metric to compute one displacement field iteration.
Precomputes the transfer functions (hidden random variables) and variances of the estimators. Also precomputes the gradient of both input images. Note that once the images are transformed to the opposite modality, the gradient of the transformed images can be used with the gradient of the corresponding modality in the same fashion as diffdemons does for monomodality images. If the flag self.use_double_gradient is True these gradients are averaged.
 use_moving_image_dynamics(original_moving_image, transformation)#
This is called by the optimizer just after setting the moving image.
EMMetric takes advantage of the image dynamics by computing the current moving image mask from the original_moving_image mask (warped by nearest neighbor interpolation)
Parameters#
 original_moving_imagearray, shape (R, C) or (S, R, C)
the original moving image from which the current moving image was generated, the current moving image is the one that was provided via ‘set_moving_image(…)’, which may not be the same as the original moving image but a warped version of it.
 transformationDiffeomorphicMap object
the transformation that was applied to the original_moving_image to generate the current moving image
 use_static_image_dynamics(original_static_image, transformation)#
This is called by the optimizer just after setting the static image.
EMMetric takes advantage of the image dynamics by computing the current static image mask from the originalstaticImage mask (warped by nearest neighbor interpolation)
Parameters#
 original_static_imagearray, shape (R, C) or (S, R, C)
the original static image from which the current static image was generated, the current static image is the one that was provided via ‘set_static_image(…)’, which may not be the same as the original static image but a warped version of it (even the static image changes during Symmetric Normalization, not only the moving one).
 transformationDiffeomorphicMap object
the transformation that was applied to the original_static_image to generate the current static image
SSDMetric
#
 class dipy.align.metrics.SSDMetric(dim, smooth=4, inner_iter=10, step_type='demons')#
Bases:
SimilarityMetric
 __init__(dim, smooth=4, inner_iter=10, step_type='demons')#
Sum of Squared Differences (SSD) Metric
Similarity metric for (monomodal) nonlinear image registration defined by the sum of squared differences (SSD)
Parameters#
 dimint (either 2 or 3)
the dimension of the image domain
 smoothfloat
smoothness parameter, the larger the value the smoother the deformation field
 inner_iterint
number of iterations to be performed at each level of the multi resolution GaussSeidel optimization algorithm (this is not the number of steps per Gaussian Pyramid level, that parameter must be set for the optimizer, not the metric)
 step_typestring
the displacement field step to be computed when ‘compute_forward’ and ‘compute_backward’ are called. Either ‘demons’ or ‘gauss_newton’
 compute_backward()#
Computes one step bringing the static image towards the moving.
Computes the updated displacement field to be used for registration of the static image towards the moving image
 compute_demons_step(forward_step=True)#
Demons step for SSD metric
Computes the demons step proposed by Vercauteren et al.[Vercauteren09] for the SSD metric.
Parameters#
 forward_stepboolean
if True, computes the Demons step in the forward direction (warping the moving towards the static image). If False, computes the backward step (warping the static image to the moving image)
Returns#
 displacementarray, shape (R, C, 2) or (S, R, C, 3)
the Demons step
References#
 [Vercauteren09] Tom Vercauteren, Xavier Pennec, Aymeric Perchant,
Nicholas Ayache, “Diffeomorphic Demons: Efficient Nonparametric Image Registration”, Neuroimage 2009
 compute_forward()#
Computes one step bringing the reference image towards the static.
Computes the update displacement field to be used for registration of the moving image towards the static image
 compute_gauss_newton_step(forward_step=True)#
Computes the GaussNewton energy minimization step
Minimizes the linearized energy function (Newton step) defined by the sum of squared differences of corresponding pixels of the input images with respect to the displacement field.
Parameters#
 forward_stepboolean
if True, computes the Newton step in the forward direction (warping the moving towards the static image). If False, computes the backward step (warping the static image to the moving image)
Returns#
 displacementarray, shape = static_image.shape + (3,)
if forward_step==True, the forward SSD GaussNewton step, else, the backward step
 free_iteration()#
Nothing to free for the SSD metric
 get_energy()#
The numerical value assigned by this metric to the current image pair
Returns the Sum of Squared Differences (data term) energy computed at the largest iteration
 initialize_iteration()#
Prepares the metric to compute one displacement field iteration.
Precomputes the gradient of the input images to be used in the computation of the forward and backward steps.
v_cycle_2d#
 dipy.align.metrics.v_cycle_2d(n, k, delta_field, sigma_sq_field, gradient_field, target, lambda_param, displacement, depth=0)#
Multiresolution GaussSeidel solver using Vtype cycles
Multiresolution GaussSeidel solver: solves the GaussNewton linear system by first filtering (GSiterate) the current level, then solves for the residual at a coarser resolution and finally refines the solution at the current resolution. This scheme corresponds to the Vcycle proposed by Bruhn and Weickert[Bruhn05].
Parameters#
 nint
number of levels of the multiresolution algorithm (it will be called recursively until level n == 0)
 kint
the number of iterations at each multiresolution level
 delta_fieldarray, shape (R, C)
the difference between the static and moving image (the ‘derivative w.r.t. time’ in the optical flow model)
 sigma_sq_fieldarray, shape (R, C)
the variance of the gray level value at each voxel, according to the EM model (for SSD, it is 1 for all voxels). Inf and 0 values are processed specially to support infinite and zero variance.
 gradient_fieldarray, shape (R, C, 2)
the gradient of the moving image
 targetarray, shape (R, C, 2)
righthand side of the linear system to be solved in the Weickert’s multiresolution algorithm
 lambda_paramfloat
smoothness parameter, the larger its value the smoother the displacement field
 displacementarray, shape (R, C, 2)
the displacement field to start the optimization from
Returns#
 energythe energy of the EM (or SSD if sigmafield[…]==1) metric at this
iteration
References#
 [Bruhn05] Andres Bruhn and Joachim Weickert, “Towards ultimate motion
estimation: combining the highest accuracy with realtime performance”, 10th IEEE International Conference on Computer Vision, 2005. ICCV 2005.
v_cycle_3d#
 dipy.align.metrics.v_cycle_3d(n, k, delta_field, sigma_sq_field, gradient_field, target, lambda_param, displacement, depth=0)#
Multiresolution GaussSeidel solver using Vtype cycles
Multiresolution GaussSeidel solver: solves the linear system by first filtering (GSiterate) the current level, then solves for the residual at a coarser resolution and finally refines the solution at the current resolution. This scheme corresponds to the Vcycle proposed by Bruhn and Weickert[1]. [1] Andres Bruhn and Joachim Weickert, “Towards ultimate motion estimation:
combining highest accuracy with realtime performance”, 10th IEEE International Conference on Computer Vision, 2005. ICCV 2005.
Parameters#
 nint
number of levels of the multiresolution algorithm (it will be called recursively until level n == 0)
 kint
the number of iterations at each multiresolution level
 delta_fieldarray, shape (S, R, C)
the difference between the static and moving image (the ‘derivative w.r.t. time’ in the optical flow model)
 sigma_sq_fieldarray, shape (S, R, C)
the variance of the gray level value at each voxel, according to the EM model (for SSD, it is 1 for all voxels). Inf and 0 values are processed specially to support infinite and zero variance.
 gradient_fieldarray, shape (S, R, C, 3)
the gradient of the moving image
 targetarray, shape (S, R, C, 3)
righthand side of the linear system to be solved in the Weickert’s multiresolution algorithm
 lambda_paramfloat
smoothness parameter, the larger its value the smoother the displacement field
 displacementarray, shape (S, R, C, 3)
the displacement field to start the optimization from
Returns#
 energythe energy of the EM (or SSD if sigmafield[…]==1) metric at this
iteration
ParzenJointHistogram
#
 class dipy.align.parzenhist.ParzenJointHistogram(nbins)#
Bases:
object
 __init__(nbins)#
Computes joint histogram and derivatives with Parzen windows
Base class to compute joint and marginal probability density functions and their derivatives with respect to a transform’s parameters. The smooth histograms are computed by using Parzen windows [Parzen62] with a cubic spline kernel, as proposed by Mattes et al. [Mattes03]. This implementation is not tied to any optimization (registration) method, the idea is that informationtheoretic matching functionals (such as Mutual Information) can inherit from this class to perform the lowlevel computations of the joint intensity distributions and its gradient w.r.t. the transform parameters. The derived class can then compute the similarity/dissimilarity measure and gradient, and finally communicate the results to the appropriate optimizer.
Parameters#
 nbinsint
the number of bins of the joint and marginal probability density functions (the actual number of bins of the joint PDF is nbins**2)
References#
 [Parzen62] E. Parzen. On the estimation of a probability density
function and the mode. Annals of Mathematical Statistics, 33(3), 10651076, 1962.
 [Mattes03] Mattes, D., Haynor, D. R., Vesselle, H., Lewellen, T. K.,
& Eubank, W. PETCT image registration in the chest using freeform deformations. IEEE Transactions on Medical Imaging, 22(1), 1208, 2003.
Notes#
We need this class in cython to allow _joint_pdf_gradient_dense_2d and _joint_pdf_gradient_dense_3d to use a nogil Jacobian function (obtained from an instance of the Transform class), which allows us to evaluate Jacobians at all the sampling points (maybe the full grid) inside a nogil loop.
The reason we need a class is to encapsulate all the parameters related to the joint and marginal distributions.
 bin_index(xnorm)#
Bin index associated with the given normalized intensity
The return value is an integer in [padding, nbins  1  padding]
Parameters#
 xnormfloat
intensity value normalized to the range covered by the histogram
Returns#
 binint
the bin index associated with the given normalized intensity
 bin_normalize_moving(x)#
Maps intensity x to the range covered by the moving histogram
If the input intensity is in [self.mmin, self.mmax] then the normalized intensity will be in [self.padding, self.nbins  self.padding]
Parameters#
 xfloat
the intensity to be normalized
Returns#
 xnormfloat
normalized intensity to the range covered by the moving histogram
 bin_normalize_static(x)#
Maps intensity x to the range covered by the static histogram
If the input intensity is in [self.smin, self.smax] then the normalized intensity will be in [self.padding, self.nbins  self.padding]
Parameters#
 xfloat
the intensity to be normalized
Returns#
 xnormfloat
normalized intensity to the range covered by the static histogram
 setup(static, moving, smask=None, mmask=None)#
Compute histogram settings to store the PDF of input images
Parameters#
 staticarray
static image
 movingarray
moving image
 smaskarray
mask of static object being registered (a binary array with 1’s inside the object of interest and 0’s along the background). If None, the behaviour is equivalent to smask=ones_like(static)
 mmaskarray
mask of moving object being registered (a binary array with 1’s inside the object of interest and 0’s along the background). If None, the behaviour is equivalent to mmask=ones_like(static)
 update_gradient_dense(theta, transform, static, moving, grid2world, mgradient, smask=None, mmask=None)#
Computes the Gradient of the joint PDF w.r.t. transform parameters
Computes the vector of partial derivatives of the joint histogram w.r.t. each transformation parameter.
The gradient is stored in self.joint_grad.
Parameters#
 thetaarray, shape (n,)
parameters of the transformation to compute the gradient from
 transforminstance of Transform
the transformation with respect to whose parameters the gradient must be computed
 staticarray, shape (S, R, C)
static image
 movingarray, shape (S, R, C)
moving image
 grid2worldarray, shape (4, 4)
we assume that both images have already been sampled at a common grid. This transform must map voxel coordinates of this common grid to physical coordinates of its corresponding voxel in the moving image. For example, if the moving image was sampled on the static image’s grid (this is the typical setting) using an aligning matrix A, then
grid2world = A.dot(static_affine)
where static_affine is the transformation mapping static image’s grid coordinates to physical space.
 mgradientarray, shape (S, R, C, 3)
the gradient of the moving image
 smaskarray, shape (S, R, C), optional
mask of static object being registered (a binary array with 1’s inside the object of interest and 0’s along the background). The default is None, indicating all voxels are considered.
 mmaskarray, shape (S, R, C), optional
mask of moving object being registered (a binary array with 1’s inside the object of interest and 0’s along the background). The default is None, indicating all voxels are considered.
 update_gradient_sparse(theta, transform, sval, mval, sample_points, mgradient)#
Computes the Gradient of the joint PDF w.r.t. transform parameters
Computes the vector of partial derivatives of the joint histogram w.r.t. each transformation parameter.
The list of intensities sval and mval are assumed to be sampled from the static and moving images, respectively, at the same physical points. Of course, the images may not be perfectly aligned at the moment the sampling was performed. The resulting gradient corresponds to the paired intensities according to the alignment at the moment the images were sampled.
The gradient is stored in self.joint_grad.
Parameters#
 thetaarray, shape (n,)
parameters to compute the gradient at
 transforminstance of Transform
the transformation with respect to whose parameters the gradient must be computed
 svalarray, shape (m,)
sampled intensities from the static image at sampled_points
 mvalarray, shape (m,)
sampled intensities from the moving image at sampled_points
 sample_pointsarray, shape (m, 3)
coordinates (in physical space) of the points the images were sampled at
 mgradientarray, shape (m, 3)
the gradient of the moving image at the sample points
 update_pdfs_dense(static, moving, smask=None, mmask=None)#
Computes the Probability Density Functions of two images
The joint PDF is stored in self.joint. The marginal distributions corresponding to the static and moving images are computed and stored in self.smarginal and self.mmarginal, respectively.
Parameters#
 staticarray, shape (S, R, C)
static image
 movingarray, shape (S, R, C)
moving image
 smaskarray, shape (S, R, C)
mask of static object being registered (a binary array with 1’s inside the object of interest and 0’s along the background). If None, ones_like(static) is used as mask.
 mmaskarray, shape (S, R, C)
mask of moving object being registered (a binary array with 1’s inside the object of interest and 0’s along the background). If None, ones_like(moving) is used as mask.
 update_pdfs_sparse(sval, mval)#
Computes the Probability Density Functions from a set of samples
The list of intensities sval and mval are assumed to be sampled from the static and moving images, respectively, at the same physical points. Of course, the images may not be perfectly aligned at the moment the sampling was performed. The resulting distributions corresponds to the paired intensities according to the alignment at the moment the images were sampled.
The joint PDF is stored in self.joint. The marginal distributions corresponding to the static and moving images are computed and stored in self.smarginal and self.mmarginal, respectively.
Parameters#
 svalarray, shape (n,)
sampled intensities from the static image at sampled_points
 mvalarray, shape (n,)
sampled intensities from the moving image at sampled_points
compute_parzen_mi#
 dipy.align.parzenhist.compute_parzen_mi(joint, joint_gradient, smarginal, mmarginal, mi_gradient)#
Computes the mutual information and its gradient (if requested)
Parameters#
 jointarray, shape (nbins, nbins)
the joint intensity distribution
 joint_gradientarray, shape (nbins, nbins, n)
the gradient of the joint distribution w.r.t. the transformation parameters
 smarginalarray, shape (nbins,)
the marginal intensity distribution of the static image
 mmarginalarray, shape (nbins,)
the marginal intensity distribution of the moving image
 mi_gradientarray, shape (n,)
the buffer in which to write the gradient of the mutual information. If None, the gradient is not computed
cubic_spline#
cubic_spline_derivative#
sample_domain_regular#
 dipy.align.parzenhist.sample_domain_regular(k, shape, grid2world, sigma=0.25, rng=None)#
Take floor(total_voxels/k) samples from a (2D or 3D) grid
The sampling is made by taking all pixels whose index (in lexicographical order) is a multiple of k. Each selected point is slightly perturbed by adding a realization of a normally distributed random variable and then mapped to physical space by the given gridtospace transform.
The lexicographical order of a pixels in a grid of shape (a, b, c) is defined by assigning to each voxel position (i, j, k) the integer index
F((i, j, k)) = i * (b * c) + j * (c) + k
and sorting increasingly by this index.
Parameters#
 kint
the sampling rate, as described before
 shapearray, shape (dim,)
the shape of the grid to be sampled
 grid2worldarray, shape (dim+1, dim+1)
the gridtospace transform
 sigmafloat
the standard deviation of the Normal random distortion to be applied to the sampled points
Returns#
 samplesarray, shape (total_pixels//k, dim)
the matrix whose rows are the sampled points
Examples#
>>> from dipy.align.parzenhist import sample_domain_regular >>> import dipy.align.vector_fields as vf >>> shape = np.array((10, 10), dtype=np.int32) >>> sigma = 0 >>> dim = len(shape) >>> grid2world = np.eye(dim+1) >>> n = shape[0]*shape[1] >>> k = 2 >>> samples = sample_domain_regular(k, shape, grid2world, sigma) >>> (samples.shape[0], samples.shape[1]) == (n//k, dim) True >>> isamples = np.array(samples, dtype=np.int32) >>> indices = (isamples[:, 0] * shape[1] + isamples[:, 1]) >>> len(set(indices)) == len(indices) True >>> (indices%k).sum() 0
reslice#
 dipy.align.reslice.reslice(data, affine, zooms, new_zooms, order=1, mode='constant', cval=0, num_processes=1)#
Reslice data with new voxel resolution defined by
new_zooms
.Parameters#
 dataarray, shape (I,J,K) or (I,J,K,N)
3d volume or 4d volume with datasets
 affinearray, shape (4,4)
mapping from voxel coordinates to world coordinates
 zoomstuple, shape (3,)
voxel size for (i,j,k) dimensions
 new_zoomstuple, shape (3,)
new voxel size for (i,j,k) after resampling
 orderint, from 0 to 5
order of interpolation for resampling/reslicing, 0 nearest interpolation, 1 trilinear etc.. if you don’t want any smoothing 0 is the option you need.
 modestring (‘constant’, ‘nearest’, ‘reflect’ or ‘wrap’)
Points outside the boundaries of the input are filled according to the given mode.
 cvalfloat
Value used for points outside the boundaries of the input if mode=’constant’.
 num_processesint, optional
Split the calculation to a pool of children processes. This only applies to 4D data arrays. Default is 1. If < 0 the maximal number of cores minus
num_processes + 1
is used (enter 1 to use as many cores as possible). 0 raises an error.
Returns#
 data2array, shape (I,J,K) or (I,J,K,N)
datasets resampled into isotropic voxel size
 affine2array, shape (4,4)
new affine for the resampled image
Examples#
>>> from dipy.io.image import load_nifti >>> from dipy.align.reslice import reslice >>> from dipy.data import get_fnames >>> f_name = get_fnames('aniso_vox') >>> data, affine, zooms = load_nifti(f_name, return_voxsize=True) >>> data.shape == (58, 58, 24) True >>> zooms (4.0, 4.0, 5.0) >>> new_zooms = (3.,3.,3.) >>> new_zooms (3.0, 3.0, 3.0) >>> data2, affine2 = reslice(data, affine, zooms, new_zooms) >>> data2.shape == (77, 77, 40) True
ScaleSpace
#
 class dipy.align.scalespace.ScaleSpace(image, num_levels, image_grid2world=None, input_spacing=None, sigma_factor=0.2, mask0=False)#
Bases:
object
 __init__(image, num_levels, image_grid2world=None, input_spacing=None, sigma_factor=0.2, mask0=False)#
ScaleSpace. Computes the Scale Space representation of an image. The scale space is simply a list of images produced by smoothing the input image with a Gaussian kernel with increasing smoothing parameter. If the image’s voxels are isotropic, the smoothing will be the same along all directions: at level L = 0, 1, …, the sigma is given by \(s * ( 2^L  1 )\). If the voxel dimensions are not isotropic, then the smoothing is weaker along low resolution directions. Parameters ——— image : array, shape (r,c) or (s, r, c) where s is the number of slices, r is the number of rows and c is the number of columns of the input image. num_levels : int the desired number of levels (resolutions) of the scale space image_grid2world : array, shape (dim + 1, dim + 1), optional the gridtospace transform of the image grid. The default is the identity matrix input_spacing : array, shape (dim,), optional the spacing (voxel size) between voxels in physical space. The default is 1.0 along all axes sigma_factor : float, optional the smoothing factor to be used in the construction of the scale space. The default is 0.2 mask0 : Boolean, optional if True, all smoothed images will be zero at all voxels that are zero in the input image. The default is False.
 get_affine(level)#
Voxeltospace transformation at a given level.
Returns the voxeltospace transformation associated with the subsampled image at a particular resolution of the scale space (note that this object does not explicitly subsample the smoothed images, but only provides the properties the subsampled images must have).
Parameters#
 levelint, 0 <= from_level < L, (L = number of resolutions)
the scale space level to get affine transform from
Returns#
the affine (voxeltospace) transform at the requested resolution or None if an invalid level was requested
 get_affine_inv(level)#
Spacetovoxel transformation at a given level.
Returns the spacetovoxel transformation associated with the subsampled image at a particular resolution of the scale space (note that this object does not explicitly subsample the smoothed images, but only provides the properties the subsampled images must have).
Parameters#
 levelint, 0 <= from_level < L, (L = number of resolutions)
the scale space level to get the inverse transform from
Returns#
the inverse (spacetovoxel) transform at the requested resolution or None if an invalid level was requested
 get_domain_shape(level)#
Shape the subsampled image must have at a particular level.
Returns the shape the subsampled image must have at a particular resolution of the scale space (note that this object does not explicitly subsample the smoothed images, but only provides the properties the subsampled images must have).
Parameters#
 levelint, 0 <= from_level < L, (L = number of resolutions)
the scale space level to get the subsampled shape from
Returns#
the subsampled shape at the requested resolution or None if an invalid level was requested
 get_expand_factors(from_level, to_level)#
Ratio of voxel size from pyramid level from_level to to_level.
Given two scale space resolutions a = from_level, b = to_level, returns the ratio of voxels size at level b to voxel size at level a (the factor that must be used to multiply voxels at level a to ‘expand’ them to level b).
Parameters#
 from_levelint, 0 <= from_level < L, (L = number of resolutions)
the resolution to expand voxels from
 to_levelint, 0 <= to_level < from_level
the resolution to expand voxels to
Returns#
 factorsarray, shape (k,), k = 2, 3
the expand factors (a scalar for each voxel dimension)
 get_image(level)#
Smoothed image at a given level.
Returns the smoothed image at the requested level in the Scale Space.
Parameters#
 levelint, 0 <= from_level < L, (L = number of resolutions)
the scale space level to get the smooth image from
Returns#
the smooth image at the requested resolution or None if an invalid level was requested
 get_scaling(level)#
Adjustment factor for inputspacing to reflect voxel sizes at level.
Returns the scaling factor that needs to be applied to the input spacing (the voxel sizes of the image at level 0 of the scale space) to transform them to voxel sizes at the requested level.
Parameters#
 levelint, 0 <= from_level < L, (L = number of resolutions)
the scale space level to get the scalings from
Returns#
the scaling factors from the original spacing to the spacings at the requested level
 get_sigmas(level)#
Smoothing parameters used at a given level.
Returns the smoothing parameters (a scalar for each axis) used at the requested level of the scale space
Parameters#
 levelint, 0 <= from_level < L, (L = number of resolutions)
the scale space level to get the smoothing parameters from
Returns#
the smoothing parameters at the requested level
 get_spacing(level)#
Spacings the subsampled image must have at a particular level.
Returns the spacings (voxel sizes) the subsampled image must have at a particular resolution of the scale space (note that this object does not explicitly subsample the smoothed images, but only provides the properties the subsampled images must have).
Parameters#
 levelint, 0 <= from_level < L, (L = number of resolutions)
the scale space level to get the subsampled shape from
Returns#
the spacings (voxel sizes) at the requested resolution or None if an invalid level was requested
IsotropicScaleSpace
#
 class dipy.align.scalespace.IsotropicScaleSpace(image, factors, sigmas, image_grid2world=None, input_spacing=None, mask0=False)#
Bases:
ScaleSpace
 __init__(image, factors, sigmas, image_grid2world=None, input_spacing=None, mask0=False)#
IsotropicScaleSpace.
Computes the Scale Space representation of an image using isotropic smoothing kernels for all scales. The scale space is simply a list of images produced by smoothing the input image with a Gaussian kernel with different smoothing parameters.
This specialization of ScaleSpace allows the user to provide custom scale and smoothing factors for all scales.
Parameters#
 imagearray, shape (r,c) or (s, r, c) where s is the number of
slices, r is the number of rows and c is the number of columns of the input image.
 factorslist of floats
custom scale factors to build the scale space (one factor for each scale).
 sigmaslist of floats
custom smoothing parameter to build the scale space (one parameter for each scale).
 image_grid2worldarray, shape (dim + 1, dim + 1), optional
the gridtospace transform of the image grid. The default is the identity matrix.
 input_spacingarray, shape (dim,), optional
the spacing (voxel size) between voxels in physical space. The default if 1.0 along all axes.
 mask0Boolean, optional
if True, all smoothed images will be zero at all voxels that are zero in the input image. The default is False.
logger#
 dipy.align.scalespace.logger()#
Instances of the Logger class represent a single logging channel. A “logging channel” indicates an area of an application. Exactly how an “area” is defined is up to the application developer. Since an application can have any number of areas, logging channels are identified by a unique string. Application areas can be nested (e.g. an area of “input processing” might include subareas “read CSV files”, “read XLS files” and “read Gnumeric files”). To cater for this natural nesting, channel names are organized into a namespace hierarchy where levels are separated by periods, much like the Java or Python package namespace. So in the instance given above, channel names might be “input” for the upper level, and “input.csv”, “input.xls” and “input.gnu” for the sublevels. There is no arbitrary limit to the depth of nesting.
StreamlineDistanceMetric
#
 class dipy.align.streamlinear.StreamlineDistanceMetric(num_threads=None)#
Bases:
object
 __init__(num_threads=None)#
An abstract class for the metric used for streamline registration.
If the two sets of streamlines match exactly then method
distance
of this object should be minimum.Parameters#
 num_threadsint, optional
Number of threads to be used for OpenMP parallelization. If None (default) the value of OMP_NUM_THREADS environment variable is used if it is set, otherwise all available threads are used. If < 0 the maximal number of threads minus num_threads + 1 is used (enter 1 to use as many threads as possible). 0 raises an error. Only metrics using OpenMP will use this variable.
 abstract distance(xopt)#
calculate distance for current set of parameters.
 abstract setup(static, moving)#
BundleMinDistanceMetric
#
 class dipy.align.streamlinear.BundleMinDistanceMetric(num_threads=None)#
Bases:
StreamlineDistanceMetric
Bundlebased Minimum Distance aka BMD
This is the cost function used by the StreamlineLinearRegistration.
Methods#
setup(static, moving) distance(xopt)
References#
[Garyfallidis14]Garyfallidis et al., “Direct nativespace fiber bundle alignment for group comparisons”, ISMRM, 2014.
 __init__(num_threads=None)#
An abstract class for the metric used for streamline registration.
If the two sets of streamlines match exactly then method
distance
of this object should be minimum.Parameters#
 num_threadsint, optional
Number of threads to be used for OpenMP parallelization. If None (default) the value of OMP_NUM_THREADS environment variable is used if it is set, otherwise all available threads are used. If < 0 the maximal number of threads minus num_threads + 1 is used (enter 1 to use as many threads as possible). 0 raises an error. Only metrics using OpenMP will use this variable.
BundleMinDistanceMatrixMetric
#
 class dipy.align.streamlinear.BundleMinDistanceMatrixMetric(num_threads=None)#
Bases:
StreamlineDistanceMetric
Bundlebased Minimum Distance aka BMD
This is the cost function used by the StreamlineLinearRegistration
Methods#
setup(static, moving) distance(xopt)
Notes#
The difference with BundleMinDistanceMetric is that this creates the entire distance matrix and therefore requires more memory.
 __init__(num_threads=None)#
An abstract class for the metric used for streamline registration.
If the two sets of streamlines match exactly then method
distance
of this object should be minimum.Parameters#
 num_threadsint, optional
Number of threads to be used for OpenMP parallelization. If None (default) the value of OMP_NUM_THREADS environment variable is used if it is set, otherwise all available threads are used. If < 0 the maximal number of threads minus num_threads + 1 is used (enter 1 to use as many threads as possible). 0 raises an error. Only metrics using OpenMP will use this variable.
 distance(xopt)#
Distance calculated from this Metric.
Parameters#
 xoptsequence
List of affine parameters as an 1D vector
 setup(static, moving)#
Setup static and moving sets of streamlines.
Parameters#
 staticstreamlines
Fixed or reference set of streamlines.
 movingstreamlines
Moving streamlines.
Notes#
Call this after the object is initiated and before distance.
Num_threads is not used in this class. Use
BundleMinDistanceMetric
for a faster, threaded and less memory hungry metric
BundleMinDistanceAsymmetricMetric
#
 class dipy.align.streamlinear.BundleMinDistanceAsymmetricMetric(num_threads=None)#
Bases:
BundleMinDistanceMetric
Asymmetric Bundlebased Minimum distance.
This is a cost function that can be used by the StreamlineLinearRegistration class.
 __init__(num_threads=None)#
An abstract class for the metric used for streamline registration.
If the two sets of streamlines match exactly then method
distance
of this object should be minimum.Parameters#
 num_threadsint, optional
Number of threads to be used for OpenMP parallelization. If None (default) the value of OMP_NUM_THREADS environment variable is used if it is set, otherwise all available threads are used. If < 0 the maximal number of threads minus num_threads + 1 is used (enter 1 to use as many threads as possible). 0 raises an error. Only metrics using OpenMP will use this variable.
BundleSumDistanceMatrixMetric
#
 class dipy.align.streamlinear.BundleSumDistanceMatrixMetric(num_threads=None)#
Bases:
BundleMinDistanceMatrixMetric
Bundlebased Sum Distance aka BMD
This is a cost function that can be used by the StreamlineLinearRegistration class.
Methods#
setup(static, moving) distance(xopt)
Notes#
The difference with BundleMinDistanceMatrixMetric is that it uses uses the sum of the distance matrix and not the sum of mins.
 __init__(num_threads=None)#
An abstract class for the metric used for streamline registration.
If the two sets of streamlines match exactly then method
distance
of this object should be minimum.Parameters#
 num_threadsint, optional
Number of threads to be used for OpenMP parallelization. If None (default) the value of OMP_NUM_THREADS environment variable is used if it is set, otherwise all available threads are used. If < 0 the maximal number of threads minus num_threads + 1 is used (enter 1 to use as many threads as possible). 0 raises an error. Only metrics using OpenMP will use this variable.
JointBundleMinDistanceMetric
#
 class dipy.align.streamlinear.JointBundleMinDistanceMetric(num_threads=None)#
Bases:
StreamlineDistanceMetric
Bundlebased Minimum Distance for joint optimization.
This cost function is used by the StreamlineLinearRegistration class when running halfway streamline linear registration for unbiased groupwise bundle registration and atlasing.
It computes the BMD distance after moving both static and moving bundles to a halfway space in between both.
Methods#
setup(static, moving) distance(xopt)
Notes#
In this metric both static and moving bundles are treated equally (i.e., there is no static reference bundle as both are intended to move). The naming convention is kept for consistency.
 __init__(num_threads=None)#
An abstract class for the metric used for streamline registration.
If the two sets of streamlines match exactly then method
distance
of this object should be minimum.Parameters#
 num_threadsint, optional
Number of threads to be used for OpenMP parallelization. If None (default) the value of OMP_NUM_THREADS environment variable is used if it is set, otherwise all available threads are used. If < 0 the maximal number of threads minus num_threads + 1 is used (enter 1 to use as many threads as possible). 0 raises an error. Only metrics using OpenMP will use this variable.
StreamlineLinearRegistration
#
 class dipy.align.streamlinear.StreamlineLinearRegistration(metric=None, x0='rigid', method='LBFGSB', bounds=None, verbose=False, options=None, evolution=False, num_threads=None)#
Bases:
object
 __init__(metric=None, x0='rigid', method='LBFGSB', bounds=None, verbose=False, options=None, evolution=False, num_threads=None)#
Linear registration of 2 sets of streamlines [Garyfallidis15].
Parameters#
 metricStreamlineDistanceMetric,
If None and fast is False then the BMD distance is used. If fast is True then a faster implementation of BMD is used. Otherwise, use the given distance metric.
 x0array or int or str
Initial parametrization for the optimization.
 If 1D array with:
a) 6 elements then only rigid registration is performed with the 3 first elements for translation and 3 for rotation. b) 7 elements also isotropic scaling is performed (similarity). c) 12 elements then translation, rotation (in degrees), scaling and shearing is performed (affine).
Here is an example of x0 with 12 elements:
x0=np.array([0, 10, 0, 40, 0, 0, 2., 1.5, 1, 0.1, 0.5, 0])
This has translation (0, 10, 0), rotation (40, 0, 0) in degrees, scaling (2., 1.5, 1) and shearing (0.1, 0.5, 0).
 If int:
 6
x0 = np.array([0, 0, 0, 0, 0, 0])
 7
x0 = np.array([0, 0, 0, 0, 0, 0, 1.])
 12
x0 = np.array([0, 0, 0, 0, 0, 0, 1., 1., 1, 0, 0, 0])
 If str:
 “rigid”
x0 = np.array([0, 0, 0, 0, 0, 0])
 “similarity”
x0 = np.array([0, 0, 0, 0, 0, 0, 1.])
 “affine”
x0 = np.array([0, 0, 0, 0, 0, 0, 1., 1., 1, 0, 0, 0])
 methodstr,
‘L_BFGS_B’ or ‘Powell’ optimizers can be used. Default is ‘L_BFGS_B’.
 boundslist of tuples or None,
If method == ‘L_BFGS_B’ then we can use bounded optimization. For example for the six parameters of rigid rotation we can set the bounds = [(30, 30), (30, 30), (30, 30),
(45, 45), (45, 45), (45, 45)]
That means that we have set the bounds for the three translations and three rotation axes (in degrees).
 verbosebool, optional.
If True, if True then information about the optimization is shown. Default: False.
 optionsNone or dict,
Extra options to be used with the selected method.
 evolutionboolean
If True save the transformation for each iteration of the optimizer. Default is False. Supported only with Scipy >= 0.11.
 num_threadsint, optional
Number of threads to be used for OpenMP parallelization. If None (default) the value of OMP_NUM_THREADS environment variable is used if it is set, otherwise all available threads are used. If < 0 the maximal number of threads minus num_threads + 1 is used (enter 1 to use as many threads as possible). 0 raises an error. Only metrics using OpenMP will use this variable.
References#
[Garyfallidis15]Garyfallidis et al. “Robust and efficient linear registration of whitematter fascicles in the space of streamlines”, NeuroImage, 117, 124–140, 2015
[Garyfallidis14]Garyfallidis et al., “Direct nativespace fiber bundle alignment for group comparisons”, ISMRM, 2014.
[Garyfallidis17]Garyfallidis et al. Recognition of white matter bundles using local and global streamlinebased registration and clustering, Neuroimage, 2017.
 optimize(static, moving, mat=None)#
Find the minimum of the provided metric.
Parameters#
 staticstreamlines
Reference or fixed set of streamlines.
 movingstreamlines
Moving set of streamlines.
 matarray
Transformation (4, 4) matrix to start the registration.
mat
is applied to moving. Default value None which means that initial transformation will be generated by shifting the centers of moving and static sets of streamlines to the origin.
Returns#
map : StreamlineRegistrationMap
StreamlineRegistrationMap
#
 class dipy.align.streamlinear.StreamlineRegistrationMap(matopt, xopt, fopt, matopt_history, funcs, iterations)#
Bases:
object
 __init__(matopt, xopt, fopt, matopt_history, funcs, iterations)#
A map holding the optimum affine matrix and some other parameters of the optimization
Parameters#
 matrixarray,
4x4 affine matrix which transforms the moving to the static streamlines
 xoptarray,
1d array with the parameters of the transformation after centering
 foptfloat,
final value of the metric
 matrix_historyarray
All transformation matrices created during the optimization
 funcsint,
Number of function evaluations of the optimizer
 iterationsint
Number of iterations of the optimizer
JointStreamlineRegistrationMap
#
 class dipy.align.streamlinear.JointStreamlineRegistrationMap(xopt, fopt, matopt_history, funcs, iterations)#
Bases:
object
 __init__(xopt, fopt, matopt_history, funcs, iterations)#
A map holding the optimum affine matrices for halfway streamline linear registration and some other parameters of the optimization.
xopt is optimized by StreamlineLinearRegistration using the JointBundleMinDistanceMetric. In that case the mat argument of the optimize method needs to be np.eye(4) to avoid streamline centering.
This constructor derives and stores the transformations to move both static and moving bundles to the halfway space.
Parameters#
 xoptarray
1d array with the parameters of the transformation.
 foptfloat
Final value of the metric.
 matopt_historyarray
All transformation matrices created during the optimization.
 funcsint
Number of function evaluations of the optimizer.
 iterationsint
Number of iterations of the optimizer.
logger#
 dipy.align.streamlinear.logger()#
Instances of the Logger class represent a single logging channel. A “logging channel” indicates an area of an application. Exactly how an “area” is defined is up to the application developer. Since an application can have any number of areas, logging channels are identified by a unique string. Application areas can be nested (e.g. an area of “input processing” might include subareas “read CSV files”, “read XLS files” and “read Gnumeric files”). To cater for this natural nesting, channel names are organized into a namespace hierarchy where levels are separated by periods, much like the Java or Python package namespace. So in the instance given above, channel names might be “input” for the upper level, and “input.csv”, “input.xls” and “input.gnu” for the sublevels. There is no arbitrary limit to the depth of nesting.
bundle_sum_distance#
 dipy.align.streamlinear.bundle_sum_distance(t, static, moving, num_threads=None)#
MDF distance optimization function (SUM).
We minimize the distance between moving streamlines as they align with the static streamlines.
Parameters#
 tndarray
t is a vector of affine transformation parameters with size at least 6. If the size is 6, t is interpreted as translation + rotation. If the size is 7, t is interpreted as translation + rotation + isotropic scaling. If size is 12, t is interpreted as translation + rotation + scaling + shearing.
 staticlist
Static streamlines
 movinglist
Moving streamlines. These will be transformed to align with the static streamlines
 num_threadsint, optional
Number of threads. If 1 then all available threads will be used.
Returns#
cost: float
bundle_min_distance#
 dipy.align.streamlinear.bundle_min_distance(t, static, moving)#
MDFbased pairwise distance optimization function (MIN).
We minimize the distance between moving streamlines as they align with the static streamlines.
Parameters#
 tndarray
t is a vector of affine transformation parameters with size at least 6. If size is 6, t is interpreted as translation + rotation. If size is 7, t is interpreted as translation + rotation + isotropic scaling. If size is 12, t is interpreted as translation + rotation + scaling + shearing.
 staticlist
Static streamlines
 movinglist
Moving streamlines.
Returns#
cost: float
bundle_min_distance_fast#
 dipy.align.streamlinear.bundle_min_distance_fast(t, static, moving, block_size, num_threads=None)#
MDFbased pairwise distance optimization function (MIN).
We minimize the distance between moving streamlines as they align with the static streamlines.
Parameters#
 tarray
1D array. t is a vector of affine transformation parameters with size at least 6. If the size is 6, t is interpreted as translation + rotation. If the size is 7, t is interpreted as translation + rotation + isotropic scaling. If size is 12, t is interpreted as translation + rotation + scaling + shearing.
 staticarray
N*M x 3 array. All the points of the static streamlines. With order of streamlines intact. Where N is the number of streamlines and M is the number of points per streamline.
 movingarray
K*M x 3 array. All the points of the moving streamlines. With order of streamlines intact. Where K is the number of streamlines and M is the number of points per streamline.
 block_sizeint
Number of points per streamline. All streamlines in static and moving should have the same number of points M.
 num_threadsint, optional
Number of threads to be used for OpenMP parallelization. If None (default) the value of OMP_NUM_THREADS environment variable is used if it is set, otherwise all available threads are used. If < 0 the maximal number of threads minus num_threads + 1 is used (enter 1 to use as many threads as possible). 0 raises an error.
Returns#
cost: float
Notes#
This is a faster implementation of
bundle_min_distance
, which requires that all the points of each streamline are allocated into an ndarray (of shape N*M by 3, with N the number of points per streamline and M the number of streamlines). This can be done by calling dipy.tracking.streamlines.unlist_streamlines.
bundle_min_distance_asymmetric_fast#
 dipy.align.streamlinear.bundle_min_distance_asymmetric_fast(t, static, moving, block_size)#
MDFbased pairwise distance optimization function (MIN).
We minimize the distance between moving streamlines as they align with the static streamlines.
Parameters#
 tarray
1D array. t is a vector of affine transformation parameters with size at least 6. If the size is 6, t is interpreted as translation + rotation. If the size is 7, t is interpreted as translation + rotation + isotropic scaling. If size is 12, t is interpreted as translation + rotation + scaling + shearing.
 staticarray
N*M x 3 array. All the points of the static streamlines. With order of streamlines intact. Where N is the number of streamlines and M is the number of points per streamline.
 movingarray
K*M x 3 array. All the points of the moving streamlines. With order of streamlines intact. Where K is the number of streamlines and M is the number of points per streamline.
 block_sizeint
Number of points per streamline. All streamlines in static and moving should have the same number of points M.
Returns#
cost: float
remove_clusters_by_size#
 dipy.align.streamlinear.remove_clusters_by_size(clusters, min_size=0)#
progressive_slr#
 dipy.align.streamlinear.progressive_slr(static, moving, metric, x0, bounds, method='LBFGSB', verbose=False, num_threads=None)#
Progressive SLR.
This is an utility function that allows for example to do affine registration using Streamlinebased Linear Registration (SLR) [Garyfallidis15] by starting with translation first, then rigid, then similarity, scaling and finally affine.
Similarly, if for example, you want to perform rigid then you start with translation first. This progressive strategy can helps with finding the optimal parameters of the final transformation.
Parameters#
static : Streamlines moving : Streamlines metric : StreamlineDistanceMetric x0 : string
Could be any of ‘translation’, ‘rigid’, ‘similarity’, ‘scaling’, ‘affine’
 boundsarray
Boundaries of registration parameters. See variable DEFAULT_BOUNDS for example.
 methodstring
L_BFGS_B’ or ‘Powell’ optimizers can be used. Default is ‘L_BFGS_B’.
 verbosebool, optional.
If True, log messages. Default:
 num_threadsint, optional
Number of threads to be used for OpenMP parallelization. If None (default) the value of OMP_NUM_THREADS environment variable is used if it is set, otherwise all available threads are used. If < 0 the maximal number of threads minus num_threads + 1 is used (enter 1 to use as many threads as possible). 0 raises an error. Only metrics using OpenMP will use this variable.
References#
[Garyfallidis15]Garyfallidis et al. “Robust and efficient linear registration of whitematter fascicles in the space of streamlines”, NeuroImage, 117, 124–140, 2015
slr_with_qbx#
 dipy.align.streamlinear.slr_with_qbx(static, moving, x0='affine', rm_small_clusters=50, maxiter=100, select_random=None, verbose=False, greater_than=50, less_than=250, qbx_thr=(40, 30, 20, 15), nb_pts=20, progressive=True, rng=None, num_threads=None)#
Utility function for registering large tractograms.
For efficiency, we apply the registration on cluster centroids and remove small clusters.
Parameters#
static : Streamlines moving : Streamlines
 x0str, optional.
rigid, similarity or affine transformation model (default affine)
 rm_small_clustersint, optional
Remove clusters that have less than rm_small_clusters
 maxiterint, optional
Maximum number of iterations to perform.
 select_randomint, optional.
If not, None selects a random number of streamlines to apply clustering Default None.
 verbosebool, optional
If True, logs information about optimization. Default: False
 greater_thanint, optional
Keep streamlines that have length greater than this value (default 50)
 less_thanint, optional
Keep streamlines have length less than this value (default 250)
 qbx_thrvariable int
Thresholds for QuickBundlesX (default [40, 30, 20, 15])
 nb_ptsint, optional
Number of points for discretizing each streamline (default 20)
 progressiveboolean, optional
(default True)
 rngnp.random.Generator
If None creates random generator in function.
 num_threadsint, optional
Number of threads to be used for OpenMP parallelization. If None (default) the value of OMP_NUM_THREADS environment variable is used if it is set, otherwise all available threads are used. If < 0 the maximal number of threads minus num_threads + 1 is used (enter 1 to use as many threads as possible). 0 raises an error. Only metrics using OpenMP will use this variable.
Notes#
The order of operations is the following. First short or long streamlines are removed. Second, the tractogram or a random selection of the tractogram is clustered with QuickBundles. Then SLR [Garyfallidis15] is applied.
References#
[Garyfallidis15]Garyfallidis et al. “Robust and efficient linear
registration of whitematter fascicles in the space of streamlines”, NeuroImage, 117, 124–140, 2015 .. [Garyfallidis14] Garyfallidis et al., “Direct nativespace fiber
bundle alignment for group comparisons”, ISMRM, 2014.
[Garyfallidis17]Garyfallidis et al. Recognition of white matter
bundles using local and global streamlinebased registration and clustering, Neuroimage, 2017.
groupwise_slr#
 dipy.align.streamlinear.groupwise_slr(bundles, x0='affine', tol=0, max_iter=20, qbx_thr=[4], nb_pts=20, select_random=10000, verbose=False, rng=None)#
Function to perform unbiased groupwise bundle registration.
All bundles are moved to the same space by iteratively applying halfway streamline linear registration in pairs. With each iteration, bundles get closer to each other until the procedure converges and there is no more improvement.
Parameters#
 bundleslist
List with streamlines of the bundles to be registered.
 x0str, optional
rigid, similarity or affine transformation model. Default: affine.
 tolfloat, optional
Tolerance value to be used to assume convergence. Default: 0.
 max_iterint, optional
Maximum number of iterations. Depending on the number of bundles to be registered this may need to be larger. Default: 20.
 qbx_thrvariable int, optional
Thresholds for Quickbundles used for clustering streamlines and reduce computational time. If None, no clustering is performed. Higher values cluster streamlines into a smaller number of centroids. Default: [4].
 nb_ptsint, optional
Number of points for discretizing each streamline. Default: 20.
 select_randomint, optional
Maximum number of streamlines for each bundle. If None, all the streamlines are used. Default: 10000.
 verbosebool, optional
If True, logs information. Default: False.
 rngnp.random.Generator
If None, creates random generator in function. Default: None.
References#
[Garyfallidis15]Garyfallidis et al. “Robust and efficient linear
registration of whitematter fascicles in the space of streamlines”, NeuroImage, 117, 124–140, 2015 .. [Garyfallidis14] Garyfallidis et al., “Direct nativespace fiber
bundle alignment for group comparisons”, ISMRM, 2014.
[Garyfallidis17]Garyfallidis et al. Recognition of white matter
bundles using local and global streamlinebased registration and clustering, Neuroimage, 2017.
get_unique_pairs#
 dipy.align.streamlinear.get_unique_pairs(n_bundle, pairs=None)#
Make unique pairs from n_bundle bundles.
The function allows to input a previous pairs assignment so that the new pairs are different.
Parameters#
 n_bundleint
Number of bundles to be matched in pairs.
 pairsarray, optional
array containing the indexes of previous pairs.
compose_matrix44#
 dipy.align.streamlinear.compose_matrix44(t, dtype=<class 'numpy.float64'>)#
Compose a 4x4 transformation matrix.
Parameters#
 tndarray
This is a 1D vector of affine transformation parameters with size at least 3. If the size is 3, t is interpreted as translation. If the size is 6, t is interpreted as translation + rotation. If the size is 7, t is interpreted as translation + rotation + isotropic scaling. If the size is 9, t is interpreted as translation + rotation + anisotropic scaling. If size is 12, t is interpreted as translation + rotation + scaling + shearing.
Returns#
 Tndarray
Homogeneous transformation matrix of size 4x4.
decompose_matrix44#
 dipy.align.streamlinear.decompose_matrix44(mat, size=12)#
Given a 4x4 homogeneous matrix return the parameter vector.
Parameters#
 matarray
Homogeneous 4x4 transformation matrix
 sizeint
Size of the output vector. 3, for translation, 6 for rigid, 7 for similarity, 9 for scaling and 12 for affine. Default is 12.
Returns#
 tndarray
One dimensional ndarray of 3, 6, 7, 9 or 12 affine parameters.
average_bundle_length#
find_missing#
 dipy.align.streamwarp.find_missing(lst, cb)#
Find unmatched streamline indices in moving bundle.
Parameters#
 lstList
List of integers containing all the streamlines indices in moving bundle.
 cbList
List of integers containing streamline indices of the moving bundle that were not matched to any streamline in static bundle.
Returns#
 list
List containing unmatched streamlines from moving bundle
bundlewarp#
 dipy.align.streamwarp.bundlewarp(static, moving, dist=None, alpha=0.3, beta=20, max_iter=15, affine=True)#
Register two bundles using nonlinear method.
Parameters#
 staticStreamlines
Reference/fixed bundle
 movingStreamlines
Target bundle that will be moved/registered to match the static bundle
 distfloat, optional.
Precomputed distance matrix (default None)
 alphafloat, optional
Represents the tradeoff between regularizing the deformation and having points match very closely. Lower value of alpha means high deformations (default 0.3)
 betaint, optional
Represents the strength of the interaction between points Gaussian kernel size (default 20)
 max_iterint, optional
Maximum number of iterations for deformation process in mlCPD method (default 15)
 affineboolean, optional
If False, use rigid registration as starting point (default True)
Returns#
 deformed_bundleStreamlines
Nonlinearly moved bundle (warped bundle)
 moving_alignedStreamlines
Linearly moved bundle (affinely moved)
 distnp.ndarray
Float array containing distance between moving and static bundle
 matched_pairsnp.ndarray
Int array containing streamline correspondences between two bundles
 warpnp.ndarray
Nonlinear warp map generated by BundleWarp
References#
[Chandio2023]Chandio et al. “BundleWarp, streamlinebased nonlinear registration of white matter tracts.” bioRxiv (2023): 202301.
bundlewarp_vector_filed#
 dipy.align.streamwarp.bundlewarp_vector_filed(moving_aligned, deformed_bundle)#
Calculate vector fields.
Vector field computation as the difference between each streamline point in the deformed and linearly aligned bundles
Parameters#
 moving_alignedStreamlines
Linearly (affinely) moved bundle
 deformed_bundleStreamlines
Nonlinearly (warped) bundle
Returns#
 offsetsList
Vector field modules
 directionsList
Unitary vector directions
colors : List
bundlewarp_shape_analysis#
 dipy.align.streamwarp.bundlewarp_shape_analysis(moving_aligned, deformed_bundle, no_disks=10, plotting=False)#
Calculate bundle shape difference profile.
Bundle shape difference analysis using magnitude from BundleWarp displacements and BUAN.
Depending on the number of points of a streamline, and the number of segments requested, multiple points may be considered for the computation of a given segment; a segment may contain information from a single point; or some segments may not contain information from any points. In the latter case, the segment will contain an
np.nan
value. The pointtosegment mapping is defined by theassignment_map()
: for each segment index, the point information of the matching index positions, as returned byassignment_map()
, are considered for the computation.Parameters#
 moving_alignedStreamlines
Linearly (affinely) moved bundle
 deformed_bundleStreamlines
Nonlinearly (warped) moved bundle
 no_disksint, optional
Number of segments to be created along the length of the bundle
 plottingBoolean, optional
Plot bundle shape profile
Returns#
 shape_profilenp.ndarray
Float array containing bundlewarp displacement magnitudes along the length of the bundle
 stdvnp.ndarray
Float array containing standard deviations
compute_energy_ssd_2d#
 dipy.align.sumsqdiff.compute_energy_ssd_2d(delta_field)#
Sum of squared differences between two 2D images
Computes the Sum of Squared Differences between the static and moving image. Those differences are given by delta_field
Parameters#
 delta_fieldarray, shape (R, C)
the difference between the static and moving image (the ‘derivative w.r.t. time’ in the optical flow model)
Returns#
 energyfloat
the SSD energy at this iteration
Notes#
The numeric value of the energy is used only to detect convergence. This function returns only the energy corresponding to the data term (excluding the energy corresponding to the regularization term) because the GreedySyN algorithm is an unconstrained gradient descent algorithm in the space of diffeomorphisms: in each iteration it makes a step along the negative smoothed gradient –of the data term– and then makes sure the resulting diffeomorphisms are invertible using an explicit inversion algorithm. Since it is not clear how to reflect the energy corresponding to this reprojection to the space of diffeomorphisms, a more precise energy computation including the regularization term is useless. Instead, convergence is checked considering the dataterm energy only and detecting oscilations in the energy profile.
compute_energy_ssd_3d#
 dipy.align.sumsqdiff.compute_energy_ssd_3d(delta_field)#
Sum of squared differences between two 3D volumes
Computes the Sum of Squared Differences between the static and moving volume Those differences are given by delta_field
Parameters#
 delta_fieldarray, shape (R, C)
the difference between the static and moving image (the ‘derivative w.r.t. time’ in the optical flow model)
Returns#
 energyfloat
the SSD energy at this iteration
Notes#
The numeric value of the energy is used only to detect convergence. This function returns only the energy corresponding to the data term (excluding the energy corresponding to the regularization term) because the GreedySyN algorithm is an unconstrained gradient descent algorithm in the space of diffeomorphisms: in each iteration it makes a step along the negative smoothed gradient –of the data term– and then makes sure the resulting diffeomorphisms are invertible using an explicit inversion algorithm. Since it is not clear how to reflect the energy corresponding to this reprojection to the space of diffeomorphisms, a more precise energy computation including the regularization term is useless. Instead, convergence is checked considering the dataterm energy only and detecting oscilations in the energy profile.
compute_residual_displacement_field_ssd_2d#
 dipy.align.sumsqdiff.compute_residual_displacement_field_ssd_2d(delta_field, sigmasq_field, gradient_field, target, lambda_param, d, residual)#
The residual displacement field to be fit on the next iteration
Computes the residual displacement field corresponding to the current displacement field in the Multiresolution GaussSeidel solver proposed by Bruhn and Weickert [Bruhn05].
Parameters#
 delta_fieldarray, shape (R, C)
the difference between the static and moving image (the ‘derivative w.r.t. time’ in the optical flow model)
 sigmasq_fieldarray, shape (R, C)
the variance of the gray level value at each voxel, according to the EM model (for SSD, it is 1 for all voxels). Inf and 0 values are processed specially to support infinite and zero variance.
 gradient_fieldarray, shape (R, C, 2)
the gradient of the moving image
 targetarray, shape (R, C, 2)
righthand side of the linear system to be solved in the Weickert’s multiresolution algorithm
 lambda_paramfloat
smoothness parameter in the objective function
 darray, shape (R, C, 2)
the current displacement field to compute the residual from
 residualarray, shape (R, C, 2)
the displacement field to put the residual to
Returns#
 residualarray, shape (R, C, 2)
the residual displacement field. If residual was None a input, then a new field is returned, otherwise the same array is returned
References#
 [Bruhn05] Andres Bruhn and Joachim Weickert, “Towards ultimate motion
estimation: combining highest accuracy with realtime performance”, 10th IEEE International Conference on Computer Vision, 2005. ICCV 2005.
compute_residual_displacement_field_ssd_3d#
 dipy.align.sumsqdiff.compute_residual_displacement_field_ssd_3d(delta_field, sigmasq_field, gradient_field, target, lambda_param, disp, residual)#
The residual displacement field to be fit on the next iteration
Computes the residual displacement field corresponding to the current displacement field (given by ‘disp’) in the Multiresolution GaussSeidel solver proposed by Bruhn and Weickert [Bruhn].
Parameters#
 delta_fieldarray, shape (S, R, C)
the difference between the static and moving image (the ‘derivative w.r.t. time’ in the optical flow model)
 sigmasq_fieldarray, shape (S, R, C)
the variance of the gray level value at each voxel, according to the EM model (for SSD, it is 1 for all voxels). Inf and 0 values are processed specially to support infinite and zero variance.
 gradient_fieldarray, shape (S, R, C, 3)
the gradient of the moving image
 targetarray, shape (S, R, C, 3)
righthand side of the linear system to be solved in the Weickert’s multiresolution algorithm
 lambda_paramfloat
smoothness parameter in the objective function
 disparray, shape (S, R, C, 3)
the current displacement field to compute the residual from
 residualarray, shape (S, R, C, 3)
the displacement field to put the residual to
Returns#
 residualarray, shape (S, R, C, 3)
the residual displacement field. If residual was None a input, then a new field is returned, otherwise the same array is returned
References#
 [Bruhn05] Andres Bruhn and Joachim Weickert, “Towards ultimate motion
estimation: combining highest accuracy with realtime performance”, 10th IEEE International Conference on Computer Vision, 2005. ICCV 2005.
compute_ssd_demons_step_2d#
 dipy.align.sumsqdiff.compute_ssd_demons_step_2d(delta_field, gradient_moving, sigma_sq_x, out)#
Demons step for 2D SSDdriven registration Computes the demons step for SSDdriven registration ( eq. 4 in [Bruhn05] ) Parameters ——— delta_field : array, shape (R, C) the difference between the static and moving image (the ‘derivative w.r.t. time’ in the optical flow model) gradient_field : array, shape (R, C, 2) the gradient of the moving image sigma_sq_x : float parameter controlling the amount of regularization. It corresponds to \(\sigma_x^2\) in algorithm 1 of Vercauteren et al.[Vercauteren09] out : array, shape (R, C, 2) if None, a new array will be created to store the demons step. Otherwise the provided array will be used. Returns —— demons_step : array, shape (R, C, 2) the demons step to be applied for updating the current displacement field energy : float the current ssd energy (before applying the returned demons_step) References ——— [Bruhn05] Andres Bruhn and Joachim Weickert, “Towards ultimate motion estimation: combining highest accuracy with realtime performance”, 10th IEEE International Conference on Computer Vision, 2005. ICCV 2005. [Vercauteren09] Vercauteren, T., Pennec, X., Perchant, A., & Ayache, N. (2009). Diffeomorphic demons: efficient nonparametric image registration. NeuroImage, 45(1 Suppl), S6172. doi:10.1016/j.neuroimage.2008.10.040
compute_ssd_demons_step_3d#
 dipy.align.sumsqdiff.compute_ssd_demons_step_3d(delta_field, gradient_moving, sigma_sq_x, out)#
Demons step for 3D SSDdriven registration Computes the demons step for SSDdriven registration ( eq. 4 in [Bruhn05] ) Parameters ——— delta_field : array, shape (S, R, C) the difference between the static and moving image (the ‘derivative w.r.t. time’ in the optical flow model) gradient_field : array, shape (S, R, C, 2) the gradient of the moving image sigma_sq_x : float parameter controlling the amount of regularization. It corresponds to \(\sigma_x^2\) in algorithm 1 of Vercauteren et al.[Vercauteren09] out : array, shape (S, R, C, 2) if None, a new array will be created to store the demons step. Otherwise the provided array will be used. Returns —— demons_step : array, shape (S, R, C, 3) the demons step to be applied for updating the current displacement field energy : float the current ssd energy (before applying the returned demons_step) References ——— [Bruhn05] Andres Bruhn and Joachim Weickert, “Towards ultimate motion estimation: combining highest accuracy with realtime performance”, 10th IEEE International Conference on Computer Vision, 2005. ICCV 2005. [Vercauteren09] Vercauteren, T., Pennec, X., Perchant, A., & Ayache, N. (2009). Diffeomorphic demons: efficient nonparametric image registration. NeuroImage, 45(1 Suppl), S6172. doi:10.1016/j.neuroimage.2008.10.040
iterate_residual_displacement_field_ssd_2d#
 dipy.align.sumsqdiff.iterate_residual_displacement_field_ssd_2d(delta_field, sigmasq_field, grad, target, lambda_param, displacement_field)#
One iteration of a large linear system solver for 2D SSD registration
Performs one iteration at one level of the Multiresolution GaussSeidel solver proposed by Bruhn and Weickert [Bruhn05].
Parameters#
 delta_fieldarray, shape (R, C)
the difference between the static and moving image (the ‘derivative w.r.t. time’ in the optical flow model)
 sigmasq_fieldarray, shape (R, C)
the variance of the gray level value at each voxel, according to the EM model (for SSD, it is 1 for all voxels). Inf and 0 values are processed specially to support infinite and zero variance.
 gradarray, shape (R, C, 2)
the gradient of the moving image
 targetarray, shape (R, C, 2)
righthand side of the linear system to be solved in the Weickert’s multiresolution algorithm
 lambda_paramfloat
smoothness parameter of the objective function
 displacement_fieldarray, shape (R, C, 2)
current displacement field to start the iteration from
Returns#
 max_displacementfloat
the norm of the maximum change in the displacement field after the iteration
References#
 [Bruhn05] Andres Bruhn and Joachim Weickert, “Towards ultimate motion
estimation: combining highest accuracy with realtime performance”, 10th IEEE International Conference on Computer Vision, 2005. ICCV 2005.
iterate_residual_displacement_field_ssd_3d#
 dipy.align.sumsqdiff.iterate_residual_displacement_field_ssd_3d(delta_field, sigmasq_field, grad, target, lambda_param, disp)#
One iteration of a large linear system solver for 3D SSD registration
Performs one iteration at one level of the Multiresolution GaussSeidel solver proposed by Bruhn and Weickert [Bruhn05].
Parameters#
 delta_fieldarray, shape (S, R, C)
the difference between the static and moving image (the ‘derivative w.r.t. time’ in the optical flow model)
 sigmasq_fieldarray, shape (S, R, C)
the variance of the gray level value at each voxel, according to the EM model (for SSD, it is 1 for all voxels). Inf and 0 values are processed specially to support infinite and zero variance.
 gradarray, shape (S, R, C, 3)
the gradient of the moving image
 targetarray, shape (S, R, C, 3)
righthand side of the linear system to be solved in the Weickert’s multiresolution algorithm
 lambda_paramfloat
smoothness parameter of the objective function
 disparray, shape (S, R, C, 3)
the displacement field to start the optimization from
Returns#
 max_displacementfloat
the norm of the maximum change in the displacement field after the iteration
References#
 [Bruhn05] Andres Bruhn and Joachim Weickert, “Towards ultimate motion
estimation: combining highest accuracy with realtime performance”, 10th IEEE International Conference on Computer Vision, 2005. ICCV 2005.
solve_2d_symmetric_positive_definite#
 dipy.align.sumsqdiff.solve_2d_symmetric_positive_definite(A, y, det)#
Solves a 2variable symmetric positivedefinite linear system Solves the symmetric positivedefinite linear system \(Mx = y\) given by:: M = [[A[0], A[1]], [A[1], A[2]]] Parameters ——— A : array, shape (3,) the array containing the entries of the symmetric 2x2 matrix y : array, shape (2,) righthand side of the system to be solved Returns —— out : array, shape (2,) the array the output will be stored in
solve_3d_symmetric_positive_definite#
 dipy.align.sumsqdiff.solve_3d_symmetric_positive_definite(g, y, tau)#
Solves a 3variable symmetric positivedefinite linear system Solves the symmetric semipositivedefinite linear system \(Mx = y\) given by \(M = (g g^{T} + \tau I)\). Parameters ——— g : array, shape (3,) the vector in the outer product above y : array, shape (3,) righthand side of the system to be solved tau : double \(\tau\) in \(M = (g g^{T} + \tau I)\) Returns —— out : array, shape (3,) the array the output will be stored in is_singular : int 1 if M is singular, otherwise 0
AffineTransform2D
#
AffineTransform3D
#
RigidIsoScalingTransform2D
#
 class dipy.align.transforms.RigidIsoScalingTransform2D#
Bases:
Transform
 __init__()#
Rigid isoscaling transform in 2D.
(rotation + translation + scaling)
The parameter vector theta of length 4 is interpreted as follows: theta[0] : rotation angle theta[1] : translation along the x axis theta[2] : translation along the y axis theta[3] : isotropic scaling
RigidIsoScalingTransform3D
#
 class dipy.align.transforms.RigidIsoScalingTransform3D#
Bases:
Transform
 __init__()#
Rigid isoscaling transform in 3D.
(rotation + translation + scaling) The parameter vector theta of length 7 is interpreted as follows: theta[0] : rotation about the x axis theta[1] : rotation about the y axis theta[2] : rotation about the z axis theta[3] : translation along the x axis theta[4] : translation along the y axis theta[5] : translation along the z axis theta[6] : isotropic scaling
RigidScalingTransform2D
#
 class dipy.align.transforms.RigidScalingTransform2D#
Bases:
Transform
 __init__()#
Rigid scaling transform in 2D.
(rotation + translation + scaling)
The parameter vector theta of length 5 is interpreted as follows: theta[0] : rotation angle theta[1] : translation along the x axis theta[2] : translation along the y axis theta[3] : scaling along the x axis theta[4] : scaling along the y axis
RigidScalingTransform3D
#
 class dipy.align.transforms.RigidScalingTransform3D#
Bases:
Transform
 __init__()#
Rigid Scaling transform in 3D (rotation + translation + scaling).
The parameter vector theta of length 9 is interpreted as follows: theta[0] : rotation about the x axis theta[1] : rotation about the y axis theta[2] : rotation about the z axis theta[3] : translation along the x axis theta[4] : translation along the y axis theta[5] : translation along the z axis theta[6] : scaling in the x axis theta[7] : scaling in the y axis theta[8] : scaling in the z axis
RigidTransform2D
#
RigidTransform3D
#
 class dipy.align.transforms.RigidTransform3D#
Bases:
Transform
 __init__()#
Rigid transform in 3D (rotation + translation) The parameter vector theta of length 6 is interpreted as follows: theta[0] : rotation about the x axis theta[1] : rotation about the y axis theta[2] : rotation about the z axis theta[3] : translation along the x axis theta[4] : translation along the y axis theta[5] : translation along the z axis
RotationTransform2D
#
RotationTransform3D
#
ScalingTransform2D
#
ScalingTransform3D
#
Transform
#
 class dipy.align.transforms.Transform#
Bases:
object
Base class (contract) for all transforms for affine image registration Each transform must define the following (fast, nogil) methods:
_jacobian(theta, x, J): receives a parameter vector theta, a point in x, and a matrix J with shape (dim, len(theta)). It must writes in J, the Jacobian of the transform with parameters theta evaluated at x.
_get_identity_parameters(theta): receives a vector theta whose length is the number of parameters of the transform and sets in theta the values that define the identity transform.
_param_to_matrix(theta, T): receives a parameter vector theta, and a matrix T of shape (dim + 1, dim + 1) and writes in T the matrix representation of the transform with parameters theta
This base class defines the (slow, convenient) python wrappers for each of the above functions, which also do parameter checking and raise a ValueError in case the provided parameters are invalid.
 __init__()#
 get_dim()#
 get_identity_parameters()#
Parameter values corresponding to the identity transform
Returns#
 thetaarray, shape (n,)
the n parameter values corresponding to the identity transform
 get_number_of_parameters()#
 jacobian(theta, x)#
Jacobian function of this transform
Parameters#
 thetaarray, shape (n,)
vector containing the n parameters of this transform
 xarray, shape (dim,)
vector containing the point where the Jacobian must be evaluated
Returns#
 Jarray, shape (dim, n)
Jacobian matrix of the transform with parameters theta at point x
TranslationTransform2D
#
TranslationTransform3D
#
compose_vector_fields_2d#
 dipy.align.vector_fields.compose_vector_fields_2d(d1, d2, premult_index, premult_disp, time_scaling, comp)#
Computes the composition of two 2D displacement fields
Computes the composition of the two 2D displacements d1 and d2. The evaluation of d2 at nonlattice points is computed using trilinear interpolation. The actual composition is computed as:
comp[i] = d1[i] + t * d2[ A * i + B * d1[i] ]
where t = time_scaling, A = premult_index and B=premult_disp and i denotes the voxel coordinates of a voxel in d1’s grid. Using this parameters it is possible to compose vector fields with arbitrary discretizations: let R and S be the voxeltospace transformation associated to d1 and d2, respectively then the composition at a voxel with coordinates i in d1’s grid is given by:
comp[i] = d1[i] + R*i + d2[Sinv*(R*i + d1[i])]  R*i
(the R*i terms cancel each other) where Sinv = S^{1} we can then define A = Sinv * R and B = Sinv to compute the composition using this function.
Parameters#
 d1array, shape (R, C, 2)
first displacement field to be applied. R, C are the number of rows and columns of the displacement field, respectively.
 d2array, shape (R’, C’, 2)
second displacement field to be applied. R’, C’ are the number of rows and columns of the displacement field, respectively.
 premult_indexarray, shape (3, 3)
the matrix A in the explanation above
 premult_disparray, shape (3, 3)
the matrix B in the explanation above
 time_scalingfloat
this corresponds to the time scaling ‘t’ in the above explanation
 comparray, shape (R, C, 2)
the buffer to write the composition to. If None, the buffer is created internally
Returns#
 comparray, shape (R, C, 2), same dimension as d1
on output, this array will contain the composition of the two fields
 statsarray, shape (3,)
on output, this array will contain three statistics of the vector norms of the composition (maximum, mean, standard_deviation)
compose_vector_fields_3d#
 dipy.align.vector_fields.compose_vector_fields_3d(d1, d2, premult_index, premult_disp, time_scaling, comp)#
Computes the composition of two 3D displacement fields
Computes the composition of the two 3D displacements d1 and d2. The evaluation of d2 at nonlattice points is computed using trilinear interpolation. The actual composition is computed as:
comp[i] = d1[i] + t * d2[ A * i + B * d1[i] ]
where t = time_scaling, A = premult_index and B=premult_disp and i denotes the voxel coordinates of a voxel in d1’s grid. Using this parameters it is possible to compose vector fields with arbitrary discretization: let R and S be the voxeltospace transformation associated to d1 and d2, respectively then the composition at a voxel with coordinates i in d1’s grid is given by:
comp[i] = d1[i] + R*i + d2[Sinv*(R*i + d1[i])]  R*i
(the R*i terms cancel each other) where Sinv = S^{1} we can then define A = Sinv * R and B = Sinv to compute the composition using this function.
Parameters#
 d1array, shape (S, R, C, 3)
first displacement field to be applied. S, R, C are the number of slices, rows and columns of the displacement field, respectively.
 d2array, shape (S’, R’, C’, 3)
second displacement field to be applied. R’, C’ are the number of rows and columns of the displacement field, respectively.
 premult_indexarray, shape (4, 4)
the matrix A in the explanation above
 premult_disparray, shape (4, 4)
the matrix B in the explanation above
 time_scalingfloat
this corresponds to the time scaling ‘t’ in the above explanation
 comparray, shape (S, R, C, 3), same dimension as d1
the buffer to write the composition to. If None, the buffer will be created internally
Returns#
 comparray, shape (S, R, C, 3), same dimension as d1
on output, this array will contain the composition of the two fields
 statsarray, shape (3,)
on output, this array will contain three statistics of the vector norms of the composition (maximum, mean, standard_deviation)
Notes#
If d1[s,r,c] lies outside the domain of d2, then comp[s,r,c] will contain a zero vector.
create_circle#
 dipy.align.vector_fields.create_circle(nrows, ncols, radius)#
Create a binary 2D image where pixel values are 1 iff their distance to the center of the image is less than or equal to radius.
Parameters#
 nrowsint
number of rows of the resulting image
 ncolsint
number of columns of the resulting image
 radiusint
the radius of the circle
Returns#
 carray, shape (nrows, ncols)
the binary image of the circle with the requested dimensions
create_harmonic_fields_2d#
 dipy.align.vector_fields.create_harmonic_fields_2d(nrows, ncols, b, m)#
Creates an invertible 2D displacement field
Creates the invertible displacement fields used in Chen et al. eqs. 9 and 10 [1]
Parameters#
 nrowsint
number of rows in the resulting harmonic field
 ncolsint
number of columns in the resulting harmonic field
 b, mfloat
parameters of the harmonic field (as in [1]). To understand the effect of these parameters, please consider plotting a deformed image (a circle or a grid) under the deformation field, or see examples in [1]
Returns#
 darray, shape (nrows, ncols, 2)
the harmonic displacement field
 invarray, shape (nrows, ncols, 2)
the analytical inverse of the harmonic displacement field
 [1] Chen, M., Lu, W., Chen, Q., Ruchala, K. J., & Olivera, G. H. (2008).
A simple fixedpoint approach to invert a deformation field. Medical Physics, 35(1), 81. doi:10.1118/1.2816107
create_harmonic_fields_3d#
 dipy.align.vector_fields.create_harmonic_fields_3d(nslices, nrows, ncols, b, m)#
Creates an invertible 3D displacement field
Creates the invertible displacement fields used in Chen et al. eqs. 9 and 10 [1] computing the angle theta along zslides.
Parameters#
 nslicesint
number of slices in the resulting harmonic field
 nrowsint
number of rows in the resulting harmonic field
 ncolsint
number of columns in the resulting harmonic field
 b, ffloat
parameters of the harmonic field (as in [1]). To understand the effect of these parameters, please consider plotting a deformed image (e.g. a circle or a grid) under the deformation field, or see examples in [1]
Returns#
 darray, shape (nslices, nrows, ncols, 3)
the harmonic displacement field
 invarray, shape (nslices, nrows, ncols, 3)
the analytical inverse of the harmonic displacement field
 [1] Chen, M., Lu, W., Chen, Q., Ruchala, K. J., & Olivera, G. H. (2008).
A simple fixedpoint approach to invert a deformation field. Medical Physics, 35(1), 81. doi:10.1118/1.2816107
create_random_displacement_2d#
 dipy.align.vector_fields.create_random_displacement_2d(from_shape, from_grid2world, to_shape, to_grid2world, rng=None)#
Creates a random 2D displacement ‘exactly’ mapping points of two grids
Creates a random 2D displacement field mapping points of an input discrete domain (with dimensions given by from_shape) to points of an output discrete domain (with shape given by to_shape). The affine matrices bringing discrete coordinates to physical space are given by from_grid2world (for the displacement field discretization) and to_grid2world (for the target discretization). Since this function is intended to be used for testing, voxels in the input domain will never be assigned to boundary voxels on the output domain.
Parameters#
 from_shapearray, shape (2,)
the grid shape where the displacement field will be defined on.
 from_grid2worldarray, shape (3,3)
the gridtospace transformation of the displacement field
 to_shapearray, shape (2,)
the grid shape where the deformation field will map the input grid to.
 to_grid2worldarray, shape (3,3)
the gridtospace transformation of the mapped grid
 rngnumpy.rnadom.Generator class, optional
Numpy’s random generator for setting seed values when needed. Default is None.
Returns#
 outputarray, shape = from_shape
the random displacement field in the physical domain
 int_fieldarray, shape = from_shape
the assignment of each point in the input grid to the target grid
create_random_displacement_3d#
 dipy.align.vector_fields.create_random_displacement_3d(from_shape, from_grid2world, to_shape, to_grid2world, rng=None)#
Creates a random 3D displacement ‘exactly’ mapping points of two grids
Creates a random 3D displacement field mapping points of an input discrete domain (with dimensions given by from_shape) to points of an output discrete domain (with shape given by to_shape). The affine matrices bringing discrete coordinates to physical space are given by from_grid2world (for the displacement field discretization) and to_grid2world (for the target discretization). Since this function is intended to be used for testing, voxels in the input domain will never be assigned to boundary voxels on the output domain.
Parameters#
 from_shapearray, shape (3,)
the grid shape where the displacement field will be defined on.
 from_grid2worldarray, shape (4,4)
the gridtospace transformation of the displacement field
 to_shapearray, shape (3,)
the grid shape where the deformation field will map the input grid to.
 to_grid2worldarray, shape (4,4)
the gridtospace transformation of the mapped grid
 rngnumpy.rnadom.Generator class, optional
Numpy’s random generator for setting seed values when needed. Default is None.
Returns#
 outputarray, shape = from_shape
the random displacement field in the physical domain
 int_fieldarray, shape = from_shape
the assignment of each point in the input grid to the target grid
create_sphere#
 dipy.align.vector_fields.create_sphere(nslices, nrows, ncols, radius)#
Create a binary 3D image where voxel values are 1 iff their distance to the center of the image is less than or equal to radius.
Parameters#
 nslicesint
number if slices of the resulting image
 nrowsint
number of rows of the resulting image
 ncolsint
number of columns of the resulting image
 radiusint
the radius of the sphere
Returns#
 carray, shape (nslices, nrows, ncols)
the binary image of the sphere with the requested dimensions
downsample_displacement_field_2d#
 dipy.align.vector_fields.downsample_displacement_field_2d(field)#
Downsamples the 2D input vector field by a factor of 2 Downsamples the input vector field by a factor of 2. The value at each pixel of the resulting field is the average of its surrounding pixels in the original field.
Parameters#
 fieldarray, shape (R, C)
the vector field to be downsampled
Returns#
 downarray, shape (R’, C’)
the downsampled displacement field, where R’= ceil(R/2), C’=ceil(C/2),
downsample_displacement_field_3d#
 dipy.align.vector_fields.downsample_displacement_field_3d(field)#
Downsamples the input 3D vector field by a factor of 2
Downsamples the input vector field by a factor of 2. This operation is equivalent to dividing the input image into 2x2x2 cubes and averaging the 8 vectors. The resulting field consists of these average vectors.
Parameters#
 fieldarray, shape (S, R, C)
the vector field to be downsampled
Returns#
 downarray, shape (S’, R’, C’)
the downsampled displacement field, where S’ = ceil(S/2), R’= ceil(R/2), C’=ceil(C/2)
downsample_scalar_field_2d#
 dipy.align.vector_fields.downsample_scalar_field_2d(field)#
Downsamples the input 2D image by a factor of 2
Downsamples the input image by a factor of 2. The value at each pixel of the resulting image is the average of its surrounding pixels in the original image.
Parameters#
 fieldarray, shape (R, C)
the image to be downsampled
Returns#
 downarray, shape (R’, C’)
the downsampled displacement field, where R’= ceil(R/2), C’=ceil(C/2)
downsample_scalar_field_3d#
 dipy.align.vector_fields.downsample_scalar_field_3d(field)#
Downsamples the input volume by a factor of 2
Downsamples the input volume by a factor of 2. The value at each voxel of the resulting volume is the average of its surrounding voxels in the original volume.
Parameters#
 fieldarray, shape (S, R, C)
the volume to be downsampled
Returns#
 downarray, shape (S’, R’, C’)
the downsampled displacement field, where S’ = ceil(S/2), R’= ceil(R/2), C’=ceil(C/2)
gradient#
 dipy.align.vector_fields.gradient(img, img_world2grid, img_spacing, out_shape, out_grid2world)#
Gradient of an image in physical space
Parameters#
 img2D or 3D array, shape (R, C) or (S, R, C)
the input image whose gradient will be computed
 img_world2gridarray, shape (dim+1, dim+1)
the spacetogrid transform matrix associated to img
 img_spacingarray, shape (dim,)
the spacing between voxels (voxel size along each axis) of the input image
 out_shapearray, shape (dim,)
the number of (slices), rows and columns of the sampling grid
 out_grid2worldarray, shape (dim+1, dim+1)
the gridtospace transform associated to the sampling grid
Returns#
 outarray, shape (R’, C’, 2) or (S’, R’, C’, 3)
the buffer in which to store the image gradient, where (S’), R’, C’ are given by out_shape
invert_vector_field_fixed_point_2d#
 dipy.align.vector_fields.invert_vector_field_fixed_point_2d(d, d_world2grid, spacing, max_iter, tolerance, start=None)#
Computes the inverse of a 2D displacement fields
Computes the inverse of the given 2D displacement field d using the fixedpoint algorithm [1].
 [1] Chen, M., Lu, W., Chen, Q., Ruchala, K. J., & Olivera, G. H. (2008).
A simple fixedpoint approach to invert a deformation field. Medical Physics, 35(1), 81. doi:10.1118/1.2816107
Parameters#
 darray, shape (R, C, 2)
the 2D displacement field to be inverted
 d_world2gridarray, shape (3, 3)
the spacetogrid transformation associated to the displacement field d (transforming physical space coordinates to voxel coordinates of the displacement field grid)
 spacing :array, shape (2,)
the spacing between voxels (voxel size along each axis)
 max_iterint
maximum number of iterations to be performed
 tolerancefloat
maximum tolerated inversion error
 startarray, shape (R, C)
an approximation to the inverse displacement field (if no approximation is available, None can be provided and the start displacement field will be zero)
Returns#
 parray, shape (R, C, 2)
the inverse displacement field
Notes#
We assume that the displacement field is an endomorphism so that the shape and voxeltospace transformation of the inverse field’s discretization is the same as those of the input displacement field. The ‘inversion error’ at iteration t is defined as the mean norm of the displacement vectors of the input displacement field composed with the inverse at iteration t.
invert_vector_field_fixed_point_3d#
 dipy.align.vector_fields.invert_vector_field_fixed_point_3d(d, d_world2grid, spacing, max_iter, tol, start=None)#
Computes the inverse of a 3D displacement fields
Computes the inverse of the given 3D displacement field d using the fixedpoint algorithm [1].
 [1] Chen, M., Lu, W., Chen, Q., Ruchala, K. J., & Olivera, G. H. (2008).
A simple fixedpoint approach to invert a deformation field. Medical Physics, 35(1), 81. doi:10.1118/1.2816107
Parameters#
 darray, shape (S, R, C, 3)
the 3D displacement field to be inverted
 d_world2gridarray, shape (4, 4)
the spacetogrid transformation associated to the displacement field d (transforming physical space coordinates to voxel coordinates of the displacement field grid)
 spacing :array, shape (3,)
the spacing between voxels (voxel size along each axis)
 max_iterint
maximum number of iterations to be performed
 tolfloat
maximum tolerated inversion error
 startarray, shape (S, R, C)
an approximation to the inverse displacement field (if no approximation is available, None can be provided and the start displacement field will be zero)
Returns#
 parray, shape (S, R, C, 3)
the inverse displacement field
Notes#
We assume that the displacement field is an endomorphism so that the shape and voxeltospace transformation of the inverse field’s discretization is the same as those of the input displacement field. The ‘inversion error’ at iteration t is defined as the mean norm of the displacement vectors of the input displacement field composed with the inverse at iteration t.
is_valid_affine#
 dipy.align.vector_fields.is_valid_affine(M, dim)#
reorient_vector_field_2d#
 dipy.align.vector_fields.reorient_vector_field_2d(d, affine)#
Linearly transforms all vectors of a 2D displacement field
Modifies the input displacement field by multiplying each displacement vector by the given matrix. Note that the elements of the displacement field are vectors, not points, so their last homogeneous coordinate is zero, not one, and therefore the translation component of the affine transform will not have any effect on them.
Parameters#
 darray, shape (R, C, 2)
the displacement field to be reoriented
 affine: array, shape (3, 3)
the matrix to be applied
reorient_vector_field_3d#
 dipy.align.vector_fields.reorient_vector_field_3d(d, affine)#
Linearly transforms all vectors of a 3D displacement field
Modifies the input displacement field by multiplying each displacement vector by the given matrix. Note that the elements of the displacement field are vectors, not points, so their last homogeneous coordinate is zero, not one, and therefore the translation component of the affine transform will not have any effect on them.
Parameters#
 darray, shape (S, R, C, 3)
the displacement field to be reoriented
 affinearray, shape (4, 4)
the matrix to be applied
resample_displacement_field_2d#
 dipy.align.vector_fields.resample_displacement_field_2d(field, factors, out_shape)#
Resamples a 2D vector field to a custom target shape
Resamples the given 2D displacement field on a grid of the requested shape, using the given scale factors. More precisely, the resulting displacement field at each grid cell i is given by
D[i] = field[Diag(factors) * i]
Parameters#
 factorsarray, shape (2,)
the scaling factors mapping (integer) grid coordinates in the resampled grid to (floating point) grid coordinates in the original grid
 out_shapearray, shape (2,)
the desired shape of the resulting grid
Returns#
 expandedarray, shape = out_shape + (2, )
the resampled displacement field
resample_displacement_field_3d#
 dipy.align.vector_fields.resample_displacement_field_3d(field, factors, out_shape)#
Resamples a 3D vector field to a custom target shape
Resamples the given 3D displacement field on a grid of the requested shape, using the given scale factors. More precisely, the resulting displacement field at each grid cell i is given by
D[i] = field[Diag(factors) * i]
Parameters#
 factorsarray, shape (3,)
the scaling factors mapping (integer) grid coordinates in the resampled grid to (floating point) grid coordinates in the original grid
 out_shapearray, shape (3,)
the desired shape of the resulting grid
Returns#
 expandedarray, shape = out_shape + (3, )
the resampled displacement field
simplify_warp_function_2d#
 dipy.align.vector_fields.simplify_warp_function_2d(d, affine_idx_in, affine_idx_out, affine_disp, out_shape)#
Simplifies a nonlinear warping function combined with an affine transform
Modifies the given deformation field by incorporating into it a an affine transformation and voxeltospace transforms associated to the discretization of its domain and codomain. The resulting transformation may be regarded as operating on the image spaces given by the domain and codomain discretization. More precisely, the resulting transform is of the form:
T[i] = W * d[U * i] + V * i
Where U = affine_idx_in, V = affine_idx_out, W = affine_disp.
Parameters#
 darray, shape (R’, C’, 2)
the nonlinear part of the transformation (displacement field)
 affine_idx_inarray, shape (3, 3)
the matrix U in eq. (1) above
 affine_idx_outarray, shape (3, 3)
the matrix V in eq. (1) above
 affine_disparray, shape (3, 3)
the matrix W in eq. (1) above
 out_shapearray, shape (2,)
the number of rows and columns of the sampling grid
Returns#
 outarray, shape = out_shape
the deformation field out associated with T in eq. (1) such that: T[i] = i + out[i]
Notes#
Both the direct and inverse transforms of a DiffeomorphicMap can be written in this form:
 Direct: Let D be the voxeltospace transform of the domain’s
discretization, P be the prealign matrix, Rinv the spacetovoxel transform of the reference grid (the grid the displacement field is defined on) and Cinv be the spacetovoxel transform of the codomain’s discretization. Then, for each i in the domain’s grid, the direct transform is given by
T[i] = Cinv * d[Rinv * P * D * i] + Cinv * P * D * i
and we identify U = Rinv * P * D, V = Cinv * P * D, W = Cinv
 Inverse: Let C be the voxeltospace transform of the codomain’s
discretization, Pinv be the inverse of the prealign matrix, Rinv the spacetovoxel transform of the reference grid (the grid the displacement field is defined on) and Dinv be the spacetovoxel transform of the domain’s discretization. Then, for each j in the codomain’s grid, the inverse transform is given by
Tinv[j] = Dinv * Pinv * d[Rinv * C * j] + Dinv * Pinv * C * j
and we identify U = Rinv * C, V = Dinv * Pinv * C, W = Dinv * Pinv
simplify_warp_function_3d#
 dipy.align.vector_fields.simplify_warp_function_3d(d, affine_idx_in, affine_idx_out, affine_disp, out_shape)#
Simplifies a nonlinear warping function combined with an affine transform
Modifies the given deformation field by incorporating into it an affine transformation and voxeltospace transforms associated with the discretization of its domain and codomain.
The resulting transformation may be regarded as operating on the image spaces given by the domain and codomain discretization. More precisely, the resulting transform is of the form:
T[i] = W * d[U * i] + V * i
Where U = affine_idx_in, V = affine_idx_out, W = affine_disp.
Parameters#
 darray, shape (S’, R’, C’, 3)
the nonlinear part of the transformation (displacement field)
 affine_idx_inarray, shape (4, 4)
the matrix U in eq. (1) above
 affine_idx_outarray, shape (4, 4)
the matrix V in eq. (1) above
 affine_disparray, shape (4, 4)
the matrix W in eq. (1) above
 out_shapearray, shape (3,)
the number of slices, rows and columns of the sampling grid
Returns#
 outarray, shape = out_shape
the deformation field out associated with T in eq. (1) such that: T[i] = i + out[i]
Notes#
Both the direct and inverse transforms of a DiffeomorphicMap can be written in this form:
 Direct: Let D be the voxeltospace transform of the domain’s
discretization, P be the prealign matrix, Rinv the spacetovoxel transform of the reference grid (the grid the displacement field is defined on) and Cinv be the spacetovoxel transform of the codomain’s discretization. Then, for each i in the domain’s grid, the direct transform is given by
T[i] = Cinv * d[Rinv * P * D * i] + Cinv * P * D * i
and we identify U = Rinv * P * D, V = Cinv * P * D, W = Cinv
 Inverse: Let C be the voxeltospace transform of the codomain’s
discretization, Pinv be the inverse of the prealign matrix, Rinv the spacetovoxel transform of the reference grid (the grid the displacement field is defined on) and Dinv be the spacetovoxel transform of the domain’s discretization. Then, for each j in the codomain’s grid, the inverse transform is given by
Tinv[j] = Dinv * Pinv * d[Rinv * C * j] + Dinv * Pinv * C * j
and we identify U = Rinv * C, V = Dinv * Pinv * C, W = Dinv * Pinv
sparse_gradient#
 dipy.align.vector_fields.sparse_gradient(img, img_world2grid, img_spacing, sample_points)#
Gradient of an image in physical space
Parameters#
 img2D or 3D array, shape (R, C) or (S, R, C)
the input image whose gradient will be computed
 img_world2gridarray, shape (dim+1, dim+1)
the spacetogrid transform matrix associated to img
 img_spacingarray, shape (dim,)
the spacing between voxels (voxel size along each axis) of the input image
 sample_points: array, shape (n, dim)
list of points where the derivative will be evaluated (one point per row)
Returns#
 outarray, shape (n, dim)
the gradient at each point stored at its corresponding row
transform_2d_affine#
 dipy.align.vector_fields.transform_2d_affine(image, ref_shape, affine=None)#
Transforms a 2D image by an affine transform with bilinear interp.
Deforms the input image under the given affine transformation using trilinear interpolation. The shape of the resulting image is given by ref_shape. If the affine matrix is None, it is taken as the identity.
Parameters#
 imagearray, shape (R, C)
the input image to be transformed
 ref_shapearray, shape (2,)
the shape of the resulting image
 affinearray, shape (3, 3)
the affine transform to be applied
Returns#
 outarray, shape (R’, C’)
the transformed image
Notes#
The reason it is necessary to provide the intended shape of the resulting image is because the affine transformation is defined on all R^{2} but we must sample a finite lattice. Also the resulting shape may not be necessarily equal to the input shape, unless we are interested on endomorphisms only and not general diffeomorphisms.
transform_2d_affine_nn#
 dipy.align.vector_fields.transform_2d_affine_nn(image, ref_shape, affine=None)#
Transforms a 2D image by an affine transform with NN interpolation
Deforms the input image under the given affine transformation using nearest neighbor interpolation. The shape of the resulting image is given by ref_shape. If the affine matrix is None, it is taken as the identity.
Parameters#
 imagearray, shape (R, C)
the input image to be transformed
 ref_shapearray, shape (2,)
the shape of the resulting image
 affinearray, shape (3, 3)
the affine transform to be applied
Returns#
 outarray, shape (R’, C’)
the transformed image
Notes#
The reason it is necessary to provide the intended shape of the resulting image is because the affine transformation is defined on all R^{2} but we must sample a finite lattice. Also the resulting shape may not be necessarily equal to the input shape, unless we are interested on endomorphisms only and not general diffeomorphisms.
transform_3d_affine#
 dipy.align.vector_fields.transform_3d_affine(volume, ref_shape, affine)#
Transforms a 3D volume by an affine transform with trilinear interp.
Deforms the input volume under the given affine transformation using trilinear interpolation. The shape of the resulting transformation is given by ref_shape. If the affine matrix is None, it is taken as the identity.
Parameters#
 volumearray, shape (S, R, C)
the input volume to be transformed
 ref_shapearray, shape (3,)
the shape of the resulting volume
 affinearray, shape (4, 4)
the affine transform to be applied
Returns#
 outarray, shape (S’, R’, C’)
the transformed volume
Notes#
The reason it is necessary to provide the intended shape of the resulting volume is because the affine transformation is defined on all R^{3} but we must sample a finite lattice. Also the resulting shape may not be necessarily equal to the input shape, unless we are interested on endomorphisms only and not general diffeomorphisms.
transform_3d_affine_nn#
 dipy.align.vector_fields.transform_3d_affine_nn(volume, ref_shape, affine=None)#
Transforms a 3D volume by an affine transform with NN interpolation
Deforms the input volume under the given affine transformation using nearest neighbor interpolation. The shape of the resulting volume is given by ref_shape. If the affine matrix is None, it is taken as the identity.
Parameters#
 volumearray, shape (S, R, C)
the input volume to be transformed
 ref_shapearray, shape (3,)
the shape of the resulting volume
 affinearray, shape (4, 4)
the affine transform to be applied
Returns#
 outarray, shape (S’, R’, C’)
the transformed volume
Notes#
The reason it is necessary to provide the intended shape of the resulting volume is because the affine transformation is defined on all R^{3} but we must sample a finite lattice. Also the resulting shape may not be necessarily equal to the input shape, unless we are interested on endomorphisms only and not general diffeomorphisms.
warp_2d#
 dipy.align.vector_fields.warp_2d(image, d1, affine_idx_in=None, affine_idx_out=None, affine_disp=None, out_shape=None)#
Warps a 2D image using bilinear interpolation
Deforms the input image under the given transformation. The warped image is computed using bilinear interpolation and is given by:
warped[i] = image[ C * d1[A*i] + B*i ]
where A = affine_idx_in, B = affine_idx_out, C = affine_disp and i denotes the discrete coordinates of a voxel in the sampling grid of shape = out_shape.
Parameters#
 imagearray, shape (R, C)
the input image to be transformed
 d1array, shape (R’, C’, 2)
the displacement field driving the transformation
 affine_idx_inarray, shape (3, 3)
the matrix A in eq. (1) above
 affine_idx_outarray, shape (3, 3)
the matrix B in eq. (1) above
 affine_disparray, shape (3, 3)
the matrix C in eq. (1) above
 out_shapearray, shape (2,)
the number of rows and columns of the sampling grid
Returns#
 warpedarray, shape = out_shape
the transformed image
Notes#
To illustrate the use of this function, consider a displacement field d1 with gridtospace transformation R, an image with gridtospace transformation T and let’s say we want to sample the warped image on a grid with gridtospace transformation S (sampling grid). For each voxel in the sampling grid with discrete coordinates i, the warped image is given by:
warped[i] = image[Tinv * ( d1[Rinv * S * i] + S * i ) ]
where Tinv = T^{1} and Rinv = R^{1}. By identifying A = Rinv * S, B = Tinv * S, C = Tinv we can use this function to efficiently warp the input image.
warp_2d_nn#
 dipy.align.vector_fields.warp_2d_nn(image, d1, affine_idx_in=None, affine_idx_out=None, affine_disp=None, out_shape=None)#
Warps a 2D image using nearest neighbor interpolation
Deforms the input image under the given transformation. The warped image is computed using nearestneighbor interpolation and is given by:
warped[i] = image[ C * d1[A*i] + B*i ]
where A = affine_idx_in, B = affine_idx_out, C = affine_disp and i denotes the discrete coordinates of a voxel in the sampling grid of shape = out_shape.
Parameters#
 imagearray, shape (R, C)
the input image to be transformed
 d1array, shape (R’, C’, 2)
the displacement field driving the transformation
 affine_idx_inarray, shape (3, 3)
the matrix A in eq. (1) above
 affine_idx_outarray, shape (3, 3)
the matrix B in eq. (1) above
 affine_disparray, shape (3, 3)
the matrix C in eq. (1) above
 out_shapearray, shape (2,)
the number of rows and columns of the sampling grid
Returns#
 warpedarray, shape = out_shape
the transformed image
Notes#
To illustrate the use of this function, consider a displacement field d1 with gridtospace transformation R, an image with gridtospace transformation T and let’s say we want to sample the warped image on a grid with gridtospace transformation S (sampling grid). For each voxel in the sampling grid with discrete coordinates i, the warped image is given by:
warped[i] = image[Tinv * ( d1[Rinv * S * i] + S * i ) ]
where Tinv = T^{1} and Rinv = R^{1}. By identifying A = Rinv * S, B = Tinv * S, C = Tinv we can use this function to efficiently warp the input image.
warp_3d#
 dipy.align.vector_fields.warp_3d(volume, d1, affine_idx_in=None, affine_idx_out=None, affine_disp=None, out_shape=None)#
Warps a 3D volume using trilinear interpolation
Deforms the input volume under the given transformation. The warped volume is computed using trilinear interpolation and is given by:
warped[i] = volume[ C * d1[A*i] + B*i ]
where A = affine_idx_in, B = affine_idx_out, C = affine_disp and i denotes the discrete coordinates of a voxel in the sampling grid of shape = out_shape.
Parameters#
 volumearray, shape (S, R, C)
the input volume to be transformed
 d1array, shape (S’, R’, C’, 3)
the displacement field driving the transformation
 affine_idx_inarray, shape (4, 4)
the matrix A in eq. (1) above
 affine_idx_outarray, shape (4, 4)
the matrix B in eq. (1) above
 affine_disparray, shape (4, 4)
the matrix C in eq. (1) above
 out_shapearray, shape (3,)
the number of slices, rows and columns of the sampling grid
Returns#
 warpedarray, shape = out_shape
the transformed volume
Notes#
To illustrate the use of this function, consider a displacement field d1 with gridtospace transformation R, a volume with gridtospace transformation T and let’s say we want to sample the warped volume on a grid with gridtospace transformation S (sampling grid). For each voxel in the sampling grid with discrete coordinates i, the warped volume is given by:
warped[i] = volume[Tinv * ( d1[Rinv * S * i] + S * i ) ]
where Tinv = T^{1} and Rinv = R^{1}. By identifying A = Rinv * S, B = Tinv * S, C = Tinv we can use this function to efficiently warp the input image.
warp_3d_nn#
 dipy.align.vector_fields.warp_3d_nn(volume, d1, affine_idx_in=None, affine_idx_out=None, affine_disp=None, out_shape=None)#
Warps a 3D volume using using nearestneighbor interpolation
Deforms the input volume under the given transformation. The warped volume is computed using nearestneighbor interpolation and is given by:
warped[i] = volume[ C * d1[A*i] + B*i ]
where A = affine_idx_in, B = affine_idx_out, C = affine_disp and i denotes the discrete coordinates of a voxel in the sampling grid of shape = out_shape.
Parameters#
 volumearray, shape (S, R, C)
the input volume to be transformed
 d1array, shape (S’, R’, C’, 3)
the displacement field driving the transformation
 affine_idx_inarray, shape (4, 4)
the matrix A in eq. (1) above
 affine_idx_outarray, shape (4, 4)
the matrix B in eq. (1) above
 affine_disparray, shape (4, 4)
the matrix C in eq. (1) above
 out_shapearray, shape (3,)
the number of slices, rows and columns of the sampling grid
Returns#
 warpedarray, shape = out_shape
the transformed volume
Notes#
To illustrate the use of this function, consider a displacement field d1 with gridtospace transformation R, a volume with gridtospace transformation T and let’s say we want to sample the warped volume on a grid with gridtospace transformation S (sampling grid). For each voxel in the sampling grid with discrete coordinates i, the warped volume is given by:
warped[i] = volume[Tinv * ( d1[Rinv * S * i] + S * i ) ]
where Tinv = T^{1} and Rinv = R^{1}. By identifying A = Rinv * S, B = Tinv * S, C = Tinv we can use this function to efficiently warp the input image.