ITK  5.0.0
Insight Segmentation and Registration Toolkit
itkMattesMutualInformationImageToImageMetricv4.h
Go to the documentation of this file.
1 /*=========================================================================
2  *
3  * Copyright Insight Software Consortium
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  * http://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 
30 namespace itk
31 {
32 
97 template <typename TFixedImage, typename TMovingImage, typename TVirtualImage = TFixedImage,
98  typename TInternalComputationValueType = double,
99  typename TMetricTraits = DefaultImageToImageMetricTraitsv4<TFixedImage,TMovingImage,TVirtualImage,TInternalComputationValueType>
100  >
102  public ImageToImageMetricv4<TFixedImage, TMovingImage, TVirtualImage, TInternalComputationValueType, TMetricTraits>
103 {
104 public:
105  ITK_DISALLOW_COPY_AND_ASSIGN(MattesMutualInformationImageToImageMetricv4);
106 
109  using Superclass = ImageToImageMetricv4<TFixedImage, TMovingImage, TVirtualImage,
110  TInternalComputationValueType,TMetricTraits>;
113 
115  itkNewMacro(Self);
116 
119 
121  using MeasureType = typename Superclass::MeasureType;
122  using DerivativeType = typename Superclass::DerivativeType;
123  using DerivativeValueType = typename DerivativeType::ValueType;
124 
125  using FixedImageType = typename Superclass::FixedImageType;
126  using FixedImagePointType = typename Superclass::FixedImagePointType;
127  using FixedImageIndexType = typename Superclass::FixedImageIndexType;
128  using FixedImagePixelType = typename Superclass::FixedImagePixelType;
129  using FixedImageGradientType = typename Superclass::FixedImageGradientType;
130 
131  using MovingImagePointType = typename Superclass::MovingImagePointType;
132  using MovingImagePixelType = typename Superclass::MovingImagePixelType;
133  using MovingImageGradientType = typename Superclass::MovingImageGradientType;
134 
135  using MovingTransformType = typename Superclass::MovingTransformType;
136  using JacobianType = typename Superclass::JacobianType;
137  using VirtualImageType = typename Superclass::VirtualImageType;
138  using VirtualIndexType = typename Superclass::VirtualIndexType;
139  using VirtualPointType = typename Superclass::VirtualPointType;
140  using VirtualPointSetType = typename Superclass::VirtualPointSetType;
141 
143  using FixedSampledPointSetPointer = typename Superclass::FixedSampledPointSetPointer;
144 
145  /* Image dimension accessors */
146  static constexpr typename TVirtualImage::ImageDimensionType VirtualImageDimension = TVirtualImage::ImageDimension;
147  static constexpr typename TFixedImage::ImageDimensionType FixedImageDimension = TFixedImage::ImageDimension;
148  static constexpr typename TMovingImage::ImageDimensionType MovingImageDimension = TMovingImage::ImageDimension;
149 
155  itkSetClampMacro( NumberOfHistogramBins, SizeValueType, 5, NumericTraits<SizeValueType>::max() );
156  itkGetConstReferenceMacro(NumberOfHistogramBins, SizeValueType);
158 
159  void Initialize() override;
160 
162  //NOTE: floating point precision is not as stable.
163  // Double precision proves faster and more robust in real-world testing.
164  using PDFValueType = TInternalComputationValueType;
165 
169 
174  const typename JointPDFType::Pointer GetJointPDF () const
175  {
176  if( this->m_ThreaderJointPDF.empty() )
177  {
178  return typename JointPDFType::Pointer(nullptr);
179  }
180  return this->m_ThreaderJointPDF[0];
181  }
183 
191  {
192  return this->m_JointPDFDerivatives;
193  }
194 
195  void FinalizeThread( const ThreadIdType threadId ) override;
196 
197 protected:
199  ~MattesMutualInformationImageToImageMetricv4() override = default;
200 
202  friend class MattesMutualInformationImageToImageMetricv4GetValueAndDerivativeThreader< ThreadedIndexedContainerPartitioner, Superclass, Self >;
204  MattesMutualInformationImageToImageMetricv4GetValueAndDerivativeThreader< ThreadedImageRegionPartitioner< Superclass::VirtualImageDimension >, Superclass, Self >;
206  MattesMutualInformationImageToImageMetricv4GetValueAndDerivativeThreader< ThreadedIndexedContainerPartitioner, Superclass, Self >;
207 
208  void PrintSelf(std::ostream& os, Indent indent) const override;
209 
218 
222 
225  virtual void GetValueCommonAfterThreadedExecution();
226 
227  OffsetValueType ComputeSingleFixedImageParzenWindowIndex( const FixedImagePixelType & value ) const;
228 
230  SizeValueType m_NumberOfHistogramBins{50};
239 
243 
246  using PRatioArrayType = std::vector<PRatioType>;
247 
249 
252  mutable std::vector<OffsetValueType> m_JointPdfIndex1DArray;
253 
255  mutable std::vector<PDFValueType> m_MovingImageMarginalPDF;
256  mutable std::vector<std::vector<PDFValueType> > m_ThreaderFixedImageMarginalPDF;
257 
259  typename std::vector<typename JointPDFType::Pointer> m_ThreaderJointPDF;
260 
261  /* \class DerivativeBufferManager
262  * A helper class to manage complexities of minimizing memory
263  * needs for mattes mutual information derivative computations
264  * per thread.
265  *
266  * Thread safety note:
267  * A seperate object is used locally per each thread. Only the members
268  * m_ParentJointPDFDerivativesLockPtr and m_ParentJointPDFDerivatives
269  * are shared between threads and access to m_ParentJointPDFDerivatives
270  * is controlled with the m_ParentJointPDFDerivativesLockPtr mutex lock.
271  * \ingroup ITKMetricsv4
272  */
274  {
276 public:
277  /* All these methods are thread safe except ReduceBuffer */
278 
279  void Initialize( size_t maxBufferLength, const size_t cachedNumberOfLocalParameters,
280  std::mutex * parentDerivativeLockPtr,
281  typename JointPDFDerivativesType::Pointer parentJointPDFDerivatives);
282 
283  void DoubleBufferSize();
284 
286  m_MemoryBlock(0)
287  {
288  }
289 
290  ~DerivativeBufferManager() = default;
291 
293  {
294  return this->m_CachedNumberOfLocalParameters;
295  }
296 
301  void CheckAndReduceIfNecessary();
302 
306  void BlockAndReduce();
307 
308  // If offset is same as previous offset, then accumulate with previous
310  {
311  m_BufferOffsetContainer[m_CurrentFillSize] = offset;
312  PDFValueType * PDFBufferForWriting = m_BufferPDFValuesContainer[m_CurrentFillSize];
313  ++m_CurrentFillSize;
314  return PDFBufferForWriting;
315  }
316 
321  void ReduceBuffer();
322 
323 private:
324  // How many AccumlatorElements used
325  size_t m_CurrentFillSize{0};
326  // Continguous chunk of memory for efficiency
327  std::vector<PDFValueType> m_MemoryBlock;
328  // The (number of lines in the buffer) * (cells per line)
330  std::vector<PDFValueType *> m_BufferPDFValuesContainer;
331  std::vector<OffsetValueType> m_BufferOffsetContainer;
334  // Pointer handle to parent version
336  // Smart pointer handle to parent version
338  };
339 
340  std::vector<DerivativeBufferManager> m_ThreaderDerivativeManager;
343 
345 
348  mutable std::vector<DerivativeType> m_LocalDerivativeByParzenBin;
349 
350 private:
352  virtual void ComputeResults() const;
353 
354 };
355 
356 } // end namespace itk
357 
358 #ifndef ITK_MANUAL_INSTANTIATION
359 #include "itkMattesMutualInformationImageToImageMetricv4.hxx"
360 #endif
361 
362 #endif
TPixel PixelType
Definition: itkImage.h:95
Light weight base class for most itk classes.
typename Superclass::MovingTransformType MovingTransformType
Define numeric traits for std::vector.
unsigned long SizeValueType
Definition: itkIntTypes.h:83
typename FixedImageType::IndexType FixedImageIndexType
Class for partitioning of an ImageRegion.
typename MetricTraits::FixedImageGradientType FixedImageGradientType
typename MetricTraits::MovingImageGradientType MovingImageGradientType
typename FixedImageType::PointType FixedImagePointType
typename Superclass::VirtualImageType VirtualImageType
typename DerivativeType::ValueType DerivativeValueType
typename MovingImageType::PixelType MovingImagePixelType
Derivative of a BSpline kernel used for density estimation and nonparameteric regression.
BSpline kernel used for density estimation and nonparameteric regression.
typename Superclass::MeasureType MeasureType
typename FixedSampledPointSetType::Pointer FixedSampledPointSetPointer
Computes the mutual information between two images to be registered using the method of Mattes et al...
typename Superclass::VirtualPointType VirtualPointType
typename MovingImageType::PointType MovingImagePointType
typename Superclass::VirtualPointSetType VirtualPointSetType
typename Superclass::RegionType RegionType
Definition: itkImage.h:139
typename Superclass::SizeType SizeType
Definition: itkImage.h:128
unsigned int ThreadIdType
Definition: itkIntTypes.h:99
typename Superclass::VirtualIndexType VirtualIndexType
typename Superclass::IndexType IndexType
Definition: itkImage.h:121
typename FixedImageType::PixelType FixedImagePixelType
typename Superclass::DerivativeType DerivativeType
typename Superclass::JacobianType JacobianType
signed long OffsetValueType
Definition: itkIntTypes.h:94
Templated n-dimensional image class.
Definition: itkImage.h:75