第5章 Python 数字图像处理(DIP) - 图像复原与重建15

标题

最小均方误差（维纳）滤波

最小均方误差（维纳）滤波

目标是求未污染图像fff的一个估计f^\hat{f}f^，使它们之间的均方误差最小。
e2=E{(f−f^)2}(5.80)e^2 = E \big\{(f - \hat{f})^2 \big\} \tag{5.80}e2=E{(f−f^)2}(5.80)

误差函数的最小值在频率域中的表达比如下：
F^(u,v)=[H∗(u,v)Sf(u,v)Sf(u,v)∣H(u,v)∣2+Sη(u,v)]G(u,v)=[H∗(u,v)∣H(u,v)∣2+Sη(u,v)/Sf(u,v)]G(u,v)=[1H(u,v)∣H(u,v)∣2∣H(u,v)2+K]G(u,v)(5.81)\begin{aligned} \hat{F}(u, v) & = \Bigg[\frac{H^*(u, v) S_{f}(u, v)}{S_{f}(u, v)|H(u, v)|^2 + S_{\eta}(u, v)} \Bigg] G(u, v) \\ & = \Bigg[\frac{H^*(u, v) }{|H(u, v)|^2 + S_{\eta}(u, v) / S_{f}(u, v)} \Bigg] G(u, v) \\ & = \Bigg[\frac{1}{H(u,v)} \frac{|H(u,v)|^2}{|H(u,v)^2 + K} \Bigg]G(u,v) \end{aligned} \tag{5.81}F^(u,v)=[Sf(u,v)∣H(u,v)∣2+Sη(u,v)H∗(u,v)Sf(u,v)]G(u,v)=[∣H(u,v)∣2+Sη(u,v)/Sf(u,v)H∗(u,v)]G(u,v)=[H(u,v)1∣H(u,v)2+K∣H(u,v)∣2]G(u,v)(5.81)

注：逆滤波与维纳滤波都要求未退化图像和噪声的功率谱是已知的。

# 运动模糊PSF与谱
fft_shift = np.fft.fftshift(PSF)
fft = np.fft.fft2(PSF)
spect = spectrum_fft(fft)
plt.figure(figsize=(8, 8))
plt.subplot(1,2,1), plt.imshow(PSF, 'gray')
plt.subplot(1,2,2), plt.imshow(spect, 'gray')
plt.show()

# 仿真运动模糊
def motion_process(image_size, motion_angle, degree=15):"""This function has some problem"""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_tanif slope_tan<=1:for i in range(degree):offset=round(i*slope_tan)    #((center_position-i)*slope_tan)PSF[int(center_position+offset),int(center_position-offset)]=1return PSF / PSF.sum()  #对点扩散函数进行归一化亮度else:for i in range(degree):offset=round(i*slope_cot)PSF[int(center_position-offset),int(center_position+offset)]=1return PSF / PSF.sum()def get_motion_dsf(image_size, motion_angle, motion_dis):"""Get motion PSFparam: image_size: input image shapeparam: motion_angle: blur motion angleparam: motion_dis: blur distant, the greater value, more blurredreturn normalize PSF"""PSF = np.zeros(image_size)  # 点扩散函数x_center = (image_size[0] - 1) / 2y_center = (image_size[1] - 1) / 2sin_val = np.sin(motion_angle * np.pi / 180)cos_val = np.cos(motion_angle * np.pi / 180)# 将对应角度上motion_dis个点置成1for i in range(motion_dis):x_offset = round(sin_val * i)y_offset = round(cos_val * i)PSF[int(x_center - x_offset), int(y_center + y_offset)] = 1return PSF / PSF.sum()    # 归一化# 对图片进行运动模糊
def make_blurred(input, PSF, eps):"""blurred image with PSFparam: input: input imageparam: PSF: input PSF maskparam: eps: epsilon, very small value, to make sure not divided or multiplied by zeroreturn blurred image"""input_fft = np.fft.fft2(input)              #  image FFTPSF_fft = np.fft.fft2(PSF)+ eps             # PSF FFT plus epsilonblurred = np.fft.ifft2(input_fft * PSF_fft) # image FFT multiply PSF FFTblurred = np.abs(np.fft.fftshift(blurred))return blurreddef inverse_filter(input, PSF, eps): """inverse filter using FFT to denoiseparam: input: input imageparam: PSF: known PSFparam: eps: epsilon"""input_fft = np.fft.fft2(input)PSF_fft = np.fft.fft2(PSF) + eps           #噪声功率，这是已知的，考虑epsilonresult = np.fft.ifft2(input_fft / PSF_fft) #计算F(u,v)的傅里叶反变换result = np.abs(np.fft.fftshift(result))return resultdef wiener_filter(input, PSF, eps, K=0.01):"""wiener filter for image denoiseparam: input: input imageparam: PSF: input the PSF maskparam: eps: epsilonparam: K=0.01: K value for wiener fuctionreturn image after wiener filter"""input_fft = np.fft.fft2(input)PSF_fft = np.fft.fft2(PSF) + epsPSF_fft_1 = np.conj(PSF_fft) / (np.abs(PSF_fft)**2 + K)# 按公式，居然得不到正确的值
#     PSF_abs = PSF * np.conj(PSF)
#     PSF_fft_1 = (1 / (PSF + eps)) * (PSF_abs / (PSF_abs + K))result = np.fft.ifft2(input_fft * PSF_fft_1)result = np.abs(np.fft.fftshift(result))return result# 要实现的功能都在这里调用
if __name__ == "__main__":image = cv2.imread('DIP_Figures/DIP3E_Original_Images_CH05/Fig0526(a)(original_DIP).tif', 0)# 显示原图像plt.figure(1, figsize=(6, 6))plt.title("Original Image"), plt.imshow(image, 'gray')plt.xticks([]), plt.yticks([])plt.figure(2, figsize=(18, 12))# 进行运动模糊处理PSF = get_motion_dsf(image.shape[:2], -50, 100)blurred = make_blurred(image, PSF, 1e-3)plt.subplot(231), plt.imshow(blurred, 'gray'), plt.title("Motion blurred")plt.xticks([]), plt.yticks([])# 逆滤波result = inverse_filter(blurred, PSF, 1e-3)   plt.subplot(232), plt.imshow(result, 'gray'), plt.title("inverse deblurred")plt.xticks([]), plt.yticks([])# 维纳滤波result = wiener_filter(blurred, PSF, 1e-3)     plt.subplot(233), plt.imshow(result, 'gray'), plt.title("wiener deblurred(k=0.01)")plt.xticks([]), plt.yticks([])# 添加噪声,standard_normal产生随机的函数blurred_noisy = blurred + 0.1 * blurred.std() * np.random.standard_normal(blurred.shape)   # 显示添加噪声且运动模糊的图像plt.subplot(234), plt.imshow(blurred_noisy, 'gray'), plt.title("motion & noisy blurred")plt.xticks([]), plt.yticks([])# 对添加噪声的图像进行逆滤波result = inverse_filter(blurred_noisy, PSF, 0.1 + 1e-3)    plt.subplot(235), plt.imshow(result, 'gray'), plt.title("inverse deblurred")plt.xticks([]), plt.yticks([])# 对添加噪声的图像进行维纳滤波result = wiener_filter(blurred_noisy, PSF, 0.1 + 1e-3)         plt.subplot(236), plt.imshow(result, 'gray'), plt.title("wiener deblurred(k=0.01)")plt.xticks([]), plt.yticks([])plt.tight_layout()plt.show()

# 要实现的功能都在这里调用
if __name__ == "__main__":image = cv2.imread('DIP_Figures/DIP3E_Original_Images_CH05/Fig0526(a)(original_DIP).tif', 0)# 显示原图像plt.figure(1, figsize=(6, 6))plt.title("Original Image"), plt.imshow(image, 'gray')plt.xticks([]), plt.yticks([])plt.figure(2, figsize=(18, 12))# 进行运动模糊处理PSF = get_motion_dsf(image.shape[:2], -45, 95)blurred = cv2.imread('DIP_Figures/DIP3E_Original_Images_CH05/Fig0529(d)(medium_noise_var_pt01).tif', 0)plt.subplot(231), plt.imshow(blurred, 'gray'), plt.title("Motion blurred")plt.xticks([]), plt.yticks([])# 逆滤波result = inverse_filter(blurred, PSF, 1e-3)   plt.subplot(232), plt.imshow(result, 'gray'), plt.title("inverse deblurred")plt.xticks([]), plt.yticks([])# 维纳滤波result = wiener_filter(blurred, PSF, 1e-3, K=0.01)     plt.subplot(233), plt.imshow(result, 'gray'), plt.title("wiener deblurred(k=0.01)")plt.xticks([]), plt.yticks([])# 添加噪声,standard_normal产生随机的函数blurred_noisy = blurred + 0.1 * blurred.std() * np.random.standard_normal(blurred.shape)   # 显示添加噪声且运动模糊的图像plt.subplot(234), plt.imshow(blurred_noisy, 'gray'), plt.title("motion & noisy blurred")plt.xticks([]), plt.yticks([])# 对添加噪声的图像进行逆滤波result = inverse_filter(blurred_noisy, PSF, 0.1 + 1e-3)    plt.subplot(235), plt.imshow(result, 'gray'), plt.title("inverse deblurred")plt.xticks([]), plt.yticks([])# 对添加噪声的图像进行维纳滤波result = wiener_filter(blurred_noisy, PSF, 0.1 + 1e-3)         plt.subplot(236), plt.imshow(result, 'gray'), plt.title("wiener deblurred(k=0.01)")plt.xticks([]), plt.yticks([])plt.tight_layout()plt.show()

第5章 Python 数字图像处理(DIP) - 图像复原与重建15 - 最小均方误差（维纳）滤波相关推荐

第5章 Python 数字图像处理(DIP) - 图像复原与重建12 - 空间滤波 - 使用频率域滤波降低周期噪声 - 陷波滤波、最优陷波滤波
标题使用频率域滤波降低周期噪声陷波滤波深入介绍最优陷波滤波本章陷波滤波器有部分得出的结果不佳,如果有更好的解决方案,请赐教,不胜感激. 使用频率域滤波降低周期噪声陷波滤波深入介绍零相移滤波 ...
第5章 Python 数字图像处理(DIP) - 图像复原与重建1 - 高斯噪声
本章主要讲图像复原与重建,首先是了解一下各种噪声的特点与模型,还有形成的方法.一些重点的噪声,如高斯噪声,均匀噪声,伽马噪声,指数噪声,还有椒盐噪声等. 本章主要的噪声研究方法主要是加性噪声. 标题 ...
第5章 Python 数字图像处理(DIP) - 图像复原与重建17 - 由投影重建图像、雷登变换、投影、反投影、反投影重建
标题由投影重建图像投影和雷登变换 Johann Radon 反投影滤波反投影重建由投影重建图像本由投影重建图像,主要是雷登变换与雷登把变换的应用,所以也没有太多的研究,只为了保持完整性,而添 ...
第5章 Python 数字图像处理(DIP) - 图像复原与重建16 - 约束最小二乘方滤波、几何均值滤波
标题约束最小二乘方滤波几何均值滤波约束最小二乘方滤波 F^(u,v)=[H∗(u,v)∣H(u,v)∣2+γ∣P(u,v)∣2]G(u,v)(5.89)\hat{F}(u,v) = \bigg[ ...
第5章 Python 数字图像处理(DIP) - 图像复原与重建14 - 逆滤波
标题逆滤波逆滤波逆滤波逆滤波图像的退化函数已知或者由前面的方法获取退化函数,则可以直接逆滤波 F^(u,v)=G(u,v)H(u,v)(5.78)\hat{F}(u,v) = \frac{G ...
第5章 Python 数字图像处理(DIP) - 图像复原与重建13 - 空间滤波 - 线性位置不变退化 - 退化函数估计、运动模糊函数
标题线性位置不变退化估计退化函数采用观察法估计退化函数采用试验法估计退化函数采用建模法估计退化函数运动模糊函数 OpenCV Motion Blur 在这一节中,得到的结果,有些不是很好, ...
第5章 Python 数字图像处理(DIP) - 图像复原与重建11 - 空间滤波 - 自适应滤波器 - 自适应局部降噪、自适应中值滤波器
标题自适应滤波器自适应局部降噪滤波器自适应中值滤波器自适应滤波器自适应局部降噪滤波器均值是计算平均值的区域上的平均灰度,方差是该区域上的图像对比度 g(x,y)g(x, y)g(x,y)噪 ...
第5章 Python 数字图像处理(DIP) - 图像复原与重建10 - 空间滤波 - 统计排序滤波器 - 中值、最大值、最小值、中点、修正阿尔法均值滤波器
标题统计排序滤波器中值.最大值.最小值.中点滤波器修正阿尔法均值滤波器统计排序滤波器中值.最大值.最小值.中点滤波器 f^(x,y)=median{g(r,c)}(5.27)\hat{f ...
第5章 Python 数字图像处理(DIP) - 图像复原与重建9 - 空间滤波 - 均值滤波器 - 算术平均、几何平均、谐波平均、反谐波平均滤波器
标题只存在噪声的复原 - 空间滤波均值滤波器算术平均滤波器几何均值滤波器谐波平均滤波器反(逆)谐波平均滤波器只存在噪声的复原 - 空间滤波仅被加性噪声退化 g(x,y)=f(x,y)+ ...

第5章 Python 数字图像处理(DIP) - 图像复原与重建15 - 最小均方误差（维纳）滤波

标题

最小均方误差（维纳）滤波

第5章 Python 数字图像处理(DIP) - 图像复原与重建15 - 最小均方误差（维纳）滤波相关推荐

最新文章

热门文章