第七讲 快速傅里叶变换 (FFT) Q&A 办公室: 51971617 手 机：13466573224 手 机：13466573224 Email: tangzhyguo@163.com

本讲在分析直接计算DFT的特点的基础上介绍DFT的快速算法-快速傅里叶变换(FFT)；同时简要介绍了FFT算法的发展历程；此外还要介绍FFT的两种最常用的算法－－基于时间抽取的FFT（DIT：库利－图基算法）和基于频率抽取的FFT（DIF：桑德－图基算法）。

一、直接计算DFT存在的问题 N点序列x(n)的DFT变换定义为： 计算X(k)的运算量：需要N2次复数乘法，N(N－1)次复数加法。在N较大时计算量很大。 例如：N＝1024时, 需要1,048,576次复数乘法, 即4,194,304次实数乘法

对于象雷达、通信、声纳等需要实时处理的信号，因为其运算量更大，所以无法满足信号处理的实时性要求。迫切需要有新的算法。 二、DFT运算的特点 实际上，DFT运算中包含有大量的重复运算。在WN 矩阵中，虽然其中有N2个元素，但由于WN的周期性，其中只有N个独立的值，即 ，且这N个值也有一些对称关系。总之，WN因子具有如下所述周期性及对称性：

1.对称性 2.周期性 由上述特性还可得出： 利用上述对称特性，可使DFT运算中有些项可以合并，这样，可使乘法次数减少大约一半；利用WN矩阵的对称性及周期性，可以将长序列的DFT分解为短序列的DFT，N越小，运算量能够减少。 例如，对于四点的DFT，直接计算需要16次复数乘法，根据上述特性可以有以下形式的算法：

则有： 第二列和第三列交换，则：

由此得出： 从上例可知，通过应用对称性和周期性，4点的DFT实际上只需要进行一次复数乘法。

三、FFT发展简介 FFT的实质：快速傅里叶变换(FFT)并不是一种新的变换，是为了改进和提高离散傅里叶变换(DFT)运算速度基于DFT运算特点而发展起来的DFT快速算法，其实质还是DFT。 FFT发展的原因：DFT是信号分析与处理中的一种重要变换，广泛应用于通信、图像处理、雷达及声纳等领域，由于其计算量与变换区间长度N的平方成正比，在N较大时，计算量很大，使得直接应用DFT进行实时处理信号是不现实的。

FFT的发展历程： 1965年，J. W. Cooley和J. W. Tukey巧妙应用DFT中W因子的周期性及对称性提出了最早的FFT，这是基于时间抽取的FFT。具有里程碑式的贡献(运算量缩短两个数量级) 1966年，G. Sand提出了基于频率抽取的FFT算法 1975年，Winogard提出WFTA法；1977年Kolha和Parks提出素因子算法（PFA） 1984年，P. Dohamel和H. Hollmann提出分裂基快速算法，进一步减少了计算量，提高了计算速度（目前最理想的算法）

FFT的各种算法 纵观FFT的发展历程，FFT算法分成两大类： (1) 针对N等于2的整数次幂的算法，如基2算法、基4算法、实因子算法和分裂基算法 (2)针对N不等于2的整数次幂的算法，其以Winograd为代表的一类算法(素因子法PFA、Winograd算法WFTA) 简要介绍库利-图基算法和桑得-图基算法

四、按时间抽取(DIT)的FFT－－库利-图基算法 1.基本原理 设序列x(n)的长度为N，且满足N＝2M，M为自然数，按n的奇偶将x(n)分解为两个N/2的子序列： x1(r)=x(2r), r=0,1,2,…,N/2-1 x2(r)=x(2r+1) r=0,1,2,…,N/2-1 则x(n)的DFT为：

因为 ，所以： 上式中X1(k)和X2(k)分别为x2(r)和x2(r)的N/2点DFT，即 由于X1(k)和X2(k)均以N/2为周期，且 ，以X(k)又可表示为：

即将一个N点的DFT分解成为两个N/2点的DFT。上述运算可用右下图来表示，称为蝶形运算符号。 从右图可知，要完成一个蝶形运算需要进行一次复数相乘和两次复数相加运算。

下图是N＝8时的一个分解运算图。 从上图可知，经过一次分解后，计算一个N点的DFT共需要计算两个N/2点FFT和N/2个蝶形运算。

计算一个N/2点DFT需要(N/2)2复数乘和N/2(N/2-1)次复数加法。所以按刚才的方法计算N点DFT总的运算量为2(N/2)2+N/2=N(N+1)/2≈N2/2(N>>1时)复数次乘法和N(N/2-1)+2N/2=N2/2次复数加法运算。 由此可见，仅仅经过一次分解就能使运算量减少近一半！ 因为N/2仍然是偶数，可以作进一步的分解： 与第一次分解相同，将x1(r)按奇偶分解成两个N/4 的子序列x3(l)和x4(l), 即：

则，X1(k)又可表示为： 同理，X3(k)和X4(k)的周期性和WN的对称性，到最后我们能够得到：

同理可得： 其中：

这样，又将N/2点的DFT分解为两个N/4点的DFT。依次类推，经过M-1次分解，最后将N点DFT分解成N/2个2点DFT。一个完整的8点DFT-FFT运算流图如下图所示。

2.运算量的比较 从上述分析过程可知，在N＝2M时，每一级都由N/2个蝶形运算构成，即每级都需要N/2次复数乘和N次复数加，所以总的复数乘的次数为： 总的复数加的次数为： 直接计算时复数乘的次数为N2，加为N(N-1)次。当N>>1时， ，使运算量大大减少。

以N＝1024为例，其运算量与直接计算的比例为： 即运算效率提高了200多倍。易知N越大，优越性越明显。另外，在N＝2048时，直接运算需要3个小时，而采用FFT则只需不到一分钟就能完成！ 3.DIT－FFT的运算规律 (1)原位计算

根据运算流图可知，DIT-FFT的运算很有规律。N＝2M点的FFT共进行M级运算，每级运算有N/2个蝶形运算构成；同一级中，每个蝶形的两个输入数据只对计算本蝶形有用，并且每个蝶形的输入、输出数据节点又同在一条水平线上，这意味着计算完一个蝶形后所得数据可立即存入原输入数据所占用的存贮单元，这样，经过M级运算后，原来存放输入序列数据的N个存贮单元中并依次存放了X(k)的N个值。这种利用同一存贮单元存贮计算输入、输出数据的方法称为原位（址）计算。

输入数据、中间运算结果和最后输出均用同一存储器。 (0)=X0(0) X1(0) X2(0) X3(0)=X(0) (4)=X0(1) X1(1) X2(1) X3(1)=X(1) (2)=X0(2) X3(2)=X(2) (6)=X0(3) X3(3)=X(3) (1)=X0(4) X1(4) X2(4) X3(4)=X(4) (5)=X0(5) X3(5)=X(5) (3)=X0(6) X3(6)=X(6) (7)=X0(7) X1(7) X2(7) X3(7)=X(7) W N -1 2 1 3 .

(2)旋转因子的变化规律 在每个蝶形的运算过程中，都要乘以因子 ，称其为旋转因子，p称为旋转因子指数。但各级的旋转因子和循环方式都有所不同。 旋转因子与运算级数有一定的关系，若用L表示运算级数，对于N＝2M的一般情况，第L级的旋转因子为：

(3)倒位序规律 从运算流图可以看出，原位计算时，FFT的输出X(k)是按正常顺序排列在存储单元中，即按X(0)，X(1),···，X(7)的顺序排列，但是这时输入x(n)都不是按自然顺序存储的，这看起来好象是“混乱无序”的，实际上是有规律的，我们称之为倒位序。 造成倒位序的原因是输入x(n)按标号n的偶奇的不断分组而造成。

x(0) x(4) x(2) x(6) x(1) x(5) x(3) x(7) 45 倒位序实现 输入序列先按自然顺序存入存储单元,然后经变址运算来实现倒位序排列，设输入序列的序号为n,二进制为(n2 n1 n0 )2 ,倒位序顺序用 表示,其倒位序二进制为(n0 n1 n2 )2 。 A(1) A(2) A(3) A(4) A(5) A(6) A(7) A(8) 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) 变址处理方法 存储单元 自然顺序 变址 倒位序

例如 ，N=8时如下表： 自然顺序n 二进制n n n 倒位序二进制n n n 倒位顺序n ^ 2 1 0 0 1 2 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

4.蝶形运算两节点的距离:2m-1 其中,m表示第m列,且m =1,… ,L 例如N=8=23 ,第一级(列)距离为21-1=1, 第二级(列)距离为22-1=2， 第三级(列)距离为23-1=4。

五、按频率抽取(DIF)的FFT－－桑德-图基算法 库利图基法是将输入序列按其顺序是奇数还是偶数来分解为越来越短的序列；桑德图基法是把输出序列X(k)按其顺序的偶奇来分解为越来越短的序列。 1.原理 设序列x(n)长度为N＝2M，首先将x(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为偶数时：

k为奇数时： 可见,上面两式均为N/2的DFT。

3.蝶形运算 -1

4.N=8时,按k的奇偶分解过程 先蝶形运算，后DFT: W N 1 2 3 -1 -1 -1 -1

再将N/2点DFT按k的奇偶分解为两个N/4点的DFT,如此进行下去,直至分解为2点DFT。 以下是8点的DIF-DFT流程： 5.仿照DIT的方法 再将N/2点DFT按k的奇偶分解为两个N/4点的DFT,如此进行下去,直至分解为2点DFT。 以下是8点的DIF-DFT流程： x(0) X(0) x(1) X(4) x(2) X(2) x(3) X(6) x(4) X(1) x(5) X(5) x(6) X(3) x(7) X(7) -1 W N 1 2 3

B.原位运算 每级(列)都是由N/2个蝶形运算构成,即 r -1 W N

C.蝶形运算两节点的距离 一般公式为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。

D.DIF法与DIT法的异同 1.相同点 (1)进行原位运算； (2)运算量相同,均为（N/2) Log2N次复乘、N Log2N次复加。 2.不同点 (1)DIT输入为倒位序,输出为自然顺序； DIF正好与此相反。但DIT也有输入为自然顺序,输出为倒位序的情况。

(2)蝶形运算不同 A.DIT 用矩阵表示 1 = 1

B.DIF 用矩阵表示 = 1

A. DIT B.DIF (3)两种蝶形运算的关系－－互为转置（矩阵）； 将流图的所有支路方向都反向,交换输入和输出，即可得到另一种蝶形。 1 1 1 1 将流图的所有支路方向都反向,交换输入和输出，即可得到另一种蝶形。

六、IDFT的快速算法--IFFT FFT算法，同样可以适用于离散博里叶反变换(IDFT)运算，并简称为IFFT，即快速博里叶反变换，从IDFT公式看： 而DFT公式 比较以上两式可知，只要把DFT运算中的每一个系数 换成 ，并且最后再乘以常数1/N, 则

以上所有时间抽取或频率抽取的FFT算法都可以拿来运算IDFT。 此外，有一种直接采用FFT程序计算IFFT的方法。 对 取共轭，则： 两边同时再取共轭，有：

上述情况表明：在进行IDFT运算的时候，可以先将X(k)取共轭，直接调用FFT子程序，然后再取共轭并乘以1/N就可得到序列x(n)。