ME227计算机视觉：还原马赛克

Posted 2021-10-27 西兰花12138

ME227计算机视觉：还原马赛克

问题重述：

​ 课上已经讲过通过平均方法马赛克很难还原真实图像，那么假如我们采用如下函数进行图像模糊化处理：

$$f(x,y)=a{f}_0(x,y)\times sin(\frac{\pi }{8}x+b)+c{f}_0(x,y)\times cos(\frac{\pi }{8}x+d)$$
其中f0为原始像素值，f为处理后像素值，a, b, c,d 为非零常数且未知。请分析通过这种方法给图片打模糊马赛克后是否能够还原?为什么?假设模糊的像素数目为10×10。

问题分析：

​ 如果采用反解方程的形式，不妨算一算可以知道未知数太多，已知量太少，不足以反解出答案，以还原图像。

​ 由二维傅里叶变换，题中对原图像的模糊化处理可以看作是在原图像上添加一个频率为$\frac{\pi }{8}$的振动。对于二维傅里叶变换的理解可以参考知乎文章《形象理解二维傅里叶变换》，那么只需要将模糊化图像的频率域中的频率为$\frac{\pi }{8}$的分量置零即可，即频谱图距离原点为$\frac{\pi }{8}$的点置为零。

结果

原图为：

模糊化后的图为：

去噪后的图为：

代码：

import cv2
import numpy as np
import matplotlib.pyplot as plt
import math
img = cv2.imread("houguoyu.png", 0)
sp = img.shape
height = sp[0]
width =sp[1]
#按照题述要求打码，随机设置：a=1,b=2,c=3,d=4
new = np.zeros((height,width,1),np.uint8)
for i in range(height):
    for j in range(width):
        new[i,j]= img[i,j]*math.sin(math.pi/8*i+2) + 3*math.cos(math.pi/8*j+4)

#下面对原图和打码后的图进行二维傅里叶变换
f = np.fft.fft2(img)
f_n = np.fft.fft2(new)
fshift = np.fft.fftshift(f)
fshift_n = np.fft.fftshift(f_n)

magnitude_spectrum = 20 * np.log(np.abs(fshift)) #原图的幅值图
magnitude_spectrum_n = 20 * np.log(np.abs(fshift_n)) #打码后的幅值图

#去码方式：将频率为pi/8的分量去掉，即将二维傅里叶变换后得到的频域图中距离原点为9的点全部置零
qumashift = np.zeros_like(fshift)
for i in range(height):
    for j in range(width):
        if (i^2+j^2 == math.pi/8):
            qumashift[i][j] = 0
        else:
            qumashift[i][j] = fshift[i][j]


magnitude_spectrum_q = 20 * np.log(np.abs(qumashift)+1) #去码后的幅值图

#下面进行傅里叶逆变换
ishift = np.fft.ifftshift(qumashift)
iimg = np.fft.ifft2(ishift)
iimg = np.abs(iimg)

#原图
plt.subplot(121)
plt.imshow(img , cmap = 'gray')
plt.title('original')
plt.axis('off')

#原图的幅值图
plt.subplot(122)
plt.imshow(magnitude_spectrum , cmap = 'gray')
plt.title("magnitude_spectrum")
plt.axis('off')

plt.show()
#打码的图
plt.subplot(121)
plt.imshow(new , cmap = 'gray')
plt.title('noise')
plt.axis('off')

#打码的幅值图
plt.subplot(122)
plt.imshow(magnitude_spectrum_n , cmap = 'gray')
plt.title("magnitude_spectrum_n")
plt.axis('off')

plt.show()
#去码后的图
plt.subplot(121)
plt.imshow(iimg , cmap = 'gray')
plt.title('quma')
plt.axis('off')

#去码后的幅值图
plt.subplot(122)
plt.imshow(magnitude_spectrum_q , cmap = 'gray')
plt.title("magnitude_spectrum_q")
plt.axis('off')

plt.show()