python 二维FFT

二维FFT常用在图像处理上，首先要能理解二维FFT的意义，否则很难明白它到底是怎么工作的。

1. 高频与低频

import numpy as np
def demo_fft():
    fig = plt.figure()
    from skimage import io
    image = io.imread(r'C:\Temp/sample2.png')
    image = np.mean(image, axis=2)
    M, N = image.shape
    print(image.shape)
    # show image
    ax = fig.add_subplot(2, 3, 1)
    ax.imshow(image, cmap='gray')

    F = np.fft.fftn(image)
    F_magnitude = np.abs(F)
    F_magnitude = np.fft.fftshift(F_magnitude)
    ax = fig.add_subplot(2, 3, 4)
    ax.imshow(np.log(1 + F_magnitude), cmap='gray', extent=(-N // 2, N // 2, -M // 2, M // 2))

    # set the low frequency section to 0
    K = 40
    F_shift = np.fft.fftshift(F)
    F_shift[M // 2 - K: M // 2 + K, N // 2 - K: N // 2 + K] = 0
    ax = fig.add_subplot(2, 3, 5)
    ax.imshow(np.log(1 + np.abs(F_shift)), cmap='gray', extent=(-N // 2, N // 2, -M // 2, M // 2))
    F_shift = np.fft.ifftshift(F_shift)

    image_filtered = np.real(np.fft.ifft2(F_shift))
    ax = fig.add_subplot(2, 3, 2)
    ax.imshow(image_filtered, cmap='gray')

    F_copy = F.copy()
    zero_array= np.zeros((M, N))
    K = 40
    zero_array[M // 2 - K: M // 2 + K, N // 2 - K: N // 2 + K] = 1
    F_shift = np.fft.fftshift(F_copy)
    F_shift = F_shift * zero_array
    ax = fig.add_subplot(2, 3, 6)
    ax.imshow(np.log(1 + np.abs(F_shift)), cmap='gray', extent=(-N // 2, N // 2, -M // 2, M // 2))
    image_filtered = np.real(np.fft.ifft2(F_shift))
    ax = fig.add_subplot(2, 3, 3)
    ax.imshow(image_filtered, cmap='gray')

    plt.show()
demo_fft()

点击查看原图

第一列是原图和对应的频率信息，第二列是去除低频部分后，FFT逆变换得到的图像。第三列是去除高频部分后FFT逆变换得到的图像。
从第二列可以看出高频贡献了图像的细节。从白到黑的边界保留了下来。而原图中大片的白与大片的黑在这个图中没什么区别。
第三列中保留了原图中的亮部与灰部，而由黑到白的临界线却很模糊。细小的白线黑线也没能显示。所以低频贡献了图像的明暗。

2.工作原理理解

import numpy as np
def demo_fft():
    fig = plt.figure()
    from skimage import io
    for i in range(1,4):
        image = io.imread(f'C:/Users/UC212310/Desktop/Temp/image{i}.png')
        image = np.mean(image, axis=2)
        M, N = image.shape
        # show image
        ax = fig.add_subplot(2, 3, i)
        ax.imshow(image, cmap='gray')

        F = np.fft.fft2(image)
        F_magnitude = np.abs(F)
        F_magnitude = np.fft.fftshift(F_magnitude)
        ax = fig.add_subplot(2, 3, i+3)
        ax.imshow(np.log(1 + F_magnitude), cmap='gray', extent=(-N // 2, N // 2, -M // 2, M // 2))
    plt.show()
demo_fft()

Figure_1.png

二维FFT就是先对行做次一维FFT，这样每个元素都是关于行频率信息了，然后再对列做一维FFT,这样每个元素都包含了行和列的频率信息。每个元素都是个复数，取绝对值可得到振幅，从实部与虚部的比值可等到相位，在二维矩阵的位置信息包含了频率大小和方向。方向在一维FFT中是不用考虑的。
FFT2的结果也是正频率从0到高然后负频率从高到0.fftshift()之后会将低频放到中间位置。
第一幅图的频谱是中间一条白线，也就是说许多个正弦波沿横向传播。纵向上没有变化。
第三幅图的频谱是十字形加一条从左下角到右上角的直线。说明原图在横向，纵向都有变化，变化的方向从左下角到右上角。
从中心到频谱图上某一点构成的向量方向就是这个波传播的方向。
正负对称才能消除虚部，这点与一维FFT原理一致。