ITK  5.0.0
Insight Segmentation and Registration Toolkit
itkJointHistogramMutualInformationImageToImageMetricv4.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 
19 #ifndef itkJointHistogramMutualInformationImageToImageMetricv4_h
20 #define itkJointHistogramMutualInformationImageToImageMetricv4_h
21 
23 #include "itkImage.h"
25 
28 
29 namespace itk
30 {
43 template<typename TFixedImage,typename TMovingImage,typename TVirtualImage = TFixedImage,
44  typename TInternalComputationValueType = double,
45  typename TMetricTraits = DefaultImageToImageMetricTraitsv4<TFixedImage,TMovingImage,TVirtualImage,TInternalComputationValueType>
46  >
48  public ImageToImageMetricv4<TFixedImage, TMovingImage, TVirtualImage, TInternalComputationValueType, TMetricTraits>
49 {
50 public:
51  ITK_DISALLOW_COPY_AND_ASSIGN(JointHistogramMutualInformationImageToImageMetricv4);
52 
55  using Superclass = ImageToImageMetricv4<TFixedImage, TMovingImage, TVirtualImage,
56  TInternalComputationValueType,TMetricTraits>;
59 
61  itkNewMacro(Self);
62 
65 
67  using CoordinateRepresentationType = typename Superclass::CoordinateRepresentationType;
68 
72  using InternalComputationValueType = TInternalComputationValueType;
73 
75  using ParametersType = typename Superclass::ParametersType;
76  using ParametersValueType = typename Superclass::ParametersValueType;
77  using NumberOfParametersType = typename Superclass::NumberOfParametersType;
78 
80  using MeasureType = typename Superclass::MeasureType;
81  using DerivativeType = typename Superclass::DerivativeType;
82  using FixedImagePointType = typename Superclass::FixedImagePointType;
83  using FixedImagePixelType = typename Superclass::FixedImagePixelType;
84  using FixedImageGradientType = typename Superclass::FixedGradientPixelType;
85  using MovingImagePointType = typename Superclass::MovingImagePointType;
86  using MovingImagePixelType = typename Superclass::MovingImagePixelType;
87  using MovingImageGradientType = typename Superclass::MovingGradientPixelType;
88 
89  using FixedTransformJacobianType = typename Superclass::FixedTransformType::JacobianType;
90  using MovingTransformJacobianType = typename Superclass::MovingTransformType::JacobianType;
91 
92  using VirtualImageType = typename Superclass::VirtualImageType;
93  using VirtualIndexType = typename Superclass::VirtualIndexType;
94  using VirtualPointType = typename Superclass::VirtualPointType;
95  using VirtualPointSetType = typename Superclass::VirtualPointSetType;
96 
97  /* Image dimension accessors */
98  static constexpr typename TVirtualImage::ImageDimensionType VirtualImageDimension = TVirtualImage::ImageDimension;
99  static constexpr typename TMovingImage::ImageDimensionType MovingImageDimension = TMovingImage::ImageDimension;
100 
102  using PDFValueType = TInternalComputationValueType;
103 
112 
114  itkGetModifiableObjectMacro(JointPDF, JointPDFType );
115 
116  // Declare the type for the derivative calculation
119  using JPDFGradientImagePointer = typename JPDFGradientImageType::Pointer;
120 
123  using MarginalGradientImagePointer = typename MarginalGradientImageType::Pointer;
124 
130 
136 
137 
139  itkSetClampMacro( NumberOfHistogramBins, SizeValueType, 5, NumericTraits< SizeValueType >::max() );
140  itkGetConstReferenceMacro(NumberOfHistogramBins, SizeValueType );
142 
144  itkSetMacro(VarianceForJointPDFSmoothing, TInternalComputationValueType);
145  itkGetMacro(VarianceForJointPDFSmoothing, TInternalComputationValueType);
147 
149  void Initialize() override;
150 
151  MeasureType GetValue() const override;
152 
153 protected:
156 
160  void InitializeForIteration() const override;
161 
163  MeasureType ComputeValue() const;
164 
167  inline void ComputeJointPDFPoint( const FixedImagePixelType fixedImageValue, const MovingImagePixelType movingImageValue, JointPDFPointType & jointPDFpoint ) const;
168 
170  friend class JointHistogramMutualInformationComputeJointPDFThreaderBase< ThreadedIndexedContainerPartitioner, Self >;
171  friend class JointHistogramMutualInformationComputeJointPDFThreader< ThreadedImageRegionPartitioner< Self::VirtualImageDimension >, Self >;
172  friend class JointHistogramMutualInformationComputeJointPDFThreader< ThreadedIndexedContainerPartitioner, Self >;
173 
175  JointHistogramMutualInformationComputeJointPDFThreader< ThreadedImageRegionPartitioner< Self::VirtualImageDimension >, Self >;
177  JointHistogramMutualInformationComputeJointPDFThreader< ThreadedIndexedContainerPartitioner, Self >;
178 
179  typename JointHistogramMutualInformationDenseComputeJointPDFThreaderType::Pointer m_JointHistogramMutualInformationDenseComputeJointPDFThreader;
180  typename JointHistogramMutualInformationSparseComputeJointPDFThreaderType::Pointer m_JointHistogramMutualInformationSparseComputeJointPDFThreader;
181 
182  friend class JointHistogramMutualInformationGetValueAndDerivativeThreader< ThreadedImageRegionPartitioner< Superclass::VirtualImageDimension >, Superclass, Self >;
183  friend class JointHistogramMutualInformationGetValueAndDerivativeThreader< ThreadedIndexedContainerPartitioner, Superclass, Self >;
184 
186  JointHistogramMutualInformationGetValueAndDerivativeThreader< ThreadedImageRegionPartitioner< Superclass::VirtualImageDimension >, Superclass, Self >;
188  JointHistogramMutualInformationGetValueAndDerivativeThreader< ThreadedIndexedContainerPartitioner, Superclass, Self >;
189 
191  void PrintSelf(std::ostream & os, Indent indent) const override;
192 
194  SizeValueType m_JointHistogramTotalCount{0};
195 
196 private:
199 
202 
205 
207  TInternalComputationValueType m_VarianceForJointPDFSmoothing;
208 
211  TInternalComputationValueType m_FixedImageTrueMin;
212  TInternalComputationValueType m_FixedImageTrueMax;
213  TInternalComputationValueType m_MovingImageTrueMin;
214  TInternalComputationValueType m_MovingImageTrueMax;
215  TInternalComputationValueType m_FixedImageBinSize;
216  TInternalComputationValueType m_MovingImageBinSize;
217 
218  TInternalComputationValueType m_JointPDFSum;
220 
221  TInternalComputationValueType m_Log2;
223 
224 };
225 
226 } // end namespace itk
227 
228 #ifndef ITK_MANUAL_INSTANTIATION
229 #include "itkJointHistogramMutualInformationImageToImageMetricv4.hxx"
230 #endif
231 
232 #endif
TPixel PixelType
Definition: itkImage.h:95
Light weight base class for most itk classes.
Define numeric traits for std::vector.
unsigned long SizeValueType
Definition: itkIntTypes.h:83
TInternalComputationValueType ParametersValueType
Class for partitioning of an ImageRegion.
typename MetricTraits::FixedImageGradientType FixedImageGradientType
typename MetricTraits::MovingImageGradientType MovingImageGradientType
typename FixedImageType::PointType FixedImagePointType
typename Superclass::VirtualImageType VirtualImageType
Computes the gradient of an image by convolution with the first derivative of a Gaussian.
typename MovingImageType::PixelType MovingImagePixelType
Provide a threaded computation of the joint PDF for JointHistogramMutualInformationImageToImageMetric...
typename Superclass::MeasureType MeasureType
TInternalComputationValueType CoordinateRepresentationType
typename Superclass::VirtualPointType VirtualPointType
typename MovingImageType::PointType MovingImagePointType
signed long IndexValueType
Definition: itkIntTypes.h:90
typename Superclass::VirtualPointSetType VirtualPointSetType
typename Superclass::SpacingType SpacingType
Definition: itkImage.h:143
Linearly interpolate an image at specified positions.
typename Superclass::VirtualIndexType VirtualIndexType
Computes the mutual information between two images to be registered using the method referenced below...
Processes points for JointHistogramMutualInformationImageToImageMetricv4 GetValueAndDerivative().
typename FixedImageType::PixelType FixedImagePixelType
typename Superclass::DerivativeType DerivativeType
TInternalComputationValueType InternalComputationValueType
typename Superclass::FixedTransformJacobianType FixedTransformJacobianType
typename Superclass::MovingTransformJacobianType MovingTransformJacobianType
Templated n-dimensional image class.
Definition: itkImage.h:75