1. 概要

周波数フィルタリング (Spectral Filtering)とは，画像を二次元フーリエ変換（2D Fourier Transform）を用いて，周波数領域側のスペクトル画像に対してフィルタリング処理をおこなう，画像フィルタリングの１種である．

この記事では，画像の周波数フィルタリングの処理手順のうち，二次元離散フーリエ変換（２節）と，その応用例として帯域通過フィルタ（3節）についてのみ述べる．

周波数フィルタリングの前には，まず前処理として，二次元フーリエ変換を用いて周波数領域（frequency domain）の周波数スペクトル（frequency spectrum）である複素値画像への変換をおこなう．二次元フーリエ変換により，画像の輝度値は，周波数ごとの正弦波成分に分解される二次元周波数スペクトル画像となる（2節）．よって，一次元信号のフーリエ変換後と同様に，フーリエ変換後は「帯域通過フィルタ処理」が可能となる．例えばローパスフィルタによるノイズ削減ができたり（3節），高周波数帯域のみ強調する「コントラスト強調」が行えるなど，「空間フィルタリングより単純計算で済み，周波数領域でのフィルタリングが実現できる」という利点が出てくる．

また，画像の周波数領域における処理は，離散フーリエ変換の派生系である「離散コサイン変換」を用いた画像の圧縮などにも応用される，「画像・映像の符号化」にも関連した技術である．他にも，劣化画像の復元や，医療画像再構成などの復元系処理でも，画像のフーリエ変換は応用されている．

ちなみに，この記事中では一次元フーリエ変換について多少は復習的な記述はおこなうものの，読者はすでに理解しているものとみなす．よって，一次元フーリエ変換について知らない方は， [References > 参考書籍] に挙げた書籍等や他Webサイト等で勉強していただきたい．

親記事：畳み込み(convolution): 画像フィルタリングの基本演算

1.1 記事の構成

2節以降は，以下の構成にて，画像のフーリエ変換と帯域通過フィルタの応用について紹介していく．

• 2節　画像の周波数フィルタリング
  • 2.1 周波数フィルタリングの手順
  • 2.2 二次元フーリエ変換
  • 2.3 逆二次元フーリエ変換
• 3節 帯域通過フィルタ
• 4節 まとめ

2. 画像の周波数フィルタリング

この2節では，基本的な「二次元デジタル画像の周波数フィルタリング」について，手順とその原理を紹介する．

2.1 周波数フィルタリングの手順

まず入力をグレー画像$f \in \mathbb{R}^{w \times h}$とする（$w$と$h$は入力画像の幅と高さ）．

周波数フィルタリングの処理は３段階で行われる：

1. 【フーリエ変換】：二次元フーリエ変換$\mathscr{F}$により，複素周波数スペクトル画像の$F \in \mathbb{C}^{w \times h}$へ変換する（図1 左部矢印）．
2. 【周波数フィルタリング】次に$F$に対して画像と同サイズの任意のフィルタ$H \in \mathbb{C}^{w \times h}$を用いて周波数フィルタリングをおこなう（図1 下部矢印）．
3. 【逆フーリエ変換】逆二次元フーリエ変換$\mathscr{F}^{-1}\{G(u,v)\}$により，$G \in \mathbb{R}^{w \times h}$をフィルタリング済みの画像$g \in \mathbb{R}^{w \times h}$へ変換する（図1 右部矢印）．

空間フィルタリングでは，元の画像空間における座標$(x,y)$に基づく，フィルタ窓内での畳み込み処理をおこなのであった （図1 上部矢印）．それに対し，周波数フィルタリングでは，二次元周波数空間$(u,v)$において，周波数スペクトル$F$のフィルタリングを行う．

一次元信号の時間軸$t$の代わりに，画像のフーリエ変換では，画像信号の$x$軸や$y$軸の空間方向を代わりに時系列方向とみなすので，$(u,v)$のことを空間周波数(spatial frequency)とよぶ．つまり，二次元周波数空間上にある信号の（通常は）グレー画像を，二次元周波数スペクトルに変換するのが，二次元画像でのフーリエ変換である．

以上より，周波数フィルタリングの3手順を，今度は次節以降の式形式と対応させやすい「画素値での表示」で再度示しておく：

1. フーリエ変換：
   • 二次元フーリエ変換$\mathscr{F}$により，入力画像の全座標の輝度を用いて，周波数ドメイン画像上の座標$(u,v)$における値$F(u,v)= \mathscr{F}\{f(x,y)\}$へ写像する(2.2節，式(2.3))．
   • $u$ は $x$ 方向の周波数で，$v$ は $y$ 方向の周波数とする．
2. 周波数フィルタリング：
   • $F(u,v)$の$H(u,v)$による周波数フィルタリングを実施
   • $G(u,v) = F(u,v) \cdot H(u,v)$を得る(3節)．
3. 逆フーリエ変換：
   • $G(u,v)$に逆二次元フーリエ変換$\mathscr{F}^{-1}\{G(u,v)\}$をおこない，元の空間ドメインの画像の輝度値 $g(x,y)$へと逆写像する(2.3節，式(2.5))．

※ この３手順中には，周波数フィルタリングの前処理と後処理として「象限の入れ替え処理（2.2.2節）」も実施するが，ここでは全体の見通しをよくするために省いた．

2.1.1 フーリエ対：画像領域の計算と，フーリエ変換の計算の対応関係

畳み込みに用いる2つの一次元信号(関数)の値 $x(t),h(t)$間には，以下の畳み込み定理が成立する：

\begin{equation}
y(t) = x(t) * h(t) \textcolor{blue}{\leftrightarrow} Y(w) = X(w) \cdot H(w) \tag{2.1}
\end{equation}

この「2つの関数の畳み込みは，フーリエ変換すると要素積の演算に変わる」関係のおかげで，手順2の「周波数フィルタリング」は，周波数スペクトル側で要素積 $G(u,v) = F(u,v) \cdot H(u,v)$を計算するだけ済む．

したがって，式(2.1)の画像上の二次元空間$(x,y)$上での空間フィルタリング(畳み込み)は，周波数スペクトルの二次元空間$(u,v)$上の周波数フィルタリング（要素積）と，お互い等価な計算であることがわかる：

\begin{align}
\mathscr{F}\{ f(x,y) * h(x,y)\}
&=\mathscr{F}\{f(x,y)\} \cdot \mathscr{F}\{h(x,y)\} \\
&=F(u,v) \cdot H(u,v)\tag{2.2}
\end{align}

したがって（前処理のフーリエ変換と後処理の逆フーリエ変換の計算コストも含めた上で），画像の空間フィルタリングよりも，周波数フィルタリングの全計算コストが少なく済むケースが多いという利点がある．詳しくはFFTの記事やCNNの畳み込み高速化の記事などで述べるつもりだが，空間フィルタのカーネルサイズが大きいほど，フーリエ変換+周波数フィルタリングで等価な処理を行った方が計算コストが抑えられて有利になる．つまり，等価なフィルタ処理同士なら，周波数フィルタリングを行った方が計算が速く済むことが多い．

また，二次元フーリエ変換後は，「周波数成分ごとに分離されている周波数スペクトル信号」になるので，ローパスフィルタなど「周波数空間であることを活用する処理」も実施できる利点も出てくる（3節）．

2.2 ニ次元フーリエ変換：Spectral Filtering の前処理

2.1節と図１で既にみたように，周波数フィルタリングでは，前処理として二次元フーリエ変換（2D Fourier Transform）を実施したのちに，その結果の$F$に対して周波数フィルタリングをおこなう．

二次元フーリエ変換は，入力の二次元信号を複素スペクトル値へ変換する．具体的には，以下の変換式により，$f(x,y)$ で重み付けされた基底関数(複素正弦波)の和を取ることで，変換後の画素値$F(u,v)$を計算する：

\begin{align}
F(u,v)=\frac{1}{w \cdot h} \sum_{x=0}^{w-1}\sum_{y=0}^{h-1} \underbrace{f(x,y)}_{\text{重み係数}} \cdot \underbrace{\exp \left[ -i \cdot 2\pi \cdot \left(\frac{xu}{w} + \frac{yv}{h}\right) \right]}_{\text{基底関数}} \tag{2.3}
\end{align}

この式(2.3)の変換を，出力画像の周波数 $u = \{0,1, \ldots , w-1 \}$と $v = \{0,1, \ldots , h-1 \}$の全座標$(u,v)$において計算することで，フーリエ変換された結果の複素スペクトル画像 $F \in \mathbb{C}^{w \times h}$ が出来上がる．

ここで，式(2.3)の解釈のために，ß任意の実数 $\alpha$に対して，以下のオイラーの公式を考えてみる：

\begin{equation}
\exp(i \cdot \alpha) = \cos \alpha + i \cdot \sin \alpha \tag{2.4}
\end{equation}

この式から，右辺の元は別々の正弦波$\cos \alpha$と余弦波$\sin \alpha$の2つの信号を，左辺の1つの複素正弦波の１点 $\exp(i \cdot \alpha)$ としてとらえることができる．すると，二次元フーリエ変換[式(2.3)]は「正弦波と余弦波の重み付き和」を計算し，それを複素数値として出力していると解釈できる．

2.2.1 二次元フーリエ変換での振幅と位相

しかし$F$が複素数のままだと，虚部があるので二次元配列全体を画像として可視化もできないし，この後のフィルタリング処理的にも扱いづらい．

そこで，一次元フーリエ変換と同様に，$F$が複素正弦波であることを生かして，元の複素数を，「振幅」と「位相」の2つの実数表現に変換する．具体的には，$F$を以下の２チャンネルの画像にそれぞれ変換したのちに，それらの画像を可視化する：

• 振幅スペクトル画像：$F$の振幅値$|F(u,v)|$を，各画素に格納した画像．
• 位相スペクトル画像：$F$の位相値$arg(F(u,v))$を，各画素に格納した画像．

また，振幅や位相に対する処理をなったりすることもある．（この記事の帯域フィルタリングでは複素数$F$のまま処理を行なうので，可視化目的でしか振幅スペクトルへの変換は行わない）．

2.2.2 象限の入れ替え処理

二次元フーリエ変換では，変換後の複素信号$F$の周波数軸$(u,v)$に，それぞれ$2\pi$の周期性を仮定している（図2-a 左図）．従って，$\alpha$が $\alpha = [0, 2\pi)$の範囲をはみ出した場合，$\alpha$が範囲内に留まるように0から再び循環した周波数の値へ平行移動させる．しかし，二次元フーリエ変換直後の複素画像$F$は，画像の4隅が低周波で，真ん中が高周波の配置で得られるが，この配置だと周波数範囲を区切るための単純なマスク図形を作りづらい．

そこで，周波数フィルタリングを実施する前に，$F$の4象限を，平行移動で入れ替える処理を通常行なう (図2-a)．また，周波数フィルタリング後は，逆フーリエ変換前に，入れ替えた象限を元に戻す処理を行ってから逆フーリエ変換を行なう．つまり，図1では描かなかったが，帯域周波数フィルタリングの前と後に，この「象限の入れ替え処理」と「象限を元に戻す処理」を実施する．具体的には，Matlabやnumpyなどでは，象限入れ替えと象限戻しの関数として，fftshift関数とifftshift関数が，それぞれ提供されている．

象限入れ替えを行った結果作られた$F^{\prime}$は，振幅スペクトル$|F^{\prime}|$を可視化すると，いつもフーリエ変換の解説でよく見るような，「中央に低周波数係数がきて，4隅に高周波数成分が来る配置」の振幅スペクトル画像になっていることがわかる（図2-a 右図）

2.3 逆ニ次元フーリエ変換：Spectral Filtering の後処理

周波数フィルタリング処理後の複素画像$G$には，画像の逆二次元フーリエ変換（Inverse 2D Fourier Transform）をおこない，空間領域（spatial domain）の画像gへと戻す．ここでは順フーリエ変換との対称性をみたいので，$F$から$f$への逆二次元フーリエ変換を行なうことを考えたい．

以下の変換式により，$F$の全画素$F(u,v)$を用いて，画素値$f(x,y)$への逆二次元フーリエ変換をおこなうことができる（2.2の逆変換）：

\begin{align}
f(x,y)= \sum_{u=0}^{v-1}\sum_{y=0}^{h-1} \underbrace{F(u,v)}_{重み係数} \underbrace{\exp \left[i \cdot 2\pi \cdot \left(\frac{xu}{w} + \frac{yv}{h}\right) \right]}_{基底関数} \tag{2.5}
\end{align}

逆フーリエ変換の式(2.5)は，基底関数の和として各画素$f(x,y)$を再構成することで，元の画像空間信号へと戻す変換である．$F(u,v)$は逆フーリエ変換時の基底関数の重み係数である．

こうして式(2.3)と式(2.7)により，元の座標空間と，周波数空間を行き来することができる．

2.3.1 逆2次元フーリエ変換の基底関数

式(2.3)は，基底関数の重み付き積として表現されているのだが，基底関数がどういう関数であり，式(2.3)は基底関数からどのように構成されている変換であるのかをこの節で抑えておきたい．

次に，式(2.4)から逆二次元フーリエ変換の基底関数（basis function）を以下のように定義できる：

\begin{equation}
\phi(x,y;u,v) = \exp[i \cdot 2\pi \cdot(ux +vy)] =\exp(i \alpha)\tag{2.6}
\end{equation}

すると，基底関数の式(2.5)を用いて，式(2.5)の逆二次元フーリエ変換を以下のように書き直すことができる：

\begin{align}
f(x,y) = \frac{1}{w \cdot h} \sum_{x=0}^{w-1}\sum_{y=0}^{h-1} F(u,v)\cdot \phi(x,y;u,v) \tag{2.7}
\end{align}

このように，二次元画像$f$の離散空間座標全体において「重み係数（元信号）と複素正弦波の積」の画像空間全体における積分（重み係数$f(x,y)$でスケールした基底関数$\phi(x,y;u,v)$の積分）が，二次元フーリエ変換だと捉えられるようになった（※ 同様に，逆一次元フーリエ変換でも，「重み係数群のベクトル」と「基底関数群のベクトル」の内積を計算することで，$t$方向に対する離散積分値を取るのであった）．

元画像を再構成するための，各基底関数の重み係数である$F(x,y)$はフーリエ係数（Fourier coefficients）と呼ばれる．

既に述べてきた通りであるが，この象限入れ替え後の$F$の振幅スペクトル$|F|$では，中央から遠くなる位置ほど，同心円の法線方向に傾いている二次元正弦波が基底関数として配置されている（図2-b）．フーリエ変換に寄って，これらの基底関数群の重み付け和として，元画像が表現されている

よって，象限入れ替え後は，$F(u,v)$に対して「円を用いた単純な帯域通過マスク」(3節)による周波数フィルタリングを使用しやすくなっている．

3. 帯域通過フィルタへの応用

バンドパスフィルタは，電気分野や，デジタル音響処理などの1次元信号波形に対して，フーリエ変換後に行う「特定周波数の信号のみをフィルタリングで通過させる」処理である．必要な周波数帯域の信号のみを残して不要な他の周波数帯域の信号を除去することで，特定帯域の信号だけを抽出するのが目的である．

画像でも同様に，二次元フーリエ変換した画像に象限移動（2.2.2節）を行った$F$に対して，帯域通過フィルタ(band pass filter)のカーネル$H$により，特定の周波数帯域（バンド）の信号だけを残した画像をつくることができる．

スペクトル画像$F(u,v)$に適用する，2次元バンドパスフィルタ$H(u,v)$には，主に以下の3つに分類できる：

• ローパスフィルタ (低域周波数のみ通過．(図3-a)
• ハイパスフィルタ (高域周波数のみ通過．(図3-b）
• バンドパスフィルタ (特定の中域周波数のみ通過．(図3-c）

これらのように，「通過帯域だけ残す二値マスク画像」を設計するだけで，バンドパスフィルタを簡単に設計することができる．目当ての通過させて残らせたい周波数帯域の信号のみ，マスクをかけて通過させ，それ以外の帯域の信号は遮蔽して0にする．画像処理で行うゆえ，カーネルの各画素の値$H(u,v)$には，1 or 0の画素値のみを持つマスク画像を用いる(図3)．図3の上半分の3つの図は，フィルタの振幅スペクトル値を$(u,v)$座標系で眺めたものである．ピンクの領域がマスクの値1に相当し，他は値は0である．下半分の3つの図は，上半分の各マスクの振幅スペクトル値を$v=0$の断面からみたものである．

よって，2.1節の手順2で既に見たように，以下の単純な「2値マスクによるマスク計算」により，周波数フィルタリング後のスペクトル $G$ を得ることができる：

\begin{equation}
G(u,v) = F(u,v) \cdot H(u,v) \tag{3.1}
\end{equation}

この$G$に象限を戻す処理を行った後，逆フーリエ変換を行なうと，周波数フィルタリングが終了した画像が手に入り，2.1節（図1）の全ての処理が終了する．

3.1 ローパスフィルタとハイパスフィルタの結果画像

画像におけるバンドパスフィルタの機能を確認したい．そこで，ローパスフィルタとハイパスフィルタを用いて，画像の周波数フィルタリングをおこなった結果を可視化した（図4）．

図3の中間部では，まず入力画像の離散フーリエ変換を行い複素スペクトル画像$F$を計算し，その$F$の振幅スペクトルを可視化した

※ この図4でも，単純化によって処理全体を見やすくしやすくするため，「象限の入れ替え処理」を図1同様に省略した．

ここで，複素スペクトル画像$F$にローパスフィルタを適用したのち，逆離散フーリエ変換をおこなうと，低周波数成分のみが残った画像が得られる(図4 – 左下)．これは，ブラー系の空間フィルタを使用した結果と同様に，輝度の変換が少ない領域や，変化しても変化周期が長い「低周波数領域」のみが残った画像が得られる．

一方，フーリエ変換 $F(u,v)$ にハイパスフィルタを適用したのち，逆離散フーリエ変換をおこなうと，高周波数成分のみが残った画像が得られる(図4 – 右下)．これはエッジ抽出空間フィルタやラプラシアンフィルタを適用した画像のような，高域の細かいテクスチャ成分とエッジ候補領域だけが残った画像を得ることができる．

4. まとめ

この記事では，二次元フーリエ変換を用いた，画像の「周波数フィルタリング」の基礎的な内容と処理手順についてまとめた．

また．周波数フィルタリングの基本的な例として帯域通過(バンドパス)フィルタについて紹介した．

参考書籍

• ディジタル画像処理 [改訂第二版] , ディジタル画像処理編集委員会, 公益財団法人画像情報教育振興協会（CG-ARTS協会）, 2020.
• 音響学入門ペディア，日本音響学会(編), コロナ社, 2017
• 道具としてのフーリエ解析, 涌井 良幸，涌井 貞美，日本実業出版社，2014
• Concise Computer Vision, Reinhard Klette, Springer, 2014.

References

参照外部リンク

• Discrete Fourier Transform – Wikipedia
• バンドパスフィルタ – Wikipedia
• OpenCV-Python tutorial – フーリエ変換
• Fourier Transform – Wikipedia
• 新潟大学 医用画像工学実験 ，5. 周波数領域における画像処理【フーリエ変換の概要】
• 離散フーリエ変換(DFT:Discrete Fourier Transform)　| イメージングソリューション

関連記事

• 親記事：畳み込み(convolution): 画像フィルタリングの基本演算
• 空間フィルタリング(Spatial Filtering)【画像フィルタリング その１】