基于故障特征系数模板的变转速滚动轴承故障诊断

程卫东， 赵德尊， 刘东东

(北京交通大学 机械与电子控制工程学院，北京 100044)

基于故障特征系数模板的变转速滚动轴承故障诊断

程卫东， 赵德尊， 刘东东

(北京交通大学 机械与电子控制工程学院，北京 100044)

与恒转速相比，机械中普遍存在的变转速工作模式使滚动轴承的故障诊断更加困难；另外变转速条件下的常规方法—阶比分析存在误差以及计算效率方面的问题，因此，提出了基于故障特征系数模板的滚动轴承故障诊断方法。该方法主要包括六部分：①根据目标轴承的几何参数计算其故障特征系数以设定模板；②利用快速谱峭度滤波算法对滚动轴承振动信号进行滤波；③根据Hilbert变换以及短时傅里叶变换计算滤波信号的包络时频图；④通过峰值搜索算法从滤波信号的包络时频图中提取瞬时故障特征频率趋势线；⑤根据转速脉冲信号计算滚动轴承的转速曲线；⑥瞬时故障特征频率与瞬时转频相比获取瞬时故障特征系数，进而通过故障特征系数模板实现滚动轴承的故障诊断。随即以变转速情况下的故障轴承仿真信号以及实测的外圈故障、内圈故障和健康轴承的振动信号为例验证了该算法的有效性。

变转速；滚动轴承；故障特征系数模板；故障诊断

在机械、石化、能源、航空等不同领域中，变转速工作模式普遍存在[1-2]。变转速主要表现在两方面，一是大范围的转速变化过程，如旋转机械的启停；另一方面表现为转速波动准平稳过程，如受载荷、工况的影响转速发生波动。变转速条件下，利用包络分析为核心的滚动轴承故障诊断方法获得的频谱图将变得模糊，无法诊断轴承故障。因此变转速滚动轴承的故障诊断引起了学者的广泛关注。

目前，阶比分析[3-5]是变转速旋转机械故障诊断的重要的技术。其原理是获得相对于参考轴的恒定角增量采样，将时域非周期信号转化为角域的周期信号。阶比分析技术主要包括硬件阶比跟踪和计算阶比跟踪两种[6]。硬件阶比跟踪是通过特定的硬件设备(采样率合成器、抗混叠跟踪滤波器、转速计和采样率可调的采样装置) 实现信号的恒定角增量采样。该算法对硬件设备的要求较高，而且只在转速变化较缓的工况下效果明显，因此使用范围受到一定限制。在上述背景下，计算阶比跟踪运用而生，该算法的特点是采用恒定采样率同步采集振动信号和转速信号，通过相关算法用软件实现非周期信号到周期信号的转换。相对于硬件阶比跟踪，计算阶比跟踪的硬件成本以及传感器安装的复杂性大幅降低，因此一经提出便受到很多学者的关注，并广泛应用于工程实践。然而，计算阶比分析过程中存在的缺陷难以忽略。SAAVEDRA等[7]人认为对振动信号进行阶比跟踪时所用到的插值算法存在计算误差，而且该误差无法避免。另外阶比跟踪过程中鉴相时标的获取需要求解大量的高次方程，计算复杂，而且方程会出现无解的情况。

鉴于如今变转速工作模式下旋转机械的诊断方法都是以阶比分析为核心，而该算法又存在一定的缺陷，一些学者尝试直接根据振动信号的特征通过短时傅里叶变化、小波变换等时频分析方法实现旋转机械的故障诊断。LI等[8]利用旋转机械振动信号的快速傅里叶变换，小波变换以及双谱构造故障特征，应用于马尔科夫模型实现变转速工况下的旋转机械故障诊断；MELTZER等[9]通过对振动信号的小波幅值进行时域到极坐标的转换以实现齿轮故障诊断，并能准确识别故障发生的位置及故障的严重程度； LI等[10]利用广义傅里叶变换对以小波变换为基础的Synchrosqueezing Wavelet transform算法进行改进，通过得到的时频分析结果实现升降速情况下的齿轮故障诊断；上述算法虽然并没有直接应用于滚动轴承的故障特征识别，但是依然为变转速工作模式下滚动轴承故障诊断提供思路上的指导。

通过以上分析，本文从故障轴承振动信号的特征出发，提出了基于故障特征系数模板的滚动轴承故障诊断方法。该方法主要包括以下六部分：①根据目标轴承的几何参数计算其故障特征系数以设定匹配模板；②利用快速谱峭度滤波算法对滚动轴承振动信号进行滤波；③根据Hilbert变换以及短时傅里叶变换计算滤波信号的包络时频图；④通过峰值搜索算法从滤波信号的包络时频图中提取瞬时故障特征频率趋势线；⑤根据转速脉冲信号计算滚动轴承的转速曲线；⑥瞬时故障特频率与瞬时转频相比获取瞬时故障特征系数，进而通过故障特征系数模板实现滚动轴承的故障诊断。随即以变转速情况下的故障轴承仿真信号以及实测的外圈故障、内圈故障和健康轴承的振动信号为例验证了该算法的有效性。

1 算法原理

1.1 瞬时故障特征频率趋势线及其提取算法

当轴承某一位置表面出现损伤时，受载运行过程中该损伤点撞击其他元件表面将产生脉冲，故障脉冲每转出现的次数即为故障特征频率。当轴承外圈、内圈以及滚珠出现故障时对应的故障特征频率计算公式分别为

(1)

(2)

(3)

式中：n为滚动体个数；fr为轴承转频；D为滚动轴承节圆直径；d为滚珠直径；α为接触角。

对式(1)～式(3)进行改写，计算轴承外圈、内圈以及滚珠的故障特征系数Fo、Fi和Fb分别为

(4)

(5)

(6)

恒转速下，故障轴承的包络谱中存在明显的表征瞬时故障特征频率及其倍频的突出峰值。根据该现象以及式(4)～式(5)不难得到以下两点：①当轴承类型确定时，其内圈、外圈以及滚珠故障时所对应的故障特征频率与轴承转频的比值，即故障特征系数唯一确定；②当轴承转频变化时，其故障特征频率也以与转频相同的变化规律同步变化，而且在故障轴承振动信号的包络时频图中表现出明显幅值优势。因此可以利用峰值搜索算法从轴承信号包络时频图中提取瞬时故障特征频率(InstantaneousFaultCharacteristicFrequency,IFCF)。本文利用短时傅里叶变换计算轴承包络信号的时频图，计算公式为

(7)

式中：x(τ)为滚动轴承的振动信号；g(τ)为高斯窗函数。

通过峰值搜索算法从式(7)计算的包络时频谱中提取IFCF，其过程总结为：①将包络信号时频图中的时间轴等分为i份(t1,t2,…,ti)，搜索每一时刻对应的时频图中的最大幅值点；②记录上述幅值点所对应的频率值，上述频率值即为相应时刻的故障特征频率；③根据时间大小对频率值进行排序，得一个IFCF的集合，即IFCF趋势线

IFCFi=arg max{IFSi} (i=1,2,…,n)

(8)

式中：IFS(InstantaneousFrequencySpectrum)为瞬时频谱；i为瞬时频谱对应的序号,取值为1～n；arg max函数表示ITSi取得最大值时对应的频域横坐标，瞬时频谱的总数将由原始信号与窗函数的长度决定。

1.2 快速谱峭度滤波

本节利用快速谱峭度滤波算法对故障轴承振动信号进行滤波以更准确的提取IFCF趋势线。针对峭度指标在衡量振动信号中冲击成分多少时容易受到噪声干扰这一问题，ANTONI等对峭度算法进行了改进提出了谱峭度指标这一概念，公式为

(9)

式中：〈·〉为均值；|·|为取模；H(t,f)为故障轴承振动信号x(t)在频率f处的时频复包络，可由STFT获取。

谱峭度的计算方法主要包括STFT、小波变换以及基于快速Kurtogram的算法。本文选用快速Kurtogram算法计算普峭度值。该算法构造一系列树状双子带带通解析滤波器组来实现的。具体的计算过程如下：

(1) 计算各层复包络。分别构造两个具有线性相位的准解析低通h0(n)和高通h1(n)滤波器

h0=h(n)ejπn/4=h(n)[cos(πn/4)+jsin(πn/4)]

(10)

h1=h(n)ej3πn/4=h(n)[cos(3πn/4)+jsin(3πn/4)]

(11)

(12)

(13)

式中，*表示卷积。

(2) 计算各级滤波包络信号的峭度值

(14)

(3)找出谱峭度图中峭度值最大的区域，确定其对应的中心频率和滤波带宽，进而通过带通滤波获取反映轴承故障的目标信号。

2 基于故障特征系数模板的滚动轴承故障诊断方法

针对变转速工作模式下，以包络分析为核心的轴承故障诊断方法出现的频谱模糊现象，本文根据IFCF趋势线与轴承转频具有固定的比例关系以及其在故障轴承振动信号的包络时频图中的可提取性等特性提出了基于故障特征系数模板的滚动轴承故障诊断方法。该方法的流程图如图1所示，具体步骤如下：

(1) 根据滚动轴承的几何参数计算其故障特征系数用以设定匹配模板；

(2) 通过快速谱峭度算法确定最优的滤波参数对故障轴承振动信号进行滤波；

(3) 利用Hilbert变换和STFT获得滤波信号的包络时频图；

(4) 通过峰值搜索算法提取包络时频图中的IFCF趋势线；

(5) 根据转速脉冲信号计算滚动轴承的转频曲线；

(6) 瞬时故障特征频率与瞬时转频相比获取瞬时故障特征系数(InstantaneousFaultCharacteristicCoefficient,IFCC)，进而通过预设故障特征系数模板实现滚动轴承的故障诊断。

图1 基于故障特征系数模板的滚动轴承故障诊断方法流程图Fig.1 Flowchart of rolling element bearing fault diagnosis based on fault characteristic coefficient matching template

3 仿真分析

本节构造了升速情况下的故障轴承仿真信号对提出的方法进行了验证。为实现变转速条件下故障轴承振动信号的仿真，本文对恒转速条件下故障轴承振动信号的仿真模型[11]进行了改进。改进后的模型为

(15)

式中：xv(t)为变转速下的故障轴承振动信号；n(t)为高斯白噪声；Ai=a·ti+b为第i个冲击对应的幅值；B为衰减系数；wr为由轴承故障激起的共振频率；τi=λi·T为第i次故障冲击时刻的滑移误差；Ti为故障冲击对应的时刻坐标，计算过程简述如下：

本文设定轴承外圈模拟故障特征系数Fo=3.5，即轴承每转一周所激起的故障冲击次数为3.5。预设轴承转频方程v(t)=4t+35；故障轴承转过的角度与故障冲击对应的时刻坐标之间的关系可由式(16)计算

(16)

式中：θi=i/Fo(i=1,2,…,K)根据以上分析即可求出冲击对应的时刻坐标Ti，如图2所示。

仿真信号模型的其他参数如表1所示。

表1 仿真模型参数

图2 升速仿真轴承信号的故障冲击示意图Fig.2 Impulses of simulated signals under speed-up condition

根据上述仿真模型及参数得到故障轴承仿真信号的时域波形图如图3所示。图4表示原始信号的快速谱峭度图，从图中可以确定尺度为1.5时仿真信号的峭度值达到最大，中心频率为5 000 Hz，与预设值完全相等。利用图4中确定的最优滤波参数对仿真信号进行带通滤波。计算滤波信号的包络时频图如图5所示。图5中含有明显的IFCF趋势线及其倍频。作为对比，直接求取仿真信号的包络时频图如图6所示，从图6中无法清楚找的IFCF趋势线及其倍频。因此说明快速谱峭度滤波算法对提取轴承信号中的故障冲击成分的有效性。

图3 升速仿真故障轴承信号的时域波形图Fig.3 Simulated faulty bearing signal in time domain

图4 仿真信号的快速谱峭度图Fig.4 Kurtogram of simulated signal

图5 滤波信号的包络时频图Fig.5 Envelope time-frequency representation of filtered signal

图6 仿真信号的包络时频图Fig.6 Envelope time-frequency representation of simulated signal

利用峰值搜索算法从图5中提取IFCF趋势线。图7上部的方块为提取到的IFCF趋势线，下部直线则为理论的转频曲线。通过对比可以看出，提取的IFCF趋势线与预设的转频曲线变化趋势相同。将IFCF趋势线与转频曲线相比，得到瞬时故障特征系数如图8中的十字符号。图8中的五角星代表预设的故障特征系数模板。从图中可以看出获取的瞬时故障特征系数与预设的模板几乎想重合，其中平均误差为0.58%。从而说明本文算法对于变转速滚动轴承故障诊断的有效性。

图7 IFCF趋势线与预设的转频曲线Fig.7 IFCF trend and designed rotational frequency curve

图8 故障轴承模板匹配结果Fig.8 Faulty bearing matching result

4 实验验证

本节利用滚动轴承实验台测取变转速工作模式下的轴承外圈故障、内圈故障振动信号以及健康轴承振动信号对本文算法的有效性进一步验证。其中目标轴承的具体参数见表2。

表2 目标轴承参数

根据表2中滚动轴承的几何参数计算得到的外圈故障、内圈故障以及滚珠故障所对应的故障特征系数分别为Fo=2.548、Fi=4.452，和Fb=1.7。

根据计算得到的故障特征系数预设故障匹配模板如图9中五角星所示，其中从上到下依次代表内圈故障、外圈故障和滚珠故障。

图9 外圈故障轴承模板匹配结果Fig.9 Faulty bearing matching result

首先利用加速度传感器测取外圈故障的滚动轴承振动信号，编码器测量对应的转速脉冲信号。采集卡的采样频率设置为24 000 Hz。振动信号的时域波形图如图10所示，对应的转速脉冲信号(局部)如图11所示，采样时长为3.5 s。转速脉冲信号转化为转速曲线的步骤主要包括以下三步：①记录转速脉冲信号中下降沿的位置并计算两个下降沿之间的时间；②通过编码器的参数确定两个脉冲间的旋转角度，进而计算两个脉冲间的平均转频，得到转频离散点；③通过多项式曲线拟合上述离散点得到转频曲线。计算得到的转频曲线如图12中的下部曲线所示，转频从30 Hz升高到49 Hz进而又降到33 Hz。

图10 外圈故障轴承振动信号的时域波形图Fig.10 Outer faulty bearing signal in time domain

图11 外圈故障轴承信号对应的局部转速脉冲信号Fig.11 Local pulse signal of outer faulty bearing signal

图12 IFCF趋势线与测量的转频曲线Fig.12 IFCF trend and measured rotational frequency curve

利用本文算法计算的振动信号的快速谱峭度图如图13所示，其中在尺度1.5，中心频率为6 000 Hz处对应着最大峭度值。利用最大峭度值对应的最优滤波参数对振动信号进行滤波，计算滤波信号的包络时频图如图14所示。图14中明显的存在IFCF趋势线及其倍频。利用峰值搜索算法提取的IFCF趋势线如图12中上部正方形所示。作为对比，图12给出了实测的轴承转频(下部实线)。不难发现提取的IFCF趋势线与实测转频具体相同的变化趋势。根据IFCF趋势线和转频计算得到的故障特征系数如图9中十字符号所示。与外圈故障特征系数相重合，平均误差为0.23%，说明该轴承存在外圈故障。根据以上分析可知本文算法可以有效的识别出滚动轴承的故障类型。

其次，利用内圈存在故障的滚动轴承振动信号进一步验证本文算法的有效性，转速变化范围为2 460～1 380 r/min。图15为内圈故障轴承振动信号的时域波形图，利用本文算法对该振动信号进行处理，得到的故障特征系数如图16中的十字符号。通过模板匹配不难发现通过本文算法得到的故障特征系数与模板中的内圈故障特征系数几乎相重合，平均误差为0.14%，从而说明该轴承存在内圈故障。

最后，利用本文算法对健康轴承振动信号进行分析，其中图17为健康轴承信号的时域波形图。图18为模板匹配结果。计算得到的故障特征系数凌乱分布，与实际的故障特征系数没有任何关系。

图13 外圈故障轴承振动信号的快速谱峭度图Fig.13 Kurtogram of the outer faulty bearing signal

图14 滤波信号的包络时频图Fig.14 Envelope time-frequency representation of filtered signal

图15 内圈故障轴承振动信号的时域波形图Fig.15 Inner faulty bearing signal in time domain

图16 内圈故障轴承模板匹配结果Fig.16 Inner faulty bearing matching result

图17 健康轴承振动信号的时域波形图Fig.17 Healthy bearing signal in time domain

图18 健康轴承模板匹配结果Fig.18 Healthy bearing matching result

5 结 论

本文提出了基于故障特征系数模板的滚动轴承故障诊断方法。该算法可以有效的识别变转速条件下滚动轴承是否发生故障以及故障类型。

(1) 滚动轴承内圈、外圈以及滚珠故障时所对应的故障特征频率与轴承转频具有固定的比例关系，而且故障轴承振动信号的包络时频图中IFCF趋势线具有明显的幅值优势，易于提取。

(2) 快速谱峭度算法可以有效的确定轴承振动信号的最优滤波参数，进而得到包含故障冲击成分的滤波信号。

(3) 本文算法无需阶比跟踪即可对变转速条件下的滚动轴承进行故障诊断，步骤更为简练，也进一步丰富了现有的故障诊断理论。

[ 1 ] CHENG W D, GAO R X, WANG J J, et al. Envelope deformation in computed order tracking and error in order analysis[J]. Mechanical System and Signal Processing, 2014, 48(1/2): 92-102.

[ 2 ] ZIMROZ R, BARTELMUS W, BARSZCZ T, et al. Diagnostics of bearings in presence of strong operating conditions non-stationarity — a procedure of load-dependent features processing with application to wind turbines bearings[J]. Mechanical System and Signal Processing, 2014, 46(1): 16-27.

[ 3 ] FYFE K R, MUNCK E D S. Analysis of computed order tracking[J]. Mechanical System and Signal Processing，1997, 11(2):187-205.

[ 4 ] RANDALL R B, ANTONI J. Rolling element bearing diagnostics — a tutorial[J].Mechanical Systems and Signal Processing, 2011, 25(2):485-520.

[ 5 ] 刘亭伟，郭瑜，李斌，等. 基于谱峭度的滚动轴承故障包络阶比跟踪分析[J].振动与冲击，2012, 31(17):149-164. LIU Tingwei, GUO Yu, LI Bin, et al. Envelope order tracking analysis for rolling element bearing faults based on spectral kurtosis[J]. Journal of Vibration and Shock, 2012, 31(17):149-164.

[ 6 ] 郭瑜，秦树人.旋转机械非稳定信号的伪转速跟踪阶比分析[J].振动与冲击,2004,23(1):61-69. GUO Yu, QIN Shuren. Pseudo-speed tracking order analysis for non-stationary vibration signals of rotating machinery[J]. Journal of Vibration and Shock, 2004,23(1):61-69.

[ 7 ] SAAVEDRA P N, RODRIGUEZ C G. Accurate assessment of computed order tracking [J]. Shock and Vibration, 2006, 13 (1):13-21.

[ 8 ] LI Z, WU Z, HE Y, et al. Hidden Markov model-based fault diagnostics method in speed-up and speed-down process for rotating machinery[J]. Mechanical Systems and Signal Processing, 2005, 19(2):329-339.

[ 9 ] MELTZER G, DIEN N P. Fault diagnosis in gears operating under non-stationary rotational speed using polar wavelet amplitude maps[J]. Mechanical Systems and Signal Processing, 2004, 18(5):985-992.

[10] LI C, LIANG M. A generalized synchrosqueezing transform for enhancing signal time-frequency representation[J]. Signal Processing, 2012, 92(9):2264-2274.

[11] ANTONI J. Cyclic spectral analysis of rolling-element bearing signals: facts and fictions[J].Journal of Sound and Vibration, 2007, 304(3/4/5):497-529.

Rolling element bearing fault diagnosis based on fault characteristic coefficienttemplate under time-varying rotating speed

CHENG Weidong, ZHAO Dezun, LIU Dongdong

(School of Mechanical, Electronic and Control Engineering, Beijng Jiaotong University, Beijing 100044, China)

Compared to constant rotating speed condition, it is more difficult to diagnose bearing faults under time-varying rotating speed. On the other hand, the common method under time-varying rotating speed, i.e., the order ratio analysis has some problems, such as, larger error and lower computing efficiency. In order to solve these problems, the method of rolling element bearing fault diagnosis based on fault characteristic coefficient template was proposed nere. It consisted of six main steps: ① a template was set up based on the fault characteristic coefficient computed according to geometric parameters of the target bearing; ② the bearing vibration signals were filtered via the fast spectral kurtosis filtering method; ③ the envelope time-frequency figures of the filtered signals were computed based on Hilbert transformation and short time Fourier transformation; ④ the instantaneous fault characteristic frequency trend curve was extracted from the envelope time-frequency figures with the spectral peak search algorithm; ⑤ the speed curve of the bearing was calculated with the rotating speed pulse signals; ⑥ the ratio of the instantaneous fault characteristic frequency to the instantaneous rotating frequency is the instantaneous fault characteristic coefficient for diagnosing the fault type of the bearing. The effectiveness of the proposed method was validated using both simulated and actual measured rolling element bearing faulty vibration signals.

time-varying rotating speed; rolling element bearing; fault characteristic coefficient template

国家自然科学基金资助项目(51275030)

2015-12-07 修改稿收到日期：2016-02-29

程卫东 男，教授，1967年生 E-mail：wdcheng@bjtu.edu.cn

TH113；TH165

A

10.13465/j.cnki.jvs.2017.07.019