快速傅里叶变换

单击此处编辑母版标题样式,单击此处编辑母版文本样式,第二级,第三级,第四级,第五级,*,单击此处编辑母版标题样式,单击此处编辑母版文本样式,第二级,第三级,第四级,第五级,*,单击此处编辑母版标题样式,*,单击此处编辑母版文本样式,第二级,第三级,第四级,第五级,单击此处编辑母版标题样式,单击此处编辑母版文本样式,第二级,第三级,第四级,第五级,单击此处编辑母版标题样式,单击此处编辑母版文本样式,第二级,第三级,第四级,第五级,*,单击此处编辑母版标题样式,单击此处编辑母版文本样式,第二级,第三级,第四级,第五级,*,单击此处编辑母版标题样式,*,*,单击此处编辑母版文本样式,第二级,第三级,第四级,第五级,单击此处编辑母版标题样式,单击此处编辑母版文本样式,第二级,第三级,第四级,第五级,*,*,单击此处编辑母版标题样式,单击此处编辑母版文本样式,第二级,第三级,第四级,第五级,*,*,单击此处编辑母版标题样式,单击此处编辑母版文本样式,第二级,第三级,第四级,第五级,*,*,快速傅里叶变换,本章目录,直接计算DFT的问题与改进的途径,按,时间抽取,的,基,2-FFT,算法,按,频率抽取,的基,2-FFT,算,法,快速,傅里叶逆变换,(IFFT),算法,Matlab,实现,2,4.1,引言,DFT在实际应用中很重要: 可以计算信号的频谱、功率谱和线性卷积等。,直接按DFT变换进展计算，当序列长度N很大时，计算量非常大，所需时间会很长,实时处理难以实现。,1965年，图基和库利发表了?机器计算快速傅立叶级数的一种算法?论文后，很快形成了快速计算DFT的计算机算法FFT。Fast Fourier Transform),FFT并不是一种与DFT不同的变换，而是DFT的一种快速计算的算法。,3,4.2 直接计算DFT的问题与改进的途径,DFT,的运算量,设复序列,x,(,n,),长度为,N,点，其,DFT,为,k,=0,，,，,N,-1,1计算一个X(k) 值的运算量,复数乘法次数：,N,复数加法次数：,N,1,2计算全部N个X(k) 值的运算量,复数乘法次数：,N,2,复数加法次数：,N,(,N,1),4,DFT,运算量的结论,N,点,DFT,的复数乘法次数举例,N,N2,N,N2,2,4,64,4049,4,16,128,16384,8,64,256,65 536,16,256,512,262 144,32,1024,1024,1 048 576,结论：当N很大时，其运算量很大，对实时性很强的信号处理来说，要求计算速度快，因此需要改进DFT的计算方法，以大大减少运算次数。,5,4.2.2,减少运算工作量的途径,主要原理是利用系数 的以下特性对DFT进展分解：,2对称性,1周期性,3可约性,另外，,6,4.3,按时间抽取的基,2-FFT,算法,算法原理,按时间抽取基-2FFT算法与直接计算DFT运算量的比较,按时间抽取的FFT算法的特点,按时间抽取FFT算法的其它形式流程图,DIT-FFT(Decimation-In-Time),7,4.3.1,算法原理,设,N,2,M,，将,x,(,n,),按,n,的奇偶分为两组：,r,=0,，,1,，,，,那么,8,式中，,X,1,(,k,),和,X,2,(,k,),分别是,x,1,(,n,),和,x,2,(,n,),的,N,/2,的,DFT,。,另外，式中,k,的取值范围是：,0,，,1,，,，,N/2,1,。,9,因此， 只能计算出,X,(k),的前一半值。,后一半,X,(k),值，,N/2,，,N/2,1,，,，,N-1,？,同理可得,10,因此可得后半局部X(k),由前半局部X(k),k=,0,，,1,，,，,N/2,1,k,=0,，,1,，,，,N/2,1,11,k=,0,，,1,，,，,N/2,1,结论：,因此，只要求出,2,个,N/2,点的,DFT,，即,X,1,(k),和,X,2,(k),，再经过这种运算就可求出全部,X,(k),的值,12,新概念：蝶形运算,蝶形运算式,蝶形运算,信号流图符号,X,1,(,k,),X,2,(,k,),蝶形运算的运算量：每次蝶形含一次复数乘和两,次复数加,13,以,8,点为例第一次按奇偶分解,以,N=8,为例，,分解为,2,个,4,点的,DFT,，然后做,8/2=4,次蝶形运算即可求出所有,8,点,X,(k),的值。,分解一次后所需的运算量,2,个,N/2,的,DFT,N/2,蝶形,14,运算量比较,复数乘法次数：,N,2,复数加法次数：,N,(,N,1),复数乘法次数：,2*(,N/,2,),2,+,N,/2=,N,2,/2+,N,/2,复数加法次数：,2*(,N/2),(,N/2,1)+2*,N,/2=,N,2,/2,N,点,DFT,的运算量,分解一次后所需的运算量,2,个,N/2,的,DFT,N/2,蝶形：,通过一次分解后，运算工作量减少了差不多一半。,每次蝶形含一次复数乘和两次复数加,15,进一步按奇偶分解,由于N2M，因而N/2仍是偶数 ，可以进一步把每个N/2点,子序列再按其奇偶局部分解为两个N/4点的子序列。,以,N,/2,点序列,x,1,(,r,),为例,那么有,k,=0,1,16,且,k,=0,1,由此可见，一个,N,/2,点,DFT,可分解成两个,N,/4,点,DFT,。,同理，也可对x2(n)进展同样的分解，求出X2(k)。,17,以,8,点为例第二次按奇偶分解,概念：信号流图,18,算法原理,对此例,N,=8,，最后剩下的是,4,个,N,/4= 2,点的,DFT,，,2,点,DFT,也可以由蝶形运算来完成。,N=2,这说明，N=2M的DFT可全部由蝶形运算来完成。,因此：N必须是2的整数次幂. “基2,蝶形运算,19,以,8,点为例第三次按奇偶分解,N=8,按时间抽取法,FFT,信号流图,x(n),抽取,X(k),顺序,DIT-FFT(Decimation-In-Time),FFT,算法,20,4.3.2 按时间抽取基2-FFT算法与直接计算DFT运算量的比较,由按时间抽取法,FFT,的信号流图可知，当,N,=2,M,时，共有,级,蝶形运算；每级都由,个蝶形运算组成，而每个蝶形有,次复乘、,次复加，因此每级运算都需,次复乘和,次复加。,M,N,/2,N,/2,1,2,N,21,FFT,算法中，这样,级运算总共需要：,M,复数乘法：,复数加法：,直接,DFT,算法运算量,复数乘法：,复数加法：,N,2,N,(,N,1),直接计算,DFT,与,FFT,算法的计算量之比为,M,N=2,M,22,FFT算法与直接DFT算法运算量的比较,N,N,2,计算量之比,M,N,N,2,计算量之比,M,2,4,1,4.0,128,16 384,448,36.6,4,16,4,4.0,256,65 536,1 024,64.0,8,64,12,5.4,512,262 144,2 304,113.8,16,256,32,8.0,1024,1 048 576,5 120,204.8,32,1028,80,12.8,2048,4 194 304,11 264,372.4,64,4049,192,21.4,23,24,N=1024; M=80;,x=1:M,zeros(1,N-M);,t=cputime;,y1=fft(x,N);,time_fft=cputime-t;,t1=cputime;,y2=dft(x,N);,time_dft=cputime-t1;,FFT算法与直接DFT算法运算量的比较,25,L=1,L=2,L=3,倒序,4.3.3,按时间抽取的,FFT,算法的特点,26,同址运算原位运算,序列的逆序排列,蝶形运算两节点间的距离,一样旋转因子节点间的距离,(旋转因子确实定,27,同址运算原位运算,某一列任何两个节点k 和j 的节点变量进展蝶形运算后，得到结果为下一列k、j两节点的节点变量，而和其他节点变量无关。这种原位运算构造可以节省存储单元，降低设备本钱。,运算前,运算后,例,同址运算原位运算,28,观察原位运算规律,29,序列的逆序排列,由于,x,(,n,),被反复地按奇、偶分组，所以流图输,入端的,排列不再是顺序的，但仍有规律可循：,因为,N,=2,M,，,对于任意 n0n N-1)，可以用M个二进制码表示为：,n 反复按奇、偶分解时，即按二进制码的“0 “1 分解。,序列的逆序排列,30,倒位序的树状图N=8,31,码位的倒位序,(N=8),自然顺序,n,二进制数,倒位序二进制数,倒位序顺序数,0,000,000,0,1,001,100,4,2,010,010,2,3,011,110,6,4,100,001,1,5,101,101,5,6,110,011,3,7,111,111,7,I,J,规律：最左端加,1,，向右进位,倒序,32,倒位序的变址处理N=8,I=J,不换,IJ,不换,IX= czt(x, M, W, A),其中，x是待变换的时域信号x(n)，其长度为N，M是变换长度，W确定变换的步长，A确定变换的起点。假设M= N，A= 1，那么CZT变成DFT。,75,4.8.1,用,FFT,计算实信号的快速线性卷积,设一离散线性移不变系统的冲激响应为,其输入信号为,.,其输出为,.,并且,的长度为,N,点,的长度为,M,点,计算其线性卷积。,76,用,FFT,算 的步骤：,77,FFT,FFT,IFFT,x,x(n),h(n),y(n),X(k),H(k),Y(k),流程图,程序参考第三章幻灯片,78,4.8.2 用FFT进展谱分析的Matlab实现,例4.8.2 设模拟信号 ，以 tn (n=0: N-1) 进展取样，试用fft函数对其做频谱分析。N分别为：(1) N=10；(2) N=50；(3) N=100；(2) N=1024。,程序清单如下,%,计算,N,=10,的,FFT,并绘出其幅频曲线,N=10;n=0:N-1;t=0.01*n;,q=n*2*pi/N;,x=2*sin(20*pi*t)+5*cos(40*pi*t);,y=fft(x,N);,figure(1),subplot(2,2,1),plot(q,abs(y),title(FFT N=10),79,例,4.8.2,程序清单,%,计算,N,=50,的,FFT,并绘出其幅频曲线,N=50;n=0:N-1;t=0.01*n;,q=n*2*pi/N;,x=2*sin(20*pi*t)+5*cos(40*pi*t);,y=fft(x,N);,figure(1),subplot(2,2,2),plot(q,abs(y),title(FFT N=50),80,%,计算,N,=100,的,FFT,并绘出其幅频曲线,N=100;n=0:N-1;t=0.01*n;,q=n*2*pi/N;,x=2*sin(20*pi*t)+5*cos(40*pi*t);,y=fft(x,N);,figure(1),subplot(2,2,3),plot(q,abs(y),title(FFT N=100),81,%,计算,N,=1024,的,FFT,并绘出其幅频曲线,N=1024;n=0:N-1;t=0.01*n;,q=n*2*pi/N;,x=2*sin(20*pi*t)+5*cos(40*pi*t);,y=fft(x,N);,figure(1),subplot(2,2,4),plot(q,abs(y),title(FFT N=1024),82,谢谢！,83,谢谢观赏,