【youcans 的 OpenCV 例程300 篇】105. 湍流模糊退化模型

5.3 模型法估计退化函数

估计图像复原中所用的退化函数，主要有三种方法：观察法、试验法和数学建模方法。

分析导致退化的原因，根据基本原理提出退化模型，如湍流导致的模糊、匀速运动导致的模糊，可以基于模型更加准确地估计退化函数。

下面以运动模糊和大气湍流模型为例，采用退化模型对图像的退化建模。

例程 9.19：湍流模糊退化模型

湍流是自然界中普遍存在的一种复杂的流动现象。物体通过湍流大气成像时，受到湍流效应的影响，出现光强闪烁、光束方向漂移、光束宽度扩展及接收面上相位的起伏，造成图像模糊和抖动，甚至扭曲变形。

Hufnagel and Stanley 根据大气湍流的物理特性提出一种退化模型：

$H(u,v) = e^{-k (u^2+v^2)^{5/6}}$

k 是湍流常数，反映湍流的强烈程度。

    # 9.19: 湍流模糊退化模型 (turbulence blur degradation model)def getDegradedImg(image, Huv):  # 根据退化模型生成退化图像rows, cols = image.shape[:2]  # 图片的高度和宽度# (1) 中心化, centralized 2d array f(x,y) * (-1)^(x+y)mask = np.ones((rows, cols))mask[1::2, ::2] = -1mask[::2, 1::2] = -1imageCen = image * mask# (2) 快速傅里叶变换dftImage = np.zeros((rows, cols, 2), np.float32)dftImage[:, :, 0] = imageCencv2.dft(dftImage, dftImage, cv2.DFT_COMPLEX_OUTPUT)  # 快速傅里叶变换 (rows, cols, 2)# (4) 构建 频域滤波器传递函数:Filter = np.zeros((rows, cols, 2), np.float32)  # (rows, cols, 2)Filter[:, :, 0], Filter[:, :, 1] = Huv, Huv# (5) 在频率域修改傅里叶变换: 傅里叶变换 点乘 滤波器传递函数dftFilter = dftImage * Filter# (6) 对修正傅里叶变换 进行傅里叶逆变换，并只取实部idft = np.ones((rows, cols), np.float32)  # 快速傅里叶变换的尺寸cv2.dft(dftFilter, idft, cv2.DFT_REAL_OUTPUT + cv2.DFT_INVERSE + cv2.DFT_SCALE)  # 只取实部# (7) 中心化, centralized 2d array g(x,y) * (-1)^(x+y)mask2 = np.ones(dftImage.shape[:2])mask2[1::2, ::2] = -1mask2[::2, 1::2] = -1idftCen = idft * mask2  # g(x,y) * (-1)^(x+y)# (8) 截取左上角，大小和输入图像相等imgDegraded = np.uint8(cv2.normalize(idftCen, None, 0, 255, cv2.NORM_MINMAX))  # 归一化为 [0,255]# print(image.shape, dftFilter.shape, imgDegraded.shape)return imgDegradeddef turbulenceBlur(img, k=0.001):  # 湍流模糊传递函数# H(u,v) = exp(-k(u^2+v^2)^5/6)M, N = img.shape[1], img.shape[0]u, v = np.meshgrid(np.arange(M), np.arange(N))radius = (u - M//2)**2 + (v - N//2)**2kernel = np.exp(-k * np.power(radius, 5/6))return kernel# 读取原始图像img = cv2.imread("../images/Fig0525a.tif", 0)  # flags=0 读取为灰度图像# 生成湍流模糊图像HBlur1 = turbulenceBlur(img, k=0.001)  # 湍流模糊传递函数imgBlur1 = getDegradedImg(img, HBlur1)  # 生成湍流模糊图像HBlur2 = turbulenceBlur(img, k=0.0025)imgBlur2 = getDegradedImg(img, HBlur2)plt.figure(figsize=(9, 6))plt.subplot(131), plt.title("origin"), plt.axis('off'), plt.imshow(img, 'gray')plt.subplot(132), plt.title("turbulence blur(k=0.001)"), plt.axis('off'), plt.imshow(imgBlur1, 'gray')plt.subplot(133), plt.title("turbulence blur(k=0.0025)"), plt.axis('off'), plt.imshow(imgBlur2, 'gray')plt.tight_layout()plt.show()