龚廷恺　袁晓辉　王细洋

1.华中科技大学，武汉，430074　　2.南昌航空大学，南昌，330063



最小熵反褶积的数学形态法在滚动轴承故障特征提取中的应用

龚廷恺1，2袁晓辉1王细洋2

1.华中科技大学，武汉，4300742.南昌航空大学，南昌，330063

针对强噪声背景下滚动轴承故障特征提取，提出了基于最小熵反褶积的数学形态法。该方法先应用最小熵反褶积算法加强信号中的冲击特性，再利用数学形态法进行故障特征提取，其中选取具有双向脉冲提取能力的DIF滤波器作为形态算子，并以峭度值作为结构元素长度选取依据。仿真信号和滚动轴承的内外故障实例分析表明该方法具有较好的特征提取效果。通过对比发现:最小熵反褶积算法能够增大信号中峭度值，有效加强信号脉冲特性。

最小熵反褶积；数学形态；故障诊断；滚动轴承

0　引言

滚动轴承是旋转机械中最常用的通用零部件之一，同时也是最易损的零件。根据统计，旋转机械的故障中有30%是由轴承故障引起的，因此它的工作状态对机器运行具有巨大的影响[1]。当故障点出现在轴承的内外圈及滚动体上时，周期性脉冲会产生。而且这些脉冲信号的频宽范围常常覆盖了轴承和机械系统固有振动频率，易引起轴承甚至机械系统的共振，导致机器振动和产生噪声，甚至损坏。因此这些冲击脉冲的成功提取成为轴承故障诊断的一种有效方法。

对于这些冲击脉冲性的振动信号分析有很多方法，如局部特征分解、经验模态分解和小波变换等[2-6]。相比前述方法，一种非线性的数学形态技术近来得到了广泛关注。文献[7]应用形态闭算子作为形态算子，以0.6倍的脉冲周期作为结构元素长度，成功提取出了轴承内外圈故障。但这种方法在选取结构元素长度时需要依赖经验知识。为了克服这种缺陷，Dong等[8]提出了一种基于信噪比的改进形态分析方法。然而，在强噪声背景下形态方法所能提取出的故障形态特征，即脉冲冲击特性是非常有限的。因此在形态方法处理前，希望通过一种方法能够简化信号，突出信号中的冲击特性，这样形态方法就能提取出更多的故障形态特征，便于故障诊断。因此，基于上述考虑本文提出了一种基于最小熵反褶积的数学形态法(mathematical morphology based on minimum entropy deconvolution, MMBMED)，并应用于故障轴承的振动信号特征提取中。该方法通过最小熵反褶积(minimum eneropy deconvolution, MED)方法进行预处理，而后应用数学形态方法提取故障特征，选取具有双向脉冲提取能力的DIF滤波器作为形态算子，结构元素长度的确定通过峭度作为选取原则。仿真信号和滚动轴承内外圈故障分析表明，基于最小熵反褶积的数学形态方法具有较好的故障特征提取效果。

1　相关理论

1.1数学形态法

数学形态学的基本思想是用结构简单、紧凑的结构元素去量度和提取图像中的对应形状，以达到对图像分析和识别的目的。形态学的非线性特性可以简化图像，在保持图像本质基础上除去噪声部分，最初是作为图像处理工具的，后来才发展到一维信号处理。

数学形态学中包含4个基本形态运算算子,即膨胀算子、腐蚀算子、开算子和关算子。它们的定义分别如下：

设f(n)为一维多值信号，s(n)为结构元素。H、K分别为f(n)和s(n)的定义域，H=(0,1，…，n-1)，K=(0，1，…，m-1)，且m

(f⊕s)(n)=min(f(n-m)+s(n))

(1)

m∈{0,1,…，M-1}

(fΘs)(n)=min(f(n+m)-s(n))

(2)

m∈{0,1,…M-1}

根据上述定义，形态开运算和形态闭运算的定义分别为

(3)

(f•s)(n)=(f⊕sΘs)(n)

(4)

Nikolaou等[7]通过分析发现4种形态算子对脉冲信号的作用不同。如腐蚀运算对负脉冲具有平滑作用，会减少脉冲信号中的正脉冲部分。膨胀运算具有和腐蚀运算相反的效果。闭算子能够保留信号中的正脉冲信号，开算子能够保留信号中的负脉冲信号。显然，实际应用中很难得到正负脉冲的先验知识。因此，常常用4个基本算子的算术组合作为新的形态算子。例如，形态闭和开的差分运算得到的DIF(differencefilter)算子，其数学定义如下式所示：

(5)

该算子具有同时提取正负脉冲信号的特性，因此具有较好的脉冲提取能力。

1.2最小熵反褶积算法

Wiggins[9]首次提出了MED的概念，并应用在地震信号处理当中。求取MED是指寻求一个反转滤波器的最优系数集合，使得所恢复的输出信号具有最大峭度。因为峭度是一个对故障敏感的因子，反映了信号中的峰值(即脉冲属性)，增大信号峭度值有助于轴承中的故障诊断。

图1说明了MED算法原理。信号e(n)通过结构滤波器后的输出和噪声v(n)共同作为一个反转滤波器 f 的输入，输出为y(n)。为了处理后有更好的效果，我们希望反转滤波器f的输出y(n)与最初的输入信号e(n)尽可能相似。虽然系统输入e(n)是未知的，但希望这些就是脉冲序列。Endo等[10]总结了MED方法具体过程：

图1　MED算法原理

(1)最小熵算法为长度L的反转滤波器 f找到一个系数集合，优化输出信号y(n)的峭度：

(6)

其中，[·]表示信号为离散信号。

(2)反转滤波器的卷积表达式为

(7)

(3)当目标函数(式(6))出现最大峭度值时，即

(8)

(4)对式(7)两边求导可得

(9)

(5)由式(7)和式(9)可将式(8)写成：

(10)

式(10)写成矩阵形式如下：

b=Af

(11)

式中，A为输出信号y(n)的自相关矩阵；b为滤波器的输出立方与输入信号的互相关加权和；f为反转滤波器系数矩阵。

由式(11)通过迭代计算可以计算出反转滤波器系数：

f=A-1b

(12)

由以上分析可以归纳出最小熵递归迭代算法步骤如下：

(1)初始化反转滤波器的系数，设定初始状态f(0)中元素全为1；

(2)通过输出信号x(0)和反转滤波器初始值f(0)(按式(7))计算y(0)；

(3)通过y(0)和x(0)(按式(11))计算列矢量b(1)(b(1)为第一次计算b的结果)；

(4)计算得到滤波器系数矩阵f(1)=A-1b(1)；

(5)计算预设定误差。如果误差小于设定值，结束计算。否则，返回步骤(3)。

2　信号仿真

仿真信号由谐波信号、脉冲信号和噪声三个信号复合构成，其数学表达式如下：

x(t)=x1+x2+x3

(13)

其中，谐波信号x1=cos(60t)+cos(130t)；x2为周期为32 Hz的指数衰减序列，该信号在一个周期内的衰减函数为e-10tsin(10t)；x3为均值为0，标准偏差为1的高斯白噪声。采样频率fs为2048 Hz,采样时间为1 s。图2为仿真信号的时域图和频谱图。在图2b中脉冲特征频率32 Hz、64 Hz等成分由于谐波信号和噪声抑制，没有被提取出来，因此采用建议方法进行分析。新方法处理过程为：①先采用MED方法预处理；②通过数学形态法提取特征，以式(5)中的DIF算子作为形态算子，选择扁平型高度为0的结构元素。对于其长度，本文以脉冲重复周期T的0.1倍为间隔，构造9个不同长度的结构元素，即 {0.1T，0.2T，0.3T，0.4T，0.5T，0.6T，0.7T，0.8T，0.9T}，分别计算9个不同长度的峭度值，最后选取最大峭度的长度作为结构元素长度。

(a)时域图

(b)频谱图图2　仿真信号

仿真信号经过MED处理后，继续进行数学形态分析。计算9个不同长度的结构元素，发现长度为0.1T时的峭度值7.7879为最大，因此长度选择0.1T，相应结果的频谱图见图3a。从图3a中发现脉冲特征频率32 Hz，以及它的谐波成分(64～160 Hz) 都被提取出来了，谐波成分和噪声得到了较好的抑制。特征提取效果说明最小熵反褶积和数学形态法相结合能够有效提取出脉冲特征成分，较好地抑制噪声。

为了说明MED方法的作用，仿真信号将直接运用新方法的第二步进行处理与对比。计算发现在0.1T时的峭度值只有4.2503，其结果频谱图见图3b。图中只发现了几个脉冲特征频率(64 Hz,96 Hz,128 Hz和160 Hz)。与图3a对比可知：仿真信号未经过MED预处理，结构元素长度在0.1T时峭度值从4.2503增大到7.7879，峭度值明显增大。图3b虽然提取出了几个脉冲特征频率，但它们周围的噪声较大，而且脉冲特征频率32 Hz没有提取出来。因此通过这个对比说明MED算法能够优化信号的峭度，加强信号的脉冲特性。

(a)有MED处理

(b)无MED处理图3　仿真信号频谱图对比

3　实例分析

3.1滚动轴承外圈故障分析

实例分析的数据来自于美国凯西西楚大学轴承试验中心。实验中外圈故障测试轴承为SKF6205-2RSJEMSKF深沟球轴承。采样频率为12kHz，采样点数为12288。转速为1750r/min，由计算可得外圈故障通过频率(theballpassingfrequencyouterrace,BPFO)约为105.1Hz。外圈故障原始振动信号如图4a所示。在图4b中BPFO等特征成分没有显现，因此采用本文方法处理。

(a)时域图

(b)频谱图图4　外圈故障原始信号

外圈故障信号经MED处理后，计算发现当结构元素长度为0.9T时峭度值达到最大，为8.2639，因此结构元素长度为0.9T，结果如图5所示。图5b中清楚提取出了BPFO和它的几个谐波成分(2BPFO，3BPFO，4BPFO，5BPFO和6BPFO)。此外转速fr也被找到。上述分析表明新方法有效诊断出了滚动轴承外圈故障。

(a)时域图

(b)频谱图图5　外圈故障的新方法分析结果

3.2滚动轴承内圈故障分析

为了进一步说明新方法对轴承内圈故障的有效性，以一个内圈故障为例进行分析。测试轴承型号为6203-2RSJEMSKF。采样频率为12kHz，采样点数为12 288，转速为1772r/min。计算可知外圈故障通过频率(theballpassingfrequencyinnerrace,BPFI)约为146Hz。内圈故障原始振动信号如图6a所示，其频普图如图6b所示。

(a)时域图

(b)频谱图图6　内圈故障原始信号

新方法处理后计算发现结构元素长度为0.6T时峭度值最大，为14.7514，因此结构元素长度为0.6T，相应结果如图7所示。在图7b频谱图中，内圈故障特征频率BPFI成为了主要成分，它的三个谐波成分(2BPFI,3BPFI和4BPFI)也能够识别出来，同时还能发现转速 fr和它的2倍频成分。此外特征频率BPFI及它的谐波两侧还清晰找到了丰富的边频带，即频率成分为146κ±30i(其中，i=1,2;κ=1,2,3,4)。因此，这些结果说明新方法对滚动轴承内圈故障的特征提取也是有效的。

(a)时域图

(b)频谱图图7　内圈故障的新方法分析结果

3.3有无最小熵反褶积方法的轴承内外圈故障对比

与仿真信号中的对比相同，轴承内外圈故障原始振动信号直接采用数学形态法处理，外圈和内圈故障时结构元素长度仍然是0.9T(峭度值2.9685)和0.6T(峭度值6.0521)，相应频谱图分别为图8a和图8b。图8a外圈故障中只有BPFO和2BPFO被提取出来。图8b内圈故障中，也只有BPFI和2BPFI被识别，而且特征频率两侧的边频带也未发现。与图5b和图7b相比：经过MED方法处理后，外圈故障时峭度值从2.9685增大到8.2639，所提取BPFO的谐波成分更多(从BPFO～6BPFO)。内圈故障时峭度值从6.0521增大到14.7514，所提取BPFI的谐波成分更多(从BPFI～4BPFO)，以及它们两侧边频带也较为明显。因此，事实说明MED算法的峭度值增大效果明显，有效加强了轴承内圈故障信号的脉冲特征，特征提取效果更好。

(a)外圈无MED分析

(b)内圈无MED分析图8　轴承外圈和内圈故障时的频谱图对比

4　结语

为了提取出强噪声背景下滚动轴承内外圈故障的特征频率，提出了基于最小熵反褶积的数学形态法。新方法利用MED方法进行预处理，以提高原始信号峭度值，加强故障的冲击特性。然后运用数学形态方法提取故障特征。选取了DIF滤波器作为形态算子，结构元素长度以9个不同长度中峭度值作为选取依据。通过仿真信号和滚动轴承内外圈故障分析说明新方法能够提取出丰富的故障特征。通过对比发现最小熵反褶积使得信号峭度值增大明显，轴承故障信号冲击特征得到了加强，故障特征提取效果明显。

[1]McFaddenPD,SmithJD.VibrationMonitoringofRollingElementBearingsbytheHigh-frequencyResonanceTechnique—AReview[J].TribologyInternational, 1984, 17:3-10.

[2]GeorgoulasG,LoutasT,StyliosCD,etal.BearingFaultDetectionBasedonHybridEnsembleDetectorandEmpiricalModeDecomposition[J].MechanicalSystemsandSignalProcessing, 2013, 41:510-525.

[3]GuoW,TsePW,ANovelSignalCompressionMethodBasedonOptimalEnsembleEmpiricalModeDecompositionforBearingVibrationSignals[J].JournalofSoundandVibration, 2013, 332:423-441.

[4]QiuH,LeeJ,LinJ,etal,WaveletFilter-basedWeakSignatureDetectionMethodandItsApplicationonRollingElementBearingPrognostics[J].JournalofSoundandVibration, 2006,289:1066-1090.

[5]杨宇, 程军圣. 基于LCD降噪和VPMCD的滚动轴承故障诊断方法[J]. 中国机械工程, 2013, 24(24):3338-3344.

YangYu,ChengJunsheng.ARollingBearingFaultDiagnosisMethodBasedonLCDDenoisingandVPMCD[J].ChinaMechanicalEngineering, 2013, 24(24):3338-3344.

[7]NikolaouNG,AntoniadisIA.ApplicationofMorphologicalOperatorsasEnvelopeExtractorsforImpulsive-typePerodicSignals[J].MechanicalSystemsandSignalProcessing, 2003,17:1147-1162.

[8]DongY,LiaoM,ZhangX,etal.FaultsDiagnosisofRollingElementBearingsBasedonModifiedMorphologicalMethod[J].MechanicalSystemsandSignalProcessing, 2011,25:1276-1286.

[9]WigginsRA.MinimumEntropyDeconvolution[J].Geoexploration, 1978,16:21-35.

[10]EndoH,RandallR.EnhancementofAutoregressiveModelBasedGearToothFaultDetectionTechniquebytheUseofMinimumEntropyDeconvolutionFilter[J].MechanicalSystemsandSignalProcessing, 2007,21:906-919.

(编辑王艳丽)

Applications of Mathematical Morphology Method to Fault Feature Extraction of Rolling Bearings Based on Minimum Entropy Deconvolution

Gong Tingkai1,2Yuan Xiaohui1Wang Xiyang2

1.Huazhong University of Science and Technology, Wuhan,430074 2.Nanchang Hangkong University, Nanchang, 330063

Aiming at the extractions of fault features of rolling bearings under the strong noise background, a novel method, called mathematical morphology based on minimum entropy deconvolution (MMBMED) was proposed herein. In this method, MED was first introduced to enhance the impact behaviors of the signals. And the fault characteristics of defective bearings were extracted by MM method, where the DIF filter with catching the bidirectional pulses was adopted as morphological operator, and kurtosis was employed as the criterion to the length selection of structural elements. The effectiveness of the new method was validated by both of simulation signals and vibration signals of rolling bearings with the outer and inner race faults. Finally, the results disclose that MED is effective to increase signal kurtosis and to strengthen the impulsive characteristics.

minimum entropy deconvolution(MED); mathematical morphology(MM); fault diagnosis; rolling bearing

2015-07-13

国家自然科学基金资助项目(51379080，AA201406168)

TH 165；TN911

10.3969/j.issn.1004-132X.2016.18.009

龚廷恺，男， 1977年生。华中科技大学数字化水电工程学院博士，南昌航空大学飞行器工程学院讲师。主要研究方向为旋转机械故障诊断及信号处理。袁晓辉，男， 1971年生。华中科技大学数字化水电工程学院教授、博士研究生导师。王细洋，男，1967年生。南昌航空大学飞行器工程学院教授。