[scikit-image] 98. 2Dフーリエ変換で窓関数を適用して不連続性を除去(skimage.filters window)

はじめに

skimage.filters の windowを画像のフーリエ変換に適用する方法について説明する。

解説

モジュールのインポートなど

	import matplotlib.pyplot as plt
	import numpy as np
	from scipy.fftpack import fft2, fftshift
	from skimage import img_as_float
	from skimage.color import rgb2gray
	from skimage.filters import window,gaussian

view raw window_func_img_fft.py hosted with ❤ by GitHub

バージョン

	#version
	import matplotlib
	import scipy
	import skimage
	print(np.__version__)
	print(matplotlib.__version__)
	print(scipy.__version__)
	print(skimage.__version__)
	1.21.2
	3.4.3
	1.7.1
	0.18.3

view raw window_func_img_fft.py hosted with ❤ by GitHub

画像の読み込み

	astrophytum = plt.imread("kabutomaru.jpg")
	image = img_as_float(rgb2gray(astrophytum))
	image.shape
	#(1419, 1428)

view raw window_func_img_fft.py hosted with ❤ by GitHub

アストロフィツム属の兜丸の写真を用いる。rgb2grayでグレースケール画像に変換する。

窓関数をかける

	gimage = gaussian(image,sigma=5)
	wimage = gimage * window('hann', image.shape)

view raw window_func_img_fft.py hosted with ❤ by GitHub

gaussianフィルターでノイズを滑らかにした後に窓関数をかける。利用できるのはscipyの窓関数であり、ここでは’hann’を用いる。そのほかに使用できる窓関数は下記ページに載っている。

Window functions (scipy.signal.windows) — SciPy v1.15.2 Manual

docs.scipy.org

フーリエ変換

	image_f = np.abs(fftshift(fft2(gimage)))
	wimage_f = np.abs(fftshift(fft2(wimage)))
	image_f_log = np.log(image_f)
	wimage_f_log = np.log(wimage_f)

view raw window_func_img_fft.py hosted with ❤ by GitHub

scipyのfft2でフーリエ変換をする。

結果の表示

	fig, axes = plt.subplots(2, 2, figsize=(8, 8),dpi=140)
	ax = axes.ravel()
	ax[0].set_title("Original image")
	ax[0].imshow(gimage, cmap='gray')
	ax[1].set_title("Windowed image")
	ax[1].imshow(wimage, cmap='gray')
	ax[2].set_title("Original FFT (frequency)")
	of = ax[2].imshow(image_f_log,cmap="turbo")
	plt.colorbar(of,ax=ax[2])
	ax[3].set_title("Window + FFT (frequency)")
	wf = ax[3].imshow(wimage_f_log, cmap='turbo')
	plt.colorbar(wf,ax=ax[3])
	plt.savefig("window_func_img_fft_1.png",dpi=100)
	plt.show()

view raw window_func_img_fft.py hosted with ❤ by GitHub

左上が元画像、右上が窓関数をかけた画像で、左下が元画像のフーリエ変換の結果で、右下が窓関数をかけた画像のフーリエ変換の結果となる。窓関数の効果により、十字模様が消えて、鮮明な結果が得られている。

縦に線分析

	pos=714
	fig = plt.figure(dpi=120,figsize=(8,6))
	ax = fig.subplot_mosaic("AB;CC")
	ax["A"].set_title("Original FFT (frequency)")
	of = ax["A"].imshow(image_f_log,cmap="turbo")
	ax["A"].plot([pos,pos],[0,1418],"ro-")
	plt.colorbar(of,ax=ax["A"])
	ax["B"].set_title("Window + FFT (frequency)")
	wf = ax["B"].imshow(wimage_f_log, cmap='turbo')
	plt.colorbar(wf,ax=ax["B"])
	ax["B"].plot([pos,pos],[0,1418],"ro-")
	line1 = image_f_log[:,pos]
	line2 = wimage_f_log[:,pos]
	ax["C"].plot(line1,label="Original",alpha=.8)
	ax["C"].plot(line2,label="Windowed",alpha=.8)
	ax["C"].legend()
	plt.savefig("window_func_img_fft_2.png",dpi=100)
	plt.show()

view raw window_func_img_fft.py hosted with ❤ by GitHub

縦方向のラインプロファイルは以下のようになる。

横に線分析

	pos=709
	fig = plt.figure(dpi=120,figsize=(8,6))
	ax = fig.subplot_mosaic("AB;CC")
	ax["A"].set_title("Original FFT (frequency)")
	of = ax["A"].imshow(image_f_log,cmap="turbo")
	ax["A"].plot([0,1427],[pos,pos],"ro-")
	plt.colorbar(of,ax=ax["A"])
	ax["B"].set_title("Window + FFT (frequency)")
	wf = ax["B"].imshow(wimage_f_log, cmap='turbo')
	plt.colorbar(wf,ax=ax["B"])
	ax["B"].plot([0,1427],[pos,pos],"ro-")
	line1 = image_f_log[pos,:]
	line2 = wimage_f_log[pos,:]
	ax["C"].plot(line1,label="Original",alpha=.8)
	ax["C"].plot(line2,label="Windowed",alpha=.8)
	ax["C"].legend()
	plt.savefig("window_func_img_fft_3.png",dpi=100)
	plt.show()

view raw window_func_img_fft.py hosted with ❤ by GitHub

横方向のラインプロファイルは以下のようになる。

コードをダウンロード(.pyファイル)

コードをダウンロード(.ipynbファイル)

参考

Using window functions with images — skimage 0.26.0rc0.dev0 documentation

scikit-image.org

fft2 — SciPy v1.15.2 Manual

docs.scipy.org