fft2.txt

本程序是关于二维傅里叶变换计算相关的

#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);
    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));
}


}