python-傅里叶变换

一般来说，频谱分析指的是将信号做傅里叶变换从而进行分析。频谱分析是包括幅频谱和相频谱两张图的。不过最常用的是幅频谱。

把时域的信息转换为频域表示，可以看出一些在时域状态下看不到的特征

滤波变得很容易，可以针对某个特征频率过滤

N个采样点，经过FFT之后，就可以得到N个点的FFT结果。为了方便进行FFT运算，通常N取2的整数次方。

假设采样频率为Fs，信号频率F，采样点数为N。那么FFT之后结果就是一个为N点的复数。每一个点就对应着一个频率点。这个点的模值，就是该频率值下的幅度特性。

第一个点对应的频率为0Hz（即直流分量），最后一个点N的下一个点对应采样频率Fs

假设原始信号的峰值为A，那么FFT的结果的每个点（除了第一个点直流分量之外）的模值就是A的N/2倍。而第一个点就是直流分量，它的模值就是直流分量
的N倍。而每个点的相位呢，就是在该频率下的信号的相位。

对一信号做傅立叶级数分解，第一项必为常数项，常数项自然对应的就是直流分量，谱为0频
FFT变换后第一条谱线是直流分量实际上就是信号的均值，也称直流分量。在许多应用中，做FFT之前，常常需要去均值，这样第一个点就是0了。

那么频域的“0”是什么呢？cos（0t）就是一个周期无限长的正弦波，也就是一条直线！所以在频域，0频率也被称为直流分量，在傅里叶级数的叠加中，它仅仅影响全部波形相对于数轴整体向上或是向下而不改变波的形状。

FFT之后得到的一系列复数，是波形对应频率下的幅度特征，注意这个是幅度特征（特征值）不是幅值。

为什么FFT后幅值要除以N/2

原因是DFT的频谱是用谱密度定义的，即它的幅值表示的是单位带宽的幅值。 [公式] 个离散点的DFT（我这里说的是实数DFT）将产生 [公式] 个频率点，频率的序号是从 [公式] ，需要注意，如果是复数，[公式] 点DFT将产生 [公式] 个频率点。所以 [公式] 点实数DFT以后，频谱带宽是 [公式] ，每个频率点占的带宽是 [公式] ，所以每个频率的实际幅值需要用DFT后的幅值乘以 [公式] ，也就是除以 [公式] 。但是注意，频率序号为 [公式] 和 [公式] 的两个点带宽只占中间频率点的一半，也就是占 [公式] 的带宽，所以首尾两个点的幅值需要乘以 [公式] ，也就是除以 [公式] 。FFT是DFT的快速算法，所以FFT也需要这样处理，才能得到每个频率点真正的幅值。

X 轴表示频率，Y轴表示振幅

假设采样频率Fs，信号频率F，信号长度L，采样点数N。那么FFT之后结果就是一个为N点的复数。每一个点就对应着一个频率点。这个点的模值，就是该频率值下的幅度特性。

偶数项的振幅都是 0，
相位单位是度，就是对应圆位置的角度。 波形可以看做圆周运动在直线上的映射

FT（Fourier Transformation）：傅里叶变换。

DTFT（Discrete-time Fourier Transform）：离散时间傅里叶变换。这里的“离散时间”指的是时域上式离散的，也就是计算机进行了采样。不过傅里叶变换后的结果依然是连续的。

DFT（Discrete Fourier Transform）：离散傅里叶变换。在DTFT之后，将傅里叶变换的结果也进行离散化，就是DFT。

也就是说：FT时域连续、频域连续；DTFT时域离散、频域连续；DFT时域离散、频域离散。

FFT（Fast Fourier Transformation）：快速傅里叶变换。就是DFT的快速算法，一般工程应用时用的都是这种算法。

FS（Fourier Series）：傅里叶级数。是针对时域连续周期信号提出的，结果是离散的频域结果。

DFS（Discrete Fourier Series）：离散傅里叶级数。是针对时域离散周期信号提出的，DFS与DFT的本质是一样的。

FFT代码

import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
from matplotlib.pylab import mpl

from scipy import fft, ifft

mpl.rcParams['font.sans-serif'] = ['SimHei']  # 显示中文
mpl.rcParams['axes.unicode_minus'] = False  # 显示负号


def load_data(path):
    file = open(path, 'r')
    Lines = file.read().splitlines()
    df = pd.DataFrame(Lines, columns=['t1'], dtype=float)
    return df
    pass

def draw_raw_figure(df):
    plt.subplot(211)
    N = len(df['t1'])
    x = np.arange(N)  # 频率个数
    y = df['t1']
    plt.grid(True, linestyle='-.')
    plt.plot(x, y, 'r', label='raw data')

    plt.show()
    pass


def harmonic(df):
    y = df['t1']
    fft_y = fft(df['t1'])  # 快速傅里叶变换
    # FFT之后得到的那一串复数是波形对应频率下的幅度特征
    # 第一个数是 0Hz频率叫直流分量，
    N = len(df)
    x = np.arange(N)  # 频率个数
    half_x = x[range(int(N / 2))]  # 取一半区间

    """
    FFT变换后的复数模 - 幅度
    真实幅值，需要把除了第1个点(i=0)以及最后一个点(i=N/2)除以N以外，其余点需要求得的模除以N/2
    """
    fft_y[0] /=N
    fft_y[len(fft_y) -1] /=N
    fft_y[1:len(fft_y)-2] /= N/2
    abs_y = np.abs(fft_y)  # 取复数的绝对值，即复数的模(双边频谱)
    angle_y = np.angle(fft_y)  # 取复数的角度
    # normalization_y = abs_y / N  # 归一化处理（双边频谱）
    normalization_half_y = abs_y[range(int(N / 2))]  # 由于对称性，只取一半区间（单边频谱）

    plt.subplot(211)
    plt.plot(x, y)
    plt.title('原始波形')
    plt.grid()

    plt.subplot(212)
    plt.bar(half_x[1:], normalization_half_y[1:])
    # plt.plot(half_x[1:], normalization_half_y[1:], 'blue')
    plt.title('单边振幅谱', fontsize=9, color='blue')
    plt.grid()


    sinwave = normalization_half_y[1] * np.sin(2 * np.pi * 200 * x)

    plt.show()
    pass


if __name__ == "__main__":
    # spc的数据
    df = load_data('d:\\1.txt')
    # draw_raw_figure(df)
    harmonic(df)
    pass