CN113567127A - 一种基于时频特征分离的滚动轴承退化指标提取方法 - Google Patents

Abstract

一种基于时频特征分离的滚动轴承退化指标提取方法，首先对轴承原始信号计算其时频谱；其次，通过对时频矩阵中同一时刻极大值点和同一频率极大值点个数的确定，结合时频矩阵行累积结果确定谐波分布位置和范围，保留谐波分布区域，滤除其它区域，将处理后的时频矩阵逆短时傅里叶变换至时域以实现谐波成分的提取；从原始信号中减去谐波成分后，再阈值降噪进行故障冲击成分的提取，最后对提取的故障冲击成分计算其相对均方根值作为滚动轴承退化指标；本发明基于时频特征预先滤除时频谱中突出的谐波成分以增强故障冲击成分，避免了阈值降噪提取冲击成分时产生严重失真，对提取的故障冲击成分计算滚动轴承退化指标，该指标可有效表征轴承退化过程。

Description

一种基于时频特征分离的滚动轴承退化指标提取方法

技术领域

本发明属于滚动轴承故障诊断领域，具体涉及一种基于时频特征分离的滚动轴承退化指标提取方法。

背景技术

滚动轴承作为旋转设备的关键零部件，其健康状态对整个设备的安全运行起着决定性的作用，因此实现滚动轴承运行状态和退化过程的监测显得十分必要。

滚动轴承故障往往表现为轴承工作表面的点蚀、裂纹和剥落等局部损伤，在工作过程中，这些局部损伤被撞击，便会产生周期性的冲击成分。然而，由于工业传动系统的复杂性和运行工况的时变性导致原始信号中包含大量的干扰信号，因此，如何从复杂分量中准确提取冲击成分成为描述滚动轴承退化过程的研究重点。原始信号受电磁频率、转频相关频率等成分的影响，信号中包含大量突出谐波成分，冲击成分被淹没，导致现有技术所提取的冲击成分存在严重失真，所提取的退化指标并不能表征轴承真实的退化过程。

发明内容

为了克服上述现有技术的缺点，本发明的目的在于提供了一种基于时频特征分离的滚动轴承退化指标提取方法，基于谐波、冲击、随机噪声信号成分在时频图中的特征，采用统计学思想分离出故障冲击成分，进而提取滚动轴承性能退化指标，该方法基于时频特征将时频谱中突出的谐波成分剔除，进而加强冲击特征，便于后续进行故障冲击成分的提取。

为了达到上述目的，本发明采取的技术方案为：

一种基于时频特征分离的滚动轴承退化指标提取方法，包括以下步骤：

步骤1：首先通过振动加速度传感器采集滚动轴承由正常到故障整个过程的原始信号，然后计算原始信号x(t)的短时傅里叶时频谱；在时频谱中，谐波成分表现为突出的横条纹，冲击成分表现为明亮的局部区域，而随机噪声表现为随机分布的噪点；

步骤2：谐波成分的识别与提取，时频谱实际为时频矩阵的可视化体现，且时频谱中亮点源自于时频矩阵中的极大值，因此对明亮横条纹的识别即转化为对时频矩阵中每一列中对应同一时刻极大值点和每一行中对应同一频率极大值点个数的确定，结合时频矩阵行累积确定谐波分布位置和范围，保留谐波分布区域，将其它区域能量尽可能降低，然后将处理后的时频矩阵逆短时傅里叶变换至时域以实现谐波成分的提取；

步骤3：从原始信号中减去步骤2提取的谐波成分后，进行冲击成分的提取，冲击成分在时频谱中表现为明亮的局部区域，因此，利用阈值降噪去除时频谱中背景噪声，保留故障冲击区域，将处理后的矩阵返回至时域即得冲击成分；

步骤4：对于步骤3分离出的故障冲击成分，分别计算其相对均方根值作为滚动轴承退化指标，并采用七点滑移平均法对相对均方根趋势进行处理。

所述步骤1中短时傅里叶变换表达式如下：

其中，t表示时间，τ表示时延，f表示频率，x(t)表示输入信号，h(t)为高斯窗函数，其表达式如下：

所述步骤2谐波提取方法具体为：

2.1)针对原始信号经过短时傅里叶变换得到的时频矩阵S1，首先采用1σ准则对时频矩阵S1每一行分别进行处理，即将每行数据中1σ区间外的数据均用该行最小值S1_min(i)替代，其表达式如下：

其中，1≤i≤P，1≤j≤Q，i，j均取整数，P、Q分别为时频矩阵S1的行数和列数；S1_min(i)为时频矩阵S1第i行绝对值的最小值，μ_i为时频矩阵S1第i行取绝对值后的均值，σ_i为时频矩阵S1第i行取绝对值后的标准差，它们的表达式如下：

S1_min(i)＝min{|S1(i,1)|,|S1(i,2)|,···,|S1(i,j)|,···,||1(i,Q)||

2.2)考虑同一时间即各列中的极大值点，并进行标记，得到与时频矩阵S11同长同宽的标记矩阵Slabel，其表达式如下：

其中S11(i,j)表示步骤2.1)中时频矩阵S11对应第i行、第j列的元素；

2.3)统计标记矩阵Slabel中同一频率下即各行中出现的极大值点个数Num(i)，i表示为时频矩阵第i行，其表达式如下：

2.4)计算时频矩阵S11的行累积结果Product(i)，并标记行累积结果中出现的极大值点为1，标记极小值点为-1，其它点为0，得到标记向量Plabel，上述两者的表达式如下：

Product(i)＝S11(i,1)×S11(i,2)×···×S11(i,j)×···×S11(i,Q1)

其中，1≤i≤P1，1≤j≤Q1，i，j均取整数，P1、Q1分别为步骤2.1)中时频矩阵S11的行数和列数；

2.5)当时频矩阵S11第i行累积结果为极大值点即Plabel(i)值为1，且Num(i)大于等于该行元素个数一半时，即认为该行对应频率为谐波成分；通过在行累积结果中搜寻谐波频率左右相邻的极小值点，将其作为谐波在时频谱内分布范围的上下限，保持时频矩阵S11中谐波分布区域的值不变，将其它区域均置为矩阵最小值，得到时频矩阵NS1；

2.6)将时频矩阵NS1逆短时傅里叶变换至时域即得滤出的谐波成分h(t)。

所述步骤3具体为：

3.1)从原始信号x(t)中减去步骤2得到的谐波成分h(t)，对剩余信号y(t)重新计算其时频矩阵S2(T,F)；

3.2)对时频矩阵S2(T,F)划分阈值区间，依次选取阈值进行阈值降噪，阈值区间[a,b]、步长Δd、阈值TH(k)的计算如下：

S_min＝min(min(|S2(T,F)|))，S_max＝max(max(|S2(T,F)|))

a＝S_min+0.1*(S_max-S_min)，b＝S_min+0.9*(S_max-S_min)

Δd＝(S_max-S_min)/100，TH(k)＝a+k*Δd

其中，S_min、S_max分别表示步骤3.1)中时频矩阵S2取绝对值后的最大值和最小值，a表示阈值区间的下限，b表示阈值区间的下限，k表示第k次阈值降噪，TH(k)为第k次进行阈值降噪所取的阈值；

3.3)将时频矩阵S2(T,F)中小于阈值TH(k)的位置用S_min替换得到新的时频矩阵NS2，并将时频矩阵NS2逆短时傅里叶变换得到滤波后的时域信号z^k(t)，计算时域信号z^k(t)与z^k-1(t)的均方误差E(k)如下：

其中N为信号长度，z^k(t_n)表示第k次阈值降噪得到的时域信号的第n个元素；

3.4)当第k次均方误差为第1次均方误差的十分之一时，满足迭代终止条件，输出z^k(t)，即为分离出的故障冲击信号。

所述步骤4具体为：

4.1)对步骤3得到的故障冲击信号z^k(t)计算其均方根值Z_rms，表达如下：

4.2)选取4.1)中计算得到的滚动轴承振动平稳期内的前M组处理后的信号的均方根值，将其均值作为标准值，随后计算原均方根值与标准值之比，进而得到相对均方根值Z_rrms，表达式如下：

其中l＝1,2,3,···,NF为组号，NF为所采集的原始信号的组数；

4.3)采用七点滑移平均对步骤4.2)得到的相对均方根趋势进行平滑处理，其表达式如下：

其中，Z_srrms为相对均方根Z_rrms经滑移平均得到的平滑相对均方根；以Z_srrms为纵坐标，l为横坐标绘制滚动轴承退化趋势。

本发明的有益效果为：

滚动轴承的振动往往具有非平稳性，通过有效的时频分析能够准确、全面挖掘蕴含其中的故障信息。本发明是基于滚动轴承原始信号中各信号成分在时频谱中的不同表现形式，依次实现各信号成分的分离；首先是谐波成分的分离，考虑其在时频谱中表现为明亮的横条纹，因此通过确定谐波成分在时频谱的中心位置和分布范围，将其它区域能量缩小，并将处理后的时频矩阵返回至时域以实现谐波成分的分离；从原始信号中分离出谐波成分后，突出了故障冲击成分，为后续阈值降噪处理提供支撑。其后为故障冲击成分的分离，考虑其在时频谱中表现为明亮的局部区域，对时频矩阵元素的分布区间进行划分，从小到大依次选取阈值进行降噪，当阈值降噪后的结果满足迭代收敛条件时输出当前阈值降噪后的时频矩阵，将其返回至时域即得分离出的故障冲击信号；最后通过对分离出的故障冲击信号计算滚动轴承退化指标，该指标可有效表征轴承退化过程。

附图说明

图1为本发明的流程图。

图2为XJTU-SY(西安交通大学雷亚国老师课题组轴承数据)滚动轴承外圈故障数据(Bearing1_1)中第10组原始信号时域波形以及短时傅里叶时频谱。

图3为XJTU-SY滚动轴承外圈故障数据(Bearing1_1)中第80组原始信号时域波形以及短时傅里叶时频谱。

图4为第10组原始信号滤除谐波前频谱图以及滤除谐波后频谱图。

图5为第80组原始信号滤除谐波前频谱图以及滤除谐波后频谱图。

图6为第10组原始信号经过滤除谐波以及阈值降噪后得到的时域波形。

图7为第80组原始信号经过滤除谐波以及阈值降噪后得到的时域波形。

图8为XJTU-SY滚动轴承外圈故障数据(Bearing1_1)经过本发明方法处理后获取的滚动轴承退化趋势，并与用原始信号获取的轴承退化趋势做对比。

具体实施方式

下面结合实施例和附图对本发明进一步详细说明。

如图1所示，一种基于时频特征分析的滚动轴承退化趋势获取方法，采用XJTU-SY轴承数据中的轴承外圈故障的原始信号，包括以下步骤：

步骤1：首先通过振动加速度传感器间隔采集滚动轴承由正常到故障全过程的原始信号，每隔一分钟采一次数据，每次采样时长为1.28s，采样频率为25.6KHz，共采集123组；其中第10组原始信号(出现故障前)和第80组原始信号(出现故障后)前0.1s时域波形以及短时傅里叶变换时频谱如图2、图3所示，在时频谱中，谐波成分表现为突出的横条纹，冲击成分表现为明亮的局部区域，而随机噪声表现为随机分布的噪点；短时傅里叶变换窗函数选择高斯窗，窗长为256，短时傅里叶变换以及高斯窗函数表达式如下：

其中，t表示时间，τ表示时延，f表示频率，x(t)表示输入信号，h(t)为高斯窗函数，其表达式如下：

参照图2、图3，发现出现故障前其原始信号时域波形无明显冲击，其时频谱上主要为谐波成分(明亮的横条纹)和随机成分(随机分布的亮点)；出现故障后，其原始信号时域波形存在明显的周期性冲击，其时频谱中也可看到与之对应的周期性分布的亮点区域，同时也存在与谐波成分对应的明亮横条纹；

步骤2：谐波成分的识别与提取，鉴于谐波成分的规律性较强，且在时频谱中较为突出，优先考虑提取谐波成分；时频谱实际为时频矩阵的可视化体现，且时频谱中亮点源自于时频矩阵中的极大值，因此对明亮横条纹的识别即转化为对时频矩阵中每一列中对应同一时刻极大值点和每一行中对应同一频率极大值点个数的确定，前者考虑谐波在时频谱中突出的特征，后者考虑了谐波在时频谱中的条纹特征，结合时频矩阵行累积结果确定谐波分布位置和范围，保留谐波分布区域，将其它区域能量尽可能降低，然后将处理后的时频矩阵逆短时傅里叶变换至时域以实现谐波成分的提取；

谐波提取方法如下：

2.1)针对原始信号经过短时傅里叶变换得到的时频矩阵S1，首先采用1σ准则对时频矩阵S1每一行分别进行处理，以剔除时频谱中能量较高的冲击亮点，避免将较密的冲击区域误认为是谐波区域，即将每行数据中1σ区间外的数据均用该行最小值S1_min(i)替代，其表达式如下：

其中，1≤i≤P，1≤j≤Q，i，j均取整数，P、Q分别为时频矩阵S1的行数和列数；S1_min(i)为时频矩阵S1第i行绝对值的最小值，μ_i为时频矩阵S1第i行取绝对值后的均值，σ_i为时频矩阵S1第i行取绝对值后的标准差，它们的表达式如下：

S1_min(i)＝min{|S1(i,1)|,|S1(i,2)|,···,|S1(i,j)|,···,|S1(i,Q)|}

2.2)考虑同一时间即各列中的极大值点，并进行标记，得到与上述时频矩阵S11同长同宽的标记矩阵Slabel，其表达式如下：

其中S11(i,j)表示步骤2.1)中时频矩阵S11对应第i行、第j列的元素；

2.3)统计标记矩阵Slabel中同一频率下即各行中出现的极大值点个数Num(i)，i表示为时频矩阵第i行，其表达式如下：

2.4)计算时频矩阵S11的各行元素累积结果Product(i)，并标记行累积结果中出现的极大值点为1，标记极小值点为-1，其它点为0，得到标记向量Plabel，上述两者的表达式如下：

Product(i)＝S11(i,1)×S11(i,2)×···×S11(i,j)×···×S11(i,Q1)

其中，1≤i≤P1，1≤j≤Q1，i，j均取整数，P1、Q1分别为步骤2.1)时频矩阵S11的行数和列数；

2.5)当时频矩阵S11第i行累积结果为极大值点即Plabel(i)值为1，且Num(i)大于等于该行元素个数一半时，即认为该行对应频率为谐波成分；通过在行累积结果中搜寻谐波频率左右相邻的极小值点，将其作为谐波在时频谱内分布范围的上下限，保持时频矩阵S11中谐波分布区域的值不变，将其它区域均置为矩阵最小值，得到时频矩阵NS1；

2.6)将时频矩阵NS1逆短时傅里叶变换至时域即得滤出的谐波成分h(t)；第10组原始信号滤除谐波前后频谱如图4所示，第80组原始信号滤除谐波前后频谱如图5所示，对比滤波前后的频谱，发现原始信号中较为突出的谐波成分在经过基于时频谱的谐波剔除后得到很好的抑制，同时也保留了故障冲击的能量；

步骤3：从原始信号中滤除步骤2提取的谐波成分后，进行冲击成分的提取，冲击成分在时频谱中表现为明亮的局部区域，因此，利用阈值降噪去除时频谱中背景噪声，保留故障冲击区域，将处理后的矩阵返回至时域即得冲击成分；具体为：

3.1)从原始信号x(t)中减去步骤2得到的谐波成分h(t)，得到剩余信号y(t)，对剩余信号y(t)重新计算其时频谱S2(T,F)；

3.2)对时频矩阵S2(T,F)划分阈值区间，依次选取阈值进行阈值降噪，阈值区间[a,b]、步长Δd、阈值TH(k)的计算如下：

S_min＝min(min(|S2(T,F)|))，S_max＝max(max(|S2(T,F)|))

a＝S_min+0.1*(S_max-S_min)，b＝S_min+0.9*(S_max-S_min)

Δd＝(S_max-S_min)/100，TH(k)＝a+k*Δd

其中，S_min、S_max分别表示步骤3.1)中时频矩阵S2取绝对值后的最大值和最小值，a表示阈值区间的下限，b表示阈值区间的下限，k表示第k次阈值降噪，TH(k)为第k次进行阈值降噪所取的阈值；

3.3)将时频矩阵S2(T,F)中小于阈值TH(k)的位置用S_min替换得到新的时频矩阵NS2，并将时频矩阵NS2逆短时傅里叶变换得到滤波后的时域信号z^k(t)，计算z^k(t)与z^k-1(t)的均方误差E(k)如下：

其中N为信号长度，z^k(t_n)表示第k次阈值降噪得到的时域信号z^k(t)的第n个元素；

3.4)当第k次均方误差为第1次均方误差的十分之一时，满足迭代终止条件，输出z^k(t)，即为分离出的故障冲击信号；第10组原始信号经谐波剔除和阈值降噪后的时域波形如图6所示，第80组原始信号经谐波剔除和阈值降噪后的时域波形如图7所示；

步骤4：对于步骤3分离出的故障冲击信号，分别计算其相对均方根值，将其作为滚动轴承退化指标，并采用七点滑移平均对相对均方根趋势进行处理以消除随机振动特性的影响，具体为：

4.1)对步骤3得到的故障冲击信号z^k(t)计算其均方根值Z_rms，计算公式如下：

4.2)选取4.1)中计算得到的滚动轴承振动平稳期内的前三十组处理后的信号的均方根值，将其均值作为标准值，随后计算原均方根值与标准值之比，进而得到相对均方根值Z_rrms，计算公式如下：

其中l＝1,2,3,···,123为组号；

4.3)采用七点滑移平均对4.2)得到的相对均方根趋势进行平滑处理，其表达式如下：

其中Z_srrms为相对均方根经滑移平均处理后得到的平滑相对均方根；以Z_srrms为纵坐标，l为横坐标绘制滚动轴承相对均方根趋势。

本实施例绘制滚动轴原始信号及经本发明处理后的信号的相对均方根趋势如图8所示，图中共123组原始信号，1-72组为滚动轴承故障前正常阶段，该阶段原始信号以及处理后的信号的相对均方根值均为1左右。从第73组开始，轴承原始信号即处理后的信号的相对均方根值出现急剧爬升，到第100组时，原始信号相对均方根值达到4.6，处理后的信号相对均方根值达到28.9，后者爬升的幅度明显原高于前者，约为前者的7.7倍。可以看出本发明方法从原始信号中分离出故障冲击成分，通过计算分离出的故障冲击成分的相对均方根值获取滚动轴承的退化趋势，该指标有效的反映了滚动轴承的退化过程，且比从原始信号中提取的退化指标更加敏感。

Claims (5)

1.一种基于时频特征分离的滚动轴承退化指标提取方法，其特征在于，包括以下步骤：

步骤1：首先通过振动加速度传感器采集滚动轴承由正常到故障整个过程的原始信号，然后计算原始信号x(t)的短时傅里叶时频谱；在时频谱中，谐波成分表现为突出的横条纹，冲击成分表现为明亮的局部区域，而随机噪声表现为随机分布的噪点；

步骤2：谐波成分的识别与提取，时频谱实际为时频矩阵的可视化体现，且时频谱中亮点源自于时频矩阵中的极大值，因此对明亮横条纹的识别即转化为对时频矩阵中每一列中对应同一时刻极大值点和每一行中对应同一频率极大值点个数的确定，结合时频矩阵行累积确定谐波分布位置和范围，保留谐波分布区域，将其它区域能量尽可能降低，然后将处理后的时频矩阵逆短时傅里叶变换至时域以实现谐波成分的提取；

步骤3：从原始信号中滤除步骤2提取的谐波成分后，进行冲击成分的提取，冲击成分在时频谱中表现为明亮的局部区域，因此，利用阈值降噪去除时频谱中背景噪声，保留故障冲击区域，将处理后的矩阵返回至时域即得冲击成分；

步骤4：对于步骤3得到的故障冲击成分，分别计算其相对均方根值作为滚动轴承退化指标，并采用七点滑移平均法对相对均方根趋势进行处理。

2.根据权利要求1所述的一种基于时频特征分离的滚动轴承退化指标提取方法，其特征在于，所述步骤1中短时傅里叶变换表达式如下：

其中，t表示时间，τ表示时延，f表示频率，x(t)表示输入信号，h(t)为高斯窗函数，其表达式如下：

3.根据权利要求2所述的一种基于时频特征分离的滚动轴承退化指标提取方法，其特征在于，所述步骤2谐波提取方法具体为：

2.1)针对原始信号经过短时傅里叶变换得到的时频矩阵S1，首先采用1σ准则对时频矩阵S1每一行分别进行处理，即将每行数据中1σ区间外的数据均用该行最小值S1_min(i)替代，其表达式如下：

其中，1≤i≤P，1≤j≤Q，i，j均取整数，P、Q分别为时频矩阵S1的行数和列数；S1_min(i)为时频矩阵S1第i行绝对值的最小值，μ_i为时频矩阵S1第i行取绝对值后的均值，σ_i为时频矩阵S1第i行取绝对值后的标准差，它们的表达式如下：

S1_min(i)＝min{|S1(i,1)|,|S1(i,2)|,···,|S1(i,j)|,···,|S1(i,Q)|}

2.2)考虑同一时间即各列中的极大值点，并进行标记，得到与时频矩阵S11同长同宽的标记矩阵Slabel，其表达式如下：

其中S11(i,j)表示时频矩阵S11对应第i行、第j列的元素；

2.3)统计标记矩阵Slabel中同一频率下即各行中出现的极大值点个数Num(i)，i表示为时频矩阵第i行，其表达式如下：

2.4)计算时频矩阵S11的各行元素累积结果Product(i)，并标记行累积结果中出现的极大值点为1，标记极小值点为-1，其它点为0，得到标记向量Plabel，上述两者的表达式如下：

Product(i)＝S11(i,1)×S11(i,2)×···×S11(i,j)×···×S11(i,Q1)

其中，1≤i≤P1，1≤j≤Q1，i，j均取整数，P1、Q1分别为步骤2.1)中时频矩阵S11的行数和列数；

2.5)当时频矩阵S11第i行累积结果为极大值点即Plabel(i)值为1，且Num(i)大于等于该行元素个数一半时，即认为该行对应频率为谐波成分；通过在行累积结果中搜寻谐波频率左右相邻的极小值点，将其作为谐波在时频谱内分布范围的上下限，保持时频矩阵S11中谐波分布区域的值不变，将其它区域均置为矩阵最小值，得到时频矩阵NS1；

2.6)将时频矩阵NS1逆短时傅里叶变换至时域即得滤出的谐波成分h(t)。

4.根据权利要求3所述的一种基于时频特征分离的滚动轴承退化指标提取方法，其特征在于，所述步骤3具体为：

3.1)从原始信号x(t)中减去步骤2得到的谐波成分h(t)，对剩余信号y(t)重新计算其时频矩阵S2(T,F)；

3.2)对时频矩阵S2(T,F)划分阈值区间，依次选取阈值进行阈值降噪，阈值区间[a,b]、步长Δd、阈值TH(k)的计算如下：

S_min＝min(min(|S2(T,F)|))，S_max＝max(max(|S2(T,F)|))

a＝S_min+0.1*(S_max-S_min)，b＝S_min+0.9*(S_max-S_min)

Δd＝(S_max-S_min)/100，TH(k)＝a+k*Δd

其中，S_min、S_max分别表示步骤3.1)中时频矩阵S2取绝对值后的最大值和最小值，a表示阈值区间的下限，b表示阈值区间的下限，k表示第k次阈值降噪，TH(k)为第k次进行阈值降噪所取的阈值；

3.3)将时频矩阵S2(T,F)中小于阈值TH(k)的位置用S_min替换得到新的时频矩阵NS2，并将时频矩阵NS2逆短时傅里叶变换得到滤波后的时域信号z^k(t)，计算时域信号z^k(t)与z^k-1(t)的均方误差E(k)如下：

其中N为信号长度，z^k(t_n)表示第k次阈值降噪得到的时域信号的第n个元素；

3.4)当第k次均方误差为第1次均方误差的十分之一时，满足迭代终止条件，输出z^k(t)，即为分离出的故障冲击信号。

5.根据权利要求4所述的一种基于时频特征分离的滚动轴承退化指标提取方法，其特征在于，所述步骤4具体为：

4.1)对步骤3得到的故障冲击信号z^k(t)计算其均方根值，表达式如下：

4.2)选取4.1)中计算得到的滚动轴承振动平稳期内的前M组处理后的信号的均方根值，将其均值作为标准值，随后计算原均方根值与标准值之比，进而得到相对均方根值Z_rrms，表达式如下：

其中l＝1,2,3,···,NF为组号，NF为数据所包含的组数；

4.3)采用七点滑移平均对步骤4.2)得到的相对均方根趋势进行平滑处理，其表达式如下：

其中，Z_srrms为相对均方根Z_rrms经滑移平均得到的平滑相对均方根，以Z_srrms为纵坐标，l为横坐标绘制滚动轴承退化趋势。

Images