青山还在，夕阳红了好几次——杨慎的《临江仙》·滚滚长江东逝水

---

前言

本文在参考互联网上许多关于傅里叶变换的信息的基础上，从如何使用程序将采样获得的值转换为傅里叶，总结了在程序中采样获得的离散值（FFT）相关流程。本文不仅总结了众所周知的时间傅里叶变换，还总结了空间傅里叶变换，还总结了如何在时间和空间上采样获得的值。本文采用平面波y=A*sin(w*t - k*x)例如，通过三个简单的程序演示上述三种情况：①固定空间是一个点，采样得到空间点随时间变化的一维数值序列，然后傅里叶变换这个一维序列，得到频谱图（A'(f)图）。②固定时间是一个时刻。此时，采样得到随空间变化的一维数值序列，然后变换一维序列，获得波矢&振幅图（A'(k)图）。③对时间范围内和空间范围内采样得到的数值矩阵进行二维傅里叶变换，以得到色散图（A'(k,f)图）。这些内容是适用于各个领域的一般基本技能。最后，以微磁学中的自旋波为例，通过分析现有频谱和色散分析程序包的源代码，演示如何从微磁模拟获得的数据中提取自旋波的频谱和色散信息。本文不涉及傅里叶变换的原理和公式，而是将其视为信号分析工具，演示普通人如何操作。

---

一、傅里叶变换的作用

1.对傅里叶变换的直观理解

假如您在网上看过以下关于傅里叶变换的信息： 2014年更新了傅里叶分析的掐死教程(完整版).06.06 Python 中 FFT 傅里叶快速分析 图像理解二维傅里叶转换 【numpy】几种fft函数的使用 所以我相信你对傅里叶变换有一个基本的了解，如果你没有看到它也没关系！就我而言，我不需要深入了解傅里叶变换的原理，但我必须学习如何使用它FFT处理采样数据并获得最终图像，如频谱图，如色散图。因此，对傅里叶变换功能的直观理解应从其结果图像中理解。例如，在频谱图中，可以看到原始信号中每个频率重量的比例，系统可以根据频率带宽的分布获得允许频率和禁止频率。 维基百科全书上的图像生动地表达了傅里叶变换的功能：时域上复杂的原始信号可以分解为时域中不同频率和振幅的简单正弦波，并绘制每个正弦波的频率和相应和相应振幅的结果，即频谱图 A ~ ( f ) \tilde{A}(f) A~(f)。 事实上，我认为在原始信号中获观的例子，即三棱镜对白光的色散：

当一束白光（输入的原始信号）通过三棱镜时，三棱镜将白光扩展到不同的颜色（输出的每个频率波）。然后傅里叶变换就像三棱镜。当原始信号变换傅里叶时，可以获得形成原始信号的频率贡献。

2.傅里叶变家

傅里叶改变家族有四种获取相应频谱的方法：

（a）傅里叶连续变换（FT）：经过时域连续、非周期波形FT之后，得到的频谱是连续的，非周期的。 （b）傅里叶级数连续（FS）：在时域域连续、周期性的波形FS之后得到的频谱是离散的，非周期的。 （c）傅里叶离散变换（DTFT）：经过时域上的非周期波形DTFT之后，获得的频谱是连续的、周期性的。 （d）傅里叶级数离散（DFS）：经过时域上的波形DFS之后，得到的频谱是离散和循环的。

就实际情况而言，我们在等待时间间间隔采样和保存原始信号获得的是时域上的离散值，计算机在傅里叶转换这些输入的离散值后获得的结果只能是离散值。因此，上述四种方法似乎只适用于计算机DFS，只有该方法的输入输出是离散值，但由于其对输入的离散值有周期性要求，因此有其改进方法DFT，这种方法对输入的离散值没有特别的限制在大多数情况下，傅里叶变换是指DFT。DFT输入的离散值如下： 假设对一时变信号进行采样后得到N个数值： x ( 0 ) ， x ( 1 ) ， x ( 2 ) . . . x ( N ? 1 ) x(0)，x(1)，x(2)...x(N-1) x(0)，x(1)，x(2)...x(N?1)，那把这N个数值作为DFT的输入，对于N个输入的数值来说： f ( k ) = ∑ n = 0 N − 1 x ( n ) ⋅ e − j ⋅ 2 π ⋅ n N ⋅ k , k = 0 ， 1 ， 2... N − 1 f(k)= \displaystyle\sum_{n=0}^{N-1}{x(n)} \sdot e^{-j \sdot \frac{2\pi \sdot n}{N}\sdot k},k=0，1，2...N-1 f(k)=n=0∑N−1​x(n)⋅e−j⋅N2π⋅n​⋅k,k=0，1，2...N−1 这个公式表示：①假设DFT有N个输入数值的话，那么经过傅里叶变换后也得到了N个输出点，因此频谱图中就有N个点。

②式中的 e − j ⋅ 2 π ⋅ n N ⋅ k e^{-j \sdot \frac{2\pi \sdot n}{N}\sdot k} e−j⋅N2π⋅n​⋅k是一个复数，因此DFT的结果也是一个复数。但是对于第一个输出点 f ( k = 0 ) f(k=0) f(k=0)来说，它的DFT的结果值为 f ( 0 ) = ∑ n = 0 N − 1 x ( n ) ⋅ e − j ⋅ 2 π ⋅ n N ⋅ 0 = x ( 0 ) + x ( 1 ) + x ( 2 ) + . . . + x ( N − 1 ) f(0)= \displaystyle\sum_{n=0}^{N-1}{x(n)} \sdot e^{-j \sdot \frac{2\pi \sdot n}{N} \sdot 0}=x(0)+x(1)+x(2)+...+x(N-1) f(0)=n=0∑N−1​x(n)⋅e−j⋅N2π⋅n​⋅0=x(0)+x(1)+x(2)+...+x(N−1) ，即第一个点的DFT结果是一个实数值。

③对于输出的每一个点来说，DFT都需要进行N次运算，对于总共有N个输入的采样值来说，那么DFT的算法运算复杂度为 O ( N 2 ) O(N^2) O(N2)。可见输入的点数越多，运算次数也越多。

囿于原始DFT算法的运算速度，于是一种用于加速DFT计算的算法横空出世，那就是快速傅里叶变换（FFT），相比DFT，它的运算复杂度为 N l o g N NlogN NlogN，当输入的采样点数越多，FFT算法的速度相比原始的DFT快得多。所以在程序的数学库中用于一维傅里叶变换的函数为fft()，该函数便实现了上式的运算流程。

二、时间、空间、时空的傅里叶变换

1.时间的傅里叶变换

首先需要明白一点，用于当做时间傅里叶变换的输入的采样值是怎么获得的？在现实世界中，或者通过示波器的探头接触某个器件的输出位置点进行等时间间隔采样，以此获得N个随时间变化的采样值；或者通过放在空间中某一位置点的传感器等时间间隔采样，也能获取N个随时间变化的采样值。。。它们都有一个特点：把空间固定为一点，只等时间间隔采样在该空间位置点上的原始输出信号。

假设原始信号的最高频率分量为f（非角频率w），时间采样频率为Fs_time，采样点数为N。 以上三个条件之间有如下关系：①根据奈奎斯特采样定理，采样频率至少是原始信号的最高频率分量的两倍，即 F s _ t i m e ≥ 2 ∗ f Fs\_time≥2*f Fs_time≥2∗f。②采样频率的倒数即采样周期Ts，采样周期乘以采样点数就是采样时间。 接下来，将从时间的FFT得到的频谱图为切入点，来讲解如何获取频谱图中每一个点，并绘制出频谱图的基本流程。 首先得到频谱图中每个点的横坐标： 对N个采样数值进行FFT后，也有N个复数结果值： f ( k ) ， k = 0 ， 1 ， 2.. N − 1 f(k)，k=0，1，2..N-1 f(k)，k=0，1，2..N−1，它们对应N个频率点，这些频率点对应的频率值为 F s _ t i m e N ⋅ k \frac{Fs\_time}{N} \sdot k NFs_time​⋅k ,其中 F s _ t i m e N \frac{Fs\_time}{N} NFs_time​称为频率分辨率，即频谱图（ A ~ ( f ) \tilde{A}(f) A~(f)）中的两个点横坐标之间的最小间隔，因此要得到更小的频率分辨率（也就是让频谱曲线看起来更加连续），就只能增加采样点数（但不能降低采样频率），也就是延长采样时间。 比如当频率点k=0时，其代表的频率值即为 f = F s _ t i m e N ⋅ 0 f=\frac{Fs\_time}{N} \sdot 0 f=NFs_time​⋅0 = 0Hz，当频率点k=1时，其代表的频率值即为 f = F s _ t i m e N ⋅ 1 f=\frac{Fs\_time}{N} \sdot 1 f=NFs_time​⋅1 = F s _ t i m e N \frac{Fs\_time}{N} NFs_time​Hz。。。于是便得到了频谱图中所有点的横坐标。显然，如此计算的话，频谱图中的点的最大横坐标为Fs_time，也就是说，频谱图的横轴范围是受时间采样频率Fs_time影响的。 接着是得到频谱图中每个点的纵坐标： 对于FFT得到的除了第一个点（0Hz）以外的所有频率点，其傅里叶振幅 A ~ \tilde{A} A~为：首先对该点的复数结果取模，接着再除以采样点数N，最后再乘以2即可。对于第一个点（即直流分量0Hz），它的傅里叶振幅为取模之后，只除以采样点数N即可。如此，便得到了频谱图（ A ~ ( f ) \tilde{A}(f) A~(f)）的N个点的纵坐标。 看起来有了频谱图中每个点的横坐标和纵坐标，按理说就可以绘制出一个完整的频谱图了。但是，实际上，对于FFT之后的前面一半的频率点k=0，，， N 2 − 1 \frac{N}{2}-1 2N​−1，它们对应的频率确实为0Hz，，， F s _ t i m e N ⋅ ( N 2 − 1 ) \frac{Fs\_time}{N} \sdot (\frac{N}{2}-1) NFs_time​⋅(2N​−1)Hz，但是对于后一半部分的频率点k= N 2 \frac{N}{2} 2N​，，，N-1，它们对应的频率为 − N 2 ⋅ F s _ t i m e N -\frac{N}{2} \sdot \frac{Fs\_time}{N} −2N​⋅NFs_time​，，， − F s _ t i m e N -\frac{Fs\_time}{N} −NFs_time​。显然对于时间的傅里叶变换得到的频谱图来说，负频率没有物理意义，所以需要舍弃负频率的点，使频谱图的宽度减少一半。

什么，你居然问我怎么晓得如此操作的，有什么依据吗？

其实，网上关于傅里叶变换的资料众多，质量也是参差不齐，让人找不到北，而恰好MATLAB官方的fft()函数的文档中有关于时间的傅里叶变换的例程，见MATLAB官方参考链接，这个示例演示了如何获取频谱。而更巧的是Python的Numpy库的fft()例程居然没有获取信号频谱的例子，那我就搬一下砖，直接将MATLAB官方的傅里叶变换的例程移植到Python程序中，如下：

#################
#author:YQYUN
#date:22/6/18
#desc:从MATLAB官方例程移植而来，演示时间的傅里叶变换
#desc:固定空间为一个点，得到在这个位置点上采样得到时间范围内的信号
#################

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

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

#时间采样频率
Fs_time =1000
#采样周期
T = 1/Fs_time
#采样点数
M = 1500
#采样时间序列
t = np.arange(0, M-1) * T

# 设置需要采样的信号
y = 0.7 * np.sin(2 * np.pi * 50 * t)+np.sin(2 * np.pi * 120 * t)
#加入噪声信号
y = y + 2 * np.random.randn(np.size(t))

#画原始波形
plt.subplot(121)
plt.plot(1000*t[0:50], y[0:50])
plt.xlabel('t(s)')
plt.ylabel('A')
plt.title('原始波形')

#快速傅里叶变换得到一个复数序列，取复数的绝对值，即复数的模(双边频谱)
abs_y = np.abs(np.fft.fft(y))
#除了第一个点（0Hz）外，对模值要先除以采样点数量，再乘以2
fft_y = abs_y / M * 2
#对于第一个点（0Hz），只需对模值除以采样点数量即可 
fft_y[0] /= 2
#幅度取一半区间（单边频谱）
fft_y = fft_y[0:int(M / 2)]
#频率取一半区间
F = np.arange(0,int(M / 2)) * Fs_time / M

#将FFT结果保存到文件
FFTResult_time = np.zeros((int(M / 2),2))
FFTResult_time[:,0] = F[:]
FFTResult_time[:,1] = fft_y[:]
np.savetxt('FFTResult_time.txt', FFTResult_time,delimiter=' ')
print('时间FFT结果已经写入FFTResult_time.txt，第一列为频率，第二列为振幅')

#画图
plt.subplot(122)
plt.plot(F, fft_y)
plt.xlabel('f(Hz)')
plt.ylabel('A')
plt.title('单边频谱')
plt.show()

首先看一下运行结果，对比MATLAB（上）和Python（下）的结果，如下：

由于加入了随机的噪声信号，所以两者的时域波形是不一样的，比如再次运行一下这个程序： 以上三个图像的时域波形看起来都是乱七八糟的，没有什么规律可言，但是从它们的频谱图中可以很明显的看出50和120H