ITK  6.0.0
Insight Toolkit
itkMattesMutualInformationImageToImageMetric.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 itkMattesMutualInformationImageToImageMetric_h
19#define itkMattesMutualInformationImageToImageMetric_h
20
22#include "itkPoint.h"
23#include "itkIndex.h"
25#include "itkArray2D.h"
26
27#include <memory> // For unique_ptr.
28#include <mutex>
29
30
31namespace itk
32{
116template <typename TFixedImage, typename TMovingImage>
118 : public ImageToImageMetric<TFixedImage, TMovingImage>
119{
120public:
121 ITK_DISALLOW_COPY_AND_MOVE(MattesMutualInformationImageToImageMetric);
122
128
130 itkNewMacro(Self);
131
133 itkOverrideGetNameOfClassMacro(MattesMutualInformationImageToImageMetric);
134
136 using typename Superclass::TransformType;
137 using typename Superclass::TransformPointer;
138 using typename Superclass::TransformJacobianType;
139 using typename Superclass::InterpolatorType;
140 using typename Superclass::MeasureType;
141 using typename Superclass::DerivativeType;
142 using typename Superclass::ParametersType;
143 using typename Superclass::FixedImageType;
144 using typename Superclass::MovingImageType;
145 using typename Superclass::MovingImagePointType;
146 using typename Superclass::FixedImageConstPointer;
147 using typename Superclass::MovingImageConstPointer;
148 using typename Superclass::BSplineTransformWeightsType;
149 using typename Superclass::BSplineTransformIndexArrayType;
150
151 using typename Superclass::CoordinateRepresentationType;
152 using typename Superclass::FixedImageSampleContainer;
153 using typename Superclass::ImageDerivativesType;
154 using typename Superclass::WeightsValueType;
155 using typename Superclass::IndexValueType;
156
158
160 static constexpr unsigned int MovingImageDimension = MovingImageType::ImageDimension;
161
169 void
170 Initialize() override;
171
174 GetValue(const ParametersType & parameters) const override;
175
177 void
178 GetDerivative(const ParametersType & parameters, DerivativeType & derivative) const override;
179
181 void
183 MeasureType & value,
184 DerivativeType & derivative) const override;
185
192 itkSetClampMacro(NumberOfHistogramBins, SizeValueType, 5, NumericTraits<SizeValueType>::max());
193 itkGetConstReferenceMacro(NumberOfHistogramBins, SizeValueType);
220 itkSetMacro(UseExplicitPDFDerivatives, bool);
221 itkGetConstReferenceMacro(UseExplicitPDFDerivatives, bool);
222 itkBooleanMacro(UseExplicitPDFDerivatives);
226 using PDFValueType = double; // NOTE: floating point precision is not as stable. Double precision proves faster and
227 // more robust in real-world testing.
228
232
237 const typename JointPDFType::Pointer
239 {
240 if (this->m_MMIMetricPerThreadVariables == nullptr)
241 {
242 return JointPDFType::Pointer(nullptr);
243 }
244 return this->m_MMIMetricPerThreadVariables[0].JointPDF;
245 }
256 {
257 if (this->m_MMIMetricPerThreadVariables == nullptr)
258 {
259 return JointPDFDerivativesType::Pointer(nullptr);
260 }
261 return this->m_MMIMetricPerThreadVariables[0].JointPDFDerivatives;
262 }
265protected:
268 void
269 PrintSelf(std::ostream & os, Indent indent) const override;
270
271private:
280
284
286 void
288
290 void
292
294 void
296 unsigned int sampleNumber,
297 int pdfMovingIndex,
298 const ImageDerivativesType & movingImageGradientValue,
299 PDFValueType cubicBSplineDerivativeValue) const;
300
301 void
302 GetValueThreadPreProcess(ThreadIdType threadId, bool withinSampleThread) const override;
303 void
304 GetValueThreadPostProcess(ThreadIdType threadId, bool withinSampleThread) const override;
305 // NOTE: The signature in base class requires that movingImageValue is of type double
306 bool
308 SizeValueType fixedImageSample,
309 const MovingImagePointType & mappedPoint,
310 double movingImageValue) const override;
311
312 void
313 GetValueAndDerivativeThreadPreProcess(ThreadIdType threadId, bool withinSampleThread) const override;
314 void
315 GetValueAndDerivativeThreadPostProcess(ThreadIdType threadId, bool withinSampleThread) const override;
316 // NOTE: The signature in base class requires that movingImageValue is of type double
317 bool
319 SizeValueType fixedImageSample,
320 const MovingImagePointType & mappedPoint,
321 double movingImageValue,
322 const ImageDerivativesType & movingImageGradientValue) const override;
323
325 SizeValueType m_NumberOfHistogramBins{ 50 };
326 PDFValueType m_MovingImageNormalizedMin{ 0.0 };
327 PDFValueType m_FixedImageNormalizedMin{ 0.0 };
328 PDFValueType m_FixedImageTrueMin{ 0.0 };
329 PDFValueType m_FixedImageTrueMax{ 0.0 };
330 PDFValueType m_MovingImageTrueMin{ 0.0 };
331 PDFValueType m_MovingImageTrueMax{ 0.0 };
332 PDFValueType m_FixedImageBinSize{ 0.0 };
333 PDFValueType m_MovingImageBinSize{ 0.0 };
334
338
339 mutable PRatioArrayType m_PRatioArray{};
340
342 using MarginalPDFType = std::vector<PDFValueType>;
343 mutable MarginalPDFType m_MovingImageMarginalPDF{};
344
346 {
349
351
354
358
360
362 };
363
364#if !defined(ITK_WRAPPING_PARSER)
365 itkPadStruct(ITK_CACHE_LINE_ALIGNMENT, MMIMetricPerThreadStruct, PaddedMMIMetricPerThreadStruct);
366 itkAlignedTypedef(ITK_CACHE_LINE_ALIGNMENT, PaddedMMIMetricPerThreadStruct, AlignedMMIMetricPerThreadStruct);
367 // Due to a bug in older version of Visual Studio where std::vector resize
368 // uses a value instead of a const reference, this must be a pointer.
369 // See
370 // https://thetweaker.wordpress.com/2010/05/05/stdvector-of-aligned-elements/
371 // https://connect.microsoft.com/VisualStudio/feedback/details/692988
372 std::unique_ptr<AlignedMMIMetricPerThreadStruct[]> m_MMIMetricPerThreadVariables;
373#endif
374
375 bool m_UseExplicitPDFDerivatives{ true };
376 mutable bool m_ImplicitDerivativesSecondPass{ false };
377};
378} // end namespace itk
379
380#ifndef ITK_MANUAL_INSTANTIATION
381# include "itkMattesMutualInformationImageToImageMetric.hxx"
382#endif
383
384#endif
Array class with size defined at construction time.
Definition: itkArray.h:48
Derivative of a BSpline kernel used for density estimation and nonparametric regression.
BSpline kernel used for density estimation and nonparametric regression.
A templated class holding a n-Dimensional covariant vector.
An image region represents a structured region of data.
Computes similarity between regions of two images.
typename TransformType::OutputPointType MovingImagePointType
std::vector< FixedImageSamplePoint > FixedImageSampleContainer
Templated n-dimensional image class.
Definition: itkImage.h:89
TPixel PixelType
Definition: itkImage.h:108
Control indentation during Print() invocation.
Definition: itkIndent.h:50
Light weight base class for most itk classes.
Computes the mutual information between two images to be registered using the method of Mattes et al.
void GetValueAndDerivative(const ParametersType &parameters, MeasureType &value, DerivativeType &derivative) const override
void GetDerivative(const ParametersType &parameters, DerivativeType &derivative) const override
void ComputeFixedImageParzenWindowIndices(FixedImageSampleContainer &samples)
void ComputePDFDerivatives(ThreadIdType threadId, unsigned int sampleNumber, int pdfMovingIndex, const ImageDerivativesType &movingImageGradientValue, PDFValueType cubicBSplineDerivativeValue) const
bool GetValueAndDerivativeThreadProcessSample(ThreadIdType threadId, SizeValueType fixedImageSample, const MovingImagePointType &mappedPoint, double movingImageValue, const ImageDerivativesType &movingImageGradientValue) const override
const JointPDFDerivativesType::Pointer GetJointPDFDerivatives() const
void GetValueAndDerivativeThreadPreProcess(ThreadIdType threadId, bool withinSampleThread) const override
std::unique_ptr< AlignedMMIMetricPerThreadStruct[]> m_MMIMetricPerThreadVariables
itkPadStruct(ITK_CACHE_LINE_ALIGNMENT, MMIMetricPerThreadStruct, PaddedMMIMetricPerThreadStruct)
~MattesMutualInformationImageToImageMetric() override=default
void GetValueThreadPostProcess(ThreadIdType threadId, bool withinSampleThread) const override
itkAlignedTypedef(ITK_CACHE_LINE_ALIGNMENT, PaddedMMIMetricPerThreadStruct, AlignedMMIMetricPerThreadStruct)
void GetValueThreadPreProcess(ThreadIdType threadId, bool withinSampleThread) const override
MeasureType GetValue(const ParametersType &parameters) const override
void PrintSelf(std::ostream &os, Indent indent) const override
bool GetValueThreadProcessSample(ThreadIdType threadId, SizeValueType fixedImageSample, const MovingImagePointType &mappedPoint, double movingImageValue) const override
void GetValueAndDerivativeThreadPostProcess(ThreadIdType threadId, bool withinSampleThread) const override
Define additional traits for native types such as int or float.
Superclass::ParametersType ParametersType
SmartPointer< Self > Pointer
The "itk" namespace contains all Insight Segmentation and Registration Toolkit (ITK) classes....
unsigned int ThreadIdType
Definition: itkIntTypes.h:102
long IndexValueType
Definition: itkIntTypes.h:93
unsigned long SizeValueType
Definition: itkIntTypes.h:86
long OffsetValueType
Definition: itkIntTypes.h:97