高阶频率加权能量算子在列车轴箱轴承故障诊断中的应用

刘泽潮,张兵,易彩,吴文逸,黄晨光

(西南交通大学牵引动力国家重点实验室,610031,成都)

轴箱轴承作为轨道交通车辆走行部中的关键旋转部件,不但承受着车辆全部的垂向载荷,同时遭受着强烈的轮对冲击。恶劣的服役环境使其容易发生突发性失效,从而威胁列车的运行安全[1],因此需要对轴承的服役状态进行实时监控,保证列车的运行安全。复杂的服役环境导致采集到的振动信号成分更为复杂,对故障冲击的识别提出了更高的要求[2]。

根据循环平稳理论可知,轴承的故障信号具有典型的二阶循环平稳特性,即信号的时域信号是非周期的,但是由故障冲击引起的瞬时能量波动具有典型的周期特性[3]。因此,可通过识别轴承振动信号中的瞬时能量变化,实现轴承故障诊断与识别。

根据振动方程可知,在计算故障冲击引起的能量波动时不仅应当包含瞬时幅值部分,而且应当包含瞬时频率部分的权重。希尔伯特变换(HT)通过对故障冲击包络信息的计算,实现信号中冲击信息的识别[4-5]。但HT在对信号进行解调时,计算的仅仅是信号的瞬时幅值,不包含瞬时频率部分的权重。当信号信噪比较低且存在其他干扰频率时,解调效果明显下降[6]。特别是当干扰频率的幅值较大时,在HT解调谱中的主要谱线是干扰频率及其调制后的频率成分[7]。为更好地提取信号中瞬时能量的波动,Teager提出了Teager能量算子(TEO)[8]。TEO通过原始信号、信号一阶导数和二阶导数的非线性组合,对信号中的瞬时能量波动进行了有效估计,实现了瞬态冲击特征的提取[9-13],但是TEO的本质是对信号中瞬时能量的估计,原则上其不应当存在无意义的负能量。为避免像TEO一样存在无意义的负值,且又同时在解调时保留信号的瞬时频率部分的权重,O′Toole提出了一种非负的能量算子——频率加权能量算子(FWEO)[14]。FWEO通过计算信号导数的HT,得到信号的瞬时能量估计,故而其不存在负能量部分[15]。众多学者通过对FWEO的研究表明,FWEO具有比HT和TEO更好的抗干扰特性[16]。但是,当信号中的噪声与干扰频率的能量较大时,FWEO依然无法实现信号瞬时冲击特征的识别。

FWEO的核心是通过求导加入瞬时频率的权重,在此基础上本文通过计算信号的高阶导数,加大瞬时频率部分的权重,提高解调的可靠性。将所提出的解调方法称为高阶频率加权能量算子(HFWEO)。同时,提出相关峭度准则,确定合适的HFWEO阶次,保证解调的可靠性,并通过仿真信号验证了所提算法在强干扰环境下依然具有较好的鲁棒性。

1 基础理论

1.1 振动方程

根据牛顿力学方程,无阻尼单自由度线性系统自由振动方程可表示为[17]

(1)

x(t)=Acos(ωt+φ)

(2)

式中:A为振动的瞬时幅值;ω为系统的共振频率;φ为初始相位。

系统振动的总能量E为

(3)

将式(2)代入式(3),则系统振动的总能量E可近似表示为

(4)

1.2 频率加权能量算子

FWEO通过对信号求导,在计算信号的瞬时能量时加入瞬时频率部分的权重,故而FWEO又叫做导数包络能量算子[12]。

频率加权能量算子的表达式为

(5)

式中:S[·]表示信号的复包络;H[·]表示希尔伯特变换;φ[x(t)]为频率加权能量算子。

根据式(2)可知

(6)

(7)

将式(6)(7)代入式(5),可得

φ[x(t)]=(-Aωsin(ωt+φ))2+

(Aωcos(ωt+φ))2=A2ω2

(8)

式(8)表明,FWEO与TEO一样计算的是信号的瞬时能量,而非HT的瞬时幅值平方[15]。由于FWEO本质上是计算信号导数的包络,所以没有像TEO一样存在无意义负能量部分。

2 高阶频率加权能量算子

2.1 基本原理

FWEO通过对信号求导,在计算瞬时能量时引入瞬时频率部分的权重,从而使得FWEO具有同Teager能量算子相同的性质。鉴于此,本文提出通过求取信号的高阶导数,加大瞬时频率部分的权重,进而提高能量算子在干扰情况下解调的鲁棒性。因此,将所提的方法称为高阶频率加权能量算子(HFWEO)。

HFWEO的表达式为

ξ[x(t),m]=S[xm(t)]=|xm(t)+

jH[xm(t)]|2=xm(t)2+H[xm(t)]2

(9)

式中:xm(t)表示x(t)的m阶导数;ξ[x(t),m]为m阶频率加权能量算子。

根据式(2),可得

xm(t)=Aωmcos(ωt+φ+mπ/2)

(10)

H[xm(t)]=Aωmsin(ωt+φ+mπ/2)

(11)

将式(10)和(11)代入式(9),可得

ξ[x(t),m]=(xm(t))2+(H[xm(t)])2=

(Aωmcos(ωt+φ+mπ/2))2+

(Aωmsin(ωt+φ+mπ/2))2=A2ω2m

(12)

根据式(12)可知,当阶数m等于1时,HFWEO即为FWEO。随着阶数m的增加,瞬时频率ω部分的权重呈指数级增长,故HFWEO可以更好地追踪信号中瞬时能量的微弱变化,从而实现信号中故障冲击的分离。

(13)

根据递推公式,离散信号x(n)的m阶导数xm(n)可表示为

xm(n)=xm-1(n+1)-xm-1(n)

(14)

将式(14)代入式(9),可得HFWEO的离散形式

ξ[x(n),m]=S[xm(n)]=

|xm(n)+jH[xm(n)]|2=

(xm-1(n+1)-xm-1(n))2+

H[xm-1(n+1)-xm-1(n)]2

(15)

式(2)的离散形式为

x(n)=Acos(ωn+φ)

(16)

将式(16)代入式(15),可得

ξ[x(n),m]=22mA2sin2m(ω/2)

(17)

当ω/2<π/4时sin(ω/2)≈ω/20,式(17)可近似为

(18)

此时,离散HFWEO与连续HFWEO解调结果近似相等。通过上述分析可知,离散化误差是由sin(ω/2)≈ω/2导致,且sin2m(ω/2)与(ω/2)2m的差值会随着阶次m的增大而增大。连续的差分运算会导致计算误差的增大,故应当选取合适的阶次m。

2.2 解调特性分析

通常单个轴承的故障信号x(t)可以简化为指数信号与正弦信号的调幅信号[3],表达式为

x(t)=Ae-β tcos(tω+φ)

(19)

式中:β为阻尼系数。

谐波干扰信号为

v(t)=Lhcosωht

(20)

式中:Lh为谐波干扰信号的幅值;ωh为谐波干扰信号的频率。

信号r(t)=x(t)+v(t)中包含轴承的故障信号x(t)和谐波干扰信号v(t)。信号r(t)的FWEO解调结果为

(21)

式中

λ(t)=-2AβLhωhe-β tsin(-tωh+tω+φ)+

2AωLhωhe-β tcos(-tωh+tω+φ)

(22)

从式(21)可以看出,第1部分A2(β2+ω2)e-2β t为轴承故障信号的平方包络A2e-2β t再乘以加权系数(β2+ω2),此部分对应的是轴承故障特征频率的低频部分。第2部分为干扰信号与轴承信号的互相调制部分,该部分位于高频。第3部分与干扰信号相关。FWEO的信号干扰比σSIR可以通过幅值调制部分(式(21)的第1部分)和干扰部分(式(21)的第3部分)进行计算[19]

σSIR(φ(r(t)))=

(23)

式中:Tp为轴承故障冲击间隔。

信号r(t)的二阶HFWEO解调结果为

式中

二阶HFWEO的σSIR为

(24)

二阶HFWEO与FWEO的σSIR之比为

(25)

信号r(t)的m阶HFWEO的σSIR为

(26)

则

(27)

通过式(27)可知,随着HFWEO阶次m的增大,HFWEO的σSIR逐渐提高。阶次m越大,HFWEO对信号中的谐波干扰抑制效果越好。

通过对比不同HFWEO解调结果的信噪比σSNR,评估HFWEO阶次对噪声的影响。HFWEO解调结果的σSNR定义为

(28)

式中:x(t)为轴承的故障信号;n(t)为白噪声。

原始信号的σSNR为5 dB,图1为根据式(28)计算的1至10阶HFWEO解调结果的σSNR,其中一阶HFWEO即FWEO。从图中可以看出,随着阶次的增大HFWEO的σSNR逐渐降低,这是由于连续的差分运算,导致解调结果中包含更多的高频噪声。图2为σSNR=5 dB轴承故障仿真信号FWEO和10阶HFWEO的解调结果,显然FWEO结果中可以发现明显的故障冲击,但是10阶HFWEO结果中的故障冲击很微弱。

图1 不同阶次HFWEO解调结果σSNR

(a)FWEO结果(σSNR=6.22 dB)

(b)10阶HFWE结果(σSNR=0.463 7 dB)图2 两种方法解调结果的对比

综上可知,随着阶次的增大,HFWEO对谐波干扰信号的抑制效果越来越好,但是随着连续的差分运算,导致HFWEO解调结果中引入了更多的高频噪声。故应当设置合适的阶次m,保证HFWEO在对谐波干扰具有较好的抑制效果的同时又避免引入更多的噪声。

2.3 阶次确定

通过上述分析可知,阶次对HFWEO的解调结果有着至关重要的影响,故在此提出相关峭度准则(CK)KCK来确定合适的阶次m。

KCK指标不仅可以对序列中的冲击特性进行表征,同时通过平移周期T,可以对序列中的周期性进行表征[20]。序列yn的KCK表达式为

(29)

式中:H为平移周期的个数;T为平移周期长度;N为序列长度。

在此,yn代表不同阶次下的HFWEO解调结果ξ(yn,m)。通过式(29)可以得到不同阶次HFWEO解调结果的KCK值,则最优的HFWEO阶次mbest为

mbest=arg maxm(KCK(ξ(yn,m),T,H))

(30)

3 仿真验证

使用仿真信号对所提算法进行验证,仿真信号的表达式为

u(t-mTp-τm)

(31)

式中:Am是冲击的幅值,Am=1;L是冲击个数;Tp是故障冲击的间隔,特征频率fc=1/Tp=104.5 Hz;ωr是故障引起共振的频率,ωr=4 000 Hz;τm为随机滑移系数,通常取0.01Tp～0.02Tp;u(t)是单位阶跃函数;β是阻尼系数,β=1 500 N·s/m。

信号的采样频率fs=12 000 Hz,信号时间长度为2 s。在信号中加入不同信号噪声比的白噪声与不同信号干扰比的正弦干扰信号[17],以研究所提方法在不同噪声与干扰情况下的解调效果。

首先,对比几种不同的解调方法它们在高σSNR与高σSIR情况下的解调效果。在信号中加入σSNR=0 dB的白噪声和σSIR=-5 dB的正弦干扰信号,干扰信号的频率分别为50、270、2 043和2 810 Hz。加噪后信号的时域波形与频谱如图3a和3b所示。虽然从仿真信号的频谱在图3b中可以发现很明显的干扰频率及噪声,但由于信号的σSNR与σSIR都较高,依然可以很明显地观察到共振频带。

(a)时域信号

(b)傅里叶谱图3 σSNR=0 dB和σSIR=-5 dB时仿真信号

(a)包络谱

(b)Teager谱

(c)FWEO谱图4 σSNR=0 dB和σSIR=-5 dB时3种解调方法结果的对比

使用HT、TEO、FWEO与HFWEO 4种方法对仿真信号进行解调,HT、TEO、FWEO解调结果如图4所示。可以发现,在高σSNR与σSIR情况下,HT、TEO、FWEO 3种方法的解调谱中都可以发现故障特征频率fc及其倍频,解调谱中存在干扰频率。这说明3种解调方法对干扰信号都较敏感,只是由于干扰频率的能量较低,所以从解调谱中依然可以发现故障特征频率及其倍频。

图5 HFWEO相关峭度

图6 σSNR=0 dB和σSIR=-5 dB时的仿真信号HFWEO谱

图5为1至10阶的HFWEO相关峭度,第2阶HFWEO的相关峭度值最大,当阶次继续增大时,由于连续差分运算引入了更多的噪声,所以相关峭度值逐渐降低。图6为归一化后1至10阶的HFWEO解调谱。从不同阶次的解调谱中都可以识别到故障特征频率fc及其倍频,同时可以发现解调谱中不存在干扰频率。但随着阶次的增大,解调谱中的背景噪声也逐渐增大,故障特征频率的倍频数量减少。综上所述,HFWEO通过增加瞬时频率部分的权重,保证在解调时对干扰频率具有更好的抑制效果。同时,通过相关峭度指标可以确定合适的HFWEO阶次,确保HFWEO在对谐波干扰进行抑制的同时引入更少的高频噪声。

保持信号的σSNR=0 dB不变,将σSIR降低为-20 dB,对比4种方法对高σSNR、低σSIR信号的解调效果。加噪信号的时域波形与频谱如图7a和7b所示。由于信号的σSIR较低,此时谐波频率的幅值较大,故从信号的频谱中很难识别出共振频带。

(a)时域信号

(b)傅里叶谱图7 σSNR=0 dB和σSIR=-20 dB时的仿真信号

(a)包络谱

(b)Teager谱

(c)FWEO谱图8 σSNR=0 dB和σSIR=-20 dB时的3种解调方法结果对比

HT、TEO和FWEO的解调谱如图8所示。由于信号的σSIR较低,所以从这3种方法的解调谱中都只能识别出干扰频率,而无法识别出故障特征频率fc及其倍频。使用HFWEO对仿真信号进行解调,1至10阶HFWEO相关峭度如图9所示,归一化后1至10阶的HFWEO谱如图10所示。在前4阶HFWEO解调谱中只能识别到干扰频率,此时相关峭度值也较低;到第5阶HFWEO解调谱中发现故障特征频率fc及其倍频,但是依然可以发现干扰频率,此时相关峭度值依然较小。随着阶次的增大,在7阶以后HFWEO解调谱中只能识别到故障特征频率fc及其倍频,而不再存在干扰频率。此时,相关峭度值较大,且在第9阶时HFWEO的相关峭度值达到最大。同时,随着阶次的增大,高频噪声也在逐渐增大,第10阶HFWEO的相关峭度值反而降低。这说明,HT、TEO与FWEO的抗谐波干扰性较差,而本文所提的HFWEO通过连续的导数运算,加大了瞬时频率部分的权重,使得其在较强干扰时仍然具有较好的鲁棒性。同时,相关峭度准则可以准确地得到合适的HFWEO阶次,保证HFWEO的可靠性。

图9 HFWEO相关峭度

图10 σSNR=0 dB和σSIR=-20 dB时的仿真信号HFWEO谱

(a)时域信号

为进一步验证4种解调算法对低σSNR和低σSIR信号的解调效果,在仿真信号中加入σSNR=-5 dB的白噪声与σSIR=-20 dB的干扰信号。加噪仿真信号的时域波形与频谱如图11a与11b所示。由于信号的σSNR与σSIR都很低,所以从频谱中同样很难识别出共振频带。

(b)傅里叶谱图11 σSNR=-5 dB和σSIR=-20 dB时的仿真信号

分别使用HT、TEO和FWEO对加噪后仿真信号进行解调,解调谱如12所示。可以看出,解调谱中的主要频率成分均是干扰频率,无法识别出故障特征频率fc及其倍频。

(a)包络谱

(b)Teager谱

(c)FWEO谱图12 σSNR=-5 dB和σSIR=-20 dB时的3种解调方法结果对比

图13和图14所示分别为1至10阶HFWEO的相关峭度值和归一化后的解调谱,显然当阶次较低时,HFWEO同样无法对信号中的瞬态冲击进行提取,解调谱中只能发现干扰频率,故此时相关峭度较低。随着阶次的增加,干扰频率的幅值逐渐减小,而故障特征频率fc及其倍频的幅值逐渐增加。当阶次大于7以后,HFWEO解调谱中不存在干扰频率,但随着阶次的增加,解调谱中的噪声也在增大,从而影响对故障特征频率的观测效果,所以第7阶时的HFWEO相关峭度值最大。

图13 HFWEO相关峭度

图14 σSNR=-5 dB和σSIR=-20 dB时的仿真信号HFWEO谱

综上所述,本文所提的HFWEO能量算子对干扰具有很好的抑制作用,可以在低σSIR情况下具有良好的解调效果,同时所提的相关峭度准则可以准确地确定最佳的HFWEO阶次。

4 试验验证

4.1 铁路货车轮对跑合试验台

为验证HFWEO对实际振动信号的解调效果,本文对货车轮对跑合试验台的振动数据进行分析,关于试验台的详细介绍见文献[18]。测试轴承型号为197726,在实际运行过程中,轴承形成内圈损伤,内圈滚道出现明显剥落。轮对的旋转速度为465 r/min,轴承的旋转频率为fr=7.75 Hz,采样频率为fs=12.8 kHz。根据轴承参数,计算得到内圈故障频率fBPFI=88.25 Hz。现对1 s的振动数据进行分析,振动信号的时域波形及其频谱如图15a和15b所示。当轴承存在内圈损伤时,振动信号的解调谱中应当出现转频fr、内圈故障频率fBPFI及其倍频和以内圈故障频率fBPFI为中心且以转频fr为间隔的边频带[3]。

使用HT对信号进行解调,包络谱如16a所示。包络谱中故障特征频率fBPFI的谱线并不明显,并且没有发现以fBPFI为中心、以转频fr为间隔的边频带。然后,使用TEO与FWEO对信号进行解调,解调谱如图16b和16c所示。同包络谱相同,从解调谱中同样只能隐约观测到fBPFI,但依然无法识别到以fBPFI为中心且以转频为fr间隔的边频带。这说明HT、TEO和FWEO都无法对信号中的轴承故障冲击进行提取,因此无法从解调谱中发现与内圈故障对应的频谱特征。

(a)时域信号

(b)傅里叶谱图15 铁路货车轮对跑合试验台振动信号及频谱

(a)包络谱

(b)Teager谱

(c)FWEO谱图16 3种解调方法结果铁路货车轮对跑合试验台振动数据解调对比

使用本文所提的HFWEO对信号进行解调,并计算各阶次HFWEO相关峭度值,其结果如图17a所示。第7阶HFWEO的相关峭度值最大,其解调谱如17b所示。从解调谱中不但可以识别到转频fr、内圈故障特征频率fBPFI,而且还可以识别到以fBPFI为中心且以转频fr为间隔的边频带。这与内圈故障对应的频谱特征相吻合,说明HFWEO可以很好地从信号中分离出故障冲击信息。

(a)HFWEO相关峭度

(b)第6阶HFWEO谱图17 铁路货车轮对跑合试验台振动数据HFWEO解调结果

在货车轮对跑合试验台的振动信号中噪声与干扰较大,轴承故障冲击的能量较微弱,因此传统解调方法无法对故障冲击进行提取。根据式(27)可知,随着阶数m的增加,HFWEO的抗干扰性越来越好,故从高阶HFWEO谱中可以识别到故障特征频率。同时,相关峭度准则可以准确地确定HFWEO的阶次。

4.2 高速列车轮对跑合试验台

为进一步验证本文所提算法在复杂环境下的有效性,使用高速列车轮对跑合试验台的振动数据对算法进行验证。图18所示为高速列车轮对跑合试验台,试验台包括电机、驱动轮对、高速动车组轮对、高速动车组轴箱和载荷加载装置。电机带动驱动轮对驱动高速动车组轮对,进而带动轴箱内的轴承旋转,载荷加载装置可以对轴箱施加垂向载荷,振动传感器位于轴箱上。

图18 高速列车轮对跑合试验台

在轴承的外滚道表面以120°夹角设置3处人工损伤,3处损伤的深度为1 mm,宽度分别为1、3和5 mm,损伤部位如图19所示。

图19 轴承外圈损伤示意图

试验时,在轴箱垂向上施加50 kN的垂向载荷,采样频率为10 kHz。驱动轮对的转速为100 km/h,则轴承的旋转频率fr=10.28 Hz,因此根据轴承的几何参数和转频可以计算得到外圈故障特征频率fBPFO=83.29 Hz。

图20a和20b所示为1 s的高速列车轮对跑合试验台振动信号及其频谱。显然,由于试验台的结构复杂,因而信号中的频率成分更加丰富。现场恶劣的采集环境导致信号中存在50、100 Hz的谐波干扰频率。故高速列车轮对跑合试验台的振动信号比货车轮对跑合试验台的振动信号更加复杂,对算法的鲁棒性要求更高。

(a)时域信号

(b)傅里叶谱图20 高速列车轮对跑合试验台振动信号及频谱

使用HT、TEO和FWEO对信号进行解调,解调谱如图21所示。由于谐波干扰和噪声的影响,从它们的解调谱中都无法识别到故障特征频率。在包络谱中可以发现很明显的干扰频率,这说明HT对谐波干扰较敏感。TEO和FWEO计算的是信号的瞬时能量,对谐波干扰的抑制效果优于HT。但是,振动信号中还存在轮对冲击的成分,且其能量大于轴承故障冲击能量,因此从TEO与FWEO解调谱使用本文所提HFWEO对振动信号进行解调,1至10阶HFWEO相关峭度如22a所示,第10阶相关峭度最大,其对应的HFWEO谱如22b所示。从HFWEO中可以发现很明显的故障特征频率fBPFO及其二倍频,说明本文算法成功诊断出轴承外圈故障。这是因为,轮对冲击的能量虽然较大,但是轴承故障冲击的瞬时频率比轮对冲击的瞬时频率更高[1],HFWEO通过不断加大瞬时频率ω的权重,从而对瞬时频率较高的冲击提取效果更好,从HFWEO谱中可以更明显识别出轴承外圈故障特征频率fBPFO及其倍频。

(a)包络谱

(b)Teager谱

中仅可以识别到转频fr及其二倍频。

(c)FWEO谱图21 高速列车车轮对跑合试验台振动数据4种解调方法结果对比

(a)HFWEO相关峭度

(b)第10阶HFWEO谱图22 高速列车车轮对跑合试验台振动数据HFWEO解调结果

5 结 论

通过增加瞬时频率部分的权重,提高了HFWEO抗干扰的能力,保证在干扰均较大时,依然具有良好的解调效果。但是,连续的差分运算导致解调时引入了高频噪声,因此需要确定合适的阶数。通过相关峭度准则可以确定合适的HFWEO阶次,保证算法在对干扰信号具有较好鲁棒性的同时引入更少的高频噪声。仿真信号与试验信号分析表明,本文算法可以在信号中存在严重干扰情况下依然可以实现轴承故障冲击的解调,从解调谱中识别故障特征频率。