Go to the documentation of this file.00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017 #ifndef __itkMultivariateLegendrePolynomial_h
00018 #define __itkMultivariateLegendrePolynomial_h
00019
00020 #include "itkIndent.h"
00021 #include <vector>
00022 #include "itkArray.h"
00023
00024 namespace itk {
00025
00071 class ITK_EXPORT MultivariateLegendrePolynomial
00072 {
00073 public:
00074 typedef MultivariateLegendrePolynomial Self;
00075
00076 typedef std::vector< double > DoubleArrayType;
00077 typedef std::vector< unsigned long > ULongArrayType;
00078 typedef std::vector< long > LongArrayType;
00079
00081 typedef DoubleArrayType CoefficientArrayType;
00082
00085 typedef Array< double > ParametersType;
00086
00088 typedef ULongArrayType DomainSizeType;
00089 typedef LongArrayType IndexType;
00090
00092 MultivariateLegendrePolynomial( unsigned int dimension,
00093 unsigned int degree,
00094 const DomainSizeType & domainSize );
00095
00097 virtual ~MultivariateLegendrePolynomial();
00098
00100 unsigned int GetDimension(void) const
00101 { return m_Dimension; }
00102
00104 unsigned int GetDegree(void) const
00105 { return m_Degree; }
00106
00113 unsigned int GetNumberOfCoefficients(void) const
00114 { return m_NumberOfCoefficients; }
00115
00117 const DomainSizeType & GetDomainSize( void ) const
00118 { return m_DomainSize; }
00119
00121 class CoefficientVectorSizeMismatch
00122 {
00123 public:
00124 CoefficientVectorSizeMismatch(int given, int required)
00125 {
00126 m_Required = required;
00127 m_Given = given;
00128 }
00129
00130 int m_Required;
00131 int m_Given;
00132 };
00133
00138 void SetCoefficients(const CoefficientArrayType& coef)
00139 throw (CoefficientVectorSizeMismatch);
00140
00141 void SetCoefficients(const ParametersType& coef)
00142 throw (CoefficientVectorSizeMismatch);
00143
00145 const CoefficientArrayType& GetCoefficients(void) const;
00146
00149 double Evaluate(IndexType& index)
00150 {
00151 if (m_Dimension == 2)
00152 {
00153 if (index[1] != m_PrevY)
00154 {
00155
00156 double norm_y = m_NormFactor[1] *
00157 static_cast<double>( index[1] - 1 );
00158 this->CalculateXCoef(norm_y, m_CoefficientArray);
00159 m_PrevY = index[1];
00160 }
00162
00163
00164 double norm_x = m_NormFactor[0] *
00165 static_cast<double>( index[0] - 1 );
00166
00167 return LegendreSum(norm_x, m_Degree, m_CachedXCoef);
00168 }
00169 else if (m_Dimension == 3)
00170 {
00171 if (index[2] != m_PrevZ )
00172 {
00173
00174 double norm_z = m_NormFactor[2] *
00175 static_cast<double>( index[2] - 1 );
00176 this->CalculateYCoef(norm_z, m_CoefficientArray);
00177 m_PrevZ = index[2];
00178 }
00179
00180 if (index[1] != m_PrevY)
00181 {
00182
00183 double norm_y = m_NormFactor[1] *
00184 static_cast<double>( index[1] - 1 );
00185 this->CalculateXCoef(norm_y, m_CachedYCoef);
00186 m_PrevY = index[1];
00187 }
00188
00189
00190 double norm_x = m_NormFactor[0] *
00191 static_cast<double>( index[0] - 1 );
00192 return this->LegendreSum(norm_x, m_Degree, m_CachedXCoef);
00193 }
00194 return 0;
00195 }
00196
00198 unsigned int GetNumberOfCoefficients();
00199
00201 unsigned int GetNumberOfCoefficients(unsigned int dimension, unsigned int degree);
00202
00208 class SimpleForwardIterator
00209 {
00210 public:
00211 SimpleForwardIterator (MultivariateLegendrePolynomial* polynomial)
00212 {
00213 m_MultivariateLegendrePolynomial = polynomial;
00214 m_Dimension = m_MultivariateLegendrePolynomial->GetDimension();
00215 m_DomainSize = m_MultivariateLegendrePolynomial->GetDomainSize();
00216 m_Index.resize(m_Dimension);
00217 std::fill(m_Index.begin(), m_Index.end(), 0);
00218 }
00219
00220 void Begin( void )
00221 {
00222 m_IsAtEnd = false;
00223 for (unsigned int dim = 0; dim < m_Dimension; dim++)
00224 {
00225 m_Index[dim] = 0;
00226 }
00227 }
00228
00229 bool IsAtEnd()
00230 { return m_IsAtEnd; }
00231
00232 SimpleForwardIterator& operator++()
00233 {
00234 for (unsigned int dim = 0; dim < m_Dimension; dim++)
00235 {
00236 if (m_Index[dim] < static_cast<int>(m_DomainSize[dim] - 1))
00237 {
00238 m_Index[dim] += 1;
00239 return *this;
00240 }
00241 else
00242 {
00243 if (dim == m_Dimension - 1 )
00244 {
00245 m_IsAtEnd = true;
00246 break;
00247 }
00248 else
00249 {
00250 m_Index[dim] = 0;
00251 }
00252 }
00253 }
00254 return *this;
00255 }
00256
00257 double Get()
00258 { return m_MultivariateLegendrePolynomial->Evaluate(m_Index); }
00259
00260 private:
00261 MultivariateLegendrePolynomial* m_MultivariateLegendrePolynomial;
00262 unsigned int m_Dimension;
00263 DomainSizeType m_DomainSize;
00264 IndexType m_Index;
00265 bool m_IsAtEnd;
00266 };
00267
00268 void Print(std::ostream& os);
00269
00270 protected:
00271 void PrintSelf(std::ostream& os, Indent indent) const;
00272 double LegendreSum(const double x, int n,
00273 const CoefficientArrayType& coef,
00274 int offset = 0);
00275 void CalculateXCoef(double norm_y, const CoefficientArrayType& coef);
00276 void CalculateYCoef(double norm_z, const CoefficientArrayType& coef);
00277
00278 private:
00279 DomainSizeType m_DomainSize;
00280 unsigned int m_Dimension;
00281 unsigned int m_Degree;
00282 unsigned int m_NumberOfCoefficients;
00283 bool m_MultiplicativeBias;
00284
00285 CoefficientArrayType m_CoefficientArray;
00286 CoefficientArrayType m_CachedXCoef;
00287 CoefficientArrayType m_CachedYCoef;
00288 CoefficientArrayType m_CachedZCoef;
00289
00290 DoubleArrayType m_NormFactor;
00291 long m_PrevY;
00292 long m_PrevZ;
00293 };
00294
00295 std::ostream& operator<< (std::ostream& os,
00296 MultivariateLegendrePolynomial& poly);
00297 }
00298 #endif
00299