fft.zip_512 fft_fft_fft编程；标准fft算法

#include "math.h" #include "stdio.h" float max=0,sum=0; struct compx { float real; float imag; } compx ; 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) //*xin为需要FFT变换的数组，N为转换的点数 { int f,m,LH,nm,i,k,j,L; double p , ps ; long int le,B,ip; float pi; struct compx w,t; LH=N>>1; f=N; for(m=1;(f=f>>1)!=1;m++){;} nm=N-2; j=N>>1; for(i=1;i<=nm;i++) // 变址运算 { if(i<j){t=xin[j];xin[j]=xin[i];xin[i]=t;} k=LH; while(j>=k){j=j-k;k=k/2;} j=j+k; } for(L=1;L<=m;L++) //实现蝶形运算 { le=(long int)pow(2,L); B=le>>1; pi=3.14159; for(j=0;j<=B-1;j++) { p=pow(2,m-L)*j; ps=2*pi/N*p; w.real=cos(ps); w.imag=-sin(ps); for(i=j;i<=N-1;i=i+le) { ip=i+B; t=EE(xin[ip],w); 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; } } } } /******************************************************************** *功能： 失真度的计算 * *****************************************************************/ float THDcount(unsigned int *result) { unsigned int i=0,n=2,m,t=0; float THD=0; max=result[1]; for(i=2;i<256;i++) { if(result[i]>max) { max=result[i]; t=i; } } m=t*n; while(m<256) { sum=result[m]*result[m]+sum; n++; m=t*n; } sum=sum/(max*max); THD=100*sqrt(sum); return THD; }