ITK  4.2.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 {
44 template<class TFixedImage,class TMovingImage,class TVirtualImage = TFixedImage>
46  public ImageToImageMetricv4<TFixedImage, TMovingImage, TVirtualImage>
47 {
48 public:
49 
55 
57  itkNewMacro(Self);
58 
62 
64  typedef typename Superclass::CoordinateRepresentationType
66 
70 
72  typedef typename Superclass::ParametersType ParametersType;
73  typedef typename Superclass::ParametersValueType ParametersValueType;
74  typedef typename Superclass::NumberOfParametersType NumberOfParametersType;
75 
77  typedef typename Superclass::MeasureType MeasureType;
78  typedef typename Superclass::DerivativeType DerivativeType;
79  typedef typename Superclass::FixedImagePointType FixedImagePointType;
80  typedef typename Superclass::FixedImagePixelType FixedImagePixelType;
81  typedef typename Superclass::FixedGradientPixelType FixedImageGradientType;
82  typedef typename Superclass::MovingImagePointType MovingImagePointType;
83  typedef typename Superclass::MovingImagePixelType MovingImagePixelType;
84  typedef typename Superclass::MovingGradientPixelType MovingImageGradientType;
85 
86  typedef typename Superclass::FixedTransformType::JacobianType FixedTransformJacobianType;
87  typedef typename Superclass::MovingTransformType::JacobianType MovingTransformJacobianType;
88 
89  typedef typename Superclass::VirtualImageType VirtualImageType;
90  typedef typename Superclass::VirtualIndexType VirtualIndexType;
91  typedef typename Superclass::VirtualPointType VirtualPointType;
92  typedef typename Superclass::VirtualPointSetType VirtualPointSetType;
93 
94  /* Image dimension accessors */
95  itkStaticConstMacro(VirtualImageDimension, ImageDimensionType,
97  itkStaticConstMacro(MovingImageDimension, ImageDimensionType,
99 
102 
111 
113  itkGetConstObjectMacro( JointPDF, JointPDFType );
114 
115  // Declare the type for the derivative calculation
118  typedef typename JPDFGradientImageType::Pointer JPDFGradientImagePointer;
119 
122  typedef typename MarginalGradientImageType::Pointer MarginalGradientImagePointer;
123 
129 
135 
136 
138  itkSetClampMacro( NumberOfHistogramBins, SizeValueType, 5, NumericTraits< SizeValueType >::max() );
139  itkGetConstReferenceMacro(NumberOfHistogramBins, SizeValueType );
141 
143  itkSetMacro(VarianceForJointPDFSmoothing, InternalComputationValueType);
144  itkGetMacro(VarianceForJointPDFSmoothing, InternalComputationValueType);
146 
148  virtual void Initialize() throw (itk::ExceptionObject);
149 
150  virtual MeasureType GetValue() const;
151 
152 protected:
154  virtual ~JointHistogramMutualInformationImageToImageMetricv4();
155 
159  virtual void InitializeForIteration() const;
160 
162  MeasureType ComputeValue() const;
163 
166  inline void ComputeJointPDFPoint( const FixedImagePixelType fixedImageValue, const MovingImagePixelType movingImageValue, JointPDFPointType & jointPDFpoint ) const;
167 
171  friend class JointHistogramMutualInformationComputeJointPDFThreader< ThreadedIndexedContainerPartitioner, Self >;
172 
173  typedef JointHistogramMutualInformationComputeJointPDFThreader< ThreadedImageRegionPartitioner< Self::VirtualImageDimension >, Self >
175  typedef JointHistogramMutualInformationComputeJointPDFThreader< ThreadedIndexedContainerPartitioner, Self >
177 
178  typename JointHistogramMutualInformationDenseComputeJointPDFThreaderType::Pointer m_JointHistogramMutualInformationDenseComputeJointPDFThreader;
179  typename JointHistogramMutualInformationSparseComputeJointPDFThreaderType::Pointer m_JointHistogramMutualInformationSparseComputeJointPDFThreader;
180 
182  friend class JointHistogramMutualInformationGetValueAndDerivativeThreader< ThreadedIndexedContainerPartitioner, Superclass, Self >;
183 
184  typedef JointHistogramMutualInformationGetValueAndDerivativeThreader< ThreadedImageRegionPartitioner< Superclass::VirtualImageDimension >, Superclass, Self >
186  typedef JointHistogramMutualInformationGetValueAndDerivativeThreader< ThreadedIndexedContainerPartitioner, Superclass, Self >
188 
190  void PrintSelf(std::ostream & os, Indent indent) const;
191 
193  SizeValueType m_JointHistogramTotalCount;
194 
195 private:
196  JointHistogramMutualInformationImageToImageMetricv4(const Self &); //purposely not implemented
197  void operator=(const Self &); //purposely not implemented
198 
200  typename MarginalPDFType::Pointer m_FixedImageMarginalPDF;
201 
203  typename MarginalPDFType::Pointer m_MovingImageMarginalPDF;
204 
206  mutable typename JointPDFType::Pointer m_JointPDF;
207 
209  InternalComputationValueType m_VarianceForJointPDFSmoothing;
210 
212  SizeValueType m_NumberOfHistogramBins;
213  InternalComputationValueType m_FixedImageTrueMin;
214  InternalComputationValueType m_FixedImageTrueMax;
215  InternalComputationValueType m_MovingImageTrueMin;
216  InternalComputationValueType m_MovingImageTrueMax;
217  InternalComputationValueType m_FixedImageBinSize;
218  InternalComputationValueType m_MovingImageBinSize;
219 
221  JointPDFSpacingType m_JointPDFSpacing;
222 
225 
226 };
227 
228 } // end namespace itk
229 
230 #ifndef ITK_MANUAL_INSTANTIATION
231 #include "itkJointHistogramMutualInformationImageToImageMetricv4.hxx"
232 #endif
233 
234 #endif
235