重新分布法

重新分布法（Reassignment Method）是一种锐化讯号之时频表示的方法，借由将资料映射至较靠近原始讯号之真实支撑区域的时频座标来实现。此方法曾被不同学者独立提出，并有重映射（remapping）、时频重新分布（time-frequency reassignment）以及修正滑动视窗法（modified moving-window method）等别称。以时频谱或短时距傅立叶变换而言，重新分布法可借由估算局部的瞬时频率以及群延迟，使模糊的时频资料点重新定位并清晰化。当讯号可借由分析视窗进行时域和频域的分离时，这项时频座标的重新映射是相当精准的。

## 概述

<image>

许多信号的能量分布都具有随着时间和频率变化的特性。举例而言，任何有起始或者结束的音频讯号之能量分布都将随着时间改变，而大多数声音的能量不论在时间上还是频率上都有相当大的变化；时频分析便经常用于寻找此类讯号的特征，这种分析法将原为一维的时域讯号转换为二维的时域–频域函数。而经分析转换后的结果便称为时频表示，它可以用于描述一个讯号随时间的频谱能量分布变化，就像乐谱描述乐曲中高低音随时间的变化一样。

在音频讯号分析中，时频谱是最常用的时频表示法；这可能是因为该方法长久以来为学界所熟知，且能够避免产生不易解读的“交叉项”（cross-terms）。不过，计算时频谱所须使用的视窗函数使得分析结果必须在时域或频域的分辨率中做一取舍，也就是时频谱在时域及频域的分辨率必定至少有其中一者不理想。这项时频表示的重新分布法便是透过先把资料映射至较接近原始讯号之真实支撑区域的座标，以将时频谱或其它时频表示下较模糊的资料点重新定位并聚焦。

## 时频谱

时频谱为最为人所熟知的时频表示法之一，其定义为短时距傅立叶变换值的平方：

$$SP_x(t,f)=|X(t,f)|^2=\left|\int_{-\infty}^{\infty} w(t-\tau)x(\tau)e^{-j2\pi f \tau}d\tau\right|^2.$$

虽然短时距的相位频谱确实含有与讯号相关的重要时态资讯，但这些资讯通常难以解读，故一般而言在短时距频谱分析中只考虑其短时距强度频谱。以时频表示法而言，时频谱具有相对较差的分辨率。它在时域与频域的分辨率取决于分析时所选择的视窗函数，且其中一者愈清晰将导致另一者愈模糊。

相对于时频谱，维格纳–魏尔分布则是另一种分辨率较高的时频表示法；本质上，维格纳–魏尔分布的做法与短时距傅立叶变换类似，惟其选择的视窗函数与讯号之间有完好的对应。维格纳–魏尔分布在时域和频域中皆可得到高分辨率，但它同时也具有高度的非线性和非局部性等特质，因此对噪声相当敏感且容易产生覆盖讯号本身的交叉项，使得含有多成分的讯号分析后的资讯难以解读。

利用双线性转换得到时频表示的科恩系列分布函数则可视为“平滑化”的维格纳–魏尔分布。它借由一平滑核降低分布对于噪声的敏感度，并抑制交叉项的产生，但此举将使时域和频域的分布模糊化；这不仅会降低频谱的分辨率，同时也可能使原始维格纳准几率分布中实际上没有能量的区域得到非零的能量值。

时频谱便属于科恩系列分布之一；它也可视为一种平滑化的维格纳–魏尔分布，惟其使用之平滑核即相等于其分析视窗函数的维格纳–魏尔分布。重新分布法便是先将维格纳–魏尔分布平滑化，再将它于原始讯号成分之真实支撑区域附近的座标点重新聚焦，进而达到锐化的目的。实验结果显示，此方法可以减少任一科恩系列分布在时域与频域的模糊现象。对时频谱而言，重新分布法便是利用讯号的短时距相位频谱来修正其时频座标点，再将其重新映射回较靠近所分析讯号之真实支撑区域的位置。

## 方法

重新分布法首先由小寺邦彦（Kunihiko Kodera）、罗杰·干德林（Roger Gendrin）和克劳地·德·凡德瑞（Claude de Villedary）三人提出，并命名为“修正滑动视窗法”（Modified Moving-Window Method）。他们的方法为每个资料点指定一组新的时频座标，使其更佳反映原始讯号的能量分布，进而改善了标准滑动视窗法（等同于时频谱）的时频分辨率。

标准的滑动视窗法乃是将一时域讯号$x(t)$分解为一组基底讯号$h_{\omega}(t)$的加权总和，其中代表每项基底的系数则为$\epsilon( t, \omega )$。基底讯号的定义如下：

$$h_{\omega}(t) = h(t) e^{j \omega t}$$

其中，$h(t)$为一实值之低通核函数，如短时距傅立叶变换中的视窗函数。此般分解下的各项系数则可定义如下：

$$\begin{align} \epsilon( t, \omega ) &= \int x(\tau) h( t - \tau ) e^{ -j \omega \left[ \tau - t \right]} d\tau \\ <pre> &= e^{ j \omega t} \int x(\tau) h( t - \tau ) e^{ -j \omega \tau } d\tau \\ &= e^{ j \omega t} X(t, \omega) \\ &= X_{t}(\omega) \\ &= M_{t}(\omega) e^{j \phi_{\tau}(\omega)}. </pre> \end{align}$$

其中，$X_{t}(\omega)$为讯号$x(t)$在时间位移为$t$、使用$h(t)$做为视窗函数下所得到的傅立叶变换；而$M_{t}(\omega)$与$\phi_{\tau}(\omega)$则分别代表$X_{t}(\omega)$的强度以及相位值。

$x(t)$可透过滑动视窗的系数以下式重建：

$$\begin{align} x(t) & = \iint X_{\tau}(\omega) h^{*}_{\omega}(\tau - t) d\omega d\tau \\ <pre> & = \iint X_{\tau}(\omega) h( \tau - t ) e^{ -j \omega \left[ \tau - t \right]} d\omega d\tau \\ &= \iint M_{\tau}(\omega) e^{j \phi_{\tau}(\omega)} h( \tau - t ) e^{ -j \omega \left[ \tau - t \right]} d\omega d\tau \\ &= \iint M_{\tau}(\omega) h( \tau - t ) e^{ j \left[ \phi_{\tau}(\omega) - \omega \tau+ \omega t \right] } d\omega d\tau. </pre> \end{align}$$

当一个讯号的强度频谱为$M(t,\omega)$，且其时间上较相位上的变化来得慢时，上述重建式中积分值最主要的贡献来源即满足下列稳相条件（phase stationarity condition）的$t,\omega$及其邻近区域：

$$\begin{align} \frac{\partial}{\partial \omega} \left[ \phi_{\tau}(\omega) - \omega \tau + \omega t\right] & = 0, \\ \frac{\partial}{\partial \tau} \left[ \phi_{\tau}(\omega) - \omega \tau + \omega t \right] & = 0. \end{align}$$

亦即满足下列定义的$\hat{t}, \hat{\omega}$及其邻近区域：

$$\begin{align} \hat{t}(\tau, \omega) & = \tau - \frac{\partial \phi_{\tau}(\omega)}{\partial \omega} = -\frac{\partial \phi(\tau, \omega)}{\partial \omega}, \\ \hat{\omega}(\tau, \omega) & = \frac{\partial \phi_{\tau}(\omega)}{\partial \tau} = \omega + \frac{\partial \phi(\tau, \omega)}{\partial \tau}. \end{align}$$

此现象在光学等领域中称为稳相近似：对周期或准周期函数而言，其傅立叶相位频谱在其振荡频率附近相对于时间的变化较慢，反之在外围区域则较快。以时域上能量较为集中的脉冲讯号为例，其相位频谱在脉冲时间点附近的变化较慢，而较远的区域则变化较快。

在讯号重建时，合成波形之正值与负值于相位快速变化的频率区域会因破坏性干涉而互相抵消。只有相位变化缓慢（即稳相）的区域会对重建的讯号有显著的影响，故贡献极大值（讯号的重心）即发生在相位随时间与频率变化最缓慢的地方。

由此知，可以直接从（一般在建构时频谱时被忽略的）短时距傅立叶转换之相位值求得时频座标，且其值即等于讯号的局部群延迟$\hat{t}_g(t,\omega)$和瞬时频率$\hat{\omega}_i(t,\omega)$。这些时频值代表的是经由视窗函数过滤后的讯号，因而属于局部（local）的量且无法表示整个讯号的全域（global）特征。

修正滑动视窗法（亦即重新分布法）便是将计算$\epsilon(t,\omega)$的时频点从原始的$t,\omega$重新移至局部极大贡献点$\hat{t}(t,\omega), \hat{\omega}(t,\omega)$。由于本质上这个局部极大点和质量分布中心点有相似之处，故也常类比成“重心”。该类比也有助于做出下述的理解：一个讯号的频谱能量分布对重心的贡献只在“有能量”的区域才具有意义，因此对时频谱上量值为零的座标点使用重新分布法是没有意义的。

## 快速时频重新分布法

在数字信号处理中，最常见的做法是先在时域与频域分别对讯号进行取样，并使用离散傅立叶变换计算时域讯号$x(n)$的频域变换$X(k)$。因离散的时频资料点无法直接计算偏微分，小寺邦彦等人所提出的重新分布法并不能直接套用于离散形式的短时距傅立叶变换；甚至有人指出，这项难处就是使得重新分布法无法得到更广泛应用的主因。

然而，偏微分的值可以利用有限差分法来估计。举例而言，在相位频谱上两个邻近的时间点（$t+\Delta t/ 2$、$t-\Delta t/ 2$）取值、并将两值相减后除以时间差$\Delta t$，便可得到该处频谱于时域的偏微分估计值（频域亦同理）：

$$\begin{align} \frac{\partial \phi(t, \omega)}{\partial t} & \approx \frac{1}{\Delta t} \left[ \phi \left (t + \frac{\Delta t}{2}, \omega \right ) - \phi \left (t - \frac{\Delta t}{2}, \omega \right ) \right] \\ \frac{\partial \phi(t, \omega)}{\partial \omega} & \approx \frac{1}{\Delta \omega}\left[ \phi \left (t, \omega+ \frac{\Delta \omega}{2} \right ) - \phi \left (t, \omega-\frac{\Delta \omega}{2} \right ) \right] \end{align}$$

在频谱上，某些区域的相位变化主要源自于附近某单一弦波成分的旋转；这些区域内的相位值为一线性函数，因此对于足够微小的$\Delta t$及$\Delta \omega$而言，上式的有限差分可为相位值的偏微分做一良好的近似。

道格拉斯·尼尔森（Douglas Nelson）也曾独立于小寺邦彦等人，开发了利用短时距相位频谱的偏微分来改善时频谱之分辨率的类似方法。稍加推导便可容易地证明他提出的“跨频谱平面”（cross-spectral surface）所计算之偏微分近似值和前述的有限差分法无异。

弗朗索瓦·俄歇（François Auger）及派翠克·弗朗德兰（Patrick Flandrin）则指出，小寺邦彦等人所提出的时频谱重新分布法可透过下述的一般化，延用于任一科恩系列的时频表示法：

$$\begin{align} \hat{t} (t,\omega) &= t - \frac{\iint \tau \cdot W_{x}(t-\tau,\omega -\nu) \cdot \Phi(\tau,\nu) d\tau d\nu} {\iint W_{x} \left (t-\tau,\omega -\nu \right ) \cdot \Phi (\tau,\nu) d\tau d\nu } \\ \hat{\omega} (t,\omega) & = \omega - \frac{\iint \nu \cdot W_{x}(t-\tau,\omega -\nu) \cdot \Phi(\tau,\nu) d\tau d\nu} {\iint W_{x}(t-\tau,\omega -\nu) \cdot \Phi(\tau,\nu) d\tau d\nu} \end{align}$$

其中，$W_{x}(t,\omega)$为$x(t)$的维格纳–魏尔分布，而$\Phi(t,\omega)$则是定义该分布的核函数。他们更提出了一项不须直接计算相位之偏微分，便可快速且准确地得到重新分布后之时频谱的方法。以时频谱而言，其重新分布法的计算如下：

$$\begin{align} <pre> \hat{t} (t,\omega) & = t - \Re \left \{ \frac{ X_{\mathcal{T}h}(t,\omega) \cdot X^*(t,\omega) }{ | X(t,\omega) |^2 } \right \} \\ \hat{\omega}(t,\omega) & = \omega + \Im \left \{ \frac{ X_{\mathcal{D}h}(t,\omega) \cdot X^*(t,\omega) }{ | X(t,\omega) |^2 } \right \} </pre> \end{align}$$

其中，$X(t,\omega)$是讯号以视窗函数$h(t)$所得到的短时距傅立叶变换、$X_{\mathcal{T}h}(t,\omega)$是讯号以时间加权之视窗函数$h_{\mathcal{T}}(t) = t \cdot h(t)$所得到的短时距傅立叶变换、$X_{\mathcal{D}h}(t,\omega)$则是讯号以对时间取导数之视窗函数$h_{\mathcal{D}}(t) = \tfrac{d}{dt}h(t)$所得到的短时距傅立叶变换；故知，利用辅助视窗函数$h_{\mathcal{T}}(t)$以及$h_{\mathcal{D}}(t)$，可得到含原始在内共三种傅立叶变换的结果。接着，将这三种变换分别在$t,\omega$处取值并做代数上的组合，便可在任一时频座标点$t,\omega$计算重新分布法的结果。由于上述的方法只是在短时距频谱上的某单一时频点取值，而不须做任何导数的计算，因此可视作离散短时距傅立叶变换之重新分布的快速算法。

上述方法有一限制，即$| X(t,\omega) |^2$之值必须为非零；不过，这项限制一般而言皆会被满足，因为只在频谱有能量的区域才有值可以“重新分布”，而方法本身在频谱值为零处也没有意义。

## 可分离性

<image>

<image>

短时距傅立叶变换经常用于估算一个多成分讯号（例如准谐波的乐器音）中个别单一成分的振幅和相位，而时频重新分布则能将变换后的时频点集中至能量分布的局部重心，使频谱更加清晰。

对单一成分的讯号而言，其瞬时频率可由该讯号经任意短时距傅立叶变换后得到的相位之偏微分估算。若讯号可拆解为多个成分：

$$x(t) = \sum_{n} A_{n}(t) e^{j \theta_{n}(t)}$$

且定义各成分的瞬时频率为其相位对时间的偏微分，即：

$$\omega_{n}(t) = \frac{d \theta_{n}(t)}{d t}$$

则可以借由设计适当的带通滤波器使其通带只包含某一特定讯号成分，再从该滤波器经变换后之频谱的相位值，计算该讯号成分的瞬时频率。

尼尔森称这项频域中的条件为“可分离性”（separability），且所有透过此方法分析的讯号皆须满足该条件。若该条件并未满足，则因不能透过短时距傅立叶变换估算个别成分的参数，而无法得到欲求的多成分讯号拆解；在这样的情况下，必须采用其它的分析视窗函数，以满足可分离性的条件。

若某特定的短时距频谱分析视窗可用于拆解讯号中不同频率的成分，则每个滤波器经短时距傅立叶变换后的脉冲响应便是（至多一项）单一主要讯号成分（即能量最大）经滤波后的结果，故$X(t,\omega_0)$之相位值对时间的偏微分值，即相等于频率$\omega_0$下主要成分之相位值对时间的偏微分值。因此，若某瞬时频率为$\omega_{n}(t)$的讯号成分$x_n(t)$ 在$\omega_0$附近为主要成分时，则该成分的瞬时频率可透过原始讯号经短时距傅立叶变换后位于$\omega_0$的相位值求得，即：

$$\begin{align} \omega_{n}(t) &= \frac{\partial}{\partial t} \arg\{ x_{n}(t) \} \\ <pre> &= \frac{\partial }{\partial t} \arg\{ X(t,\omega_{0}) \} </pre> \end{align}$$

频域的可分离性要求每个带通滤波器经短时距傅立叶变换后，最多只可通过一项讯号成分，而在时域同样也有类似的限制：讯号中两个不同时间的成分（即事件）彼此间分隔须够大，使得任一时域视窗下的分段，都不得包含超过一个事件。此即时域下的可分离性条件，且等义于任两事件间的间隔须大于短时距傅立叶变换滤波器之脉冲响应长度，也就是$h(t)$中非零取样点所涵盖的范围。

一般而言，任一多成分讯号有无限多种相同合理的拆解方式；不过，可分离性必须是在所欲求得的拆解下进行考虑。举例而言，在分析一个语音讯号时，若使用时间范围较声门脉冲间隔长的分析视窗，则足以将不同谐波分离；但个别的声门脉冲则会因单一视窗包含了多个脉冲而被模糊，也就是所选的分析视窗使得讯号在时域中不具可分离性。反之，若使用时间范围远小于声门脉冲间隔的分析视窗，则能够成功分离每个脉冲，因为每个视窗中最多只包含一个脉冲；但不同谐波便会彼此模糊交叠，因为视窗在频域中的范围比单一谐波还要宽，也就是所选的分析视窗使得讯号在频域中不具可分离性。

## 外部链接

- TFTB — Time-Frequency ToolBox
- SPEAR - Sinusoidal Partial Editing Analysis and Resynthesis
- Loris - Open-source software for sound modeling and morphing
- SRA - A web-based research tool for spectral and roughness analysis of sound signals
- Sparse time-frequency representations - PNAS