现代信号处理教程 - 胡广书(清华)

81为了看清图3.3.4中交叉项的行为，我们将该图作了旋转，因此，水平方向为频率，垂直方向为时间。

图3.3.3 例3.3.3的WVD 图3.3.4 例3.3.4的WVD例3.3.5 令()2142t x t e ααπ-⎛⎫= ⎪⎝⎭（3.3.5）可求出其WVD 为()22,2exp[]x W t t ααΩ=--Ω （3.3.6）这是一个二维的高斯函数，，且()Ω,t W x 是恒正的，如图3.3.5所示。

由该图可以看出，该高斯信号的WVD 的中心在()()0,0,=Ωt 处，峰值为2。

参数α控制了WVD 在时间和频率方向上的扩展。

α越大，在时域扩展越小，而在频域扩展越大，反之亦然。

其WVD 的等高线为一椭圆。

当WVD 由峰值降到1-e 时，该椭圆的面积π=A 。

它反映了时－频平面上的分辨率。

如果令 ()2142t h t e ααπ-⎛⎫=⎪⎝⎭，()2142t x t e ββπ-⎛⎫= ⎪⎝⎭，则()t x 的谱图 ()⎥⎦⎤⎢⎣⎡Ω+-+-+=Ω2221exp 2,βαβααββααβt t STFT x （3.3.7）82图3.3.5 例3.3.5的WVD,（a ）高斯信号，（b ）高斯信号的WVD它也是时－频平面上的高斯函数。

当其峰值降到1-e 时，椭圆面积π2=A 。

这一结果说明，WVD 比STFT 有着更好的时－频分辨率。

如果令 ()()tj et t x t x 001Ω-= （3.3.8）式中()t x 是（3.3.5）式的高斯函数。

()t x 1是()t x 的时移加调制，其WVD 是：()12200,2exp[()()/]x W t t t ααΩ=---Ω-Ω （3.3.9）它将（3.3.6）式的()Ω,t W x 由()()0,0,=Ωt 移至()()00,,Ω=Ωt t 处。

其WVD 图形请读者自己画出。

例3.3.6 令 ()2201422j t t j t z t ee e αβαπΩ-⎛⎫= ⎪⎝⎭（3.3.10） 它是由（3.3.5）式的()t x 与83()202j tj t y t Aee βΩ= （3.3.11）相乘而得到的（在（3.3.9）式中，A=1）。

我们已熟知()t y 为线性调频Chirp 信号，其WVD 是：()()20,2y W t A t πδβΩ=Ω-Ω- （3.3.12）其图形我们已在图1.1.2中给出。

（3.3.9）式()t z 的WVD 是()220,2exp[()/]z W t t t αβαΩ=--Ω-Ω- （3.3.13）它是恒正的。

显然，高斯信号、Chirp 信号都是（3.3.9）式()t z 的特例。

()Ω,t W z 如图3.3.6所示。

图3.3.6 例3.3.6的WVD ，（a ）Chirp 信号，（b ）Chirp 信号的WVD例3．3．7 令()x t 为一多普勒信号。

所谓多普勒信号指的是一个物体相对一个位置不变的“观察者（如雷达）”运动时，“观察者”所听到或所记录到的该物体运动的信号，如其运动84的速度或发出的声音。

众所周知，当该运动物体接近和远离“观察者”时，其信号当频率会发生变化。

图3.3.7 给出了该信号当时域波形、频谱及时－频分布。

由该图可看出信号的能量随时间和频率的分布。

图3.3.6 例3.3.6的WVD3.4 Wigner 分布的实现如同许许多多的其它信号处理的算法一样，我们最终的目的是要将它们应用于科研或工程的实际。

这时所遇到的问题同样是信号的离散化及数据的有限长问题。

在（3.1.2）式中，若令对信号()t x 的抽样间隔为s T ，即s nT t =，并令s kT =2τ，则s kT 2=τ，这样，（3.1.2）式对τ的积分变成对k 的求和，即 ()()()∑∞-∞=Ω-*-+=Ωk T k j sssssx sekT nT x kT nT x T t W 22, （3.4.1）若将s T 归一化为1，并考虑到相对离散信号的频率s T Ω=ω［19］，则上式变为：()()()∑∞-∞=-*-+=k k j x ek n x k n x t W ωω22, （3.4.2）85我们知道，将()t x 变成()n x ，则()t x 的频谱()Ωj X 将变成周期的频谱()ωj e X ，周期为π2，且π2对应的抽样频率为s f 。

与此同样的是，()t x 的WVD ()Ω,t W x 也变成周期的()ω,n W x ，但是，由（3.4.2）式，()ω,n W x 的周期为π，即： ()()()()∑∞-∞=+-*-+=+k k j x ek n x k n x n W πωπω22, （3.4.3）众所周知，若()t x 的最高频率为max f ，那么，抽样频率至少应满足max 2f f s ≥。

这是由抽样定理所决定的。

如若按max 2f f s =对()t x 抽样，那么用抽样后的()n x 做WVD ，由于其周期变为π，因此在WVD 中必将产生严重的混迭。

解决这一问题的直接方案是提高抽样频率，要求s f 至少要满足max 4f f s ≥ （3.4.4） 但是，一旦信号()t x 由max 2f f s =抽样变成()n x 后，要想对()t x 重新抽样是困难的。

解决该问题的较为简便的方法有两个：１、采用解析信号，由解析信号的性质可知，将()t x 作Hilbert 变换得到()t x ˆ，按()()()t xj t x t z ˆ+=构成解析信号。

()t z 只包含()t x 的正频率部分。

这样，既可减轻由正、负频率分量所引起的交叉项干扰，又可在保持原有抽样频率max 2f f s =的情况下，避免了频域的混迭；２、对()n x 作插值，人为地将其抽样频率s f 提高。

具体办法是：若想将抽样频率s f 提高一倍，则可将()n x 每两点之间插入一个零，然后再让该信号通过一低通数字滤波器，从而将插入的这些零值变成原信号相应点的插值，有关插值的原理，详见本书第5章，此处不在讨论。

对（3.4..3）式，现余下两个问题要解决。

一是频率ω仍需离散化，二是式中对k 的求和需要取有限长。

式中k 是信号x 的时间序号，n 代表时移。

现令()()()k n x k n x k n r x -+=*, （3.4.5）86并假定()k x 的长度为N ，即1,,1,0-=N k ，现分析一下()k n r x ,的取值情况。

如图3.4.1（a ）所示，()1,,1,0-=N k x ，将()k x 翻转得()k x -，现将()k x 、()k x -分别向左和向右移动n 个时刻，如取2=n ，6=N ，则如图3.4.1（b ）和（c ）所示。

（a ）、()k x（b ）、()k x +2（c ）、()k x -2图3.4.1 ()k n r x ,的解释当6=N 时，读者不难写出：(){}(){}(){}(){}(){}(){}⎪⎪⎪⎪⎭⎪⎪⎪⎪⎬⎫===-===-===-===-======******************０ 时， ，，时， ，，，，３时， ，，，，时， ，，时， 时，１５２４３３４２５１k x x k r n k x x x x x x k r n k x x x x x x x x x x k r n k x x x x x x x x x x k r n k x x x x x x k r n k x x k r n x x x x x x 55354453041322314002112000,551~1,442~2,32~2,221~1,110,00（3.4.6） 假定将()k n r x ,都扩充成N 点序列，即在其后补零，那么，（3.4.2）式可写成 ()()∑-=kkl Njx x ek n r l n W π4,2, （3.4.7）式中，()k n r x ,在任一n 时刻的长度都变成了N ，k 的实际取值情况由（3.4.4）决定。

这样，对（3.1.6）式稍作变化后即可用DFT 来实现。

上述过程即是离散WVD 的思路。

具体实现方k⋅⋅⋅⋅)(k x 012345k⋅⋅⋅⋅)2(k x +01232-1-k⋅⋅⋅⋅)2(k x -0121-2-3-87案见MA TLAB 中High －Order Spectral Analysis Toolbox 中的Wig2.m 文件。

注意：()l n W x ,在频域是周期的，周期为π（归一化频率为0.5）。

该方法有明显的缺点，即在不同的n 下，计算()k n r x ,时所利用的()k x 的点数有着明显的不同（见（3.4.6）式）。

此外，由于WVD 是二次函数的分布，有交叉项存在，因此，希望能对这种交叉项有所抑制。

针对这两个原因，人们自然提出了“加窗WVD ”，即“伪WVD （Pseudo WVD ，PWVD ）”。

取窗函数()n w ，()n w 应是实对称的函数，假定其宽度为14-L ，即当L k 2≥时，()0=k w ，用()k w 乘()k n r x ,，由（3.4.5）式得：()()()() k n x k n x k w k n pr x -+=*, (3.4.8) 将其代入（3.4.2）式，有()()()()()()()()()()()()()()()()()()() n x n x e k n x k n x k w n x n x e k n x k n x k w e k n x k n x k w e k n x k n x k w k n PW L k k j L k k j L k kj L L k kj x *-=-**-=-*--=-*---=-*-⎭⎬⎫⎩⎨⎧-+=--++-+=-+=∑∑∑∑22Re 42222,12021202012212122ωωωω现将ω离散化，可将π2分成L 2等份，即l Ll 22πω=，则（3.4.9）式变为： ()()()()()() n x n x e k n x k n x k w l n PW L k L kl j x *-=⎪⎭⎫⎝⎛-*-⎪⎭⎪⎬⎫⎪⎩⎪⎨⎧-+=∑2Re 4,12024π (3.4.10) 式中12~0-=L l 。

注意, 由该式，有()() iL l n Pw l n PW x x +=,, (3.4.11) 即()l n PW x ,以L 为周期。

这样，若按（3.4.10）式计算L 2点FFT ，则求出的()l n PW x ,将有一半的冗余。

通常，我们假定：()() ，当L k k n x k n x >=-+*0 (3.4.12)这样，（3.4.10）式可变成88()()()()()() n x n x e k n x k n x k w l n PW L k kl L j x *-=⎪⎭⎫⎝⎛-*-⎪⎭⎪⎬⎫⎪⎩⎪⎨⎧-+=∑2Re 4,102π (3.4.13) 这是一个标准的L 点的FFT 。