第三章_FFT_2014S

1、1第三章第三章 离散傅里叶变换离散傅里叶变换DFTDFT及快速算法及快速算法nFFTFFTu 基基2 2时域抽取的时域抽取的FFTFFT算法算法 u 基基2 2频域抽取的频域抽取的FFTFFT算法算法nIFFTIFFTu 基基2 2时域抽取的时域抽取的IFFTIFFT算法算法u 基基2 2频域抽取的频域抽取的IFFTIFFT算法算法2直接计算直接计算DFTDFT的运算量分析的运算量分析N N 点有限长序列点有限长序列 x(n)x(n)的的 DFT DFT 变换对的定义为：变换对的定义为：10( )( )NnkNnX kx n W 0,1kN 2jNNWe 101( )( )NnkNkx nX

2、k WN 0,1nN 其中其中假设假设 x(n) 是复序列，是复序列， 同时同时 X(k) 一般也是复数一般也是复数3如如 N=512、1024 和和 8192 时，时，DFT 的乘法运算的乘法运算 5122 = 218 = 262144（26万次）万次） 10242 = 220 = 1048576（105万次）万次） 81922 = 226 = 67108864（6千千7百万次）百万次）对于大对于大 N，在实际中是不能接受的，无法，在实际中是不能接受的，无法“实时实时”应用应用 DFTCooley 与与 Turkey 提出的提出的 FFT 算法，大大减少算法，大大减少了计算次数。如了计算次数

3、。如 N =512 时，时，FFT 的乘法次数约的乘法次数约为为 2000 次，提高了约次，提高了约 128 倍，而且倍，而且简化随简化随 N 的增加而巨增的增加而巨增，因而，用数值方法计算频谱得到因而，用数值方法计算频谱得到实际应用实际应用直接计算直接计算DFTDFT的运算量分析的运算量分析4 nkmnkNmNWW 可可约约性性/nknk mNN mWW / 2(/ 2)01 1kNkNNNNNWWWW 特特殊殊点点：2jnknkNNWe 2jmnkmNe *()() ()nknkN n kn N kNNNNWWWW 对对称称性性NknkNNWW nNnkNNWWW WN Nknkn的性质的

4、性质()() nkNn kn NkNNNWWW周周期期性性5u以四点DFT为例： 直接计算需要直接计算需要：42 = 16 次复数乘次复数乘) 3 () 2() 1 () 0(111111111111) 3 () 2() 1 () 0(14141414xxxxWWWWXXXX) 3() 2() 1 () 0() 3() 2() 1 () 0(94643404644424043424140404040404xxxxWWWWWWWWWWWWWWWWXXXX 1414) 3 () 1 () 2() 0() 3 () 3 () 1 () 2() 0() 2() 3 () 1 () 2() 0() 1

5、() 3 () 1 () 2() 0() 0(WxxxxXxxxxXWxxxxXxxxxX只需要 1 次复数乘)3() 1 ()2()0(111111111111)3()2() 1 ()0(14141414xxxxWWWWXXXX6基本思路基本思路：通过迭代运算，利用通过迭代运算，利用 的周期的周期性和对称性，用低点数的性和对称性，用低点数的DFT完成高点数的完成高点数的DFT计算计算 knNjknNeW)2(DFTDFT的改进途径的改进途径1. 时间抽取（时间抽取（Decimation-in-Time） FFT算法算法将序列将序列x(n) 逐次分解为较短的子序列进行逐次分解为较短的子序列进行

6、DFT计算计算2. 频率抽取（频率抽取（Decimation-in-Frequency）FFT算法算法将离散傅里叶变换系数序列将离散傅里叶变换系数序列X(k)分解为较短的子序列分解为较短的子序列7基基2 2时间抽选法（时间抽选法（R2-DITR2-DIT）1, 1 , 0,)()(10NkWnxkXNnknN设设x(n)序列点数序列点数 N = 2L，L 为整数。若不满足，则补零为整数。若不满足，则补零N 为为 2 的整数幂的的整数幂的 FFT 算法称基算法称基-2 FFT算法算法 将序列将序列 x(n) 按按 n 的的奇偶奇偶分成两组：分成两组：即一组由偶数序号组成，另一组由奇数序号组成即一

7、组由偶数序号组成，另一组由奇数序号组成12( )(2 )0,1,12( )(21)x rxrNrx rxr 偶偶数数组组奇奇数数组组8120) 12(120)2(10) 12()2()()(NrrkNNrrkNNnknNWrxWrxWnxkX11222200( )( )( )( )kNNNkrkkrNNNrra r WWb r WWA kB k 22k rkrNNWW12021202) 12()2(NrkrNkNNrkrNWrxWWrx9分解定理分解定理 、 是是N/2点的点的DFT 而而X(k)中中k的范围为的范围为 ，还应进一步考虑还应进一步考虑 的情况，即的情况，即 121 , 0Nk)

8、(kA)(kB10N12NN2()()2)(2(2)NkNWNX kNA kkNB 10分解定理分解定理nkNW)(222NkrNrkNWW利用利用 的周期性：的周期性：)()()()()2(12021202120)2(222kAWraWWraWrakNANrrkNNrrkNrNrkNrNNN)()2(kBkNBkNkNNNkNNWWWW2)2(2()()()()22(2)NkNkNNNNXA kkA kWkWB kB 121 , 0Nk同理同理又又所以所以11分解定理分解定理( )( )( )()( )( )2kkNN X kA kB kNX kA kB kWW 12, 1 , 0Nk因此：

9、因此：X(k)与与A(k)、B(k)的关系可用的关系可用蝶形信号流图蝶形信号流图表示表示1)(kA)(kBkNW)()(kBWkAkN)()(kBWkAkN整个整个 的计算，可以分解为的计算，可以分解为前、后半部分的运算前、后半部分的运算。而只要。而只要求出前一半，就可以求出整个序列。求出前一半，就可以求出整个序列。)(kX12分解定理分解定理 例例 N=23=8； k=0,1,2,313分解定理分解定理114442/400(4 )(42)NNklkklNNNllx l WWx lW 1202)2()(NrkrNWrxkA 我们可以再按奇、偶数把我们可以再按奇、偶数把A(k), B(k)进一步

10、分解为进一步分解为 两个两个N/4点的点的DFT 令令 则则 ) 12 (2l lr及LN21144(2 )(21)2200( )(4 )(42)NNklklNNllA kx l Wx lW 14, 1 , 0NkN2()kWCD kk14分解定理分解定理1404/221404140)(4/)(2140)(4)24()4()24()4()4(4444NlklNNkNNlklNNllkNkNNllkNWlxWWWlxWlxWWlxNkANNNN14, 1 , 0Nk2( )( )kNkWC kD2N( )( )kCkWkD 142/242jNNjNNeeW15分解定理分解定理 2 点点 DFTx

11、(0)x(2)x(4)x(6) 2 点点 DFTx(6)x(4)x(2)x(0)A(0)A(1)A(2)A(3)C(0)C(1)D(0)D(1)W80W82-1-1由两个由两个N/4点点DFT组合成一个组合成一个N/2点点DFT2N2N( )( )( )()( )( )4kkA kC kD kNA kC kD kWW14, 1 , 0Nk16分解定理分解定理因为因为N N为为2 2的正整数次幂的正整数次幂, ,可以继续分解，直到最后可以继续分解，直到最后分解为分解为2 2点点DFTDFT 2 点点 DFTx(0)x(2)x(4)x(6) 2 点点 DFTx(6)x(4)x(2)x(0)A(0)

12、A(1)A(2)A(3)C(0)C(1)D(0)D(1)W80W82-1-1x(0)x(2)x(4)x(6)x(6)x(4)x(2)x(0)A(0)A(1)A(2)A(3)C(0)C(1)D(0)D(1)W80W82-1-1-1-117分解定理分解定理对于分解到最后的对于分解到最后的2 2点点DFTDFT的信号流图的信号流图蝶形图蝶形图 -1) 0( x) 4( x01NW (0)(0)(4)Cxx(1)(0)(4)Cxx（只需要一次加法，一次减法）（只需要一次加法，一次减法） 18分解定理分解定理N=8时，基时，基2时间抽选法时间抽选法FFT的信号流图的信号流图)0(X)7(X)6(X)5(

13、X)4(X)3(X)2(X)1 (X)0(x)7(x)6(x)5(x)4(x)3(x)2(x)1 (x)0(x)1 (x)2(x)3(x)4(x)5(x)6(x)7(x点24点42点8111111111111108W28W08W28W08W18W38W28W1m2m3m08W08W08W08W19例例 使用基使用基 2 时间抽选法时间抽选法 FFT 流图，计算下列流图，计算下列数据的数据的 8点点 DFT：x=4,-3,2,0,-1,-2,3,112()/j4-320-1-23142-13-30-214-123-3-201355-1-5-11-1355j-5-11j85+j-25-j-4-1+

14、j-6-1-j85+j-25-j-4j26jj245+j+ j2-2+6j5-j+ j2125+j- j2-2-6j5-j- j2j12() /j128Wj 081W28Wj 081W1 1 1 1 1 1 1 1 1 1 1 1 20蝶形计算蝶形计算 任何一个任何一个N=2M的的DFT，都可以通过都可以通过M次分解，最后成为次分解，最后成为2点点的的 DFT来计算。来计算。M次分解构成了从次分解构成了从x(n)到到X(k)的的M级迭代计级迭代计算，算，每级由每级由N/2个蝶形组成个蝶形组成。计算时总是上节点的值和下。计算时总是上节点的值和下节点的值乘以节点的值乘以 然后相加减，形成下一级两个

15、节点值，这然后相加减，形成下一级两个节点值，这种计算的基本关系是种计算的基本关系是：rNW)()()(1qxWpxpxmrNmm)()()(1qxWpxqxmrNmm)()(1上节点pxm)()(1下节点qxm)(pxmrNWmqx)(1m m级级m+1m+1级级（总的蝶形数目总的蝶形数目= = MN/2MN/2个个）21蝶形运算自成独立单元蝶形运算自成独立单元 蝶形运算自成独立单元，即蝶形运算自成独立单元，即 ， 只只 与与 ， 有关，而与其他节点的值无关，同时有关，而与其他节点的值无关，同时 也不参与另外的蝶形运算。因此，就可把计算结也不参与另外的蝶形运算。因此，就可把计算结果果 ， 放入

16、计算前放入计算前 的存储单元中，的存储单元中，而不用再建新的存储单元而不用再建新的存储单元 同址计算同址计算(原位计算原位计算)(1pxm)(1qxm)(pxm)(qxm)(1pxm)(1qxm)(),(qxpxmm同址计算同址计算同址计算所需存储量仅同址计算所需存储量仅等于等于给定数据所需的存储量，给定数据所需的存储量，这可大大节省存储单元这可大大节省存储单元22运算量的减少运算量的减少（1）每个蝶形运算包含两次复数加法，一次复数乘每个蝶形运算包含两次复数加法，一次复数乘 法法（第一级的第一级的2 2点点DFTDFT可除外可除外）（2）N点点FFT需进行需进行 次迭代运算次迭代运算（3）每次

17、迭代运算包括每次迭代运算包括N/2个蝶形个蝶形NM2log 复数加法的次数为：复数加法的次数为： 复数乘法的次数为：复数乘法的次数为：222loglog2FNANNN 221log Nlog22FNNMN -1rNW23运算量的减少运算量的减少N点点DFT： 次复数乘法和次复数乘法和N(N-1)次复数加法次复数加法N点点FFT： 次复数乘法，次复数乘法， 次复数加法次复数加法2log2NNNN2log2N当当N N越大，越大，FFTFFT优势越明显优势越明显24FDMM/FFTFFT算法与算法与DFTDFT算法复数乘法次数比较算法复数乘法次数比较25倒置重排倒置重排输入输入x(n)是是“混序混

18、序”排列排列的，的，“混混 序序”并不是说输入是并不是说输入是杂乱无章的，实际上它是有规律的。如果输入杂乱无章的，实际上它是有规律的。如果输入x(n)的序号的序号用二进制码来用二进制码来 表示，就可以发现输入的顺序恰好是正序表示，就可以发现输入的顺序恰好是正序输入的输入的“码位倒置码位倒置”26倒置重排倒置重排FFT的时间抽选法按的时间抽选法按n的奇偶分离而形成倒置重排的原理如图：的奇偶分离而形成倒置重排的原理如图：2 100 12()()x n nnx n nn输入序列重排实际上就是输入序列重排实际上就是 相互交换相互交换在在C语言中，位操作中循环左移（右移）语言中，位操作中循环左移（右移）

19、2 1 0210()()(000)(0)(001)(1)(010)(2)(011)(3)(100)(4)(101)(5)(110)(6)(111)(7)n n nxxxxxxxxxxxxxxxxxn n nx 十十进进制制01010101(000)(0)(100)(4)(010)(2)(110)(6)(001)(1)(101)(5)(011)(3)(111)(7)xxxxxxxxxxxxxxxx010011n0n1n2271)流程图的每一级的基本计算单元都是一个蝶形；流程图的每一级的基本计算单元都是一个蝶形； 2)输入输入x(n)不按自然顺序排列，称为不按自然顺序排列，称为“混序混序”排列，而

20、排列，而输输 出出 X(k)按自然顺序排列，称为按自然顺序排列，称为“正序正序”排列，因而排列，因而要要 对输入进行对输入进行“变址变址”；3)由于流程图的基本运算单元为蝶形，所以可以进行由于流程图的基本运算单元为蝶形，所以可以进行 “同址同址”计算计算R2-DITR2-DITFFTFFT总结总结28基基 2 时间抽选法信号流图（时间抽选法信号流图（N=8）x-1W8218点点24点点42点点x(0)x(1)x(2)x(3)x(4)x(5)x(6)x(7)x(0)x(4)x(2)x(6)x(1)x(5)x(3)x(7)X(0)X(1)X(2)X(3)X(4)X(5)X(6)X(7)-1-1-1

21、-1-1-1-1-1-1-1-1W80W80W82W80W81W82W83(4)8F(8)8F(2)8F8PC(0)C(1)D(0)D(1)E(0)E(1)F(0)F(1)A(0)A(1)A(2)A(3)B(0)B(1)B(2)B(3)m = 0 级m = 1级m = 2 级位序重排XDIT-FFT重重排排序序29Matlab Matlab 的实现的实现 Matlab 提供了内建的提供了内建的 X = fft(x, N) 函函数来计算矢量数来计算矢量 x 的的 DFT30算法原理：算法原理：时间抽选法时间抽选法：是把输入是把输入 x(n) 按奇偶分按奇偶分解解成两个子序列，而对应的输出是按成两

22、个子序列，而对应的输出是按自自然顺序排列的然顺序排列的X(k)的前、后两半的前、后两半频率抽选法频率抽选法：是把输入是把输入 x(n) 按照自然按照自然顺序前后对半分开顺序前后对半分开，而输出，而输出 X(k) 逐项逐项分解成偶数点子序列和奇数点子序列分解成偶数点子序列和奇数点子序列3.4.3 3.4.3 基基2 2频率抽选法频率抽选法 (DIF-FFT)(DIF-FFT)31基基2 2频率抽选法频率抽选法FFTFFTDIFDIF2 1102( )( )( )NnkkNnnnNNNX kx n Wx n W 2 12 1()0022( )()NNnkkNNnnnNx n Wx nWN ()2

23、120 ( )()2NkNnkNnNNx nx nWW 0 (1)2Nn DFT： 如果将求和分成两部分，一部分是如果将求和分成两部分，一部分是 ，另一，另一部分是部分是 ，则可得：，则可得： 10)()(NnnkNWnxkX (1)2NnN32基基2 2频率抽选法频率抽选法因为因为按按k的奇偶性的奇偶性将将 分解为分解为偶数组偶数组和和奇数组奇数组：12NNW为奇数为偶数kkWkkNN, 1, 1) 1()2(2 1012NnkNnkNX kx nx nW( ) ( ) ()() )(kX1202)2()()2(NnnrNWNnxnxrX1202)2()(NnnrNWNnxnx12()021

24、(21) ( )()2NnrNnNX rx nx nW 1202)2()(NnnrNnNWWNnxnx所以：所以：0,1,.12Nr 33基基2 2频率抽选法频率抽选法若令若令则下式表示两个则下式表示两个N/2点点DFT运算：运算：( )( )()2Ng nx nx n ( )( )()2Nhnx nx n 12, 0Nn12201220(2 )(2( )(1)NnrNnNnNnrnNg nhXrWXrWn W -1)(nx)2(NnxnNW( )( )()2Ng nx nx n( ) ( )()2nnNNNh nWx nx nW0,1,.12Nr 34N/2-1)(nx)2(NnxnNW(

25、)( )()2Ng nx nx n( ) ( )()2nnNNNh nWx nx nW-11( )X k2( )XknNW12(2 )( )( )XkX kXk12(21) ( )( )nNXkX kX k W求求N/2点点DFT35基基2 2频率抽选法频率抽选法将将N点点DFT分解为两个分解为两个N/2点点DFT（N8）()()()2Ng nx nx n ()()()2Nh nx nx n 4 4点点DFTDFT 4 4点点DFTDFT1)0( x) 1 ( x)2( x) 3( x)4( x)5( x)6( x)7( x)0(X) 1 (X(4)X(5)X(2)X(3)X)6(X)7(X0

26、8W18W28W38W111(0)g(1)g(2)g(3)g(0)h(1)h(2)h(3)h36继续分解继续分解37m=0m=1m=2) 0( x) 1 ( x) 2( x) 3 ( x) 4( x) 5 ( x) 6( x) 7( x) 7(X) 6(X) 5 (X) 4(X) 3 (X) 2(X) 1 (X) 0(X) 0(X) 4(X) 2(X) 6(X) 1 (X) 5 (X) 3 (X) 7(X11111111111108W38W28W18W08W28W28W08WN=8N=8时，基时，基2 2频率抽选法的整个信号流图如下：频率抽选法的整个信号流图如下：直到分解为直到分解为2点点DF

27、T为止为止重重排排序序38时间抽选法和频率抽选法的异同时间抽选法和频率抽选法的异同 蝶形计算蝶形计算 时间抽选法、频率抽选法都是把一个N（ ）点DFT先分解成两个N/2（ ）点DFT，再分解成四个N/4点DFT，直到最后分解为N/2个2点DFT 最基本的运算是蝶形运算。时间抽选法和频率时间抽选法和频率抽选法的蝶形运算不同：抽选法的蝶形运算不同：m212m39DIT: 先复乘后加减，先复乘后加减，W 因子在上下节点都有体现因子在上下节点都有体现DIF: 先减后复乘，先减后复乘，W因子仅体现在下节点因子仅体现在下节点DIT-FFT和和DIF-FFT区别区别11( )( )( )( )( )( )m

28、mmmkNkNmmxpxpxqxqxpxWqW11( )( )( ) ( )( )( )mnmmmNmmxpWxpxpxqxxqq)()(1上节点pxm)()(1下节点qxm)(pxmrNmWqx)(1m级m+1级()mxqnNW4022logFNmN乘乘法法：2logFaNN加加法法：运算量相同运算量相同时间抽选法和频率抽选法的异同时间抽选法和频率抽选法的异同DIT 输入重排，输入重排，DIF 输出重排，即输出重排，即1210011()()倒倒置置MMX nn nnX n nn41时间抽选法和频率抽选法的异同时间抽选法和频率抽选法的异同频率抽选法实现的另一方法频率抽选法实现的另一方法信号流图

29、转置定理信号流图转置定理：若流图全部支路方向反向，支：若流图全部支路方向反向，支 路增益不变，输入输出变量交路增益不变，输入输出变量交 换位置，则传输函数保持不变换位置，则传输函数保持不变)(nxca1z)(ny11)(aczczH)(ny)(nxca1z11)(aczczH42DIT 和和 DIF 的基本蝶形互为的基本蝶形互为信号流图转置信号流图转置DITDIF时间抽选法和频率抽选法的异同时间抽选法和频率抽选法的异同(0)x(2)x(1)x(3)x1 1 1 1 (0)X(1)X(2)X(3)X04W14W1 1 1 1 (0)x(1)x(2)x(3)x(0)X(2)X(1)X(3)X04W

30、14W输入倒位序，输出自然序（输入倒位序，输出自然序（DITDIT）输入自然序，输出倒位序（输入自然序，输出倒位序（DIFDIF）43除了基除了基 2 的的 FFT 算法外，还有基算法外，还有基 R （R=3，4，）的快速算法。原理和基）的快速算法。原理和基 2 的类似，可以进一步降低复数乘次数的类似，可以进一步降低复数乘次数基基4算法原理（自学）算法原理（自学）基基4 4时间抽选法时间抽选法N = RVN = 2V44IDFTIDFT的快速算法的快速算法IFFTIFFT由于由于DFT变换对的对称性，变换对的对称性，IFFT的算法与前面所的算法与前面所求的求的FFT算法几乎完全相同：算法几乎完

31、全相同： 10:()()n kNNnD F TXkx n W101:( )( )NNnkkIDFTx nXk WN FFT 算法中的分组方式、排序方式以及蝶形运算结构算法中的分组方式、排序方式以及蝶形运算结构都可用于都可用于IFFT 算法的设计，而这就是可算法的设计，而这就是可依据现有的依据现有的 FFT 算法直接得出算法直接得出 IFFT 算法算法的原因。的原因。 DFT 和和 IDFT 的的区别区别： 因子因子 W 的指数相差一个负号；的指数相差一个负号； 相差一个因子相差一个因子 1/N。45设有序列 x(n)，其 DFT 为 X(k)，则 IDFT 为:12122( )( )( )()

32、( )( )kNkN X kX kWXkNX kX kWXk 在在 FFT FFT 的的时间抽选时间抽选法中法中: : DIT-IFFT 对于对于 IFFT 算法，输入是算法，输入是 X(k) 和和 X(k+N/2)，输出是，输出是 X1(k) 和和 X2(k)。解上式可以得到解上式可以得到 X1(K)、X2(K) ：12122212( )( )()( )( )()kNN X kX kX kNXkWX kX kX(k)X(k+N/2)X1(k)X2(k)-11/2kN1W2101( )( )NnkNkx nX k WN468 点点 DIT-IFFT 算法算法-1x(0)x(1)x(2)x(3)

33、x(4)x(5)x(6)x(7)x(0)x(4)x(2)x(6)x(1)x(5)x(3)x(7)X(0)X(1)X(2)X(3)X(4)X(5)X(6)X(7)-1-1-1-1-1-1-1-1-10.5W8-00.5W8-20.5W8-00.5W8-10.5W8-20.5W8-3-1-10.5W8-00.5W8-2121 21 21 说明：说明：1）分组过程是按时间序列分组过程是按时间序列 x(n) 的奇的奇偶性在时域上展开的，故称此法为时间抽选偶性在时域上展开的，故称此法为时间抽选算法算法 DIT-IFFT；2）1/N 的分解，的分