基于二次EEMD的Wigner-Ville分布旋转机械故障信号分析及试验研究

郭 奇，刘卜瑜，史立波，李 波

(1.燕山大学 环境与化学工程学院，秦皇岛;2.中国石油吉林石化公司精细化学品厂，吉林)

经过100多年的发展，Fourier变换不但已成为一个重要的数学分支，而且还是信号分析与处理的重要工具。然而，Fourier变换也存在不足，即缺乏时间和频率的定位功能，在分辨率上的局限性以及对非平稳信号分析的局限性［1］。工程上旋转设备故障信号大多是非平稳信号，故障信号分析的关键便成为如何对采集到的非平稳信号进行分析。因此大量的时频分析方法应运而生。其中，由Ville在1947年提出的Wigner-Ville分布凭借其优良的时变特性和较高的时频分辨率一直被应用于信号的分析和处理领域。然而，由于其不满足可加性而引起的交叉项却成了不可忽视的缺陷，并使真实频率成分之间产生虚假频率成分，严重干扰对信号的准确分析。在一定程度上限制了Wigner-Ville分布在实际中的发展应用。

从交叉项产生的原因出发，可得到如下启发，首先将叠加的信号分解成单个分量，分别对各分量进行WVD计算，然后对计算结果求和，即可避免WVD在计算过程中不满足可加性的缺陷。对信号进行单分量分解有多种方法，其中近年来发展起来的经验模态分解(Empirical Mode Decomposition，EMD)和聚合经验模态分解(Ensemble EmpiricalMode Decomposition，EEMD)方法很受青睐。然而，这两种方法也并非完美，EMD方法在分解过程中容易出现频率混叠，影响分解的效果，不能得到真正的单一频率成分，而EEMD分解虽然较EMD方法在抗混方面有所改善，但在消除Wigner-Ville分布交叉项方面还是不太理想。对此，本文提出利用二次EEMD消除Wigner-Ville分布交叉项方法，先用EEMD将信号分解成若干近似单一频率的分量，对这些分量再次用EEMD方法进行二次分解，得到效果较好的单一频率成分，最后对新得到的单一分量分别进行Wigner－Vill时频分析，并对分析的结果求和，这样便很好地消除了交叉项的影响。

1 Wigner-Ville分布

Wigner-Ville分布是一种二次型分布，由信号x(t)本身或其频谱X(ω)定义为如下两种等价的方式:

其中:* 表示复数共轭。

式(1)中的x*(t－τ/2)x(t+τ/2)称为信号的瞬时相关函数，因此Wigner-Ville分布实质上是对信号的瞬时相关函数的Fourier变换，其结果能够反映信号的时频特征。

WVD有很多优点，但却不满足可加性，考虑如下信号:

将其代入式(1)，信号x(t)的WVD可写为:

其中:

WVDx1x2(t，ω)与 WVDx2x1(t，ω)称为互 Wigner-Ville分布，是复值的，并且可看出:

因此，WVDx1x2(t，ω)+WVDx2x1(t，ω)是实值的。式(4)可简写为:

由此可看出，两个信号和的Wigner-Ville分布并不是简单的两个信号各自的Wigner-Ville分布之和，附加项2Re［WVDx1x2(t，ω)］通常称为交叉项［1］。

2 EEMD基本原理

经验模态分解(Empirical Mode Decomposition，EMD)方法是由Huang等［3］提出的一种直观的、直接的、后验的和自适应的分解方法。该方法通过一定的筛选条件将信号分解成为若干个相对平稳的、互不相关的本征模函数(Intrinsic Mode Function，IMF)。但当信号中存在异常事件(如间歇性成分、脉冲干扰等)时，分解得到的IMF就会存在模态混叠现象。由文献［3］的研究成果知，模态混叠的产生与极值点选择有关，利用极值点信息通过三次样条插值生成的包络为异常事件的局部包络与真实信号包络的组合，经该包络计算出均值并筛选出的IMF包含了信号固有模式和与异常事件模式，即模态混叠现象［2］。为改善该现象，Wu等［4］提出了 EEMD方法，通过加入高斯白噪声解决EMD分解过程中出现的模态混叠问题。EEMD的计算方法为:

(1)初始化所要分解的IMF个数N和总的分解次数M(第一轮分解取m=1);

(2)对所输入的信号加入随机正态分布的白噪声序列，并对被加入白噪声的信号做归一化处理;

(3)应用EMD对归一化的信号进行分解，得到N个IMF;

(4)当m＜M时，重复(2)、(3)步骤，并使m=m+1;

(5)将分解得到的M组IMFs对应集成平均，当M足够大时，对应剔除添加的白噪声序列，而得到输入信号的本征模函数组合，见式(9)、式(10):

其中:r(t)为余项［3－5］。

3 利用二次EEMD消除交叉项

WVD的交叉项是因被分析信号中含有多个频率成分在计算过程中不满足可加性而造成的。如果将被分析信号分解成为一系列单一频率成分的信号，对其信号分别进行WVD计算，再将计算结果求和，便可消除交叉项的影响［6］。

3.1 二次EEMD思想基础

EMD分解IMF的过程是基于Huang的理想假设，但在实际信号分解过程中还存在比较棘手的问题。试验证明，实际多分量信号经EMD分解时容易产生虚假的IMF或者分解不出IMF(即模态混叠现象)［7］，因而分解出的IMF只能被称为近似单频率成分，对此文献［7－8］均给出了详细的论证。因模态混叠的存在，Wu等［5］提出了EEMD的思想以改善其分解效果。由于受求取信号包络和筛选条件的限制，EEMD方法只是起到了改善分解效果的作用。即，经过EEMD分解出的“单频”成分还可能存在模态混叠(“单频”成分中仍然可能存在两个或两个以上的频率成分，在3.2节仿真结果中验证)。因此本文认为有必要对首次分解出的IMFs进行第二轮的EEMD分解。由于EEMD分解过程依靠随机加入均匀分布的白噪声来改善其分解效果，并且信噪比可调，使二次分解后得到效果较好的IMFs成为可能。另外，由EEMD计算程序可知EEMD分解过程中强制规定了所要分解IMFs的数量，即TNM=fix(log2(N))－1(TNM为IMF数量，N为采样点数)，这对有些信号成分可能会造成过分解或者分解不完全。基于此，在第二轮分解结束以后，本文对分解得到的IMFs以原信号为基准进行相关系数的筛选［1，9］，该方法可以有效地剔除由采样率不足、样条插值以及过分解等引起的伪分量，从而得到质量较高的单频率成分。至此可得出二次EEMD消除交叉项算法如下:

(1)用EEMD方法将被分析信号分解成为一系列满足式(10)的单一频率成分;

(2)利用EEMD方法对步骤(1)中所得到的IMF进行二次分解，从而获得新的IMFs;

(3)对步骤(2)中所得到的IMFs以原信号为基准进行相关系数的筛选以获得质量较高的单一频率成分;

(4)对步骤(3)中所得到的一系列IMF分别进行WVD计算:

(5)将各分量的WVD计算结果求和即得计算结果:

3.2 仿真信号模拟验证

合成仿真信号:

可以看出信号中包含频率分别为10 Hz的正弦信号，20 Hz的余弦信号，40 Hz的正弦信号。信号时域波形如图1所示。图2是经过EEMD分解得到的IMFs，其数量由EEMD的计算程序强制规定(TNM=fix(log2(N))－1，N为采样点数)，第一条曲线和最后一条曲线分别为原始信号和余量。图3为再次利用EEMD方法对图2中第一个IMF分量分解后经过相关系数法筛选出的IMFs，第一条和最后一条曲线分别为原始信号和余量。图3的结果证实了经过EEMD方法分解后的部分IMF是可以再分的，同时也证明了二次EEMD方法的可行性。图4是对信号直接做WVD分析后的结果，从图中看出，除在10 Hz、20 Hz及40 Hz处有真实的频率成分外，在15 Hz、25 Hz及30 Hz附近也有频率成分。而这些频率成分在信号中并不存在，此即为WVD计算过程中出现的交叉项。为消除交叉项，采用3.1节中方法对信号进行处理，图5为一次EEMD处理后的结果，可以看出25 Hz处的虚假成分已被消除，但在15 Hz与30 Hz处虚假成分依然存在;图6图7为二次EEMD处理后的灰度图和三维图，在图中可清楚地看到大部分虚假成分已被消除，用二次EEMD方法消除WVD交叉项效果很好。

图4 未经EEMD分解的Wigner-Ville分布Fig.4 The WVD of the simulated signal not based on EEMD

图5 经过EEMD分解的Wigner-Ville分布Fig.5 The WVD of the simulated signal based on EEMD

图6 经过二次EEMD分解的Wigner-Ville分布Fig.6 The 2 － D WVD of the simulated signal based on dual EEMD

图7 经过二次EEMD分解的Wigner-Ville分布三维图Fig.7 The 3 － D WVD of the simulated signal based on dual EEMD

4 碰摩故障试验研究

本次试验采用DGZH－A型多功能转子实验台模拟汽轮机转子碰摩故障。将试验台左侧第一个轮盘施加适当的载荷，使转子旋转过程中产生不平衡力，从而在轴承处产生摩擦，如图8，转子转速设定为3000 r/min，采样频率为1280 Hz，采样点数为1024，使用加速度型传感器采集数据。

图8 DGZH－A型多功能转子实验台Fig.8 DGZH-A Multi-function rotor experiment bench

由动静件局部碰摩的故障特征可知，局部摩擦引起的振动频率中包含有不平衡引起的转速频率，即ω=3000/60=50 Hz，同时由于摩擦振动是非线性振动，所以还会包含有2ω、3ω、…一些高次谐波［10］。

为验证二次EEMD方法分析的准确度，采用京航公司HG8904C多通道数据采集故障诊断系统进行故障数据采集与分析，图9为HG8904C分析结果的时域图和频域图的截图。可以看出专家系统提供的频域图中出现了 50 Hz、100 Hz、150 Hz、200 Hz 的频率，此为动静件碰摩故障特征频率，可确定试验台此时存在动静件碰摩故障。

图9 HG8904C分析结果的时域图和频域图Fig.9 The time-domain chart and frequency-domain chart of the fault signal givened by HG8904C

图10为用Matlab还原的采集信号的时域图，由图11～图14可以看出采用二次EEMD方法处理过的实际信号的Wigner-Ville分布中的交叉项已基本被消除，体现故障特征的转速频率及其倍频(50 Hz、100 Hz、150 Hz、200 Hz)清晰可见，此结果与专家系统分析结果基本一致。由此证明了利用二次EEMD消除Wigner-Ville分布交叉项的方法是可行的。

图10 信号的时域图Fig.10 The time-domain chart of the fault signal

图11 原始故障信号的Wgner-Ville分布灰度图Fig.11 The WVD of the original fault signal

图12 经过EEMD分解的Wgner-Ville分布灰度图Fig.12 The WVD of the fault signal based on EEMD

图13 经过二次EEMD分解的Wgner-Ville分布灰度图Fig.13 The 2 － D WVD of the fault signal based on dual EEMD

图14 经过二次EEMD分解的Wgner-Ville分布的三维图Fig.14 The 3 － D WVD of the fault signal based on dual EEMD

5 与EMD方法消除WVD交叉项比较

为克服WVD交叉项的干扰，对第4节中采集的数据利用EMD消除交叉项的方法进行分析，见图15，分解后经过相关系数法筛选得到的IMFs(该方法可以得到质量较好的IMFs，见文献［9］)，第一条曲线和最后一条曲线分别为原始信号和余量。从图中可以看出EMD方法成功分解出四个单一频率成分，与故障特征频率的数量相吻合。图16、图17，分别为经EMD分解后WVD的灰度图和三维图，由图看出虽然50 Hz、100 Hz、150 Hz、200 Hz的频率成分清晰可见，但在每单一频率成分之间出现许多不连续的干扰，由文献［7－8］可知，此为EMD方法分解多分量信号时产生的虚假频率成分。通过图15、图16与图13、图14对比，可以发现利用二次EEMD方法分解比利用EMD方法分解效果更好。

图15 EMD分解后经过筛选的IMFsFig.15 The IMFs given by EMD and correlation coefficient

6 结论

本文针对WVD分析时产生的交叉项问题，提出利用二次EEMD消除Wigner-Ville分布交叉项的方法。该方法在有效地消除WVD分布交叉项的同时保留了WVD分布的各种优良特性。随后结合汽轮机动静件摩擦故障的实例，用二次EEMD消除Wigner-Ville分布交叉项的方法对故障的振动信号进行分析，获得了故障信号真实的频率成分并对故障进行诊断，所得结果与HG8904C多通道数据采集故障诊断系统所得到的分析结果一致。与利用EMD消除WVD交叉项的方法对比二次EEMD消除WVD交叉项的方法可获得更好的效果，从而证明了该方法在实际应用中的可行性。

［1］褚福磊，彭志科，冯志鹏，等.机械故障诊断中的现代信号处理方法［M］.北京:科学出版社，2009.

［2］曹冲锋，杨世锡，杨将新.大型旋转机械非平稳振动信号的EEMD降噪方法［J］.振动与冲击，2009，28(9):33－38.

［3］Huang N E，Shen Z，Long S R，et al.The empiricalmode decomposition and the hilbert spectrum for nonlinear and nonstartionary time series analysis［J］. Proc. Roy. Soc.London，1998，1971(454):903 －995.

［4］Wu Z H，Huang N E.A study of the characteris－tics of white noise using the empiricalmode decomposition method［J］.Proc.Roy.Soc.London，2005，l460A:1579 －1611.

［ 5］Wu Z H，Huang N E.Ensemble empirical mode decomposition:a noise-assisted data analysismethod［J］.Advances in Adaptive Data Analysis，2008，11:1 －4.

［6］韩中合，朱霄珣，李文华，等.基于EMD消除Wigner－vill分布交叉项的研究［J］.汽轮机技术，2010，52(3):211－214.

［7］沈国际，陶利民，陈仲生.多频信号经验模态分解的理论研究及应用［J］.振动工程学报，2005，18(1):91－94.

［8］全学海，丁宣浩，蒋英春.基于经验模态分解的筛选条件研究［J］.南阳师范学院学报，2010，9(3):22－26.

［9］林 丽，余 轮.基于相关系数的EMD改进算法［J］.计算机与数字工程，2008，36(12):28－29.

［10］盛兆顺，尹琦岭.设备状态监测与故障诊断技术及应用［M］.北京:化学工业出版社，2009.

［11］齐 天，裘 焱，吴亚锋.利用聚合经验模态分解抑制振动信号中的模态混叠［J］.噪声与振动控制，2010，2:103－106.

［12］张 超，陈建军.EEMD方法和EMD方法抗模态混叠对比研究［J］.振动与冲击，2010，29(S):87－90.