人工智能-数学基础-傅里叶变换与卷积

秒速五厘米 2023-01-14 11:59 136阅读 0赞

## 一.周期函数 ##

**周期函数可以进行傅里叶展开:**  
 f ( x ) = a 0 + ∑ ( a n sin ⁡ ( 2 π n x ) + b n cos ⁡ ( 2 π n x ) ) f(x) = a\_0 + \\sum(a\_n\\sin(2\\pi nx) + b\_n\\cos(2\\pi nx)) f(x)=a0\+∑(ansin(2πnx)\+bncos(2πnx))

""" 1.周期函数 """
    import numpy as np
    import matplotlib.pyplot as plt
    
    plt.switch_backend("TkAgg")
    
    x = np.linspace(0, 2 * np.pi, 10000)
    y = np.sin(x * 2) + np.sin(100 * x) * 0.3
    plt.plot(x, y)
    plt.show()

![在这里插入图片描述][watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3h5eTEwMjg_size_16_color_FFFFFF_t_70_pic_center]  
**图中函数包含有高频信息和低频信息, 整个大的周期曲线称为低频, 小的周期曲线称为高频**

## 二.傅里叶变换 ##

### 2.1 连续数据傅里叶变换 ###

**傅⾥叶变换假设被变换的函数（或称为信号） x(t)是连续的，并且是非周期的，此时定义⼀个变换：**  
 X ( ω ) = F ( x ( t ) ) = ∫ − ∞ + ∞ x ( t ) e − i ω t d t X(\\omega) = F(x(t)) = \\int\_\{-\\infty\}^\{+\\infty\}x(t)e^\{-i \\omega t\}dt X(ω)=F(x(t))=∫−∞\+∞x(t)e−iωtdt  
 此 时 称 X ( ω ) 为 频 谱 ， 此 公 式 将 函 数 x ( t ) 由 时 间 域 （ 时 间 t 为 ⾃ 变 量 ） 变 换 到 了 频 率 域 （ 以 ω 为 ⾃ 变 量 ） 此时称X(\\omega) 为频谱，此公式将函数x(t) 由时间域（时间t为⾃变量）变换到了频率域（以\\omega为⾃变量） 此时称X(ω)为频谱，此公式将函数x(t)由时间域（时间t为⾃变量）变换到了频率域（以ω为⾃变量）

**2.1.1** **下面看下对函数 f(x) = sin(2x) + 0.3sin(100x)+ 0.6sin(200x)进行傅里叶变换:**  
**原始图像:**

""" 原始函数图像 """
    import numpy as np
    import matplotlib.pyplot as plt
    
    plt.switch_backend("TkAgg")
    
    x = np.linspace(0, 2 * np.pi, 10000)
    y = np.sin(x * 2) + np.sin(100 * x) * 0.3 + np.sin(200 * x) * 0.6
    
    plt.plot(x, y)
    plt.show()

![在这里插入图片描述][watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3h5eTEwMjg_size_16_color_FFFFFF_t_70_pic_center 1]  
**此函数包含一个低频,一个次高频,一个高频信号**

**进行傅里叶变换:**

""" 傅里叶变换 """
    import numpy as np
    import matplotlib.pyplot as plt
    
    plt.switch_backend("TkAgg")
    
    x = np.linspace(0, 2 * np.pi, 10000)
    y = np.sin(x * 2) + np.sin(100 * x) * 0.3 + np.sin(200 * x) * 0.6
    
    # 进行傅里叶变换(此时f为频谱)
    f = np.fft.fft(y)
    # 傅里叶变换后此时f为复数域域,将之变成实数域
    plt.plot(np.abs(f))
    plt.show()

![在这里插入图片描述][watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3h5eTEwMjg_size_16_color_FFFFFF_t_70_pic_center 2]  
**放大一下:**  
![在这里插入图片描述][watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3h5eTEwMjg_size_16_color_FFFFFF_t_70_pic_center 3]  
**我们会发现,在2,100,200的值上,分别得出了低频,次高频,高频的最大值**

### 2.2.连续数据傅里叶反变换与滤波 ###

**也可以通过反变换将频率域的函数反变换为时间域：**  
 x ( t ) = F − 1 ( X ( ω ) ) = 1 2 π ∫ − ∞ + ∞ X ( ω ) e i ω t d ω x(t)=F^\{-1\}(X(\\omega)) = \\frac\{1\}\{2\\pi\} \\int\_\{-\\infty\}^\{+\\infty\}X(\\omega)e^\{i \\omega t\}d\\omega x(t)=F−1(X(ω))=2π1∫−∞\+∞X(ω)eiωtdω  
**二中的函数存在低频、次高频和高频信号, 如果我们想去掉高频信号, 或者次高频或高频信号的话,可以利用傅里叶反变换**

### 2.2.1 过滤低频信号 ###

""" 傅里叶反变换 """
    import numpy as np
    import matplotlib.pyplot as plt
    
    plt.switch_backend("TkAgg")
    
    x = np.linspace(0, 2 * np.pi, 10000)
    y = np.sin(x * 2) + np.sin(100 * x) * 0.3 + np.sin(200 * x) * 0.6
    
    # 进行傅里叶变换(此时f为频谱)
    f = np.fft.fft(y)
    # 将低频信息附近的值设为0(3与-3是根据二中频谱图像得出)
    f[:3] = 0
    f[-3:] = 0
    
    # 进行傅里叶的反变换
    y = np.fft.ifft(f)
    # 取实数域部分
    y = np.real(y)
    
    plt.plot(x, y)
    plt.show()

![在这里插入图片描述][watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3h5eTEwMjg_size_16_color_FFFFFF_t_70_pic_center 4]  
**可以看到, 此时低频信号已经被过滤掉了. 这里叫做高通滤波, 相当于傅里叶展开,只保留高阶系数**

### 2.2.2 过滤高频信号 ###

""" 傅里叶反变换 """
    import numpy as np
    import matplotlib.pyplot as plt
    
    plt.switch_backend("TkAgg")
    
    x = np.linspace(0, 2 * np.pi, 10000)
    y = np.sin(x * 2) + np.sin(100 * x) * 0.3 + np.sin(200 * x) * 0.6
    
    # 进行傅里叶变换(此时f为频谱)
    f = np.fft.fft(y)
    # 将低频信息附近的值设为0(3与-3是根据二中频谱图像得出)
    f[4:250] = 0
    f[-250:-4] = 0
    
    # 进行傅里叶的反变换
    y = np.fft.ifft(f)
    # 取实数域部分
    y = np.real(y)
    
    plt.plot(x, y)
    plt.show()

![在这里插入图片描述][watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3h5eTEwMjg_size_16_color_FFFFFF_t_70_pic_center 5]  
**此时高频与次高频的信号被过滤了, 只剩下低频信号,叫做低通滤波. 相当于傅里叶展开,只保留低阶系数**

### 2.3.离散数据傅里叶反变换与滤波器 ###

**在⼯程实践中绝⼤部分使⽤的均是离散傅⾥叶变换，即现在计算机处理中所谓的傅⾥叶变换信号和频谱都是周期性的，均是离散的。**

> **1.连续傅⾥叶变换信号和频谱都是非周期性的；  
> 2.周期信号的频谱是离散的；  
> 3.离散信号的频谱是周期性的；  
> 4.离散傅⾥叶变换（Discrete Fourier Transform，DFT）的信号和频谱都是周期性的；**

**离散傅⾥叶变换计算公式:**  
 X ( k ) = ∑ n = 0 N − 1 x ( n ) e − i 2 π n N k X(k) = \\sum\_\{n=0\}^\{N-1\}x(n)e^\{-i\\frac\{2\\pi n\}\{N\}k\} X(k)=n=0∑N−1x(n)e−iN2πnk  
**反变换公式:**  
 x ( k ) = 1 N ∑ n = 0 N − 1 X ( n ) e i 2 π n N k x(k) = \\frac\{1\}\{N\} \\sum\_\{n=0\}^\{N-1\}X(n)e^\{i\\frac\{2\\pi n\}\{N\}k\} x(k)=N1n=0∑N−1X(n)eiN2πnk

## 二.傅里叶变换与卷积 ##

**傅⾥叶变换的另外⼀个作⽤就是滤波**  
 假 设 两 个 信 号 x 1 ( t ) , x 2 ( t ) 以 及 其 傅 ⾥ 叶 频 谱 X 1 ( ω ) , X 2 ( ω ) 假设两个信号x\_1(t),x\_2(t) 以及其傅⾥叶频谱X\_1(\\omega),X\_2(\\omega) 假设两个信号x1(t),x2(t)以及其傅⾥叶频谱X1(ω),X2(ω)  
 两 个 频 谱 进 行 相 乘 ： X 3 = X 1 ( ω ) X 2 ( ω ) 两个频谱进行相乘：X\_3=X\_1(\\omega)X\_2(\\omega) 两个频谱进行相乘：X3=X1(ω)X2(ω)  
 X 3 为 滤 波 后 的 信 号 , 频 谱 相 乘 等 于 时 间 阈 上 的 卷 积 ( 也 是 折 积 ) X\_3为滤波后的信号, 频谱相乘等于时间阈上的卷积(也是折积) X3为滤波后的信号,频谱相乘等于时间阈上的卷积(也是折积)  
 等 价 于 ： x 3 ( τ ) = x 1 ( t ) ∗ x 2 ( t ) = = ∫ − ∞ + ∞ x 1 ( τ ) x 2 ( t − τ ) d τ 等价于： x\_3(\\tau)=x\_1(t)\*x\_2(t)== \\int\_\{-\\infty\}^\{+\\infty\} x\_1(\\tau)x\_2(t-\\tau) d \\tau 等价于：x3(τ)=x1(t)∗x2(t)==∫−∞\+∞x1(τ)x2(t−τ)dτ  
 时 间 域 上 的 卷 积 可 以 代 替 频 谱 相 乘 方 式 的 滤 波 时间域上的卷积可以代替频谱相乘方式的滤波 时间域上的卷积可以代替频谱相乘方式的滤波 时 间 域 卷 积 计 算 即 是 频 率 域 乘 法 计 算 时间域卷积计算即是频率域乘法计算 时间域卷积计算即是频率域乘法计算

![在这里插入图片描述][watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3h5eTEwMjg_size_16_color_FFFFFF_t_70_pic_center 6]  
**进行卷积后得出:**

""" 卷积 """
    import matplotlib.pyplot as plt
    import numpy as np
    
    plt.switch_backend("TkAgg")
    
    # 设置plt可以写中文
    plt.rcParams['font.sans-serif'] = ['SimHei']
    plt.rcParams['axes.unicode_minus'] = False
    
    
    # x1称为信号, x2称为滤波器(也称卷积核心), 卷积也称为时间域滤波
    def conv(x1, x2):
        # 根据卷积公式, 将x2翻转
        x2 = x2[::-1]
        l1 = len(x1)
        l2 = len(x2)
        x3 = np.zeros([l1])
        for i in range(l1 - l2 + 1):
            x3[i] = np.sum(x1[i:i + l2] * x2)
        return x3
    
    
    x1 = np.linspace(-5, 5, 1000)
    y1 = np.sin(2 * x1) + 0.2 * np.sin(60 * x1)
    
    # 绘制函数
    plt.subplot(121)
    plt.plot(x1, y1)
    plt.title("滤波前")
    
    # 设置滤波器
    g = np.ones([32]) / 32
    # 进行卷积计算
    y1 = conv(y1, g)
    
    plt.subplot(122)
    plt.title("滤波后")
    plt.plot(x1, y1)
    plt.show()

![在这里插入图片描述][watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3h5eTEwMjg_size_16_color_FFFFFF_t_70_pic_center 7]  
**以上为一维数据处理**

## 三.图像数据处理 ##

**图像数据可以看成是二维数据.我们先看一下图像数据是什么样的:**

利用opencv库进行图像处理.安装库:

pip install opencv-python

**用我们喜爱的姚明举下例子 ^ . ^：**![在这里插入图片描述][watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3h5eTEwMjg_size_16_color_FFFFFF_t_70_pic_center 8]

import cv2
    
    # 读取图片
    img = cv2.imread("ym.png")
    # 打印图像信息
    print(img.shape, img.dtype)
    # 展示图片
    cv2.imshow("img", img)
    # 等待输入字符才会终止程序
    cv2.waitKey(0)

**打印的图片信息和类型:**

(574, 766, 3) uint8

**其中小括号里的574、766、3分别代表高、宽和颜色的三原色(顺序为蓝，绿，红)，uint8表示8位的无符号整形, 表示该矩阵中每一个数字都是8位的。**

### 3.1 过滤颜色 ###

import cv2
    
    # 读取图片
    img = cv2.imread("ym.png")
    # 打印图像信息
    print(img.shape, img.dtype)
    
    # 将蓝和红设置为0，只留绿色
    img[:, :, 0] *= 0
    img[:, :, 2] *= 0
    
    # 展示图片
    cv2.imshow("img", img)
    # 等待输入字符才会终止程序
    cv2.waitKey(0)

![在这里插入图片描述][watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3h5eTEwMjg_size_16_color_FFFFFF_t_70_pic_center 9]

### 3.2 图像的卷积 ###

![在这里插入图片描述][watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3h5eTEwMjg_size_16_color_FFFFFF_t_70_pic_center 10]  
**3.2.1 平均模糊**

import cv2 
    import numpy as np 
    
    # 滤波器，平均模糊
    kernel = np.ones([16, 16]) / (16*16)
    
    img = cv2.imread("ym.png")
    # 进行滤波（卷积运算）， -1 代表对每一个颜色单独的进行处理
    img = cv2.filter2D(img, -1, kernel)
    cv2.imshow("img", img)
    cv2.waitKey(0)

![在这里插入图片描述][watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3h5eTEwMjg_size_16_color_FFFFFF_t_70_pic_center 11]  
**3.2.2 边界滤波**  
一般图像中可以通过物体的颜色，来区分不同的物体，那么我们在设置滤波器的时候，可以将相同颜色的像素点进行过滤，只保留物体边缘处颜色不一样的地方。

**提取纵向边界，可以看到纵向部分更突出了:**

import cv2 
    import numpy as np 
    
    # 提取纵向边界
    kernel = np.ones([4, 4])
    kernel[:, :2] = -1
    """ [[-1. -1. 1. 1.] [-1. -1. 1. 1.] [-1. -1. 1. 1.] [-1. -1. 1. 1.]] """
    
    img = cv2.imread("ym.png")
    img = cv2.filter2D(img, -1, kernel)
    cv2.imshow("img", img)
    cv2.waitKey(0)

![在这里插入图片描述][watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3h5eTEwMjg_size_16_color_FFFFFF_t_70_pic_center 12]  
**提取横向边界，可以看到横向边界更突出了:**

import cv2 
    import numpy as np 
    
    # 提取横向边界
    kernel = np.ones([4, 4])
    kernel[:2, :] = -1
    """ [[-1. -1. -1. -1.] [-1. -1. -1. -1.] [ 1. 1. 1. 1.] [ 1. 1. 1. 1.]] """
    
    img = cv2.imread("ym.png")
    img = cv2.filter2D(img, -1, kernel)
    cv2.imshow("img", img)
    cv2.waitKey(0)

![在这里插入图片描述][watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3h5eTEwMjg_size_16_color_FFFFFF_t_70_pic_center 13]  
**3.2.3 更复杂的滤波，例如浮雕滤波**

import cv2 
    import numpy as np 
    
    # 浮雕滤波
    kernel = np.array([[-2, -1, 0], 
                       [-1, 1, 1], 
                       [0, 1, 2]])
    img = cv2.imread("ym.png")
    img = cv2.filter2D(img, -1, kernel)
    cv2.imshow("img", img)
    cv2.waitKey(0)

![在这里插入图片描述][watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3h5eTEwMjg_size_16_color_FFFFFF_t_70_pic_center 14]

### 3.2 视频处理 ###

**我们同样可以对视频进行卷积过滤**

import cv2 
    import numpy as np 
    
    # 滤波器，平均模糊
    kernel = np.ones([16, 16]) / (16*16)
    
    # 读取摄像的句柄,调去第一个摄像头
    cap = cv2.VideoCapture(0) 
    # 构建读取循环
    while True:
        ret, img = cap.read()
        img = cv2.filter2D(img, -1, kernel)
        cv2.imshow("img", img)
        # 每隔100ms没响应则继续读取
        ret = cv2.waitKey(100)
        if ret == 97:
            break
    cv2.destroyAllWindows()

![在这里插入图片描述][watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3h5eTEwMjg_size_16_color_FFFFFF_t_70_pic_center 15]

[watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3h5eTEwMjg_size_16_color_FFFFFF_t_70_pic_center]: /images/20221022/3d6e7c841800472c9992c2bfe330453b.png
[watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3h5eTEwMjg_size_16_color_FFFFFF_t_70_pic_center 1]: /images/20221022/f0acb0361cde4cb0981fc23db730feaf.png
[watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3h5eTEwMjg_size_16_color_FFFFFF_t_70_pic_center 2]: /images/20221022/8ec8ec184ceb41eab9c7041b993ad67a.png
[watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3h5eTEwMjg_size_16_color_FFFFFF_t_70_pic_center 3]: /images/20221022/d91922c817cb4835af492164e83da649.png
[watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3h5eTEwMjg_size_16_color_FFFFFF_t_70_pic_center 4]: /images/20221022/50612bc0174a46b4a66b15bbb511984c.png
[watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3h5eTEwMjg_size_16_color_FFFFFF_t_70_pic_center 5]: /images/20221022/90ac0d0e954b417a82a825bdc37dbf76.png
[watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3h5eTEwMjg_size_16_color_FFFFFF_t_70_pic_center 6]: /images/20221022/5127b05bacfa4272a5d307a11a2803f5.png
[watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3h5eTEwMjg_size_16_color_FFFFFF_t_70_pic_center 7]: /images/20221022/b262a395aedf41159751d31cf15dd4d4.png
[watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3h5eTEwMjg_size_16_color_FFFFFF_t_70_pic_center 8]: /images/20221022/297147b431b843faae95e80362d52ba4.png
[watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3h5eTEwMjg_size_16_color_FFFFFF_t_70_pic_center 9]: /images/20221022/4e3cc4c5bdc0453cbbcea6ef29a2d25c.png
[watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3h5eTEwMjg_size_16_color_FFFFFF_t_70_pic_center 10]: /images/20221022/737bd3cfc74d426688c87215826414d6.png
[watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3h5eTEwMjg_size_16_color_FFFFFF_t_70_pic_center 11]: /images/20221022/3e070c0c9b3a49e7a445f3a3d4c59364.png
[watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3h5eTEwMjg_size_16_color_FFFFFF_t_70_pic_center 12]: /images/20221022/e54cb8a9388b40e5825b1a9f41f2f587.png
[watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3h5eTEwMjg_size_16_color_FFFFFF_t_70_pic_center 13]: /images/20221022/bc67845a8e7844788c696d43ed8c04bc.png
[watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3h5eTEwMjg_size_16_color_FFFFFF_t_70_pic_center 14]: /images/20221022/8fab3dba2a4147dc95751d5b2d2f0cd6.png
[watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3h5eTEwMjg_size_16_color_FFFFFF_t_70_pic_center 15]: /images/20221022/86439ec4c1f84bb3835400577ccd35f5.png