DSP实现FFT的代码

《DSP实现FFT的代码》由会员上传分享，免费在线阅读，更多相关内容在行业资料-天天文库。

1、DSP实现的方法//#include"DSP281x_Device.h"//DSP281xHeaderfileIncludeFile//#include"DSP281x_Examples.h"//DSP281xExamplesIncludeFile#include#include#include#definepiatan(1.0)*4//计算PI的值typedefstruct//定义表示复数的结构体{floatreal;floatimag;}FUSHU;FUSHUx

2、[1024];intw[1024];FUSHUFUCHENG(FUSHUb1,FUSHUb2)//复数相乘函数{FUSHUb3;b3.real=b1.real*b2.real-b1.imag*b2.imag;b3.imag=b1.real*b2.imag+b1.imag*b2.real;return(b3);}voidFFT(FUSHU*xin,intN)//FFT运算{intm,LH,I,k,J,M,K;floatp,ps;intB,N1;FUSHUw,T;M=10;//对应N值为//下面是倒序的程序LH=N/2;J

3、=LH;N1=N-2;/*变址运算*/for(I=1;I<=N1;I++){if(I=K){J=J-K;K=K/2;}J=J+K;}//下面是DIT-FFT运算程序for(m=1;m<=M;m++){B=pow(2.0,m-1);for(J=0;J<=B-1;J++){p=pow(2.0,M-m)*J;ps=2*pi/N*p;w.real=cos(ps);w.imag=-sin(ps);for(k=J;k<=N-1;k=

4、k+pow(2.0,m)){T=FUCHENG(xin[k+B],w);xin[k+B].real=xin[k].real-T.real;xin[k+B].imag=xin[k].imag-T.imag;xin[k].real=xin[k].real+T.real;xin[k].imag=xin[k].imag+T.imag;}}}}voidIFFT(FUSHU*xin,intN)//IFFT运算{intm,LH,I,k,J,M,K;floatp,ps;intB,N1;FUSHUw,T;M=10;//对应N值为//下面

5、是倒序的程序LH=N/2;J=LH;N1=N-2;/*变址运算*/for(I=1;I<=N1;I++){if(I=K){J=J-K;K=K/2;}J=J+K;}for(m=1;m<=M;m++){B=pow(2.0,m-1);for(J=0;J<=B-1;J++){p=pow(2.0,M-m)*J;ps=2*pi/N*p;w.real=cos(ps);w.imag=sin(ps);//与FFT符号相反for(k=J;k<

6、=N-1;k=k+pow(2.0,m)){T=FUCHENG(xin[k+B],w);xin[k+B].real=xin[k].real-T.real;xin[k+B].imag=xin[k].imag-T.imag;xin[k].real=xin[k].real+T.real;xin[k].imag=xin[k].imag+T.imag;}}}for(k=0;k

7、intmain(){intn;for(n=0;n<1024;n++)//输入波形x(t)=cos(2*pi*100*t){x[n].real=cos(2*pi*100*n/1024);//采样频率为kHzx[n].imag=0;}FFT(x,1024);//进行FFT运算for(n=0;n<1024;n++){w[n]=sqrt(x[n].real*x[n].real+x[n].imag*x[n].imag);//幅度谱}IFFT(x,1024);//进行IFFT运算for(n=0;n<1024;n++)//比较逆变换

8、之后的数据与原来的数据差{if(fabs(cos(2*pi*100*n/1024)-x[n].real)>0.001)break;}if(n<1024)puts("FFT不正确");if(n==1024)puts("FFT正确");return0;}