8点DIF-FFT课程设计(共20页).docx

《8点DIF-FFT课程设计(共20页).docx》由会员分享，可在线阅读，更多相关《8点DIF-FFT课程设计(共20页).docx（20页珍藏版）》请在淘文阁 - 分享文档赚钱的网站上搜索。

1、精选优质文档-倾情为你奉上摘要快速傅里叶FFT是将一个大点数N的DFT分解为若干小点的DFT的组合。工作量会明显降低，从而大大提高离散傅里叶变换DFT的运算速度。因而各个学科技术领域广泛的使用了FFT技术，她大大推动了信号处理技术的进步，现已成为数字信号处理强有力的工具，而DIF是FFT算法中按频率抽取的方法。本课设比较全面的叙述快速傅里叶变换中DIF的算法、原理、特点，并完成了基于MATLAB的实现。关键词：快速傅里叶变换FFT、MATLAB、DIF1 FFT算法简介及原理特点1.1 FFT算法简介快速傅立叶变换（FFT）并不是一种新的变换，而是离散傅立叶变换（DFT）的一种快速算法。DFT

2、的计算在数字信号处理中非常有用。例如在FIR滤波器设计中会遇到从h(n)求H(k)或由H(k)计算h(n)，这就要计算DFT；信号的谱分析对通信、图像传输、雷达等都是很重要的，也要计算DFT。因直接计算DFT的计算量与变换区间长度N的平方成正比，当N较大时，计算量太大。自从1965年图基（J. W. Tukey）和库利（T. W. Coody）在计算数学（Math. Computation , Vol. 19, 1965）杂志上发表了著名的机器计算傅立叶级数的一种算法论文后，桑德(G. Sand)-图基等快速算法相继出现，又经人们进行改进，很快形成一套高效运算方法，这就是快速傅立叶变换简称FF

3、T（Fast Fourier Transform）。这种算法使DFT的运算效率提高12个数量级。1.2 基2 FFT算法原理特点1.2.1直接计算DFT的问题设x(n)为N点有限长序列，其DFT正变换为= ， k=0，1，N-1其反变换(IDFT) x(n)= ，n=0，1，N-1二者的差别只在于的指数符号不同，以及差一个常数乘因子1/N，因而下面我们只讨论DFT正变换的运算量，反变换的运算量是完全相同的。考虑x(n)为复数序列的一般情况，每计算一个X(k)，需要N次复数乘法以及（N-1）次复数加法。因此，对所有N个k值，共需N2次复数乘法及N（N-1）次复数加法运算。所以直接计算DFT，乘法

4、次数和加法次数都是和N2成正比的，当N很大时，运算量是很可观的，因而需要改进对DFT的计算方法，以减少运算次数。1.2.2改进的途径及DIF的引出 下面讨论减少运算工作量的途径。仔细观察DFT的运算就可看出，利用系数以下固有特性，就可减小DFT的运算量：（1）的对称性 （）*=（2）的周期性 =（3）的可约性 =由此可得：=，=-1，=-。这样利用这些特性，可以将长序列的DFT分解为短序列的DFT。快速傅立叶变换算法正是基于这样的思路而发展起来的。它的算法基本上可以分成两大类：时域抽取法FFT（DIT-FFT）和频域抽取法FFT（DIF-FFT）。2 按频率抽选DIF的基2 FFT算法2.1

5、DIF的原理设序列点数为N=2M ，M为整数。 ，k=0，1，N-1先把输入按n 的顺序分成前后两半： ，k=0，1，N-1=+=+=，k=0，1，N-1 1，k为偶数=(-1)k= -1，k为奇数将X（K）分解成偶数组与奇数组，当k取偶数即k=2r时（r=0，1，N/2-1），=当k取奇数即k=2r+1时（r=0，1，N/2-1），=令 n=0，1，N/2-1 =，r=0，1，N/2-1= x1(n)、x2(n)和x(n)之间可用下图所示的蝶形运算符号表示： -1 图2.1用下图表示，N=8的一次分解流图：图2.2由于N=2M，N/2仍是偶数，继续将N/2点DFT分成偶数组和奇数组。图2.3

6、表示N=8时二次分解运算流图：图2.3最后完整的分解流图为下图：图2.4这种算法是对X（K）进行奇偶抽取分解的结果，所以称之为频域抽取法FFT。DIF-FFT算法与DIT-FFT算法类似，不同的是DIF-FFT算法输入序列为自然顺序，而输出为倒序排列。另外，蝶形运算略有不同，DIT-FFT蝶形先乘后加，而DIF-FFT蝶形先加后乘。上述两种FFT的算法流图形式不是唯一的。只要保证各节点所连支路及其传输系数不变，改变输入与输出点以及中间结点的排列顺序，就可以得到其他变形的FFT运算流图。2.2 DIF的特点1. 原位运算由图2.4可以看出，DIF-FFT的运算过程很有规律。N=2M点的FFT共进

7、行M级运算，每级由N/2个蝶形运算组成。同一级中，每个蝶形的两个输入数据只对计算本蝶形有用，而且每个蝶形的输入、输出数据结点又同在一条水平线上，这就意味着计算完一个蝶形后，所得输出数据可立即存入原输入数据所占用的存储单元。这样，经过M级运算后，原来存放输入序列数据的N个存储单元中便依次存放X（K）的N个值。这种利用同一存储单元存储蝶形计算输入、输出数据的方法称为原位计算。原位计算可节省大量内存，从而使设备成本降低。2. 旋转因子的变化规律以8点的FFT为例，第一级蝶形，r=0,1，2；第二级蝶形，r=0,1；第三级的蝶形，r=0。依次类推，对于M级蝶形，旋转因子的指数为 r=J2M-1 ，J=

8、0，1，2，3，这样就可以算出每一级的旋转因子。对于M级的任一蝶形运算所对应的旋转因子的指数，可以 如下方法得到：1将待求的蝶形输入节点中上面节点的行标号值k写成L位二进制数；2将此二进制数乘以2的M-1次方，即将L位二进制数左移M-1位，右边的空位补零，然后从低位到高位截取L位，即所得指数r所对应的二进制数。3. 蝶形运算两节点之间的“距离”第一级蝶形每个蝶形运算量节点的“距离”为4，第二级每个蝶形运算另节点的“距离”为2，第三级蝶形每个蝶形运算量节点的“距离”为1。依次类推：对于N等于2的L次方的DIF-FFT，可以得到第M级蝶形每个蝶形运算量节点的“距离”为2的L-M次方。3 DIF的编

9、程思想及MATLAB实现3.1 DIF的编程思想DIF-FFT程序包括变址（倒位序）和L级递推计算（N=，L为正整数）两部分。1变址DIF-FFT是输出倒位序的变址处理，设x (i) 表示存放自然顺序输入数据的内存单元，x（j）表示存放倒位序序数的内存单元，I、J=0，1，N-1，当I=J时，不用变址；当I J时，需要变址；但是当ij时，进行变址在先，故在IJ时，就不需要再变址了，否则变址两次等于不变。其中本程序使用的“反向进位加法”。也可以用bin2dec函数可以实现倒位序。2L级递推计算整个L级递推过程由三个嵌套循环构成。外层的一个循环控制L（L=）级的顺序运算；内层的两个循环控制同一级（

10、M相同）各蝶形结的运算，其中最内层循环控制同一种（即中的r相同）蝶形结的运算。其循环变量为I，I用来控制同一种蝶形结运算。其步进值为蝶形结的间距值LE=，同一种蝶形结中参加运算的两节点的间距为LE1=点。第二层循环，其循环变量J用来控制计算不同种（系数不同）的碟形结，J的步进值为1。也可以看出，最内层循环完成每级的蝶形结运算，第二层循环则完成系数的运算。最外层循环，用循环变量M来控制运算的级数，M为1到L，步进值为1，当M改变时，则LE1，LE和系数U都会改变。3.2 MATLAB实现通过上面的编程思想分析，可得出如下的MATLAB实现的代码：function Xk=DPS-DIF(xn,N)

11、%本程序对输入序列实现DIF基2算法，点数取小于等于长度的2的幂次 N=8;n=log2(2nextpow2(length(xn); %求得x长度对应的2的最低幂次nN=2n; if length(xn)N xn=xn,zeros(1,N-length(xn); %若长度不是2的幂，补0到2的整数幂 end %蝶形运算开始M=log2(N); %“级”的数量for m=0:M-1 %“级”循环开始 Num =2m; %每一级中组的个数 Interval_G=N/2m; %每一级中组与组之间的间距 Interval_U=N/2(m+1); %每一组中相关运算单元之间的间距 Cycle_Count

12、= Interval_U -1; %每一组内的循环次数 Wn=exp(-j*2*pi/Interval_G); %旋转因子 for g=1:Num %“组”循环开始 Interval_one=(g-1)*Interval_G; %第g组中蝶形运算变量1的偏移量 Interval_two=(g-1)*Interval_Gp+Interval_U; %第g组中蝶形运算变量2的偏移量 for r=0:Cycle_Count; %“组内”循环开始 k=r+1; %“组内”序列的下标 xn(k+Interval_one)=xn(k+Interval_one)+xn(k+Interval_two); %第

13、m级，第g组的蝶形运算式1 xn(k+Interval_two)=xn(k+Interval_one)-xn(k+Interval_two)-xn(k+Interval_one)*Wnr; %第m级，第g组的蝶形运算式2 end endend%序列排序开始n1=fliplr(dec2bin(0:N-1);%码位倒置步骤1：将码位转换为二进制，再进行倒序n2=bin2dec(n1); %码位倒置步骤2：将码位转换为十进制后翻转for i=1:N Xk(i)=xn(n2(i)+1); %根据码位倒置的结果，将xn重新排序,存入Xk中end4 程序验证4.1 8点的FFT运算 编写主函数，在主函数中

14、输入一个8点序列分别调用自己编写的FFT函数，和MATLAB本身系统的FFT函数并比较两个结果是否相等，以判断自己编写的FFT程序是否正确。xn=0:7;m=1:8;N=8x1=DPS-DIF (xn,N)x2=fft(xn)x3=abs(x1);x4=abs(x2);x5=angle(x1);x6=angle(x2);subplot(3,1,1)stem(m,xn);title(输入的离散序列)subplot(3,1,2)stem(m,x3);title(经过DIF_FFT后得到的频谱的幅度)subplot(3,1,3)stem(m,x5);title(经过DIF_FFT后得到的频谱的相位)

15、figure (2)subplot(3,1,1)stem(m,xn);title(输入的离散序列)subplot(3,1,2)stem(m,x4);title(经过库函数fft后得到的频谱的幅度)subplot(3,1,3)stem(m,x6);title(经过库函数fft后得到的频谱的相位)当使用编写的程序进行8点的DIF-FFT计算时结果如下：N =8x1 =Columns 1 through 428.0000 -4.0000 + 9.6569i -4.0000 + 4.0000i -4.0000 + 1.6569iColumns 5 through 8-4.0000 -4.0000 -

16、1.6569i -4.0000 - 4.0000i -4.0000 - 9.6569ix2 =Columns 1 through 428.0000 -4.0000 + 9.6569i -4.0000 + 4.0000i -4.0000 + 1.6569iColumns 5 through 8-4.0000 -4.0000 - 1.6569i -4.0000 - 4.0000i -4.0000 - 9.6569i并得到如下仿真图：图4.1 调用自编程序的8点FFT图4.2 调用库函数的8点FFT4.2 16点的FFT运算跟8点的相同，在主函数中输入一个16点序列分别调用自己编写的FFT函数，和M

17、ATLAB本身系统的FFT函数并比较两个结果是否相等，以判断自己编写的FFT程序是否正确。xn=0:15;m=1:16;N=16x1=DPS-DIF (xn,N)x2=fft(xn)x3=abs(x1);x4=abs(x2);x5=angle(x1);x6=angle(x2);subplot(3,1,1)stem(m,xn);title(输入的离散序列)subplot(3,1,2)stem(m,x3);title(经过DIF_FFT后得到的频谱的幅度)subplot(3,1,3)stem(m,x5);title(经过DIF_FFT后得到的频谱的相位)figure (2)subplot(3,1,

18、1)stem(m,xn);title(输入的离散序列)subplot(3,1,2)stem(m,x4);title(经过库函数fft后得到的频谱的幅度)subplot(3,1,3)stem(m,x6);title(经过库函数fft后得到的频谱的相位)当使用编写的程序进行8点的DIF-FFT计算时结果如下：N =16x1 =1.0e+002 *Columns 1 through 4 1.2000 -0.0800 + 0.4022i -0.0800 + 0.1931i -0.0800 + 0.1197iColumns 5 through 8 -0.0800 + 0.0800i -0.0800 +

19、0.0535i -0.0800 + 0.0331i -0.0800 + 0.0159iColumns 9 through 12 -0.0800 -0.0800 - 0.0159i -0.0800 - 0.0331i -0.0800 - 0.0535iColumns 13 through 16 -0.0800 - 0.0800i -0.0800 - 0.1197i -0.0800 - 0.1931i -0.0800 - 0.4022ix2 =1.0e+002 *Columns 1 through 4 1.2000 -0.0800 + 0.4022i -0.0800 + 0.1931i -0.08

20、00 + 0.1197iColumns 5 through 8 -0.0800 + 0.0800i -0.0800 + 0.0535i -0.0800 + 0.0331i -0.0800 + 0.0159iColumns 9 through 12 -0.0800 -0.0800 - 0.0159i -0.0800 - 0.0331i -0.0800 - 0.0535iColumns 13 through 16 -0.0800 - 0.0800i -0.0800 - 0.1197i -0.0800 - 0.1931i -0.0800 - 0.4022i并得到如下仿真图：图4.3 调用自编程序的1

21、6点FFT图4.4 调用库函数的16点FFT通过观察比较，得到的序列各点的值以及直观的通过图形，可以得到自己编写的DIF_FFT函数实现对序列进行FFT变换得到的结果与库函数FFT得到的结果是一样的。说明DIF_FFT子程序是正确的。从图中也可以看出有限长序列通过FFT后得到的频域为离散的。从理论讲，有限长序列经过离散傅里叶变换后，得到的频谱为离散的，从而也说明了FFT是DFT的优化方法，也属于DFT。这个程序可以实现基2FFT,但是如果想在运行时直接输入要变换的点数就不行，必须在调用FFT函数前现将要算的序列定义好，这是这个程序的不足之处。但是该程序可以计算不是2的整数次幂的序列。所以在主程

22、序中，输入序列必须给出才能进行FFT变换。5心得与体会通过做这次程序设计，我对MATLAB编程有了进一步的掌握，对数字信号处理在MATLAB中的实现有了更深的体会，同时也对数字信号处理的的理论知识有了更深刻的认识，在学习基2FFT算法时有许多地方不理解比如如何进行倒位序，蝶形运算是怎么回事等，通过重新复习相关知识，编写程序对以前一些迷惑的地方理解的更透彻，学会了理论与实践相结合的方法。快速傅里叶FFT是将一个大点数N的DFT分解为若干小点的DFT的组合。工作量会明显降低，从而大大提高离散傅里叶变换DFT的运算速度。因而各个学科技术领域广泛的使用了FFT技术，她大大推动了信号处理技术的进步，现已

23、成为数字信号处理强有力的工具，而DIF是FFT算法中按频率抽取的方法。本课设比较全面的叙述快速傅里叶变换中DIF的算法、原理、特点，并完成了基于MATLAB的实现。理论是实践的基础，只有掌握了相关的理论知识才能更好更轻松的实践。当理论某些细节不是很理解时，可以通过编程仿真来实现，将仿真结果与理论结合起来进行对比理解这样会容易点。不过这限于一些小的地方，当很多都不懂时，就不行了因为你很多地方不懂时就不能实现编程了。如果自己根据FFT的计算方法，计算一个序列经FFT变换后的结果，会花很长的实践但是直接调用已编好的程序，会轻松的得到结果。可见，使用程序的好处。在以后的学习中，可以多根据理论，在编程上去实现，这样也有助于学习。参考文献1 刘泉. 数字信号处理原理与实现. 电子工业出版社. 2005. 62 程佩青. 数字信号处理教程. 清华大学出版社. 2008. 5 3 刘卫国. MATLAB程序设计教程 . 中国水利水电出版社 . 2008. 64 余成波. 数字信号处理及其MATLAB实现. 清华大学出版社. 2008. 25 陈怀琛. 数字信号处理教程：MATLAB释义及实现. 电子出版社 . 2006.6专心-专注-专业