ITK
4.1.0
Insight Segmentation and Registration Toolkit
|
00001 /*========================================================================= 00002 * 00003 * Copyright Insight Software Consortium 00004 * 00005 * Licensed under the Apache License, Version 2.0 (the "License"); 00006 * you may not use this file except in compliance with the License. 00007 * You may obtain a copy of the License at 00008 * 00009 * http://www.apache.org/licenses/LICENSE-2.0.txt 00010 * 00011 * Unless required by applicable law or agreed to in writing, software 00012 * distributed under the License is distributed on an "AS IS" BASIS, 00013 * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. 00014 * See the License for the specific language governing permissions and 00015 * limitations under the License. 00016 * 00017 *=========================================================================*/ 00018 #ifndef __itkGradientDifferenceImageToImageMetric_h 00019 #define __itkGradientDifferenceImageToImageMetric_h 00020 00021 #include "itkImageToImageMetric.h" 00022 00023 #include "itkSobelOperator.h" 00024 #include "itkNeighborhoodOperatorImageFilter.h" 00025 #include "itkPoint.h" 00026 #include "itkCastImageFilter.h" 00027 #include "itkResampleImageFilter.h" 00028 00029 namespace itk 00030 { 00057 template< class TFixedImage, class TMovingImage > 00058 class ITK_EXPORT GradientDifferenceImageToImageMetric: 00059 public ImageToImageMetric< TFixedImage, TMovingImage > 00060 { 00061 public: 00062 00064 typedef GradientDifferenceImageToImageMetric Self; 00065 typedef ImageToImageMetric< TFixedImage, TMovingImage > Superclass; 00066 00067 typedef SmartPointer< Self > Pointer; 00068 typedef SmartPointer< const Self > ConstPointer; 00069 00071 itkNewMacro(Self); 00072 00074 itkTypeMacro(GradientDifferenceImageToImageMetric, ImageToImageMetric); 00075 00077 typedef typename Superclass::RealType RealType; 00078 typedef typename Superclass::TransformType TransformType; 00079 typedef typename Superclass::TransformPointer TransformPointer; 00080 typedef typename Superclass::TransformParametersType TransformParametersType; 00081 typedef typename Superclass::TransformJacobianType TransformJacobianType; 00082 00083 typedef typename Superclass::MeasureType MeasureType; 00084 typedef typename Superclass::DerivativeType DerivativeType; 00085 typedef typename Superclass::FixedImageType FixedImageType; 00086 typedef typename Superclass::MovingImageType MovingImageType; 00087 typedef typename Superclass::FixedImageConstPointer FixedImageConstPointer; 00088 typedef typename Superclass::MovingImageConstPointer MovingImageConstPointer; 00089 00090 typedef typename TFixedImage::PixelType FixedImagePixelType; 00091 typedef typename TMovingImage::PixelType MovedImagePixelType; 00092 00093 itkStaticConstMacro(FixedImageDimension, unsigned int, TFixedImage::ImageDimension); 00095 typedef itk::Image< FixedImagePixelType, 00096 itkGetStaticConstMacro(FixedImageDimension) > 00097 TransformedMovingImageType; 00098 00099 typedef itk::ResampleImageFilter< MovingImageType, TransformedMovingImageType > 00100 TransformMovingImageFilterType; 00101 00104 typedef itk::Image< RealType, itkGetStaticConstMacro(FixedImageDimension) > FixedGradientImageType; 00105 00106 typedef itk::CastImageFilter< FixedImageType, FixedGradientImageType > CastFixedImageFilterType; 00107 typedef typename CastFixedImageFilterType::Pointer CastFixedImageFilterPointer; 00108 00109 typedef typename FixedGradientImageType::PixelType FixedGradientPixelType; 00110 00113 itkStaticConstMacro(MovedImageDimension, unsigned int, MovingImageType::ImageDimension); 00114 00115 typedef itk::Image< RealType, itkGetStaticConstMacro(MovedImageDimension) > MovedGradientImageType; 00116 00117 typedef itk::CastImageFilter< TransformedMovingImageType, MovedGradientImageType > CastMovedImageFilterType; 00118 typedef typename CastMovedImageFilterType::Pointer CastMovedImageFilterPointer; 00119 00120 typedef typename MovedGradientImageType::PixelType MovedGradientPixelType; 00121 00123 void GetDerivative(const TransformParametersType & parameters, 00124 DerivativeType & derivative) const; 00125 00127 MeasureType GetValue(const TransformParametersType & parameters) const; 00128 00130 void GetValueAndDerivative(const TransformParametersType & parameters, 00131 MeasureType & Value, DerivativeType & derivative) const; 00132 00135 virtual void Initialize(void) 00136 throw ( ExceptionObject ); 00137 00139 void WriteGradientImagesToFiles(void) const; 00140 00143 itkSetMacro(DerivativeDelta, double); 00144 itkGetConstReferenceMacro(DerivativeDelta, double); 00145 protected: 00146 GradientDifferenceImageToImageMetric(); 00147 virtual ~GradientDifferenceImageToImageMetric() {} 00148 void PrintSelf(std::ostream & os, Indent indent) const; 00150 00152 void ComputeMovedGradientRange(void) const; 00153 00155 void ComputeVariance(void) const; 00156 00158 MeasureType ComputeMeasure(const TransformParametersType & parameters, 00159 const double *subtractionFactor) const; 00160 00161 typedef NeighborhoodOperatorImageFilter< 00162 FixedGradientImageType, FixedGradientImageType > FixedSobelFilter; 00163 00164 typedef NeighborhoodOperatorImageFilter< 00165 MovedGradientImageType, MovedGradientImageType > MovedSobelFilter; 00166 private: 00167 GradientDifferenceImageToImageMetric(const Self &); //purposely not 00168 // implemented 00169 void operator=(const Self &); //purposely not 00170 00171 // implemented 00172 00174 mutable MovedGradientPixelType m_Variance[FixedImageDimension]; 00175 00177 mutable MovedGradientPixelType m_MinMovedGradient[MovedImageDimension]; 00178 mutable MovedGradientPixelType m_MaxMovedGradient[MovedImageDimension]; 00179 00181 mutable FixedGradientPixelType m_MinFixedGradient[FixedImageDimension]; 00182 mutable FixedGradientPixelType m_MaxFixedGradient[FixedImageDimension]; 00183 00185 typename TransformMovingImageFilterType::Pointer m_TransformMovingImageFilter; 00186 00188 CastFixedImageFilterPointer m_CastFixedImageFilter; 00189 00190 SobelOperator< FixedGradientPixelType, 00191 itkGetStaticConstMacro(FixedImageDimension) > 00192 m_FixedSobelOperators[FixedImageDimension]; 00193 00194 typename FixedSobelFilter::Pointer m_FixedSobelFilters[itkGetStaticConstMacro(FixedImageDimension)]; 00195 00196 ZeroFluxNeumannBoundaryCondition< MovedGradientImageType > m_MovedBoundCond; 00197 ZeroFluxNeumannBoundaryCondition< FixedGradientImageType > m_FixedBoundCond; 00198 00200 CastMovedImageFilterPointer m_CastMovedImageFilter; 00201 00202 SobelOperator< MovedGradientPixelType, 00203 itkGetStaticConstMacro(MovedImageDimension) > 00204 m_MovedSobelOperators[MovedImageDimension]; 00205 00206 typename MovedSobelFilter::Pointer m_MovedSobelFilters[itkGetStaticConstMacro(MovedImageDimension)]; 00207 00208 double m_DerivativeDelta; 00209 }; 00210 } // end namespace itk 00211 00212 #ifndef ITK_MANUAL_INSTANTIATION 00213 #include "itkGradientDifferenceImageToImageMetric.hxx" 00214 #endif 00215 00216 #endif 00217