FFT

定義

• そのまま計算すると、$N$ 個の出力それぞれに $N$ 回の掛け算。合わせて $N^2$ 回[1]
• このページは速く計算する仕組みだけを扱う

半分に分ける

時間間引き（偶数番と奇数番に分ける）[1][5]:

$E[k]$

偶数番 $x[2n]$ の $N/2$ 点の DFT

$O[k]$

奇数番 $x[2n+1]$ の $N/2$ 点の DFT

$W_N^k$

回転因子（大きさ 1 の複素数）

周波数間引き（前半と後半に分ける）[2][6]:

• $N$ 点の DFT 1 つが $N/2$ 点の DFT 2 つになる。2 の累乗なら 1 点まで ${\log }_{2}N$ 回分けられる
• 1 点の DFT は値そのもの

バタフライ

• 1 つのバタフライは複素数の掛け算 1 回と足し算 2 回[5]
• 1 段に $N/2$ 個、${\log }_{2}N$ 段
• 出力を入力と同じ場所に書ける（その場の計算）[6]

ビット反転

• 偶数と奇数に分けることを繰り返すと、入力は番号の 2 進数を逆に読んだ順に並ぶ[1][5]
• 時間間引きは並べ替えが先、周波数間引きは後[6]
• 2 回反転すると戻る。入れ替えは対ごとに 1 回

途中の値

• 時間間引きの段 $s$ の後: 並んだ $2^s$ 行ずつが、$N/2^s$ おきに取った入力の $2^s$ 点の DFT
• 周波数間引きの段 $s$ の後: 並んだ $N/2^s$ 行ずつの DFT が、$2^s$ おきの出力 $X[k_0+2^{s}j]$
• どちらも組の先頭の番号は、組の番号をビット反転したもの

演算の数

• 実数の演算の数は、複素数の掛け算を実数の掛け算 4 回と足し算 2 回とし、$W_N^0=1$ と $W_N^{N/4}=-i$ は掛けず、$W_N^{N/8}$ の向きは 2 回で済ませた数。$2N{\log }_{2}N-7N+12$ と $3N{\log }_{2}N-3N+4$（表 8.5.1）[6]
• 複素数の掛け算を実数の掛け算 3 回と足し算 3 回で数えると $\frac{3}{2}N{\log }_{2}N-5N+8$ と $\frac{7}{2}N{\log }_{2}N-5N+8$（式 (25)、表 1・2）[4]
• すべて掛けると $5N{\log }_{2}N$。FFTW は速さの目安にこの数を使う[8]
• 100 MHz の Pentium で 1,024 点の DFT は約 25 s、FFT は約 70 ms[5]

ほかの分け方

• $N=r_1{r}_2\cdots$ と分ければ 2 の累乗でなくてよい。基数 $r$ の手間は $rN{\log }_{r}N$ に比例し、3 が最も少ないが 2 や 4 との差は約 6 %[1]

逆変換

共役（${}^{\ast }$）をとって同じ FFT に通し、共役をとって $N$ で割ります。

歴史

• 1805 年ごろ、Gauss が小惑星の軌道の計算で同じ分け方を使った。死後 1866 年に出版[3]
• 1965 年、Cooley と Tukey が計算機のアルゴリズムとして発表[1]
• 1966 年、Gentleman と Sande が回転因子を後で掛ける形を示した[2]

関連

• FFT で畳み込む: 畳み込み。2 次元: 画像の畳み込み
• 窓関数と音のスペクトル、2 本の実数を 1 回の FFT で: スペクトラムアナライザ
• FFT で求める連続ウェーブレット変換: ウェーブレット変換

このツールで決めたこと

• FFT はスペクトラムアナライザと同じ実装（並べ替えの表と回転因子の表）。図はそれを段ごとに止めて記録したもの
• 節点の矢印は、列ごとに最大の大きさで割って描く
• 拡大の図は今のバタフライ 1 つを複素平面に描く
• 時間間引きは b を回してから足し引き、周波数間引きは足し引きしてから差を回す
• 縮尺は 5 つの点（a・b・回した値・和・差）が入る大きさ
• 「複素数の掛け算」は回転因子が 1 でも数える。実数の演算の数は自明な回転因子を除く

出典

1. 直接の計算は N² 回、このアルゴリズムは 2N log₂ N 回未満（297 頁）、基数 r の手間（式 (12)・299 頁の表）、ビットの反転（300 頁）J. W. Cooley and J. W. Tukey, “An algorithm for the machine calculation of complex Fourier series,” Mathematics of Computation, vol. 19, no. 90, pp. 297–301, 1965, doi:10.1090/S0025-5718-1965-0178586-1
2. 回転因子（565 頁）、Cooley の版と Sande の版（566〜567 頁）W. M. Gentleman and G. Sande, “Fast Fourier transforms: for fun and profit,” Proc. AFIPS Fall Joint Computer Conference, pp. 563–578, 1966, doi:10.1145/1464291.1464352
3. Gauss の 1805 年の計算と 1866 年の出版（再録版 266・269・271 頁）M. T. Heideman, D. H. Johnson, and C. S. Burrus, “Gauss and the history of the fast Fourier transform,” IEEE ASSP Magazine, vol. 1, no. 4, pp. 14–21, 1984, doi:10.1109/MASSP.1984.1162257。再録 Archive for History of Exact Sciences, vol. 34, pp. 265–277, 1985
4. 演算の数え方（270 頁）、基数 2 と分割基数の式（式 (25)・(31)）、表 1・2P. Duhamel and M. Vetterli, “Fast Fourier transforms: a tutorial review and a state of the art,” Signal Processing, vol. 19, no. 4, pp. 259–299, 1990, doi:10.1016/0165-1684(90)90158-U
5. 分解（図 12-2）、ビット反転（図 12-3）、バタフライ（図 12-6）、流れ図（図 12-7）、速さの比べ（図 12-8）S. W. Smith, The Scientist and Engineer’s Guide to Digital Signal Processing, California Technical Publishing, 1997, ch. 12 The Fast Fourier Transform
6. その場の計算と並べ替え、時間間引きと周波数間引き（2.3 節）、演算の数（8.5 節、表 8.5.1）C. S. Burrus, Fast Fourier Transforms, OpenStax CNX
7. 分割基数の 4N log₂ N − 6N + 8、改良した分割基数（式 (1)）S. G. Johnson and M. Frigo, “A modified split-radix FFT with fewer arithmetic operations,” IEEE Trans. Signal Processing, vol. 55, no. 1, pp. 111–119, 2007, doi:10.1109/TSP.2006.882087
8. 速さの目安 5N log₂ NM. Frigo and S. G. Johnson, FFT Benchmark Methodology, FFTW