(完整word版)快速傅里叶变换(FFT)原理及源程序(20210413200648)

1、测试信号分析及处理课程作业 快速傅里叶变换 一、程序设计思路 快速傅里叶变换的目的是减少运算量，其用到的方法是分级进行运算。全部 计算分解为M级，其中M =|og2N ；在输入序列Xi中是按码位倒序排列的， 输出序列X k是按顺序排列；每级包含 N个蝶形单元，第i级有个群，每个 2 2 群有21 个蝶形单元；每个蝶形单元都包含乘 WN和-wN系数的运算，每个蝶形 单元数据的间隔为2iJ ,i为第i级；同一级中各个群的系数 W分布规律完全相 同。 将输入序列Xi按码位倒序排列时，用到的是倒序算法一一雷德算法。自然 序排列的二进制数，其下面一个数总比上面的数大 1,而倒序二进制数的下面一 个数是上

2、面一个数在最高位加1并由高位向低位仅为而得到的。 若已知某数的倒序数是J ,求下一个倒序数，应先判断J的最高位是否为0， 与k = N进行比较即可得到结果。如果k J，说明最高位为0,应把其变成1, 2 即J -，这样就得到倒序数了。如果J，即J的最高位为1，将最高位化为 2 0,即J - N ,再判断次高位；与k =进行比较，若为0,将其变位1,即J , 244 即得到倒序数，如果次高位为1,将其化为0，再判断下一位即从高位到低 位依次判断其是否为1,为1将其变位0,若这一位为0,将其变位1,即可得到 倒序数。若倒序数小于顺序数，进行换位，否则不变，防治重复交换，变回原数。 注：因为0的倒序

3、数为0,所以可从1开始进行求解。 二、程序设计框图 (1)倒序算法一一雷德算法流程图 三、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 n/2;i+) m=i; /得到流程图的共几级 /如果ij,即进行变址 j=2*j; if(j=n )break; n仁n-1; for(j=0,i=0;i n1;i+) if(ivj) tr=xj; xj=xi; xi=tr; k=n/2; while(k(j+1) /求j的下一个倒位序 /如果k(j+1),

4、表示j的最高 位为1 j=j-k;/把最高位变成0 k=k/2;/k/2，比较次高位，依次类推，逐个比较，直到某个位为0 j=j+k;/ 把 0 改为 1 for(i=0;i n ;i+=2) tr=xi; xi=t 叶xi+1; xi+1=tr-xi+1; n2=1; for(l=1;l=m;l+)/ 控制蝶形结级数 n 4=n2; n 2=2* n4; n1=2* n2; e=6.28318530718/n1; for(i=0;in;i+=n1)/控制同一蝶形结运算，即计算系数相同蝶形结 tr=xi; xi=t 叶xi+n2; xi+n 2=tr-xi+n2; xi+n2+n4=-xi+n

5、2+n4; a=e; for(j=2;j=(n4-1);j+)/控制计算不同种蝶形结，即计算系数不同 的蝶形结 i仁 i+j; i2=i-j+n2; i3=i+j+n2; i4=i-j+n1; cc=cos(a); ss=s in( a); a=a+e; t1=cc*xi3+ss*xi4; t2=ss*xi3-cc*xi4; xi4=xi2-t2; xi3=-xi2-t2; xi2=xi1-t1; xi1=xi1+t1; 四、计算实例及运行结果 设输入序列x(i)为 x(i) =sin200二i：t,(i =0,1,2,，n1) 其离散傅里叶变换为 N 二 X(k)八 x(i)wNk,(k =

6、0,1,2,n-1) i ：5 .j 这里W =e n。选n=512,计算离散傅里叶变换X(k) 所用软件为Turbo c 2.0,操作界面如图1所示 图1 Turbo c 2.0操作界面 程序运行结束后的界面如图2所示 回|2 国3 C:Wi ndc32c - . exe 5.1911638+J 3-6933914+J 2.B497513+J 2.3099713+J 1-5360773+J 1 .6627248+J 1-4549735+J 1.2924530+J L.1 double x; int i,j,k,l,i1,i2,i3,i4,n4,m,n1,n2; double a,e,cc,s

7、s,tr,t1,t2; for(j=1,i=1;i n/2;i+) m=i; j=2*j; if(j=n )break; n仁n-1; for(j=0,i=0;i n1;i+) if(ivj) tr=xj; xj=xi; xi=tr; k=n/2; while(k(j+1) j=j-k; k=k/2; j=j+k; for(i=0;i n ;i+=2) tr=xi; xi=t 叶xi+1; xi+1=tr-xi+1; n2=1; for(l=1;l=m;l+) n 4=n2; n 2=2* n4; n1=2* n2; e=6.28318530718/n1; for(i=0;i n ;i+=n1

8、) tr=xi; xi=t 叶xi+n2; xi+n 2=tr-xi+n2; xi+n2+n4=-xi+n2+n4; a=e; for(j=2;j=( n4-1);j+) i仁 i+j; i2=i-j+n2; i3=i+j+n2; i4=i-j+n1; cc=cos(a); ss=s in( a); a=a+e; t1=cc*xi3+ss*xi4; t2=ss*xi3-cc*xi4; xi4=xi2-t2; xi3=-xi2-t2; xi2=xi1-t1; xi1=xi1+t1; mai n() FILE *p; int i,j, n; double dt=0.001; double x512

9、; p=fope n(d:123.c,w); n=512; for(i=0;i n;i+) xi=si n(200*pi*i*dt); for(i=0;i n;i+) fprin tf(p,%10.7f,xi); fprin tf(p,n); prin tf(%10.7f,xi); prin tf(n); fft(x, n); fprin tf(p,n DISCRETE FOURIER TRANSFORMS); prin tf(n DISCRETE FOURIER TRANSFORMS); fprin tf(p,%10.7f,x0); prin tf(%10.7f,x0); fprin tf(

10、p,%10.7f+J%10.7fn,x1,x n-1); prin tf(%10.7f+J%10.7fn,x1,x n-1); for(i=2;i n/2;i+=2) fprin tf(p,%10.7f+J%10.7f,xi,x n-i); fprin tf(p,%10.7f+J%10.7f,xi+1,x n-i-1); fprin tf(p,n); prin tf(%10.7f+J%10.7f,xi,x n-i); prin tf(%10.7f+J%10.7f,xi+1,x n-i-1); prin tf(n); fprin tf(p,%10.7f,x n/2); prin tf(%10.7f,x n/2); fprin tf(p,%10.7f+J%10.7fn,x n/2-1,-x n/2+1); for(i=2;i n/2;i+=2) fprin tf(p,%10.7f+J%10.7f,x n/2-i,-x n/2+i); fprin tf(p,%10.7f+J%10.7f,x n/2-i-1,-x n/2+i+1); fprin tf(p,n); prin tf(%10.7