itk/clitkVectorBSplineInterpolateImageFunction.h

   1 /*=========================================================================
   2   Program:   vv                     http://www.creatis.insa-lyon.fr/rio/vv
   3
   4   Authors belong to:
   5   - University of LYON              http://www.universite-lyon.fr/
   6   - Léon Bérard cancer center       http://oncora1.lyon.fnclcc.fr
   7   - CREATIS CNRS laboratory         http://www.creatis.insa-lyon.fr
   8
   9   This software is distributed WITHOUT ANY WARRANTY; without even
  10   the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
  11   PURPOSE.  See the copyright notices for more information.
  12
  13   It is distributed under dual licence
  14
  15   - BSD        See included LICENSE.txt file
  16   - CeCILL-B   http://www.cecill.info/licences/Licence_CeCILL-B_V1-en.html
  17 ======================================================================-====*/
  18 #ifndef __clitkVectorBSplineInterpolateImageFunction_h
  19 #define __clitkVectorBSplineInterpolateImageFunction_h
  20 #include "clitkVectorBSplineDecompositionImageFilter.h"
  21
  22 // First make sure that the configuration is available.
  23 // This line can be removed once the optimized versions
  24 // gets integrated into the main directories.
  25 #include "itkConfigure.h"
  26
  27 // #ifdef ITK_USE_OPTIMIZED_REGISTRATION_METHODS
  28 // #include "itkOptBSplineInterpolateImageFunction.h"
  29 // #else
  30
  31 #include <vector>
  32
  33 #include "itkImageLinearIteratorWithIndex.h"
  34 #include "itkVectorInterpolateImageFunction.h"
  35 #include "vnl/vnl_matrix.h"
  36
  37 #include "itkBSplineDecompositionImageFilter.h"
  38 #include "itkConceptChecking.h"
  39 #include "itkCovariantVector.h"
  40
  41
  42 namespace clitk
  43 {
  44
  45 template <
  46   class TImageType,
  47   class TCoordRep = double,
  48   class TCoefficientType = double >
  49 class ITK_EXPORT VectorBSplineInterpolateImageFunction :
  50     public itk::VectorInterpolateImageFunction<TImageType,TCoordRep>
  51 {
  52 public:
  53   /** Standard class typedefs. */
  54   typedef VectorBSplineInterpolateImageFunction       Self;
  55   typedef itk::VectorInterpolateImageFunction<TImageType,TCoordRep>  Superclass;
  56   typedef itk::SmartPointer<Self>                    Pointer;
  57   typedef itk::SmartPointer<const Self>              ConstPointer;
  58
  59   /** Run-time type information (and related methods). */
  60   itkTypeMacro(VectorBSplineInterpolateImageFunction, InterpolateImageFunction);
  61
  62
  63   /** New macro for creation of through a Smart Pointer */
  64   itkNewMacro( Self );
  65
  66   /** OutputType typedef support. */
  67   typedef typename Superclass::OutputType OutputType;
  68
  69   /** InputImageType typedef support. */
  70   typedef typename Superclass::InputImageType InputImageType;
  71
  72   /** Dimension underlying input image. */
  73   itkStaticConstMacro(ImageDimension, unsigned int,Superclass::ImageDimension);
  74
  75   /** Index typedef support. */
  76   typedef typename Superclass::IndexType IndexType;
  77
  78   /** ContinuousIndex typedef support. */
  79   typedef typename Superclass::ContinuousIndexType ContinuousIndexType;
  80
  81   /** PointType typedef support */
  82   typedef typename Superclass::PointType PointType;
  83
  84   //JV the vector dimension
  85   itkStaticConstMacro(VectorDimension, unsigned int,Superclass::Dimension);
  86
  87
  88   /** Iterator typedef support */
  89   typedef itk::ImageLinearIteratorWithIndex<TImageType> Iterator;
  90
  91   /** Internal Coefficient typedef support */
  92   typedef TCoefficientType CoefficientDataType;
  93
  94   //JV
  95   typedef itk::Vector<CoefficientDataType, itkGetStaticConstMacro(VectorDimension)> CoefficientImagePixelType;
  96   typedef  itk::Image<CoefficientImagePixelType , itkGetStaticConstMacro(ImageDimension) > CoefficientImageType;
  97
  98   /** Define filter for calculating the BSpline coefficients */
  99   //JV make vectorial
 100   typedef clitk::VectorBSplineDecompositionImageFilter<TImageType, CoefficientImageType>
 101   CoefficientFilter;
 102   typedef typename CoefficientFilter::Pointer CoefficientFilterPointer;
 103
 104   /** Evaluate the function at a ContinuousIndex position.
 105    *
 106    * Returns the B-Spline interpolated image intensity at a
 107    * specified point position. No bounds checking is done.
 108    * The point is assume to lie within the image buffer.
 109    *
 110    * ImageFunction::IsInsideBuffer() can be used to check bounds before
 111    * calling the method. */
 112   virtual OutputType EvaluateAtContinuousIndex(
 113     const ContinuousIndexType & index ) const;
 114
 115   /** Derivative typedef support */
 116   typedef itk::CovariantVector<OutputType,
 117                           itkGetStaticConstMacro(ImageDimension)
 118     > CovariantVectorType;
 119
 120   CovariantVectorType EvaluateDerivative( const PointType & point ) const
 121   {
 122     ContinuousIndexType index;
 123     this->GetInputImage()->TransformPhysicalPointToContinuousIndex( point, index );
 124     return ( this->EvaluateDerivativeAtContinuousIndex( index ) );
 125   }
 126
 127   CovariantVectorType EvaluateDerivativeAtContinuousIndex(
 128     const ContinuousIndexType & x ) const;
 129
 130
 131   /** Get/Sets the Spline Order, supports 0th - 5th order splines. The default
 132    *  is a 3rd order spline. */
 133   void SetSplineOrder(unsigned int SplineOrder);
 134   itkGetMacro(SplineOrder, int);
 135
 136
 137   /** Set the input image.  This must be set by the user. */
 138   virtual void SetInputImage(const TImageType * inputData);
 139
 140
 141   /** The UseImageDirection flag determines whether image derivatives are
 142    * computed with respect to the image grid or with respect to the physical
 143    * space. When this flag is ON the derivatives are computed with respect to
 144    * the coodinate system of physical space. The difference is whether we take
 145    * into account the image Direction or not. The flag ON will take into
 146    * account the image direction and will result in an extra matrix
 147    * multiplication compared to the amount of computation performed when the
 148    * flag is OFF.  This flag is OFF by default.*/
 149   itkSetMacro( UseImageDirection, bool );
 150   itkGetMacro( UseImageDirection, bool );
 151   itkBooleanMacro( UseImageDirection );
 152
 153
 154 protected:
 155   VectorBSplineInterpolateImageFunction();
 156   virtual ~VectorBSplineInterpolateImageFunction() {};
 157   void operator=( const Self& ); //purposely not implemented
 158   void PrintSelf(std::ostream& os, itk::Indent indent) const;
 159
 160   // These are needed by the smoothing spline routine.
 161   std::vector<CoefficientDataType>    m_Scratch;        // temp storage for processing of Coefficients
 162   typename TImageType::SizeType       m_DataLength;  // Image size
 163   unsigned int                        m_SplineOrder; // User specified spline order (3rd or cubic is the default)
 164
 165   typename CoefficientImageType::ConstPointer       m_Coefficients; // Spline coefficients
 166
 167 private:
 168   VectorBSplineInterpolateImageFunction( const Self& ); //purposely not implemented
 169   /** Determines the weights for interpolation of the value x */
 170   void SetInterpolationWeights( const ContinuousIndexType & x,
 171                                 const vnl_matrix<long> & EvaluateIndex,
 172                                 vnl_matrix<double> & weights,
 173                                 unsigned int splineOrder ) const;
 174
 175   /** Determines the weights for the derivative portion of the value x */
 176   void SetDerivativeWeights( const ContinuousIndexType & x,
 177                              const vnl_matrix<long> & EvaluateIndex,
 178                              vnl_matrix<double> & weights,
 179                              unsigned int splineOrder ) const;
 180
 181   /** Precomputation for converting the 1D index of the interpolation neighborhood
 182     * to an N-dimensional index. */
 183   void GeneratePointsToIndex(  );
 184
 185   /** Determines the indicies to use give the splines region of support */
 186   void DetermineRegionOfSupport( vnl_matrix<long> & evaluateIndex,
 187                                  const ContinuousIndexType & x,
 188                                  unsigned int splineOrder ) const;
 189
 190   /** Set the indicies in evaluateIndex at the boundaries based on mirror
 191     * boundary conditions. */
 192   void ApplyMirrorBoundaryConditions(vnl_matrix<long> & evaluateIndex,
 193                                      unsigned int splineOrder) const;
 194
 195
 196   Iterator                  m_CIterator;    // Iterator for traversing spline coefficients.
 197   unsigned long             m_MaxNumberInterpolationPoints; // number of neighborhood points used for interpolation
 198   std::vector<IndexType>    m_PointsToIndex;  // Preallocation of interpolation neighborhood indicies
 199
 200   CoefficientFilterPointer     m_CoefficientFilter;
 201
 202   // flag to take or not the image direction into account when computing the
 203   // derivatives.
 204   bool m_UseImageDirection;
 205
 206
 207 };
 208
 209 } // namespace itk
 210
 211 #ifndef ITK_MANUAL_INSTANTIATION
 212 #include "clitkVectorBSplineInterpolateImageFunction.txx"
 213 #endif
 214
 215 #endif
 216
 217 //#endif