【傅里叶变换】DFT详细分析与Python实现

uestc_hjw

已于 2022-10-24 20:29:15 修改

阅读量4.2k

点赞数

文章标签： 1024程序员节 python 傅里叶分析

于 2022-10-24 20:27:47 首次发布

本文链接：https://blog.csdn.net/zhengyuyin/article/details/127499584

版权

一. 思考几个问题

我们知道无论是Matlab还是Python实现FFT都是十分容易的，但是需要思考几个问题：

狄拉克函数？sinc函数？
DFT的计算公式是什么？DFT后长度不变？
频谱的横坐标有何物理意义？横轴范围是多少？分辨率又是多少？
采样点数、采样频率、采样持续时间是否理解？更改采样点数或采样频率，对频谱有何影响？
频率谱幅值和时域信号幅值的关系？
频谱中F[0]的物理意义是什么？
为什么得到的频谱（除去第一个点）是中心对称的？
对图像进行FFT，频谱原点在哪里？如何进行中心化？
对图像进行FFT，得到的图像物理意义是什么？

二. 狄拉克函数和sinc函数

三. DFT的计算公式

四. 频谱幅值和时域信号幅值的关系

注意：这里所说的频谱，指的是FFT后取模得到的曲线。

五. 频谱的中心对称性

六. 取样和频率间隔间的关系

七. Python实现一维FFT

''' Author: Jiuwu Hao '''
import numpy as np
from scipy.fftpack import fft,fft2,fftshift
from math import pi
import matplotlib.pyplot as plt

def hjw_dft(fs,signal):
    # 信号持续时间T,即样本点数*采样间隔
    T = len(signal)*(1/fs)
    # x轴整体范围时采样持续时间，分辨率是采样间隔
    x = np.arange(0,T,1/fs)
    plt.plot(x,signal)
    plt.title('Original')
    plt.show()

    # dft后，长度不变，是复数表示，想要频谱图需要取模
    dft = fft(signal)
    dft = np.abs(dft)/(len(dft)/2)
    dft[0] /= 2
    # 归一化
    # dft = (dft-np.min(dft))/(np.max(dft)-np.min(dft))
    # 在频率域中，横轴范围为采样频率，分辨率是1/T，即持续时间的倒数
    dft_x = np.arange(0,fs,1/T)
    plt.plot(dft_x,dft)
    plt.title('DFT')
    plt.show()

    # 由于Nyquist采样定律和信号的对称性都可以证明，只要取一半的频率区间即可
    dft_half = dft[0:int(len(dft)/2)]
    dft_x_half = dft_x[0:int(len(dft)/2)]
    plt.plot(dft_x_half,dft_half)
    plt.title('DFT_Half')
    plt.show()


fs = 2000 # 采样率为2000Hz
n = 10000 # 采样点数为10000
T = n*(1/fs)  # 采样持续时间
# 注意：在构造测试函数时，不能x = np.arange(n)，要符合实际物理意义的采样频率
x = np.arange(0,T,1/fs)  
# 信号频率为100，200Hz,直流分量为4
y = 3*np.sin(2*pi*100*x) + 5*np.cos(2*pi*200*x) +4  

if __name__ == '__main__':
    hjw_dft(fs,y)

代码效果：