【OpenCV 完整例程】87. 频率域钝化掩蔽-CFANZ编程社区

【OpenCV 完整例程】87. 频率域钝化掩蔽

4.2 频率域钝化掩蔽

简单地，从原始图像中减去一幅平滑处理的钝化图像，也可以实现图像锐化效果，称为钝化掩蔽。

令 $f_{LP}(x,y)$ 表示低通滤波的平滑图像，则：
$g_{mask} (x,y) = f(x,y) - f_{LP}(x,y) \\ g(x,y) = f(x,y) + k * g_{mask}(x,y)， k>0$
当 k>1 时，实现高提升滤波；当 k=1 时，实现钝化掩蔽；k<1时，可以减弱钝化掩蔽的强度。

原图减去模糊图的结果为模板，输出图像等于原图加上加权后的模板，当权重为1得到非锐化掩蔽，当权重大于1成为高提升滤波。

在频率域实现钝化掩蔽，高频强调滤波器传递函数为：
$g(x,y)= J^{-1} \{ [1+k H_{HP}(u,v)]F(u,v) \}$
高频强调滤波的通用公式是：
$g(x,y)= J^{-1} \{ [k_1 + k_2 H_{HP}(u,v)]F(u,v) \}$
式中， $k_1 \geq 0$ 偏移传递函数的值， $k_2 \geq 0$ 控制高频的贡献。

例程 8.26：频率域钝化掩蔽

# OpenCVdemo08.py
# Demo08 of OpenCV
# 8. 图像的频率域滤波
# Copyright 2021 Youcans, XUPT
# Crated：2021-12-30

    # 8.26：频率域钝化掩蔽
    def gaussHighPassFilter(shape, radius=10):  # 高斯高通滤波器
        # 高斯滤波器：# Gauss = 1/(2*pi*s2) * exp(-(x**2+y**2)/(2*s2))
        u, v = np.mgrid[-1:1:2.0/shape[0], -1:1:2.0/shape[1]]
        D = np.sqrt(u**2 + v**2)
        D0 = radius / shape[0]
        kernel = 1 - np.exp(- (D ** 2) / (2 *D0**2))
        return kernel

    def dft2Image(image):  #　最优扩充的快速傅立叶变换
        # 中心化, centralized 2d array f(x,y) * (-1)^(x+y)
        mask = np.ones(image.shape)
        mask[1::2, ::2] = -1
        mask[::2, 1::2] = -1
        fImage = image * mask  # f(x,y) * (-1)^(x+y)

        # 最优 DFT 扩充尺寸
        rows, cols = image.shape[:2]  # 原始图片的高度和宽度
        rPadded = cv2.getOptimalDFTSize(rows)  # 最优 DFT 扩充尺寸
        cPadded = cv2.getOptimalDFTSize(cols)  # 用于快速傅里叶变换

        # 边缘扩充(补0), 快速傅里叶变换
        dftImage = np.zeros((rPadded, cPadded, 2), np.float32)  # 对原始图像进行边缘扩充
        dftImage[:rows, :cols, 0] = fImage  # 边缘扩充，下侧和右侧补0
        cv2.dft(dftImage, dftImage, cv2.DFT_COMPLEX_OUTPUT)  # 快速傅里叶变换
        return dftImage


    # 高频强调滤波 + 直方图均衡化
    image = cv2.imread("../images/Fig0459a.tif", flags=0)  # flags=0 读取为灰度图像
    rows, cols = image.shape[:2]  # 图片的高度和宽度
    print(rows, cols)

    # 快速傅里叶变换
    dftImage = dft2Image(image)  # 快速傅里叶变换 (rPad, cPad, 2)
    rPadded, cPadded = dftImage.shape[:2]  # 快速傅里叶变换的尺寸, 原始图像尺寸优化
    # 构建 高斯高通滤波器 (Gauss low pass filter)
    hpFilter = gaussHighPassFilter((rPadded, cPadded), radius=10)  # 高斯高通滤波器
    # 在频率域修改傅里叶变换: 傅里叶变换 点乘 低通滤波器
    dftHPfilter = np.zeros(dftImage.shape, dftImage.dtype)  # 快速傅里叶变换的尺寸(优化尺寸)
    for j in range(2):
        dftHPfilter[:rPadded, :cPadded, j] = dftImage[:rPadded, :cPadded, j] * hpFilter
    # 对高通傅里叶变换 执行傅里叶逆变换，并只取实部
    idft = np.zeros(dftImage.shape[:2], np.float32)  # 快速傅里叶变换的尺寸(优化尺寸)
    cv2.dft(dftHPfilter, idft, cv2.DFT_REAL_OUTPUT + cv2.DFT_INVERSE + cv2.DFT_SCALE)
    # 中心化, centralized 2d array g(x,y) * (-1)^(x+y)
    mask2 = np.ones(dftImage.shape[:2])
    mask2[1::2, ::2] = -1
    mask2[::2, 1::2] = -1
    idftCen = idft * mask2  # g(x,y) * (-1)^(x+y)
    # 截取左上角，大小和输入图像相等
    result = np.clip(idftCen, 0, 255)  # 截断函数，将数值限制在 [0,255]
    imgHPF = result.astype(np.uint8)
    imgHPF = imgHPF[:rows, :cols]

    # # =======高频增强滤波===================
    k1 = 0.5
    k2 = 0.75
    # 在频率域修改傅里叶变换: 傅里叶变换 点乘 低通滤波器
    hpEnhance = np.zeros(dftImage.shape, dftImage.dtype)  # 快速傅里叶变换的尺寸(优化尺寸)
    for j in range(2):
        hpEnhance[:rPadded, :cPadded, j] = dftImage[:rPadded, :cPadded, j] * (k1 + k2*hpFilter)
    # 对高通傅里叶变换 执行傅里叶逆变换，并只取实部
    idft = np.zeros(dftImage.shape[:2], np.float32)  # 快速傅里叶变换的尺寸(优化尺寸)
    cv2.dft(hpEnhance, idft, cv2.DFT_REAL_OUTPUT + cv2.DFT_INVERSE + cv2.DFT_SCALE)
    # 中心化, centralized 2d array g(x,y) * (-1)^(x+y)
    mask2 = np.ones(dftImage.shape[:2])
    mask2[1::2, ::2] = -1
    mask2[::2, 1::2] = -1
    idftCen = idft * mask2  # g(x,y) * (-1)^(x+y)
    # 截取左上角，大小和输入图像相等
    result = np.clip(idftCen, 0, 255)  # 截断函数，将数值限制在 [0,255]
    imgHPE= result.astype(np.uint8)
    imgHPE = imgHPE[:rows, :cols]

    # =======直方图均衡===================
    imgEqu = cv2.equalizeHist(imgHPE)  # 使用 cv2.qualizeHist 完成直方图均衡化变换

    plt.figure(figsize=(9, 6))
    plt.subplot(221), plt.imshow(image, 'gray'), plt.title("Origin"), plt.xticks([]), plt.yticks([])
    plt.subplot(222), plt.imshow(imgHPF, 'gray'), plt.title("Gauss high-pass filter"), plt.xticks([]), plt.yticks([])
    plt.subplot(223), plt.imshow(imgHPE, 'gray'), plt.title("High frequency emphasis"), plt.xticks([]), plt.yticks([])
    plt.subplot(224), plt.imshow(imgEqu, 'gray'), plt.title("Histogram of equalized"), plt.xticks([]), plt.yticks([])
    plt.tight_layout()
    plt.show()

在这里插入图片描述

（本节完）

youcans@xupt 原创作品，转载必须标注原文链接

Crated：2022-2-1