20127548 

 没调试过，大虾的，参考一下吧
/*****************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    <=nm1;i++)
{
 if(i    <j){t=xin[j];xin[j]=xin[i];xin[i]=t;}
 k=nv2;
 while(k    <j){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********************/

#include     <math.h>    
#include     <stdio.h>    
#include     <stdlib.h>    
#include "typedef.h"
   

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

main()
{

int i=1;
for(;i    <0x101;i++)
{ 
s[i].real=sin(pp*i/32);
s[i].imag=0;
}
 
FFT(s,Num);

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

}

发表时间:2008年10月15日15:19:59