【数字信号处理】【傅里叶分析】【FFT】快速傅里叶变换的完整公式推导

给出离散傅里叶变换DFT的公式：
$$\begin{array}{rl}DFT[x(n)]=X(m)& =\sum _{n=0}^{N-1}x(n)exp(-j\frac{2\pi mn}{N}),m\in [0,N-1]\\ & =x(0)exp(-j\frac{2\pi m0}{N})+x(1)exp(-j\frac{2\pi m1}{N})+x(2)exp(-j\frac{2\pi m2}{N})+x(3)exp(-j\frac{2\pi m3}{N})+...\end{array}$$

FFT将DFT公式分成奇和偶两部分，即，将含$x(0),x(2),...$的项取出，作为偶数项，将含$x(1),x(3),...$的项取出，作为奇数项，有
$$\begin{array}{rl}X(m)& =\sum _{n=0}^{\frac{N}{2}-1}x(2n)exp(-j\frac{2\pi m2n}{N})+\sum _{n=0}^{\frac{N}{2}-1}x(2n+1)exp(-j\frac{2\pi m(2n+1)}{N}),m\in [0,N-1]\\ & =\sum _{n=0}^{\frac{N}{2}-1}x(2n)exp(-j\frac{2\pi mn}{\frac{N}{2}})+\exp (-j\frac{2\pi m}{N})\sum _{n=0}^{\frac{N}{2}-1}x(2n+1)\exp (-j\frac{2\pi mn}{\frac{N}{2}}),m\in [0,\frac{N}{2}-1]\end{array}$$
可以看到经过变换后，$m$的取值范围变成了一半。现在考虑另一半的m，有
$$\begin{array}{rl}X(m+\frac{N}{2})& =\sum _{n=0}^{\frac{N}{2}-1}x(2n)\exp (-j\frac{2\pi n(m+\frac{N}{2})}{\frac{N}{2}})+\exp (-j\frac{2\pi (m+\frac{N}{2})}{N})\sum _{n=0}^{\frac{N}{2}-1}x(2n+1)\exp (-j\frac{2\pi n(m+\frac{N}{2})}{\frac{N}{2}}),m\in [0,\frac{N}{2}-1]\\ & =\sum _{n=0}^{\frac{N}{2}-1}x(2n)\exp (-j\frac{2\pi mn}{\frac{N}{2}})\exp (-j2\pi n)+\exp (-j\frac{2\pi m}{N})\exp (-j\pi )\sum _{n=0}^{\frac{N}{2}-1}x(2n+1)\exp (-j\frac{2\pi nm}{\frac{N}{2}})\exp (-j2\pi n),m\in [0,\frac{N}{2}-1]\\ & =\sum _{n=0}^{\frac{N}{2}-1}x(2n)exp(-j\frac{2\pi mn}{\frac{N}{2}})-\exp (-j\frac{2\pi m}{N})\sum _{n=0}^{\frac{N}{2}-1}x(2n+1)\exp (-j\frac{2\pi mn}{\frac{N}{2}}),m\in [0,\frac{N}{2}-1]\end{array}$$
令
$$W_N^{\theta }=\exp (-j\frac{2\pi \theta }{N})$$则
$$\begin{gathered} X(m)=\sum _{n=0}^{\frac{N}{2}-1}x(2n)W_{\frac{N}{2}}^{mn}+W_N^m\sum _{n=0}^{\frac{N}{2}-1}x(2n+1)W_{\frac{N}{2}}^{mn},m\in [0,\frac{N}{2}-1] \\ X(m+\frac{N}{2})=\sum _{n=0}^{\frac{N}{2}-1}x(2n)W_{\frac{N}{2}}^{mn}-W_N^m\sum _{n=0}^{\frac{N}{2}-1}x(2n+1)W_{\frac{N}{2}}^{mn},m\in [0,\frac{N}{2}-1] \end{gathered}$$
注意到
$$\begin{gathered} \sum _{n=0}^{\frac{N}{2}-1}x(2n)W_{\frac{N}{2}}^{mn}① \\ \sum _{n=0}^{\frac{N}{2}-1}x(2n+1)W_{\frac{N}{2}}^{mn}② \end{gathered}$$
①和②分别是$N$点序列$x(n)$中对偶数项与奇数项的$\frac{N}{2}$点的DFT。

上述两个式子$X(m)$和$X(m+\frac{N}{2})$为FFT的核心公式，即一个$N$点DFT可以分成2个$\frac{N}{2}$点的DFT，并且前$\frac{N}{2}$点DFT与后$\frac{N}{2}$点DFT之间，只有正负号的不同，所以可以使用前$\frac{N}{2}$点DFT公式中的元素，计算后$\frac{N}{2}$点DFT中的元素，使计算量降低。

根据该思路，可以将①和②继续分解，分别求①中奇数项与偶数项的$\frac{N}{4}$点DFT，和②中奇数项和偶数项的$\frac{N}{4}$点DFT。这样一直分解下去，直到分解成最简单的2点DFT为止。注意到，在FFT中，N点DFT会一分为二，二分为四，四分为八，最后会分解成多个2点DFT，所以N必须是2的幂次。

接下来，本文用一个$N=8$点FFT进行说明。

首先写出8点DFT原始公式
$$DFT[x(n)]=X(m)=\sum _{n=0}^{7}x(n)exp(-j\frac{2\pi mn}{N}),m\in [0,7]$$
根据上文的FFT核心公式，8点DFT可以分为两个4点DFT，即
$$\begin{gathered} X(m)=\sum _{n=0}^{3}x(2n)W_4^{mn}+W_8^m\sum _{n=0}^{3}x(2n+1)W_4^{mn},m\in [0,3] \\ X(m+4)=\sum _{n=0}^{3}x(2n)W_4^{mn}-W_8^m\sum _{n=0}^3(2n+1)W_4^{mn},m\in [0,3] \end{gathered}$$
令
$$\begin{gathered} X_1(m)=\sum _{n=0}^{3}x(2n)W_4^{mn} \\ {X}_2(m)=\sum _{n=0}^{3}x(2n+1)W_4^{mn} \end{gathered}$$
则
$$\begin{gathered} X(m)=X_1(m)+W_8^m{X}_2(m),m\in [0,3] \\ X(m+4)=X_1(m)-W_8^m{X}_2(m),m\in [0,3] \end{gathered}$$
1、对$X_1(m)$进行奇偶分解，有
$$\begin{array}{rl}{X}_1(m)& =\sum _{n=0}^{1}x(2(2n))W_4^{m2n}+\sum _{n=0}^{1}x(2(2n+1))W_4^{m(2n+1)},m\in [0,3]\\ & =\sum _{n=0}^{1}x(4n)\exp (-\frac{2\pi m2n}{4})+W_4^m\sum _{n=0}^{1}x(4n+2)\exp (-j\frac{2\pi m2n}{4}),m\in [0,3]\\ & =\sum _{n=0}^{1}x(4n)W_2^{mn}+W_4^m\sum _{n=0}^{1}x(4n+2)W_2^{mn},m\in [0,1]\\ & =x(0)W_2^0+x(4)W_2^m+W_4^m[x(2)W_2^0+x(6)W_2^m],m\in [0,1]\\ & =x(0)+x(4)W_2^m+W_4^m[x(2)+x(6)W_2^m],m\in [0,1]\end{array}$$
即
$$\begin{gathered} X_1(0)=x(0)+x(4)+x(2)+x(6) \\ {X}_1(1)=x(0)-x(4)+W_4^1[x(2)-x(6)] \end{gathered}$$
考虑$m+2\in [2,3]$的情况
$$\begin{array}{rl}{X}_1(m+2)& =\sum _{n=0}^{1}x(4n)W_2^{(m+2)n}+W_4^{m+2}\sum _{n=0}^{1}x(4n+2)W_2^{(m+2)n},m\in [0,1]\\ & =\sum _{n=0}^{1}x(4n)W_2^{mn}\exp (-j2\pi n)+W_4^m\exp (-j\pi )\sum _{n=0}^{1}x(4n+2)W_2^{mn}\exp (-j2\pi n),m\in [0,1]\\ & =\sum _{n=0}^{1}x(4n)W_2^{mn}-W_4^m\sum _{n=0}^{1}x(4n+2)W_2^{(m+2)n},m\in [0,1]\\ & =x(0)W_2^0+x(4)W_2^m-W_4^m[x(2)W_2^0+x(6)W_2^m],m\in [0,1]\\ & =x(0)+x(4)W_2^m-W_4^m[x(2)+x(6)W_2^m],m\in [0,1]\end{array}$$
即
$$\begin{gathered} X_1(2)=x(0)+x(4)-[x(2)+x(6)] \\ {X}_1(3)=x(0)-x(4)-W_4^1[x(2)-x(6)] \end{gathered}$$
2、对$X_2(m)$进行奇偶分解，有
$$\begin{array}{rl}{X}_2(m)& =\sum _{n=0}^{1}x(2(2n)+1)W_4^{m2n}+\sum _{n=0}^{1}x(2(2n+1)+1)W_4^{m(2n+1)},m\in [0,3]\\ & =\sum _{n=0}^{1}x(4n+1)\exp (-\frac{2\pi m2n}{4})+W_4^m\sum _{n=0}^{1}x(4n+3)\exp (-j\frac{2\pi m2n}{4}),m\in [0,3]\\ & =\sum _{n=0}^{1}x(4n+1)W_2^{mn}+W_4^m\sum _{n=0}^{1}x(4n+3)W_2^{mn},m\in [0,1]\\ & =x(1)W_2^0+x(5)W_2^m+W_4^m[x(3)W_2^0+x(7)W_2^m],m\in [0,1]\\ & =x(1)+x(5)W_2^m+W_4^m[x(3)+x(7)W_2^m],m\in [0,1]\end{array}$$
即
$$\begin{gathered} X_2(0)=x(1)+x(5)+x(3)+x(7) \\ {X}_2(1)=x(1)-x(5)+W_4^1[x(3)-x(7)] \end{gathered}$$
考虑$m+2\in [2,3]$的情况
$$\begin{array}{rl}{X}_2(m+2)& =\sum _{n=0}^{1}x(4n+1)W_2^{(m+2)n}+W_4^{m+2}\sum _{n=0}^{1}x(4n+3)W_2^{(m+2)n},m\in [0,1]\\ & =\sum _{n=0}^{1}x(4n+1)W_2^{mn}\exp (-j2\pi n)+W_4^m\exp (-j\pi )\sum _{n=0}^{1}x(4n+3)W_2^{mn}\exp (-j2\pi n),m\in [0,1]\\ & =\sum _{n=0}^{1}x(4n+1)W_2^{mn}-W_4^m\sum _{n=0}^{1}x(4n+3)W_2^{mn},m\in [0,1]\\ & =x(1)W_2^0+x(5)W_2^m-W_4^m[x(3)W_2^0+x(7)W_2^m],m\in [0,1]\\ & =x(1)+x(5)W_2^m-W_4^m[x(3)+x(7)W_2^m],m\in [0,1]\end{array}$$
即
$$\begin{gathered} X_2(2)=x(1)+x(5)-[x(3)+x(7)] \\ {X}_2(3)=x(1)-x(5)-W_4^1[x(3)-x(7)] \end{gathered}$$
将$X_1(m)$和$X_2(m)$分解后的表达式代回$X(m)$和$X(m+4)$中，有
$$\begin{array}{rl}X(m)& =X_1(m)+W_8^m{X}_2(m)\\ X(m+2)& =X_1(m+2)+W_8^{m+2}{X}_2(m+2)\\ X(m+4)& =X_1(m)-W_8^m{X}_2(m)\\ X((m+2)+4)& =X_1(m+2)-W_8^{m+2}{X}_2(m+2)\\ m& \in [0,1]\end{array}$$
将$m\in [0,1]$代入，得到$X(0),X(1),...,X(7)$的表达式
X ( 0 ) = x ( 0 ) + x ( 4 ) + x ( 2 ) + x ( 6 ) + x ( 1 ) + x ( 5 ) + x ( 3 ) + x ( 7 ) X ( 1 ) = x ( 0 ) − x ( 4 ) + W 4 1 [ x ( 2 ) − x ( 6 ) ] + W 8 1 { x ( 1 ) − x ( 5 ) + W 4 1 [ x ( 3 ) − x ( 7 ) ] } X ( 2 ) = x ( 0 ) + x ( 4 ) − [ x ( 2 ) + x ( 6 ) ] + W 8 2 { x ( 1 ) + x ( 5 ) − [ x ( 3 ) + x ( 7 ) ] } X ( 3 ) = x ( 0 ) − x ( 4 ) − W 4 1 [ x ( 2 ) − x ( 6 ) ] + W 8 3 { x ( 1 ) − x ( 5 ) − W 4 1 [ x ( 3 ) − x ( 7 ) ] } X ( 4 ) = x ( 0 ) + x ( 4 ) + x ( 2 ) + x ( 6 ) − [ x ( 1 ) + x ( 5 ) + x ( 3 ) + x ( 7 ) ] X ( 5 ) = x ( 0 ) − x ( 4 ) + W 4 1 [ x ( 2 ) − x ( 6 ) ] − W 8 1 { x ( 1 ) − x ( 5 ) + W 4 1 [ x ( 3 ) − x ( 7 ) ] } X ( 6 ) = x ( 0 ) + x ( 4 ) − [ x ( 2 ) + x ( 6 ) ] − W 8 2 { x ( 1 ) + x ( 5 ) − [ x ( 3 ) + x ( 7 ) ] } X ( 7 ) = x ( 0 ) − x ( 4 ) − W 4 1 [ x ( 2 ) − x ( 6 ) ] − W 8 3 { x ( 1 ) − x ( 5 ) − W 4 1 [ x ( 3 ) − x ( 7 ) ] } \begin{aligned} X(0)&=x(0)+x(4)+x(2)+x(6)+x(1)+x(5)+x(3)+x(7)\\ X(1)&=x(0)-x(4)+W_4^1[x(2)-x(6)]+W_8^1\{x(1)-x(5)+W_4^1[x(3)-x(7)]\}\\ X(2)&=x(0)+x(4)-[x(2)+x(6)]+W_8^2\{x(1)+x(5)-[x(3)+x(7)]\}\\ X(3)&=x(0)-x(4)-W_4^1[x(2)-x(6)]+W_8^3\{x(1)-x(5)-W_4^1[x(3)-x(7)]\}\\ X(4)&=x(0)+x(4)+x(2)+x(6)-[x(1)+x(5)+x(3)+x(7)]\\ X(5)&=x(0)-x(4)+W_4^1[x(2)-x(6)]-W_8^1\{x(1)-x(5)+W_4^1[x(3)-x(7)]\}\\ X(6)&=x(0)+x(4)-[x(2)+x(6)]-W_8^2\{x(1)+x(5)-[x(3)+x(7)]\}\\ X(7)&=x(0)-x(4)-W_4^1[x(2)-x(6)]-W_8^3\{x(1)-x(5)-W_4^1[x(3)-x(7)]\}\\ \end{aligned} 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(4)+W41​[x(2)−x(6)]+W81​{x(1)−x(5)+W41​[x(3)−x(7)]}=x(0)+x(4)−[x(2)+x(6)]+W82​{x(1)+x(5)−[x(3)+x(7)]}=x(0)−x(4)−W41​[x(2)−x(6)]+W83​{x(1)−x(5)−W41​[x(3)−x(7)]}=x(0)+x(4)+x(2)+x(6)−[x(1)+x(5)+x(3)+x(7)]=x(0)−x(4)+W41​[x(2)−x(6)]−W81​{x(1)−x(5)+W41​[x(3)−x(7)]}=x(0)+x(4)−[x(2)+x(6)]−W82​{x(1)+x(5)−[x(3)+x(7)]}=x(0)−x(4)−W41​[x(2)−x(6)]−W83​{x(1)−x(5)−W41​[x(3)−x(7)]}​
即
X ( 0 ) = x ( 0 ) + x ( 4 ) + x ( 2 ) + x ( 6 ) + x ( 1 ) + x ( 5 ) + x ( 3 ) + x ( 7 ) X ( 1 ) = x ( 0 ) − x ( 4 ) − j [ x ( 2 ) − x ( 6 ) ] + exp ⁡ ( − j π 4 ) { x ( 1 ) − x ( 5 ) − j [ x ( 3 ) − x ( 7 ) ] } X ( 2 ) = x ( 0 ) + x ( 4 ) − [ x ( 2 ) + x ( 6 ) ] − j { x ( 1 ) + x ( 5 ) − [ x ( 3 ) + x ( 7 ) ] } X ( 3 ) = x ( 0 ) − x ( 4 ) + j [ x ( 2 ) − x ( 6 ) ] + exp ⁡ ( − j 3 π 4 ) { x ( 1 ) − x ( 5 ) + j [ x ( 3 ) − x ( 7 ) ] } X ( 4 ) = x ( 0 ) + x ( 4 ) + x ( 2 ) + x ( 6 ) − [ x ( 1 ) + x ( 5 ) + x ( 3 ) + x ( 7 ) ] X ( 5 ) = x ( 0 ) − x ( 4 ) − j [ x ( 2 ) − x ( 6 ) ] − exp ⁡ ( − j π 4 ) { x ( 1 ) − x ( 5 ) − j [ x ( 3 ) − x ( 7 ) ] } X ( 6 ) = x ( 0 ) + x ( 4 ) − [ x ( 2 ) + x ( 6 ) ] + j { x ( 1 ) + x ( 5 ) − [ x ( 3 ) + x ( 7 ) ] } X ( 7 ) = x ( 0 ) − x ( 4 ) + j [ x ( 2 ) − x ( 6 ) ] − exp ⁡ ( − j 3 π 4 ) { x ( 1 ) − x ( 5 ) + j [ x ( 3 ) − x ( 7 ) ] } \begin{aligned} X(0)&=x(0)+x(4)+x(2)+x(6)+x(1)+x(5)+x(3)+x(7)\\ X(1)&=x(0)-x(4)-j[x(2)-x(6)]+\exp(-j\frac{\pi}4)\{x(1)-x(5)-j[x(3)-x(7)]\}\\ X(2)&=x(0)+x(4)-[x(2)+x(6)]-j\{x(1)+x(5)-[x(3)+x(7)]\}\\ X(3)&=x(0)-x(4)+j[x(2)-x(6)]+\exp(-j\frac{3\pi}4)\{x(1)-x(5)+j[x(3)-x(7)]\}\\ X(4)&=x(0)+x(4)+x(2)+x(6)-[x(1)+x(5)+x(3)+x(7)]\\ X(5)&=x(0)-x(4)-j[x(2)-x(6)]-\exp(-j\frac{\pi}4)\{x(1)-x(5)-j[x(3)-x(7)]\}\\ X(6)&=x(0)+x(4)-[x(2)+x(6)]+j\{x(1)+x(5)-[x(3)+x(7)]\}\\ X(7)&=x(0)-x(4)+j[x(2)-x(6)]-\exp(-j\frac{3\pi}4)\{x(1)-x(5)+j[x(3)-x(7)]\}\\ \end{aligned} 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(4)−j[x(2)−x(6)]+exp(−j4π​){x(1)−x(5)−j[x(3)−x(7)]}=x(0)+x(4)−[x(2)+x(6)]−j{x(1)+x(5)−[x(3)+x(7)]}=x(0)−x(4)+j[x(2)−x(6)]+exp(−j43π​){x(1)−x(5)+j[x(3)−x(7)]}=x(0)+x(4)+x(2)+x(6)−[x(1)+x(5)+x(3)+x(7)]=x(0)−x(4)−j[x(2)−x(6)]−exp(−j4π​){x(1)−x(5)−j[x(3)−x(7)]}=x(0)+x(4)−[x(2)+x(6)]+j{x(1)+x(5)−[x(3)+x(7)]}=x(0)−x(4)+j[x(2)−x(6)]−exp(−j43π​){x(1)−x(5)+j[x(3)−x(7)]}​