FFT (Chirp z-変換アルゴリズム)

長さが素数の場合はこの間書いたので，今回はどんな長さの場合にも適用できる FFT のアルゴリズムをまとめておく．これは Chirp z-変換アルゴリズムあるいは Bluestein の FFT と呼ばれるアルゴリズムで，Rader の FFT と同じく畳み込み演算を用いるが，使用方法が異なる．

準備

$x(n)\;\;n=0,\cdots,N-1$ のスペクトルを $X(k)\;\;k=0,\cdots,N-1$ とする．つまり
$\;\;\;\;X(k)=\sum_{n=0}^{N-1} x(n) W_N^{kn}\;\;{\rm where}\;W_N=\exp \left( -\frac{2 \pi j}{N} \right)$

変形

回転子 $W_N$ の指数 $kn$ を次のように変形する．
$\;\;\;\;kn=-\frac{(k-n)^2}{2}+\frac{k^2}{2}+\frac{n^2}{2}$
これを元の式に代入してやると，
$\;\;\;\;\begin{array} X(k) &=& \sum_{n=0}^{N-1} x(n) W_N^{-\frac{(k-n)^2}{2}+\frac{k^2}{2}+\frac{n^2}{2}} \\ &=& W_N^{\frac{k^2}{2}} \sum_{n=0}^{N-1} \left\{ x(n) W_N^{\frac{n^2}{2}} \right\} W_N^{-\frac{(k-n)^2}{2}} \\ &=& W_{2N}^{k^2} \sum_{n=0}^{N-1} \left\{ x(n) W_{2N}^{n^2} \right\} W_{2N}^{-(k-n)^2} \\ \end{array}$
となる．ここで， $\{\alpha_n\},\{\beta_n\}$ を次のように定める．
$\;\;\;\;\begin{array} \alpha_n &=& x(n)W_{2N}^{n^2} \\ \beta_n &=& W_{2N}^{-n^2} \end{array}$
この $\{\alpha_n\},\{\beta_n\}$ を使うと，スペクトル $X(k)$ の計算式が畳み込みの形に置き換えられることがわかる．
$\;\;\;\;X(k)=\beta_k^\ast \sum_{n=0}^{N-1} \alpha_n \beta_{k-n}=\beta_k^{\ast} (\alpha_k \ast \beta_k)$
あとはこの畳み込み演算を FFT で行えばよい．

FFT による畳み込み演算

使用する FFT の長さ

FFT で畳み込み演算を行う場合，一般には巡回畳み込みになってしまう．しかし，ここでは線形畳み込みを行いたいので，FFT の長さ $L$ に条件が付く．畳み込み演算を行う 2 つの数列の長さをそれぞれ $M,N$ とすれば， $L \geq M+N-1$ であれば線形畳み込みを行うことができる．

今回の場合， $M=N$ なので， $L$ の条件は $L \geq 2N-1$ である．この条件を満たしていれば $L$ はどんな値でも構わないが，Bluestein の FFT では $L$ に 2 のべき乗の値を選ぶことが多い．これは長さが 2 のべき乗の FFT は，Cooley-Tukey の手法などを用いれば簡単に構成でき， $O(N \log N)$ であることが保証されているためである．

パディング

$\{\alpha_n\},\{\beta_n\}$ はともに長さ $N$ であるため，長さ $L$ の FFT を使用するには，パディング(いわゆるゼロ・パディング)を行い，長さを $L$ にしなければならない．

$\alpha_n$ についてはあまり問題はない．新しい数列を $\{\alpha'_l\}$ とすれば，単純に
$\;\;\;\;\alpha'_l =\left\{ \begin{array} \alpha_l \; & (0 \leq l < N) \\ 0 & (N \leq l < L) \end{array} \right.$
とすればよい．

ところが， $\{\beta_n\}$ の場合は畳み込みの式の中で $k-n$ が負になることがある．そのため，畳み込み演算を考慮して，新しい数列 $\{\beta'_l\}$ を次のように定める．
$\;\;\;\; \left\{ \begin{array}{rcl} \beta'_0 &=& \beta_0 & (l=0) \\ \beta'_l = \beta'_{L-l} &=& \beta_l & (0 < l < N) \\ \beta'_l &=& 0 & (N \leq l \leq L-N) \end{array} \right.$

畳み込み

FFTを用い，これら 2 つの数列 $\{\alpha'_l\},\{\beta'_l\}$ から，それぞれスペクトル $A'_l,B'_l$ を得る．それぞれ掛け合わせた $\Gamma'_l=A'_l B'_l$ を要素とする数列から IFFT で $\{\gamma'_l\}$ を得る．これで $\alpha_k$ と $\beta_k$ の畳み込み演算を行ったことになる．
$\;\;\;\;\alpha_k \ast \beta_k = \gamma'_k$