第4章快速傅里叶变换课件

,单击此处编辑母版标题样式,单击此处编辑母版文本样式,第二级,第三级,第四级,第五级,#,第,4,章 快速傅里叶变换,2,第,4,章,快速傅里叶变换,4.1,引言,4.2,直接计算,DFT,的问题及改进的途径,4.3,按时间抽取（,DIT,）的基,-2FFT,算法,4.4,按频率抽取（,DIF,）的基,-2FFT,算法,4.5 IDFT,的高效算法,4.6,实序列的,FFT,算法,4.7,线性调频,Z,变换（,CZT,）,3,4.1,引言,有限长序列通过离散傅里叶变换（,DFT,）将其频域离,散化成有限长序列，但其计算量太大，很难实时处理，因此,引出了快速傅里叶变换（,FFT,）。,FFT,并不是一种新的变换形式，它只是,DFT,的一种,快,速算法,，并且根据对序列分解与选取方法的不同产生了多种,算法。,FFT,在离散傅里叶反变换、线性卷积和线性相关等方,面有重要应用。,4,4.2,直接计算,DFT,的问题及改进的途径,4.2.1,直接计算,DFT,的运算量问题,若,x,(,n,),是,N,点序列，实现,x,(,n,),的,DFT,，即求出,N,个,X,(,k,),，,需要,N,2,次复数乘法，,N,(,N,-1),次复数加法。,例如：,x,(,n,),：,N,=1024,，,N,2,=1048576,式中旋转因子,(,Twiddle Factor,),：,一般情况下，信号序列,x(n),及其频谱序列,X(k),都是用复数来表示的，,W,N,当然也是复数。因此，计算,DFT,的一个值,X(k),需要进行,N,次复数乘法（与,1,相乘也包括在内）和,N-1,次复数加法；所以，直接计算,N,点的,DFT,需要进行,N,2,次复数乘法和,N(N-1),复数加法,。,显然，直接计算,N,点的,IDFT,所需的复乘和复加的次数也是这么多。当,N,足够大时，,N,2,N(N-1),因此，,DFT,与,IDFT,的运算次数与,N,2,成正比，随着,N,的增加，运算量将急剧增加，而在实际问题中，,N,往往是较大的，因此有必要对,DFT,与,IDFT,的计算方法予以改进。,8,解决耗时的,乘法问题,是将数字信号处理理论用于实际的,关键问题。,特别是,40,年前，计算机的速度相当慢。因此，很,多学者对解决,DFT,的快速计算问题产生了极大的兴趣。,Cooley J W,Tukey J W.An algorithm for the machine,computation of complex Fourier series.,Mathematics of Computation,1965,pp297301,的正式开端！,DSP,9,4.2.2,改善途径,FFT,的思路：,如何充分利用这些关系,？,对称性,周期性,减少运算量的基本途径,旋转因子的主要性质为：,周期性：,对称性：,可约性：,(4),特殊点：,DFT,和,IDFT,的快速算法的导出主要是根据 因子的特性。,减少运算量的基本途径,FFT,算法核心思想：,利用,DFT,系数 的,特性；合并,DFT,运算中的某些,项；,把长,序列,DFT,转换成短序列,DFT,，从而减少运算量。,把长序列变短序列原因：,DFT,运算量和,N,2,成正比，,N,降低，运算量就大大下降。,12,利用,W,N,nk,的对称性和周期性，将大点数的,DFT,分解成若,干个小点数的,DFT,，,FFT,正是基于这个基本思路发展起来的。,N,点,DFT,N/2,点,DFT,N/4,点,DFT,2,点,DFT,1,个,2,个,4,个,N/2,个,13,以四点,DFT,为例,14,几个乘法？,分为：,按频率抽取（,DIFDecimation-In-Frequency,）法；,按时间抽取（,DITDecimation-In-Time,）法 ：,其中：基,-2 FFT,是最基本的,FFT,算法，它要求,FFT,的点数,N,=2,L,（,L,为正整数）；,如果序列长度不满足这一条件，需要用后补零值点的方法来补齐。,FFT,是在时域将序列逐次分解为长度减半的奇序号子序列和偶序号子序列，用子序列的,DFT,来实现整个序列,DFT,；,问题是,如何分最有效？,16,4.3,按时间抽取（,DIT,）的基,-2FFT,算法,4.3.1,算法原理,设序列,x,(n),长度为,N,2,M,，,M,为整数,(,若不满足，则在序列尾部补零,),。将,x,(,n,),（,n,=0,1,N,-1,）按,n,的奇偶分解成两组,N,/2,点的子序列,x,1,(,r,)=,x,(2,r,),x,2,(,r,)=,x,(2,r,+1),r,=0,1,N,/2-1,（长度为,N,/2,）,则,20,蝶形运算,X,1,(k),、,X,2,(k),都是,N/2,点的,DFT,，它们各自又可分成,N/4,点的,DFT,，如此继续分下去，直至两点,DFT,。两点,DFT,不需要乘法运算：,2,点,DFT,的运算流图,28,由于每一步分解都是按输入序列在时间上的奇偶次序，,分解成两个半长的子序列，所以称,“,按时间抽取算法,”,。,第,一级,L=1,第二,级,L=2,第三,级,L=3,蝶形运算单元,对于,N=2,M,点,FFT,，一共有,M,级（,L=1,2,M,）；每级有,N/2,个蝶,形；每级有,N/2,L,组。,时域奇偶分,频域对半分,29,4.3.2,运算量比较,对,N,=2,M,，共可分,M,次,，即,m,=0,1,M,-1,，每一级有,N/2,个如下的蝶形单元,一个蝶形单元只需,一,次复数乘法，,两,次复数加法。,总共复数乘法：,复数加法：,可以共享,31,直接计算,DFT,与,FFT,算法的计算量之比为,如：,N=1024,：,运算效率大大提高,.N,越大，,优势越,明显,!,32,4.3.3,算法特点,1,）,蝶形运算,2,）,原位运算,两点构成一个蝶形单元，并且这两点不再参与别的蝶形,单元的运算。,所谓原位运算，就是当数据输入到存储器后，每一级运,算的结果仍然贮存在这同一组存储器中，直到最后输出，中,间无需其它存储器。,33,由于每一步分解都是按输入序列在时间上的奇偶次序，,分解成两个半长的子序列，所以称,“,按时间抽取算法,”,。,第,一级,L=1,第二,级,L=2,第三,级,L=3,蝶形运算单元,对于,N=2,M,点,FFT,，一共有,M,级（,L=1,2,M,）；每级有,N/2,个蝶,形；每级有,N/2,L,组。,时域奇偶分,频域对半分,从图可以看出，整个算法流图可以分为四段，（,0,）段为倒序重排，后面,3,段为,3(,log,2,8=3),次迭代运算：首先计算,2,点,DFT,，,然后将,2,点,DFT,的结果组合成,4,点,DFT,，,最后将,4,点,DFT,组合为,8,点,DFT,。,因此，对于,N,点,FFT,，,只需要一列存储,N,个复数的存储器。,开始时，输入序列的,N,个数放于此,N,个存储器内，倒序重排后仍存于这,N,个存储器中，每一次迭代运算后的结果也仍然存于这,N,个存储器中，整个运算完成后，这,N,个存储器中即为所求的频谱序列,X(k),（,k=0,、,1,、,.,、,N-1,）。,这就是所谓的同址计算，这样可以大大节约存储元件。,36,3,）倒位序,输入序列,x,(n),的排列次序不符合自然顺序。此现象是由,于按,n,的奇偶分组进行,DFT,运算而造成的，这种排列方式称,为“倒位序”。,所谓,“,倒位序,”,，是指按二进制表示的数字首尾位置颠,倒，重新按十进制读数。,倒位序：,输入序列,经过,M,级分解后，,其排列顺序变为,x(0)x(4)x(2)x(6)x(1)x(5)x(3)x(7),，服从所谓的“,码位倒置,”的规律，称之为倒位序。,造成倒位序的原因,：,将其按标号的,偶奇,的不断分组，,每次分解总是将,偶序列放在上面,，,把奇序列放在下面,。,首先最低位按,0,、,1,分为偶、奇两组，,接着次低位也按,0,、,1,分组，,依此类推,右图为描述倒位序的树状图,(N=8),可以看出，将输入序号的二进制表示（,n,2,n,1,n,0,）,位置颠倒，得到（,n,0,n,1,n,2,），,就是相应的倒序的二进制序号。因此，输入序列按倒序重排，实际上就是将序号为（,n,2,n,1,n,0,）,的元素与序号为（,n,0,n,1,n,2,）,的元素的位置相互交换。,39,倒位序的实现,在实际运算时，先按自然顺序将输入序列存入存储单元，再通过,变址运算,将自然顺序变换成按,DIT-FFT,算法要求的顺序。,变址的过程可以用程序安排加以实现,-,称为,“整序”或“重排”,（采用码位倒读）,注意,:,(1),当,I=J,时，数据不必调换；,(2),当,IJ,时，必须将原来存放数据x(,I,)送入暂存器R，再将x(,J,)送入x(,I,),R中内容送x(,J,).进行数据对调。,(3)为了避免再次调换已调换过的数据，保证调换只进行一次，否则又变回原状,。,J,I,（或,I1,向内盘旋（内缩）；,W,0,1,向外盘旋（外伸）,表示相邻分析点间的夹角,。,这些,z,平面上的抽样点所沿的周线是一条螺旋线，参数,W,0,控制周线盘旋的倾斜率。若,W,0,大于,1,，则随着,k,的增加，周线向内盘旋，趋向原点；若,W,0,小于,1,，则随着,k,的增加，周线向外盘旋。若,W,0,=1,，,螺旋线实际上是一段圆弧，而如果又有,A,0,=1,，,则这段圆弧是单位圆的一部分。,A,0,和,0,分别为第一个抽样点（,k=0,）,的模和幅角，其余抽样点沿螺旋周线按角度间隔,0,分布。,要计算：,式中,N,为序列,x(n),的长度。由恒等式：,并且令：,，,就可以得到：,2.CZT,的实现,改变变量，将,n,换为,r,，,k,换为,n,，,则有：,n=0,1,M-1,其中，可以看作线性时不变系统,h(n),对输入信号,y(n),的响应。,图,：,CZT,算法的计算过程,75,说明：,76,例：用函数,chirp,产生线性调频信号，并画出波形。,t=0:0.001:2;,%2 secs 1kHz sample rate,y=chirp(t,0,1,150);%Start DC,cross 150Hz at t=1s,plot(t,y),axis(0 0.5-1 1),图,4-20,用函数,chrip,产生的线性调频信号波形,77,3.CZT,的实现步骤,（,1,）形成,h,L,(n),序列,（,2,）,（,3,）,78,(4),(5),(6),(7),4.CZT,的特点,与标准的,FFT,算法相比较，,CZT,算法有以下特点：,(1)输入和输出序列长度不需要相等，即不需要,M=N。,(2)N,与,M,都不需要是2的正整数幂。,(3),z,k,间的间隔可以任意选择，这样就可得到任意的分辨率可达到频域,“,细化,”,的目的；而且当所需间隔不同时，可以分段计算。,(4),z,k,点所沿周线不必须是圆弧。,(5)起始点的选择可以是任意的，因此便于从任意的频率上开始对输入数据进行频谱分析。,(6)当 表示式中,A=1，W=e,-j2/,N,，,并且,N=M,时，,X(z,k,),即为,x(n),的,DFT，,因此利用,CZT,算法也能快速计算,x(n),的,DFT，,而且不要求,N,为2的正整数幂。,图,4-21 Chirp-Z,变换的圆周卷积示意图（教材,P102,）,例4.1,已知序列,x(n)=(0.5),n,0n12，,试计算序列,x(n),在单位圆上的,CZT,并与该序列的,DFT,进行比较。,解：,运行结果如图4.,22,所示，从图中可以看出，此时序列在单位圆上的,CZT,就等于该序列的,DFT。,图,4.,22,序列 在单位圆上的,CZT,及其,DFT,例4.2,假设序列,x(n),由4个频率分别为6,Hz、6.3 Hz、9 Hz,和8,Hz,的正弦序列组合而成，抽样频率为40,Hz，,时域抽样200点。,（1）,用,CZT,计算,DFT；,（2）,直接计算,DFT；,（3）,在510,Hz,的频段范围求,CZT。,(,a)CZT,(,b)DFT,图,4.,23 CZT,的应用,(,c)510Hz CZT,在图,4.,23,中，（,a,）,图是利用,CZT,求取,DFT,，（,b,）,图是直接求,DFT,，,两者结果相同；图中频率为,6,Hz,和,6.3,Hz,的两个正弦信号的频谱不易分辩出。,图（,c,）,是在,5,10,Hz,这个频率范围内求出的,CZT,，,它的起始点不在,0,Hz,处，而且，由于它的分点比较细，所以图中四个正弦信号的频谱都可以分辨出来。,本章小结,1,，本章基本没有引入新的物理概念，而是从数学方法着手寻找,DFT,的快速算法。所谓快速算法，是指根据原始变换定义算法的运算规律及其中某些算子的特殊性质（如 的周期性和对称性等），找出减少乘法和加法运算次数的有效途径，实现原始变换的各种高效算法。,2,，一种 好的快速算法可使变换速度提高几个数量级，例如当 时，,DIT-DIF,算法可以使,N,点的,DFT,的乘法运算量减少到原来的,3,，快速算法很多，而且还在不断研究发展中，较成熟的算法已都有现成的程序。通过对,DIT-DIF,两种算法的学习，不仅要熟悉快速算法的特点、运算效率和适用情况，更要注重研究快速算法的基本思想和实现途径，为进一步学习掌握别的快速算法以及开发新的快速算法奠定基础。,第四章习题,1,3,4(1),6,*,复习思考题与习题（第四章）,P125,