瞬时频率

在信号处理中，观察信号的瞬时频率是很重要的课题。假设一实信号 $x(t)\,$ 可写成指数信号的N项相加(有无穷多种表示法，以 $N\,$ 小的为宜)，即
$x(t)=\sum _{k=1}^{N}a_{k}\cdot e^{j\cdot \phi _{k}(t)}$， 其中 $a_{k}\,$ 为虚常数。
则瞬时频率(以频率表示)$f_{k}(t)={\frac {\phi _{k}^{\prime }(t)}{2\pi }}$， k=1,...,N

瞬时频率 (Hz)

以频率来表示(单位为赫兹)： $f_{k}(t)={\frac {\phi _{k}^{\prime }(t)}{2\pi }}$， k=1,...,N

瞬时频率 (rad/s)

以角频率来表示(单位为弧度每秒)： $f_{k}(t)=\phi _{k}^{\prime }(t)$， k=1,...,N

以解析讯号法定义瞬时频率

直观上，瞬时频率为相位的微分。对于自然界中的实数讯号，如何定义讯号的相位。Gabor提出解析讯号法(Analytic Signal Method)，将实数讯号表示为对应的复数讯号，即可定义复数讯号的大小与相位，将实数讯号的瞬时频率求出。 实数讯号$s\left(t\right)$的解析讯号(Analytic Signal)$z\left(t\right)$定义为

解析函数的极座标表示

瞬时相位

瞬时频率

$z(t)=s(t)+{\frac {j}{\pi }}\int _{-\infty }^{\infty }{\frac {s(\tau )}{t-\tau }}\,d\tau .\,$

其中虚数项为实数讯号$s\left(t\right)$的希尔伯特转换(Hilbert Transform)，将它定义为${\widehat {s}}\left(t\right)$。称作解析函数的理由是，此型式的复数函数满足柯西-里曼(Cauchy-Riemann)的可微分条件，称之为解析函数(Analytic Function)。因此，解析讯号$z\left(t\right)$可以表示为

$z(t)=s(t)+j{\widehat {s}}\left(t\right)=m(t)\cdot e^{j\theta (t)}\,$

其中

$m(t)={\sqrt {s^{2}(t)+{\widehat {s}}^{2}(t)}}$ ; $\theta (t)=\arctan \left({{\widehat {s}}(t) \over s(t)}\right)\,$

如果$s\left(t\right)$是没有任何限制条件的时间讯号，计算出来的瞬时频率可能不是正确的结果。对于平均值为零的局部对称讯号而言，前述定义的瞬时频率才具有物理意义。在1998年，黄锷(Norden E. Huang)博士提出一个有效的算法，将讯号先行分解成具有局部对称之分量，以正确地求得资料的瞬时频率。这个方法称为希尔伯特-黄转换(Hilbert Huang Transform, HHT)。

以下简单的例子来说明，对于平均值为零的讯号，此瞬时频率的定义才具有物理意义。对于一个弦波讯号$s\left(t\right)$，

$s\left(t\right)=\beta +\cos {(t)}$

考虑三种情况: (1) $\beta =0$ (2) $0<\beta <1$ (3) $\beta >1$

(1) $\beta =0$: 当弦波讯号平均值为零时，$z\left(t\right)$在复数平面上的描述是以座标原点为中心的单位圆，它的相位角$\theta (t)$则是以座标原点为中心，反时针方向呈线性递增，其图形为斜率1的直线，而瞬时频率是一个常数值。

(2) $0<\beta <1$: $z\left(t\right)$在复数平面上仍然是一个单位圆，但其圆心从原点偏移了$\beta$个单位，其相角$\theta (t)$不再呈现线性递增，瞬时频率出现震荡的现象，而不是应有的常数值。

(3) $\beta >1$: 因为$\beta$值超过单位圆的半径，因此$z\left(t\right)$的圆心在单位圆之外。如此相位角$\theta (t)$在[$-\pi$/$2$, $\pi$/$2$]震荡，瞬时频率出现负值，与原讯号的特性有极大的差别。

观察瞬时频率的重要性

因为在目前许多数位信号处理的应用上都与信号的频谱或信号的带宽有很大的关系。
若能确实地侦测信号的瞬时频率，则通道带宽可以被可适性(adaptive)的决定，如此一来能更有效地利用系统资源，提高系统效能。

相关应用

• 调变(modulation)
• 多工方式(multiplexing)
• 滤波器的设计(filter design)
• 信号压缩(data compression)
• 信号分析(signal analysis)
• 信号辨识(signal identification)
• 语音信号处理(acoustical signal processing)
• 制作系统的模型(system modling)
• 雷达系统的分析(rader system analysis)
• 取样(sampling)

如何观察信号的瞬时频率

瞬时频率为常数－使用傅里叶变换

当瞬时频率为常数即 $\phi (t)\,$ 为一阶时间函数，使用傅里叶变换做信号分析。
由于从傅里叶变换中是无法观察出信号频谱随着时间改变的变化。
故只有当瞬时频率为常数，不是时间的函数时，便可使用傅里叶变换做信号分析。

瞬时频率不为常数－使用时频分析

当瞬时频率不为常数即 $\phi (t)\,$ 为高阶时间函数，使用时频分析做信号分析。
从时频分析可观察出信号频率随着时间变化的改，这是傅里叶变换无法做到的。
因此当瞬时频率为时间的函数，使用时频分析做信号分析，如下图，可以确切地观察到信号瞬时频率的变化。
$x(t)={\begin{cases}cos(\pi t)\ \ \ t\leq 10\\cos(3\pi t)\ \ \ 10<t\leq 20\ \ \ \\cos(2\pi t)\ \ \ t>20\end{cases}}$

解析讯号法

信号为单一的sinusoid曲线－使用希尔伯特转换 (Hilbert Transform)

当信号只有一个成分，且曲线是周期性的，具有固定的振幅、频率和相，可以使用希尔伯特转换，计算信号的相位求瞬时频率。

瞬时频率$\,={\frac {1}{2\pi }}{\frac {d}{dt}}\theta ,\quad \theta =tan^{-1}{\frac {x_{H}(t)}{x_{(}t)}}$

$x_{H}(t)={\frac {1}{\pi }}\int _{-\infty }^{\infty }{\frac {x(\tau )}{t-\tau }}d\tau \,$

举例:

$cos(2\pi ft)\;\xrightarrow {Hilbert} \;sin(2\pi ft)$

$\theta =tan^{-1}(tan(2\pi ft))=2\pi ft$

${\frac {1}{2\pi }}{\frac {d}{dt}}\theta =f$

信号为多个成分的非sinusoid曲线－使用希尔伯特-黄转换 (Hilbert Huang Transform)

当信号为复数函数、非sinusoid曲线、有多个成分或$\theta$有多个解，此时可以先将信号分成sinusoid-like成分和趋势，再运用Hilbert transform或短时距傅里叶变换 (STFT) 求零交叉的数量。

1. 使用经验模态分解 (EMD) 将信号分解为一系列本征模态函数 (IMFs）的振荡模态和趋势

2. 对每个IMF算瞬时频率:

• 使用Hilbert transform
• 计算STFT，瞬时频率 = ${\frac {{\text{number of zero-crossings between}}\;t-B\;{\text{and}}\;t+B}{4B}}$
• 直接计算零交叉，零交叉 = ${\frac {\text{number of periods}}{2}}$

同时参阅

• 信号处理
• 傅里叶变换
• 时频分析

参考资料

• Jian-Jiun Ding, class lecture of Time Frequency Analysis and Wavelet transform, Graduate Institute of Communication Engineering, National Taiwan University, Taipei, Taiwan, 2007.
• Leon Coen, Time-Frequency Analysis, Prentice Hall, 1995.
• 陈韦佑, "以希尔伯特-黄转换法为GPS接收机抑制调频干扰", 国立台湾大学电机工程研究所硕士论文, 2007.