ITK  5.4.0
Insight Toolkit
itkMutualInformationImageToImageMetric.h
Go to the documentation of this file.
1 /*=========================================================================
2  *
3  * Copyright NumFOCUS
4  *
5  * Licensed under the Apache License, Version 2.0 (the "License");
6  * you may not use this file except in compliance with the License.
7  * You may obtain a copy of the License at
8  *
9  * https://www.apache.org/licenses/LICENSE-2.0.txt
10  *
11  * Unless required by applicable law or agreed to in writing, software
12  * distributed under the License is distributed on an "AS IS" BASIS,
13  * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
14  * See the License for the specific language governing permissions and
15  * limitations under the License.
16  *
17  *=========================================================================*/
18 #ifndef itkMutualInformationImageToImageMetric_h
19 #define itkMutualInformationImageToImageMetric_h
20 
21 #include "itkImageToImageMetric.h"
22 #include "itkPoint.h"
23 
24 #include "itkIndex.h"
25 #include "itkKernelFunctionBase.h"
26 
27 namespace itk
28 {
93 template <typename TFixedImage, typename TMovingImage>
94 class ITK_TEMPLATE_EXPORT MutualInformationImageToImageMetric : public ImageToImageMetric<TFixedImage, TMovingImage>
95 {
96 public:
97  ITK_DISALLOW_COPY_AND_MOVE(MutualInformationImageToImageMetric);
98 
104 
106  itkNewMacro(Self);
107 
109  itkOverrideGetNameOfClassMacro(MutualInformationImageToImageMetric);
110 
112  using typename Superclass::TransformType;
113  using typename Superclass::TransformPointer;
114  using typename Superclass::TransformJacobianType;
115  using typename Superclass::InterpolatorType;
116  using typename Superclass::MeasureType;
117  using typename Superclass::DerivativeType;
118  using typename Superclass::ParametersType;
119  using typename Superclass::FixedImageType;
120  using typename Superclass::MovingImageType;
121  using typename Superclass::FixedImageConstPointer;
122  using MovingImageCosntPointer = typename Superclass::MovingImageConstPointer;
123 
130 
132 
134  static constexpr unsigned int MovingImageDimension = MovingImageType::ImageDimension;
135 
137  void
138  GetDerivative(const ParametersType & parameters, DerivativeType & derivative) const override;
139 
142  GetValue(const ParametersType & parameters) const override;
143 
145  void
146  GetValueAndDerivative(const ParametersType & parameters,
147  MeasureType & value,
148  DerivativeType & derivative) const override;
149 
154  void
155  SetNumberOfSpatialSamples(unsigned int num);
156 
158  itkGetConstReferenceMacro(NumberOfSpatialSamples, unsigned int);
159 
165  itkSetClampMacro(MovingImageStandardDeviation,
166  double,
169  itkGetConstReferenceMacro(MovingImageStandardDeviation, double);
177  itkSetClampMacro(FixedImageStandardDeviation,
178  double,
181  itkGetConstMacro(FixedImageStandardDeviation, double);
186  itkSetObjectMacro(KernelFunction, KernelFunctionType);
187  itkGetModifiableObjectMacro(KernelFunction, KernelFunctionType);
191 protected:
193  ~MutualInformationImageToImageMetric() override = default;
194  void
195  PrintSelf(std::ostream & os, Indent indent) const override;
196 
197 private:
204  {
205  public:
206  SpatialSample() { FixedImagePointValue.Fill(0.0); }
207  ~SpatialSample() = default;
208 
210  double FixedImageValue{ 0.0 };
211  double MovingImageValue{ 0.0 };
212  };
213 
215  using SpatialSampleContainer = std::vector<SpatialSample>;
216 
219  mutable SpatialSampleContainer m_SampleA{};
220 
223  mutable SpatialSampleContainer m_SampleB{};
224 
225  unsigned int m_NumberOfSpatialSamples{};
226  double m_MovingImageStandardDeviation{};
227  double m_FixedImageStandardDeviation{};
228  double m_MinProbability{};
229 
230  typename KernelFunctionType::Pointer m_KernelFunction{};
231 
243  virtual void
244  SampleFixedImageDomain(SpatialSampleContainer & samples) const;
245 
246  /*
247  * Calculate derivatives of the image intensity at the specified point with respect to the transform parameters.
248  *
249  * \todo This should really be done by the mapper.
250  *
251  * \todo This is a temporary solution until this feature is implemented
252  * in the mapper. This solution only works for any transform
253  * that support ComputeJacobianWithRespectToParameters()
254  */
255  void
256  CalculateDerivatives(const FixedImagePointType &, DerivativeType &, TransformJacobianType &) const;
257 
258  using typename Superclass::CoordinateRepresentationType;
260 
261  typename DerivativeFunctionType::Pointer m_DerivativeCalculator{};
262 };
263 } // end namespace itk
264 
265 #ifndef ITK_MANUAL_INSTANTIATION
266 # include "itkMutualInformationImageToImageMetric.hxx"
267 #endif
268 
269 #endif
itk::ImageToImageMetric
Computes similarity between regions of two images.
Definition: itkImageToImageMetric.h:54
itk::SingleValuedCostFunction::MeasureType
double MeasureType
Definition: itkSingleValuedCostFunction.h:50
itk::OptimizerParameters< TInternalComputationValueType >
itk::ImageToImageMetric::MovingImageIndexType
typename MovingImageType::IndexType MovingImageIndexType
Definition: itkImageToImageMetric.h:98
itk::ImageToImageMetric::FixedImagePointType
typename TransformType::InputPointType FixedImagePointType
Definition: itkImageToImageMetric.h:99
itk::ImageToImageMetric::FixedImageIndexValueType
typename FixedImageIndexType::IndexValueType FixedImageIndexValueType
Definition: itkImageToImageMetric.h:97
itkPoint.h
itk::CentralDifferenceImageFunction
Calculate the derivative by central differencing.
Definition: itkCentralDifferenceImageFunction.h:76
itk::KernelFunctionBase< double >
itk::ImageToImageMetric::MovingImagePointType
typename TransformType::OutputPointType MovingImagePointType
Definition: itkImageToImageMetric.h:100
itk::MutualInformationImageToImageMetric::SpatialSample::SpatialSample
SpatialSample()
Definition: itkMutualInformationImageToImageMetric.h:206
itkKernelFunctionBase.h
itk::SmartPointer< Self >
itk::Indent
Control indentation during Print() invocation.
Definition: itkIndent.h:49
itk::MutualInformationImageToImageMetric::SpatialSample::FixedImagePointValue
FixedImagePointType FixedImagePointValue
Definition: itkMutualInformationImageToImageMetric.h:209
itk::IndexValueType
long IndexValueType
Definition: itkIntTypes.h:90
itk::MutualInformationImageToImageMetric
Computes the mutual information between two images to be registered.
Definition: itkMutualInformationImageToImageMetric.h:94
itk::MutualInformationImageToImageMetric::SpatialSample
Definition: itkMutualInformationImageToImageMetric.h:203
itk::GTest::TypedefsAndConstructors::Dimension2::IndexType
ImageBaseType::IndexType IndexType
Definition: itkGTestTypedefsAndConstructors.h:50
itk::LightObject
Light weight base class for most itk classes.
Definition: itkLightObject.h:55
itkImageToImageMetric.h
itkIndex.h
itk::NumericTraits
Define additional traits for native types such as int or float.
Definition: itkNumericTraits.h:59
itk
The "itk" namespace contains all Insight Segmentation and Registration Toolkit (ITK) classes....
Definition: itkAnnulusOperator.h:24
itk::Array
Array class with size defined at construction time.
Definition: itkArray.h:47
itk::Point
A templated class holding a geometric point in n-Dimensional space.
Definition: itkPoint.h:53
itk::ImageToImageMetric::FixedImageIndexType
typename FixedImageType::IndexType FixedImageIndexType
Definition: itkImageToImageMetric.h:96
itk::MutualInformationImageToImageMetric::MovingImageCosntPointer
typename Superclass::MovingImageConstPointer MovingImageCosntPointer
Definition: itkMutualInformationImageToImageMetric.h:122
itk::MutualInformationImageToImageMetric::SpatialSampleContainer
std::vector< SpatialSample > SpatialSampleContainer
Definition: itkMutualInformationImageToImageMetric.h:215