上次我解释了DFT (离散傅里叶变换)的导出，但这个DFT是一个基本且非常有用的表达式，用于数值计算中的傅里叶变换。然而，当试图根据公式计算时，需要大量的计算量，并且没有实际使用，但是Cooley和Turkey在1965年开发了快速傅里叶变换(Fast Fourier Transform, FFT)算法。自那时起，它变得非常流行，现在它是数字信号处理最有用的技术。最近，电视广播也是数字化的，FFT技术用于称为OFDM的数字调制技术。
这次，简单说明一下FFT (快速傅里叶变换)的计算方法。
如前所述，DFT从N个时间信号xn (n=0、.....、N-1)产生N个频谱X k (k=0、.....、N-1)。也就是说;
.................................（1）
其中:
.................................（2）
通常，公式 (1) 中的DFT计算为复数，需要N2次复数乘法才能得到N点的Xk。现在N可以用多个整数的乘积表示的情况下，把公式 (1) 的演算顺利地分解成小组，并且使用W的周期性 (周期N)，尽量减少乘法的次数，高速进行DFT演算的算法是FFT。
对于FFT计算点N (时间函数的样本点)，通常根据算术分区的简单性和数字计算机程序实现的兼容性选择2的幂 (N=2m) 。例如，N=1024 (m=10)，N=2048 (m=11) 。这称为两个基本FFT。
在下面的说明中，N=8 (m=3) 。将N=8代入公式 (1);
.................................（3）
.................................（4）
按钮，将选定控件在Tab键次序中下移一个位置。W表示旋转因子，因为它能够以负角 (顺时针) 将单位圆旋转π/4 (45°) 。旋转因子W在旋转一周时具有相同的值，而在旋转半周时具有相反的符号。当N=8时，单位圆上的具体点如图1所示。
-
图1单位圆上N=8时的旋转因子W
接着，将时间信号x n分割为偶数列和奇数列，则式 (3);
.................................（5）
其中r=0、...和3。
式 (5) 表示偶数项4点DFT和奇数项4点DFT乘以旋转因子W k之和。注意W k =-W k-4 (图1)，此操作的流程图如图2所示。
-
图2偶数项和奇数项分别分解为2个4点DFT时的式 (5) 的信号流程图。
以完全相同的方式，通过将四点DFT分解为两个两点DFT,最终八点DFT的FFT为如图3所示的信号流程图。请注意，此方法是分割时间信号的方法，因此称为时间拉长间隔法和。
-
图3基于时间拉长间隔法的8点FFT运算的信号流程图
如果仔细查看图3中的流程图，您会发现每个计算列都只从与前一列相同的节点对进行计算。即，FFT运算以公式 (6) 形式的积和对为基本运算。
.................................（6）
该操作如图4所示，它被称为蝴蝶操作，因为它的形状像蝴蝶的翅膀。
-
图4基于时间拉长间隔法的蝴蝶演算
式 (6) 的意思是，将一对节点数据设为xα、xβ，从时间信号 (复数) 的存储器中调出该数据，进行右边的计算，将其结果存储在相同的xα、xβ存储器中。这样，作为FFT基础运算的蝶形运算，由于只需一对节点数据即可进行，因此被称为In-Place存储器运算，由于存储量少且存储器访问次数也少，因此具有效率非常高的优点。
关于图3中的流程图，还有一点需要注意。时间信号xn的序列不再按照正常采样顺序。因此，在执行FFT操作之前，必须按图3所示重新排列。如果您尝试将数字序列的顺序数字设置为二进制数字，则可以轻松重新排列它们。
表1显示了N=8的示例。这样，通过将原始序列的序号重写为二进制数字并翻转位序列，可以导出如图3所示的序列。将此操作位倒序我们称之为(位反转，BTR)。此操作也可以在普通计算机上轻松完成。
表1基于N=8时的比特逆序的排序
| 原顺序 | 位元逆序 | ||
| 10进制数 | 二进制数 | 二进制数 | 10进制数 |
| 0 | 0 | 0 | 0 |
| 1 | 1 | 100 | 4 |
| 2 | 10 | 10 | 2 |
| 3 | 11 | 110 | 6 |
| 4 | 100 | 1 | 1 |
| 5 | 101 | 101 | 5 |
| 6 | 110 | 11 | 3 |
| 7 | 111 | 111 | 7 |
这次，通过将DFT运算结果即X k (k=0、.....、N-1)分割为偶数和奇数，与时间疏除法同样，能够导出频率疏除法。图5是在N=8时使用频率稀疏方法的8点FFT的信号流程图。
-
图5基于频率间隔法的8点FFT运算的信号流程图
在这种情况下，蝶形运算如公式 (7) 所示，其信号流程图如图6所示。
.................................（7）
-
图6基于频率间隔法的蝴蝶演算
如图5所示，在频率间隔法中，通常运算结果的X k不按频率顺序排列，需要进行BTR处理。通过与表1相反的处理，可以获得正确的顺序列。
还有许多其他计算方法，例如使用时间稀疏方法对结果进行BTR,或首先使用频率稀疏方法进行BTR,但这些方法很少使用，因为访问旋转因子W时需要进行BTR处理。
与普通DFT相比，此FFT计算可节省多少计算时间?式 (6) 的蝶形运算中，作为进行1次复数乘法，由于在N点的FFT中进行NN log 2 2的蝶形运算，与DFT的乘法次数N2相比，N=1024 (=2 10) 的情况下约为1/200, N=16384 (=2 14) 的情况下，约为1/2300,可惊人地缩短。
另外，通过减少乘法次数，可以减少舍入误差引起的计算误差，被称为FFT算法的一大优点。
最后，总结一下。
- FFT是一种快速计算DFT的算法。
- FFT计算点N通常采用2的幂，例如1024或2048,称为2基FFT。
- FFT算法大致有时间疏除法和频率疏除法。
- FFT的基本运算被称为蝴蝶运算，因为是In-Place内存运算，所以可以节约内存。
- 与通常的DFT相比，FFT可以大幅减少乘法次数，因此不仅可以提高运算速度，还可以降低计算误差。
【关键词】
DFT (离散傅立叶变换)、FFT (快速傅立叶变换)、2基FFT、旋转因子、时间拉长间隔法、蝴蝶运算、位倒序 (BTR)、就地存储运算、频率拉长间隔法
【参考资料】
- “快速傅立叶变换”E.ORAN BRIGHAM著科学技术出版社
- “数字傅里叶分析 (I) -基础-”由Kenichi Shirodo 新冠公司撰写
- “信号处理”森下岩/小畑秀文共著朝仓书店
(摘自2013年5月23日发行的电子邮件杂志)
