基-2FFT算法的软件实现(共13页)

1、实验二 基-2FFT算法的软件实现一、实验目的1、 加深对DFT算法原理和基本性质的理解；2、 熟悉FFT算法的流程；3、 了解FFT算法的应用。二、基本原理1、 DFT算法原理（见教材第三章）2、 按时间抽取（DIT）的-2FFT算法（1）算法原理序列x(n)的N(N=2-M)点DFT为，k=0, 1, , N-1 （2.1）将式（2.1）按n的奇偶性分解为 （2.2）令, ，因为, 所以式（2.2）可写成 （2.3）式（2.3）说明，按n的奇偶性将x(n)分解为两个N/2长的序列x1(l)和x2(l)，则N点DFT可分解为两个N/2点DFT来计算。用X1(k)和X2(k)分别表示 （2.4

2、） （2.5）将(2.4)式和(2.5)式代入(2.31)式，并利用和X1(k)、 X2(k)的隐含周期性可得到： （2.6）这样，将N点DFT的计算分解为计算两个N/2点离散傅立叶变换X1(k)和X2(k)，再计算(2.6)式。为了将如上分解过程用运算流图表示，以便估计其运算量，观察运算规律，总结编程方法，先介绍一种表示(2.6)式的蝶形图。图2.1蝶形运算图图2.2 8点DFT一次时域抽取分解运算流图根据图2.2可以求得第一次分解后的运算量。图2.2包括两个N/2点DFT和N/2个蝶形，每个N/2点DFT需要次复数乘法和次复数加法运算，每个蝶形只有一次复数乘法运算和两次复数加法运算。所以，

3、总的复数乘法次数为总的复数加法次数为 （2.7） （2.8）N=8点DIT-FFT的运算流图如图2.3(a)所示。根据WkN/m=WkmN，将图2.3(a)转换成如图2.3(b)所示的标准形式的运算流图图2.3 N=8点DIT-FFT的运算流图（2）算法效率由图2.3可见，N=2M时，DIT-FFT运算流图由M级蝶形构成，每级有N/2个蝶形。因此，每级需要N/2次复数乘法运算和N次复数加法运算，M级形共需复数乘法次数CM(2)和复数加法次数CA(2)分别为 （2.9）CA(2) =N·M=N lb N （2.10）式中，lb N=log2 N。直接计算N点DFT的复数乘法次数为N2,

4、复数加法次数为(N-1)N。当N1时, N2/CM(2)>>1，所以N越大，DIT-FFT运算效率越高。DIT-FFT算法与DFT所需乘法次数与N的关系曲线如图2.4所示。例如，N=210=1024时，DIT-FFT的运算效率为 （2.11）而当N=211=2048时， （2.12）图2.4 DIT-FFT与DFT所需乘法次数比较曲线（3）运算规律的确定第m级运算，一个DIT蝶型运算的的两接点“距离”为，所以 （2.13）r的求解：（1）将式（2.13）中第一个节点的标号k表示成L（）位二进制数；（2）把此二进制数乘上，即将L位二进制数左移位（注意m是第m级运算），把右边空出的位置

5、补0，此数即为所求r的二进制数。序列的倒序DIT-FFT算法的输入序列的排序看起来似乎很乱，但仔细分析就会发现这种倒序是很有规律的。由于，因此顺序数可用L位二进制数()表示。M次偶奇时域抽选过程如图2.5所示。第一次按最低位的0和1将x(n)分解为偶奇两组，第二次又按次低位的0、 1值分别对偶奇组分解； 依次类推，第L次按位分解，最后所得二进制倒序数如图2.5所示。图2.5 形成例序的树状图(N=23)形成倒序J后，将原存储器中存放的输入序列重新按倒序排列。设原输入序列x(n)先按自然顺序存入数组A中。例如，对N=8, A(0)，A(1)，A(2)，A(7)中依次存放着x(0),x(1), ,

6、 x(7)。对x(n)的重新排序(倒序)规律如图2.6所示。图2.6倒序规律三、实验仪器计算机四、实验要求及内容用所学过的编程语言，自行设计出一个按时间抽取的、输入倒位序、输出顺序的基-2 FFT算法程序。五、实验报告（1）简述实验目的及原理；（2）画出程序流程框图；（3）主要给出实验内容的程序（要求程序模块化并加注释）。 程序流程框图实验的程序#include<iostream.h>#include<math.h>#include<string.h>#define PI 3.1415926class Plural/复数类float real;/实部floa

7、t img;/虚部public:Plural(float r,float i)real=r;img=i;Plural(float r)/通过构造函数重载使用数与复数乘法real=r;img=0;Plural()real=0;img=0; friend Plural* fft(float X,int n); /fft（）；friend Plural operator *(Plural p1,Plural p2);/重载乘法运算符friend Plural operator -(Plural p1,Plural p2);/重载减法运算符friend Plural operator +(Plural

8、 p1,Plural p2);/重载加法运算符friend Plural* daoxu(Plural* x,int n); /倒序 Plural operator =(Plural p) ; /重载赋值符号 friend Plural W(int N,int i) ; /生成旋转因子 void show() ; /输出real+imgi;Plural operator *(Plural p1,Plural p2)/复数乘法return Plural(p1.real *p2.real -p1.img *p2.img ,p1.real *p2.img +p1.img *p2.real );Plur

9、al operator +(Plural p1,Plural p2)/复数加法return Plural(p1.real+p2.real,p1.img +p2.img );Plural operator -(Plural p1,Plural p2)/复数加法return Plural(p1.real-p2.real,p1.img -p2.img );Plural Plural:operator =(Plural p) /重载赋值符号 real=p.real ; img=p.img ; return *this;/*生成旋转因子*Plural W(int N,int i) float r;flo

10、at im;r=(float)cos(2*PI*float(i)/float(N);im=(float)sin(-1)*2*PI*float(i)/float(N);return Plural(r,im);/*输出函数show（）*void Plural:show() /输出real+imgi if(real=0) if(img=0) cout<<0<<" " else cout<<img<<"i " else if(img<0) cout<<real<<img<<

11、"i " else if(img=0) cout<<real<<" " else cout<<real<<"+"<<img<<"i " /*2的n次方*/int mi2(int n)int m=1;for(int i=0;i<n;i+)m*=2;/m=2nreturn m;/*由二进制数转化为十进制数*/int twoten(int xu,int i)int m=0;for(int j=0;j<i;j+) m+=xuj*mi2(i-

12、j-1);/m=xu0*1+xu1*2+xu2*4+xu3*8+xu4*16+.return m;/*对二进制数倒序加一*/void add(int xu,int i)xu0+;for(int j=0;j<i;j+)if(xuj=2) xuj=0; xuj+1+;/例如xu=0000,1000,0100,1100，0010./*倒序*/void daoxu(Plural x,int n)float m=float (n);int i=0; for(i=0;m>1;i+)/得到n是2的多少次方m/=2;/m=log(n)/log(2) int M=mi2(i); int* xu=ne

13、w inti;/定义一个长度为i的数组，当做2进制数for(int j=0;j<i;j+)xuj=0; Plural* jia=new PluralM;/定义一个长度为M的数组jia，以保存倒序 for(j=0;j<M;j+) int mm; mm=twoten(xu,i);/得到数组xu所对应的十进制数 jiaj=xmm; add(xu,i);/二进制数组最左端加1 for(j=0;j<n;j+)/撤销动态生成的数组xj=jiaj;if(jia) delete jia;/*快速傅里叶变换FFT()*Plural* fft(Plural X,int m)int n=1;for

14、(;m>mi2(n);n+);/得到m等于n次方Plural* A=new Pluralm;for(int i=0;i<n;i+)/共有n级for(int j=0;j<mi2(n-1-i);j+)/第i级有2(n-i)族for(int k=0;k<mi2(i);k+)/每族有2(m-1)个蝶形运算int m,n;m=k+j*mi2(i+1);n=k+mi2(i)+j*mi2(i+1); Plural a1,a2;if(i=0)a1=Xm+W(mi2(i+1),k)*Xn; /蝶形运算 a2=Xm-W(mi2(i+1),k)*Xn;/蝶形运算 Am=a1; An=a2;e

15、lse a1=Am+W(mi2(i+1),k)*An; /蝶形运算 a2=Am-W(mi2(i+1),k)*An;/蝶形运算 Am=a1; An=a2;for(i=0;i<m;i+)/撤销动态生成的数组Xi=Ai;delete A;return X;/*主函数*/void main()int N=0; cout<<"输入变换的点数N："<<endl;cin>>N;Plural *x=new PluralN ;cout<<"数据 实部：虚部（少于N点时以#结束）:"<<endl; for(i

16、nt i=0;i<N;i+)char s;float x1,x2;cin>>x1;s=(char)x1;if(s='#') break;elsecin>>x2;xi=Plural(x1,x2); x1=0;x2=0; cout<< "您输入的数组："<<endl; for( i=0;i<N;i+) xi.show(); if(i+1)%4=0) cout<<endl;cout<<endl;daoxu(x,N);/调用倒序daoxu()函数cout<< "倒序后的数组："<<endl; for( i=0;i<N;i+) xi.s