Python으로 Fourier filter씌워 이미지 재구성해보기

지난 나의 글 - https://physikk.tistory.com/35 에서 Fouirer filter를 소개하면서 "이미지에서 명암의 변화를 기준으로 Fourier transform을 수행하여 다시 이미지를 재구성 할 수 있는데, 이 때 고주파 성분을 저주파수 성분 대비 상대적으로 증폭시키면 원래 이미지의 가장자리를 보다 선명하게 만들 수 있다. 반대로 고주파 성분을 제거하면 그림이 뿌옇게 되며 blur 효과를 얻을 수 있다." 고 설명했다.

이번 게시글에서는 Google colab. 환경에서 python을 이용해 iu사진을 직접 푸리에 변환한 뒤, 필터링해서 이미지를 재구성 해 보겠다.

컴퓨터 비전 처리에서는 다음 사진의 주인공 - 레나가 많이 쓰인다.

위키피디아에 따르면, 세밀함, 평면, 그림자, 질감이 적절하게 조화되어 아주 적합한 테스트용 샘플이라고 하는데, 나는 그냥 아이유 사진을 가지고 할 예정이다.

먼저 다음의 모듈을 import 해주자

import glob
from skimage import io
import matplotlib.pyplot as plt
import numpy as np

우리는 Numpy module에 있는 FFT2를 이용해서 Fourier transformation 할 것이다. 이에 관련된 기술적이고 상세한 설명은 다른 블로그에도 잘 작성되어있으니, 우리는 우선 써보겠다. 먼저 이미지를 불러와준다.

#image import  
our_file=glob.glob('/content/drive/MyDrive/iu.png')[0]
#showing original image
io.imshow(our_file)
plt.show()
img = io.imread(our_file, as_gray=True)

여기서 포인트가 두 가지 정도 있다.

1. 이미지를 흑백으로 전처리하는 이유는, 이미지의 Noise를 감소하기 위함이다. 노이즈가 심하다고 판단되는 경우에는 가우시안 블러처리까지 해주는 경우도 있다.

2. 그럼 왜 흑백처리하기 전에 사진을 불러왔나? - 앞으로 계속 흑백으로 볼꺼라서, 마지막으로 컬러로 한 번 감상하기 위함이다.

이제 Fourier Transformation해주자.

fig=plt.figure(figsize=(12,8))
f = np.fft.fft2(img)
fshift = np.fft.fftshift(f)
magnitude_spectrum = 20*np.log(np.abs(fshift))
phase = np.angle(fshift)

plt.subplot(131),plt.imshow(img, cmap = 'gray')
plt.title('Input Image',fontsize = 14, fontweight='bold'), plt.xticks([]), plt.yticks([])
plt.subplot(132),plt.imshow(magnitude_spectrum, cmap = 'gray')
plt.title('Magnitude Spectrum',fontsize = 14, fontweight='bold'), plt.xticks([]), plt.yticks([])
plt.subplot(133),plt.imshow(phase, cmap = 'gray')
plt.title('Phase Spectrum',fontsize = 14, fontweight='bold'), plt.xticks([]), plt.yticks([])
plt.savefig('iu',dpi=300,bbox_inches='tight')

plt.show()

numpy모듈의 fft2를 이용해서 푸리에 변환해주고, 같은 모듈의 fftshift를 이용해서 원점을 왼쪽 위에서 센터로 가져왔다. 그리고 magnitude spectrum, phase spectrum, original을 subplot으로 한 눈에 비교할 수 있게 그려보았다. Phase spectrum은 별 의미가 없고, Magnitude spectrum을 이용해서 filtering을 해보자.

원점을 센터로 가져왔으니, 가운데로 갈 수록 저주파수이고, 가장자리로 갈 수록 고주파수이다. 먼저 다음과 같이 가운데 (저주파수)를 가려보자.

fig=plt.figure(figsize=(8,8))
(w, h) = fshift.shape
half_w, half_h = int(w/2), int(h/2)
# high pass filter
n = 15
fshift[half_w-n:half_w+n+1,half_h-n:half_h+n+1] = 0 # select first 30x30 (low) frequencies
plt.imshow( (20*np.log10( 0.1 + fshift)).astype(int),cmap='gray')
plt.axis('off')
plt.savefig('highfrequencyfilter',dpi=300,bbox_inches='tight')
plt.show()

w by h 의 matrix shape에서 w/2, h/2로 원점을 찾아주고, 거기서부터 -15, +15하여 첫 30번째까지의 주파수들을 0으로 세팅해준다. 그러면 이렇게 된다.

이를 이제 이미지로 재조합(?)할 차례이다. inverse Fourier transformation해주면 된다. 고주파수만 남았으니까, hiu라고 칭하겠다. imaginary part는 쓸데없으니 버려준다.

hiu = np.fft.ifft2(fshift).real
plt.figure(figsize=(10,10))
plt.subplot(111),plt.imshow(hiu, cmap = 'gray')
plt.title('Blocked low frequency image',fontsize = 14, fontweight='bold'), plt.xticks([]), plt.yticks([])
plt.axis('off')
plt.savefig('hiu',dpi=300,bbox_inches='tight')
plt.show()

보다싶이 테두리가 진한(?) 부분들만 남았다! 원래사진과 비교하면 다음과 같다.

이번에는 고주파수를 없에보자. 아까 했던 방식대로인데, 이번에는 뒤에서부터 30번째 주파수를 없엘것이다.

f = np.fft.fft2(img)
fshift = np.fft.fftshift(f)
magnitude_spectrum = 20*np.log(np.abs(fshift))
phase = np.angle(fshift)
fig=plt.figure(figsize=(8,8))
(w, h) = fshift.shape
# high pass filter
n = 30
fshift[0:n,:] = 0 
fshift[:,0:n] = 0 
fshift[w-n:w,:] = 0 
fshift[:,h-n:h] = 0 

plt.imshow( (20*np.log10( 0.1 + fshift)).astype(int),cmap='gray')
plt.axis('off')
plt.savefig('lowfrequencyfilter',dpi=300,bbox_inches='tight')
plt.show()

이를 이미지로 재조합하자. 저주파수만 남겼으니까, liu라고 하겠다.

liu = np.fft.ifft2(fshift).real
plt.figure(figsize=(10,10))
plt.subplot(111),plt.imshow(liu, cmap = 'gray')
plt.title('Blocked high frequency image',fontsize = 14, fontweight='bold'), plt.xticks([]), plt.yticks([])
plt.axis('off')
plt.savefig('liu',dpi=300,bbox_inches='tight')
plt.show()

hiu와 달리 고주파수가 날아간 liu는 테두리부분이 날아간 것을 확인할 수 있다. 그런데 위에서 블러처리가 되어야 한다고 했는데, 그림이 그렇게 블러가 들어간거같지는 않다고 생각할 수 있다. 화끈하게(?) (500x500)에서 440x440을 0으로 처리해보면 다음과 같은 결과가 나온다.