ITK  4.6.0
Insight Segmentation and Registration Toolkit
itkMattesMutualInformationImageToImageMetric.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 __itkMattesMutualInformationImageToImageMetric_h
19 #define __itkMattesMutualInformationImageToImageMetric_h
20 
21 #include "itkImageToImageMetric.h"
22 #include "itkPoint.h"
23 #include "itkIndex.h"
25 #include "itkArray2D.h"
26 
27 #include "itkSimpleFastMutexLock.h"
28 
29 
30 namespace itk
31 {
111 template <typename TFixedImage, typename TMovingImage>
113  public ImageToImageMetric<TFixedImage, TMovingImage>
114 {
115 public:
116 
122 
124  itkNewMacro(Self);
125 
129 
131  typedef typename Superclass::TransformType TransformType;
145 
151 
153 
155  itkStaticConstMacro(MovingImageDimension, unsigned int,
156  MovingImageType::ImageDimension);
157 
165  virtual void Initialize(void)
166  throw ( ExceptionObject );
167 
169  MeasureType GetValue(const ParametersType & parameters) const;
170 
172  void GetDerivative(const ParametersType & parameters, DerivativeType & Derivative) const;
173 
175  void GetValueAndDerivative(const ParametersType & parameters, MeasureType & Value, DerivativeType & Derivative) const;
176 
183  itkSetClampMacro( NumberOfHistogramBins, SizeValueType,
184  5, NumericTraits<SizeValueType>::max() );
185  itkGetConstReferenceMacro(NumberOfHistogramBins, SizeValueType);
187 
212  itkSetMacro(UseExplicitPDFDerivatives, bool);
213  itkGetConstReferenceMacro(UseExplicitPDFDerivatives, bool);
214  itkBooleanMacro(UseExplicitPDFDerivatives);
216 
218  typedef double PDFValueType; //NOTE: floating point precision is not as stable. Double precision proves faster and more robust in real-world testing.
219 
221  typedef Image<PDFValueType, 2> JointPDFType;
222  typedef Image<PDFValueType, 3> JointPDFDerivativesType;
223 
228  const typename JointPDFType::Pointer GetJointPDF () const
229  {
230  if( this->m_MMIMetricPerThreadVariables == ITK_NULLPTR )
231  {
232  return JointPDFType::Pointer(ITK_NULLPTR);
233  }
234  return this->m_MMIMetricPerThreadVariables[0].JointPDF;
235  }
237 
245  {
246  if( this->m_MMIMetricPerThreadVariables == ITK_NULLPTR )
247  {
248  return JointPDFDerivativesType::Pointer(ITK_NULLPTR);
249  }
250  return this->m_MMIMetricPerThreadVariables[0].JointPDFDerivatives;
251  }
253 
254 protected:
255 
258  void PrintSelf(std::ostream & os, Indent indent) const;
259 
260 private:
261 
262  // purposely not implemented
264  // purposely not implemented
265  void operator=(const Self &);
266 
275 
279 
282 
284  void ComputePDFDerivatives(ThreadIdType threadId, unsigned int sampleNumber, int movingImageParzenWindowIndex,
286  & movingImageGradientValue,
287  PDFValueType cubicBSplineDerivativeValue) const;
288 
289  virtual void GetValueThreadPreProcess(ThreadIdType threadId, bool withinSampleThread) const;
290  virtual void GetValueThreadPostProcess(ThreadIdType threadId, bool withinSampleThread) const;
291  //NOTE: The signature in base class requires that movingImageValue is of type double
292  virtual bool GetValueThreadProcessSample(ThreadIdType threadId, SizeValueType fixedImageSample,
293  const MovingImagePointType & mappedPoint,
294  double movingImageValue) const;
295 
296  virtual void GetValueAndDerivativeThreadPreProcess( ThreadIdType threadId, bool withinSampleThread) const;
297  virtual void GetValueAndDerivativeThreadPostProcess( ThreadIdType threadId, bool withinSampleThread) const;
298  //NOTE: The signature in base class requires that movingImageValue is of type double
299  virtual bool GetValueAndDerivativeThreadProcessSample(ThreadIdType threadId, SizeValueType fixedImageSample,
300  const MovingImagePointType & mappedPoint,
301  double movingImageValue, const ImageDerivativesType &
302  movingImageGradientValue) const;
303 
314 
318 
322 
324 
326  typedef std::vector< PDFValueType > MarginalPDFType;
328 
330  {
333 
335 
338 
342 
344 
346  };
347 
348 #if !defined(__GCCXML__)
349  itkPadStruct( ITK_CACHE_LINE_ALIGNMENT, MMIMetricPerThreadStruct,
350  PaddedMMIMetricPerThreadStruct);
351  itkAlignedTypedef( ITK_CACHE_LINE_ALIGNMENT, PaddedMMIMetricPerThreadStruct,
352  AlignedMMIMetricPerThreadStruct );
353  // Due to a bug in older version of Visual Studio where std::vector resize
354  // uses a value instead of a const reference, this must be a pointer.
355  // See
356  // http://thetweaker.wordpress.com/2010/05/05/stdvector-of-aligned-elements/
357  // http://connect.microsoft.com/VisualStudio/feedback/details/692988
358  mutable AlignedMMIMetricPerThreadStruct * m_MMIMetricPerThreadVariables;
359 #endif
360 
363 };
364 } // end namespace itk
365 
366 #ifndef ITK_MANUAL_INSTANTIATION
367 #include "itkMattesMutualInformationImageToImageMetric.hxx"
368 #endif
369 
370 #endif
Array class with size defined at construction time.
Definition: itkArray.h:50
virtual bool GetValueThreadProcessSample(ThreadIdType threadId, SizeValueType fixedImageSample, const MovingImagePointType &mappedPoint, double movingImageValue) const
virtual bool GetValueAndDerivativeThreadProcessSample(ThreadIdType threadId, SizeValueType fixedImageSample, const MovingImagePointType &mappedPoint, double movingImageValue, const ImageDerivativesType &movingImageGradientValue) const
Superclass::RegionType RegionType
Definition: itkImage.h:140
itkAlignedTypedef(ITK_CACHE_LINE_ALIGNMENT, PaddedMMIMetricPerThreadStruct, AlignedMMIMetricPerThreadStruct)
Light weight base class for most itk classes.
const JointPDFDerivativesType::Pointer GetJointPDFDerivatives() const
BSplineTransformType::WeightsType BSplineTransformWeightsType
SmartPointer< Self > Pointer
Definition: itkImage.h:81
void ComputePDFDerivatives(ThreadIdType threadId, unsigned int sampleNumber, int movingImageParzenWindowIndex, const ImageDerivativesType &movingImageGradientValue, PDFValueType cubicBSplineDerivativeValue) const
Superclass::ParametersValueType CoordinateRepresentationType
Represent the size (bounds) of a n-dimensional image.
Definition: itkSize.h:52
signed long OffsetValueType
Definition: itkIntTypes.h:154
MeasureType GetValue(const ParametersType &parameters) const
TransformType::Pointer TransformPointer
CovariantVector< double, itkGetStaticConstMacro(MovingImageDimension) > ImageDerivativesType
virtual void GetValueAndDerivativeThreadPreProcess(ThreadIdType threadId, bool withinSampleThread) const
Superclass::DerivativeType DerivativeType
MovingImageType::ConstPointer MovingImageConstPointer
TransformType::OutputPointType MovingImagePointType
CubicBSplineDerivativeFunctionType::Pointer m_CubicBSplineDerivativeKernel
unsigned long SizeValueType
Definition: itkIntTypes.h:143
TPixel PixelType
Definition: itkImage.h:89
Derivative of a BSpline kernel used for density estimation and nonparameteric regression.
BSpline kernel used for density estimation and nonparameteric regression.
virtual void GetValueThreadPreProcess(ThreadIdType threadId, bool withinSampleThread) const
Transform points and vectors from an input space to an output space.
Definition: itkTransform.h:82
Superclass::ParametersType ParametersType
Superclass::IndexType IndexType
Definition: itkImage.h:122
Standard exception handling object.
void GetDerivative(const ParametersType &parameters, DerivativeType &Derivative) const
BSplineTransformIndexArrayType::ValueType IndexValueType
void ComputeFixedImageParzenWindowIndices(FixedImageSampleContainer &samples)
std::vector< FixedImageSamplePoint > FixedImageSampleContainer
TransformType::JacobianType TransformJacobianType
virtual void GetValueThreadPostProcess(ThreadIdType threadId, bool withinSampleThread) const
Computes the mutual information between two images to be registered using the method of Mattes et al...
BSplineTransformWeightsType::ValueType WeightsValueType
InterpolateImageFunction< MovingImageType, CoordinateRepresentationType > InterpolatorType
Control indentation during Print() invocation.
Definition: itkIndent.h:49
FixedImageType::ConstPointer FixedImageConstPointer
BSplineDerivativeKernelFunction< 3, PDFValueType > CubicBSplineDerivativeFunctionType
virtual void GetValueAndDerivativeThreadPostProcess(ThreadIdType threadId, bool withinSampleThread) const
Define additional traits for native types such as int or float.
Superclass::MeasureType MeasureType
Computes similarity between regions of two images.
A templated class holding a geometric point in n-Dimensional space.
Definition: itkPoint.h:51
A templated class holding a n-Dimensional covariant vector.
Templated n-dimensional image class.
Definition: itkImage.h:75
itkPadStruct(ITK_CACHE_LINE_ALIGNMENT, MMIMetricPerThreadStruct, PaddedMMIMetricPerThreadStruct)
void PrintSelf(std::ostream &os, Indent indent) const
void GetValueAndDerivative(const ParametersType &parameters, MeasureType &Value, DerivativeType &Derivative) const
unsigned int ThreadIdType
Definition: itkIntTypes.h:159
BSplineTransformType::ParameterIndexArrayType BSplineTransformIndexArrayType