ITK  6.0.0
Insight Toolkit
itkJointHistogramMutualInformationImageToImageMetricv4.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
19#ifndef itkJointHistogramMutualInformationImageToImageMetricv4_h
20#define itkJointHistogramMutualInformationImageToImageMetricv4_h
21
23#include "itkImage.h"
25
28
29namespace itk
30{
43template <typename TFixedImage,
44 typename TMovingImage,
45 typename TVirtualImage = TFixedImage,
46 typename TInternalComputationValueType = double,
47 typename TMetricTraits =
48 DefaultImageToImageMetricTraitsv4<TFixedImage, TMovingImage, TVirtualImage, TInternalComputationValueType>>
50 : public ImageToImageMetricv4<TFixedImage, TMovingImage, TVirtualImage, TInternalComputationValueType, TMetricTraits>
51{
52public:
54
57 using Superclass =
61
63 itkNewMacro(Self);
64
67
69 using typename Superclass::CoordinateRepresentationType;
70
74 using InternalComputationValueType = TInternalComputationValueType;
75
77 using typename Superclass::ParametersType;
78 using typename Superclass::ParametersValueType;
79 using typename Superclass::NumberOfParametersType;
80
82 using typename Superclass::MeasureType;
83 using typename Superclass::DerivativeType;
84 using typename Superclass::FixedImagePointType;
85 using typename Superclass::FixedImagePixelType;
86 using FixedImageGradientType = typename Superclass::FixedGradientPixelType;
87 using typename Superclass::MovingImagePointType;
88 using typename Superclass::MovingImagePixelType;
89 using MovingImageGradientType = typename Superclass::MovingGradientPixelType;
90
91 using FixedTransformJacobianType = typename Superclass::FixedTransformType::JacobianType;
92 using MovingTransformJacobianType = typename Superclass::MovingTransformType::JacobianType;
93
94 using VirtualImageType = typename Superclass::VirtualImageType;
95 using typename Superclass::VirtualIndexType;
96 using typename Superclass::VirtualPointType;
97 using typename Superclass::VirtualPointSetType;
98
99 /* Image dimension accessors */
100 static constexpr typename TVirtualImage::ImageDimensionType VirtualImageDimension = TVirtualImage::ImageDimension;
101 static constexpr typename TMovingImage::ImageDimensionType MovingImageDimension = TMovingImage::ImageDimension;
102
104 using PDFValueType = TInternalComputationValueType;
105
114
116 itkGetModifiableObjectMacro(JointPDF, JointPDFType);
117
118 // Declare the type for the derivative calculation
122
126
132
138
139
141 itkSetClampMacro(NumberOfHistogramBins, SizeValueType, 5, NumericTraits<SizeValueType>::max());
142 itkGetConstReferenceMacro(NumberOfHistogramBins, SizeValueType);
146 itkSetMacro(VarianceForJointPDFSmoothing, TInternalComputationValueType);
147 itkGetMacro(VarianceForJointPDFSmoothing, TInternalComputationValueType);
151 void
152 Initialize() override;
153
155 GetValue() const override;
156
157protected:
160
164 void
165 InitializeForIteration() const override;
166
170
173 inline void
175 const MovingImagePixelType movingImageValue,
176 JointPDFPointType & jointPDFpoint) const;
177
179 ThreadedImageRegionPartitioner<Self::VirtualImageDimension>,
180 Self>;
181 friend class JointHistogramMutualInformationComputeJointPDFThreaderBase<ThreadedIndexedContainerPartitioner, Self>;
183 ThreadedImageRegionPartitioner<Self::VirtualImageDimension>,
184 Self>;
185 friend class JointHistogramMutualInformationComputeJointPDFThreader<ThreadedIndexedContainerPartitioner, Self>;
186
188 JointHistogramMutualInformationComputeJointPDFThreader<ThreadedImageRegionPartitioner<Self::VirtualImageDimension>,
189 Self>;
191 JointHistogramMutualInformationComputeJointPDFThreader<ThreadedIndexedContainerPartitioner, Self>;
192
194 m_JointHistogramMutualInformationDenseComputeJointPDFThreader;
196 m_JointHistogramMutualInformationSparseComputeJointPDFThreader;
197
199 ThreadedImageRegionPartitioner<Superclass::VirtualImageDimension>,
200 Superclass,
201 Self>;
202 friend class JointHistogramMutualInformationGetValueAndDerivativeThreader<ThreadedIndexedContainerPartitioner,
203 Superclass,
204 Self>;
205
208 ThreadedImageRegionPartitioner<Superclass::VirtualImageDimension>,
209 Superclass,
210 Self>;
212 JointHistogramMutualInformationGetValueAndDerivativeThreader<ThreadedIndexedContainerPartitioner, Superclass, Self>;
213
215 void
216 PrintSelf(std::ostream & os, Indent indent) const override;
217
219 SizeValueType m_JointHistogramTotalCount{ 0 };
220
221private:
223 typename MarginalPDFType::Pointer m_FixedImageMarginalPDF{};
224
226 typename MarginalPDFType::Pointer m_MovingImageMarginalPDF{};
227
229 mutable typename JointPDFType::Pointer m_JointPDF{};
230
232 TInternalComputationValueType m_VarianceForJointPDFSmoothing{};
233
235 SizeValueType m_NumberOfHistogramBins{};
236 TInternalComputationValueType m_FixedImageTrueMin{};
237 TInternalComputationValueType m_FixedImageTrueMax{};
238 TInternalComputationValueType m_MovingImageTrueMin{};
239 TInternalComputationValueType m_MovingImageTrueMax{};
240 TInternalComputationValueType m_FixedImageBinSize{};
241 TInternalComputationValueType m_MovingImageBinSize{};
242
243 TInternalComputationValueType m_JointPDFSum{};
244 JointPDFSpacingType m_JointPDFSpacing{};
245
246 TInternalComputationValueType m_Log2{};
248};
249
250} // end namespace itk
251
252#ifndef ITK_MANUAL_INSTANTIATION
253# include "itkJointHistogramMutualInformationImageToImageMetricv4.hxx"
254#endif
255
256#endif
Computes the gradient of an image by convolution with the first derivative of a Gaussian.
Templated n-dimensional image class.
Definition: itkImage.h:89
TPixel PixelType
Definition: itkImage.h:108
Provide a threaded computation of the joint PDF for JointHistogramMutualInformationImageToImageMetric...
Processes points for JointHistogramMutualInformationImageToImageMetricv4 GetValueAndDerivative().
Computes the mutual information between two images to be registered using the method referenced below...
void ComputeJointPDFPoint(const FixedImagePixelType fixedImageValue, const MovingImagePixelType movingImageValue, JointPDFPointType &jointPDFpoint) const
Light weight base class for most itk classes.
Linearly interpolate an image at specified positions.
Define additional traits for native types such as int or float.
Class for partitioning of an ImageRegion.
SmartPointer< Self > Pointer
The "itk" namespace contains all Insight Segmentation and Registration Toolkit (ITK) classes....
long IndexValueType
Definition: itkIntTypes.h:93
unsigned long SizeValueType
Definition: itkIntTypes.h:86