用C語言實現FFT算法

/*****************fft programe*********************/
#include typedef.h
#include math.h

本文引用地址：http://www.104case.com/article/150637.htm

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=nm1;i++)
{
if(ij)
{
t=xin[j];
xin[j]=xin[i];
xin[i]=t ;
}
k=nv2 ;
while(kj)
{
j=j-k ;
k=k/2 ;
}
j=j+k ;
}

{
int le,lei,ip ;
float pi ;
for(l=1;l=m;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=lei;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=N;i=i+le)
{
/* 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********************/

#includemath.h>
#includestdio.h>
#includestdlib.h>
#include typedef.h

float result[257];
struct compx s[257];
int Num=256 ;
const float pp=3.14159 ;

main()
{
int i=1 ;
for(;i0x101;i++)
{
s[i].real=sin(pp*i/32);
s[i].imag=0 ;
}

FFT(s,Num);

for(i=1;i0x101;i++)
{
result[i]=sqrt(pow(s[i].real,2)+pow(s[i].imag,2));
}
}

c語言相關文章:c語言教程