用C语言实现FFT算法
用C语言实现FFT算法
/*****************fft programe*********************/#include "typedef.h" #include "math.h"
struct compx EE(struct compx b1,struct compx b2){ struct compx b3 ; b3.real=b1.real*b2.real-b1.imag*b2.imag ; b3.imag=b1.real*b2.imag+b1.imag*b2.real ; return(b3);}
void FFT(struct compx*xin,int N){ int f,m,nv2,nm1,i,k,j=1,l ; /*int f,m,nv2,nm1,i,k,j=N/2,l;*/ struct compx v,w,t ; nv2=N/2 ; f=N ; for(m=1;(f=f/2)!=1;m++) { ; } nm1=N-1 ; /*变址运算*/ for(i=1;i { if(i { t=xin[j]; xin[j]=xin[i]; xin[i]=t ; } k=nv2 ; while(k { j=j-k ; k=k/2 ; } j=j+k ; } { int le,lei,ip ; float pi ; for(l=1;l { le=pow(2,l); // 这里用的是L而不是1 !!!! lei=le/2 ; pi=3.14159 ; v.real=1.0 ; v.imag=0.0 ; w.real=cos(pi/lei); w.imag=-sin(pi/lei); for(j=1;j { /*double p=pow(2,m-l)*j; double ps=2*pi/N*p; w.real=cos(ps); w.imag=-sin(ps);*/ for(i=j;i { /* w.real=cos(ps); w.imag=-sin(ps);*/ ip=i+lei ; t=EE(xin[ip],v); xin[ip].real=xin[i].real-t.real ; xin[ip].imag=xin[i].imag-t.imag ; xin[i].real=xin[i].real+t.real ; xin[i].imag=xin[i].imag+t.imag ; } v=EE(v,w); } } } return ;}
/*****************main programe********************/
#include#include#include#include "typedef.h"
float result[257];struct compx s[257];int Num=256 ;const float pp=3.14159 ;
main(){ int i=1 ; for(;i { s[i].real=sin(pp*i/32); s[i].imag=0 ; } FFT(s,Num); for(i=1;i { result[i]=sqrt(pow(s[i].real,2)+pow(s[i].imag,2)); }}