#include <itkMutualInformationImageToImageMetric.h>
Computes the mutual information between two images to be registered.
MutualInformationImageToImageMetric computes the mutual information between a fixed and moving image to be registered.
This class is templated over the FixedImage type and the MovingImage type.
The fixed and moving images are set via methods SetFixedImage() and SetMovingImage(). This metric makes use of user specified Transform and Interpolator. The Transform is used to map points from the fixed image to the moving image domain. The Interpolator is used to evaluate the image intensity at user specified geometric points in the moving image. The Transform and Interpolator are set via methods SetTransform() and SetInterpolator().
The method GetValue() computes of the mutual information while method GetValueAndDerivative() computes both the mutual information and its derivatives with respect to the transform parameters.
The calculations are based on the method of Viola and Wells where the probability density distributions are estimated using Parzen windows.
By default a Gaussian kernel is used in the density estimation. Other option include Cauchy and spline-based. A user can specify the kernel passing in a pointer a KernelFunctionBase using the SetKernelFunction() method.
Mutual information is estimated using two sample sets: one to calculate the singular and joint pdf's and one to calculate the entropy integral. By default 50 samples points are used in each set. Other values can be set via the SetNumberOfSpatialSamples() method.
Quality of the density estimate depends on the choice of the kernel's standard deviation. Optimal choice will depend on the images. It is can be shown that around the optimal variance, the mutual information estimate is relatively insensitive to small changes of the standard deviation. In our experiments, we have found that a standard deviation of 0.4 works well for images normalized to have a mean of zero and standard deviation of 1.0. The variance can be set via methods SetFixedImageStandardDeviation() and SetMovingImageStandardDeviation().
Implementation of this class is based on [134].
Definition at line 91 of file itkMutualInformationImageToImageMetric.h.
Classes | |
class | SpatialSample |
Static Public Member Functions | |
static Pointer | New () |
![]() | |
static bool | GetGlobalWarningDisplay () |
static void | GlobalWarningDisplayOff () |
static void | GlobalWarningDisplayOn () |
static Pointer | New () |
static void | SetGlobalWarningDisplay (bool val) |
![]() | |
static void | BreakOnError () |
static Pointer | New () |
Static Public Attributes | |
static constexpr unsigned int | MovingImageDimension = MovingImageType::ImageDimension |
![]() | |
static constexpr unsigned int | FixedImageDimension = TFixedImage::ImageDimension |
static constexpr unsigned int | MovingImageDimension = TMovingImage::ImageDimension |
Protected Member Functions | |
MutualInformationImageToImageMetric () | |
void | PrintSelf (std::ostream &os, Indent indent) const override |
~MutualInformationImageToImageMetric () override=default | |
![]() | |
virtual void | ComputeImageDerivatives (const MovingImagePointType &mappedPoint, ImageDerivativesType &gradient, ThreadIdType threadId) const |
void | GetValueAndDerivativeMultiThreadedInitiate () const |
void | GetValueAndDerivativeMultiThreadedPostProcessInitiate () const |
virtual void | GetValueAndDerivativeThread (ThreadIdType threadId) const |
virtual void | GetValueAndDerivativeThreadPostProcess (ThreadIdType threadId, bool withinSampleThread) const |
virtual void | GetValueAndDerivativeThreadPreProcess (ThreadIdType threadId, bool withinSampleThread) const |
virtual bool | GetValueAndDerivativeThreadProcessSample (ThreadIdType threadId, SizeValueType fixedImageSample, const MovingImagePointType &mappedPoint, double movingImageValue, const ImageDerivativesType &movingImageGradientValue) const |
ImageToImageMetric () | |
virtual void | PreComputeTransformValues () |
virtual void | SampleFixedImageIndexes (FixedImageSampleContainer &samples) const |
virtual void | SampleFixedImageRegion (FixedImageSampleContainer &samples) const |
virtual void | SampleFullFixedImageRegion (FixedImageSampleContainer &samples) const |
virtual void | SynchronizeTransforms () const |
virtual void | TransformPoint (unsigned int sampleNumber, MovingImagePointType &mappedPoint, bool &sampleOk, double &movingImageValue, ThreadIdType threadId) const |
virtual void | TransformPointWithDerivatives (unsigned int sampleNumber, MovingImagePointType &mappedPoint, bool &sampleOk, double &movingImageValue, ImageDerivativesType &movingImageGradient, ThreadIdType threadId) const |
~ImageToImageMetric () override=default | |
void | GetValueMultiThreadedInitiate () const |
void | GetValueMultiThreadedPostProcessInitiate () const |
virtual void | GetValueThread (ThreadIdType threadId) const |
virtual void | GetValueThreadPreProcess (ThreadIdType threadId, bool withinSampleThread) const |
virtual bool | GetValueThreadProcessSample (ThreadIdType threadId, SizeValueType fixedImageSample, const MovingImagePointType &mappedPoint, double movingImageValue) const |
virtual void | GetValueThreadPostProcess (ThreadIdType threadId, bool withinSampleThread) const |
![]() | |
SingleValuedCostFunction ()=default | |
~SingleValuedCostFunction () override | |
![]() | |
CostFunctionTemplate ()=default | |
CostFunctionTemplate ()=default | |
~CostFunctionTemplate () override=default | |
~CostFunctionTemplate () override=default | |
![]() | |
Object () | |
bool | PrintObservers (std::ostream &os, Indent indent) const |
void | PrintSelf (std::ostream &os, Indent indent) const override |
virtual void | SetTimeStamp (const TimeStamp &timeStamp) |
~Object () override | |
![]() | |
virtual LightObject::Pointer | InternalClone () const |
LightObject () | |
virtual void | PrintHeader (std::ostream &os, Indent indent) const |
virtual void | PrintTrailer (std::ostream &os, Indent indent) const |
virtual | ~LightObject () |
Private Types | |
using | CoordinateRepresentationType |
using | DerivativeFunctionType = CentralDifferenceImageFunction<MovingImageType, CoordinateRepresentationType> |
using | SpatialSampleContainer = std::vector<SpatialSample> |
Private Member Functions | |
void | CalculateDerivatives (const FixedImagePointType &, DerivativeType &, TransformJacobianType &) const |
virtual void | SampleFixedImageDomain (SpatialSampleContainer &samples) const |
Private Attributes | |
DerivativeFunctionType::Pointer | m_DerivativeCalculator {} |
double | m_FixedImageStandardDeviation {} |
KernelFunctionType::Pointer | m_KernelFunction {} |
double | m_MinProbability {} |
double | m_MovingImageStandardDeviation {} |
unsigned int | m_NumberOfSpatialSamples {} |
SpatialSampleContainer | m_SampleA {} |
SpatialSampleContainer | m_SampleB {} |
using itk::MutualInformationImageToImageMetric< TFixedImage, TMovingImage >::ConstPointer = SmartPointer<const Self> |
Definition at line 100 of file itkMutualInformationImageToImageMetric.h.
|
private |
Type used for representing point components
Definition at line 66 of file itkImageToImageMetric.h.
|
private |
Definition at line 256 of file itkMutualInformationImageToImageMetric.h.
using itk::SingleValuedCostFunction::DerivativeType |
DerivativeType type alias. It defines a type used to return the cost function derivative.
Definition at line 133 of file itkSingleValuedCostFunction.h.
using itk::ImageToImageMetric< TFixedImage, TMovingImage >::FixedImageConstPointer |
Definition at line 79 of file itkImageToImageMetric.h.
using itk::MutualInformationImageToImageMetric< TFixedImage, TMovingImage >::FixedImageIndexType = typename FixedImageType::IndexType |
Index and Point type alias support
Definition at line 122 of file itkMutualInformationImageToImageMetric.h.
using itk::MutualInformationImageToImageMetric< TFixedImage, TMovingImage >::FixedImageIndexValueType = typename FixedImageIndexType::IndexValueType |
Definition at line 123 of file itkMutualInformationImageToImageMetric.h.
using itk::MutualInformationImageToImageMetric< TFixedImage, TMovingImage >::FixedImagePointType = typename TransformType::InputPointType |
Definition at line 125 of file itkMutualInformationImageToImageMetric.h.
using itk::ImageToImageMetric< TFixedImage, TMovingImage >::FixedImageType |
Type of the fixed Image.
Definition at line 77 of file itkImageToImageMetric.h.
using itk::ImageToImageMetric< TFixedImage, TMovingImage >::InterpolatorType |
Type of the Interpolator Base class
Definition at line 105 of file itkImageToImageMetric.h.
using itk::MutualInformationImageToImageMetric< TFixedImage, TMovingImage >::KernelFunctionType = KernelFunctionBase<double> |
Definition at line 128 of file itkMutualInformationImageToImageMetric.h.
using itk::SingleValuedCostFunction::MeasureType |
MeasureType type alias. It defines a type used to return the cost function value.
Definition at line 130 of file itkSingleValuedCostFunction.h.
using itk::MutualInformationImageToImageMetric< TFixedImage, TMovingImage >::MovingImageCosntPointer = typename Superclass::MovingImageConstPointer |
Definition at line 119 of file itkMutualInformationImageToImageMetric.h.
using itk::MutualInformationImageToImageMetric< TFixedImage, TMovingImage >::MovingImageIndexType = typename MovingImageType::IndexType |
Definition at line 124 of file itkMutualInformationImageToImageMetric.h.
using itk::MutualInformationImageToImageMetric< TFixedImage, TMovingImage >::MovingImagePointType = typename TransformType::OutputPointType |
Definition at line 126 of file itkMutualInformationImageToImageMetric.h.
using itk::ImageToImageMetric< TFixedImage, TMovingImage >::MovingImageType |
Type of the moving Image.
Definition at line 72 of file itkImageToImageMetric.h.
using itk::SingleValuedCostFunction::ParametersType |
ParametersType type alias. It defines a position in the optimization search space.
Definition at line 136 of file itkSingleValuedCostFunction.h.
using itk::MutualInformationImageToImageMetric< TFixedImage, TMovingImage >::Pointer = SmartPointer<Self> |
Definition at line 99 of file itkMutualInformationImageToImageMetric.h.
using itk::MutualInformationImageToImageMetric< TFixedImage, TMovingImage >::Self = MutualInformationImageToImageMetric |
Standard class type aliases.
Definition at line 97 of file itkMutualInformationImageToImageMetric.h.
|
private |
SpatialSampleContainer type alias support
Definition at line 212 of file itkMutualInformationImageToImageMetric.h.
using itk::MutualInformationImageToImageMetric< TFixedImage, TMovingImage >::Superclass = ImageToImageMetric<TFixedImage, TMovingImage> |
Definition at line 98 of file itkMutualInformationImageToImageMetric.h.
using itk::ImageToImageMetric< TFixedImage, TMovingImage >::TransformJacobianType |
Definition at line 93 of file itkImageToImageMetric.h.
using itk::ImageToImageMetric< TFixedImage, TMovingImage >::TransformPointer |
Definition at line 89 of file itkImageToImageMetric.h.
using itk::ImageToImageMetric< TFixedImage, TMovingImage >::TransformType |
Type of the Transform Base class
Definition at line 87 of file itkImageToImageMetric.h.
|
protected |
Referenced by GetNameOfClass().
|
overrideprotecteddefault |
|
private |
|
virtual |
Create an object from an instance, potentially deferring to a factory. This method allows you to create an instance of an object that is exactly the same type as the referring object. This is useful in cases where an object has been cast back to a base class.
Reimplemented from itk::LightObject.
|
overridevirtual |
Get the derivatives of the match measure.
Implements itk::SingleValuedCostFunction.
|
virtual |
Set/Get the fixed image intensity standard deviation. This defines the kernel bandwidth used in the joint probability distribution calculation. Default value is 0.4 which works well for image intensities normalized to a mean of 0 and standard deviation of 1.0. Value is clamped to be always greater than zero.
|
virtual |
Set/Get the kernel function. This is used to calculate the joint probability distribution. Default is the GaussianKernelFunction.
|
virtual |
Set/Get the kernel function. This is used to calculate the joint probability distribution. Default is the GaussianKernelFunction.
|
virtual |
Set/Get the moving image intensity standard deviation. This defines the kernel bandwidth used in the joint probability distribution calculation. Default value is 0.4 which works well for image intensities normalized to a mean of 0 and standard deviation of 1.0. Value is clamped to be always greater than zero.
References itk::NumericTraits< T >::NonpositiveMin().
|
overridevirtual |
Reimplemented from itk::ImageToImageMetric< TFixedImage, TMovingImage >.
References MutualInformationImageToImageMetric().
|
virtual |
Get the number of spatial samples.
References itk::NumericTraits< T >::NonpositiveMin().
|
overridevirtual |
Get the value.
Implements itk::SingleValuedCostFunction.
|
overridevirtual |
Get the value and derivatives for single valued optimizers.
Reimplemented from itk::SingleValuedCostFunction.
|
static |
Method for creation through the object factory.
|
overrideprotectedvirtual |
Methods invoked by Print() to print information about the object including superclasses. Typically not called by the user (use Print() instead) but used in the hierarchical print process to combine the output of several classes.
Reimplemented from itk::ImageToImageMetric< TFixedImage, TMovingImage >.
|
privatevirtual |
Uniformly select samples from the fixed image buffer.
Each sample consists of:
|
virtual |
Set/Get the fixed image intensity standard deviation. This defines the kernel bandwidth used in the joint probability distribution calculation. Default value is 0.4 which works well for image intensities normalized to a mean of 0 and standard deviation of 1.0. Value is clamped to be always greater than zero.
References itk::NumericTraits< T >::max().
|
virtual |
Set/Get the kernel function. This is used to calculate the joint probability distribution. Default is the GaussianKernelFunction.
|
virtual |
Set/Get the moving image intensity standard deviation. This defines the kernel bandwidth used in the joint probability distribution calculation. Default value is 0.4 which works well for image intensities normalized to a mean of 0 and standard deviation of 1.0. Value is clamped to be always greater than zero.
References itk::NumericTraits< T >::max().
void itk::MutualInformationImageToImageMetric< TFixedImage, TMovingImage >::SetNumberOfSpatialSamples | ( | unsigned int | num | ) |
Set the number of spatial samples. This is the number of image samples used to calculate the joint probability distribution. The number of spatial samples is clamped to be a minimum of 1. Default value is 50.
|
private |
Definition at line 258 of file itkMutualInformationImageToImageMetric.h.
|
private |
Definition at line 224 of file itkMutualInformationImageToImageMetric.h.
|
private |
Definition at line 227 of file itkMutualInformationImageToImageMetric.h.
|
private |
Definition at line 225 of file itkMutualInformationImageToImageMetric.h.
|
private |
Definition at line 223 of file itkMutualInformationImageToImageMetric.h.
|
private |
Definition at line 222 of file itkMutualInformationImageToImageMetric.h.
|
mutableprivate |
Container to store sample set A - used to approximate the probability density function (pdf).
Definition at line 216 of file itkMutualInformationImageToImageMetric.h.
|
mutableprivate |
Container to store sample set B - used to approximate the mutual information value.
Definition at line 220 of file itkMutualInformationImageToImageMetric.h.
|
staticconstexpr |
Enum of the moving image dimension.
Definition at line 131 of file itkMutualInformationImageToImageMetric.h.