diff options
author | Mark Borgerding <mark@borgerding.net> | 2010-01-22 00:35:03 -0500 |
---|---|---|
committer | Mark Borgerding <mark@borgerding.net> | 2010-01-22 00:35:03 -0500 |
commit | cd7912313dc2477283de767029462d7d0e6ee8ab (patch) | |
tree | bfd5ec7f5a448e3220e592d120f6efd44034bada /unsupported/Eigen/src/FFT | |
parent | a30d42354f06b86e35838ff9e8c14b524bf1c8aa (diff) |
changed FFT function vector and Matrix args to pointer as Benoit suggested
implemented 2D Complex FFT for FFTW impl
Diffstat (limited to 'unsupported/Eigen/src/FFT')
-rw-r--r-- | unsupported/Eigen/src/FFT/ei_fftw_impl.h | 61 |
1 files changed, 59 insertions, 2 deletions
diff --git a/unsupported/Eigen/src/FFT/ei_fftw_impl.h b/unsupported/Eigen/src/FFT/ei_fftw_impl.h index a66b7398c..411ff7425 100644 --- a/unsupported/Eigen/src/FFT/ei_fftw_impl.h +++ b/unsupported/Eigen/src/FFT/ei_fftw_impl.h @@ -90,6 +90,18 @@ m_plan = fftwf_plan_dft_c2r_1d(nfft,src,dst,FFTW_ESTIMATE); fftwf_execute_dft_c2r( m_plan, src,dst); } + + inline + void fwd2( complex_type * dst,complex_type * src,int nrows,int ncols) { + if (m_plan==NULL) m_plan = fftwf_plan_dft_2d(ncols,nrows,src,dst,FFTW_FORWARD,FFTW_ESTIMATE); + fftwf_execute_dft( m_plan, src,dst); + } + inline + void inv2( complex_type * dst,complex_type * src,int nrows,int ncols) { + if (m_plan==NULL) m_plan = fftwf_plan_dft_2d(ncols,nrows,src,dst,FFTW_BACKWARD,FFTW_ESTIMATE); + fftwf_execute_dft( m_plan, src,dst); + } + }; template <> struct ei_fftw_plan<double> @@ -121,6 +133,16 @@ m_plan = fftw_plan_dft_c2r_1d(nfft,src,dst,FFTW_ESTIMATE); fftw_execute_dft_c2r( m_plan, src,dst); } + inline + void fwd2( complex_type * dst,complex_type * src,int nrows,int ncols) { + if (m_plan==NULL) m_plan = fftw_plan_dft_2d(ncols,nrows,src,dst,FFTW_FORWARD,FFTW_ESTIMATE); + fftw_execute_dft( m_plan, src,dst); + } + inline + void inv2( complex_type * dst,complex_type * src,int nrows,int ncols) { + if (m_plan==NULL) m_plan = fftw_plan_dft_2d(ncols,nrows,src,dst,FFTW_BACKWARD,FFTW_ESTIMATE); + fftw_execute_dft( m_plan, src,dst); + } }; template <> struct ei_fftw_plan<long double> @@ -152,6 +174,16 @@ m_plan = fftwl_plan_dft_c2r_1d(nfft,src,dst,FFTW_ESTIMATE); fftwl_execute_dft_c2r( m_plan, src,dst); } + inline + void fwd2( complex_type * dst,complex_type * src,int nrows,int ncols) { + if (m_plan==NULL) m_plan = fftwl_plan_dft_2d(ncols,nrows,src,dst,FFTW_FORWARD,FFTW_ESTIMATE); + fftwl_execute_dft( m_plan, src,dst); + } + inline + void inv2( complex_type * dst,complex_type * src,int nrows,int ncols) { + if (m_plan==NULL) m_plan = fftwl_plan_dft_2d(ncols,nrows,src,dst,FFTW_BACKWARD,FFTW_ESTIMATE); + fftwl_execute_dft( m_plan, src,dst); + } }; template <typename _Scalar> @@ -180,6 +212,13 @@ get_plan(nfft,false,dst,src).fwd(ei_fftw_cast(dst), ei_fftw_cast(src) ,nfft); } + // 2-d complex-to-complex + inline + void fwd2(Complex * dst, const Complex * src, int nrows,int ncols) + { + get_plan(nrows,ncols,false,dst,src).fwd2(ei_fftw_cast(dst), ei_fftw_cast(src) ,nrows,ncols); + } + // inverse complex-to-complex inline void inv(Complex * dst,const Complex *src,int nfft) @@ -194,9 +233,18 @@ get_plan(nfft,true,dst,src).inv(ei_fftw_cast(dst), ei_fftw_cast(src),nfft ); } + // 2-d complex-to-complex + inline + void inv2(Complex * dst, const Complex * src, int nrows,int ncols) + { + get_plan(nrows,ncols,true,dst,src).inv2(ei_fftw_cast(dst), ei_fftw_cast(src) ,nrows,ncols); + } + + protected: typedef ei_fftw_plan<Scalar> PlanData; - typedef std::map<int,PlanData> PlanMap; + + typedef std::map<int64_t,PlanData> PlanMap; PlanMap m_plans; @@ -205,7 +253,16 @@ { bool inplace = (dst==src); bool aligned = ( (reinterpret_cast<size_t>(src)&15) | (reinterpret_cast<size_t>(dst)&15) ) == 0; - int key = (nfft<<3 ) | (inverse<<2) | (inplace<<1) | aligned; + int64_t key = ( (nfft<<3 ) | (inverse<<2) | (inplace<<1) | aligned ) << 1; + return m_plans[key]; + } + + inline + PlanData & get_plan(int nrows,int ncols,bool inverse,void * dst,const void * src) + { + bool inplace = (dst==src); + bool aligned = ( (reinterpret_cast<size_t>(src)&15) | (reinterpret_cast<size_t>(dst)&15) ) == 0; + int64_t key = ( ( (((int64_t)ncols) << 30)|(nrows<<3 ) | (inverse<<2) | (inplace<<1) | aligned ) << 1 ) + 1; return m_plans[key]; } }; |