1c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath// This file is part of Eigen, a lightweight C++ template library 2c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath// for linear algebra. 3c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath// 4c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath// Copyright (C) 2009 Mark Borgerding mark a borgerding net 5c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath// 6c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath// This Source Code Form is subject to the terms of the Mozilla 7c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath// Public License v. 2.0. If a copy of the MPL was not distributed 8c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath// with this file, You can obtain one at http://mozilla.org/MPL/2.0/. 9c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 10c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamathnamespace Eigen { 11c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 12c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamathnamespace internal { 13c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 14c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath // FFTW uses non-const arguments 15c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath // so we must use ugly const_cast calls for all the args it uses 16c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath // 17c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath // This should be safe as long as 18c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath // 1. we use FFTW_ESTIMATE for all our planning 19c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath // see the FFTW docs section 4.3.2 "Planner Flags" 20c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath // 2. fftw_complex is compatible with std::complex 21c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath // This assumes std::complex<T> layout is array of size 2 with real,imag 22c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath template <typename T> 23c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath inline 24c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath T * fftw_cast(const T* p) 25c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath { 26c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath return const_cast<T*>( p); 27c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 28c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 29c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath inline 30c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath fftw_complex * fftw_cast( const std::complex<double> * p) 31c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath { 32c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath return const_cast<fftw_complex*>( reinterpret_cast<const fftw_complex*>(p) ); 33c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 34c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 35c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath inline 36c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath fftwf_complex * fftw_cast( const std::complex<float> * p) 37c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath { 38c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath return const_cast<fftwf_complex*>( reinterpret_cast<const fftwf_complex*>(p) ); 39c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 40c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 41c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath inline 42c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath fftwl_complex * fftw_cast( const std::complex<long double> * p) 43c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath { 44c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath return const_cast<fftwl_complex*>( reinterpret_cast<const fftwl_complex*>(p) ); 45c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 46c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 47c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath template <typename T> 48c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath struct fftw_plan {}; 49c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 50c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath template <> 51c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath struct fftw_plan<float> 52c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath { 53c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath typedef float scalar_type; 54c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath typedef fftwf_complex complex_type; 55c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath fftwf_plan m_plan; 56c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath fftw_plan() :m_plan(NULL) {} 57c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath ~fftw_plan() {if (m_plan) fftwf_destroy_plan(m_plan);} 58c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 59c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath inline 60c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath void fwd(complex_type * dst,complex_type * src,int nfft) { 61c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath if (m_plan==NULL) m_plan = fftwf_plan_dft_1d(nfft,src,dst, FFTW_FORWARD, FFTW_ESTIMATE|FFTW_PRESERVE_INPUT); 62c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath fftwf_execute_dft( m_plan, src,dst); 63c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 64c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath inline 65c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath void inv(complex_type * dst,complex_type * src,int nfft) { 66c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath if (m_plan==NULL) m_plan = fftwf_plan_dft_1d(nfft,src,dst, FFTW_BACKWARD , FFTW_ESTIMATE|FFTW_PRESERVE_INPUT); 67c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath fftwf_execute_dft( m_plan, src,dst); 68c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 69c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath inline 70c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath void fwd(complex_type * dst,scalar_type * src,int nfft) { 71c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath if (m_plan==NULL) m_plan = fftwf_plan_dft_r2c_1d(nfft,src,dst,FFTW_ESTIMATE|FFTW_PRESERVE_INPUT); 72c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath fftwf_execute_dft_r2c( m_plan,src,dst); 73c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 74c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath inline 75c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath void inv(scalar_type * dst,complex_type * src,int nfft) { 76c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath if (m_plan==NULL) 77c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath m_plan = fftwf_plan_dft_c2r_1d(nfft,src,dst,FFTW_ESTIMATE|FFTW_PRESERVE_INPUT); 78c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath fftwf_execute_dft_c2r( m_plan, src,dst); 79c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 80c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 81c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath inline 82c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath void fwd2( complex_type * dst,complex_type * src,int n0,int n1) { 83c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath if (m_plan==NULL) m_plan = fftwf_plan_dft_2d(n0,n1,src,dst,FFTW_FORWARD,FFTW_ESTIMATE|FFTW_PRESERVE_INPUT); 84c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath fftwf_execute_dft( m_plan, src,dst); 85c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 86c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath inline 87c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath void inv2( complex_type * dst,complex_type * src,int n0,int n1) { 88c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath if (m_plan==NULL) m_plan = fftwf_plan_dft_2d(n0,n1,src,dst,FFTW_BACKWARD,FFTW_ESTIMATE|FFTW_PRESERVE_INPUT); 89c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath fftwf_execute_dft( m_plan, src,dst); 90c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 91c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 92c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath }; 93c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath template <> 94c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath struct fftw_plan<double> 95c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath { 96c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath typedef double scalar_type; 97c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath typedef fftw_complex complex_type; 98c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath ::fftw_plan m_plan; 99c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath fftw_plan() :m_plan(NULL) {} 100c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath ~fftw_plan() {if (m_plan) fftw_destroy_plan(m_plan);} 101c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 102c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath inline 103c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath void fwd(complex_type * dst,complex_type * src,int nfft) { 104c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath if (m_plan==NULL) m_plan = fftw_plan_dft_1d(nfft,src,dst, FFTW_FORWARD, FFTW_ESTIMATE|FFTW_PRESERVE_INPUT); 105c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath fftw_execute_dft( m_plan, src,dst); 106c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 107c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath inline 108c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath void inv(complex_type * dst,complex_type * src,int nfft) { 109c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath if (m_plan==NULL) m_plan = fftw_plan_dft_1d(nfft,src,dst, FFTW_BACKWARD , FFTW_ESTIMATE|FFTW_PRESERVE_INPUT); 110c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath fftw_execute_dft( m_plan, src,dst); 111c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 112c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath inline 113c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath void fwd(complex_type * dst,scalar_type * src,int nfft) { 114c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath if (m_plan==NULL) m_plan = fftw_plan_dft_r2c_1d(nfft,src,dst,FFTW_ESTIMATE|FFTW_PRESERVE_INPUT); 115c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath fftw_execute_dft_r2c( m_plan,src,dst); 116c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 117c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath inline 118c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath void inv(scalar_type * dst,complex_type * src,int nfft) { 119c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath if (m_plan==NULL) 120c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath m_plan = fftw_plan_dft_c2r_1d(nfft,src,dst,FFTW_ESTIMATE|FFTW_PRESERVE_INPUT); 121c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath fftw_execute_dft_c2r( m_plan, src,dst); 122c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 123c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath inline 124c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath void fwd2( complex_type * dst,complex_type * src,int n0,int n1) { 125c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath if (m_plan==NULL) m_plan = fftw_plan_dft_2d(n0,n1,src,dst,FFTW_FORWARD,FFTW_ESTIMATE|FFTW_PRESERVE_INPUT); 126c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath fftw_execute_dft( m_plan, src,dst); 127c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 128c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath inline 129c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath void inv2( complex_type * dst,complex_type * src,int n0,int n1) { 130c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath if (m_plan==NULL) m_plan = fftw_plan_dft_2d(n0,n1,src,dst,FFTW_BACKWARD,FFTW_ESTIMATE|FFTW_PRESERVE_INPUT); 131c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath fftw_execute_dft( m_plan, src,dst); 132c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 133c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath }; 134c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath template <> 135c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath struct fftw_plan<long double> 136c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath { 137c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath typedef long double scalar_type; 138c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath typedef fftwl_complex complex_type; 139c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath fftwl_plan m_plan; 140c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath fftw_plan() :m_plan(NULL) {} 141c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath ~fftw_plan() {if (m_plan) fftwl_destroy_plan(m_plan);} 142c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 143c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath inline 144c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath void fwd(complex_type * dst,complex_type * src,int nfft) { 145c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath if (m_plan==NULL) m_plan = fftwl_plan_dft_1d(nfft,src,dst, FFTW_FORWARD, FFTW_ESTIMATE|FFTW_PRESERVE_INPUT); 146c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath fftwl_execute_dft( m_plan, src,dst); 147c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 148c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath inline 149c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath void inv(complex_type * dst,complex_type * src,int nfft) { 150c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath if (m_plan==NULL) m_plan = fftwl_plan_dft_1d(nfft,src,dst, FFTW_BACKWARD , FFTW_ESTIMATE|FFTW_PRESERVE_INPUT); 151c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath fftwl_execute_dft( m_plan, src,dst); 152c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 153c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath inline 154c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath void fwd(complex_type * dst,scalar_type * src,int nfft) { 155c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath if (m_plan==NULL) m_plan = fftwl_plan_dft_r2c_1d(nfft,src,dst,FFTW_ESTIMATE|FFTW_PRESERVE_INPUT); 156c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath fftwl_execute_dft_r2c( m_plan,src,dst); 157c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 158c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath inline 159c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath void inv(scalar_type * dst,complex_type * src,int nfft) { 160c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath if (m_plan==NULL) 161c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath m_plan = fftwl_plan_dft_c2r_1d(nfft,src,dst,FFTW_ESTIMATE|FFTW_PRESERVE_INPUT); 162c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath fftwl_execute_dft_c2r( m_plan, src,dst); 163c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 164c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath inline 165c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath void fwd2( complex_type * dst,complex_type * src,int n0,int n1) { 166c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath if (m_plan==NULL) m_plan = fftwl_plan_dft_2d(n0,n1,src,dst,FFTW_FORWARD,FFTW_ESTIMATE|FFTW_PRESERVE_INPUT); 167c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath fftwl_execute_dft( m_plan, src,dst); 168c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 169c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath inline 170c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath void inv2( complex_type * dst,complex_type * src,int n0,int n1) { 171c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath if (m_plan==NULL) m_plan = fftwl_plan_dft_2d(n0,n1,src,dst,FFTW_BACKWARD,FFTW_ESTIMATE|FFTW_PRESERVE_INPUT); 172c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath fftwl_execute_dft( m_plan, src,dst); 173c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 174c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath }; 175c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 176c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath template <typename _Scalar> 177c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath struct fftw_impl 178c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath { 179c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath typedef _Scalar Scalar; 180c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath typedef std::complex<Scalar> Complex; 181c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 182c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath inline 183c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath void clear() 184c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath { 185c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath m_plans.clear(); 186c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 187c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 188c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath // complex-to-complex forward FFT 189c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath inline 190c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath void fwd( Complex * dst,const Complex *src,int nfft) 191c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath { 192c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath get_plan(nfft,false,dst,src).fwd(fftw_cast(dst), fftw_cast(src),nfft ); 193c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 194c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 195c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath // real-to-complex forward FFT 196c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath inline 197c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath void fwd( Complex * dst,const Scalar * src,int nfft) 198c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath { 199c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath get_plan(nfft,false,dst,src).fwd(fftw_cast(dst), fftw_cast(src) ,nfft); 200c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 201c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 202c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath // 2-d complex-to-complex 203c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath inline 204c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath void fwd2(Complex * dst, const Complex * src, int n0,int n1) 205c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath { 206c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath get_plan(n0,n1,false,dst,src).fwd2(fftw_cast(dst), fftw_cast(src) ,n0,n1); 207c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 208c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 209c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath // inverse complex-to-complex 210c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath inline 211c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath void inv(Complex * dst,const Complex *src,int nfft) 212c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath { 213c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath get_plan(nfft,true,dst,src).inv(fftw_cast(dst), fftw_cast(src),nfft ); 214c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 215c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 216c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath // half-complex to scalar 217c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath inline 218c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath void inv( Scalar * dst,const Complex * src,int nfft) 219c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath { 220c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath get_plan(nfft,true,dst,src).inv(fftw_cast(dst), fftw_cast(src),nfft ); 221c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 222c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 223c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath // 2-d complex-to-complex 224c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath inline 225c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath void inv2(Complex * dst, const Complex * src, int n0,int n1) 226c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath { 227c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath get_plan(n0,n1,true,dst,src).inv2(fftw_cast(dst), fftw_cast(src) ,n0,n1); 228c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 229c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 230c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 231c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath protected: 232c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath typedef fftw_plan<Scalar> PlanData; 233c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 234c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath typedef std::map<int64_t,PlanData> PlanMap; 235c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 236c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath PlanMap m_plans; 237c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 238c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath inline 239c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath PlanData & get_plan(int nfft,bool inverse,void * dst,const void * src) 240c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath { 241c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath bool inplace = (dst==src); 242c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath bool aligned = ( (reinterpret_cast<size_t>(src)&15) | (reinterpret_cast<size_t>(dst)&15) ) == 0; 243c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath int64_t key = ( (nfft<<3 ) | (inverse<<2) | (inplace<<1) | aligned ) << 1; 244c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath return m_plans[key]; 245c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 246c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 247c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath inline 248c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath PlanData & get_plan(int n0,int n1,bool inverse,void * dst,const void * src) 249c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath { 250c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath bool inplace = (dst==src); 251c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath bool aligned = ( (reinterpret_cast<size_t>(src)&15) | (reinterpret_cast<size_t>(dst)&15) ) == 0; 252c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath int64_t key = ( ( (((int64_t)n0) << 30)|(n1<<3 ) | (inverse<<2) | (inplace<<1) | aligned ) << 1 ) + 1; 253c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath return m_plans[key]; 254c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 255c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath }; 256c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 257c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath} // end namespace internal 258c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 259c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath} // end namespace Eigen 260c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 261c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath/* vim: set filetype=cpp et sw=2 ts=2 ai: */ 262