フーリエ変換とDFTのつながり

POINT

  • フーリエ変換〜DTFT〜DFTのつながりを整理する.
  • サンプリングデータから,元の連続信号のスペクトルを得るにはどうすればよいか,という視点で考える.

  • DFTと元の連続信号のフーリエ変換はどういう関係にあるのか?
  • サンプリング信号は,どうしてインパルス列で定義するのか?
など,個人的に悩んだ点を整理しました.

【関連記事】

アナログ信号のフーリエ変換(スペクトル)を求めたい!

我々が解析したい,実世界のデータはアナログ(連続)信号です.しかし,全ての時間で信号のデータを取得するのは不可能で,一定時間ごとにデータをサンプリングすることになります.

こうして得た離散的なデータから,「元のアナログ信号のフーリエ変換(スペクトル)」を求めるにはどうしたら良いでしょうか?

サンプリングデータが無限個の場合

元のアナログ信号$g(t)$をサンプリング周期$T_{\mathrm{s}}$で「無限個」サンプリングしたデータ列$\{g(k T_{\mathrm{s}})\}_{k=-\infty}^{\infty}$を考えます.

ここで,サンプリングデータを基につくった「連続」信号

\begin{aligned}
g_{\mathrm{s}} (t)
&= \sum_{k=-\infty}^{\infty} g(k T_{\mathrm{s}}) \delta(t - kT_{\mathrm{s}}) \\
&=g(t) \sum_{k=-\infty}^{\infty} \delta(t - kT_{\mathrm{s}})
\end{aligned}
を考えると,そのフーリエ変換は
\begin{aligned}
& G_{\mathrm{s}} (f) = \mathcal{F}[g_{\mathrm{s}}] (f) \\
&= \int_{-\infty}^{\infty} \underbrace{g_{\mathrm{s}}(t)}_{\mathrlap{\displaystyle = \sum_{k=-\infty}^{\infty} g(k T_{\mathrm{s}}) \delta(t - kT_{\mathrm{s}})}} e^{-i 2\pi f t} \,\mathrm{d} t \\
&= \sum_{k=-\infty}^{\infty} g(k T_{\mathrm{s}}) e^{-i 2\pi f k T_{\mathrm{s}}}
\end{aligned}
あるいは,
\begin{aligned}
&G_{\mathrm{s}} (f) = \mathcal{F}[g_{\mathrm{s}}] (f) \\
&=G(f) \otimes \underbrace{\mathcal{F}\Biggl[\sum_{k=-\infty}^{\infty} \delta(t - kT_{\mathrm{s}}) \Biggr] (f)}_{\!\!\! \mathrlap{\displaystyle = \frac{1}{T_{\mathrm{s}}} \sum_{n=-\infty}^{\infty} \delta (f - n/T_{\mathrm{s}}) }} \\
&=\int_{-\infty}^{\infty} g(\phi) \cdot \frac{1}{T_{\mathrm{s}}} \sum_{n=-\infty}^{\infty} \delta (f - n/T_{\mathrm{s}} - \phi) \, \mathrm{d}\phi\\
&= \frac{1}{T_{\mathrm{s}}} \sum_{n=-\infty}^{\infty} G (f - n f_{\mathrm{s}} )
\end{aligned}
となります(離散時間フーリエ変換(DTFT)と呼ぶ).ここで,$f_{\mathrm{s}} = 1/T_{\mathrm{s}}$はサンプリング周波数であり,$G = \mathcal{F}[g]$は「元のアナログ信号$g$のフーリエ変換」です.

よって,サンプリングデータから$g_{\mathrm{s}}$をつくってフーリエ変換すれば,もとのアナログ信号$g$のフーリエ変換$G$を周期$f_{s}$で平行移動させながら足し合わせた$G_{\mathrm{s}}$を得られます.これを上手く使って,$G(f)$を求められないでしょうか?

もし,$G(f)$が限られた範囲でだけゼロでない値をとる場合には,$f_{s}$を十分大きく取ることで上式の各項の和が重ならないようにでき,$G_{\mathrm{s}}$から$G$を求められます.具体的には,

\begin{aligned}
\begin{cases}
\, G(f) \neq 0 & (|f| < w / 2) \\
\, G(f) = 0 & (|f| \geq w / 2)
\end{cases}
\end{aligned}
のときに,$f_{s} \geq w$となるようにサンプリング周波数を選べば,和は重なりません.$w$はナイキスト周波数と呼ばれます.

サンプリングデータが有限個の場合

現実的には,無限回サンプリングすることはできず,有限個($N$個)のサンプリングデータ$\{g(k T_{\mathrm{s}})\}_{k=0}^{N - 1}$を扱う必要があります.そこで,サンプリングデータが有限個の場合の議論を,無限個の場合の議論に帰着させることを考えます.

アナログ信号$g$に矩形窓

\begin{aligned}
w (t) &=
\begin{cases}
\, 1 & (0 \leq t \leq (N - 1) T_{\mathrm{s}} ) \\
\, 0 & (\mathrm{otherwise})
\end{cases} \\
&=\mathrm{rect} \bigl(t - (N - 1) T_{\mathrm{s}} / 2 \bigr)
\end{aligned}
を掛た信号$h = w \cdot g$をサンプリングすれば,上の$g_{\mathrm{s}}$に対応する信号は
\begin{aligned}
h_{\mathrm{s}} (t)
&= \sum_{k=-\infty}^{\infty} w(k T_{\mathrm{s}}) g(k T_{\mathrm{s}}) \delta(t - kT_{\mathrm{s}}) = w(t) g_{\mathrm{s}} (t) \\
&= \sum_{k=0}^{N - 1} g(k T_{\mathrm{s}}) \delta(t - kT_{\mathrm{s}})
\end{aligned}
となります.これは,$g_{\mathrm{s}}$を有限個($N$個)のサンプリングデータ$\{g(k T_{\mathrm{s}})\}_{k=0}^{N - 1}$で打ち切った形になっています.つまり,我々が扱える表式になっています

そこで,$h_{\mathrm{s}}$のフーリエ変換$H_{\mathrm{s}}$がどのような式で表されるのか,$g$の連続フーリエ変換とはどのような関係にあるのかを調べます.

DFT

まず,$H_{\mathrm{s}}$がどのような式で表されるのかを調べます.

$h_{\mathrm{s}}$の連続フーリエ変換は

\begin{aligned}
& H_{\mathrm{s}} (f) = \mathcal{F}[h_{\mathrm{s}}] (f) \\
&= \int_{-\infty}^{\infty} h_{\mathrm{s}}(t) e^{-i 2\pi f t} \,\mathrm{d} t \\
&= \sum_{k=0}^{N - 1} g(k T_{\mathrm{s}}) e^{-i 2\pi f k T_{\mathrm{s}}}
\end{aligned}
となり,これも$G_{\mathrm{s}} (f)$を有限個($N$個)のサンプリングデータ$\{g(k T_{\mathrm{s}})\}_{k=0}^{N - 1}$で打ち切った形になります.

$H_{\mathrm{s}} (f)$は$f$に関する連続信号なので,やはり現実では扱うことができません.計算できるのは有限個の$f$の値だけです.そこで,$t$と同様に$f$も離散化することを考えます.上の表式から,$H_{\mathrm{s}} (f)$は$f_{\mathrm{s}}$を基本周期とする周期関数なので,$[0, f_{\mathrm{s}}]$を$N$等分して

\begin{aligned}
H_{\mathrm{s}} (n f_{\mathrm{s}} / N) = \sum_{k=0}^{N - 1} g(k T_{\mathrm{s}}) e^{-i 2\pi n k / N} \\
(n = 0,...,N-1)
\end{aligned}
を考えます.これを離散フーリエ変換(DFT)と呼びます.

連続フーリエ変換との関係

上では,連続フーリエ変換$G (f) = \mathcal{F}[g] (f) $とDTFT $G_{\mathrm{s}} (f) = \mathcal{F}[g_{\mathrm{s}}] (f) $との間の関係を調べました.

よって,$H_{\mathrm{s}} (f) = \mathcal{F}[h_{\mathrm{s}}] (f)$と$G_{\mathrm{s}} (f) = \mathcal{F}[g_{\mathrm{s}}] (f) $との関係がわかれば,$H_{\mathrm{s}} (f)$と$G (f) $との関係もわかります.

$h_{\mathrm{s}} (t)$は

\begin{aligned}
h_{\mathrm{s}} (t)
= w(t) g_{\mathrm{s}} (t)
\end{aligned}
なので,その連続フーリエ変換
\begin{aligned}
& H_{\mathrm{s}} (f) = \mathcal{F}[h_{\mathrm{s}}] (f) \\
&= \mathcal{F}[w] (f) \otimes \mathcal{F}[g_{\mathrm{s}}] (f)
\end{aligned}
となります.そして,これをサンプリングしたものがDFTでした.つまり,DFTで得られるのは$G_{\mathrm{s}} (f) = \mathcal{F}[g_{\mathrm{s}}] (f) $そのものではなく,矩形窓関数のフーリエ変換
\begin{aligned}
\mathcal{F}[w] (f)
&=\mathcal{F}[\mathrm{rect} \bigl(t - (N - 1) T_{\mathrm{s}} / 2 \bigr)] (f) \\
&=\mathcal{F}[\mathrm{rect} (t)] (f)
\cdot \exp[-i2\pi f (N - 1) T_{\mathrm{s}} / 2] \\
&= \mathrm{sinc}(f) \exp[-i\pi f (N - 1) T_{\mathrm{s}}]
\end{aligned}
を畳み込んだものになります.

窓関数

ここまでは,矩形窓
\begin{aligned}
w (t) &=
\begin{cases}
\, 1 & (0 \leq t \leq (N - 1) T_{\mathrm{s}} ) \\
\, 0 & (\mathrm{otherwise})
\end{cases} \\
&=\mathrm{rect} \bigl(t - (N - 1) T_{\mathrm{s}} / 2 \bigr)
\end{aligned}
を扱ってきました.

しかし,よく考えてみれば,有限個のデータを扱うために必要な条件は,DTFTが有限和で切れることです.つまり,

\begin{aligned}
\begin{cases}
\, w (t) \neq 1 & (0 \leq t \leq (N - 1) T_{\mathrm{s}} ) \\
\, w (t) = 0 & (\mathrm{otherwise})
\end{cases}
\end{aligned}
さえ満たしていればよいわけです.したがって,この条件を満たす関数(窓関数)の中で「DFTで連続スペクトルを知る」という目的を果たすためにより適したものを選ぶことになります.
窓関数 - Wikipedia


DFTで得られるのは,サンプリングデータに窓関数を掛けた

\begin{aligned}
h_{\mathrm{s}} (t)
&= w(t) g_{\mathrm{s}} (t)
\end{aligned}
連続フーリエ変換
\begin{aligned}
& H_{\mathrm{s}} (f) = \mathcal{F}[h_{\mathrm{s}}] (f) \\
&= \mathcal{F}[w] (f) \otimes \mathcal{F}[g_{\mathrm{s}}] (f)
\end{aligned}
の$f = nf_{\mathrm{s}} / N\, (n = 0, 1, ..., N - 1)$における値でした.ここで,$\otimes$は畳込み
\begin{aligned}
(G_{1} \otimes G_{2}) (f)
&= \int_{-\infty}^{\infty} G_{1}(\phi) G_{2}(f-\phi) \,\mathrm{d} \phi
\end{aligned}
を意味します.


この式は,DFTの枠組みで表現することもできます.まず,

\begin{aligned}
h_{\mathrm{s}} (t)
&= \sum_{k=-\infty}^{\infty} w(k T_{\mathrm{s}}) g(k T_{\mathrm{s}}) \delta(t - kT_{\mathrm{s}}) = w(t) g_{\mathrm{s}} (t) \\
&= \sum_{k=0}^{N - 1} w(k T_{\mathrm{s}}) g(k T_{\mathrm{s}}) \delta(t - kT_{\mathrm{s}})
\end{aligned}
であり,DFTは
\begin{aligned}
H_{\mathrm{s}} (n f_{\mathrm{s}} / N) = \sum_{k=0}^{N - 1} w(k T_{\mathrm{s}}) g(k T_{\mathrm{s}}) e^{-i 2\pi n k / N} \\
(n = 0,...,N-1)
\end{aligned}
と表されます.したがって,$w,g$のDFTをそれぞれ$\mathcal{F}_{\mathrm{DFT}}[w],\mathcal{F}_{\mathrm{DFT}}[g]$で表せば,積のDFTの表式(離散フーリエ変換(DFT) - Notes_JP)から
\begin{aligned}
& H_{\mathrm{s}} (n f_{\mathrm{s}} / N) \\
& = \frac{1}{N} \bigl(\mathcal{F}_{\mathrm{DFT}}[w] \otimes \mathcal{F}_{\mathrm{DFT}}[g] \bigr) (n f_{\mathrm{s}} / N)
\end{aligned}
が成り立ちます.ただし,$\otimes$は巡回畳込み
\begin{aligned}
(G_{1} \otimes G_{2}) (n)
&=\sum_{m = 0}^{N-1} G_{1}(m) G_{2}(n - m)
\end{aligned}
を意味します.

まとめ

アナログ信号$g(t)$の連続フーリエ変換$G(f)$を知りたいとします.ただし,フィルタをかけるなどして
\begin{aligned}
\begin{cases}
\, G(f) \neq 0 & (|f| < w / 2) \\
\, G(f) = 0 & (|f| \geq w / 2)
\end{cases}
\end{aligned}
になっているとします(帯域制限).

仮に,サンプリング周期$T_{\mathrm{s}}$(サンプリング周波数$f_{\mathrm{s}} = 1/T_{\mathrm{s}}$)の無限個のデータ$\{g(k T_{\mathrm{s}})\}_{k=-\infty}^{\infty}$からつくった信号

\begin{aligned}
g_{\mathrm{s}} (t)
&= \sum_{k=-\infty}^{\infty} g(k T_{\mathrm{s}}) \delta(t - kT_{\mathrm{s}})
\end{aligned}
の連続フーリエ変換(DTFT)
\begin{aligned}
&G_{\mathrm{s}} (f)
&= \frac{1}{T_{\mathrm{s}}} \sum_{n=-\infty}^{\infty} G (f - n f_{\mathrm{s}} )
\end{aligned}
を求められるとします(実際は計算できない!).この場合,$f_{\mathrm{s}} \geq w$となるようにサンプリング周波数を選べば,右辺の和は重ならないので,$G_{\mathrm{s}}$から目的の$G (f)$を求められます

しかし,現実で得られるのは有限個のデータ$\{g(k T_{\mathrm{s}})\}_{k=0}^{N - 1}$だけです.そこで,

\begin{aligned}
\begin{cases}
\, w (t) \neq 1 & (0 \leq t \leq (N - 1) T_{\mathrm{s}} ) \\
\, w (t) = 0 & (\mathrm{otherwise})
\end{cases}
\end{aligned}
を満たす関数(窓関数)をサンプリングデータに掛けた関数
\begin{aligned}
h_{\mathrm{s}} (t)
& =w(t) g_{\mathrm{s}} (t) \\
&= \sum_{k=0}^{N - 1} w(k T_{\mathrm{s}}) g(k T_{\mathrm{s}}) \delta(t - kT_{\mathrm{s}})
\end{aligned}
を考えます.この関数は有限個のデータ$\{g(k T_{\mathrm{s}})\}_{k=0}^{N - 1}$だけしか必要ない上,この関数の連続フーリエ変換の$f = nf_{\mathrm{s}} / N\, (n = 0, 1, ..., N - 1)$における値(DFT)
\begin{aligned}
H_{\mathrm{s}} (n f_{\mathrm{s}} / N)
= \sum_{k=0}^{N - 1} w(k T_{\mathrm{s}}) g(k T_{\mathrm{s}}) e^{-i 2\pi n k / N} \\
(n = 0,...,N-1)
\end{aligned}
はコンピュータで計算できます(FFT).

さらに,

\begin{aligned}
& H_{\mathrm{s}} (f) = \mathcal{F}[h_{\mathrm{s}}] (f) \\
&= \mathcal{F}[w] (f) \otimes G_{\mathrm{s}} (f)
\end{aligned}
という関係があるので,性質のよい窓関数を選べば,$G_{\mathrm{s}} (f) $のスペクトルを推測できます.したがって,目的の$G (f)$を推測できます.


参考文献