8点基于DIF-FFT课程设计

武汉理工大学《信号分析与处理》课程设计11数字信号处理简介1.1数字信号处理概念相对于模拟信号，数字信号是时间、幅值都离散的信号，数字信号处理就是通过计算机或专用处理设备，用数值计算等数字方式对数字信号进行各种处理，将数字信号变换成符合要求的某种形式。数字信号处理主要包括数字滤波和数字频谱分析两大部分。例如，对数字信号进行滤波，限制其频带或滤除噪声和干扰，以提取和增强信号的有用分量；对信号进行频谱分析或功率分析，了解信号的频谱组成，以对信号进行识别。当然，凡是用数字方式对信号进行滤波，变换，增强，压缩，估计和识别等都是数字信号处理的研究范畴。数字信号处理在理论上所涉及的范围极其广泛。数字领域中的微积分，概率统计，随机过程，高等代数，数值分析，复变函数和各种变换（如傅里叶变换，Z变换，离散傅里叶变换，小波变换等）都是它的基本工具，网络理论，信号与系统等则是它的理论基础。在科学发展上，数字信号处理又和最优控制，通信理论等紧密相连，目前已成为人工智能，模式识别，神经网络等新兴学科的重要理论基础，其实现技术又和计算机科学和微电子技术密不可分。因此，数字信号处理是把经典的理论基础体系作为自身的理论基础，同时又使自己成为一系列新兴学科的理论基础。1.2数字信号处理的特点及实现方法与模拟信号处理相比，数字信号处理具有高精度、高稳定性、灵活性好、易于大规模集成等显着的优点。数字信号处理的主要研究对象是数字信号，且采用数值运算的方法达到处理的目的。数字信号处理的实现方法基本上可以分为软件实现方法、硬件实现方法和软硬件相结合的实现方法。数字信号处理的理论、算法和实现方法三者是密不可分的。武汉理工大学《信号分析与处理》课程设计22FFT算法介绍2.1DFT的定义对于有限长离散数字信号{x[n]}，0nN-1，其离散谱{X[k]}可以由离散傅氏变换（DFT）求得。DFT的定义为：210[][]NjnkNnXkxne，k=0,1,…N-1通常令2jNNeW，称为旋转因子。2.2FFT基本思想由DFT的定义可以看出，在x[n]为复数序列的情况下，直接运算N点DFT需要N的二次方次复数乘法和N（N-1）次加法。因此，对于一些相当大的N值（如1024）来说，直接计算它的DFT所作的计算量是很大的，对实现快速计算是很不利的。FFT的基本思想在于，将原有的N点序列分成两个较短的序列，这些序列的DFT可以很简单的组合起来得到原序列的DFT。例如，若N为偶数，将原有的N点序列分成两个（N/2）点序列，那么计算N点DFT将只需要约[(N/2)2·2]=N2/2次复数乘法。即比直接计算少作一半乘法。因子（N/2）2表示直接计算（N/2）点DFT所需要的乘法次数，而乘数2代表必须完成两个DFT。上述处理方法可以反复使用，即（N/2）点的DFT计算也可以化成两个（N/4）点的DFT（假定N/2为偶数），从而又少作一半的乘法。这样一级一级的划分下去一直到最后就划分成两点的FFT运算的情况。虽然这样做，要求序列的长度为n2（n=0,1,…）,但在原序列后面补零使其满足条件长度，并不影响对原序列的频域分析，仅仅只是抽样点更密了一些。DIT-FFT和DIF-FFT正是基于这种思想。武汉理工大学《信号分析与处理》课程设计32DIF-FFT2.1DIF-FFT过程推导设序列长度为2LN，L为整数（如果序列长度不满足此条件，通过在后面补零让其满足）。在把[]Xk按k的奇偶分组之前，把输入按n的顺序分成前后两半：1/21100/2/21/21()200/2120[]{[]}[][][][][]2[[][]],0,1,...,12NNNnknknkNNNnnnNNNNnknkNNnnNNknkNNnXkDFTxnxnWxnWxnWNxnWxnWNxnxnWWkN因为21NNW，则有2(1)NkkNW，所以：/210[][[](1)[]],0,1,...,12NknkNnNXkxnxnWkN按k的奇偶来讨论，k为偶数时：/2120[2][[][]],0,1,...,12NrnNnNXrxnxnWkNk为奇数时：/21(21)0[21][[][]],0,1,...,12NrnNnNXrxnxnWkN前面已经推导过2/2rkrkNNWW，所以上面的两个等式可以写为：/21/20[2][[][]],0,1,...,/212NrnNnNXrxnxnWrN/21/20[21]{[[][]]},0,1,...,/212NnnrNNnNXrxnxnWWrN通过上面的推导，[]Xk的偶数点值[2]Xr和奇数点值[21]Xr分别可以由组合而成的N/2点的序列来求得，其中偶数点值[2]Xr为输入x[n]的前半段和后半段之和序列的N/2点DFT值，奇数点值[21]Xr为输入x[n]的前半段和后半段之差再与nNW相乘序列的N/2点DFT值。武汉理工大学《信号分析与处理》课程设计4令1[][][]2Nxnxnxn，2[][[][]]2nNNxnxnxnW，则有：/21/211/22/200[2][],[21][],0,1,...,12NNrnrnNNnnNXrxnWXrxnWr这样，也可以用两个N/2点DFT来组合成一个N点DFT，组合过程如图2-1-1所示：[]xn[][]2Nxnxn[]2Nxn-1nNW[[][]]2nNNxnxnW图2-1-1由于2LN，N/2还是偶数，所以可以把N/2点的DFT奇、偶分组，这样就把一个N/2点的DFT分解成两个N/4点的DFT。当N=8时，采用DIF-FFT运算二次分解的运算流图如图2-1-2所示：图2-1-28点DIF-FFT流程图武汉理工大学《信号分析与处理》课程设计52.5DIF-FFT的特点2.5.1原位运算在DIF-FFT蝶形图中，取第m级且两输入节点分别在第k，j行的蝶形为例，讨论DIF-FFT的原位运算规律。由图可得蝶形运算的关系式可表示为()mXk=11()()mmXkXj，()mXj=[11()()mmXkXj]rNW。有上式可得的m-1级的第k行与第j行的输出1()mXk，1()mXj在运算流图中的作用是，用来计算第m级的第k行和第j行的输出()mXk，()mXj，这样当计算完()mXk，()mXj后，1()mXk，1()mXj在运算流图中就不在起作用，因此可以采用原位运算，把()mXk，()mXj直接存入原来存放1()mXk，1()mXj的存储单元。同理可以把第m级蝶形的N个输出值直接存放在第m-1级蝶形输出的N个存储单元中，这样从第一级的输入x（n）开始到最后一级的输出X(k)，只需要N个存储单元。2.5.2蝶形运算两节点之间的“距离”第一级蝶形每个蝶形运算量节点的“距离”为4，第二级每个蝶形运算另节点的“距离”为2，第三级蝶形每个蝶形运算量节点的“距离”为1。依次类推：对于N等于2的L次方的DIF-FFT，可以得到第M级蝶形每个蝶形运算量节点的“距离”为2的L-M次方。2.5.3旋转因子的变化规律以8点的FFT为例，第一级蝶形，r=0,1，2；第二级蝶形，r=0,1；第三级的蝶形，r=0。依次类推，对于M级蝶形，旋转因子的指数为r=12MJ，J=0,1,2,3，……，121M这样就可以算出每一级的旋转因子。对于M级的任一蝶形运算所对应的旋转因子的指数，可以如下方法得到：1将待求的蝶形输入节点中上面节点的行标号值k写成L位二进制数；2将此二进制数乘以2的M-1次方，即将L位二进制数左移M-1位，右边的空位补零，然后从低位到高位截取L位，即所得指数r所对应的二进制数。武汉理工大学《信号分析与处理》课程设计63DIF程序编写与验证3.1DIF程序编写FFT程序包括变址（倒位序）和L级递推计算（N=L2，L为正整数）两部分。3.1.1变址DIF-FFT是输出倒位序的变址处理，设x(i)表示存放自然顺序输入数据的内存单元，x（j）表示存放倒位序序数的内存单元，I、J=0,1,…,N-1,当I=J时，不用变址；当IJ时，需要变址；但是当ij时，进行变址在先，故在IJ时，就不需要再变址了，否则变址两次等于不变。其中本程序使用的“反向进位加法”。也可以用bin2dec函数可以实现倒位序。3.1.2L级递推计算整个L级递推过程由三个嵌套循环构成。外层的一个循环控制L（L=N2log）级的顺序运算；内层的两个循环控制同一级（M相同）各蝶形结的运算，其中最内层循环控制同一种（即rNW中的r相同）蝶形结的运算。其循环变量为I，I用来控制同一种蝶形结运算。其步进值为蝶形结的间距值LE=M2，同一种蝶形结中参加运算的两节点的间距为LE1=2ML点。第二层循环，其循环变量J用来控制计算不同种（系数不同）的碟形结，J的步进值为1。也可以看出，最内层循环完成每级的蝶形结运算，第二层循环则完成系数kNW的运算。最外层循环，用循环变量M来控制运算的级数，M为1到L，步进值为1，当M改变时，则LE1，LE和系数U都会改变。3.2程序编写流图可按如下程序流图编写程序：武汉理工大学《信号分析与处理》课程设计7开始原序列补零至L2计算序列级数M=N2logm=0否mM-1是本级中蝶形图个数group=m2每个蝶形图中元素个数uint=8/group旋转因子Wn=exp（-j*2*pi/group）利用/2120[2][[][]],0,1,...,12NrnNnNXrxnxnWkN/21(21)0[21][[][]],0,1,...,12NrnNnNXrxnxnWkN组内计算计算完毕武汉理工大学《信号分析与处理》课程设计8MATLAB实现的代码：1、在M序列中编写DIF-FFT函数：function[Xk]=DIF_FFT(xn,N);%本程序对输入序列实现DIF-FFT基2算法，%点数取小于等于长度的2的幂次n=log2(2^nextpow2(length(xn)));%求的x长度对应的2的最低幂次nN=2^n;iflength(xn)Nxn=[xn,zeros(1,N-length(xn))];%若长度不是2的幂，补0到2的整数幂end%蝶形运算开始M=log2(N);%“级”的数量form=0:M-1%“级”循环开始Num_of_Group=2^m;%每一级中组的个数Interval_of_Group=N/2^m;%每一级中组与组之间的间距Interval_of_Unit=N/2^(m+1);%每一组中相关运算单元之间的间距Cycle_Count=Interval_of_Unit-1;%每一组内的循环次数Wn=exp(-j*2*pi/Interval_of_Group);%旋转因子forg=1:Num_of_Group%“组”循环开始Interval_1=(g-1)*Interval_of_Group;%第g组中蝶形运算变量1的偏移量Interval_2=(g-1)*Interval_of_Group+Interval_of_Unit;%第g组中蝶形运算变量2的偏移量forr=0:Cycle_Count;%“组内”循环开始k=r+1;%“组内”序列的下标xn(k+Interval_1)=xn(k+Interval_1)+xn(k+Interval_2);%第m级，第g组的蝶形运算式1xn(k+Interval_2)=[xn(k+Interval_1)-xn(k+Interval_2)-xn(k+Interval_2)]*Wn^r;%第m级，第g组的蝶形运算式2武汉理工大学《信号分析与处理》课程设计9endendend%序列排序开始n1=fliplr(dec2bin([0:N-1]));%码位倒置步骤1：将码位转换为二进制，再进行倒序n2=[bin2dec(n1)];%