離散フーリエ変換
離散フーリエ変換（りさんフーリエへんかん、英語: discrete Fourier transform、DFT）とは離散化されたフーリエ変換であり、信号処理などで離散化されたデジタル信号の周波数解析などによく使われる。また偏微分方程式や畳み込み積分を効率的に計算するためにも使われる。離散フーリエ変換は（計算機上で）高速フーリエ変換(FFT)を使って高速に計算することができる。
離散フーリエ変換とは、複素関数 [math]f(x)[/math]を複素関数[math]F(t)[/math]に写す写像であって、次の式で定義されるものを言う。
[math]F(t)= \sum_{x=0}^{N-1} f(x) e^{-i\frac{2 \pi t x}{N}} \quad \quad [/math]
ここで、Nは任意の自然数、 [math]e[/math] はネイピア数、[math]i[/math] は虚数単位 ([math]i^2 = -1[/math])[1]で、[math]\pi[/math]は円周率である。このとき、{[math]x = 0,\dots,N-1[/math]}を標本点という。また、この変換を [math]\mathcal{F}_N[/math] という記号で表し、
[math]F = \mathcal{F}_N(f)[/math]
のように略記することが多い。
この逆変換にあたる逆離散フーリエ変換（英語: inverse discrete Fourier transform、IDFT）は
[math]f(x) = \frac{1}{N} \sum_{t=0}^{N-1} F(t) e^{i\frac{2\pi x t}{N}} \quad \quad[/math]
正規化係数（DFT は 1, IDFT は 1/N）や指数の符号は単なる慣習的なものであり、上式とは異なる式を扱うことがある。DFT と IDFT の差について、それぞれの正規化係数を掛けると 1 / N になることと、指数の符号が異符号であるということがだけが重要であり、根本的には同一の変換作用素と考えられる。DFT と IDFT の正規化係数を両方とも [math]\frac{1}{\sqrt{N}}[/math] にすると、両方ともユニタリ作用素（ユニタリ変換）になる。理論的にはユニタリ作用素にするのが好ましいが、実用上数値計算を行うときは上式のように正規化係数を1つにまとめて、スケーリングを一度に行うことが多い。
Contents
性質
離散フーリエ変換はフーリエ変換を離散化したものであるので、フーリエ変換と同様の性質を持つ。
完全性
離散フーリエ変換においては、有限個の標本点しか使わないため、ある関数を離散フーリエ変換し、それを逆変換した場合に、標本点以外で元の関数と一致するとは限らない。
すなわち、複素関数ｆに対して、
[math]F(t)= \sum_{x=0}^{N-1} f(x) e^{-i\frac{2 \pi t x}{N}} \quad \quad [/math]
により離散フーリエ変換を行い、それを逆変換したものをgとすると
[math]f(x)=g(x) \quad \quad x = 0,\dots,N-1[/math]
は言えるが、その他の点で{{[math]f(x)=g(x) \quad \quad [/math]}}が言えるとは限らない。これを高周波の問題、あるいはエイリアシング(aliasing)という。
選点直交性
[math]h(t,x)=e^{i2\pi\frac{t x}{N}}[/math] は内積
[math](f,g)=\sum_{n=0}^{N-1}f(n)g(n)[/math]
に関し、[math]t,t'[/math]が整数のとき直交基底である。
[math](h(t,x),h(t',x))=\sum_{n=0}^{N-1} \left(e^{ i2\pi\frac{tn}{N}}\right) \left(e^{-i2\pi\frac{t'n}{N}}\right) =N~\delta_{tt'}[/math]
[math]\delta_{tt'}[/math]はクロネッカーのデルタ
畳み込み定理と相互相関定理
[math]f(x)[/math] と [math]g(x)[/math] の畳み込み[math]f*g[/math]は次のように定義される。
[math]f*g(x) = \sum_{n=0}^{N-1} f(n) g(x-n) \quad \quad [/math]
ここで[math]f,g[/math]は次のような周期性を持つとする。
[math]f(x+N)=f(x)\quad\quad[/math]
[math]g(x+N)=g(x)\quad\quad[/math]
周期関数の畳み込みを離散フーリエ変換したものは、それぞれの離散フーリエ変換の積になる(畳み込み定理)。 つまり
[math]\mathcal{F}_N(f*g)=\mathcal{F}_N(f)\mathcal{F}_N(g)[/math]
畳み込み和を直接定義式を用いて計算すると O(N²) の計算量が掛かる。しかし、上式より一旦 DFT をしてから掛算をして、そして IDFT で戻す方法もある。DFT の高速アルゴリズムである FFT を介してこのように計算すると O(N log N) の計算量で済む。
[math]f(x)[/math]と[math]g(x)[/math]の相互相関[math]h(x)[/math]は以下のように定義される。
[math]h(x)=(f\star g)(x)= \sum_{n=0}^{N-1}f(n)g(x+n)[/math]
[math]f,g[/math]が上記の周期性を持てば、
[math](\mathcal{F}_Nh)(u)=(\mathcal{F}_N(f\star g))(u)=\mathcal{F}_N(f)(-u)\mathcal{F}_N(g)(u)[/math]
補間三角多項式との関係
実数値関数の離散フーリエ変換
応用上は、実数値関数を対象とすることが多いが、このとき、
[math](\mathcal{F}_Nf)(-t)=\overline{(\mathcal{F}_Nf)(t)}[/math]
([math]\bar{z}[/math]は[math]z[/math]の複素共役)。そのため出力[math](\mathcal{F}_Nf)(t)[/math]の半分([math]t\gt 0[/math])で全ての情報を持っていることになる。
一般化
2次元での変換
デジタル画像処理では2次元変換が画像の周波数成分を解析するのに使われる。
変換は以下のように定義される。
[math]F(u, v) = \sum_{x=0}^{M-1} \sum_{y=0}^{N-1} f(x,y) e^{- 2 \pi i \left(\frac{u x}{M} +\frac{v y}{N} \right) } \quad \quad[/math]
そして逆変換は次のようになる。
[math]f(x, y) = \frac{1}{M N} \sum_{u=0}^{M-1} \sum_{v=0}^{N-1} F(u, v) e^{2 \pi i \left(\frac{u x}{M} +\frac{v y}{N} \right) } \quad \quad[/math]
但し
- [math]f(x, y)[/math]は2次元信号(例えば画像)であり、fのx列y行成分のことである。
- [math]F(u, v)[/math]は[math]f(x, y)[/math]の2次元周波数スペクトラムである。ここでuはx成分の周波数、vはy成分の周波数である。
2次元DFT は行列を用いて簡単に記述できる。
[math]F = W f^T W[/math]
ここで
- [math]F = \begin{bmatrix} F(0, 0) & F(1, 0) & \ldots & F(M-1, 0) \\ F(0, 1) & F(1, 1) & \ldots & F(M-1, 1) \\ \vdots & \vdots & \ddots & \vdots \\ F(0, N-1) & F(1, N-1) & \ldots & F(M-1, N-1) \end{bmatrix}[/math]
- [math]W = \begin{bmatrix} \omega_n^{0 \cdot 0} & \omega_n^{0 \cdot 1} & \ldots & \omega_n^{0 \cdot (N-1)} \\ \omega_n^{1 \cdot 0} & \omega_n^{1 \cdot 1} & \ldots & \omega_n^{1 \cdot (N-1)} \\ \vdots & \vdots & \ddots & \vdots \\ \omega_n^{(N-1) \cdot 0} & \omega_n^{(N-1) \cdot 1} & \ldots & \omega_n^{(N-1) \cdot (N-1)} \\ \end{bmatrix}[/math] (注:これはユニタリ行列)
- [math]f = \begin{bmatrix} f(0, 0) & f(1, 0) & \ldots & f(M-1, 0) \\ f(0, 1) & f(1, 1) & \ldots & f(M-1, 1) \\ \vdots & \vdots & \ddots & \vdots \\ f(0, N-1) & f(1, N-1) & \ldots & f(M-1, N-1) \end{bmatrix}[/math]
行列の表記による変形
2次元DFT を行列計算によって以下のように変形できる。
- [math]F(u, v) = \sum_{x=0}^{M-1} \sum_{y=0}^{N-1} f(x,y) e^{ -2 \pi i ( \frac{u x}{M} + \frac{v y}{N} ) }[/math]
- [math]F(u, v) = \sum_{x=0}^{M-1} \left[ \sum_{y=0}^{N-1} f(x,y) e^{ -2 \pi i \frac{v y}{N} } \right] e^{ -2 \pi i \frac{u x}{M} }[/math]
- [math]F(u, v) = \sum_{x=0}^{M-1} \underbrace{ \left[ \sum_{y=0}^{N-1} f(x,y) e^{ -2 \pi i \frac{v y}{N} } \right] }_{ \mathcal{F}_y f(x, y) } e^{ -2 \pi i \frac{u x}{M} } [/math]
- [math]F_v = W f[/math]
- [math]F(u, v) = \sum_{x=0}^{M-1} F_v(x, v) e^{- 2 \pi i \frac{u x}{M} }[/math]
- [math]F = W F_v^{T}[/math]
- [math]F = W f^T W[/math]
以下上式の 1 - 7 を解説すると、
- 2次元DFTの定義
- 式1から e-2πiux/M を内側σから出した
- 内側のσは、信号[math]f(x, y)[/math]のyの次元([math]\mathcal{F}_y\{f(x, y)\}[/math]と書いた行)の1次元DFTであることを示している
- 式3で示した、[math]\mathcal{F}_y\{f(x, y)\}[/math]の行列表現
- 式3の注目箇所を[math]F_v(x, v)[/math]で書き変えた - [math]F_v(x, v)[/math]は[math]f(x, y)[/math]のx行目の行のv番目の周波数。[math]f(x, y)[/math]の1次元DFTの行を集めたものがFvである。
- 式4の1次元DFTの行列表現と同様に、FvTを使って式5を表現した
- 式4を式6に代入した結果。ここでWは対称行列であるのでW=WTとした。
[math]F_v[/math]の行は[math]f(x, y)[/math]のx行目の行を1次元DFTしたものである。ゆえに[math]F_v[/math]は[math]f(x, y)[/math]の各行の1次元DFTした結果の行ベクトルを集めたものになる。F=WFvTにおける、FvTを後から掛ける作用素は[math]F_v[/math]の列の1次元DFTしたものと同じと考えて良い。
つまり、2次元DFT（2次元フーリエ変換も同様だが）は[math]f(x, y)[/math]を、各行ごとに1次元DFTし、その結果をさらに各列ごとに1次元DFTする事と等価である。ここで、1次元DFTの計算はFFTのアルゴリズムで高速に計算できる。そのため実用上は2次元DFTも、2次元FFTとして計算される。
離散フーリエ変換表
表中で、
[math]W=e^{-i\frac{2\pi}{N}}[/math]
| 時間領域 | 周波数領域 | 備考 |
|---|---|---|
| [math]f(x)=\frac{1}{N}\sum_{k=0}^{N-1}F(k)W^{-kx}[/math] | [math]F(t)=\sum_{n=0}^{N-1}f(n) W^{tn}[/math] | IDFT,DFTのWを使った定義 |
| [math]f(x)[/math] | [math]F(t)[/math] | 定義 |
| [math]g(x)[/math] | [math]G(t)[/math] | 定義 |
| [math]af(x)+bg(x)[/math] | [math]aF(t)+bG(t)[/math] | 線形性 |
| [math]f(x){e^{i\frac{2\pi{mx}}{N}}}[/math] | [math]F(t-m)[/math] | 時間軸変調、周波数軸移動 |
| [math]f(x-a)[/math] | [math]W^{at}F(t)[/math] | 時間軸移動（正） |
| [math]f(x+a)[/math] | [math]W^{-at}F(t)[/math] | 時間軸移動（負） |
| [math]\sum_{m=0}^{N-1}f(m)g(x-m)[/math] | [math]F(t)G(t)[/math] | 時間軸畳み込み、周波数軸積 |
| [math]f(x)g(x)[/math] | [math]{\frac{1}{N}}\sum_{m=0}^{N-1}F(m)G(t-m)[/math] | 時間軸積、周波数軸畳み込み |
| [math]\overline{f(x)}[/math] | [math]\overline{F(N-t)}[/math] | 時間軸共役、周波数軸反転 |
| [math]\overline{f(N-x)}[/math] | [math]\overline{F(t)}[/math] | 時間軸反転、周波数軸共役 |
| [math]Re(f(x))[/math] | [math]\frac{1}{2}(F(t)+\overline{F(N-t)})[/math] | 時間軸実部、周波数軸偶関数 |
| [math]Im(f(x))[/math] | [math]\frac{1}{2i}(F(t)-\overline{F(N-t)})[/math] | 時間軸虚部、周波数軸奇関数 |
| [math]\frac{1}{2}(f(x)+\overline{f(N-x)})[/math] | [math]Re(F(t))[/math] | 時間軸偶関数、周波数軸実部 |
| [math]\frac{1}{2i}(f(x)-\overline{f(N-x)})[/math] | [math]Im(F(t))[/math] | 時間軸奇関数、周波数軸虚部 |
| [math]a^x\,[/math] | [math]\frac{1-a^N{W}^{tx}}{1-{a}W^t}[/math] | べき乗 |
応用
DFTは多くの広い分野で利用されている。ここでは、その中の一部を示しているだけに過ぎない。これらの応用は高速フーリエ変換(FFT)とその逆変換(IFFT)で高速に計算できることを前提としていて、定義通りにDFTを計算しているのではない。
信号解析
信号x(t)を解析するのに使われる。ここでtは時間で[0,T]の範囲をとるものとする。例えば、音声信号の場合は、x(t)は時刻tでの空気の圧力を表現することになる。
この信号はN個の等間隔の点で標本化されて、x0, x1, x2, ... , xN-1 の実数列になる。但し標本化の間隔を Δx(=T/N) とすると xk=x(k Δx) である。これのDFTである f0, ..., fN−1 をFFTで計算できる。ただし標本化定理からこれの半分(Nが偶数とすると、fN/2 + 1, ..., fN−1)は冗長であるので捨てるか無視する。
DFT から得られる |fk|/N は信号の周波数 j/T 成分の振幅の半分の値であり、振幅スペクトルと呼ばれる。また、この偏角である arg(fj) はこの周波数成分の位相を表す。また、|fk|2はパワースペクトルと呼ばれ、この周波数成分のエネルギーを表している。
fkは信号x(t)のフーリエ級数の係数に相当するものと考えることができる。そのために、無限に広がるフーリエ級数計算を、有限のサンプル点に対しての DFT を使って近似するという形になる。連続信号の場合はこれをスペクトル推定(spectral estimation)と呼び、色々な推定法がある。(例えば、DFT が有限サンプル点を扱うことに起因するスペクトル漏れの弊害を少しでも和らげるために窓関数(窓))を使うことがよく行われる。)
データ圧縮
信号処理の分野では周波数領域（フーリエ変換）で処理をすることも少なくない。例えば、画像の非可逆圧縮や音声圧縮技術などでは離散フーリエ変換の考えが用いられている。信号に対して DFT (実装上ではFFT) を行い、人間が知覚しづらい周波数成分の情報を取り除くことで、正味の情報量を削減（圧縮）する。最も単純な方法では、IDFTの際にその取り除いた周波数情報（フーリエ係数）を0にすることにより、通常のIDFTを実現する。実装の単純化（計算の効率化）のために、実数演算のみで可能な離散コサイン変換を使って2次元情報を圧縮した例がJPEGである。
偏微分方程式
大きな数の掛け算の計算
大きな数や多項式の乗算アルゴリズムにおいて、DFTを使う高速なアルゴリズムとして1971年にショーンハーゲ・ストラッセン法が発見された。数字や多項式の係数の列は畳み込み演算の対象のベクトルと見なす。するとそれぞれのベクトルを DFT して、その結果同士を掛けて、そして逆変換することで掛算の結果が得られる（つまり畳み込み定理を使う）。2007年にショーンハーゲ・ストラッセン法より理論的に高速なアルゴリズム（en:Fürer's algorithm）が見付かった。
注
- ↑ ディジタル信号処理分野の古典である、OppenheimとSchaferの著書『ディジタル信号処理』では[math]i[/math]の代わりに[math]j^2 = -1[/math]を使用している。
参考文献
- E. O. Brigham, The Fast Fourier Transform and Its Applications (Prentice-Hall, Englewood Cliffs, NJ, 1988).
- A. V. Oppenheim, R. W. Schafer, and J. R. Buck, Discrete-Time Signal Processing (Prentice-Hall, 1999).
- S. W. Smith, The Scientist and Engineer's Guide to Digital Signal Processing, (California Technical Publishing, San Diego, 2nd edition, 1999)
- A. V. Oppenheim, R. W. Schafer, 伊達玄訳, ディジタル信号処理 (コロナ社, 1978).