通过示例了解傅立叶变换
2017年4月23日，作者Ritchie Vink

在过去的几周中，我一直在研究傅立叶变换的结果，它具有许多有趣的特性，
这些特性最初对我来说并不明确。在这篇文章中，我总结了一些我发现有趣
的东西以及我对傅立叶变换学到的东西。

应用
傅立叶变换用于工程中，以确定振动信号中的主频率。当信号的主频率与结
构的固有频率相对应时，由于共振，发生的振动会被放大。发生这种情况的
程度可能会导致结构崩溃。

现在说我买了一个新的音响系统，我的客厅窗户的自然频率约为100 Hz。让
我们使用傅立叶变换（Fourier Transform），并检查将Kendrick Lamar的歌
曲“ Alright”全部调大声是否安全。傅里叶变换，通常用于将时间谱中的信
号转换为频谱。时间频谱的示例包括声波，电，机械振动等。下图显示了
Kendrick曲调的0.25秒。可以清楚地看到，它看起来像是具有不同频率的波。
实际上，它看起来像是许多波浪。

歌曲《Kendrick Lamar- Alright》时间谱

傅里叶变换
这就是傅里叶变换的用处。这种方法利用了每个非线性函数都可以表示为
（无限多）正弦波之和的事实。在下面的图中对此进行了说明，因为阶跃函数
（矩形波）由许多个正弦波模拟综合而成。

正弦波模拟的阶跃函数（矩形波）
傅立叶变换将分解一个时间信号，并将返回有关模拟该时间信号所需的所有正
弦波频率的信息。对于均匀间隔的值序列，离散傅里叶变换（DFT）定义为：

式中:

N = 样本数
n = 当前样本
xn = 时间 n 的样本值
k = 当前频率 (0 Hz to N-1 Hz)
Xk =  DFT 的结果(振幅与相位)
注意点积定义为:

 DFT 算法 可写成python:

/////////////////////////////////////////
import numpy as np
def DFT(x):
    """
    Compute the discrete Fourier Transform of the 1D array x
    :param x: (array)
    """

    N = x.size
    n = np.arange(N)
    k = n.reshape((N, 1))
    e = np.exp(-2j * np.pi * k * n / N)
    return np.dot(e, x)

////////////////////////////////////////

但是，如果我们在时间信号上运行此代码，其中包含大约10,000个值，
则计算将花费10秒钟以上！哇...太慢了。
幸运的是，一些聪明的人（Cooley和Tukey）提出了快速傅立叶变换
（FFT）算法，该算法将DFT递归地划分为较小的DFT，从而大大减少了
所需的计算时间。A standard DFT scales O(N2) while the FFT
scales O(N log(N))

探索FFT
让我们写一些代码来找出FFT实际在做什么。
首先，我们定义一个简单信号，其中包含两个正弦波。一种频率为40 Hz，
另一种频率为90 Hz。

///////////////////////////////////

t = np.linspace(0, 0.5, 500)
s = np.sin(40 * 2 * np.pi * t) + 0.5 * np.sin(90 * 2 * np.pi * t)

plt.ylabel("Amplitude")
plt.xlabel("Time [s]")
plt.plot(t, s)
plt.show()
///////////////////////////////////

复数
为了检索上述时间信号频率的频谱，我们必须对该序列进行FFT。

/////////////////////////////////
fft = np.fft.fft(s)

for i in range(2):
    print("Value at index {}:\t{}".format(i, fft[i + 1]), "\nValue at
 index {}:\t{}".format(fft.size -1 - i, fft[-1 - i]))

>>>
Value at index 0: (0.0003804834928402556-0.060555031761900024j)
Value at index 499: (0.0003804834928403944+0.060555031761903175j)
Value at index 1: (0.0015317714831371565-0.12188808528069561j)
Value at index 498: (0.0015317714831373785+0.1218880852806919j)

////////////////////////

在上述代码段中，确定了两个正弦波的FFT结果。将FFT序列的前两个值和后两
个值打印到stdout。如我们所见，结果得到了复数。如果将序列的第一个值
（指标 index 0）与序列的最后一个值（指标 499）进行比较，我们可以看到两个数的
实部是相等的，而虚数的值在大小上也相等，只有一个是正的，另一个是负的。
数字彼此复共轭。对于序列中的所有数字都是如此。
对于实数输入是 n ，它是复数 N - n的共轭.

由于序列的后半部分没有提供任何新信息，因此我们已经可以得出结论，FFT序
列的后半部分是我们需要的输出。
FFT的复数输出编号包含以下信息：

• 振幅一定频率的正弦波（能量）的。
• 某个频率正弦波的相位偏移。

通过获取数字的绝对值来检索幅度，并通过计算数字的角度来获得相位偏移。

频谱

我们对每个频率的能量都感兴趣，因此我们可以确定FFT输出的绝对值。为了更好
地了解频谱，应针对频率绘制能量。FFT的每个离散数输出对应于某个频率。频率
分辨率取决于：

综上所述，我们可以绘制简单正弦波函数的频谱。我们只绘制频谱的一半，因为
那是唯一能提供真实信息的频谱。

/////////////////////////////
fft = np.fft.fft(s)
T = t[1] - t[0]  # sampling interval
N = s.size

# 1/T = frequency
f = np.linspace(0, 1 / T, N)

plt.ylabel("Amplitude")
plt.xlabel("Frequency [Hz]")
plt.bar(f[:N // 2], np.abs(fft)[:N // 2] * 1 / N, width=1.5)  # 1 / N
is a normalization factor

/////////////////////////////////

如我们所见，FFT起作用了！它为我们提供了有关时间信号中电波频率的信息。
FFT是时间信息和频率信息之间的权衡。通过对时间信号进行FFT，所有时间信息
都将丢失，以换取频率信息。为了将有关时间和频率的信息保持在一个频谱中，
我们必须制作一个频谱图。这些是离散时间窗口上的DFT。

《Alright》歌曲
通过对肯德里克•拉马尔（Kendrick Lamar）歌曲的时间信号进行FFT结果，我们
得到了如下所示的频谱。频率刻度绘制在对数刻度上。正如我们之前假设的，窗
户的固有频率约为100 Hz。从图中可以看到，最主要的频率出现在10 1.5 -
10 2.2 Hz（30-158 Hz）之间。我的窗户固有频率恰好在歌曲的主要频率的中间，
因此可能会因音量高而产生共鸣。
现在断言以全音量收听这首歌还为时尚早。但是，如果我真的想确定自己的窗户，
则应该检查另一首歌曲的频率。