import matplotlib.pyplot as graph
import numpy as np
from numpy import fft
import math
import cv2
 
# 仿真運(yùn)動(dòng)模糊
def motion_process(image_size,motion_angle):
 PSF = np.zeros(image_size)
 print(image_size)
 center_position=(image_size[0]-1)/2
 print(center_position)
 
 slope_tan=math.tan(motion_angle*math.pi/180)
 slope_cot=1/slope_tan
 if slope_tan<=1:
  for i in range(15):
   offset=round(i*slope_tan) #((center_position-i)*slope_tan)
   PSF[int(center_position+offset),int(center_position-offset)]=1
  return PSF / PSF.sum() #對(duì)點(diǎn)擴(kuò)散函數(shù)進(jìn)行歸一化亮度
 else:
  for i in range(15):
   offset=round(i*slope_cot)
   PSF[int(center_position-offset),int(center_position+offset)]=1
  return PSF / PSF.sum()
 
#對(duì)圖片進(jìn)行運(yùn)動(dòng)模糊
def make_blurred(input, PSF, eps):
 input_fft = fft.fft2(input)# 進(jìn)行二維數(shù)組的傅里葉變換
 PSF_fft = fft.fft2(PSF)+ eps
 blurred = fft.ifft2(input_fft * PSF_fft)
 blurred = np.abs(fft.fftshift(blurred))
 return blurred
 
def inverse(input, PSF, eps):  # 逆濾波
 input_fft = fft.fft2(input)
 PSF_fft = fft.fft2(PSF) + eps #噪聲功率，這是已知的，考慮epsilon
 result = fft.ifft2(input_fft / PSF_fft) #計(jì)算F(u,v)的傅里葉反變換
 result = np.abs(fft.fftshift(result))
 return result
 
def wiener(input,PSF,eps,K=0.01):  #維納濾波，K=0.01
 input_fft=fft.fft2(input)
 PSF_fft=fft.fft2(PSF) +eps
 PSF_fft_1=np.conj(PSF_fft) /(np.abs(PSF_fft)**2 + K)
 result=fft.ifft2(input_fft * PSF_fft_1)
 result=np.abs(fft.fftshift(result))
 return result
 
image = cv2.imread('you.jpg')
image = cv2.cvtColor(image,cv2.COLOR_BGR2GRAY)
img_h=image.shape[0]
img_w=image.shape[1]
graph.figure(1)
graph.xlabel("Original Image")
graph.gray()
graph.imshow(image)  #顯示原圖像
 
graph.figure(2)
graph.gray()
#進(jìn)行運(yùn)動(dòng)模糊處理
PSF = motion_process((img_h,img_w), 60)
blurred = np.abs(make_blurred(image, PSF, 1e-3))
 
graph.subplot(231)
graph.xlabel("Motion blurred")
graph.imshow(blurred)
 
result = inverse(blurred, PSF, 1e-3) #逆濾波
graph.subplot(232)
graph.xlabel("inverse deblurred")
graph.imshow(result)
 
result=wiener(blurred,PSF,1e-3)  #維納濾波
graph.subplot(233)
graph.xlabel("wiener deblurred(k=0.01)")
graph.imshow(result)
 
blurred_noisy=blurred + 0.1 * blurred.std() * \
   np.random.standard_normal(blurred.shape) #添加噪聲,standard_normal產(chǎn)生隨機(jī)的函數(shù)
 
graph.subplot(234)
graph.xlabel("motion & noisy blurred")
graph.imshow(blurred_noisy)  #顯示添加噪聲且運(yùn)動(dòng)模糊的圖像
 
result = inverse(blurred_noisy, PSF, 0.1+1e-3) #對(duì)添加噪聲的圖像進(jìn)行逆濾波
graph.subplot(235)
graph.xlabel("inverse deblurred")
graph.imshow(result)
 
result=wiener(blurred_noisy,PSF,0.1+1e-3)   #對(duì)添加噪聲的圖像進(jìn)行維納濾波
graph.subplot(236)
graph.xlabel("wiener deblurred(k=0.01)")
graph.imshow(result)
 
graph.show()