FFT快速傅里叶变换(蝶形算法)详解-ppt课课件.ppt

1、n直接计算直接计算DFT的问题及改进的途径的问题及改进的途径 2ppt课件课件 3ppt课件课件 设复序列设复序列x(n) 长度为长度为N点，其点，其DFT为为10( )( )NnkNnX kx n Wk=0，N-1 （1）计算一个）计算一个X(k) 值的运算量值的运算量复数乘法次数：复数乘法次数： N复数加法次数：复数加法次数： N14ppt课件课件（2）计算全部）计算全部N个个X(k) 值的运算量值的运算量复数乘法次数：复数乘法次数： N2复数加法次数：复数加法次数： N(N1)（3）对应的实数运算量）对应的实数运算量1100( )( )Re ( )Im ( )ReImNNnknknkNN

2、NnnX kx n Wx njx nWjW10Re ( ) ReIm ( ) ImNnknkNNnx nWx nWRe ( ) ImIm ( ) RenknkNNjx nWx nW5ppt课件课件一次复数乘法一次复数乘法： 4次实数乘法次实数乘法 2次实数加法次实数加法 一个一个X(k) ：4N次实数乘法次实数乘法 2N+2(N-1)= 2(2N-1)次实数加法次实数加法 所以所以 整个整个N点点DFT运算共需要：运算共需要：N2(2N-1)= 2N(2N-1)实数乘法次数：实数乘法次数：4 N2实数加法次数：实数加法次数：6ppt课件课件N点点DFT的复数乘法次数举例的复数乘法次数举例NN2

3、NN22464404941612816384864256 65 536 16256512 262 144 3210281024 1 048 576 结论结论：当：当N很大时，其运算量很大，对实时性很强的信号很大时，其运算量很大，对实时性很强的信号处理来说，要求计算速度快，因此需要改进处理来说，要求计算速度快，因此需要改进DFT的计算的计算方法，以大大减少运算次数。方法，以大大减少运算次数。 7ppt课件课件 nkNW主要原理是利用系数主要原理是利用系数 的以下特性对的以下特性对DFT进行分解：进行分解： nkNW（1）对称性）对称性 ()nkNW()k N nNW（2）周期性）周期性 ()()

4、n N kn k NnkNNNWWW（3）可约性）可约性 mnknkmNNWW/nknk mNN mWW另外，另外，12/NNWkNNkNWW)2/(8ppt课件课件9ppt课件课件1(2 )( )xrx r设设N2L，将，将x(n)按按 n 的奇偶分为两组：的奇偶分为两组： 2(21)( )xrx rr =0，1， 12N10( ) ( )( )NnkNnX kDFT x nx n W则则1010)()(NnnnkNNnnnkNWnxWnx为奇数为偶数10ppt课件课件120)12(1202) 12()2(NrkrNNrrkNWrxWrx1202212021)()(NrrkNkNNrrkNW

5、rxWWrx)()(21kXWkXkN1010)()(NnnnkNNnnnkNWnxWnx为奇数为偶数式中，式中，X1(k)和和X2(k)分别是分别是x1(n)和和x2(n)的的N/2的的DFT。另外，式中另外，式中k的取值范围是：的取值范围是：0，1， ，N/21 。11ppt课件课件因此，因此， 只能计算出只能计算出X(k)的前一半值。的前一半值。12( )( )( )kNX kX kW Xk后一半后一半X(k) 值，值， N/2 ， N/2 1， ，N ？rkNW2(2)2r NkNW利用利用可得到可得到 1()2NXk2 1(2)120( )Nr NkNrx r W2 1120( )N

6、rkNrx r W)(1kX同理可得同理可得22()( )2NXkXk12ppt课件课件考虑到考虑到 kNkNNNkNNWWWW2)2(因此可得后半部分因此可得后半部分X(k) )2()2()2(221NkXWNkXNkXNkN12( )( )( )kNX kX kW Xk及前半部分及前半部分X(k) )()(21kXWkXkNk=0，1， ，N/21k=0，1， ，N/2113ppt课件课件12( )( )( )kNX kX kW Xk12( )( )( )kNX kX kW Xk蝶形运算式蝶形运算式蝶形运算信蝶形运算信号流图符号号流图符号 因此，只要求出因此，只要求出2个个N/2点的点的D

7、FT，即，即X1(k)和和X2(k)，再，再经过蝶形运算就可求出全部经过蝶形运算就可求出全部X(k)的值，运算量大大减少。的值，运算量大大减少。14ppt课件课件0NW1NW2NW3NW以以N=8为例，为例，分解为分解为2个个4点点的的DFT，然后，然后做做8/2=4次蝶形次蝶形运算即可求出运算即可求出所有所有8点点X(k)的的值。值。15ppt课件课件复数乘法次数： N2复数加法次数： N(N1)复数乘法次数： 2*(N/2)2+N/2=N2/2+N/2复数加法次数： 2*(N/2)(N/21)+2*N/2=N2/2nN点 16ppt课件课件 由于N2L，因而N/2仍是偶数 ，可以进一步把每

8、个N/2点子序列再按其奇偶部分分解为两个N/4点的子序列。 以N/2点序列x1(r)为例 1314(2 )( )0,1,1(21)( )4xlx lNlxlx l则有 rkNNrWrxkX212011)()(klNNllkNNlWlxWlx)12(21401221401) 12()2(lkNNlkNlkNNlWlxWWlx41404241403)()()()(42/3kXWkXkNk=0,1, 14N17ppt课件课件且且13/24( )( )4kNNXkXkWXkk=0,1, 14N由此可见，一个由此可见，一个N/2点点DFT可分解成两个可分解成两个N/4点点DFT。 同理，也可对同理，也可

9、对x2(n)进行同样的分解，求出进行同样的分解，求出X2(k)。 18ppt课件课件19ppt课件课件13/40( )lkNlx l W02(0)(4)xW x0(0)(4)NxW x 对此例对此例N=8，最后剩下的是，最后剩下的是4个个N/4= 2点的点的DFT，2点点DFT也可以由蝶形运算来完成。以也可以由蝶形运算来完成。以X3(k)为例。为例。 /4 133/40( )( )NlkNlXkx l Wk=0, 1即即03323(0)(0)(1)XxW x13323(1)(0)(1)XxW x12(0)(4)xW x0(0)(4)NxW x这说明，这说明，N=2M的的DFT可全部由蝶形运算来

10、完成。可全部由蝶形运算来完成。20ppt课件课件N=8按时间抽取法按时间抽取法FFT信号流图信号流图 21ppt课件课件由按时间抽取法FFT的信号流图可知，当N=2L时，共有 级蝶形运算；每级都由 个蝶形运算组成，而每个蝶形有 次复乘、 次复加，因此每级运算都需 次复乘和 次复加。 LN/2 N/2 12N22ppt课件课件这样这样 级运算总共需要：级运算总共需要： L复数乘法： NNLN2log22复数加法： NNLN2log直接直接DFT算法运算量算法运算量 复数乘法： 复数加法： N2N(N1)直接计算直接计算DFT与与FFT算法的计算量之比为算法的计算量之比为MNNNNNM222log

11、2log223ppt课件课件NN2计算量之比M NN2计算量之比M 2414.012816 38444836.641644.025665 5361 02464.0864125.4512262 1442 304113.816256328.010241 048 5765 120204.83210288012.820484 194 30411 264372.464404919221.4NN2log2NN2log224ppt课件课件n序列的逆序排列n同址运算（原位运算）n蝶形运算两节点间的距离n 的确定rNW25ppt课件课件)(01221)()(BINMMDECnnnnnn 由于由于 x(n) 被反

12、复地按奇、偶分组，所以流图输被反复地按奇、偶分组，所以流图输入端的入端的排列不再是顺序的，但仍有规律可循：排列不再是顺序的，但仍有规律可循： 因为因为 N=2M ， 对于任意对于任意 n（0n N-1)，可以用，可以用M个个二进制码表示为：二进制码表示为：10,01221nnnnnMM n 反复按奇、偶分解时，即按二进制码的反复按奇、偶分解时，即按二进制码的“0” “1” 分解。分解。n序列的逆序排列26ppt课件课件27ppt课件课件自然顺序自然顺序 n二进制数二进制数倒位序二进制数倒位序二进制数倒位序顺序数倒位序顺序数000000001001100420100102301111064100

13、0011510110156110011371111117 n28ppt课件课件29ppt课件课件 某一列任何两个节点某一列任何两个节点k 和和j 的节点变量进行蝶形运算的节点变量进行蝶形运算后，得到结果为下一列后，得到结果为下一列k、j两节点的节点变量，而和其他两节点的节点变量，而和其他节点变量无关。这种原位运算结构可以节省存储单元，节点变量无关。这种原位运算结构可以节省存储单元，降低设备成本。降低设备成本。)(kX)2(NkX)(1kX)(2kXkNW运算前运算前运算后运算后)2(NkA)(kA)2(NkA)(kA例例n同址运算（原位运算）30ppt课件课件31ppt课件课件以以N=8为例：

14、为例：第一级蝶形，距离为：第一级蝶形，距离为：第二级蝶形，距离为：第二级蝶形，距离为：第三级蝶形，距离为：第三级蝶形，距离为：规律规律：对于共：对于共L级的蝶形而言，其级的蝶形而言，其m级蝶形运算的节级蝶形运算的节 点间的距离为点间的距离为12412mn蝶形运算两节点间的距离 32ppt课件课件 rNW以N=8为例：0,10224/jWWWWmjjNrNm时，1 , 0,2422/jWWWWmjjjNrNm时，3 , 2 , 1 , 0,382jWWWWmjjjNrNm时，级：第LNM,212 , 2 , 1 , 0,12LjrNjWWLMLMLMLN2222LMLMMLMLjNjNjjNjj

15、NrNWeeWW222222rNWn 的确定 33ppt课件课件n算法原理算法原理 再把输出再把输出X(k)按按k的奇偶分组的奇偶分组先把输入按先把输入按n的顺序分成前后两半的顺序分成前后两半设序列长度为设序列长度为N=2L，L为整数为整数 前半子序列前半子序列x(n) 后半子序列后半子序列)2(Nnx 0n12N0n12N34ppt课件课件10( )( )NnkNnX kx n W由由DFT定义得定义得/2 110/2( )( )NNnknkNNnn Nx n Wx n W12/0)2(12/0)2()(NnkNnNNnnkNWNnxWnx12/02)2()(NnnkNkNNWWNnxnxk

16、=0，1， ，N35ppt课件课件由于由于 1222jNNjNNeeW/2 120( )( )()2NNknkNNnNX kx nx nWW所以所以kkNNW) 1(2则则 12/0)2() 1()()(NnnkNkWNnxnxkXk=0，1， ，N36ppt课件课件然后按然后按k的奇偶可将的奇偶可将X(k)分为两部分分为两部分 221krkrr=0，1， ，12N则式则式 12/0)2() 1()()(NnnkNkWNnxnxkX可转化为可转化为 nrNNnWNnxnxrX212/0)2()()2(12/02/)2()(NnnrNWNnxnx)12(12/0)2()() 12(rnNNnWN

17、nxnxrXnrNnNNnWWNnxnx212/0)2()(37ppt课件课件/2 1/20(2 )( )()2NnrNnNXrx nx nW令令 nNWNnxnxnxNnxnxnx)2()()()2()()(21n=0，1， ，12N代入代入 /2 120(21) ( )()2NnnrNNnNXrx nx nWWnrNNnnrNNnWnxrWnxr2120221201)() 12()()2(r=0，1， ，12N可得可得为为2个个N/2点的点的DFT，合起来正好是，合起来正好是N点点X(k)的值。的值。38ppt课件课件nNWNnxnxnxNnxnxnx)2()()()2()()(21将将称

18、为蝶形运算称为蝶形运算与时间抽选基与时间抽选基2FFT算法中的蝶形运算符号略有不同。算法中的蝶形运算符号略有不同。39ppt课件课件例例 按频率抽取，将按频率抽取，将N点点DFT分解为两个分解为两个N/2点点DFT的组合的组合(N=8) 40ppt课件课件 与时间抽取法的推导过程一样，由于与时间抽取法的推导过程一样，由于 N=2L，N/2仍然是仍然是一个偶数，因而可以将每个一个偶数，因而可以将每个N/2点点DFT的输出再分解为偶数组的输出再分解为偶数组与奇数组，这就将与奇数组，这就将N/2点点DFT进一步分解为两个进一步分解为两个N/4点点DFT。 N=841ppt课件课件n频率抽取法输入是自

19、然顺序，输出是倒位序的；时间抽取法正好相反。n频率抽取法的基本蝶形与时间抽取法的基本蝶形有所不同。n频率抽取法运算量与时间抽取法相同。n频率抽取法与时间抽取法的基本蝶形是互为转置的。 42ppt课件课件MN 2IDFT公式公式 10)(1NknkNWkXNkXIDFTnxDFT公式公式 nkNNnWnxnxDFTkX10)()()(比较可以看出，比较可以看出， nkNWnkNWMN211IDFT多出多出M个个1/2可分解到可分解到M级蝶形运算中。级蝶形运算中。43ppt课件课件44ppt课件课件10)(1)(NknkNWkXNkXIDFT10)(1NknkNWkXN10)(1NknkNWkXN

20、1( )( )IFFT X kFFT XkN( )X k 求共轭( )XkFFT 求( )FFT Xk( )FFT XkN 除以( )x n 求共轭45ppt课件课件n用用FFT进行谱分析的进行谱分析的Matlab实现实现n用用CZT进行谱分析的进行谱分析的Matlab实现实现n在在Matlab中使用的线性调频中使用的线性调频z变换函数为变换函数为czt，其调用格式为其调用格式为nX= czt(x, M, W, A)n其中，其中，x是待变换的时域信号是待变换的时域信号x(n)，其长度为，其长度为N，M是变换的长度，是变换的长度，W确定变换的步长，确定变换的步长，A确定变换确定变换的起点。若的起

21、点。若M= N，A= 1，则，则CZT变成变成DFT。46ppt课件课件例例5.1 设模拟信号设模拟信号 ，以，以 t= 0.01n (n=0: N-1) 进行取样，试用进行取样，试用fft函数对其做频谱分析。函数对其做频谱分析。N分别分别为：为：(1) N=45；(2) N=50；(3) N=55；(2) N=60。 ( )2sin(4)5cos(8)x ttt程序清单如下程序清单如下 %计算计算N=45的的FFT并绘出其幅频曲线并绘出其幅频曲线N=45;n=0:N-1;t=0.01*n;q=n*2*pi/N;x=2*sin(4*pi*t)+5*cos(8*pi*t);y=fft(x,N);

22、figure(1)subplot(2,2,1)plot(q,abs(y)title(FFT N=45)47ppt课件课件%计算计算N=50的的FFT并绘出其幅频曲线并绘出其幅频曲线N=50;n=0:N-1;t=0.01*n;q=n*2*pi/N;x=2*sin(4*pi*t)+5*cos(8*pi*t);y=fft(x,N);figure(1)subplot(2,2,2)plot(q,abs(y)title(FFT N=50)48ppt课件课件%计算计算N=55的的FFT并绘出其幅频曲线并绘出其幅频曲线N=55;n=0:N-1;t=0.01*n;q=n*2*pi/N;x=2*sin(4*pi*

23、t)+5*cos(8*pi*t);y=fft(x,N);figure(1)subplot(2,2,3)plot(q,abs(y)title(FFT N=55)49ppt课件课件%计算计算N=60的的FFT并绘出其幅频曲线并绘出其幅频曲线N=60;n=0:N-1;t=0.01*n;q=n*2*pi/N;x=2*sin(4*pi*t)+5*cos(8*pi*t);y=fft(x,N);figure(1)subplot(2,2,4)plot(q,abs(y)title(FFT N=60)50ppt课件课件02468050100150FFT N=4502468050100150FFT N=500246

24、8050100150FFT N=5502468050100150FFT N=60从图中可以看出，这几种情况下均有较好的精度。从图中可以看出，这几种情况下均有较好的精度。 51ppt课件课件分析：由分析：由t=0.01n进行取样可得，采样频率进行取样可得，采样频率fs=100Hz。而。而连续信号的最高模拟角频率为连续信号的最高模拟角频率为8 ，由，由2 f可得，可得，最高频率为最高频率为8 /2 =4Hz。因此，满足采样定理的要求。因此，满足采样定理的要求。 采样序列为采样序列为( )2cos(4)5cos(8)x nTnTn48( )2cos5cos100100 x nnn即即为周期序列，周期为周期序列，周期N=50。将程序中将程序中plot改为改为stem函数，则可以更清楚地看出频谱。函数，则可以更清楚地看出频谱。52ppt课件课件02468050100150FFT N=4502468050100150FFT N=5002468050100150FFT N=5502468050100150FFT N=6053ppt课件课件