OpenCVの基底を使った画像の変換

ここではフーリエ変換だけを扱う

フーリエ変換

目的

このチュートリアルでは

  • OpenCVを使って画像のフーリエ変換を計算する方法を学ぶ.
  • NumpyのFFTを使う方法を学ぶ.
  • フーリエ変換を使ったアプリケーションについて学ぶ.
  • 次の関数の使い方を学ぶ : cv2.dft(), cv2.idft()

理論

フーリエ変換は種々のフィルタの周波数特性を解析するために使われる.画像に対しては 2次元離散フーリエ変換 (DFT) を使って周波数領域に変換する.高速化されたアルゴリズムである 高速フーリエ変換 (FFT) はDFTの計算に使いる.これらのアルゴリズムの詳細については信号処理や画像処理の教科書を参照すること. 補足資料 の章に幾つか参考文献を挙げている.

正弦波を $x(t) = A \sin(2 \pi ft)$ と書く.ここで $$f は信号の周波数を表す.この信号を周波数領域で観測すれば,周波数 f の点にspikeが見られる.離散信号を形成するために信号を標本化すると,同じ周波数領域での信号を得られるが,周波数領域での信号は $[-\pi, \pi]$ の範囲もしくは $[0,2\pi]$ の範囲での周期性を持つ信号とみなされる(N点DFTであれば $[0,N]$ の範囲). 画像を2方向に標本化された信号とみなすことができる.横方向と縦方向にフーリエ変換を計算すれば,画像を周波数領域で表現できる.

より直観的に言うと,短時間にある正弦波信号の振幅変化が速く起こるとき、それを高周波信号と言い、振幅変化が遅ければ低周波信号と呼ぶ.全く同じ考えを画像に対して拡張する.画像中で振幅変化が急激に生じる場所はどこだろうか? それはエッジやノイズである.つまり,エッジやノイズは画像の高周波成分に対応している.振幅の変化が大きくなければ低周波成分になる( 周波数変換の直観的な説明をしている資料を 補足資料 の欄に幾つか挙げておく).

それではフーリエ変換の計算方法を見ていくことにする.

Numpyを使ったフーリエ変換

まず初めにNumpyを使ったフーリエ変換の計算方法を見ていこう.NumpyはFFTを計算するための関数 np.fft.fft2() を用意している.この関数は複素数型の配列を出力する.第1引数は入力画像をグレースケール画像として与える.第2引数は出力配列のサイズを指定するが,オプションである.指定するサイズが入力画像のサイズより大きければ入力画像はFFTの計算をする前にゼロパディングされる.入力画像のサイズより小さければ入力画像を切り取る.何も指定されなければ出力配列のサイズは入力画像のサイズと同じになる.

実際に計算をしてみると,周波数領域の原点(直流成分)が画像の左上の角に位置するようになる.直流成分を画像中心に移動させたければ,スペクトル全体を $\displaystyle \frac{N}{2}$ 両方向のずらす必要がある.この移動には np.fft.fftshift() 関数を使いる(これでより解析がしやすくなる).一度フーリエ変換を計算すれば,スペクトルの大きさが分かる (画像):

In [5]:
%matplotlib inline
import cv2
import numpy as np
from matplotlib import pyplot as plt

img = cv2.imread('messi5.jpg',0)
f = np.fft.fft2(img)
fshift = np.fft.fftshift(f)
magnitude_spectrum = 20*np.log(np.abs(fshift))

plt.figure(figsize=(15,5))
plt.subplot(121),plt.imshow(img, cmap = 'gray')
plt.title('Input Image'), plt.xticks([]), plt.yticks([])
plt.subplot(122),plt.imshow(magnitude_spectrum, cmap = 'gray')
plt.title('Magnitude Spectrum'), plt.xticks([]), plt.yticks([])
plt.show()

中心に白い領域が集中している事が分かる.これは画像が低周波成分をより多く含んでいることを意味する.

これでフーリエ変換を見ることができた.これによりハイパスフィルタといった周波数領域での処理ができるようになる.低周波成分に対して矩形windowを使ったマスク処理をすることによって低周波成分を取り除く事ができる.それから np.fft.ifftshift() 関数を使って直流成分の位置を画像の左上に戻し、 np.ifft2() 関数を使って逆フーリエ変換を適用する.最終的な結果は複素数型の配列になるので,その絶対値をとる(画像):

In [8]:
%matplotlib inline
import cv2
import numpy as np
from matplotlib import pyplot as plt

img = cv2.imread('messi5.jpg',0)
f = np.fft.fft2(img)
fshift = np.fft.fftshift(f)
magnitude_spectrum = 20*np.log(np.abs(fshift))

rows, cols = img.shape
crow,ccol = rows//2 , cols//2
fshift[crow-30:crow+30, ccol-30:ccol+30] = 0
f_ishift = np.fft.ifftshift(fshift)
img_back = np.fft.ifft2(f_ishift)
img_back = np.abs(img_back)

plt.figure(figsize=(15,10))
plt.subplot(221),plt.imshow(img, cmap = 'gray')
plt.title('Input Image'), plt.xticks([]), plt.yticks([])
plt.subplot(223),plt.imshow(img_back, cmap = 'gray')
plt.title('Image after HPF'), plt.xticks([]), plt.yticks([])

plt.subplot(224),plt.imshow(img_back)
plt.title('Result in JET'), plt.xticks([]), plt.yticks([])

plt.show()