《快速傅里叶变换(FFT)原理及源程序》由会员分享，可在线阅读，更多相关《快速傅里叶变换(FFT)原理及源程序（9页珍藏版）》请在装配图网上搜索。
1、《测试信号分析及处理》课程作业
快速傅里叶变换
1、 程序设计思路
快速傅里叶变换的目的是减少运算量，其用到的方法是分级进行运算。全部计算分解为级，其中；在输入序列中是按码位倒序排列的，输出序列是按顺序排列；每级包含个蝶形单元，第级有个群，每个群有个蝶形单元； 每个蝶形单元都包含乘和系数的运算，每个蝶形单元数据的间隔为，i为第i级； 同一级中各个群的系数分布规律完全相同。
将输入序列按码位倒序排列时，用到的是倒序算法——雷德算法。 自然序排列的二进制数，其下面一个数总比上面的数大1，而倒序二进制数的下面一个数是上面一个数在最高位加1并由高位向低位仅为而得到的。
若已知某数的倒
2、序数是，求下一个倒序数，应先判断的最高位是否为0，与进行比较即可得到结果。如果，说明最高位为0，应把其变成1，即，这样就得到倒序数了。如果，即的最高位为1，将最高位化为0，即，再判断次高位；与进行比较，若为0，将其变位1，即，即得到倒序数，如果次高位为1，将其化为0，再判断下一位……即从高位到低位依次判断其是否为1，为1将其变位0，若这一位为0，将其变位1，即可得到倒序数。若倒序数小于顺序数，进行换位，否则不变，防治重复交换，变回原数。
注：因为0的倒序数为0，所以可从1开始进行求解。
2、 程序设计框图
（1）倒序算法——雷德算法流程图
(2)FFT算法流程
3、
3、 FFT源程序
void fft(x,n)
int n;
double x[];
{int i,j,k,l,m,n1,n2;
double c,c1,e,s,s1,t,tr;
for(j=1,i=1;i
4、 //如果i
5、 //把0改为1
}
for(i=0;i
6、tr=x[i];
x[i]=tr+x[i+n2];
x[i+n2]=tr-x[i+n2];
x[i+n2+n4]=-x[i+n2+n4];
a=e;
for(j=2;j<=(n4-1);j++) //控制计算不同种蝶形结，即计算系数不同的蝶形结
{i1=i+j;
i2=i-j+n2;
i3=i+j+n2;
i4=i-j+n1;
cc=cos(a);
ss=sin(a);
a=a+e;
t1=cc*x[i3]+ss*x[i4];
t2=ss*x[i3]-cc*x[i4
7、];
x[i4]=x[i2]-t2;
x[i3]=-x[i2]-t2;
x[i2]=x[i1]-t1;
x[i1]=x[i1]+t1;
}
}
}
}
4、 计算实例及运行结果
设输入序列为
其离散傅里叶变换为
这里。选n=512,计算离散傅里叶变换。
所用软件为Turbo c 2.0，操作界面如图1所示
图1 Turbo c 2.0操作界面
程序运行结束后的界面如图2所示
图2 程序运行后的界面
例子的具体程序如下：
#include
#include
8、
#include
#define pi 3.14159265359
void fft(x,n)
int n;
double x[];
{int i,j,k,l,i1,i2,i3,i4,n4,m,n1,n2;
double a,e,cc,ss,tr,t1,t2;
for(j=1,i=1;i
9、tr;
}
k=n/2;
while(k<(j+1))
{j=j-k;
k=k/2;
}
j=j+k;
}
for(i=0;i
10、
x[i+n2]=tr-x[i+n2];
x[i+n2+n4]=-x[i+n2+n4];
a=e;
for(j=2;j<=(n4-1);j++)
{i1=i+j;
i2=i-j+n2;
i3=i+j+n2;
i4=i-j+n1;
cc=cos(a);
ss=sin(a);
a=a+e;
t1=cc*x[i3]+ss*x[i4];
t2=ss*x[i3]-cc*x[i4];
x[i4]=x[i2]-t2;
x[i3]=-x[i2]-t2;
x[i2
11、]=x[i1]-t1;
x[i1]=x[i1]+t1;
}
}
}
}
main()
{FILE *p;
int i,j,n;
double dt=0.001;
double x[512];
p=fopen("d:\123.c","w");
n=512;
for(i=0;i
12、
printf("\n");
}
fft(x,n);
fprintf(p,"\n DISCRETE FOURIER TRANSFORM\n");
printf("\n DISCRETE FOURIER TRANSFORM\n");
fprintf(p,"%10.7f",x[0]);
printf("%10.7f",x[0]);
fprintf(p,"%10.7f+J%10.7f\n",x[1],x[n-1]);
printf("%10.7f+J%10.7f\n",x[1],x[n-1]);
for(i=2;i
13、tf(p,"%10.7f+J%10.7f",x[i],x[n-i]);
fprintf(p,"%10.7f+J%10.7f",x[i+1],x[n-i-1]);
fprintf(p,"\n");
printf("%10.7f+J%10.7f",x[i],x[n-i]);
printf("%10.7f+J%10.7f",x[i+1],x[n-i-1]);
printf("\n");
}
fprintf(p,"%10.7f",x[n/2]);
printf("%10.7f",x[n/2]);
fprintf(p,"%10.7f+J%10.7f\n",x
14、[n/2-1],-x[n/2+1]);
for(i=2;i
快速傅里叶变换(FFT)原理及源程序
