00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017 #ifndef __itkCoxDeBoorBSplineKernelFunction_h
00018 #define __itkCoxDeBoorBSplineKernelFunction_h
00019
00020 #include "itkKernelFunction.h"
00021 #include "vnl/vnl_math.h"
00022 #include "vnl/vnl_real_polynomial.h"
00023 #include "vnl/vnl_matrix.h"
00024
00025 namespace itk
00026 {
00027
00055 template <unsigned int VSplineOrder = 3>
00056 class ITK_EXPORT CoxDeBoorBSplineKernelFunction
00057 : public KernelFunction
00058 {
00059 public:
00061 typedef CoxDeBoorBSplineKernelFunction Self;
00062 typedef KernelFunction Superclass;
00063 typedef SmartPointer<Self> Pointer;
00064 typedef SmartPointer<const Self> ConstPointer;
00065
00067 itkNewMacro(Self);
00068
00070 itkTypeMacro( CoxDeBoorBSplineKernelFunction, KernelFunction );
00071
00072 typedef double RealType;
00073 typedef vnl_vector<RealType> VectorType;
00074 typedef vnl_real_polynomial PolynomialType;
00075 typedef vnl_matrix<RealType> MatrixType;
00076
00078 void SetSplineOrder( unsigned int );
00079 itkGetMacro( SplineOrder, unsigned int );
00081
00083 inline RealType Evaluate( const RealType & u ) const
00084 {
00085 RealType absValue = vnl_math_abs( u );
00086 unsigned int which;
00087 if ( this->m_SplineOrder % 2 == 0 )
00088 {
00089 which = static_cast<unsigned int>( absValue+0.5 );
00090 }
00091 else
00092 {
00093 which = static_cast<unsigned int>( absValue );
00094 }
00095 if ( which < this->m_BSplineShapeFunctions.rows() )
00096 {
00097 return PolynomialType(
00098 this->m_BSplineShapeFunctions.get_row( which ) ).evaluate( absValue );
00099 }
00100 else
00101 {
00102 return NumericTraits<RealType>::Zero;
00103 }
00104 }
00106
00108 inline RealType EvaluateDerivative( const double & u ) const
00109 {
00110 RealType absValue = vnl_math_abs( u );
00111 unsigned int which;
00112 if ( this->m_SplineOrder % 2 == 0 )
00113 {
00114 which = static_cast<unsigned int>( absValue+0.5 );
00115 }
00116 else
00117 {
00118 which = static_cast<unsigned int>( absValue );
00119 }
00120 if( which < static_cast<unsigned int>( this->m_BSplineShapeFunctions.rows() ) )
00121 {
00122 RealType der = PolynomialType(
00123 this->m_BSplineShapeFunctions.get_row( which ) ).devaluate( absValue );
00124 if ( u < NumericTraits<RealType>::Zero )
00125 {
00126 return -der;
00127 }
00128 else
00129 {
00130 return der;
00131 }
00132 }
00133 else
00134 {
00135 return NumericTraits<RealType>::Zero;
00136 }
00137 }
00139
00144 MatrixType GetShapeFunctions();
00145
00150 MatrixType GetShapeFunctionsInZeroToOneInterval();
00151
00152 protected:
00153 CoxDeBoorBSplineKernelFunction();
00154 ~CoxDeBoorBSplineKernelFunction();
00155 void PrintSelf( std::ostream& os, Indent indent ) const;
00156
00157 private:
00158 CoxDeBoorBSplineKernelFunction( const Self& );
00159 void operator=( const Self& );
00160
00165 void GenerateBSplineShapeFunctions( unsigned int );
00166
00173 PolynomialType CoxDeBoor( unsigned short, VectorType, unsigned int, unsigned int );
00174
00175 MatrixType m_BSplineShapeFunctions;
00176 unsigned int m_SplineOrder;
00177
00178 };
00179
00180 }
00181
00182 #ifndef ITK_MANUAL_INSTANTIATION
00183 #include "itkCoxDeBoorBSplineKernelFunction.txx"
00184 #endif
00185
00186 #endif
00187