【OpenCV】Chapter6.频率域图像滤波_opencv 频域滤波

最近想对OpenCV进行系统学习，看到网上这份教程写得不错，于是跟着来学习实践一下。
【youcans@qq.com, youcans 的 OpenCV 例程, https://youcans.blog.csdn.net/article/details/125112487】
程序仓库：https://github.com/zstar1003/OpenCV-Learning

离散傅里叶变换

OpenCV 中的cv.dft()函数也可以实现图像的傅里叶变换，cv.idft()函数实现图像傅里叶逆变换。

cv.dft(src[, dst[, flags[, nonzeroRows]]]) → dst
cv.idft(src[, dst[, flags[, nonzeroRows]]]) → dst

参数说明：

• src：输入图像，单通道灰度图像，使用 np.float32 格式
• dst：输出图像，图像大小与 src 相同，数据类型由 flag 决定
• flag：转换标识符
  • cv.DFT_INVERSE：用一维或二维逆变换取代默认的正向变换
  • cv.DFT_SCALE：缩放比例标识，根据元素数量求出缩放结果，常与DFT_INVERSE搭配使用
  • cv.DFT_ROWS: 对输入矩阵的每行进行正向或反向的傅里叶变换，常用于三维或高维变换等复杂操作
  • cv.DFT_COMPLEX_OUTPUT：对一维或二维实数数组进行正向变换，默认方法，结果是由 2个通道表示的复数阵列，第一通道是实数部分，第二通道是虚数部分
  • cv.DFT_REAL_OUTPUT：对一维或二维复数数组进行逆变换，结果通常是一个尺寸相同的复数矩阵

注：输入图像 src 是 np.float32 格式，如图像使用 np.uint8 格式则必须先转换 np.float32 格式。

使用cv.magnitude()函数可以实现计算二维矢量的幅值

cv.magnitude(x, y[, magnitude]) → dst

参数说明：

• x：一维或多维数组，也表示复数的实部，浮点型
• y：一维或多维数组，也表示复数的虚部，浮点型，数组大小必须与 x 相同

示例程序：

"""
离散傅里叶变换
"""
import cv2
import matplotlib.pyplot as plt
import numpy as np

imgGray = cv2.imread("../img/lena.jpg", flags=0)

# cv2.dft 实现图像的傅里叶变换
imgFloat32 = np.float32(imgGray)  # 将图像转换成 float32
dft = cv2.dft(imgFloat32, flags=cv2.DFT_COMPLEX_OUTPUT)  # 傅里叶变换
dftShift = np.fft.fftshift(dft)  # 将低频分量移动到频域图像的中心

# 相位谱
phase = np.arctan2(dftShift[:, :, 1], dftShift[:, :, 0])  # 计算相位角(弧度制)
dftPhi = phase / np.pi * 180  # 将相位角转换为 [-180, 180]

# cv2.idft 实现图像的逆傅里叶变换
invShift = np.fft.ifftshift(dftShift)  # 将低频逆转换回图像四角
imgIdft = cv2.idft(invShift)  # 逆傅里叶变换
imgRebuild = cv2.magnitude(imgIdft[:, :, 0], imgIdft[:, :, 1])  # 重建图像

plt.figure(figsize=(9, 6))
plt.subplot(131), plt.title("Original image"), plt.axis('off')
plt.imshow(imgGray, cmap='gray')
plt.subplot(132), plt.title("DFT Phase"), plt.axis('off')
plt.imshow(dftPhi, cmap='gray')
plt.subplot(133), plt.title("Rebuild image with IDFT"), plt.axis('off')
plt.imshow(imgRebuild, cmap='gray')
plt.tight_layout()
plt.show()
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32

快速傅里叶变换

OpenCV 中的cv.getOptimalDFTSize()函数可以实现图像的最优 DFT 尺寸扩充。

参数说明：

• versize：数组大小
• retval：DFT 扩充的最优数组大小

示例程序：

"""
快速傅里叶变换
"""
import cv2
import matplotlib.pyplot as plt
import numpy as np

imgGray = cv2.imread("../img/lena.jpg", flags=0)
rows, cols = imgGray.shape[:2]  # 图像的行(高度)/列(宽度)

# 快速傅里叶变换(要对原始图像进行矩阵扩充)
rPad = cv2.getOptimalDFTSize(rows)  # 最优 DFT 扩充尺寸
cPad = cv2.getOptimalDFTSize(cols)  # 用于快速傅里叶变换
imgEx = np.zeros((rPad, cPad, 2), np.float32)  # 对原始图像进行边缘扩充
imgEx[:rows, :cols, 0] = imgGray  # 边缘扩充，下侧和右侧补0
dftImgEx = cv2.dft(imgEx, cv2.DFT_COMPLEX_OUTPUT)  # 快速傅里叶变换

# 傅里叶逆变换
idftImg = cv2.idft(dftImgEx)  # 逆傅里叶变换
idftMag = cv2.magnitude(idftImg[:, :, 0], idftImg[:, :, 1])  # 逆傅里叶变换幅值

# 矩阵裁剪，得到恢复图像
idftMagNorm = np.uint8(cv2.normalize(idftMag, None, 0, 255, cv2.NORM_MINMAX))  # 归一化为 [0,255]
imgRebuild = np.copy(idftMagNorm[:rows, :cols])


plt.figure(figsize=(9, 6))
plt.subplot(131), plt.title("Original image"), plt.axis('off')
plt.imshow(imgGray, cmap='gray')
plt.subplot(132), plt.title("Log-trans of DFT amp"), plt.axis('off')
dftAmp = cv2.magnitude(dftImgEx[:, :, 0], dftImgEx[:, :, 1])  # 幅度谱，中心化
dftAmpLog = np.log(1 + dftAmp)  # 幅度谱对数变换，以便于显示
plt.imshow(cv2.normalize(dftAmpLog, None, 0, 255, cv2.NORM_MINMAX), cmap='gray')
plt.subplot(133), plt.title("Rebuild image with IDFT"), plt.axis('off')
plt.imshow(imgRebuild, cmap='gray')
plt.tight_layout()
plt.show()
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37

频率域低通滤波器

频率域高斯低通滤波器

频率域高斯低通滤波器是一个掩模蒙板：

$D_0=\sigma$是截止频率，$\sigma$越小，高斯函数越狭窄，滤除的高频成分（图像细节）越多，图像越模糊。

示例程序：

"""
频率域高斯低通滤波器
"""
import cv2
import matplotlib.pyplot as plt
import numpy as np

imgGray = cv2.imread("../img/lena.jpg", flags=0)

# （1）首先对图像进行傅里叶变换
imgFloat32 = np.float32(imgGray)  # 将图像转换成 float32
dft = cv2.dft(imgFloat32, flags=cv2.DFT_COMPLEX_OUTPUT)  # 傅里叶变换
dftShift = np.fft.fftshift(dft)  # 将低频分量移动到频域图像的中心

plt.figure(figsize=(9, 6))
rows, cols = imgGray.shape[:2]  # 图片的高度和宽度
sigma2 = [0.5, 0.09, 0.01]  # square of sigma
for i in range(3):
    # 构造高斯滤波器遮罩
    x, y = np.mgrid[-1:1:2.0 / rows, -1:1:2.0 / cols]
    z = 1 / (2 * np.pi * sigma2[i]) * np.exp(-(x ** 2 + y ** 2) / (2 * sigma2[i]))
    zNorm = np.uint8(cv2.normalize(z, None, 0, 255, cv2.NORM_MINMAX))  # 归一化为 [0,255]
    maskGauss = np.zeros((rows, cols, 2), np.uint8)
    maskGauss[:, :, 0] = zNorm
    maskGauss[:, :, 1] = zNorm
    # （2）然后在频率域修改傅里叶变换
    dftTrans = dftShift * maskGauss  # 修改傅里叶变换实现滤波
    # （3）最后通过傅里叶逆变换返回空间域
    ishift = np.fft.ifftshift(dftTrans)  # 将低频逆转换回图像四角
    idft = cv2.idft(ishift)  # 逆傅里叶变换
    imgRebuild = cv2.magnitude(idft[:, :, 0], idft[:, :, 1])  # 重建图像

    plt.subplot(2, 3, i + 1), plt.title("Mask (s^2={})".format(sigma2[i])), plt.axis('off')
    plt.imshow(maskGauss[:, :, 0], cmap='gray')
    plt.subplot(2, 3, i + 4), plt.title("DFT GLPF (s^2={})".format(sigma2[i])), plt.axis('off')
    plt.imshow(imgRebuild, cmap='gray')

plt.tight_layout()
plt.show()
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39

图中s^2为${\sigma }^2$

频率域巴特沃斯低通滤波器

n阶巴特沃斯(Butterworth)低通滤波器的传递函数为：

当n较大时，巴特沃斯低通滤波器BLPF可以逼近理想低通滤波器ILPF的特性；而当n较小时，巴特沃斯低通滤波器 BLPF 可以逼近高斯低通滤波器 GLPF 的特性，同时提供从低频到高频的平滑过渡。

示例程序：

"""
频率域巴特沃斯低通滤波器
"""
import cv2
import matplotlib.pyplot as plt
import numpy as np

imgGray = cv2.imread("../img/lena.jpg", flags=0)

imgFloat32 = np.float32(imgGray)  # 将图像转换成 float32
rows, cols = imgGray.shape[:2]  # 图片的高度和宽度

# (2) 中心化, centralized 2d array f(x,y) * (-1)^(x+y)
mask = np.ones(imgGray.shape)
mask[1::2, ::2] = -1
mask[::2, 1::2] = -1
fImage = imgFloat32 * mask  # f(x,y) * (-1)^(x+y)

# (3) 快速傅里叶变换
# dftImage = fft2Image(fImage)  # 快速傅里叶变换 (rPad, cPad, 2)
rPadded = cv2.getOptimalDFTSize(rows)  # 最优 DFT 扩充尺寸
cPadded = cv2.getOptimalDFTSize(cols)  # 用于快速傅里叶变换
dftImage = np.zeros((rPadded, cPadded, 2), np.float32)  # 对原始图像进行边缘扩充
dftImage[:rows, :cols, 0] = fImage  # 边缘扩充，下侧和右侧补0
cv2.dft(dftImage, dftImage, cv2.DFT_COMPLEX_OUTPUT)  # 快速傅里叶变换

dftAmp = cv2.magnitude(dftImage[:, :, 0], dftImage[:, :, 1])  # 傅里叶变换的幅度谱 (rPad, cPad)
dftAmpLog = np.log(1.0 + dftAmp)  # 幅度谱对数变换，以便于显示
dftAmpNorm = np.uint8(cv2.normalize(dftAmpLog, None, 0, 255, cv2.NORM_MINMAX))  # 归一化为 [0,255]
minValue, maxValue, minLoc, maxLoc = cv2.minMaxLoc(dftAmp)  # 找到傅里叶谱最大值的位置

plt.figure(figsize=(9, 6))
# rows, cols = imgGray.shape[:2]  # 图片的高度和宽度
u, v = np.mgrid[0:rPadded:1, 0:cPadded:1]
D = np.sqrt(np.power((u - maxLoc[1]), 2) + np.power((v - maxLoc[0]), 2))
D0 = [20, 40, 80]  # cut-off frequency
n = 2
for k in range(3):
    # (4) 构建低通滤波器 传递函数
    # 巴特沃斯低通滤波 (Butterworth low pass filter)
    epsilon = 1e-8  # 防止被 0 除
    lpFilter = 1.0 / (1.0 + np.power(D / (D0[k] + epsilon), 2 * n))

    # (5) 在频率域修改傅里叶变换: 傅里叶变换 点乘 低通滤波器
    dftLPfilter = np.zeros(dftImage.shape, dftImage.dtype)  # 快速傅里叶变换的尺寸(优化尺寸)
    for j in range(2):
        dftLPfilter[:rPadded, :cPadded, j] = dftImage[:rPadded, :cPadded, j] * lpFilter

    # (6) 对低通傅里叶变换 执行傅里叶逆变换，并只取实部
    idft = np.zeros(dftAmp.shape, np.float32)  # 快速傅里叶变换的尺寸(优化尺寸)
    cv2.dft(dftLPfilter, idft, cv2.DFT_REAL_OUTPUT + cv2.DFT_INVERSE + cv2.DFT_SCALE)

    # (7) 中心化, centralized 2d array g(x,y) * (-1)^(x+y)
    mask2 = np.ones(dftAmp.shape)
    mask2[1::2, ::2] = -1
    mask2[::2, 1::2] = -1
    idftCen = idft * mask2  # g(x,y) * (-1)^(x+y)

    # (8) 截取左上角，大小和输入图像相等
    result = np.clip(idftCen, 0, 255)  # 截断函数，将数值限制在 [0,255]
    imgBLPF = result.astype(np.uint8)
    imgBLPF = imgBLPF[:rows, :cols]

    plt.subplot(2, 3, k + 1), plt.title("BLPF mask(D0={})".format(D0[k])), plt.axis('off')
    plt.imshow(lpFilter[:, :], cmap='gray')
    plt.subplot(2, 3, k + 4), plt.title("BLPF rebuild(D0={})".format(D0[k])), plt.axis('off')
    plt.imshow(imgBLPF, cmap='gray')

plt.tight_layout()
plt.show()
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70

频率域高通滤波器

传递函数可视化

在频率域中用 1 减去低通滤波器的传递函数，就可以得到相应的高通滤波器传递函数：

理想高通滤波器（IHPF）的传递函数为：

高斯高通滤波器（GHPF）的传递函数为：

巴特沃斯高通滤波器（BHPF）的传递函数为：

示例程序：

"""
频率域高通滤波器
"""
import cv2
import matplotlib.pyplot as plt
import numpy as np


def ideaHighPassFilter(shape, radius=10):  # 理想高通滤波器
    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 = np.ones(shape)
    kernel[D <= D0] = 0  # 理想低通滤波 (Idea low pass filter)
    return kernel


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 butterworthHighPassFilter(shape, radius=10, n=2):  # 巴特沃斯高通滤波
    u, v = np.mgrid[-1:1:2.0 / shape[0], -1:1:2.0 / shape[1]]
    epsilon = 1e-8
    D = np.sqrt(u ** 2 + v ** 2)
    D0 = radius / shape[0]
    kernel = 1.0 / (1.0 + np.power(D0 / (D + epsilon), 2 * n))

    return kernel


# 理想、高斯、巴特沃斯高通传递函数
shape = [128, 128]
radius = 32
IHPF = ideaHighPassFilter(shape, radius=radius)
GHPF = gaussHighPassFilter(shape, radius=radius)
BHPF = butterworthHighPassFilter(shape, radius=radius)

filters = ['IHPF', 'GHPF', 'BHPF']
u, v = np.mgrid[-1:1:2.0 / shape[0], -1:1:2.0 / shape[1]]
fig = plt.figure(figsize=(10, 8))
for i in range(3):
    hpFilter = eval(filters[i]).copy()

    ax1 = fig.add_subplot(3, 3, 3 * i + 1)
    ax1.imshow(hpFilter, 'gray')
    ax1.set_title(filters[i]), ax1.set_xticks([]), ax1.set_yticks([])

    ax2 = plt.subplot(3, 3, 3 * i + 2, projection='3d')
    ax2.set_title("transfer function")
    ax2.plot_wireframe(u, v, hpFilter, rstride=2, linewidth=0.5, color='c')
    ax2.set_xticks([]), ax2.set_yticks([]), ax2.set_zticks([])

    ax3 = plt.subplot(3, 3, 3 * i + 3)
    profile = hpFilter[shape[0] // 2:, shape[1] // 2]
    ax3.plot(profile), ax3.set_title("profile"), ax3.set_xticks([]), ax3.set_yticks([])

plt.show()
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63

应用二维图像

示例程序：

"""
频率域高通滤波器-图片
"""
import cv2
import matplotlib.pyplot as plt
import numpy as np


def ideaHighPassFilter(shape, radius=10):  # 理想高通滤波器
    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 = np.ones(shape)
    kernel[D <= D0] = 0  # 理想低通滤波 (Idea low pass filter)
    return kernel


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 butterworthHighPassFilter(shape, radius=10, n=2):  # 巴特沃斯高通滤波
    u, v = np.mgrid[-1:1:2.0 / shape[0], -1:1:2.0 / shape[1]]
    epsilon = 1e-8
    D = np.sqrt(u ** 2 + v ** 2)
    D0 = radius / shape[0]
    kernel = 1.0 / (1.0 + np.power(D0 / (D + epsilon), 2 * n))

    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


# (1) 读取原始图像
imgGray = cv2.imread("../img/lena.jpg", flags=0)
rows, cols = imgGray.shape[:2]  # 图片的高度和宽度

plt.figure(figsize=(10, 6))
plt.subplot(2, 3, 1), plt.title("Original"), plt.axis('off'), plt.imshow(imgGray, cmap='gray')

# (2) 快速傅里叶变换
dftImage = dft2Image(imgGray)  # 快速傅里叶变换 (rPad, cPad, 2)
rPadded, cPadded = dftImage.shape[:2]  # 快速傅里叶变换的尺寸, 原始图像尺寸优化

D0 = [20, 40, 80, 120, 160]  # radius
for k in range(5):
    # (3) 构建高通滤波器
    # hpFilter = ideaHighPassFilter((rPadded, cPadded), radius=D0[k])  # 理想高通滤波器
    # hpFilter = gaussHighPassFilter((rPadded, cPadded), radius=D0[k])  # 高斯高通滤波器
    hpFilter = butterworthHighPassFilter((rPadded, cPadded), radius=D0[k])  # 巴特沃斯高通滤波器

    # (5) 在频率域修改傅里叶变换: 傅里叶变换 点乘 低通滤波器
    dftHPfilter = np.zeros(dftImage.shape, dftImage.dtype)  # 快速傅里叶变换的尺寸(优化尺寸)
    for j in range(2):
        dftHPfilter[:rPadded, :cPadded, j] = dftImage[:rPadded, :cPadded, j] * hpFilter

    # (6) 对高通傅里叶变换 执行傅里叶逆变换，并只取实部
    idft = np.zeros(dftImage.shape[:2], np.float32)  # 快速傅里叶变换的尺寸(优化尺寸)
    cv2.dft(dftHPfilter, 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] = -1
    mask2[::2, 1::2] = -1
    idftCen = idft * mask2  # g(x,y) * (-1)^(x+y)

    # (8) 截取左上角，大小和输入图像相等
    result = np.clip(idftCen, 0, 255)  # 截断函数，将数值限制在 [0,255]
    imgHPF = result.astype(np.uint8)
    imgHPF = imgHPF[:rows, :cols]

    plt.subplot(2, 3, k + 2), plt.title("HPFilter rebuild(n={})".format(D0[k])), plt.axis('off')
    plt.imshow(imgHPF, cmap='gray')

plt.tight_layout()
plt.show()
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98

注：在第三步可选择滤波器种类

频率域图像锐化-Laplacian算子

拉普拉斯算子(Laplace)是导数算子，会突出图像中的急剧灰度变化，抑制灰度缓慢变化区域，往往会产生暗色背景下的灰色边缘和不连续图像。将拉普拉斯图像与原图叠加，可以得到保留锐化效果的图像。

下面的示例程序将空间域拉普拉斯算子锐化和频率域拉普拉斯算子进行比较

"""
频率域图像锐化
"""
import cv2
import matplotlib.pyplot as plt
import numpy as np


def LaplacianFilter(shape):  # 频域 Laplacian 滤波器
    u, v = np.mgrid[-1:1:2.0 / shape[0], -1:1:2.0 / shape[1]]
    D = np.sqrt(u ** 2 + v ** 2)
    kernel = -4 * np.pi ** 2 * D ** 2
    return kernel


def imgHPfilter(image, lpTyper="Laplacian"):  # 频域高通滤波
    # (1) 中心化, 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)

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

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

    # (4) 构建 频域滤波器传递函数: 以 Laplacian 为例
    LapFilter = LaplacianFilter((rPadded, cPadded))  # 拉普拉斯滤波器

    # (5) 在频率域修改傅里叶变换: 傅里叶变换 点乘 滤波器传递函数
    dftFilter = np.zeros(dftImage.shape, dftImage.dtype)  # 快速傅里叶变换的尺寸(优化尺寸)
    for j in range(2):
        dftFilter[:rPadded, :cPadded, j] = dftImage[:rPadded, :cPadded, j] * LapFilter

    # (6) 对高通傅里叶变换 执行傅里叶逆变换，并只取实部
    idft = np.zeros(dftImage.shape[:2], 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] = -1
    mask2[::2, 1::2] = -1
    idftCen = idft * mask2  # g(x,y) * (-1)^(x+y)

    # (8) 截取左上角，大小和输入图像相等
    result = np.clip(idftCen, 0, 255)  # 截断函数，将数值限制在 [0,255]
    imgFilter = result.astype(np.uint8)
    imgFilter = imgFilter[:rows, :cols]
    return imgFilter


# (1) 读取原始图像
img = cv2.imread("../img/lena.jpg", flags=0)
rows, cols = img.shape[:2]  # 图片的高度和宽度

# (2) 空间域 拉普拉斯算子 (Laplacian)s
# 使用 cv2.Laplacian 实现 Laplace 卷积算子
imgLaplace2 = cv2.Laplacian(img, -1, ksize=3)
imgLapReSpace = cv2.add(img, imgLaplace2)  # 恢复原图像

# (3) 频率域 拉普拉斯算子 (Laplacian)
imgLaplace = imgHPfilter(img, "Laplacian")  # 调用自定义函数 imgHPfilter()
imgLapRe = cv2.add(img, imgLaplace)  # 恢复原图像

plt.figure(figsize=(10, 6))
plt.subplot(131), plt.imshow(img, 'gray'), plt.title("Origin from NASA"), plt.xticks([]), plt.yticks([])
plt.subplot(132), plt.imshow(imgLapReSpace, 'gray'), plt.title("Spatial Lapalacian"), plt.xticks([]), plt.yticks([])
plt.subplot(133), plt.imshow(imgLapRe, 'gray'), plt.title("Freauency Lapalacian"), plt.xticks([]), plt.yticks([])
plt.tight_layout()
plt.show()
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76