ITK  6.0.0
Insight Toolkit
itkMattesMutualInformationImageToImageMetricv4.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 itkMattesMutualInformationImageToImageMetricv4_h
19#define itkMattesMutualInformationImageToImageMetricv4_h
20
23#include "itkPoint.h"
24#include "itkIndex.h"
26#include "itkArray2D.h"
28#include <mutex>
29
30namespace itk
31{
32
97template <typename TFixedImage,
98 typename TMovingImage,
99 typename TVirtualImage = TFixedImage,
100 typename TInternalComputationValueType = double,
101 typename TMetricTraits =
102 DefaultImageToImageMetricTraitsv4<TFixedImage, TMovingImage, TVirtualImage, TInternalComputationValueType>>
104 : public ImageToImageMetricv4<TFixedImage, TMovingImage, TVirtualImage, TInternalComputationValueType, TMetricTraits>
105{
106public:
107 ITK_DISALLOW_COPY_AND_MOVE(MattesMutualInformationImageToImageMetricv4);
108
115
117 itkNewMacro(Self);
118
120 itkOverrideGetNameOfClassMacro(MattesMutualInformationImageToImageMetricv4);
121
123 using typename Superclass::MeasureType;
124 using typename Superclass::DerivativeType;
126
127 using typename Superclass::FixedImageType;
128 using typename Superclass::FixedImagePointType;
129 using typename Superclass::FixedImageIndexType;
130 using typename Superclass::FixedImagePixelType;
131 using typename Superclass::FixedImageGradientType;
132
133 using typename Superclass::MovingImagePointType;
134 using typename Superclass::MovingImagePixelType;
135 using typename Superclass::MovingImageGradientType;
136
137 using typename Superclass::MovingTransformType;
138 using typename Superclass::JacobianType;
139 using VirtualImageType = typename Superclass::VirtualImageType;
140 using typename Superclass::VirtualIndexType;
141 using typename Superclass::VirtualPointType;
142 using typename Superclass::VirtualPointSetType;
143
145 using typename Superclass::FixedSampledPointSetPointer;
146
147 /* Image dimension accessors */
148 static constexpr typename TVirtualImage::ImageDimensionType VirtualImageDimension = TVirtualImage::ImageDimension;
149 static constexpr typename TFixedImage::ImageDimensionType FixedImageDimension = TFixedImage::ImageDimension;
150 static constexpr typename TMovingImage::ImageDimensionType MovingImageDimension = TMovingImage::ImageDimension;
151
157 itkSetClampMacro(NumberOfHistogramBins, SizeValueType, 5, NumericTraits<SizeValueType>::max());
158 itkGetConstReferenceMacro(NumberOfHistogramBins, SizeValueType);
161 void
162 Initialize() override;
163
165 // NOTE: floating point precision is not as stable.
166 // Double precision proves faster and more robust in real-world testing.
167 using PDFValueType = TInternalComputationValueType;
168
172
177 const typename JointPDFType::Pointer
179 {
180 if (this->m_ThreaderJointPDF.empty())
181 {
182 return typename JointPDFType::Pointer(nullptr);
183 }
184 return this->m_ThreaderJointPDF[0];
185 }
196 {
197 return this->m_JointPDFDerivatives;
198 }
199
200 void
201 FinalizeThread(const ThreadIdType threadId) override;
202
203protected:
206
208 ThreadedImageRegionPartitioner<Superclass::VirtualImageDimension>,
210 Self>;
212 ThreadedIndexedContainerPartitioner,
213 Superclass,
214 Self>;
217 ThreadedImageRegionPartitioner<Superclass::VirtualImageDimension>,
218 Superclass,
219 Self>;
222 Superclass,
223 Self>;
224
225 void
226 PrintSelf(std::ostream & os, Indent indent) const override;
227
236
240
243 virtual void
244 GetValueCommonAfterThreadedExecution();
245
247 ComputeSingleFixedImageParzenWindowIndex(const FixedImagePixelType & value) const;
248
250 SizeValueType m_NumberOfHistogramBins{ 50 };
251 PDFValueType m_MovingImageNormalizedMin{};
252 PDFValueType m_FixedImageNormalizedMin{};
253 PDFValueType m_FixedImageTrueMin{};
254 PDFValueType m_FixedImageTrueMax{};
255 PDFValueType m_MovingImageTrueMin{};
256 PDFValueType m_MovingImageTrueMax{};
257 PDFValueType m_FixedImageBinSize{};
258 PDFValueType m_MovingImageBinSize{};
259
262 using PRatioArrayType = std::vector<PRatioType>;
263
264 mutable PRatioArrayType m_PRatioArray{};
265
268 mutable std::vector<OffsetValueType> m_JointPdfIndex1DArray{};
269
271 mutable std::vector<PDFValueType> m_MovingImageMarginalPDF{};
272 mutable std::vector<std::vector<PDFValueType>> m_ThreaderFixedImageMarginalPDF{};
273
275 typename std::vector<typename JointPDFType::Pointer> m_ThreaderJointPDF{};
276
277 /* \class DerivativeBufferManager
278 * A helper class to manage complexities of minimizing memory
279 * needs for mattes mutual information derivative computations
280 * per thread.
281 *
282 * Thread safety note:
283 * A separate object is used locally per each thread. Only the members
284 * m_ParentJointPDFDerivativesMutexPtr and m_ParentJointPDFDerivatives
285 * are shared between threads and access to m_ParentJointPDFDerivatives
286 * is controlled with the m_ParentJointPDFDerivativesMutexPtr mutex lock.
287 * \ingroup ITKMetricsv4
288 */
290 {
292
293 public:
294 /* All these methods are thread safe except ReduceBuffer */
295
296 void
297 Initialize(size_t maxBufferLength,
298 const size_t cachedNumberOfLocalParameters,
299 std::mutex * parentDerivativeMutexPtr,
300 typename JointPDFDerivativesType::Pointer parentJointPDFDerivatives);
301
302 void
304
306 : m_MemoryBlock(0)
307 {}
308
310
311 size_t
313 {
314 return this->m_CachedNumberOfLocalParameters;
315 }
316
321 void
323
327 void
329
330 // If offset is same as previous offset, then accumulate with previous
333 {
334 m_BufferOffsetContainer[m_CurrentFillSize] = offset;
335 PDFValueType * PDFBufferForWriting = m_BufferPDFValuesContainer[m_CurrentFillSize];
336 ++m_CurrentFillSize;
337 return PDFBufferForWriting;
338 }
339
344 void
346
347 private:
348 // How many AccumulatorElements used
349 size_t m_CurrentFillSize{ 0 };
350 // Continguous chunk of memory for efficiency
351 std::vector<PDFValueType> m_MemoryBlock;
352 // The (number of lines in the buffer) * (cells per line)
354 std::vector<PDFValueType *> m_BufferPDFValuesContainer;
355 std::vector<OffsetValueType> m_BufferOffsetContainer;
358 // Pointer handle to parent version
360 // Smart pointer handle to parent version
362 };
363
364 std::vector<DerivativeBufferManager> m_ThreaderDerivativeManager{};
365 std::mutex m_JointPDFDerivativesLock{};
366 typename JointPDFDerivativesType::Pointer m_JointPDFDerivatives{};
367
368 PDFValueType m_JointPDFSum{};
369
372 mutable std::vector<DerivativeType> m_LocalDerivativeByParzenBin{};
373
374private:
376 virtual void
378};
379
380} // end namespace itk
381
382#ifndef ITK_MANUAL_INSTANTIATION
383# include "itkMattesMutualInformationImageToImageMetricv4.hxx"
384#endif
385
386#endif
Derivative of a BSpline kernel used for density estimation and nonparametric regression.
BSpline kernel used for density estimation and nonparametric regression.
An image region represents a structured region of data.
Templated n-dimensional image class.
Definition: itkImage.h:89
TPixel PixelType
Definition: itkImage.h:108
Light weight base class for most itk classes.
void Initialize(vcl_size_t maxBufferLength, const vcl_size_t cachedNumberOfLocalParameters, std::mutex *parentDerivativeMutexPtr, typename JointPDFDerivativesType::Pointer parentJointPDFDerivatives)
Computes the mutual information between two images to be registered using the method of Mattes et al.
void FinalizeThread(const ThreadIdType threadId) override
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....
unsigned int ThreadIdType
Definition: itkIntTypes.h:102
unsigned long SizeValueType
Definition: itkIntTypes.h:86
long OffsetValueType
Definition: itkIntTypes.h:97