快速傅里叶变换的基2FFT算法的C实现

1、快速傅里叶变换的基2fft算法的c 实现快速傅里叶变换的基2fft算法的c+实现2011-01-19 05：26快速傅里叶变换的基本原理由于公式不好显示请读者参考其它文章或书籍，本文重点给出了时域抽取法fft的c+实现。下面是算法的流程图倒序的流程图c+实现代码：1.fft.h#pragma once#ifndef fft_h#define fft_h#include cstring#include cmath using namespace std；#define pi 3.141592 class fftprivate：void changeorder(double*xr,double*x

2、i,int n)；public：fft(void)；void fft_1d(double*ctxr,double*ctxi,double*cfxr,double*cfxi,int len)；void rfft_1d(double*ctxr,double*cfxr,double*cfxi,int len)；void ifft_1d(double*cfxr,double*cfxi,double*ctxr,double*ctxi,int len)；void rifft_1d(double*cfxr,double*cfxi,double*ctxr,int len)；fft(void)；#endif 2

3、.fft.cpp#include.fft.hfft：fft()/倒序实现/xr实部，xi虚部，n为2的幂/void fft：changeorder(double*xr,double*xi,int n)double temp；int k；int lh=n/2；int j=lh；int n1=n-2；for(int i=1；i=n1；i+)if(i j)temp=xri；xri=xrj；xrj=temp；temp=xii；xii=xij；xij=temp；k=lh；while(j=k)j=j-k；k=(int)(k/2+0.5)；j=j+k；/复数fft/ctxr和ctxi的长度为len/cfxr

4、和cfxi的长度为2的幂/void fft：fft_1d(double*ctxr,double*ctxi,double*cfxr,double*cfxi,int len)int m=ceil(log(double)len)/log(2.0)；int l,b,j,p,k；double rkb,ikb；int n=1 m；double*rcos=new doublen/2；double*isin=new doublen/2；for(l=0；l n/2；l+)rcosl=cos(l*pi*2/n)；isinl=sin(l*pi*2/n)；memcpy(cfxr,ctxr,sizeof(double)

5、*len)；memcpy(cfxi,ctxi,sizeof(double)*len)；for(l=len；l n；l+)cfxrl=0；cfxil=0；changeorder(cfxr,cfxi,n)；/倒序for(l=1；l=m；l+)b=(int)(pow(2,l-1)+0.5)；for(j=0；j b；j+)p=j*(int)(pow(2,m-l)+0.5)；for(k=j；k n；k+=(int)(pow(2,l)+0.5)rkb=cfxrk+b*rcos+cfxik+b*isin；ikb=cfxik+b*rcos-cfxrk+b*isin；cfxrk+b=cfxrk-rkb；cfxi

6、k+b=cfxik-ikb；cfxrk=cfxrk+rkb；cfxik=cfxik+ikb；delete rcos；delete isin；/实数fft/ctxr的长度为len/cfxr和cfxi的长度为2的幂/void fft：rfft_1d(double*ctxr,double*cfxr,double*cfxi,int len)int m=ceil(log(double)len)/log(2.0)；int n=1 m；double*rcos=new doublen/2；double*isin=new doublen/2；for(int l=0；l n/2；l+)rcosl=cos(l*pi

7、*2/n)；isinl=sin(l*pi*2/n)；double*txr=new double(len+1)/2；double*txi=new double(len+1)/2；for(int i=0；i len/2；i+)txri=ctxri*2；txii=ctxri*2+1；if(len%2=1)txr(len-1)/2=ctxrlen-1；txi(len-1)/2=0；fft_1d(txr,txi,cfxr,cfxi,(len+1)/2)；double*x1r=new doublen/2；double*x1i=new doublen/2；double*x2r=new doublen/2；d

8、ouble*x2i=new doublen/2；x1r0=cfxr0；x1i0=0；x2r0=cfxi0；x2i0=0；for(int k=1；k n/2；k+)x1rk=(cfxrk+cfxrn/2-k)/2.0；x1ik=(cfxik-cfxin/2-k)/2.0；x2rk=(cfxik+cfxin/2-k)/2.0；x2ik=(-cfxrk+cfxrn/2-k)/2.0；double rkb,ikb；for(i=0；i n/2；i+)rkb=x2ri*rcosi+x2ii*isini；ikb=x2ii*rcosi-x2ri*isini；cfxri+n/2=x1ri-rkb；cfxii+n

9、/2=x1ii-ikb；cfxri=x1ri+rkb；cfxii=x1ii+ikb；delete txr；delete txi；delete x1r；delete x1i；delete x2r；delete x2i；delete rcos；delete isin；/复数ifft/cfxr和cfxi的长度为n(2的幂)/ctxr和ctxi的长度为len/void fft：ifft_1d(double*cfxr,double*cfxi,double*ctxr,double*ctxi,int len)int m=ceil(log(double)len)/log(2.0)；int n=1 m；doub

10、le*txr=new doublen；double*txi=new doublen；for(int i=0；i n；i+)cfxii=-cfxii；fft_1d(cfxr,cfxi,txr,txi,n)；for(i=0；i len；i+)ctxri=txri/n；ctxii=-txii/n；delete txr；delete txi；/实数ifft/cfxr和cfxi的长度为n(2的幂)/ctxr的长度为len/void fft：rifft_1d(double*cfxr,double*cfxi,double*ctxr,int len)int m=ceil(log(double)len)/log

11、(2.0)；int n=1 m；double*txr=new doublen；double*txi=new doublen；for(int i=0；i n；i+)cfxii=-cfxii；fft_1d(cfxr,cfxi,txr,txi,n)；for(int i=0；i len；i+)ctxri=txri/n；delete txr；delete txi；fft：fft(void)3.test.cpp#include iostream#includefft.husing namespace std；void main()double xr10=1,2,3,4,5,6,7,8,9,10；/实部double xi10=0,0,0,0,0,0,0,0,0,0；/虚部double*cxr,*cxi；cxr=new double16；cxi=new double16；fft f；f.rfft_1d(xr,cxr,cxi,10)；for(int i=0；i 16；i+)cout cxri+jc