Go to the documentation of this file.00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017 #ifndef __itkGradientDifferenceImageToImageMetric_h
00018 #define __itkGradientDifferenceImageToImageMetric_h
00019
00020 #include "itkImageToImageMetric.h"
00021
00022 #include "itkSobelOperator.h"
00023 #include "itkNeighborhoodOperatorImageFilter.h"
00024 #include "itkPoint.h"
00025 #include "itkCastImageFilter.h"
00026 #include "itkResampleImageFilter.h"
00027
00028 namespace itk
00029 {
00055 template < class TFixedImage, class TMovingImage >
00056 class ITK_EXPORT GradientDifferenceImageToImageMetric :
00057 public ImageToImageMetric< TFixedImage, TMovingImage>
00058 {
00059 public:
00060
00062 typedef GradientDifferenceImageToImageMetric Self;
00063 typedef ImageToImageMetric<TFixedImage, TMovingImage > Superclass;
00064
00065 typedef SmartPointer<Self> Pointer;
00066 typedef SmartPointer<const Self> ConstPointer;
00067
00069 itkNewMacro(Self);
00070
00072 itkTypeMacro(GradientDifferenceImageToImageMetric, ImageToImageMetric);
00073
00074
00076
00077 #if defined(_MSC_VER) && (_MSC_VER == 1300)
00078 typedef double RealType;
00079 #else
00080 typedef typename Superclass::RealType RealType;
00081 #endif
00082 typedef typename Superclass::TransformType TransformType;
00083 typedef typename Superclass::TransformPointer TransformPointer;
00084 typedef typename Superclass::TransformParametersType TransformParametersType;
00085 typedef typename Superclass::TransformJacobianType TransformJacobianType;
00086
00087 typedef typename Superclass::MeasureType MeasureType;
00088 typedef typename Superclass::DerivativeType DerivativeType;
00089 typedef typename Superclass::FixedImageType FixedImageType;
00090 typedef typename Superclass::MovingImageType MovingImageType;
00091 typedef typename Superclass::FixedImageConstPointer FixedImageConstPointer;
00092 typedef typename Superclass::MovingImageConstPointer MovingImageConstPointer;
00093
00094 typedef typename TFixedImage::PixelType FixedImagePixelType;
00095 typedef typename TMovingImage::PixelType MovedImagePixelType;
00096
00097 itkStaticConstMacro(FixedImageDimension, unsigned int, TFixedImage::ImageDimension);
00099 typedef itk::Image< FixedImagePixelType,
00100 itkGetStaticConstMacro( FixedImageDimension ) >
00101 TransformedMovingImageType;
00102
00103 typedef itk::ResampleImageFilter< MovingImageType, TransformedMovingImageType >
00104 TransformMovingImageFilterType;
00105
00108 typedef itk::Image< RealType,
00109 itkGetStaticConstMacro( FixedImageDimension ) >
00110 FixedGradientImageType;
00111
00112 typedef itk::CastImageFilter< FixedImageType, FixedGradientImageType >
00113 CastFixedImageFilterType;
00114 typedef typename CastFixedImageFilterType::Pointer CastFixedImageFilterPointer;
00115
00116 typedef typename FixedGradientImageType::PixelType FixedGradientPixelType;
00117
00118
00121 itkStaticConstMacro( MovedImageDimension, unsigned int,
00122 MovingImageType::ImageDimension );
00123
00124 typedef itk::Image< RealType,
00125 itkGetStaticConstMacro( MovedImageDimension ) >
00126 MovedGradientImageType;
00127
00128 typedef itk::CastImageFilter< TransformedMovingImageType, MovedGradientImageType >
00129 CastMovedImageFilterType;
00130 typedef typename CastMovedImageFilterType::Pointer CastMovedImageFilterPointer;
00131
00132 typedef typename MovedGradientImageType::PixelType MovedGradientPixelType;
00133
00134
00136 void GetDerivative( const TransformParametersType & parameters,
00137 DerivativeType & derivative ) const;
00138
00140 MeasureType GetValue( const TransformParametersType & parameters ) const;
00141
00143 void GetValueAndDerivative( const TransformParametersType & parameters,
00144 MeasureType& Value, DerivativeType& derivative ) const;
00145
00148 virtual void Initialize(void) throw ( ExceptionObject );
00149
00151 void WriteGradientImagesToFiles(void) const;
00152
00155 itkSetMacro( DerivativeDelta, double );
00156 itkGetConstReferenceMacro( DerivativeDelta, double );
00158
00159 protected:
00160 GradientDifferenceImageToImageMetric();
00161 virtual ~GradientDifferenceImageToImageMetric() {};
00162 void PrintSelf(std::ostream& os, Indent indent) const;
00163
00165 void ComputeMovedGradientRange( void ) const;
00166
00168 void ComputeVariance( void ) const;
00169
00171 MeasureType ComputeMeasure( const TransformParametersType ¶meters,
00172 const double *subtractionFactor ) const;
00173
00174 typedef NeighborhoodOperatorImageFilter<
00175 FixedGradientImageType, FixedGradientImageType > FixedSobelFilter;
00176
00177 typedef NeighborhoodOperatorImageFilter<
00178 MovedGradientImageType, MovedGradientImageType > MovedSobelFilter;
00179
00180 private:
00181 GradientDifferenceImageToImageMetric(const Self&);
00182 void operator=(const Self&);
00183
00185 mutable MovedGradientPixelType m_Variance[FixedImageDimension];
00186
00188 mutable MovedGradientPixelType m_MinMovedGradient[MovedImageDimension];
00189 mutable MovedGradientPixelType m_MaxMovedGradient[MovedImageDimension];
00190
00192 mutable FixedGradientPixelType m_MinFixedGradient[FixedImageDimension];
00193 mutable FixedGradientPixelType m_MaxFixedGradient[FixedImageDimension];
00194
00196 typename TransformMovingImageFilterType::Pointer m_TransformMovingImageFilter;
00197
00199 CastFixedImageFilterPointer m_CastFixedImageFilter;
00200
00201 SobelOperator< FixedGradientPixelType,
00202 itkGetStaticConstMacro(FixedImageDimension) >
00203 m_FixedSobelOperators[FixedImageDimension];
00204
00205 typename FixedSobelFilter::Pointer m_FixedSobelFilters[itkGetStaticConstMacro( FixedImageDimension )];
00206
00207 ZeroFluxNeumannBoundaryCondition< MovedGradientImageType > m_MovedBoundCond;
00208 ZeroFluxNeumannBoundaryCondition< FixedGradientImageType > m_FixedBoundCond;
00209
00211 CastMovedImageFilterPointer m_CastMovedImageFilter;
00212
00213 SobelOperator< MovedGradientPixelType,
00214 itkGetStaticConstMacro(MovedImageDimension) >
00215 m_MovedSobelOperators[MovedImageDimension];
00216
00217 typename MovedSobelFilter::Pointer m_MovedSobelFilters[itkGetStaticConstMacro( MovedImageDimension )];
00218
00219 double m_DerivativeDelta;
00220 };
00221
00222 }
00223
00224 #ifndef ITK_MANUAL_INSTANTIATION
00225 #include "itkGradientDifferenceImageToImageMetric.txx"
00226 #endif
00227
00228 #endif
00229