1 /* =======================================================================
2 * @author: Leonardo Florez-Valencia
3 * @email: florez-l@javeriana.edu.co
4 * =======================================================================
6 #ifndef __FourierSeries__h__
7 #define __FourierSeries__h__
14 #include <itkMatrix.h>
16 #include <itkVector.h>
19 * Planar series defined by its harmonic ellipses.
21 * The tuple [z(l),z(-l)], l != 0 defines an ellipse. The z(0)
22 * harmonic defines the series's center.
25 template< class S, unsigned int D >
27 : public std::deque< std::complex< S > >
30 typedef FourierSeries Self;
32 itkStaticConstMacro( Dimension, unsigned int, D );
34 typedef itk::Matrix< S, D, D > TMatrix;
35 typedef itk::Point< S, D > TPoint;
36 typedef typename TPoint::VectorType TVector;
37 typedef std::complex< S > TComplex;
38 typedef std::deque< TComplex > Superclass;
41 FourierSeries( unsigned int q = 1 );
42 template< class S2, unsigned int D2 >
43 FourierSeries( const FourierSeries< S2, D2 >& o );
44 template< class _TIt >
45 FourierSeries( const _TIt& b, const _TIt& e, unsigned int q );
46 virtual ~FourierSeries( );
48 template< class S2, unsigned int D2 >
49 Self& operator=( const FourierSeries< S2, D2 >& o );
51 bool operator==( const Self& o ) const { return( false ); }
52 bool operator!=( const Self& o ) const { return( true ); }
54 const unsigned int& GetRealAxis( ) const;
55 const unsigned int& GetImagAxis( ) const;
56 void SetRealAxis( const unsigned int& a );
57 void SetImagAxis( const unsigned int& a );
60 * $A=\pi\sum_{l=-q}^{q}l|z_{l}|^{2}$
65 * $\kappa(\omega)=\frac{|c'(\omega)\times c''(\omega)|}{|c'(\omega)|^{3}}$
67 S GetCurvature( const S& w ) const;
70 * Sets the order of the ellipse. This order is defined as the
71 * counterclock/clock wiseness of the series. Calculations are
72 * done by computing the series's area (\see Area) and comparing
75 * \param cw counter-clock
77 void SetOrdering( bool cw );
78 void InvertOrdering( );
79 void SetOrderingToClockWise( );
80 void SetOrderingToCounterClockWise( );
83 * Series manipulators.
85 int GetNumberOfHarmonics( ) const;
86 void SetNumberOfHarmonics( const int& q );
87 TComplex& operator[]( const int& l );
88 const TComplex& operator[]( const int& l ) const;
91 * $c(\omega)=\sum_{l=-q}^{q}z_{l}e^{jl\omega}$
93 * where $0\le\omega<2\pi$ is the angular curvilinear parameter, q is
94 * the number of harmonics defining the series, $z_{l}\in\mathbb{C}$
95 * is the l-th harmonic and $j$ is the imaginary unity.
97 TPoint operator( )( const S& w ) const;
100 * $\frac{d^{(n)}c(\omega)}{d\omega^{(n)}}=\sum_{l=-q}^{q}z_{l}e^{jl\omega}(jl)^{n}$
102 * where $n$ is the derivative number, $0\le\omega<2\pi$ is the
103 * angular curvilinear parameter, q is the number of harmonics
104 * defining the series, $z_{l}\in\mathbb{C}$ is the l-th harmonic
105 * and $j$ is the imaginary unity.
107 TVector operator( )( const S& w, const unsigned int& n ) const;
110 * Coefficient-wise addition operators
112 Self operator+( const Self& o ) const;
113 Self& operator+=( const Self& o );
114 Self operator-( const Self& o ) const;
115 Self& operator-=( const Self& o );
120 Self operator+( const TComplex& translation ) const;
121 Self& operator+=( const TComplex& translation );
122 Self operator-( const TComplex& translation ) const;
123 Self& operator-=( const TComplex& translation );
124 Self operator+( const TVector& translation ) const;
125 Self& operator+=( const TVector& translation );
126 Self operator-( const TVector& translation ) const;
127 Self& operator-=( const TVector& translation );
128 TPoint GetCenter( ) const;
131 * Scale and/or rotate
133 Self operator*( const TComplex& scale_or_rotation ) const;
134 Self& operator*=( const TComplex& scale_or_rotation );
135 Self operator/( const S& scale ) const;
136 Self& operator/=( const S& scale );
141 Self operator&( const S& phase ) const;
142 Self& operator&=( const S& phase );
143 S GetPhase( const S& w = S( 0 ) ) const;
148 void Sample( std::vector< TPoint >& p, const unsigned long& s ) const;
153 void GetEllipse( int l, S& a, S& b, S& t, S& p ) const;
154 void SetEllipse( int l, S& a, S& b, S& t, S& p );
156 template< class _TOtherScalar >
157 void SetEllipses( const std::vector< _TOtherScalar >& ellipses );
160 // Main method to compute all of series's values.
161 TComplex _Z( const S& w, const unsigned int& n ) const;
164 void _DFT( const std::vector< TComplex >& p, unsigned int q );
167 unsigned int m_RealAxis;
168 unsigned int m_ImagAxis;
171 friend std::ostream& operator<<( std::ostream& out, const Self& fs )
173 int q = fs.GetNumberOfHarmonics( );
175 out << "{" << l << ":" << fs[ l ] << "}";
176 for( ++l; l <= q; l++ )
177 out << ", {" << l << ":" << fs[ l ] << "}";
182 #ifndef ITK_MANUAL_INSTANTIATION
183 #include "FourierSeries.hxx"
184 #endif // ITK_MANUAL_INSTANTIATION
185 #endif // __FourierSeries__h__