ITK  4.1.0
Insight Segmentation and Registration Toolkit
itkFEMImageMetricLoad.h
Go to the documentation of this file.
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 __itkFEMImageMetricLoad_h
00019 #define __itkFEMImageMetricLoad_h
00020 
00021 #include "itkFEMLoadElementBase.h"
00022 
00023 #include "itkImage.h"
00024 #include "itkTranslationTransform.h"
00025 
00026 #include "itkImageRegionIteratorWithIndex.h"
00027 #include "itkNeighborhoodIterator.h"
00028 #include "itkNeighborhoodIterator.h"
00029 #include "itkNeighborhoodInnerProduct.h"
00030 #include "itkDerivativeOperator.h"
00031 #include "itkForwardDifferenceOperator.h"
00032 #include "itkLinearInterpolateImageFunction.h"
00033 #include "vnl/vnl_math.h"
00034 
00035 #include <itkMutualInformationImageToImageMetric.h>
00036 #include <itkMattesMutualInformationImageToImageMetric.h>
00037 #include <itkMeanSquaresImageToImageMetric.h>
00038 #include <itkNormalizedCorrelationImageToImageMetric.h>
00039 
00040 namespace itk
00041 {
00042 namespace fem
00043 {
00068 template <class TMoving, class TFixed>
00069 class ImageMetricLoad : public LoadElement
00070 {
00071 public:
00072 
00074   typedef ImageMetricLoad          Self;
00075   typedef LoadElement              Superclass;
00076   typedef SmartPointer<Self>       Pointer;
00077   typedef SmartPointer<const Self> ConstPointer;
00078 
00080   itkSimpleNewMacro(Self);
00081 
00083   itkTypeMacro(ImageMetricLoad, LoadElement);
00084 
00087   virtual::itk::LightObject::Pointer CreateAnother(void) const;
00088 
00089   // Necessary typedefs for dealing with images BEGIN
00090   typedef typename LoadElement::Float Float;
00091 
00092   typedef TMoving                           MovingType;
00093   typedef typename MovingType::ConstPointer MovingConstPointer;
00094   typedef MovingType *                      MovingPointer;
00095   typedef TFixed                            FixedType;
00096   typedef FixedType *                       FixedPointer;
00097   typedef typename FixedType::ConstPointer  FixedConstPointer;
00098 
00100   itkStaticConstMacro(ImageDimension, unsigned int,
00101                       MovingType::ImageDimension);
00102 
00103   typedef ImageRegionIteratorWithIndex<MovingType> RefRegionIteratorType;
00104   typedef ImageRegionIteratorWithIndex<FixedType>  TarRegionIteratorType;
00105 
00106   typedef NeighborhoodIterator<MovingType>
00107   MovingNeighborhoodIteratorType;
00108   typedef typename MovingNeighborhoodIteratorType::IndexType
00109   MovingNeighborhoodIndexType;
00110   typedef typename MovingNeighborhoodIteratorType::RadiusType
00111   MovingRadiusType;
00112   typedef NeighborhoodIterator<FixedType>
00113   FixedNeighborhoodIteratorType;
00114   typedef typename FixedNeighborhoodIteratorType::IndexType
00115   FixedNeighborhoodIndexType;
00116   typedef typename FixedNeighborhoodIteratorType::RadiusType
00117   FixedRadiusType;
00118 
00119 // IMAGE DATA
00120   typedef   typename  MovingType::PixelType                             RefPixelType;
00121   typedef   typename  FixedType::PixelType                              TarPixelType;
00122   typedef   Float                                                       PixelType;
00123   typedef   Float                                                       ComputationType;
00124   typedef   Image<RefPixelType, itkGetStaticConstMacro(ImageDimension)> RefImageType;
00125   typedef   Image<TarPixelType, itkGetStaticConstMacro(ImageDimension)> TarImageType;
00126   typedef   Image<PixelType, itkGetStaticConstMacro(ImageDimension)>    ImageType;
00127   typedef   vnl_vector<Float>                                           VectorType;
00128 
00129 // Necessary typedefs for dealing with images END
00130 
00131 // ------------------------------------------------------------
00132 // Set up the metrics
00133 // ------------------------------------------------------------
00134   typedef double
00135   CoordinateRepresentationType;
00136   typedef Transform<CoordinateRepresentationType, itkGetStaticConstMacro(ImageDimension),
00137                     itkGetStaticConstMacro(ImageDimension)>            TransformBaseType;
00138   typedef TranslationTransform<CoordinateRepresentationType,
00139                                itkGetStaticConstMacro(ImageDimension)> DefaultTransformType;
00140 
00142   typedef   ImageToImageMetric<FixedType, MovingType> MetricBaseType;
00143   typedef typename MetricBaseType::Pointer            MetricBaseTypePointer;
00144 
00145   typedef   MutualInformationImageToImageMetric<MovingType, FixedType> MutualInformationMetricType;
00146 
00147   typedef   MeanSquaresImageToImageMetric<MovingType, FixedType> MeanSquaresMetricType;
00148 
00149   typedef   NormalizedCorrelationImageToImageMetric<MovingType, FixedType> NormalizedCorrelationMetricType;
00150 
00151   typedef  MeanSquaresMetricType                        DefaultMetricType;
00152   typedef typename DefaultTransformType::ParametersType ParametersType;
00153   typedef typename DefaultTransformType::JacobianType   JacobianType;
00154 
00155   typedef unsigned long                                        ElementIdentifier;
00156   typedef VectorContainer<ElementIdentifier, Element::Pointer> ElementContainerType;
00157 // ------------------------------------------------------------
00158 // Set up an Interpolator
00159 // ------------------------------------------------------------
00160   typedef LinearInterpolateImageFunction<MovingType, double> InterpolatorType;
00161 
00163   typedef float RealType;
00164   typedef CovariantVector<RealType,
00165                           itkGetStaticConstMacro(ImageDimension)> GradientPixelType;
00166   typedef Image<GradientPixelType,
00167                 itkGetStaticConstMacro(ImageDimension)> GradientImageType;
00168   typedef SmartPointer<GradientImageType> GradientImagePointer;
00169   typedef GradientRecursiveGaussianImageFilter<ImageType,
00170                                                GradientImageType>
00171   GradientImageFilterType;
00172   //  typedef typename GradientImageFilterType::Pointer
00173   // GradientImageFilterPointer;
00175 
00176 // FUNCTIONS
00177 
00179   void SetMetric(MetricBaseTypePointer MP)
00180   {
00181     m_Metric = MP;
00182   }
00183 
00185   void SetMovingImage(MovingType *R)
00186   {
00187     m_RefImage = R;
00188     m_RefSize = m_RefImage->GetLargestPossibleRegion().GetSize();
00189   }
00191 
00192   void SetMetricMovingImage(MovingType *R)
00193   {
00194     m_Metric->SetMovingImage(R);
00195     m_RefSize = R->GetLargestPossibleRegion().GetSize();
00196   }
00197 
00199   void SetFixedImage(FixedType *T)
00200   {
00201     m_TarImage = T;
00202     m_TarSize = T->GetLargestPossibleRegion().GetSize();
00203   }
00205 
00206   void SetMetricFixedImage(FixedType *T)
00207   {
00208     m_Metric->SetFixedImage(T);
00209     m_TarSize = T->GetLargestPossibleRegion().GetSize();
00210   }
00211 
00212   MovingPointer GetMovingImage()
00213   {
00214     return m_RefImage;
00215   }
00216   FixedPointer GetFixedImage()
00217   {
00218     return m_TarImage;
00219   }
00220 
00222   void SetMetricRadius(MovingRadiusType T)
00223   {
00224     m_MetricRadius  = T;
00225   }
00226 
00228   MovingRadiusType GetMetricRadius()
00229   {
00230     return m_MetricRadius;
00231   }
00232 
00237   void SetNumberOfIntegrationPoints(unsigned int i)
00238   {
00239     m_NumberOfIntegrationPoints = i;
00240   }
00241   unsigned int GetNumberOfIntegrationPoints()
00242   {
00243     return m_NumberOfIntegrationPoints;
00244   }
00246 
00250   void SetSign(Float s)
00251   {
00252     m_Sign = s;
00253   }
00254 
00256   void SetTemp(Float s)
00257   {
00258     m_Temp = s;
00259   }
00260 
00262   void SetGamma(Float s)
00263   {
00264     m_Gamma = s;
00265   }
00266 
00270   void SetSolution(Solution::ConstPointer ptr)
00271   {
00272     m_Solution = ptr;
00273   }
00274 
00278   Solution::ConstPointer GetSolution()
00279   {
00280     return m_Solution;
00281   }
00282 
00286   Float GetMetric(VectorType InVec);
00287 
00288   VectorType GetPolynomialFitToMetric(VectorType PositionInElement, VectorType SolutionAtPosition);
00289 
00290   VectorType MetricFiniteDiff(VectorType PositionInElement, VectorType SolutionAtPosition);
00291 
00292   // FIXME - WE ASSUME THE 2ND VECTOR (INDEX 1) HAS THE INFORMATION WE WANT
00293   Float GetSolution(unsigned int i, unsigned int which = 0)
00294   {
00295     return m_Solution->GetSolutionValue(i, which);
00296   }
00297 
00298 // define the copy constructor
00299 //  ImageMetricLoad(const ImageMetricLoad& LMS);
00300 
00301   void InitializeMetric(void);
00302 
00303   ImageMetricLoad(); // cannot be private until we always use smart pointers
00304   Float EvaluateMetricGivenSolution(Element::ArrayType *el, Float step = 1.0);
00305 
00306   Float EvaluateMetricGivenSolution1(Element::ArrayType *el, Float step = 1.0);
00307 
00311   VectorType Fe(VectorType, VectorType);
00312 
00313   static Baseclass * NewImageMetricLoad(void)
00314   {
00315     return new ImageMetricLoad;
00316   }
00317 
00319   // void InitializeGradientImage();
00320   void SetMetricGradientImage(GradientImageType *g)
00321   {
00322     m_MetricGradientImage = g;
00323   }
00324   GradientImageType * GetMetricGradientImage()
00325   {
00326     return m_MetricGradientImage;
00327   }
00329 
00330   void PrintCurrentEnergy()
00331   {
00332     std::cout << " energy " << m_Energy << std::endl;
00333   }
00334   double GetCurrentEnergy()
00335   {
00336     return m_Energy;
00337   }
00338   void  SetCurrentEnergy(double e)
00339   {
00340     m_Energy = e;
00341   }
00342 
00343   // FIXME - Documentation
00344   virtual void ApplyLoad(Element::ConstPointer element, Element::VectorType & Fe);
00345 
00346 protected:
00347   virtual void PrintSelf(std::ostream& os, Indent indent) const;
00348 
00349 private:
00350   GradientImageType *m_MetricGradientImage;
00351   MovingPointer      m_RefImage;
00352   FixedPointer       m_TarImage;
00355   typename MovingType::SizeType m_RefSize;
00356   typename FixedType::SizeType  m_TarSize;
00357   unsigned int                  m_NumberOfIntegrationPoints;
00358   unsigned int                  m_SolutionIndex;
00359   unsigned int                  m_SolutionIndex2;
00360   Float                         m_Sign;
00361   Float                         m_Temp;
00362   Float                         m_Gamma;
00363 
00364   typename Solution::ConstPointer     m_Solution;
00365   MetricBaseTypePointer               m_Metric;
00366   typename TransformBaseType::Pointer m_Transform;
00367   typename InterpolatorType::Pointer  m_Interpolator;
00368 
00369   mutable double m_Energy;
00370 private:
00371 
00372 };
00373 }
00374 }  // end namespace fem/itk
00375 
00376 #ifndef ITK_MANUAL_INSTANTIATION
00377 #include "itkFEMImageMetricLoad.hxx"
00378 #endif
00379 
00380 #endif
00381