Main Page   Groups   Namespace List   Class Hierarchy   Alphabetical List   Compound List   File List   Namespace Members   Compound Members   File Members   Concepts

itkBSplineInterpolateImageFunction.h

Go to the documentation of this file.
00001 /*=========================================================================
00002 
00003   Program:   Insight Segmentation & Registration Toolkit
00004   Module:    $RCSfile: itkBSplineInterpolateImageFunction.h,v $
00005   Language:  C++
00006   Date:      $Date: 2009-04-25 12:27:05 $
00007   Version:   $Revision: 1.24 $
00008 
00009   Copyright (c) Insight Software Consortium. All rights reserved.
00010   See ITKCopyright.txt or http://www.itk.org/HTML/Copyright.htm for details.
00011 
00012   Portions of this code are covered under the VTK copyright.
00013   See VTKCopyright.txt or http://www.kitware.com/VTKCopyright.htm for details.
00014 
00015      This software is distributed WITHOUT ANY WARRANTY; without even 
00016      the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR 
00017      PURPOSE.  See the above copyright notices for more information.
00018 
00019 =========================================================================*/
00020 #ifndef __itkBSplineInterpolateImageFunction_h
00021 #define __itkBSplineInterpolateImageFunction_h
00022 
00023 // First make sure that the configuration is available.
00024 // This line can be removed once the optimized versions
00025 // gets integrated into the main directories.
00026 #include "itkConfigure.h"
00027 
00028 #ifdef ITK_USE_OPTIMIZED_REGISTRATION_METHODS
00029 #include "itkOptBSplineInterpolateImageFunction.h"
00030 #else
00031 
00032 #include <vector>
00033 
00034 #include "itkImageLinearIteratorWithIndex.h"
00035 #include "itkInterpolateImageFunction.h"
00036 #include "vnl/vnl_matrix.h"
00037 
00038 #include "itkBSplineDecompositionImageFilter.h"
00039 #include "itkConceptChecking.h"
00040 #include "itkCovariantVector.h"
00041 
00042 namespace itk
00043 {
00077 template <
00078   class TImageType, 
00079   class TCoordRep = double,
00080   class TCoefficientType = double >
00081 class ITK_EXPORT BSplineInterpolateImageFunction : 
00082     public InterpolateImageFunction<TImageType,TCoordRep> 
00083 {
00084 public:
00086   typedef BSplineInterpolateImageFunction                 Self;
00087   typedef InterpolateImageFunction<TImageType,TCoordRep>  Superclass;
00088   typedef SmartPointer<Self>                              Pointer;
00089   typedef SmartPointer<const Self>                        ConstPointer;
00090 
00092   itkTypeMacro(BSplineInterpolateImageFunction, InterpolateImageFunction);
00093 
00094 
00096   itkNewMacro( Self );
00097 
00099   typedef typename Superclass::OutputType OutputType;
00100 
00102   typedef typename Superclass::InputImageType InputImageType;
00103 
00105   itkStaticConstMacro(ImageDimension, unsigned int,Superclass::ImageDimension);
00106 
00108   typedef typename Superclass::IndexType IndexType;
00109 
00111   typedef typename Superclass::ContinuousIndexType ContinuousIndexType;
00112 
00114   typedef typename Superclass::PointType PointType;
00115 
00117   typedef ImageLinearIteratorWithIndex<TImageType> Iterator;
00118 
00120   typedef TCoefficientType CoefficientDataType;
00121   typedef Image<CoefficientDataType, 
00122                      itkGetStaticConstMacro(ImageDimension)
00123     >                      CoefficientImageType;
00124 
00126   typedef BSplineDecompositionImageFilter<TImageType, CoefficientImageType> 
00127   CoefficientFilter;
00128 
00129   typedef typename CoefficientFilter::Pointer CoefficientFilterPointer;
00130 
00139   virtual OutputType EvaluateAtContinuousIndex( 
00140     const ContinuousIndexType & index ) const; 
00141 
00143   typedef CovariantVector<OutputType,
00144                           itkGetStaticConstMacro(ImageDimension)
00145     > CovariantVectorType;
00146 
00147   CovariantVectorType EvaluateDerivative( const PointType & point ) const
00148     {
00149     ContinuousIndexType index;
00150     this->GetInputImage()->TransformPhysicalPointToContinuousIndex( point, index );
00151     return ( this->EvaluateDerivativeAtContinuousIndex( index ) );
00152     }
00153 
00154   CovariantVectorType EvaluateDerivativeAtContinuousIndex( 
00155     const ContinuousIndexType & x ) const;
00156 
00157 
00160   void SetSplineOrder(unsigned int SplineOrder);
00161   itkGetConstMacro(SplineOrder, int);
00163 
00164 
00166   virtual void SetInputImage(const TImageType * inputData);
00167 
00168 
00181   itkSetMacro( UseImageDirection, bool );
00182   itkGetConstMacro( UseImageDirection, bool );
00183   itkBooleanMacro( UseImageDirection );
00185 
00186 
00187 protected:
00188   BSplineInterpolateImageFunction();
00189   virtual ~BSplineInterpolateImageFunction() {};
00190   void PrintSelf(std::ostream& os, Indent indent) const;
00191 
00192   // These are needed by the smoothing spline routine.
00193   std::vector<CoefficientDataType>    m_Scratch;        // temp storage for processing of Coefficients
00194   typename TImageType::SizeType       m_DataLength;  // Image size
00195   unsigned int                        m_SplineOrder; // User specified spline order (3rd or cubic is the default)
00196 
00197   typename CoefficientImageType::ConstPointer       m_Coefficients; // Spline coefficients  
00198 
00199 private:
00200   BSplineInterpolateImageFunction( const Self& ); //purposely not implemented
00201   void operator=( const Self& ); //purposely not implemented
00202 
00204   void SetInterpolationWeights( const ContinuousIndexType & x, 
00205                                 const vnl_matrix<long> & EvaluateIndex, 
00206                                 vnl_matrix<double> & weights, 
00207                                 unsigned int splineOrder ) const;
00208 
00210   void SetDerivativeWeights( const ContinuousIndexType & x, 
00211                              const vnl_matrix<long> & EvaluateIndex, 
00212                              vnl_matrix<double> & weights, 
00213                              unsigned int splineOrder ) const;
00214 
00217   void GeneratePointsToIndex(  );
00218 
00220   void DetermineRegionOfSupport( vnl_matrix<long> & evaluateIndex, 
00221                                  const ContinuousIndexType & x, 
00222                                  unsigned int splineOrder ) const;
00223 
00226   void ApplyMirrorBoundaryConditions(vnl_matrix<long> & evaluateIndex, 
00227                                      unsigned int splineOrder) const;
00228 
00229 
00230   Iterator                  m_CIterator;    // Iterator for traversing spline coefficients.
00231   unsigned long             m_MaxNumberInterpolationPoints; // number of neighborhood points used for interpolation
00232   std::vector<IndexType>    m_PointsToIndex;  // Preallocation of interpolation neighborhood indicies
00233 
00234   CoefficientFilterPointer     m_CoefficientFilter;
00235   
00236   // flag to take or not the image direction into account when computing the
00237   // derivatives.
00238   bool m_UseImageDirection;
00239 
00240 };
00241 
00242 } // namespace itk
00243 
00244 #ifndef ITK_MANUAL_INSTANTIATION
00245 #include "itkBSplineInterpolateImageFunction.txx"
00246 #endif
00247 
00248 #endif
00249 
00250 #endif
00251 

Generated at Tue Sep 15 02:28:03 2009 for ITK by doxygen 1.5.8 written by Dimitri van Heesch, © 1997-2000