快速傅里叶变换

快速傅里叶变换（Fast Fourier Transform, FFT），是快速计算序列的离散傅里叶变换（DFT）或其逆变换的方法。傅里叶分析将信号从原始域（通常是时间或空间）转换到频域的表示或者逆过来转换。FFT会通过把DFT矩阵分解为稀疏（大多为零）因子之积来快速计算此类变换。因此，它能够将计算DFT的复杂度从只用DFT定义计算需要的 $O(n^2)$，降低到 $O(n \log n)$，其中 $n$ 为数据大小。

快速傅里叶变换广泛的应用于工程、科学和数学领域。这里的基本思想在1965年才得到普及，但早在1805年就已推导出来。1994年美国数学家吉尔伯特·斯特朗把FFT描述为“我们一生中最重要的数值算法”，它还被IEEE科学与工程计算期刊列入20世纪十大算法。

## 定义和速度

用FFT计算DFT会得到与直接用DFT定义计算相同的结果；最重要的区别是FFT更快。（由于舍入误差的存在，许多FFT算法还会比直接运用定义求值精确很多，后面会讨论到这一点。）

令 $X_0, \cdots , X_{N-1}$ 为复数。DFT由下式定义

$$X_k = \sum_{n=0}^{N-1} x_n e^{-{i 2\pi k \frac{n}{N}}} \qquad k = 0,\dots,N-1.$$

直接按这个定义求值需要 $O(N^2)$ 次运算：$X^k$ 共有 $N$ 个输出，每个输出需要 $N$ 项求和。直接使用DFT运算需使用 $N$ 个复数乘法（$4N$ 个实数乘法）与 $N-1$ 个复数加法（$2N-2$ 个实数加法），因此，计算使用DFT所有 $N$ 点的值需要 $N^2$ 复数乘法与 $N^2-N$ 个复数加法。FFT则是能够在 $O(N\log{}N)$ 次操作计算出相同结果的任何方法。更准确的说，所有已知的FFT算法都需要 $O(N\log{}N)$ 次运算（技术上 $O$ 只标记上界），虽然还没有已知的证据证明更低的复杂度是不可能的。

要说明FFT节省时间的方式，就得考虑复数相乘和相加的次数。直接计算DFT的值涉及到 $N^2$ 次复数相乘和 $N(N-1)$ 次复数相加（可以通过削去琐碎运算（如乘以1）来节省 $O(N)$ 次运算）。众所周知的基2库利-图基算法，$N$ 为2的幂，可以只用 $(N/2)\log_{2} N$ 次复数乘法（再次忽略乘以1的简化）和 $N\log_{2} N$ 次加法就可以得到相同结果。在实际中，现代计算机通常的实际性能通常不受算术运算的速度和对复杂主体的分析主导，但是从$O(N^2)$到$O(N\log{}N)$的总体改进仍然能够体现出来。

## 一般的简化理论

假设一个 $m\times n$ 型矩阵 $\mathbf{S}$ 可分解成列向量以及行向量相乘：

$$\mathbf{S}=\begin{bmatrix} a_1 \\ a_2 \\ \vdots \\ a_m\end{bmatrix}\begin{bmatrix} b_1 & b_2 & \cdots & b_n\end{bmatrix}$$

若 $\begin{bmatrix} a_1 & a_2 & \cdots & a_m\end{bmatrix}^T$ 有 $M_0$ 个相异的非平凡值（$a_m\ne\pm2^k,a_m\ne\pm2^ka_n$ 其中 $m\ne n$）

$\begin{bmatrix} b_1 & b_2 & \cdots & b_n\end{bmatrix}$ 有 $N_0$ 个相异的非平凡值

则 $\mathbf{S}$ 共需要 $M_0\times N_0$ 个乘法。

$$\begin{bmatrix} Z[1] \\ Z[2] \\ \vdots \\ Z[N] \end{bmatrix}= \mathbf{S}\begin{bmatrix} X[1] \\ X[2] \\ \vdots \\ X[N]\end{bmatrix}=\begin{bmatrix} a_1 \\ a_2 \\ \vdots \\ a_m\end{bmatrix}\begin{bmatrix} b_1 & b_2 & \cdots & b_n\end{bmatrix}\begin{bmatrix} X[1] \\ X[2] \\ \vdots \\ X[N]\end{bmatrix}$$

Step 1：$Z_a=b_1X[1]+b_2X[2]+\cdots+b_nX[N]$

Step 2：$Z[1]=a_1Z_a,Z[2]=a_2Z_a,\cdots,Z[N]=a_mZ_a$

简化理论的变型：

$$\mathbf{S}=\begin{bmatrix} a_1 \\ a_2 \\ \vdots \\ a_m\end{bmatrix}\begin{bmatrix} b_1 & b_2 & \cdots & b_n\end{bmatrix}+\mathbf{S}_1$$

$\mathbf{S}_1$ 也是一个 $M\times N$ 的矩阵。

若 $\mathbf{S}_1$ 有 $P_1$ 个值不等于 $0$，则 $\mathbf{S}$ 的乘法量上限为 $M_0+N_0+P_1$。

## 快速傅立叶变换乘法量的计算

假设 $N= P_1 \times P_2 \times \cdots \times P_k$，其中 $P_1,P_2, \cdots , P_k$ 彼此互质。

$P_k$ 点DFT的乘法量为 $B_k$，则 $N$ 点DFT的乘法量为：

$$\frac{N}{P_1}B_1+\frac{N}{P_2}B_2+\cdots\cdots+\frac{N}{P_k}B_k$$

假设 $N=P^c$，$P$ 是一个质数。

若 $N_1=P^a$ 点的DFT需要的乘法量为 $B_1$

且 $n_1\times n_2$ 当中（$n_1=0,1,\cdots,N_1-1, \quad n_2=0,1, \cdots , N_2-1$）

有 $D_1$ 个值不为 $\frac{N}{12}$ 及 $\frac{N}{8}$ 的倍数。

有 $D_2$ 个值为 $\frac{N}{12}$ 及 $\frac{N}{8}$ 的倍数，但不为 $\frac{N}{4}$ 的倍数。

则 $N$ 点DFT的乘法量为：

$$N_2B_1+N_1B_2+3D_1+2D_2$$

## 库利-图基算法

库利-图基算法是最常见的FFT算法。这一方法以分治法为策略递归地将长度为$N=N_1 N_2$的离散傅里叶变换分解为长度为$N_1$的$N_2$个较短序列的离散傅里叶变换，以及与$\Omicron (N)$个旋转因子的复数乘法。

这种方法以及FFT的基本思路在1965年J. W. Cooley和J. W. Tukey合作发表An algorithm for the machine calculation of complex Fourier series之后开始为人所知。但后来发现，实际上这两位作者只是重新发明了高斯在1805年就已经提出的算法（此算法在历史上数次以各种形式被再次提出）。

库利-图基算法最有名的应用，是将序列长为 $N$ 的DFT分割为两个长为 $N/2$ 的子序列的DFT，因此这一应用只适用于序列长度为2的幂的DFT计算，即基2-FFT。实际上，如同高斯和Cooley与Tukey都指出的那样，Cooley-Tukey算法也可以用于序列长度 $N$ 为任意因数分解形式的DFT，即混合基FFT，而且还可以应用于其他诸如分裂基FFT等变种。尽管Cooley-Tukey算法的基本思路是采用递归的方法进行计算，大多数传统的算法实现都将显式的递归算法改写为非递归的形式。另外，因为Cooley-Tukey算法是将DFT分解为较小长度的多个DFT，因此它可以同任一种其他的DFT算法联合使用。

### 设计思想

下面，我们用 $N$ 次单位根 $W_{N}$ 来表示 $e^{-j\frac{2\pi}{N}}$。

$W_{N}$ 的性质：

- 周期性，$W_{N}$具有周期$N$，即$W_{N}^{k+N}=W_{N}^{k}$
- 对称性：$W_{N}^{k+N/2}=-W_{N}^{k}$。
- 若$m$是$N$的约数，$W_{N}^{mkn}=W_{N/m}^{kn}$

为了简单起见，我们下面设待变换序列长度 $n=2^r$。根据上面单位根的对称性，求级数 $y_k=\sum_{n=0}^{N-1} W_{N}^{kn}x_n$ 时，可以将求和区间分为两部分：

$$\begin{matrix}y_k=\sum_{n=2t} W_{N}^{kn} x_n + \sum_{n=2t+1} W_{N}^{kn}x_n\\= \sum_{t} W_{\frac{N}{2}}^{kt}x_{2t} + W_{N}^{k}\sum_{t} W_{\frac{N}{2}}^{kt}x_{2t+1}\\= F_{even}(k) + W_{N}^{k}F_{odd}(k)&&&&&&(i\in\mathbb{Z})\end{matrix}$$

$F_{odd}(k)$ 和 $F_{even}(k)$ 是两个分别关于序列 $\left\{x_n\right\}_0^{N-1}$ 奇数号和偶数号序列N/2点变换。由此式只需计算出 $y_k$ 的前 $N/2$ 个点，对于后N/2个点可以通过复指数函数的对称性快速求解。即，注意 $F_{odd}(k)$ 和 $F_{even}(k)$ 都是周期为 $N/2$ 的函数，由单位根的对称性，有以下变换公式：

- $$y_{k+\frac{N}{2}} = F_{even}(k) - W_{N}^{k}F_{odd}(k)$$
- $y_k = F_{even}(k) + W_{N}^{k}F_{odd}(k)$。

这样，一个 $N$ 点变换就分解成了两个 $N/2$ 点变换。照这样可继续分解下去。这就是库利-图基快速傅里叶变换算法的基本原理。根据主定理不难分析出此时算法的时间复杂度为$\Omicron (N\log N)$.

### 算法实现

首先将$n=2^N$个输入点列按二进制进行编号，然后对各个编号按位倒置并按此重新排序。例如，对于一个8点变换。

001倒置以后变成 100

000 → 000

001 → 100

010 → 010

011 → 110

100 → 001

101 → 101

110 → 011

111 → 111

倒置后的编号为{0,4,2,6,1,5,3,7}。

然后将这n个点列作为输入传送到蝶形结网络中，注意将因子 $W_{N}^{k}$ 逐层加入到蝶形网络中。

### 算法复杂度

由于按蝶形结网络计算 $n$ 点变换要进行 $\log{} n$ 层计算，每层计算n个点的变换，故算法的时间复杂度为$\Omicron (n \log n)$。

## 其他算法

在Cooley-Tukey算法之外还有其他DFT的快速算法。对于长度 $N = N_1N_2$ 且 $N_1$ 与 $N_2$ 互质的序列，可以采用基于中国剩余定理的互质因子算法将N长序列的DFT分解为两个子序列的DFT。与Cooley-Tukey算法不同的是，互质因子算法不需要旋转因子。

Rader-Brenner算法是类似于Cooley-Tukey算法，但是采用的旋转因子都是纯虚数，以增加加法运算和降低了数值稳定性为代价减少了乘法运算。这方法之后被split-radix variant of Cooley-Tukey所取代，与Rader-Brenner算法相比，有一样多的乘法量，却有较少的加法量，且不牺牲数值的准确性。

Bruun以及QFT算法是不断的把DFT分成许多较小的DFT运算。（Rader-Brenner以及QFT算法是为了2的指数所设计的算法，但依然可以适用在可分解的整数上。Bruun算法则可以运用在可被分成偶数个运算的数字）。尤其是Bruun算法，把FFT看成是 $z^N-1$，并把它分解成 $z^{M-1}$与$z^{2M}+az^M+1$ 的形式。

另一个从多项式观点的快速傅立叶变换法是Winograd算法。此算法把 $z^N-1$ 分解成cyclotomic多项式，而这些多项式的系数通常为1，0，-1。这样只需要很少的乘法量（如果有需要的话），所以winograd是可以得到最少乘法量的快速傅立叶算法，对于较小的数字，可以找出有效率的算方式。更精确地说，winograd算法让DFT可以用 $2^k$ 点的DFT来简化，但减少乘法量的同时，也增加了非常多的加法量。Winograd也可以利用剩余值定理来简化DFT。

Rader算法提出了利用点数为N（N为质数）的DFT进行长度为N-1的回旋折积来表示原本的DFT，如此就可利用折积用一对基本的FFT来计算DFT。另一个prime-size的FFT算法为chirp-Z算法。此法也是将DFT用折积来表示，此法与Rader算法相比，能运用在更一般的转换上，其转换的基础为Z转换（Rabiner et al., 1969）。

## 实数或对称资料专用的算法

在许多的运用当中，要进行DFT的资料是纯实数，如此一来经过DFT的结果会满足对称性：

$\mathbf{X}_{N-k}$=$\mathbf{X}_k^*$

而有一些算法是专门为这种情形设计的（e.g. Sorensen, 1987）。另一些则是利用旧有的算法（e.g. Cooley-Tukey），删去一些不必要的演算步骤，如此省下了记忆体的使用，也省下了时间。另一方面，也可以把一个偶数长度且纯实数的DFT，用长度为原本一半的复数型态DFT来表示（实数项为原本纯实数资料的偶数项，虚数项则为奇数项）。

在Matlab上用一次N点FFT计算两个N点实数信号x,y的FFT:

```matlab
function [X,Y] = Real2FFT( x,y )

% N-point FFT is used for 2 real signals both with length N

N = length(x);

z = x+y*i ;

Z = fft(z);

X = 0.5*(Z+conj(Z(1+mod(0:-1:1-N,N))));

Y = 0.5*(Z-conj(Z(1+mod(0:-1:1-N,N))))/i;

end
```

一度人们认为，用离散哈特利转换（Discrete Hartley Transform）来处理纯实数的DFT会更有效率，但接着人们发现，对于同样点数的纯实数DFT，经过设计的FFT，可以比DHT省下更多的运算。Bruun算法是第一个试着从减少实数输入的DFT运算量的算法，但此法并没有成为人们普遍使用的方法。

对于实数输入，且输入为偶对称或奇对称的情形，可以更进一步的省下时间以及记忆体，此时DFT可以用离散余弦转换或离散正弦转换来代替（Discrete cosine/sine transforms）。由于DCT/DST也可以设计出FFT的算法，故在此种情形下，此方法取代了对DFT设计的FFT算法。

DFT可以应用在频谱分析以及做折积的运算，而在此处，不同应用可以用不同的算法来取代，列表如下：

用来做频谱分析的情况下，DFT可用下列的算法代替：

- DCT
- DST
- DHT
- 正交基底的扩展（orthogonal basis expantion）包括正交多项式（orthogonal polynomials）以及CDMA.
- Walsh（Hadamard）转换
- Haar转换
- 小波（wavelet）转换
- 时频分布（time-frequency distribution）

用来做折积的情况下，DFT可用下列的算法代替：

- DCT
- DST
- DHT
- 直接做折积（direct computing）
- 分段式DFT折积（sectioned DFT convolution）
- 威诺格拉德快速傅立叶变换算法
- 沃尔什（Walsh、Hadamard）转换
- 数论转换

## 单一路径延迟回授系统

单一路径延迟回授系统(英语：Single-path delay feedback，缩写为SDF)是一种用于实现快速傅立叶转换(Fast Fourier Transform)中运算的管线式架构(pipeline architecture)。此架构的特点是资料从输入一直到最终的输出只会通过单一的路径，并使用蝶形结(butterfly diagram)做为资料的运算单元。经过蝶形结后的输出资料会先暂存在系统的移位暂存器(shift registers)或是随机存取内存(RAM)当中，而由于只有单一路径的关系，对于每一级的蝶形结架构我们只需要一个复数乘法器来进行和FFT当中的旋转因子(twiddle factor)的乘积运算。

SDF的优点在于其架构简单并且利用硬件实现的成本较低，原因在于其单一路径的设计，使得资料会依序进行处理并输出，并不会在系统中加入多个复杂的运算单元以进行平行化的运算，事实上SDF的架构在管线式快速傅立叶转换处理器(pipelined FFT processors)中具有最高的内存使用效率，而由于内存单元数量会随快速傅立叶转换的级数而指数增长，内存使用效率进而成为影响电路面积和功耗的主要因素之一。因此在实作大型快速傅立叶转换处理器的情况下，SDF架构一直是被优先考量采用的选择。然而其缺点也同样来自于单一路径的设计，导致其吞吐量(throughput)和多路径延迟交换线路系统(英语：Multiple-path delay commutator，缩写为MDC)相比较低。

接着我们将以8个点的基-2单一路径延迟回授系统(8-point R2SDF)为例进行其算法细节的介绍，并且展示在硬件实作中每一个cycle的运算和中间结果的存储位置。

<image>

首先考虑系统输入为x0, x1, ..., x7，并在前8个cycle依序输入。

- Cycle 0~3: x0~x3依序存入第一级BF-2上方的暂存器中。
- Cycle 4: x0和x4输入第一级BF-2后输出x0'和x4'，x0'存入第二级BF-2上方的暂存器中，x4'存入第一级BF-2上方的暂存器中。
- Cycle 5: x1和x5输入第一级BF-2后输出x1'和x5'，x1'存入第二级BF-2上方的暂存器中，x5'存入第一级BF-2上方的暂存器中。
- Cycle 6: x2和x6输入第一级BF-2后输出x2'和x6'，x6'存入第一级BF-2上方的暂存器中，x0'和x2'输入第二级BF-2后输出x0`''`和x2`''`，x0`''`存入第三级BF-2上方的暂存器中，x2`''`存入第二级BF-2上方的暂存器中。
- Cycle 7: x3和x7输入第一级BF-2后输出x3'和x7'，x7'存入第一级BF-2上方的暂存器中，x1'和x3'输入第二级BF-2后输出x1`''`和x3`''`，x3`''`存入第二级BF-2上方的暂存器中，x0`''`和x1`''`输入第三级BF-2后输出x0`'''`和x1`'''`，x0`'''`直接输出，x1`'''`存入第三级BF-2上方的暂存器中。(此为第一笔输出，此后每一个cycle皆会输出一笔资料。)
- Cycle 8: x4'乘上W0后存入第二级BF-2上方的暂存器中，x2`''`乘上W0后存入第三级BF-2上方的暂存器中，x1`'''`直接输出。
- Cycle 9: x5'乘上W1后存入第二级BF-2上方的暂存器中，x3`''`乘上W2后和x2`''`输入第三级BF-2后输出x2`'''`和x3`'''`， x2`'''`直接输出，x3`'''`存入第三级BF-2上方的暂存器中。
- Cycle 10: x6'乘上W2后和x4'输入第二级BF-2后输出x4`''`和x6`''`，x6`''`存入第二级BF-2上方的暂存器中，x4`''`存入第三级BF-2上方的暂存器中，x3`'''`直接输出。
- Cycle 11: x7'乘上W3后和x5'输入第二级BF-2后输出x5`''`和x7`''`，x7`''`存入第二级BF-2上方的暂存器中，x4`''`和x5`''`输入第三级BF-2后输出x4`'''`和x5`'''`，x4`'''`直接输出，x5`'''`存入第三级BF-2上方的暂存器中。
- Cycle 12: x6`''`乘上W0后存入第三级BF-2上方的暂存器中，x5`'''`直接输出。
- Cycle 13: x7`''`乘上W2后和x6`''`输入第三级BF-2后输出x6`'''`和x7`'''`，x6`'''`直接输出，x7`'''`存入第三级BF-2上方的暂存器中。
- Cycle 14: x7`'''`直接输出。(此为最后一笔输出。)

此外，我们可以在蝶形结中的加法器和乘法器加上额外的暂存器(register)来降低系统的关键路径(critical path)，来避免某些时刻不同级之间的繁重运算必须挤在同一个cycle中运算完成，进而严重影响整体电路在timing上的表现。但如此的管线式设计(pipeline design)也会增加整体电路的延迟(latency)和暂存器使用量。

### 相似架构之硬件实作比较

- Radix-2 SDF: 2(log4N-1) 个乘法器、4log4N 个加法器、内存大小为 N-1。
- Radix-4 SDF: log4N-1 个乘法器、8log4N 个加法器、内存大小为 N-1。
- Radix-22 SDF: log4N-1 个乘法器、4log4N 个加法器、内存大小为 N-1。
- Radix-2 MDC: 2(log4N-1) 个乘法器、4log4N 个加法器、内存大小为 1.5N-2。
- Radix-4 MDC: 3(log4N-1) 个乘法器、8log4N 个加法器、内存大小为 2.5N-4。
- Radix-4 SDC: log4N-1 个乘法器、3log4N 个加法器、内存大小为 2N-2。

## 复杂度以及运算量的极限

长久以来，人们对于求出快速傅立叶变换的复杂度下限以及需要多少的运算量感到很有兴趣，而实际上也还有许多问题需要解决。即使是用较简单的情形，即 $2^k$ 点的DFT，也还没能够严谨的证明出FFT至少需要 $\Omega (N\log N)$ （不比 $N\log N$ 小）的运算量，目前也没有发现复杂度更低的算法。通常数学运算量的多寡会是运算效率好坏最主要的因素，但在现实中，有许多因素也会有很大的影响，如快取内存以及CPU均有很大的影响。

在1978年，Winograd率先导出一个较严谨的FFT所需乘法量的下限： $\Theta (N)$。当 $N=2^k$ 时，DFT只需要 $4N-2\log_{2}^2N-2\log_{2}N-4$ 次无理实数的乘法即可以计算出来。更详尽，且也能趋近此下限的算法也一一被提出（Heideman & Burrus, 1986; Duhamel, 1990）。很可惜的是，这些算法，都需要很大量的加法计算，目前的硬件无法克服这个问题。

对于所需加法量的数目，虽然我们可以在某些受限制的假设下，推得其下限，但目前并没有一个精确的下限被推导出来。1973年，Morgenstern在乘法常数趋近巨大的情形下（对大部分的FFT算法为真，但不是全部）推导出加法量的下限： $\Omega \left(N \log N \right)$。Pan（1986）在假设FFT算法的不同步的情形有其极限下证明出加法量的下限 $\Omega (NlogN)$，但一般来说，此假设相当的不明确。长度为 $N=2^k$ 的情形下，在某些假设下，Papadimitriou（1979）提出使用Cooley-Tukey算法所需的复数加法量 $N\log_{2}N$ 是最少的。到目前为止，在长度为 $N=2^k$ 情况，还没有任何FFT的算法可以让复数的加法量 比$N\log_{2}N$ 还少。

还有一个问题是如何把乘法量与加法量的总和最小化，有时候称作"演算复杂度"（在这里考虑的是实际的运算量，而不是渐近复杂度）。同样的，没有一个严谨下限被证明出来。从1968年开始，$N=2^k$点DFT而言，split-radix FFT算法需要最少的运算量，在 $N>1$ 的情形下，其需要 $4N\log_{2}N-6N+8$ 个乘法运算以及加法运算。最近有人导出更低的运算量： $\frac{34}{9}N\log_{2}N$。（Johnson and Frigo, 2007; Lundy and Van Buskirk, 2007）

大多数尝试要降低或者证明FFT复杂度下限的人都把焦点放在复数资料输入的情况，因其为最简单的情形。但是，复数资料输入的FFT算法，与实数资料输入的FFT算法，离散余弦转换（DCT），离散哈特列转换（DHT），以及其他的算法，均有很大的关连性。故任何一个算法，在复杂度上有任何的改善的话，其他的算法复杂度也会马上获得改善（Duhamel & Vetterli, 1990）。

## 快速算法查表

对于 $N$ 点的DFT:

$N$ = 1~60

N｜乘法数｜加法数｜N｜乘法数｜N｜乘法数｜N｜乘法数

1｜0｜0｜11｜46｜24｜28｜39｜182

2｜0｜4｜12｜8｜25｜148｜40｜100

3｜2｜12｜13｜52｜26｜104｜42｜124

4｜0｜16｜14｜32｜27｜114｜44｜160

5｜10｜34｜15｜40｜28｜64｜45｜170

6｜4｜36｜16｜20｜30｜80｜48｜92

7｜16｜72｜18｜32｜32｜72｜52｜208

8｜4｜52｜20｜40｜33｜160｜54｜228

9｜16｜72｜21｜62｜35｜150｜56｜156

10｜20｜88｜22｜80｜36｜64｜60｜160

$N$ < 1000

N｜乘法数｜N｜乘法数｜N｜乘法数｜N｜乘法数

63｜256｜96｜280｜192｜752｜360｜1540

64｜204｜104｜468｜204｜976｜420｜2080

66｜284｜108｜456｜216｜1020｜480｜2360

70｜300｜112｜396｜224｜1016｜504｜2300

72｜164｜120｜380｜240｜940｜512｜3180

80｜260｜128｜560｜252｜1024｜560｜3100

81｜480｜144｜436｜256｜1308｜672｜3496

84｜248｜160｜680｜288｜1160｜720｜3620

88｜412｜168｜580｜312｜1324｜784｜4412

90｜340｜180｜680｜336｜1412｜840｜4580

$N$ > 1000

N｜乘法数｜N｜乘法数｜N｜乘法数｜N｜乘法数

1008｜5356｜1440｜8680｜2520｜16540｜4032｜29488

1024｜7436｜1680｜10420｜2688｜19108｜4096｜37516

1152｜7088｜2016｜12728｜2880｜20060｜4368｜35828

1260｜7640｜2048｜16836｜3369｜24200｜4608｜36812

1344｜8252｜2304｜15868｜3920｜29900｜5040｜36860

其他 $N$ 的算法:

(1)若 $N$ 可被两个以上互质的数相乘得到

设 $N = P_1\times P_2$，且 $P_1$ 与 $P_2$ 互质(无须 $P_1$ 与 $P_2$ 为质数):

则当输入信号长度为N所需之的乘法器 $=P_1 \times B_2+P_2 \times B_1$。

其中 $B_1$ 为 $P_1$点DFT乘法量， $B_2$ 为 $P_2$ 点DFT乘法量。

例如: $N=65$ 可以拆成 $65=5\times 13$，其中5与13互质。

借由查表可得知5点 DFT乘法量为10，13点DFT乘法量为52。

所以65点DFT乘法量为 $5\times52 + 13\times10 = 390$。

扩展:

设 $N = P_1\times P_2\times\cdot\cdot\cdot\times P_K$，且 $P_1\times P_2\times\cdot\cdot\cdot\times P_K$ 互质。

$P_k$ 点DFT乘法量为 $B_k$。

则 $N$ 点DFT可分解为 $N/P_1$ 个$P_1$-point DFT+$N/P_2$ 个 $P_2$-point DFT+$\cdot\cdot\cdot$+$N/P_K$ 个$P_K$-point DFT。

总乘法量为 $\frac{N}{P_1}B_1+\frac{N}{P_2}B_2+ \cdots +\frac{N}{P_K}B_K$

(2)若N不可被两个互质的数相乘得到

设 $N = P_1\times P_2$，且 $P_1$ 与 $P_2$ 不互质:

且 $n \times m$ 当中( $n = 1,2,...,P_1-1$，$m = 1,2,...,P_2-1$ )。

有 $D_1$ 个值不为 $N/12$ 及 $N/8$ 的倍数。

有 $D_2$ 个值为 $N/12$ 或 $N/8$ 的倍数，但不为 $N/4$ 的倍数。

则 $N$ 点DFT所需之的乘法量 $=P_1 \times B_2+P_2 \times B_1+3D_1+2D_2$。

其中 $B_1$ 为 $P_1$ 点DFT乘法量， $B_2$ 为 $P_2$ 点DFT乘法量。

例如:16点DFT， $16=8 \times 2$

乘法量 $=2 \times 4 + 8 \times 0+3 \times 4+2 \times 2=24$

或 $16=4 \times 4$

乘法量 $=4 \times 0 + 4 \times 0 +3 \times 4+2 \times 4=20$

## FFT的应用

快速傅立叶变换可将卷积化乘法并以对数级数降低运算量，已成数字信号与科学计算核心工具。

其应用涵盖：

### 1.音讯与影像压缩

MP3、AAC 与 JPEG 均把讯号分块后先作 FFT／MDCT，以把时域样本转到频域；高能量系数集中于低频，小幅子载波便能重现重要讯息。接着结合人耳／视觉掩蔽模型及 Zig-zag 扫描，将低感知值系数量化为零，大幅降低位元率。最后用霍夫曼或算术编码存储系数与量化步阶，达到有损压缩与可调画质／音质的目的。若不经 FFT，逐点卷积将耗费 O(N²) 运算，无法在即时装置上完成。

### 2.无线通讯

现代蜂巢式网络与 Wi-Fi 系统皆采正交分频多工 (OFDM)。基站先将并行调变符号经 IFFT 合成等间隔子载波；终端则以 FFT 还原频域资料。此架构可在多径通道中保持子载波互不干扰，同时透过可变子载波间距支撑 5G NR 的 eMBB 与 mMTC 场景。运用循环字首可把时域折积化为频域乘法，等效频率选择通道简化为每载波一次复数乘法，计算量对比传统单载波大幅减少。

### 3.频谱分析与测试仪器

FFT 让示波器、功率分析仪与机房电力品质监控能在毫秒级显示功率频谱。将长度 N 的窗函数样本取 FFT 后，平方其复数幅度并适当标准化，即得单／双边功率谱密度 (PSD)。工程师可由此快速辨识谐波、噪声底或机械转速边带；科学领域亦用来分析地震波与星体光谱。无须昂贵的扫频式分析器，即可精确估测频率分量与相位。

### 4.图像处理与快速卷积

空间滤波、去模糊及模板匹配常利用“乘频域、回时域”策略：对影像与滤波器 PSF 取 FFT，相乘后作 IFFT，即得卷积结果，计算量 $O(N^2)$ 降为 $O(N\log{} N)$。在盲去模糊或维纳解卷积中，亦可用频域除法抑制模糊并保留高频细节。此方法支援任意遮罩大小，远优于滑动视窗直接卷积，适用医学影像与天文摄影等高分辨率场景。

### 5.雷达／声纳脉冲压缩

雷达为兼顾高解析与大探测距离，常发射长码调变脉冲并于接收端做匹配滤波。将回波与已知码序列分别取 FFT，相乘后 IFFT，即可在频域完成相关；此“脉冲压缩”技术能在固定峰值功率下获得窄等效脉冲宽，提升距离分辨率与信杂比，并大幅降低即时硬件乘加数。

### 6.多项式乘法与同态加密

在多项式乘法或卷积为瓶颈的应用（例如格基同态加密、FFT-based 大数计算）中，先将系数以 NTT/FFT 转到频域，再做逐点乘法，最后逆变换还原，即可把 $O(N^2)$ 乘法降至 $O(N\log{} N)$。加密函式库 SEAL、HElib 均内建此优化，使百万阶多项式运算在一般 CPU 上于毫秒内完成，直接左右现代隐私计算性能。

由此可见，FFT已由数学算法演变为支撑资讯社会不可或缺的基础技术。

## 参阅

- 离散傅里叶变换
- 并行快速傅里叶变换
- 快速数论变换