数字信号处理：第六章 快速傅里叶变换1

1、第六章:快速傅里叶变换The Fast Fourier Transform （FFT）FFT：快速傅里叶变换（ Fast Fourier Transform ）快速：用较少的基本运算完成傅里叶变换。意义：节约计算资源；提高实时性标志性文章： 1965年 Cooley 和 Tukey 发表的 “机器计算傅里叶级数的一种算法” 6-1 引言 采用这种算法能使计算机计算离散傅里叶变换所需要的乘法次数大为减少，特别是被变换的抽样点数N越多，FFT算法计算量的节省就越显著适用领域范围:有限长序列6-2 直接计算DFT的问题及改进的途径一.直接计算DFT的计算量： 正变换反变换二者的差别仅在指数的符号和因

2、子 ；二者有相同的计算量。 乘法：N 次 加法：N-1 次 注：通常x(n)和 都是复数，所以上述运算量为复数乘法、复数加法运算量。 复数乘法：N2 次； 复数加法： N(N-1) 次 当N很大时运算量是惊人的，如N=1024，则要完成1048576 次(一百多万次)运算，难以实现实时处理。 一个X(k)的值的计算量：例如X(1)： 全部N 个X(k)的值的计算量：实数乘法实数加法一次复乘42一次复加2一个X (k)4N2N +2 (N 1)=2(2N 1)N个X (k)(N点DFT)4N 22N (2N 1) 实数运算量：计算量大处理实时性差二. 改进的途径 系数 的特性：可推得:对称性周期

3、性复共轭对称性减少运算量的思路： 1） 利用 的上述性质合并一些乘法运算； 2） 将长序列的DFT运算分解为短序列的DFT运算。1965年，Cooley和Tukey首先提出一种FFT算法。对于N点DFT,仅需(N/2)log2N 次复数乘法运算. 例如 N=1024=210 时，需要复乘 (1024/2)log2210 =512*10=5120次。 5120/1048576=4.88% ,速度提高20倍。根据算法推导过程中数据抽选分组所在的域，可分为时间抽选法（DIT）: Decimation-In-Time频率抽选法（DIF）: Decimation-In-FrequencyFFT算法分类：

4、另外，也可按分解过程中的比值分为基2、基4等FFT算法。Radix-2 、Radix-4 6-3 按时间抽取(DIT)的基-2 FFT算法 库利-图基算法一.算法原理： 设 N=2L ，不足时可补零。(一) N点DFT化为N/2点DFT：1.将 按 n 的奇偶分为两组作DFT： n为偶数时: n为奇数时:将 按照n为奇偶分开Radix-2 其中,2.两点结论: (1) 均为N/2点的DFT。 (2) 只能确定出 时的X(k)，即前一半的结果。 同理, 这就是说，X1(k)、X2(k) 对应的后一半值，为其前一半值的重复。3. X(k)的后一半的确定：由于 （周期性）,所以又由于 ,所以 可见,

5、X(k)的后一半，也完全由X1(k)、X2(k)的 个系数确定。 结论：N点的DFT可由两个N/2点的 FT来计算实现上式运算的流图称作蝶形运算：前一半4.蝶形运算后一半（N/2个) 1 1 11-1由X1(k)、X 2(k) 表示X(k)： 对于传输系数为1的支路，可不标系数，这样就可以将蝶形画成每个蝶形的运算量： 1次复乘，2次复加。（乘 -1不算）例：N=85.计算量分析直接计算N点DFT的复乘为N2 ，复加N(N-1)；分解后与直接计算相比,计算量差不多减少 一半。两个 点DFT 个蝶形总计21一个蝶形一个 点DFT复数加法复数乘法 由于N=2L ,所以 N/2仍为偶数,可以进一步把每

6、个N/2点的序列再按其奇偶部分分解为两个N/4的子序列。例如，n为偶数时的 N/2点分解为:(二) N/2点DFT化为N/4点DFT进行N/4点的DFT,得到(偶中偶)(偶中奇)从而可得到前N/4的 X1(k)：类似可得后N/4的 X1(k)为(奇中偶)(奇中奇) 同样，n为奇数的N/2点DFT也可以分解为两个N/4点： 如此分解下去，直到分解成2点DFT。 对于N=8，N/4点即为两点DFT： 同样， 亦即 每个2点DFT也是一个蝶形运算。 (0) (4) (2) (6) (1) (5) (3) (7)WN0WN0WN0W0N-1-1-1-1X (0)X (1)X (0)X (1)X (0)

7、X (1)X (0)X (1)33445566WN0WN2WN0WN2-1-1-1-1X (0)X (1)X (2)X (3)X (0)X (1)X (2)X (3)11121222WWWWN0N1N2N3-1-1-1-1X(0)X(1)X(2)X(3)X(4)X(5)X(6)X(7)8点DFT的FFT的运算流图如下：二.运算量 当N = 2L时，共有L级蝶形，每级N / 2个蝶形， 每个蝶形有1次复数乘法、2次复数加法。复数乘法：复数加法：与直接计算DFT 比较乘法次数：点数DFT复数乘法FFT复数乘法NN 2N L/24 =2 216 48 =2 3641216 =2 425632256

8、=2 8655361024512=2 926214423041024 =2 10104857651202048 =2 114194304112644096 =2 1216,777,21624,5762 1367,108,8642 14268,435,4562 151,073,741,8242 164,294,967,296=8192* 524,288524,288三. DIT的FFT算法的特点1.原位运算(同址计算)输入数据、中间运算结果和最后输出均用同一存储器。m表示第m级迭代，k，j表示数据所在的行数。 2.倒位序规律 由图可知,输出X(k) 按正常顺序排列在存储单元,而输入是按顺序:这种

9、顺序称作倒位序, 即二进制数倒位。000， 100， 010， 110， 001， 101， 011， 111n =00n =10n =01n =11n =01n =1101010101 n2x(000) 0 x(100) 4 x(010) 2 x(110) 6 x(001) 1x(101) 5 x(011) 3 x(111) 7 （偶）（奇） 这是由奇偶分组造成的,每次都是偶数号分为一组，奇数号分为一组，对应于二进制的0、1。以N=8为例说明如下: 3.倒位序实现 输入序列先按自然顺序存入存储单元,然后经变址运算来实现倒位序排列。 设输入序列的序号为n ,二进制为(n2 n1 n0 )2 ,

10、倒位序顺序用 表示,其二进制为(n0 n1 n2 )2 , 例如 ，N=8时如下表： 0 0 0 0 0 0 0 0 1 0 0 1 1 0 0 4 2 0 1 0 0 1 0 2 3 0 1 1 1 1 0 6 4 1 0 0 0 0 1 1 5 1 0 1 1 0 1 5 6 1 1 0 0 1 1 3 7 1 1 1 1 1 1 7 自然顺序n 二进制n n n 倒位序二进制n n n 倒位顺序 n 2 1 0 0 1 2注：序号总是互为倒位序A(0) A(1) A(2) A(3) A(4) A(5) A(6) A(7)x(0) x(1) x(2) x(3) x(4) x(5) x(6)

11、 x(7)x(0) x(4) x(2) x(6) x(1) x(5) x(3) x(7)变址处理方法存储单元自然顺序变址倒位序(1) 从000开始，逐一与其倒位序对应的序号进行数据互换。(2) 为保证只交换一次，应先比较n与n的大小，再决定是否 互换 4.蝶形运算两节点的距离: 2m-1 其中m表示第m列（级）, m =1, ,L 例如 N = 8 =23 ,第一级(列)距离为21-1=1, 第二级(列)距离为22-1=2， 第三级(列)距离为23-1=4。第m级运算：每个蝶形的两节点距离为 2m1 m=1m=2m=3级号m：5. 考虑蝶形运算两节点的距离为2m-1 ,蝶形运算可表示为结论：

12、蝶形运算两节点的第一个节点为k值，表示成L位二进制数，左移L m位，把右边空出的位置补零，结果为r的二进制数。0 0 00 0 0 0 0 0 0 00 0 0L=3, m=1, k=0r=01 0 11 0 1 1 0 1 00 1 0L=3, m=2, k=5r=20 1 10 1 1 L=3, m=3, k=3r=30 1 10 1 1 例： 6.存储单元 输入序列x(n) : N个存储单元系数 ：N / 2个存储单元 共计需要 ( N+N/2 ) 个存储单元。四. DIT算法的其他形式流图输入倒位序、输出自然序输入自然序、输出倒位序输入输出均自然序相同几何形状输入倒位序输出自然序输入自

13、然序输出倒位序对于任何流图，只要保持各节点所连的支路及其传输系数不变，不论节点位置怎么排列所得流图总是等效的。 蝶形相同，运算量一样；取用系数的顺序不同。6-4 按频率抽选（DIF）的基-2 FFT算法 (桑德图基算法) 一.算法原理： 1.N点DFT的另一种表达形式 由于 故 因此 X(k)可表示为 2.N点DFT按k的奇偶分组可分为两个N/2点的DFT： 当k为偶数,即 k=2r 时,(-1)k =1； 当k为奇数,即 k=2r+1 时 (-1)k =-1 。 这时 X(k) 可分为两部分： k为偶数时： 可见,上面两式均为N/2点的DFT。k为奇数时：-13.蝶形运算4. N=8时,按k

14、的奇偶分解过程： 先蝶形运算，后DFT:x1(0)x1(1)-1x1(2)x1(3)-1x2(0)x2(1)-1x2(2)x2(3)-1N/2点DFTN/2点DFTx(0)x(7)x(1)x(2)x(3)x(4)x(5)x(6)X1(0)=X(0)X2(0)=X(1)X1(1)=X(2)X1(2)=X(4)X1(3)=X(6)X2(1)=X(3)X2(2)=X(5)X2(3)=X(7) 5.进一步的分解： 再将N/2点DFT按k的奇偶分解为两个N/4点的DFT：x3(0)x3(1)-1-1x4(0)x4(1)N/4点DFTN/4点DFTx1(0)x1(1)x1(2)x1(3)X3(0)=X1(

15、0)=X(0)X4(0)=X1(1)=X(2)X3(1)=X1(2)=X(4)X4(1)=X1(3)=X(6)如此进行下去,直至分解为2点DFT： 二. 算法特点： 1）原位运算 每级(列)都由N/2个蝶形运算构成,每个蝶形的运算为： -12）蝶形运算两节点的距离： 距离为： 2L-m=N / 2m 例如 N =23 = 8 ： (1) m=1 时的距离为 8/2=4； (2) m=2 时的距离为 8/4=2; (3) m=3 时的距离为 8/8=1。m=1m=2m=3 r 的求法： 蝶形运算两节点的第一个节点为k值，表示成L位二进制数，左移m-1位，把右边空出的位置补零，结果为 r 的二进制

16、数。 3） 的计算 由于DIF蝶形运算的两节点的距离为 N/2m , 所以蝶形运算可表为：0 1 10 1 1 0 1 10 1 1L=3, m=1, k=3r=31 0 11 0 1 1 0 1 00 1 0L=3, m=2, k=5r=21 0 01 0 0 L=3, m=3, k=4r=01 0 0 0 00 0 0 例： 1.相同点 (1)进行原位运算； (2)运算量相同,均为复乘 复加 2.不同点 (1) DIT输入为倒位序,输出为自然顺序； DIF正好与此相反。但DIT也有输入为自然顺序,输出为倒位序的情况。三. DIF法与DIT法的比较：(2) 蝶形运算不同DIT的蝶形：用矩阵表

17、示=11DIF的蝶形：用矩阵表示=11(3) 两种蝶形运算的关系： 互为转置（矩阵） 转置定理：将流图的所有支路方向倒转,并交换输入与输出，所得系统的输入输出关系不变。（图4.5的DIT与图4.17的DIF基2FFT流图实质为一个图） DIT：DIF:1111转置后即为： 6-5 IDFT的快速计算方法一. 稍微变动FFT程序和参数可实现IFFT： 比较两式的区别：比较可知,对DFT作如下修改便可实现IDFT: (1) (2) 输出前乘以常数1/N。 另外,可以将常数1/N分配到每级运算中: , 每级蝶形运算均乘 1/2。 二. 利用FFT程序直接实现IFFT: 利用FFT子程序实现IDFT的

18、步骤： （1）取共轭：即将 的虚部乘 -1 ； （2）调用FFT子程序：计算DFT； （3）取共轭：对FFT的输出再取一次共轭； （4）乘1/N ： 即得 。共轭FFT共轭乘1/ N6-6 基-4 FFT 算法基-4FFT条件：基本思想：将N点DFT化为4个 点DFT。一. 按时间抽选（DIT）：Radix-4 N/4DFTN/4DFTN/4DFTN/4DFT-j-1j-1-1-1j-jN/4DFTN/4DFTN/4DFTN/4DFT- j-1-1-1-1二. 按频率抽选（DIF）： 其中这里的系数矩阵与DIT中的矩阵相同。N/4DFTN/4DFTN/4DFTN/4DFT-1-j-1j-1-1

19、j-jN/4DFTN/4DFTN/4DFTN/4DFT- j-1-1-1-1三. 推广： 仍为4的整数倍，可再分解为4个 点 DFT的组合。 4点DFT:即所以也可以用前面的蝶形实现。频率四.计算量 一个基4蝶形需 3 次乘旋转因子（ 除外）； 每级N/4个基4蝶形； 最后分解的4点FFT不需要乘旋转因子。 所以，基4蝶形共有 (L-1) 级。基 -2FFT：6-7 N为复合数的FFT算法混合基算法基-4FFT条件：要求基-2FFT算法：要求方法一：采用其它基，如基-3、基-5、基-7等等。方法二：采用混合基，即 每级采用不同的基。例：N=120 第一级：采用基-4，分解为4个30点DFT；（基-4蝶形） 第二级：采用基-3，将30点DFT分解为3个10点DFT； （基-3蝶形） 第三级：采用基-5，将10点DFT分解为5个2点DFT； （基-5蝶形） 第四级：采用基-2，实现2点DFT。 （基-2蝶形）优点：不需要补0，不增加数据量；缺点：各级不统一，程序复杂。6-8 分裂基FFT算法分裂基：对偶序号使用基-2FFT算法， 对奇序号使用基-4FFT算法。 特点：针对 的算法中具有最少乘法次数， 且