频率域滤波学习笔记1——频率域滤波基本操作

1.基础知识

  在频率域中，图像的低频分量对应着图像中较为缓慢变化的部分，而高频分量则对应着图像中较为快速变化的部分。举例来说，对于一幅房间图像，墙壁和地板通常具有较为平滑的灰度变化，这些平滑的变化可以被认为是低频分量。而在图像中的物体边缘、纹理或其他细节部分，则表现出较快速的灰度变化，这些快速变化可以被认为是高频分量。因此，通过傅里叶变换，我们可以将图像分解为不同频率的分量，从而更好地理解图像中的灰度变化模式，并针对不同频率的分量进行处理，以实现各种图像处理任务，如平滑、锐化、边缘检测等。

  频率域中的滤波过程如下：

  中心化操作是将图像的低频分量（图像的平均亮度或颜色渐变）位于频谱的中心，而高频分量（图像的边缘和细节）分布在边缘，这使得频率分布更加直观和易于分析。

  滤波过程用公式表达为：

$$g(x,y)=Real\left\{{\mathrm{\Im }}^{-1}[H(u,v)F(u,v)]\right\}$$

  其中${\mathrm{\Im }}^{-1}$是IDFT(傅里叶逆变换), F(u,v)是输出图像f(x,y)的DFT, H(u,v)是滤波器传递函数，g(x,y)是滤波后的图像。函数F,H和g是大小为PxQ的阵列，即经过填充的输入图像的大小。

  变换中的低频与图像中缓慢变化的灰度分量有关，如房间的墙壁和室外少云的天空等。另一方面，高频由灰度的急剧过渡造成，如边缘和噪声。因此，我们可以认为衰减高频而通过低频的函数H(u,v)(称为低通滤波器)将图像模糊；具有相反性质的滤波器（称为高通滤波器）将使图像的细节更加清晰，但会降低图像的对比度。

2.低通滤波

2.1 理想低通滤波器（ILPF）

  以原点为中心的一个圆无衰减的通过所有频率，而圆外频率被截止

$$H(u,v)=\{\begin{array}{l}1,D(u,v)\le D_0\\ 0,D(u,v)>D_0\end{array}$$

式中，$D_0$是一个正常数（即圆的半径），$D(u,v)$是频率中点$（u,v）$到$P\times Q$频率域中心的距离，即：

$$D(u,v)=[(u-P/2{)}^2+(v-Q/2{)}^2{]}^{1/2}$$

注意$P$与$Q$都是填充零后的大小。

示意图：

2.2 高斯低通滤波器（GLPF）

$$H(u,v)=e^{-D(u,v{)}^2/{2D_0}^2}$$

示意图：

2.3 巴特沃斯低通滤波器（BLPF）

$$H(u,v)=\frac{1}{1+[D(u,v)/D_0{]}^{2n}}$$

示意图:

3.高通滤波

高通滤波可以由低通滤波简单变换得到，如下表所示：

理想高通滤波	高斯高通滤波	巴特沃斯高通滤波
$H(u,v)=\{\begin{array}{l}0,D(u,v)\le D_0\\ 1,D(u,v)>D_0\end{array}$	$H(u,v)=1-e^{-D(u,v{)}^2/{2D_0}^2}$	$H(u,v)=\frac{1}{1+[D_0/D(u,v){]}^{2n}}$

示意图：

4.实例

下面我们对一张飞机照片进行频率域的变换来体会频率域滤波的特点

效果图：

  可以看到，上面三幅图对图片进行了低通滤波，截止频率$D_0$设置为20,下面三幅图对图片进行了高通滤波，截止频率设置为2。我们可以发现高斯滤波与巴特沃斯滤波的效果几乎相同，这是由于这里我们对巴特沃斯滤波的阶数$n$设置为1所造成的，较高的n值可以使BLPF函数逼近ILPF的特性，而使用较低的n值可以使BLPF逼近GLPF函数的特性。
  其次还需要注意的是使用理想低通滤波器会产生振铃效应（由于理想滤波器在频率域中的截止边界是突变的，即从完全通过到完全阻止的过渡是突然的），振铃效应会在滤波后的图像中产生明显的振铃或震荡现象，表现为图像中的明暗交替区域出现周期性的光暗环或边缘，这是由于频率域中的截止边界导致的频率分量的泄漏。这种振铃效应会使得图像的细节部分产生不自然的震荡，降低图像的质量。而观察GLPF的示意图可以其截止边界较为平滑，因此不会产生振铃效应。
  ILPF振铃效应示意图如下：

 使用BLPF可以用截止频率来控制低频与高频的过渡情况（比如医学图像），其在低阶时只有较小的振铃效应和较小的负值，但远没有ILPF那样的明显，但是高阶振铃效应非常明显。

 BLPF振铃效应示意图如下：

5.总结

 在实际应用中，理想低通滤波器常用于频域图像平滑，高斯低通滤波器用于频域图像平滑并缓解振铃效应，巴特沃斯低通滤波器可调节截止频率及斜率以更精确地控制频域图像平滑效果。我们需要根据具体的问题来选择方法。

6.代码

 对三种滤波器进行基本操作的代码如下：

import numpy as np
import matplotlib.pyplot as plt
import cv2

def make_transform_fun(image,d,n,T):
    """
    :param image: 原始图像
    :param d:  截止频率
    :param n:  阶数
    :param T:  高通、低通或高斯
    :return:   传递函数
    """
    # 获取图像尺寸
    rows, cols = image.shape

    # 生成中心在原点的 x 和 y 坐标网格
    x = np.arange(cols) - cols // 2
    y = np.arange(rows) - rows // 2
    x, y = np.meshgrid(x, y)

    # 计算每个点到中心点的距离 D(u,v)
    # 创建两个网格 x 和 y 来实现，它们分别包含了每个点相对于中心点的横纵坐标值
    D = np.sqrt(x ** 2 + y ** 2)

    # 根据滤波器类型应用不同的转换,利用广播机制
    if T == 'L':  # 低通滤波器
        transform_fun = 1 / (1 + (D / d) ** (2 * n))
    elif T == 'H':  # 高通滤波器
        transform_fun = 1 / (1 + (d / (D + 1e-9)) ** (2 * n))  # 添加小值避免除以零
    elif T == 'GL':  # 高斯低通
        transform_fun = np.exp(-(D ** 2) / (2 * d ** 2))
    elif T == 'GH':  # 高斯高通
        transform_fun = 1 - np.exp(-(D ** 2) / (2 * d ** 2))

    return transform_fun

def IdealLowPassFilter(image):
    """
    理想低通滤波器
    :param img: 频域图
    :return: 幅度图像
    """
    f = np.fft.fft2(image)
    fshift = np.fft.fftshift(f)

    rows,cols = image.shape
    crow,ccol = int(rows/2),int(cols/2)     # 中心位置
    mask = np.zeros((rows,cols),np.uint8)     # 2个通道，分别存储实部和虚部
    mask[crow - 20:crow + 20,ccol - 20:ccol + 20] = 1

    # 掩膜图像和频谱图像乘积
    f = fshift * mask

    # 傅里叶逆变换
    ishift = np.fft.ifftshift(f)
    iimg = np.fft.ifft2(ishift)
    iimg = np.abs(iimg)     # 返回幅度值(将逆变换后的复数图像转换为幅度图像)

    return iimg

def IdealHighPassFilter(image):
    """
    理想高通滤波器
    :param image: 原图
    :return: 幅度图像
    """
    f = np.fft.fft2(image)
    fshift = np.fft.fftshift(f)

    rows,cols = image.shape
    crow,ccol = int(rows / 2),int(cols / 2)
    fshift[crow - 2:crow + 2,ccol - 2:ccol + 2] = 0

    # 傅里叶逆变换
    ishift = np.fft.ifftshift(fshift)
    iimg = np.fft.ifft2(ishift)
    iimg = np.abs(iimg)  # 返回幅度值(将逆变换后的复数图像转换为幅度图像)

    return iimg

def GaussianPassfilter(image,d,T):
    """
    高斯滤波器
    :param image: 原图像
    :param d:  截止频率
    """
    f = np.fft.fft2(image)
    fshift = np.fft.fftshift(f)
    T = 'G' + T
    d_fun = make_transform_fun(image,d,1,T)
    new_img = np.abs(np.fft.ifft2(np.fft.ifftshift(fshift * d_fun)))
    return new_img

def ButterworthPassfilter(image,d,n,T):
    """
    巴特沃尔滤波器
    :param image: 原图像
    :param d: 截止频率
    :param n: 阶数
    :param T: 高通或低通
    :return: 幅度图
    """
    f = np.fft.fft2(image)
    fshift = np.fft.fftshift(f)

    d_fun = make_transform_fun(image,d,n,T)
    new_img = np.abs(np.fft.ifft2(np.fft.ifftshift(fshift * d_fun)))
    return new_img

if __name__ == "__main__":
    image = cv2.imread("Airport.jpg",cv2.COLOR_BGR2GRAY)

    # 这里设置低通半径为20，高通为2
    result1 = IdealLowPassFilter(image)
    result2 = IdealHighPassFilter(image)
    result3 = GaussianPassfilter(image,20,'L')     # D0方差为20，越小波形越尖
    result4 = GaussianPassfilter(image,2,'H')
    result5 = ButterworthPassfilter(image,20,1,'L')   # n为1类似于高斯,D0方差为20
    result6 = ButterworthPassfilter(image,2,1,'H')

    plt.figure(figsize=(10,5))
    plt.subplot(231), plt.axis('off'), plt.title('ILPF Image D20'), plt.imshow(result1,cmap='gray')
    plt.subplot(234), plt.axis('off'), plt.title('IHPF Image D2'), plt.imshow(result2, cmap='gray')
    plt.subplot(232), plt.axis('off'), plt.title('GLPF Image D20'), plt.imshow(result3, cmap='gray')
    plt.subplot(235), plt.axis('off'), plt.title('GHPF Image D2'), plt.imshow(result4, cmap='gray')
    plt.subplot(233), plt.axis('off'), plt.title('BLPF Image D20'), plt.imshow(result5, cmap='gray')
    plt.subplot(236), plt.axis('off'), plt.title('BLPF Image D2'), plt.imshow(result6, cmap='gray')
    plt.savefig("result.jpg")
    plt.show()

7.参考

1.https://blog.csdn.net/qq_38463737/article/details/118682500
2.数字图像处理第四版（冈萨雷斯）
3.https://zh.wikipedia.org/wiki/%E4%BD%8E%E9%80%9A%E6%BB%A4%E6%B3%A2%E5%99%A8