库利-图基快速傅里叶变换算法 (英语:Cooley–Tukey FFT algorithm )[ 1] 是最常见的快速傅里叶变换 算法。这一方法以分治法 为策略递归 地将长度为N = N 1 N 2 的DFT分解为长度分别为N 1 和N 2 的两个较短序列的DFT,以及与旋转因子的复数乘法。这种方法以及FFT的基本思路在1965年詹姆斯·库利 和约翰·图基 合作发表《An algorithm for the machine calculation of complex Fourier series 》之后开始为人所知。但后来发现,实际上这两位作者只是重新发明了高斯 在1805年就已经提出的算法(此算法在历史上数次以各种形式被再次提出)。
库利-图基算法最有名的应用,是将序列长为N的DFT分割为两个长为N/2的子序列的DFT,因此这一应用只适用于序列长度为2的幂的DFT计算,即基2-FFT。实际上,如同高斯和库利与图基都指出的那样,库利-图基算法也可以用于序列长度N为任意因数分解形式的DFT,即混合基FFT,而且还可以应用于其他诸如分裂基FFT等变种。尽管库利-图基算法的基本思路是采用递归的方法进行计算,大多数传统的算法实现都将显示的递归算法改写为非递归的形式。另外,因为库利-图基算法是将DFT分解为较小长度的多个DFT,因此它可以同任一种其他的DFT算法联合使用。
在FFT演算法中,针对输入做不同方式的分组会造成输出顺序上的不同。如果我们使用时域抽取(Decimation-in-time),那么输入的顺序将会是位元反转排列(bit-reversed order),输出将会依序排列。但若我们采取的是频域抽取(Decimation-in-frequency),那么输出与输出顺序的情况将会完全相反,变为依序排列的输入与位元反转排列的输出。
我们可以将DFT公式中的项目在时域上重新分组,这样就叫做时域抽取(Decimation-in-time),在这里
e
−
j
(
2
π
n
k
N
)
{\displaystyle e^{-j(2\pi {\frac {nk}{N}})}}
将会被代换为旋转因子 (twiddle factor)
W
N
n
k
{\displaystyle W_{N}^{nk}}
。
X
[
k
]
=
∑
n
=
0
N
−
1
x
[
n
]
e
−
j
(
2
π
n
k
N
)
k
=
0
,
1
,
…
,
N
−
1
{\displaystyle X[k]=\sum _{n=0}^{N-1}x[n]e^{-j(2\pi {\frac {nk}{N}})}\qquad k=0,1,\ldots ,N-1}
=
∑
n
e
v
e
n
x
[
n
]
W
N
n
k
+
∑
n
o
d
d
x
[
n
]
W
N
n
k
{\displaystyle =\sum _{n\ even}x[n]W_{N}^{nk}+\sum _{n\ odd}x[n]W_{N}^{nk}}
=
∑
r
=
0
(
N
/
2
)
−
1
x
[
2
r
]
W
N
2
r
k
+
∑
r
=
0
(
N
/
2
)
−
1
x
[
2
r
+
1
]
W
N
(
2
r
+
1
)
k
{\displaystyle =\sum _{r=0}^{(N/2)-1}x[2r]W_{N}^{2rk}+\sum _{r=0}^{(N/2)-1}x[2r+1]W_{N}^{(2r+1)k}}
=
∑
r
=
0
(
N
/
2
)
−
1
x
[
2
r
]
W
N
2
r
k
+
W
N
k
∑
r
=
0
(
N
/
2
)
−
1
x
[
2
r
+
1
]
W
N
2
r
k
{\displaystyle =\sum _{r=0}^{(N/2)-1}x[2r]W_{N}^{2rk}+W_{N}^{k}\sum _{r=0}^{(N/2)-1}x[2r+1]W_{N}^{2rk}}
=
∑
r
=
0
(
N
/
2
)
−
1
x
[
2
r
]
W
N
/
2
r
k
+
W
N
k
∑
r
=
0
(
N
/
2
)
−
1
x
[
2
r
+
1
]
W
N
/
2
r
k
{\displaystyle =\sum _{r=0}^{(N/2)-1}x[2r]W_{N/2}^{rk}+W_{N}^{k}\sum _{r=0}^{(N/2)-1}x[2r+1]W_{N/2}^{rk}}
=
G
[
k
]
+
W
N
k
H
[
k
]
{\displaystyle =G[k]+W_{N}^{k}H[k]}
在这边我们要注意的是,我们所替换的G[k]与H[k]具有周期性:
{
G
[
k
+
N
2
]
=
G
[
k
]
H
[
k
+
N
2
]
=
H
[
k
]
{\displaystyle {\begin{cases}G[k+{\frac {N}{2}}]=G[k]\\H[k+{\frac {N}{2}}]=H[k]\end{cases}}}
还注意到系数具有对称性:
W
N
k
+
N
/
2
=
−
W
N
k
{\displaystyle W_{N}^{k+N/2}=-W_{N}^{k}}
上述的推导可以划成下面的图:
划红框处所示的
N
2
{\displaystyle {\frac {N}{2}}}
点DFT架构如下图所示:
划红框处所示的
N
4
{\displaystyle {\frac {N}{4}}}
点DFT架构如下图所示:
下图是一个8点DIT FFT的完整架构图。
我们可以将DFT公式中的项目在频域上重新分组,这样就叫做频域抽取(Decimation-in-frequency)。
首先先观察在频域上偶数项的部分:
X
[
2
r
]
=
∑
n
=
0
N
−
1
x
[
n
]
W
N
n
(
2
r
)
r
=
0
,
1
,
⋯
,
N
2
−
1
{\displaystyle X[2r]=\sum _{n=0}^{N-1}x[n]W_{N}^{n(2r)}\ r=0,1,\cdots ,{\frac {N}{2}}-1}
=
∑
n
=
0
N
2
−
1
x
[
n
]
W
N
2
n
r
+
∑
n
=
N
2
N
−
1
x
[
n
]
W
N
2
n
r
{\displaystyle =\sum _{n=0}^{{\frac {N}{2}}-1}x[n]W_{N}^{2nr}+\sum _{n={\frac {N}{2}}}^{N-1}x[n]W_{N}^{2nr}}
=
∑
n
=
0
N
2
−
1
x
[
n
]
W
N
2
n
r
+
∑
n
=
0
N
2
−
1
x
[
n
+
N
2
]
W
N
2
r
[
n
+
N
2
]
{\displaystyle =\sum _{n=0}^{{\frac {N}{2}}-1}x[n]W_{N}^{2nr}+\sum _{n=0}^{{\frac {N}{2}}-1}x[n+{\frac {N}{2}}]W_{N}^{2r[n+{\frac {N}{2}}]}}
∵
W
N
2
r
[
n
+
N
2
]
=
W
N
2
r
n
W
N
r
N
=
W
N
2
r
n
{\displaystyle {\color {Gray}\because W_{N}^{2r[n+{\frac {N}{2}}]}=W_{N}^{2rn}W_{N}^{rN}=W_{N}^{2rn}}}
=
∑
n
=
0
N
2
−
1
x
[
n
]
W
N
2
r
n
+
∑
n
=
0
N
2
−
1
x
[
n
+
N
2
]
W
N
2
r
n
{\displaystyle =\sum _{n=0}^{{\frac {N}{2}}-1}x[n]W_{N}^{2rn}+\sum _{n=0}^{{\frac {N}{2}}-1}x[n+{\frac {N}{2}}]W_{N}^{2rn}}
=
∑
n
=
0
N
2
−
1
(
x
[
n
]
+
x
[
n
+
N
2
]
)
W
N
2
r
n
{\displaystyle =\sum _{n=0}^{{\frac {N}{2}}-1}(x[n]+x[n+{\frac {N}{2}}])W_{\frac {N}{2}}^{rn}}
=
∑
n
=
0
N
2
−
1
g
[
n
]
W
N
2
r
n
{\displaystyle =\sum _{n=0}^{{\frac {N}{2}}-1}g[n]W_{\frac {N}{2}}^{rn}}
再观察在频域上奇数项的部分:
X
[
2
r
+
1
]
=
∑
n
=
0
N
−
1
x
[
n
]
W
N
n
(
2
r
+
1
)
r
=
0
,
1
,
⋯
,
N
2
−
1
{\displaystyle X[2r+1]=\sum _{n=0}^{N-1}x[n]W_{N}^{n(2r+1)}\ r=0,1,\cdots ,{\frac {N}{2}}-1}
=
∑
n
=
0
N
2
−
1
x
[
n
]
W
N
n
(
2
r
+
1
)
+
∑
n
=
N
2
N
−
1
x
[
n
]
W
N
n
(
2
r
+
1
)
{\displaystyle =\sum _{n=0}^{{\frac {N}{2}}-1}x[n]W_{N}^{n(2r+1)}+\sum _{n={\frac {N}{2}}}^{N-1}x[n]W_{N}^{n(2r+1)}}
=
∑
n
=
0
N
2
−
1
x
[
n
]
W
N
n
(
2
r
+
1
)
+
∑
n
=
0
N
2
−
1
x
[
n
+
N
2
]
W
N
(
2
r
+
1
)
[
n
+
N
2
]
{\displaystyle =\sum _{n=0}^{{\frac {N}{2}}-1}x[n]W_{N}^{n(2r+1)}+\sum _{n=0}^{{\frac {N}{2}}-1}x[n+{\frac {N}{2}}]W_{N}^{(2r+1)[n+{\frac {N}{2}}]}}
∵
W
N
(
2
r
+
1
)
[
n
+
N
2
]
=
W
N
(
2
r
+
1
)
n
W
N
(
2
r
+
1
)
N
2
=
−
W
N
(
2
r
+
1
)
n
{\displaystyle {\color {Gray}\because W_{N}^{(2r+1)[n+{\frac {N}{2}}]}=W_{N}^{(2r+1)n}W_{N}^{(2r+1){\frac {N}{2}}}=-W_{N}^{(2r+1)n}}}
=
∑
n
=
0
N
2
−
1
x
[
n
]
W
N
(
2
r
+
1
)
n
−
∑
n
=
0
N
2
−
1
x
[
n
+
N
2
]
W
N
(
2
r
+
1
)
n
{\displaystyle =\sum _{n=0}^{{\frac {N}{2}}-1}x[n]W_{N}^{(2r+1)n}-\sum _{n=0}^{{\frac {N}{2}}-1}x[n+{\frac {N}{2}}]W_{N}^{(2r+1)n}}
=
∑
n
=
0
N
2
−
1
(
x
[
n
]
−
x
[
n
+
N
2
]
)
W
N
n
(
2
r
+
1
)
{\displaystyle =\sum _{n=0}^{{\frac {N}{2}}-1}(x[n]-x[n+{\frac {N}{2}}])W_{N}^{n(2r+1)}}
=
∑
n
=
0
N
2
−
1
(
x
[
n
]
−
x
[
n
+
N
2
]
)
W
N
n
W
N
2
n
r
{\displaystyle =\sum _{n=0}^{{\frac {N}{2}}-1}(x[n]-x[n+{\frac {N}{2}}])W_{N}^{n}W_{\frac {N}{2}}^{nr}}
=
∑
n
=
0
N
2
−
1
(
h
[
n
]
W
N
n
)
W
N
2
r
n
{\displaystyle =\sum _{n=0}^{{\frac {N}{2}}-1}(h[n]W_{N}^{n})W_{\frac {N}{2}}^{rn}}
上述的推导可以画成下面的图:
更进一步的拆解
N
2
{\displaystyle {\frac {N}{2}}}
-point DFT的架构
下图为8点FFT下
N
4
{\displaystyle {\frac {N}{4}}}
-point DFT的架构
总结上述架构,完整的8点DIF FFT架构图为
利用库利-图基演算法进行离散傅立叶 拆解时,能够依需求而以2, 4, 8…等2的幂次方为单位进行拆解,不同的拆解方法会产生不同层数快速傅里叶变换 的架构,基底越大则层数越少,复数乘法器也越少,但是每级的蝴蝶形架构则会越复杂,因此常见的架构为2基底、4基底与8基底这三种设计。
2基底-快速傅立叶演算法(Radix-2 FFT)是最广为人知的一种库利-图基快速傅立叶演算法分支。此方法不断地将N点的FFT拆解成两个'N/2'点的FFT,利用旋转因子
W
n
k
,
W
N
2
{\displaystyle W^{nk},W^{\frac {N}{2}}}
的对称性借此来降低DFT的计算复杂度,达到加速的功效。
而其实前述有关时域抽取或是频域抽取的方法介绍,即为2基底的快速傅立叶转换法。以下展示其他种2基底快速傅立叶演算法的连线方法,此种不规则的连线方法可以让输出与输入都为循序排列,但是连线的复杂度却大大的增加。
4基底快速傅立叶变换演算法则是承接2基底的概念,在此里用时域 抽取的技巧,将原本的DFT公式拆解为四个一组的形式:
X
[
k
]
=
∑
n
=
0
N
−
1
x
[
n
]
e
−
j
(
2
π
n
k
N
)
k
=
0
,
1
,
…
,
N
−
1
{\displaystyle X[k]=\sum _{n=0}^{N-1}x[n]e^{-j(2\pi {\frac {nk}{N}})}\qquad k=0,1,\ldots ,N-1}
=
∑
n
=
0
N
4
−
1
x
[
4
n
+
0
]
W
N
(
4
n
+
0
)
k
+
∑
n
=
0
N
4
−
1
x
[
4
n
+
1
]
W
N
(
4
n
+
1
)
k
{\displaystyle =\sum _{n=0}^{{\frac {N}{4}}-1}x[4n+0]W_{N}^{(4n+0)k}+\sum _{n=0}^{{\frac {N}{4}}-1}x[4n+1]W_{N}^{(4n+1)k}}
+
∑
n
=
0
N
4
−
1
x
[
4
n
+
2
]
W
N
(
4
n
+
2
)
k
+
∑
n
=
0
N
4
−
1
x
[
4
n
+
3
]
W
N
(
4
n
+
3
)
k
{\displaystyle +\sum _{n=0}^{{\frac {N}{4}}-1}x[4n+2]W_{N}^{(4n+2)k}+\sum _{n=0}^{{\frac {N}{4}}-1}x[4n+3]W_{N}^{(4n+3)k}}
=
W
N
0
∑
n
=
0
N
4
−
1
x
[
4
n
+
0
]
W
N
4
n
k
+
W
N
1
k
∑
n
=
0
N
4
−
1
x
[
4
n
+
1
]
W
N
4
n
k
{\displaystyle =W_{N}^{0}\sum _{n=0}^{{\frac {N}{4}}-1}x[4n+0]W_{\frac {N}{4}}^{nk}+W_{N}^{1k}\sum _{n=0}^{{\frac {N}{4}}-1}x[4n+1]W_{\frac {N}{4}}^{nk}}
+
W
N
2
k
∑
n
=
0
N
4
−
1
x
[
4
n
+
2
]
W
N
4
n
k
+
W
N
3
k
∑
n
=
0
N
4
−
1
x
[
4
n
+
3
]
W
N
4
n
k
{\displaystyle +W_{N}^{2k}\sum _{n=0}^{{\frac {N}{4}}-1}x[4n+2]W_{\frac {N}{4}}^{nk}+W_{N}^{3k}\sum _{n=0}^{{\frac {N}{4}}-1}x[4n+3]W_{\frac {N}{4}}^{nk}}
W
N
0
F
0
(
k
)
+
W
N
k
F
1
(
k
)
+
W
N
2
k
F
2
(
k
)
+
W
N
3
k
F
3
(
k
)
{\displaystyle W_{N}^{0}F_{0}(k)+W_{N}^{k}F_{1}(k)+W_{N}^{2k}F_{2}(k)+W_{N}^{3k}F_{3}(k)}
在这里再利用
W
n
k
+
N
4
=
−
W
n
k
+
3
N
4
=
−
j
W
n
k
{\displaystyle W^{nk+{\frac {N}{4}}}=-W^{nk+{\frac {3N}{4}}}=-jW^{nk}}
的特性来进行与2基数FFT类似的化减方法,以降低演算法复杂度。
在库利-图基算法里,使用的基底(radix)越大,复数的乘法与记忆体的存取就越少,所带来的好处不言而喻。但是随之而来的就是实数的乘法与实数的加法也会增加,尤其计算单元的设计也就越复杂,对于可应用FFT之点数限制也就越严格。在表中我们可以看到在不同基底下所需的计算成本。
使用4096点FFT在不同基底下的计算量
动作
2基底
4基底
8基底
复数乘法
22528
15360
10752
实数乘法
0
0
8192
复数加法
49152
49152
49152
实数加法
0
0
8192
记忆体存取
49152
24576
16384
在DFT的公式中,我们重新定义x=[x(0),x(1),…,x(N-1)]T , X=[X(0),X(1),…,X(N-1)]T ,则DFT可重写为X=FN ‧x。FN 是一个N×N的DFT矩阵,其元素定义为[FN ]ij =WNij (i,j ∈ [0,N-1]),当N=8时,我们可以得到以下的F8 矩阵并且进一步将其拆解。
在拆解成三个矩阵相乘之后,我们可以将复数运算的数量从56个降至24个复数的加法。底下是8基底的SFG。要注意的是所有的输出与输入都是复数的形式,而输出与输入的排序也并非规律排列,此种排列方式是为了要达到接线的最佳化。
在2/8基底的演算法中,我们可以看到我们对于偶数项的输出会使用2基底的分解法,对于奇数项的输出采用8基底的分解法。这个做法充分利用了2基底与4基底拥有最少乘法数与加法数的特性。当使用了2基底的分解法后,偶数项的输出如下所示。
C
2
k
=
∑
n
=
0
N
2
−
1
(
x
2
n
+
x
N
2
+
n
)
W
2
n
k
{\displaystyle C_{2k}=\sum _{n=0}^{{\frac {N}{2}}-1}(x_{2n}+x_{{\frac {N}{2}}+n})W^{2nk}}
奇数项的输出则交由8基底分解来处理,如下四式所述。
C
8
k
+
1
=
∑
n
=
0
N
8
−
1
{
[
(
x
n
−
x
n
+
N
2
)
−
j
(
x
n
+
N
4
−
x
n
+
3
N
4
)
]
+
W
N
8
[
(
x
n
+
N
8
−
x
n
+
5
N
8
)
−
j
(
x
n
+
3
N
8
−
x
n
+
7
N
8
)
]
}
W
8
n
k
W
n
{\displaystyle C_{8k+1}=\sum _{n=0}^{{\frac {N}{8}}-1}{\begin{Bmatrix}[(x_{n}-x_{n+{\frac {N}{2}}})-j(x_{n+{\frac {N}{4}}}-x_{n+{\frac {3N}{4}}})]+W^{\frac {N}{8}}[(x_{n+{\frac {N}{8}}}-x_{n+{\frac {5N}{8}}})-j(x_{n+{\frac {3N}{8}}}-x_{n+{\frac {7N}{8}}})]\end{Bmatrix}}W^{8nk}W^{n}}
C
8
k
+
5
=
∑
n
=
0
N
8
−
1
{
[
(
x
n
−
x
n
+
N
2
)
−
j
(
x
n
+
N
4
−
x
n
+
3
N
4
)
]
−
W
N
8
[
(
x
n
+
N
8
−
x
n
+
5
N
8
)
−
j
(
x
n
+
3
N
8
−
x
n
+
7
N
8
)
]
}
W
8
n
k
W
5
n
{\displaystyle C_{8k+5}=\sum _{n=0}^{{\frac {N}{8}}-1}{\begin{Bmatrix}[(x_{n}-x_{n+{\frac {N}{2}}})-j(x_{n+{\frac {N}{4}}}-x_{n+{\frac {3N}{4}}})]-W^{\frac {N}{8}}[(x_{n+{\frac {N}{8}}}-x_{n+{\frac {5N}{8}}})-j(x_{n+{\frac {3N}{8}}}-x_{n+{\frac {7N}{8}}})]\end{Bmatrix}}W^{8nk}W^{5n}}
C
8
k
+
3
=
∑
n
=
0
N
8
−
1
{
[
(
x
n
−
x
n
+
N
2
)
+
j
(
x
n
+
N
4
−
x
n
+
3
N
4
)
]
+
W
3
N
8
[
(
x
n
+
N
8
−
x
n
+
5
N
8
)
+
j
(
x
n
+
3
N
8
−
x
n
+
7
N
8
)
]
}
W
8
n
k
W
3
n
{\displaystyle C_{8k+3}=\sum _{n=0}^{{\frac {N}{8}}-1}{\begin{Bmatrix}[(x_{n}-x_{n+{\frac {N}{2}}})+j(x_{n+{\frac {N}{4}}}-x_{n+{\frac {3N}{4}}})]+W^{\frac {3N}{8}}[(x_{n+{\frac {N}{8}}}-x_{n+{\frac {5N}{8}}})+j(x_{n+{\frac {3N}{8}}}-x_{n+{\frac {7N}{8}}})]\end{Bmatrix}}W^{8nk}W^{3n}}
C
8
k
+
7
=
∑
n
=
0
N
8
−
1
{
[
(
x
n
−
x
n
+
N
2
)
+
j
(
x
n
+
N
4
−
x
n
+
3
N
4
)
]
−
W
3
N
8
[
(
x
n
+
N
8
−
x
n
+
5
N
8
)
+
j
(
x
n
+
3
N
8
−
x
n
+
7
N
8
)
]
}
W
8
n
k
W
7
n
{\displaystyle C_{8k+7}=\sum _{n=0}^{{\frac {N}{8}}-1}{\begin{Bmatrix}[(x_{n}-x_{n+{\frac {N}{2}}})+j(x_{n+{\frac {N}{4}}}-x_{n+{\frac {3N}{4}}})]-W^{\frac {3N}{8}}[(x_{n+{\frac {N}{8}}}-x_{n+{\frac {5N}{8}}})+j(x_{n+{\frac {3N}{8}}}-x_{n+{\frac {7N}{8}}})]\end{Bmatrix}}W^{8nk}W^{7n}}
以上式子就是2/8基底的FFT快速演算法。在架构图上可化为L型的蝴蝶运算架构,如下图所示。
而下图表示的是一个64点的FFT使用2/8基底的架构图。虽然2/8基底的演算法缩减了
W
N
8
,
W
3
N
8
{\displaystyle W^{\frac {N}{8}},W^{\frac {3N}{8}}}
的乘法量,但是这种演算法最大的缺点是跟其他固定基底或是混合基底比较起来,他的架构较为不规则。所以在硬体上比4基底或是2基底更难实现。
为了改进Radix 2/8在架构上的不规则性,我们在这里做了一些修改,如下图4.。此修改可让架构更加规则,且所使用的加法器与乘法器数量更加减少,如下表所示。
8n 点FFT在不同演算法下所需复数乘法量
N=8n
2基底
4基底
2/4混合基底
2/4/8基底
8
2
-
2
0
64
98
76
72
48
512
1538
-
1082
824
4096
18434
13996
12774
10168
在这里我们最小的运算单元称为PE(Process Element),PE内部包含了2/8基底、2/4基底、2基底的运算,简化过的信号处理流程与蝴蝶型架构图可见下图
基底的选择越大会造成蝴蝶形架构更加复杂,控制电路也会复杂化。因此Shousheng He和Mats Torkelson在1996提出了一个2^2基底的FFT演算法,利用旋转因子的特性:
W
N
4
=
−
j
{\displaystyle W_{\frac {N}{4}}=-j}
。而–j的乘法基本上只需要将被乘数的实部虚部对调,然后将虚部加上负号即可,这样的负数乘法被定义为'简单乘法',因此可以用很简单的硬体架构来实现。利用上面的特性,22 基底FFT能用串接的方式将旋转因子以4为单位对DFT公式进行拆解,将蝴蝶形架构层数降到log4N,大幅减少复数乘法器的用量,而同时却维持了2基底蝴蝶形架构的简单性,在效能上获得改进。22 基底DIF FFT演算法的拆解方法如下列公式所述:
N点DFT公式:
X
(
k
)
=
∑
n
=
0
N
−
1
x
(
n
)
W
N
n
k
,
0
≦
k
<
N
{\displaystyle X(k)=\sum _{n=0}^{N-1}x(n)W_{N}^{nk},0\leqq k<N}
利用线性映射将n与k映射到三个维度上面
{
n
=<
N
2
n
1
+
N
4
n
2
+
n
3
>
N
k
=<
k
1
+
2
k
2
+
4
k
3
>
N
{\displaystyle {\begin{cases}n=<{\frac {N}{2}}n_{1}+{\frac {N}{4}}n_{2}+n_{3}>_{N}\\k=<k_{1}+2k_{2}+4k_{3}>_{N}\end{cases}}}
然后套用Common Factor Algorithm(CFA)
X
(
k
1
+
2
k
2
+
4
k
3
)
=
∑
n
3
=
0
N
4
−
1
∑
n
2
=
0
1
∑
n
1
=
0
1
x
(
N
2
n
1
+
N
4
n
2
+
n
3
)
W
N
(
N
2
n
1
+
N
4
n
2
+
n
3
)
(
k
1
+
2
k
2
+
4
k
3
)
{\displaystyle X(k_{1}+2k_{2}+4k_{3})=\sum _{n_{3}=0}^{{\frac {N}{4}}-1}\sum _{n_{2}=0}^{1}\sum _{n_{1}=0}^{1}x({\frac {N}{2}}n_{1}+{\frac {N}{4}}n_{2}+n_{3})W_{N}^{({\frac {N}{2}}n_{1}+{\frac {N}{4}}n_{2}+n_{3})(k_{1}+2k_{2}+4k_{3})}}
=
∑
n
3
=
0
N
4
−
1
∑
n
2
=
0
1
{
B
N
2
k
1
W
N
(
N
4
n
2
+
n
3
)
k
1
}
W
N
(
N
4
n
2
+
n
3
)
(
2
k
2
+
4
k
3
)
{\displaystyle =\sum _{n_{3}=0}^{{\frac {N}{4}}-1}\sum _{n_{2}=0}^{1}{\begin{Bmatrix}B_{\frac {N}{2}}^{k_{1}}W_{N}^{({\frac {N}{4}}n_{2}+n_{3})k_{1}}\end{Bmatrix}}W_{N}^{({\frac {N}{4}}n_{2}+n_{3})(2k_{2}+4k_{3})}}
而蝴蝶型架构会变成以下形式
B
N
2
k
1
=
x
(
N
4
n
2
+
n
3
)
+
(
−
1
)
k
1
x
(
N
4
n
2
+
n
3
+
N
2
)
{\displaystyle B_{\frac {N}{2}}^{k_{1}}=x({\frac {N}{4}}n_{2}+n_{3})+(-1)^{k_{1}}x({\frac {N}{4}}n_{2}+n_{3}+{\frac {N}{2}})}
利用旋转因子
W
N
4
=
−
j
{\displaystyle W_{\frac {N}{4}}=-j}
的特性,可以观察出
W
N
(
N
4
n
2
+
n
3
)
(
k
1
+
2
k
2
+
4
k
3
)
=
W
N
N
n
2
n
3
W
N
N
4
n
2
(
k
1
+
2
k
2
)
W
N
n
3
(
k
1
+
2
k
2
)
W
N
4
n
3
k
3
{\displaystyle W_{N}^{({\frac {N}{4}}n_{2}+n_{3})(k_{1}+2k_{2}+4k_{3})}=W_{N}^{Nn_{2}n_{3}}W_{N}^{{\frac {N}{4}}n_{2}(k_{1}+2k_{2})}W_{N}^{n_{3}(k_{1}+2k_{2})}W_{N}^{4n_{3}k_{3}}}
=
(
−
j
)
n
2
(
k
1
+
2
k
2
)
W
N
n
3
(
k
1
+
2
k
2
)
W
N
4
n
3
k
3
{\displaystyle =(-j)^{n_{2}(k_{1}+2k_{2})}W_{N}^{n_{3}(k_{1}+2k_{2})}W_{N}^{4n_{3}k_{3}}}
再将此公式带入原式中可以得到
X
(
k
1
+
2
k
2
+
4
k
3
)
=
∑
n
3
=
0
N
4
−
1
[
H
(
k
1
,
k
2
,
n
3
)
W
N
n
3
(
k
1
+
2
k
2
)
]
W
N
4
n
3
k
3
{\displaystyle X(k_{1}+2k_{2}+4k_{3})=\sum _{n_{3}=0}^{{\frac {N}{4}}-1}[H(k_{1},k_{2},n_{3})W_{N}^{n_{3}(k_{1}+2k_{2})}]W_{N}^{4n_{3}k_{3}}}
H
(
k
1
,
k
2
,
n
3
)
=
B
F
2
I
[
x
(
n
3
)
+
(
−
1
)
k
1
(
n
3
+
N
2
)
]
⏞
+
(
−
j
)
k
1
+
2
k
2
B
F
2
I
[
x
(
n
3
+
N
4
)
+
(
−
1
)
k
1
(
n
3
+
3
N
4
)
]
⏞
⏟
B
F
2
I
I
{\displaystyle H(k_{1},k_{2},n_{3})={\begin{matrix}\underbrace {{\begin{matrix}BF2I\\\overbrace {[x(n_{3})+(-1)^{k_{1}}(n_{3}+{\frac {N}{2}})]} \end{matrix}}{\begin{matrix}\\\\+(-j)^{k_{1}+2k_{2}}\end{matrix}}{\begin{matrix}BF2I\\\overbrace {[x(n_{3}+{\frac {N}{4}})+(-1)^{k_{1}}(n_{3}+{\frac {3N}{4}})]} \end{matrix}}} \\BF2II\end{matrix}}}
如上述公式所示,原本的DFT公式会被拆解成多个
H
(
k
1
,
k
2
,
n
3
)
{\displaystyle H(k_{1},k_{2},n_{3})}
,而
H
(
k
1
,
k
2
,
n
3
)
{\displaystyle H(k_{1},k_{2},n_{3})}
又可分为BF2I与BF2II两个阶层,分别会对应到之后所介绍的两种硬体架构。
一个16点的DFT公式经过一次上面所述之拆解之后可得下面的FFT架构
可以看出上图的架构保留了2基底的简单架构,然而复数乘法却降到每两级才出现一次,也就是
l
o
g
4
N
{\displaystyle log_{4}N}
次。而BF2I以及BF2II所对应的硬体架构下图:
其中BF2II硬体单元中左下角的交叉电路就是用来处理-j的乘法。
一个256点的FFT架构可以由下面的硬体来实现:
其中图下方的为一7位元宽的计数器,而此架构的控制电路相当单纯,只要将计数器的各个位元分别接上BF2I与BF2II单元即可。
下表将2基底、4基底与22 基底演算法做一比较,可以发现22 基底演算法所需要的乘法气数量为2基底的一半,加法弃用量是4基底的一半,并维持一样的记忆体用量和控制电路的简单性。
乘法器与加法器数量比较
乘法数
加法数
记忆体大小
控制电路
R2SDF
2(log4 N-1)
4log4 N
N-1
简单
R4SDF
log4 N -1
8log4 N
N-1
中等
R22 SDF
log4 N -1
4log4 N
N-1
简单
如上所述,22 演算法是将旋转因子
W
N
4
=
−
j
{\displaystyle W^{\frac {N}{4}}=-j}
视为一个简单乘法,进而由公式以及硬体上的化简获得硬体需求上的改进。而借由相同的概念,Shousheng He和Mats Torkelson进一步将旋转因子
W
N
8
=
2
2
(
1
−
j
)
{\displaystyle W^{\frac {N}{8}}={\frac {\sqrt {2}}{2}}(1-j)}
的乘法化简成一个简单乘法,而化简的方法将会在下面讲解。
2
2
{\displaystyle {\frac {\sqrt {2}}{2}}}
乘法化简
编辑
在2基数FFT演算法中的基本概念是利用旋转因子
W
n
k
,
W
N
2
{\displaystyle W^{nk},W^{\frac {N}{2}}}
的对称性,4基数演算法则是利用
W
n
k
+
N
4
=
−
W
n
k
+
3
N
4
=
−
j
W
n
k
{\displaystyle W^{nk+{\frac {N}{4}}}=-W^{nk+{\frac {3N}{4}}}=-jW^{nk}}
的特性。但是我们会发现在这些旋转因子 的对称特性中─
W
N
8
=
−
W
5
N
8
=
2
2
(
1
−
j
)
,
W
3
N
8
=
−
W
7
N
8
=
−
2
2
(
1
+
j
)
{\displaystyle W^{\frac {N}{8}}=-W^{\frac {5N}{8}}={\frac {\sqrt {2}}{2}}(1-j),W^{\frac {3N}{8}}=-W^{\frac {7N}{8}}=-{\frac {\sqrt {2}}{2}}(1+j)}
─并没有被利用到。主要是因为
2
2
(
1
−
j
)
{\displaystyle {\frac {\sqrt {2}}{2}}(1-j)}
的乘法运算会让整个FFT变得复杂,但是如果借由近似的方法,我们便可以将此一运算化简为12个加法。
2
2
=
0.70710678
=
2
−
1
+
2
−
3
+
2
−
4
+
2
−
6
+
2
−
8
+
2
−
9
{\displaystyle {\frac {\sqrt {2}}{2}}=0.70710678=2^{-1}+2^{-3}+2^{-4}+2^{-6}+2^{-8}+2^{-9}}
我们可以从上式注意到,
2
2
{\displaystyle {\frac {\sqrt {2}}{2}}}
可以被近似为五个加法的结果,所以
2
2
(
1
+
j
)
{\displaystyle {\frac {\sqrt {2}}{2}}(1+j)}
就可以被简化为只有六个加法,再算入实部与虚部的计算,总共只需12个加法器就可以运用到此一简化特性。
经由与22 基底类似的推导,并用串接的方式将旋转因子 以8为单位对DFT公式进行拆解,23 基底FFT演算法进一步将复数乘法器的用量缩减到log8 N,并同时维持硬体架构的简单性。
推导的方法与22 基底相当类似。借由这样的方法,23 基底能将乘法器的用量缩减到2基底的1/3,并同时维持一样的记忆体用量以及控制电路的简单性。
除了常在应用中见到与
2
{\displaystyle 2}
相关基底的拆解法,对于更加一般性的
N
=
N
1
N
2
{\displaystyle N=N_{1}N_{2}}
点离散傅立叶变换 问题,
我们也有办法在理论上进行拆解,将问题化为数个
N
1
{\displaystyle N_{1}}
与
N
2
{\displaystyle N_{2}}
点离散傅立叶变换 问题,并可对计算量进行估计。
而特别的是,透过互质 在数论 上的特性,对于
N
1
{\displaystyle N_{1}}
与
N
2
{\displaystyle N_{2}}
互质的情况,可以进一步节省一些运算,
在下面会特别分开讨论。
N
1
N
2
{\displaystyle N_{1}N_{2}}
非互质
编辑
为了避免之后的符号混淆,我们先将
N
1
N
2
{\displaystyle N_{1}N_{2}}
置换为
P
1
P
2
{\displaystyle P_{1}P_{2}}
,也就是说接著要将
N
=
P
1
P
2
{\displaystyle N=P_{1}P_{2}}
点离散傅立叶变换 ,
想办法拆解为数个
P
1
{\displaystyle P_{1}}
与
P
2
{\displaystyle P_{2}}
点离散傅立叶变换 问题。
接著定义要拆分的问题,要拆分的问题为
N
{\displaystyle N}
点离散傅立叶变换 ,将
f
[
n
]
{\displaystyle f[n]}
转换至
F
[
m
]
{\displaystyle F[m]}
:
F
[
m
]
=
∑
n
=
0
N
−
1
e
−
i
2
π
N
m
n
f
[
n
]
m
,
n
=
0
,
1
,
…
,
N
−
1
{\displaystyle F[m]=\sum _{n=0}^{N-1}e^{-i{\frac {2\pi }{N}}mn}f[n]\qquad m,n=0,1,\ldots ,N-1}
直观地说,这个
N
{\displaystyle N}
点离散傅立叶变换 ,将由
n
{\displaystyle n}
作为参数的函数
f
[
n
]
{\displaystyle f[n]}
,转换成由
m
{\displaystyle m}
作为参数的函数
F
[
m
]
{\displaystyle F[m]}
,
并且
m
,
n
{\displaystyle m,n}
都有
N
{\displaystyle N}
个可能的数值。
待定义好要拆分的问题,便可以开始讨论如何进行拆分,基本概念是将有
N
{\displaystyle N}
个可能的数值的
m
,
n
{\displaystyle m,n}
,
分别化为个使用两个参数进行描述的函数
m
=
[
m
1
,
m
2
]
,
n
=
[
n
1
,
n
2
]
{\displaystyle m=[m_{1},m_{2}],n=[n_{1},n_{2}]}
,并借此将原问题化为二维度离散傅立叶变换 ,
便可使用数个较小的离散傅立叶变换 问题描述整个过程。
而一种很直觉的转换方式,便是透过将
m
,
n
{\displaystyle m,n}
分别除以
P
2
,
P
1
{\displaystyle P_{2},P_{1}}
,
以商数 与馀数 来做为参数描述
m
,
n
{\displaystyle m,n}
的值:
n
=
n
1
P
1
+
n
2
m
=
m
1
+
m
2
P
2
{\displaystyle n=n_{1}P_{1}+n_{2}\qquad m=m_{1}+m_{2}P_{2}}
n
1
,
m
1
=
0
,
1
,
…
,
P
2
−
1
n
2
,
m
2
=
0
,
1
,
…
,
P
1
−
1
{\displaystyle n_{1},m_{1}=0,1,\ldots ,P_{2}-1\qquad n_{2},m_{2}=0,1,\ldots ,P_{1}-1}
其中
n
1
{\displaystyle n_{1}}
作为将
n
{\displaystyle n}
除以
P
1
{\displaystyle P_{1}}
的商数,与作为
m
{\displaystyle m}
除以
P
2
{\displaystyle P_{2}}
的馀数的
m
1
{\displaystyle m_{1}}
相同,
具有
P
2
{\displaystyle P_{2}}
个可能数值,同理
n
2
{\displaystyle n_{2}}
与
m
2
{\displaystyle m_{2}}
有
P
1
{\displaystyle P_{1}}
个可能数值。
将上述的参数代换及
N
=
P
1
P
2
{\displaystyle N=P_{1}P_{2}}
带入原式,可以得到:
F
[
m
1
+
m
2
P
2
]
=
∑
n
1
=
0
P
2
−
1
∑
n
2
=
0
P
1
−
1
e
−
i
2
π
P
1
P
2
(
m
1
+
m
2
P
2
)
(
n
1
P
1
+
n
2
)
f
[
n
1
P
1
+
n
2
]
{\displaystyle F[m_{1}+m_{2}P_{2}]=\sum _{n_{1}=0}^{P_{2}-1}\sum _{n_{2}=0}^{P_{1}-1}e^{-i{\frac {2\pi }{P_{1}P_{2}}}(m_{1}+m_{2}P_{2})(n_{1}P_{1}+n_{2})}f[n_{1}P_{1}+n_{2}]}
将右式的指数部分乘开并分项化简可以得到:
F
[
m
1
+
m
2
P
2
]
=
∑
n
1
=
0
P
2
−
1
∑
n
2
=
0
P
1
−
1
f
[
n
1
P
1
+
n
2
]
e
−
i
2
π
P
2
m
1
n
1
e
−
i
2
π
P
1
m
2
n
2
e
−
i
2
π
P
1
P
2
m
1
n
2
e
−
i
2
π
m
2
n
1
{\displaystyle F[m_{1}+m_{2}P_{2}]=\sum _{n_{1}=0}^{P_{2}-1}\sum _{n_{2}=0}^{P_{1}-1}f[n_{1}P_{1}+n_{2}]e^{-i{\frac {2\pi }{P_{2}}}m_{1}n_{1}}e^{-i{\frac {2\pi }{P_{1}}}m_{2}n_{2}}e^{-i{\frac {2\pi }{P_{1}P_{2}}}m_{1}n_{2}}e^{-i2\pi m_{2}n_{1}}}
最后透过
e
−
i
2
π
=
1
{\displaystyle e^{-i2\pi }=1}
与
P
1
P
2
=
N
{\displaystyle P_{1}P_{2}=N}
,可以得到:
F
[
m
1
+
m
2
P
2
]
=
∑
n
2
=
0
P
1
−
1
∑
n
1
=
0
P
2
−
1
f
[
n
1
P
1
+
n
2
]
e
−
i
2
π
P
2
m
1
n
1
e
−
i
2
π
N
m
1
n
2
e
−
i
2
π
P
1
m
2
n
2
{\displaystyle F[m_{1}+m_{2}P_{2}]=\sum _{n_{2}=0}^{P_{1}-1}\sum _{n_{1}=0}^{P_{2}-1}f[n_{1}P_{1}+n_{2}]e^{-i{\frac {2\pi }{P_{2}}}m_{1}n_{1}}e^{-i{\frac {2\pi }{N}}m_{1}n_{2}}e^{-i{\frac {2\pi }{P_{1}}}m_{2}n_{2}}}
观察上式,并加上括号辅助厘清运算顺序我们可以得到:
F
[
m
1
+
m
2
P
2
]
=
∑
n
2
=
0
P
1
−
1
{
{
∑
n
1
=
0
P
2
−
1
f
[
n
1
P
1
+
n
2
]
e
−
i
2
π
P
2
m
1
n
1
}
e
−
i
2
π
N
m
1
n
2
}
e
−
i
2
π
P
1
m
2
n
2
{\displaystyle F[m_{1}+m_{2}P_{2}]=\sum _{n_{2}=0}^{P_{1}-1}\left\{\left\{\sum _{n_{1}=0}^{P_{2}-1}f[n_{1}P_{1}+n_{2}]e^{-i{\frac {2\pi }{P_{2}}}m_{1}n_{1}}\right\}e^{-i{\frac {2\pi }{N}}m_{1}n_{2}}\right\}e^{-i{\frac {2\pi }{P_{1}}}m_{2}n_{2}}}
最内层的运算可以视为将双参数函数
f
[
n
1
,
n
2
]
{\displaystyle f[n_{1},n_{2}]}
中的一个参数
n
1
{\displaystyle n_{1}}
,透过离散傅立叶变换 取代为由
m
1
{\displaystyle m_{1}}
描述,
得到一个新函数
G
1
[
m
1
,
n
2
]
{\displaystyle G_{1}[m_{1},n_{2}]}
(这步因为对每个不同
n
2
{\displaystyle n_{2}}
,都需要做一次将
n
1
{\displaystyle n_{1}}
取代为
m
1
{\displaystyle m_{1}}
的转换,
共需要
P
1
{\displaystyle P_{1}}
个
P
2
{\displaystyle P_{2}}
点离散傅立叶变换 ):
F
[
m
1
+
m
2
P
2
]
=
∑
n
2
=
0
P
1
−
1
{
G
1
[
m
1
,
n
2
]
e
−
i
2
π
N
m
1
n
2
}
e
−
i
2
π
P
1
m
2
n
2
{\displaystyle F[m_{1}+m_{2}P_{2}]=\sum _{n_{2}=0}^{P_{1}-1}\left\{G_{1}[m_{1},n_{2}]e^{-i{\frac {2\pi }{N}}m_{1}n_{2}}\right\}e^{-i{\frac {2\pi }{P_{1}}}m_{2}n_{2}}}
下一层的运算则可视为单纯的乘法,将
G
1
[
m
1
,
n
2
]
{\displaystyle G_{1}[m_{1},n_{2}]}
与
e
−
i
2
π
N
m
1
n
2
{\displaystyle e^{-i{\frac {2\pi }{N}}m_{1}n_{2}}}
相乘,得到
G
2
[
m
1
,
n
2
]
{\displaystyle G_{2}[m_{1},n_{2}]}
(这步需要的计算量视
m
1
n
2
N
{\displaystyle {\frac {m_{1}n_{2}}{N}}}
的特殊性而会有变化):
F
[
m
1
+
m
2
P
2
]
=
∑
n
2
=
0
P
1
−
1
G
2
[
m
1
,
n
2
]
e
−
i
2
π
P
1
m
2
n
2
{\displaystyle F[m_{1}+m_{2}P_{2}]=\sum _{n_{2}=0}^{P_{1}-1}G_{2}[m_{1},n_{2}]e^{-i{\frac {2\pi }{P_{1}}}m_{2}n_{2}}}
最后的运算可以视为将函数
G
2
[
m
1
,
n
2
]
{\displaystyle G_{2}[m_{1},n_{2}]}
中
n
2
{\displaystyle n_{2}}
,透过离散傅立叶变换 取代为由
m
2
{\displaystyle m_{2}}
描述,
得到一个新函数
G
3
[
m
1
,
m
2
]
{\displaystyle G_{3}[m_{1},m_{2}]}
(这步因为对每个不同
m
1
{\displaystyle m_{1}}
,都需要做一次将
n
2
{\displaystyle n_{2}}
取代为
m
2
{\displaystyle m_{2}}
的转换,
共需要
P
2
{\displaystyle P_{2}}
个
P
1
{\displaystyle P_{1}}
点离散傅立叶变换 ):
F
[
m
(
=
m
1
+
m
2
P
2
)
]
=
G
3
[
m
1
,
m
2
]
{\displaystyle F[m(=m_{1}+m_{2}P_{2})]=G_{3}[m_{1},m_{2}]}
就成功仅使用
P
1
{\displaystyle P_{1}}
与
P
2
{\displaystyle P_{2}}
点离散傅立叶变换 ,描述了原先的
N
{\displaystyle N}
点离散傅立叶变换 。
而在这样的分解下,我们使用了
P
1
{\displaystyle P_{1}}
个
P
2
{\displaystyle P_{2}}
点离散傅立叶变换 与
P
2
{\displaystyle P_{2}}
个
P
1
{\displaystyle P_{1}}
点离散傅立叶变换 与一些额外的乘法,
并且这些额外使用的复数 乘法
G
1
[
m
1
,
n
2
]
×
e
−
i
2
π
N
m
1
n
2
{\displaystyle G_{1}[m_{1},n_{2}]\times e^{-i{\frac {2\pi }{N}}m_{1}n_{2}}}
,
在电脑的运算架构中,如果
m
1
n
2
N
{\displaystyle {\frac {m_{1}n_{2}}{N}}}
是
1
4
{\displaystyle {\frac {1}{4}}}
的倍数则不需要使用乘法,
如果
m
1
n
2
N
{\displaystyle {\frac {m_{1}n_{2}}{N}}}
是
1
8
,
1
12
{\displaystyle {\frac {1}{8}},{\frac {1}{12}}}
的倍数则仅需两个实数 乘法,
其他则需三个实数乘法,所以总运算量可以如下方式表示:
P
2
B
1
+
P
1
B
2
+
3
D
3
+
2
D
2
{\displaystyle P_{2}B_{1}+P_{1}B_{2}+3D_{3}+2D_{2}}
其中
B
1
{\displaystyle B_{1}}
是
P
1
{\displaystyle P_{1}}
傅立叶所需乘法数,
B
2
{\displaystyle B_{2}}
是
P
2
{\displaystyle P_{2}}
傅立叶所需乘法数,
D
3
{\displaystyle D_{3}}
是需三个实数乘法
m
1
n
2
{\displaystyle m_{1}n_{2}}
组合个数,
D
2
{\displaystyle D_{2}}
是需两个实数乘法
m
1
n
2
{\displaystyle m_{1}n_{2}}
组合个数。
而常见以
2
{\displaystyle 2}
为基底的分解则是为了使离散傅立叶变换 所需乘法数为零,这样就仅需考虑上面提到的额外乘法,可以提高效率也有较简单的结构。
N
1
N
2
{\displaystyle N_{1}N_{2}}
互质
编辑
在
N
1
N
2
{\displaystyle N_{1}N_{2}}
互质 的情况下,仍是采取和上面相近的思路来将问题进行拆分,首先,为了避免之后的符号混淆,我们同样将
N
1
N
2
{\displaystyle N_{1}N_{2}}
置换为
P
1
P
2
{\displaystyle P_{1}P_{2}}
。
接著同样定义要拆分的问题:
F
[
m
]
=
∑
n
=
0
N
−
1
e
−
i
2
π
N
m
n
f
[
n
]
m
,
n
=
0
,
1
,
…
,
N
−
1
{\displaystyle F[m]=\sum _{n=0}^{N-1}e^{-i{\frac {2\pi }{N}}mn}f[n]\qquad m,n=0,1,\ldots ,N-1}
接著就是和上面的算法有最大差异的部分,在将
m
,
n
{\displaystyle m,n}
化为个使用两个参数进行描述的函数
m
=
[
m
1
,
m
2
]
,
n
=
[
n
1
,
n
2
]
{\displaystyle m=[m_{1},m_{2}],n=[n_{1},n_{2}]}
时,
最直觉的作法便是使用商数和馀数,但在
P
1
,
P
2
{\displaystyle P_{1},P_{2}}
互质的情况下,可以有一些其他更具技巧性的选择。
当
P
1
,
P
2
{\displaystyle P_{1},P_{2}}
互质,对所有
n
=
0
,
1
,
…
,
N
−
1
{\displaystyle n=0,1,\ldots ,N-1}
我们可以找到唯一且不重复的一组
n
1
,
n
2
{\displaystyle n_{1},n_{2}}
使得:
n
=
(
(
n
1
P
1
+
n
2
P
2
)
)
N
=
n
1
P
1
+
n
2
P
2
+
c
1
N
0
≤
n
1
<
P
2
0
≤
n
2
<
P
1
{\displaystyle n=((n_{1}P_{1}+n_{2}P_{2}))_{N}=n_{1}P_{1}+n_{2}P_{2}+c_{1}N\qquad 0\leq n_{1}<P_{2}\qquad 0\leq n_{2}<P_{1}}
其中
(
(
a
)
)
N
=
a
mod
N
{\displaystyle ((a))_{N}=a\mod N}
,代表取馀数的意思,
c
1
{\displaystyle c_{1}}
是一个整数。
例如假设
N
=
15
,
P
1
=
3
,
P
2
=
5
{\displaystyle N=15,P_{1}=3,P_{2}=5}
,则
n
=
1
{\displaystyle n=1}
对应到的
(
n
1
,
n
2
)
{\displaystyle (n_{1},n_{2})}
就是
(
2
,
2
)
{\displaystyle (2,2)}
,
有
2
×
3
+
2
×
5
mod
15
=
16
mod
15
=
1
{\displaystyle 2\times 3+2\times 5\mod 15=16\mod 15=1}
。
并且对所有
n
1
,
n
2
{\displaystyle n_{1},n_{2}}
的组合(有
P
1
×
P
2
=
N
{\displaystyle P_{1}\times P_{2}=N}
组),都对应到一个特定不重复的
n
{\displaystyle n}
。
同理我们可以把
m
=
0
,
1
,
…
,
N
−
1
{\displaystyle m=0,1,\ldots ,N-1}
表示为
m
1
,
m
2
{\displaystyle m_{1},m_{2}}
的双参数函数:
m
=
(
(
m
1
P
1
+
m
2
P
2
)
)
N
=
m
1
P
1
+
m
2
P
2
+
c
2
N
0
≤
m
1
<
P
2
0
≤
m
2
<
P
1
{\displaystyle m=((m_{1}P_{1}+m_{2}P_{2}))_{N}=m_{1}P_{1}+m_{2}P_{2}+c_{2}N\qquad 0\leq m_{1}<P_{2}\qquad 0\leq m_{2}<P_{1}}
将上述的参数代换及
N
=
P
1
P
2
{\displaystyle N=P_{1}P_{2}}
带入原式,可以得到:
F
[
(
(
m
1
P
1
+
m
2
P
2
)
)
N
]
=
∑
n
1
=
0
P
2
−
1
∑
n
2
=
0
P
1
−
1
e
−
i
2
π
N
(
m
1
P
1
+
m
2
P
2
+
c
2
N
)
(
n
1
P
1
+
n
2
P
2
+
c
1
N
)
f
[
(
(
n
1
P
1
+
n
2
P
2
)
)
N
]
{\displaystyle F[((m_{1}P_{1}+m_{2}P_{2}))_{N}]=\sum _{n_{1}=0}^{P_{2}-1}\sum _{n_{2}=0}^{P_{1}-1}e^{-i{\frac {2\pi }{N}}(m_{1}P_{1}+m_{2}P_{2}+c_{2}N)(n_{1}P_{1}+n_{2}P_{2}+c_{1}N)}f[((n_{1}P_{1}+n_{2}P_{2}))_{N}]}
接著透过一次的展开化简及应用
e
−
i
2
π
=
1
{\displaystyle e^{-i2\pi }=1}
可得:
F
[
(
(
m
1
P
1
+
m
2
P
2
)
)
N
]
=
∑
n
1
=
0
P
2
−
1
∑
n
2
=
0
P
1
−
1
e
−
i
2
π
N
(
m
1
P
1
+
m
2
P
2
)
(
n
1
P
1
+
n
2
P
2
)
e
−
i
2
π
c
2
(
n
1
P
1
+
n
2
P
2
+
c
1
N
)
e
−
i
2
π
c
1
(
m
1
P
1
+
m
2
P
2
+
c
2
N
)
e
−
i
2
π
c
1
c
2
N
f
[
(
(
n
1
P
1
+
n
2
P
2
)
)
N
]
=
∑
n
1
=
0
P
2
−
1
∑
n
2
=
0
P
1
−
1
e
−
i
2
π
N
(
m
1
P
1
+
m
2
P
2
)
(
n
1
P
1
+
n
2
P
2
)
f
[
(
(
n
1
P
1
+
n
2
P
2
)
)
N
]
{\displaystyle {\begin{aligned}F[((m_{1}P_{1}+m_{2}P_{2}))_{N}]&=\sum _{n_{1}=0}^{P_{2}-1}\sum _{n_{2}=0}^{P_{1}-1}e^{-i{\frac {2\pi }{N}}(m_{1}P_{1}+m_{2}P_{2})(n_{1}P_{1}+n_{2}P_{2})}e^{-i2\pi c_{2}(n_{1}P_{1}+n_{2}P_{2}+c_{1}N)}e^{-i2\pi c_{1}(m_{1}P_{1}+m_{2}P_{2}+c_{2}N)}e^{-i2\pi c_{1}c_{2}N}f[((n_{1}P_{1}+n_{2}P_{2}))_{N}]\\&=\sum _{n_{1}=0}^{P_{2}-1}\sum _{n_{2}=0}^{P_{1}-1}e^{-i{\frac {2\pi }{N}}(m_{1}P_{1}+m_{2}P_{2})(n_{1}P_{1}+n_{2}P_{2})}f[((n_{1}P_{1}+n_{2}P_{2}))_{N}]\end{aligned}}}
再将
N
=
P
1
P
2
{\displaystyle N=P_{1}P_{2}}
带入并再透过一次的展开化简及应用
e
−
i
2
π
=
1
{\displaystyle e^{-i2\pi }=1}
可得:
F
[
(
(
m
1
P
1
+
m
2
P
2
)
)
N
]
=
∑
n
1
=
0
P
2
−
1
∑
n
2
=
0
P
1
−
1
f
[
(
(
n
1
P
1
+
n
2
P
2
)
)
N
]
e
−
i
2
π
P
2
m
1
P
1
n
1
e
−
i
2
π
P
1
m
2
P
2
n
2
e
−
i
2
π
m
1
n
2
e
−
i
2
π
m
2
n
1
=
∑
n
1
=
0
P
2
−
1
∑
n
2
=
0
P
1
−
1
f
[
(
(
n
1
P
1
+
n
2
P
2
)
)
N
]
e
−
i
2
π
P
2
m
1
P
1
n
1
e
−
i
2
π
P
1
m
2
P
2
n
2
{\displaystyle {\begin{aligned}F[((m_{1}P_{1}+m_{2}P_{2}))_{N}]&=\sum _{n_{1}=0}^{P_{2}-1}\sum _{n_{2}=0}^{P_{1}-1}f[((n_{1}P_{1}+n_{2}P_{2}))_{N}]e^{-i{\frac {2\pi }{P_{2}}}m_{1}P_{1}n_{1}}e^{-i{\frac {2\pi }{P_{1}}}m_{2}P_{2}n_{2}}e^{-i2\pi m_{1}n_{2}}e^{-i2\pi m_{2}n_{1}}\\&=\sum _{n_{1}=0}^{P_{2}-1}\sum _{n_{2}=0}^{P_{1}-1}f[((n_{1}P_{1}+n_{2}P_{2}))_{N}]e^{-i{\frac {2\pi }{P_{2}}}m_{1}P_{1}n_{1}}e^{-i{\frac {2\pi }{P_{1}}}m_{2}P_{2}n_{2}}\end{aligned}}}
观察上式,并加上括号辅助厘清运算顺序我们可以得到:
F
[
(
(
m
1
P
1
+
m
2
P
2
)
)
N
]
=
∑
n
2
=
0
P
1
−
1
{
∑
n
1
=
0
P
2
−
1
f
[
(
(
n
1
P
1
+
n
2
P
2
)
)
N
]
e
−
i
2
π
P
2
m
1
P
1
n
1
}
e
−
i
2
π
P
1
m
2
P
2
n
2
{\displaystyle F[((m_{1}P_{1}+m_{2}P_{2}))_{N}]=\sum _{n_{2}=0}^{P_{1}-1}\left\{\sum _{n_{1}=0}^{P_{2}-1}f[((n_{1}P_{1}+n_{2}P_{2}))_{N}]e^{-i{\frac {2\pi }{P_{2}}}m_{1}P_{1}n_{1}}\right\}e^{-i{\frac {2\pi }{P_{1}}}m_{2}P_{2}n_{2}}}
内层的运算可以视为将双参数函数
f
[
(
(
n
1
P
1
+
n
2
P
2
)
)
N
]
{\displaystyle f[((n_{1}P_{1}+n_{2}P_{2}))_{N}]}
中的一个参数
n
1
{\displaystyle n_{1}}
,
透过离散傅立叶变换 取代为由一个与
m
1
{\displaystyle m_{1}}
有关的变数
m
3
=
(
(
m
1
P
1
)
)
P
2
{\displaystyle m_{3}=((m_{1}P_{1}))_{P_{2}}}
描述,
得到一个新函数
G
1
[
m
3
,
n
2
]
{\displaystyle G_{1}[m_{3},n_{2}]}
(这步因为对每个不同
n
2
{\displaystyle n_{2}}
,都需要做一次将
n
1
{\displaystyle n_{1}}
取代为
m
3
{\displaystyle m_{3}}
的转换,
共需要
P
1
{\displaystyle P_{1}}
个
P
2
{\displaystyle P_{2}}
点离散傅立叶变换 ):
F
[
(
(
m
1
P
1
+
m
2
P
2
)
)
N
]
=
∑
n
2
=
0
P
1
−
1
G
1
[
m
3
,
n
2
]
e
−
i
2
π
P
1
m
2
P
2
n
2
{\displaystyle F[((m_{1}P_{1}+m_{2}P_{2}))_{N}]=\sum _{n_{2}=0}^{P_{1}-1}G_{1}[m_{3},n_{2}]e^{-i{\frac {2\pi }{P_{1}}}m_{2}P_{2}n_{2}}}
外层的运算可以视为将函数
G
1
[
m
3
,
n
2
]
{\displaystyle G_{1}[m_{3},n_{2}]}
中的参数
n
2
{\displaystyle n_{2}}
,
透过离散傅立叶变换 取代为由一个与
m
2
{\displaystyle m_{2}}
有关的变数
m
4
=
(
(
m
2
P
2
)
)
P
1
{\displaystyle m_{4}=((m_{2}P_{2}))_{P_{1}}}
描述,
得到一个新函数
G
2
[
m
3
,
m
4
]
{\displaystyle G_{2}[m_{3},m_{4}]}
(这步因为对每个不同
m
3
{\displaystyle m_{3}}
,都需要做一次将
n
2
{\displaystyle n_{2}}
取代为
m
4
{\displaystyle m_{4}}
的转换,
共需要
P
2
{\displaystyle P_{2}}
个
P
1
{\displaystyle P_{1}}
点离散傅立叶变换 ):
F
[
m
]
=
F
[
(
(
m
1
P
1
+
m
2
P
2
)
)
N
]
=
G
2
[
m
3
,
m
4
]
=
G
2
[
(
(
m
1
P
1
)
)
P
2
,
(
(
m
2
P
2
)
)
P
1
]
{\displaystyle F[m]=F[((m_{1}P_{1}+m_{2}P_{2}))_{N}]=G_{2}[m_{3},m_{4}]=G_{2}[((m_{1}P_{1}))_{P_{2}},((m_{2}P_{2}))_{P_{1}}]}
最后透过
F
{\displaystyle F}
与
G
2
{\displaystyle G_{2}}
在不同
m
1
,
m
2
{\displaystyle m_{1},m_{2}}
时的点点数值对应关系,
就成功仅使用
P
1
{\displaystyle P_{1}}
与
P
2
{\displaystyle P_{2}}
点离散傅立叶变换 ,描述了原先的
N
{\displaystyle N}
点离散傅立叶变换 。
而这个方法透过聪明的选取表达
m
,
n
{\displaystyle m,n}
的方式,使得拆解的过程中完全不需要多馀的乘法运算,
总运算量可以简单表示为:
P
2
B
1
+
P
1
B
2
+
3
D
3
+
2
D
2
{\displaystyle P_{2}B_{1}+P_{1}B_{2}+3D_{3}+2D_{2}}
其中
B
1
{\displaystyle B_{1}}
是
P
1
{\displaystyle P_{1}}
傅立叶所需乘法数,
B
2
{\displaystyle B_{2}}
是
P
2
{\displaystyle P_{2}}
傅立叶所需乘法数。
虽然这个方法可以较上面的方法节省运算量,
但这个方法也牵涉较为复杂的
m
,
n
{\displaystyle m,n}
与
m
1
,
m
2
,
n
1
,
n
2
{\displaystyle m_{1},m_{2},n_{1},n_{2}}
转换,较为不直觉且不易理解,
也会遇到许多需要取馀数的运算,可能会需要较大的空间建表进行查表法。
最后关于实际上要如何求得
n
{\displaystyle n}
与
n
1
,
n
2
{\displaystyle n_{1},n_{2}}
的转换关系,
可以先透过辗转相除法 取得一对特定的
n
11
,
n
21
{\displaystyle n_{11},n_{21}}
使得:
(
(
n
11
P
1
+
n
21
P
2
)
)
N
=
1
{\displaystyle ((n_{11}P_{1}+n_{21}P_{2}))_{N}=1}
然后我们可以知道对于任意整数
0
≤
k
<
N
{\displaystyle 0\leq k<N}
有:
(
(
k
n
11
P
1
+
k
n
21
P
2
)
)
N
=
k
=
(
(
k
n
11
P
1
−
c
1
N
+
k
n
21
P
2
−
c
2
N
)
)
N
=
(
(
(
k
n
11
−
c
1
P
2
)
P
1
+
(
k
n
21
−
c
2
P
1
)
P
2
)
)
N
{\displaystyle ((kn_{11}P_{1}+kn_{21}P_{2}))_{N}=k=((kn_{11}P_{1}-c_{1}N+kn_{21}P_{2}-c_{2}N))_{N}=(((kn_{11}-c_{1}P_{2})P_{1}+(kn_{21}-c_{2}P_{1})P_{2}))_{N}}
然后就可以得到:
n
1
k
=
(
(
k
n
11
)
)
P
2
n
2
k
=
(
(
k
n
21
)
)
P
1
{\displaystyle n_{1k}=((kn_{11}))_{P_{2}}\qquad n_{2k}=((kn_{21}))_{P_{1}}}
满足:
(
(
n
1
k
P
1
+
n
2
k
P
2
)
)
N
=
k
0
≤
n
1
k
<
P
2
0
≤
n
2
k
<
P
1
{\displaystyle ((n_{1k}P_{1}+n_{2k}P_{2}))_{N}=k\qquad 0\leq n_{1k}<P_{2}\qquad 0\leq n_{2k}<P_{1}}
Widhe, T., J. Melander, et al. (1997). Design of efficient radix-8 butterfly PEs for VLSI. Circuits and Systems, 1997. ISCAS '97., Proceedings of 1997 IEEE International Symposium on.
Lihong, J., G. Yonghong, et al. (1998). A new VLSI-oriented FFT algorithm and implementation. ASIC Conference 1998. Proceedings. Eleventh Annual IEEE International.
Duhamel, P. and H. Hollmann (1984). "Split-radix FFT algorithm." Electronics Letters 20(1): 14-16.
Vetterli, M. and P. Duhamel (1989). "Split-radix algorithms for length-pm DFT's." Acoustics, Speech and Signal Processing, IEEE Transactions on 37(1): 57-64.
Richards, M. A. (1988). "On hardware implementation of the split-radix FFT." Acoustics, Speech and Signal Processing, IEEE Transactions on 36(10): 1575-1581.
Shousheng, H. and M. Torkelson (1996). A new approach to pipeline FFT processor. Parallel Processing Symposium, 1996., Proceedings of IPPS '96, The 10th International.
Shousheng, H. and M. Torkelson (1998). Designing pipeline FFT processor for OFDM (de)modulation. Signals, Systems, and Electronics, 1998. ISSSE 98. 1998 URSI International Symposium on.