Fourier analysis

傅里叶分析的核心很简洁:任何复杂的信号，都可以拆成一系列不同频率的正弦波叠加。它包括了两个工具:

• 傅里叶级数(Fourier Series).如果数据具有周期性（每隔T秒重复一次），通过:

\begin{equation} f(t) = \frac{a_0}{2} + \sum_{n=1}^{\infty} \left( a_n \cos(n\omega_0 t) + b_n \sin(n\omega_0 t) \right) \end{equation} 其中$\omega_0=\frac{2\pi}{T}$表示基频,$f(t)$是随时间$t$ 变化的函数。

• 傅里叶变换(Fourier Transform)。现实中的数据大多不是严格周期的。傅里叶变换把周期推广到无穷大:

\begin{equation} F(\omega) = \int_{-\infty}^{\infty} f(t) \, e^{-j\omega t} \, dt \end{equation} 其中,$F(\omega)$频域函数，表示信号在频率$\omega$处的成分

计算机处理的是离散采样点，因此，在实际用的是离散傅里叶变换（Discrete Fourier Transform,DFT）: \begin{equation} X[k] = \sum_{n=0}^{N-1} x[n] \, e^{-j2\pi kn/N} \end{equation} 其中,$N$表示采样点总数，$n$表示时域采样索引，$n=0,1,2,....,N-1$,$k$表示频域索引,$k=1,2,3...,N-1$

假设在一个河口附近的监测站，每小时记录一次水位（单位：米），连续记录 8 小时。由于上游来水和下游潮汐的共同作用，水位变化看起来“忽高忽低”，没有明显规律:$h=[1.8, 0.3, −0.5, −1.2, −0.3, 0.9, 1.1, −0.1]$。现在我们想知道这段水位变化里包含哪些周期成分？各占多大？

根据DFT公式,设$f_s$表示采样率，取值为1 次/小时；$\Delta f=\frac{f_s}{N}$表示频率分辨率,取值为0.125 次/小时；$f_k=k \cdot \Delta f$表示第$k$个频点对应的实际频率:

1. $k=0$时，频率$f_0=0$次/小时，周期 $T_0$趋于无穷大，那么$X[0]=\sum_{n=0}^7h(n)=2.0$
2. $k=1$时,频率$f_1=0.125$次/小时，周期 $T_1$=8小时，那么$X[1]=\sum_{n=0}^7h(n)e^{-j2\pi n/8}=1.8 \times 1+ 0.3 \times (0.707-j0.707)+ (-0.5)\times (-j)+...+(-0.1)\times (0.707+j0.707)= 3.726+j2.944$,幅值=$\sqrt{3.726^2+2.944^2} \approx 4.75$
3. $k=2$时,频率$f_2=0.25$次/小时，周期 $T_2$=4小时，那么$X[2]=-1−j2.414$
4. $k=3$时,频率$f_3=0.375$次/小时，周期 $T_3$=2.667小时，那么$X[3]=0.900−j2.500$
5. $k=4$时,频率$f_4=0.5$次/小时，周期 $T_4$=2小时，那么$X[4]=2.2$
6. $k=5,6,7$时，它们镜像频点，实信号 DFT 满足共轭对称，也即当输入信号是实数时，DFT 结果中正频率和负频率的复数值互为共轭，所以$|X[5]|=∣X[3]∣\approx 2.66,∣X[6]∣=∣X[2]∣\approx 2.61,∣X[7]∣=∣X[1]∣\approx 4.75$

因此，根据上面的计算，我们可以看到，

1. 平均水位：0.25 米（$k=0$直流分量）。
2. 主要周期成分：
   • 8 小时周期（$k=1$）：幅度约1.19米，最强。
   • 2.67 小时周期（$k=3$）：幅度约 0.66 米。
   • 4 小时周期（$k=2$）：幅度约 0.65 米
   • 2 小时周期（$k=4$）：幅度约 0.28 米。
   • 镜像频点（$k=5,6,7$）不提供新信息，是实信号 DFT 的对称性造成的。

因此，我们可以说，8 小时和 4 小时的强周期成分，可能对应上游水库调度或局部风生波动的周期。真实潮汐周期（12.42 小时、24.84 小时）超出了 8 小时窗口的分辨能力。