滚动轴承自适应特征提取的包络谱多点峭度多级降噪方法

时间：2024-07-28

张龙蔡秉桓熊国良王朝兵,2 胡俊锋

1.华东交通大学机电与车辆工程学院，南昌，3300132.中车戚墅堰机车有限公司，常州，2130113.中国铁路南昌局集团有限公司科学技术研究所，南昌，330002

0 引言

滚动轴承广泛应用于机械、交通、航空航天等领域，是旋转机械的重要零部件之一[1]，但其恶劣的工作环境易导致故障发生，若未及时发现则可能引起严重后果，因此，准确判断滚动轴承健康状态对提高机械设备的可靠性、可用性和保障设备安全运行至关重要[2-3]。

能够表征轴承局部故障的周期性瞬态脉冲在信号采集过程中会被传递路径、强背景噪声以及高幅值偶然性冲击所干扰，导致周期性瞬态冲击特征提取困难[4]，因此，有效的特征提取方法对轴承故障检测至关重要。经验模态分解(empirical mode decomposition, EMD)作为一种频带分解方法可以有效地将信号分解为若干近似信号和细节信号。苏文胜等[5]利用EMD对原始信号进行预处理降噪后，采用谱峭度方法进行共振滤波，然而谱峭度方法所采用的峭度指标不能考虑故障特征的周期特性，同时EMD也存在分解过程缺乏理论基础、端点效应和模态混叠等问题。SELESNICK[6]提出了一种信号共振稀疏分解(resonanced-based signal sparse decomposition,RSSD)方法，该方法不同于基于频带或尺度分解的线性方法，通过信号共振属性不同(即品质因子大小差异)将故障信号分解成由持续振荡成分组成的高品质因子共振分量和由瞬态冲击成分组成的低品质因子共振分量。该方法因其在脉冲提取方面的独特优势，迅速被用于旋转机械故障诊断领域。WANG等[7]以峭度最大为指标对原始信号进行总体经验模态分解，然后采用RSSD对选取的最优固有模态分量(intrinsic mode functions， IMF)进行分解，完成故障诊断。LI等[8]首先以峭度为优化指标，利用本征特征尺度分解(intrinsic characteristic-scale decomposition，ICD)对原始信号进行预处理，进一步根据特征频域比选择RSSD最佳分量进行分析。在上述研究中，当涉及信号处理参数优化或分量选择时，大多研究采用的优化指标是易受偶然性干扰冲击影响的峭度，该指标未考虑滚动轴承故障冲击特征周期性发生的特点，从而诊断效果不佳。

传感器采集到的振动信号通常是故障瞬态脉冲和背景噪声与轴承系统脉冲响应函数(传递路径)之间的卷积结果。为了消除传递路径对轴承故障特征的干扰，一系列解卷积方法应运而生[9-11]。BARSZCZ等[12]采用最小熵解卷积(minimum entropy deconvolution, MED)方法消除传递路径影响，并通过包络谱分析实现故障诊断。RICCI等[13]考虑单一信号处理方法难以有效去除强干扰噪声的影响，提出将EMD与MED复合使用，以提高滚动轴承故障诊断效果。上述方法虽然能够起到一定的效果，但MED以最小熵(实际是最大峭度)作为衡量指标，存在与谱峭度方法同样的缺陷，即没有考虑故障冲击特征的周期特性，因此更倾向于提取单瞬态脉冲。针对上述问题，MCDONALD等[11]提出了多点最优最小熵解卷积(multipoint optimal minimum entropy deconvolution adjusted，MOMEDA)方法。SHANG等[14]将MOMEDA与改进经验小波变换(improved empirical wavelet transform，IEWT)相结合对振动信号进行特征提取，进一步强化了单一MOMEDA特征提取效果。与MED、最大相关峭度解卷积(maximum correlated kurtosis deconvolution, MCKD)方法不同的是，MOMEDA方法通过一个目标向量定义脉冲所在位置及权重，由此提取旋转机械故障振动信号中每一个周期内出现的冲击脉冲，因此可以有效考虑故障冲击特征的周期性，无需迭代计算且可以使用非整数故障周期，这样可避免重采样,同时也节省了运算时间。

本文从轴承故障信号包络谱中故障特征频率及其倍频成分表现的典型周期性特点出发，提出一种能够免疫偶然性冲击和强背景噪声的周期性冲击特征度量新指标——包络谱多点峭度(envelope spectral multipoint kurtosis，ESMK)。首先采用MOMEDA削弱信号中传递路径的影响，然后针对现有RSSD方法中品质因子依赖人工经验选择的不足，以标准化的ESMK为适应度函数，采用粒子群算法(PSO)对品质因子进行优化，进而得到一种自适应稀疏分解方法PSO-RSSD用于瞬态冲击特征提取，有效消除偶然性冲击等噪声干扰，进一步通过包络谱实现滚动轴承的故障诊断。

1 基础理论介绍

1.1 MOMEDA方法

MOMEDA方法[11]作为一种无需迭代处理的解卷积算法，能够有效识别连续脉冲，实现故障瞬态冲击成分增强的目的。假设采集到的振动信号如下：

x(n)=h*y+e

(1)

式中,e为高斯白噪声；y为故障脉冲序列；h为轴承系统传输路径的响应；*表示卷积；x为传感器采集到的振动信号。

MOMEDA算法的本质是寻找一个最优的FIR(finite impulse response)滤波器，通过对输出信号x进行解卷积处理，从而恢复故障瞬态冲击信号y，即

(2)

式中，f为滤波器系数，f=(f1,f2,…，fk,…，fL)；N为采样个数。

MOMEDA求解最优滤波器的过程是通过对多重D范数(multi D-norm，MDN)M(y,t)求最大值，定义如下：

(3)

(4)

其中,t为确定脉冲位置和权重的常数矢量。在不同采样频率下，不需要重新采样就可以提取不同的故障周期，因此，目标向量t可以用来确定脉冲信号的位置。

计算滤波器系数f的导数：

(5)

由于

Mk=(xk+L-1,xk+L-2，…，xk)T

(6)

则令式(5)等于0，进一步化简得

(7)

(8)

MOMEDA的最优滤波器和输出结果可简化为

(9)

1.2 包络谱多点峭度

由循环平稳理论可知，滚动轴承故障信号的瞬态冲击具有典型的二阶循环平稳特征，即其时域信号并非严格的周期信号，其频谱中的瞬时能量成分才是周期出现的[15]。基于此，可通过计算信号包络谱的多点峭度考察瞬时能量的波动情况。滤波后的信号可表示为

(10)

可以根据下式获得滤波信号的包络谱S(n)：

(11)

滤波后信号的包络谱表示如下：

(12)

m=0,1,…,N-1

将ESMK作为轴承局部故障周期性冲击特征的度量指标，经过标准化的ESMK可以定义为

(13)

当信号中只包含一个瞬态冲击成分时，信号具有较大的峭度值；当信号中包含多个周期分布的瞬态冲击序列时，信号的峭度值反而更小。而经过标准化的ESMK有效改善了峭度指标易受高幅值偶然性冲击所干扰的缺陷，可以在强干扰下有效识别故障周期性瞬态冲击。

1.3 共振稀疏分解(RSSD)

SELESNICK[6]提出了以可调品质因子小波变换(TQWT)[16]为基础的信号共振稀疏非线性分解方法，不再利用传统方法将信号划分到不同频带，而是根据信号中谐波信号和冲击信号所对应的品质因子的差异，采用形态分量分析(morphological component analysis，MCA)[17]方法将拥有不同振荡特性的成分进行分离，得到含有冲击信号的低品质共振分量以及含有谐波分量的高品质共振分量。品质因子Q定义为中心频率fw与带宽BW的比值：

(14)

在共振稀疏分解中，品质因子Q定义了信号的共振属性。瞬态冲击为单振荡信号时，其品质因子越小，信号就拥有越少的共振属性；谐波部分为持续振荡信号时，其品质因子越大，信号就拥有越多的共振属性。共振稀疏分解方法依据信号的共振属性，使用一个双带通滤波器组实现分解过程，如图1所示。通过TQWT分别获得高品质因子与低品质因子的基函数库，并通过迭代计算得到相应的变换系数。低通滤波器组尺度参数和高通滤波器组尺度参数的表达式分别如下：

(15)

式中，r为冗余度。

由式(15)可知，低通滤波器和高通滤波器是由尺度因子决定的，尺度因子又是由品质因子和冗余度决定的，因此不同的品质因子和冗余度决定了不同的小波滤波器。

图1 两通道滤波器组Fig.1 The two-channel filter bank

x=x1+x2+e

(16)

式中，x1主要为谐波信号；x2主要为故障冲击信号；e为背景噪声。

利用x1和x2具有不同的品质因子，采用形态学分量分析对原始信号x进行非线性分解，其稀疏表示过程实际上是一个最小化问题，假设S1和S2分别表示高品质因子和低品质因子的滤波器组,根据形态分量分析构建稀疏分解目标函数如下：

argmin(W1,W2)

(17)

式中,W1、W2分别为x1和x2在S1和S2基函数库下的变换系数；γ1、γ2为正则化参数；argmin(W1,W2)表示式(17)达到最小时W1和W2的取值。

稀疏分解可以看成寻找最优的变换系数W1和W2的过程，采用分裂增广拉格朗日收缩算法进行优化，使得目标函数最小。分解得到的高共振分量与低共振分量分别为

(18)

1.4 粒子群优化算法(PSO)

粒子群算法(PSO)是KENNEDY[18]于1995年提出的一种优化算法，具有较好的全局寻优能力。PSO算法首先在解空间中对粒子进行初始化，每个粒子都表示优化问题的一个潜在最优解，粒子的特征由位置、速度和适应度三个指标来表示。根据每个粒子的位置和速度计算当前粒子的适应度值，由适应度值的大小来判断粒子的优劣。假设在一个D维空间中，由M个粒子组成的种群中的第i个粒子的速度和位置分别为Vi=(Vi,1,Vi,2,…,Vi,D)，Xi=(Xi,1,Xi,2,…,Xi,D)，个体极值Pi=(Pi,1,Pi,2,…,Pi,D)，群体极值Pg=(Pg,1,Pg,2,…,Pg,D)。粒子每次迭代更新其速度和位置变化公式分别为

(19)

(20)

式中，ω为惯性权重，用来衡量粒子在迭代过程中继承上一个粒子速度的能力；c1、c2为加速度因子；u为代数；r1、r2为随机分布在区间[0,1]的常数。

2 本文提出的轴承故障诊断算法

能够表征轴承故障的周期性瞬态脉冲在信号采集过程中经常被传递路径、强背景噪声以及高幅值偶然性冲击所干扰，导致周期性瞬态冲击特征难以提取，同时共振稀疏分解参数品质因子的选择严重影响分解结果的好坏。

品质因子Q越大，对应的共振属性越多，反之对应的共振属性则越少，Q值过大或过小都会影响分解结果。在传统的RSSD方法中，Q值的选择严重依赖先验知识，而且目前大多数故障度量指标(以峭度为例)易受高幅值偶然性冲击等干扰的影响，难以有效度量循环瞬态冲击的周期性特征。因此本文针对上述问题，考虑到故障轴承包络谱中故障特征频率及其倍频成分表现的典型周期性特点以及PSO全局寻优的优异性能，以ESMK作为度量指标，提出了一种滚动轴承自适应特征提取的ESMK多级降噪方法。实现过程如下：①对原始信号进行MOMEDA解卷积预处理，以消除传递路径影响并初步突出冲击故障；②设置粒子群算法初始条件(加速度因子c1=c2=1.494，权重ω=1，种群规模M=20，最大迭代次数为50)，然后进一步设置高低品质因子Q1、Q2的变化范围分别为[1.0，3.0]、[4.0，12.0]，冗余因子r为3.0;③以低共振分量的ESMK最大为度量指标，采用PSO寻找最优高低品质因子Q1、Q2进行RSSD分解，得到包含瞬态冲击特征的最优低共振分量；④求最佳低共振瞬态分量的包络谱并与轴承的理论故障频率进行比较，完成故障诊断。具体流程如图2所示。

3 案例分析

3.1 内圈故障仿真信号分析

在轴承实际运行过程中，一旦轴承局部出现故障，将会产生周期性的瞬态冲击，然而，传感器所采集的轴承振动信号通常是传递路径对应的脉冲响应函数与轴承局部故障所产生的脉冲、背景噪声之间的卷积结果[19]。除此之外，振动信号还可能被外界人为因素或机械设备中其他部件导致的偶然性冲击所干扰，此类偶然性冲击在振动信号中往往表现为幅值大且不具有周期性，冲击幅值可以达到轴承故障冲击幅值的几倍，因此，偶然性冲击的峭度值会远大于轴承故障循环冲击的峭度值，从而影响最终的解调分析结果。传感器采集的振动信号x的具体产生原理如图3所示，本小节以仿真信号为例对进行分析验证。纯内圈故障瞬态冲击仿真信号如图3a所示，其中内圈故障特征频率为90 Hz，信号采样频率为20 480 Hz。

图2 本文所提方法流程图Fig.2 Flow chart of the proposed method

(a)故障瞬态冲击(b)噪声成分(d)加入噪声和干扰冲击经过传递路径后的内圈故障仿真信号(c)高幅值偶然性冲击 (e)包络谱图3 内圈故障仿真信号Fig.3 Simulated signals of inner race fault

为了使仿真信号更接近轴承实际运转时所产生的振动信号，加入幅值为0.4 m/s2的高斯随机噪声成分，结果如图3b所示。为了证明峭度易受高幅值偶然性冲击的影响，进一步在时域信号图3b中采样点1000～1060范围内人为仿真一段幅值为10 m/s2、频率为1500 Hz的正弦信号，其时域波形如图3c所示。最终与传递路径的脉冲响应函数卷积得到的信号如图3d所示。可见内圈故障瞬态冲击在高斯噪声和正弦冲击干扰下已无法从时域辨识。图3e为图3d的包络谱，从包络谱中不能找到有效的故障特征频率成分。

(a)MOMEDA预处理后的信号

(b)PSO-RSSD收敛曲线

(c)本文方法处理后的包络谱图4 内圈故障仿真信号分析结果Fig.4 Results on inner race fault simulated signals

内圈故障仿真信号的分析结果如图4所示。根据图2所示的本文方法流程，首先采用MOMEDA对原始信号进行预处理以消除传递路径的影响，结果如图4a所示。可见故障冲击成分得到初步增强，但冲击成分的周期性仍然不够明显，由于受噪声干扰严重，不足以判断轴承故障发生情况。进一步采用本文提出的PSO-RSSD方法进行瞬态冲击特征提取，设置PSO算法参数、高低品质因子的变化范围如上节所述。以提出的ESMK新指标作为PSO优化算法目标函数，寻优过程收敛曲线如图4b所示，可见适应度函数在第21次迭代后达到最大值。得到的最优高低品质因子分别为Q1=11.3488，Q2=1.0367。根据最优品质因子对MOMEDA预处理信号进行RSSD分解，得到的高低共振分量和残余分量如图5所示。图5a为信号的最优高共振分量，主要包含谐波成分；图5b为信号的最优低共振分量，其中冲击成分相对明显，主要包含瞬态冲击信息，可见瞬态冲击成分得到了明显增强。其包络谱如图4c所示，从包络谱中可以看到90.7 Hz的频率成分与内圈故障特征频率90 Hz非常接近，幅值明显且存在边频带，同时存在181.4 Hz、272.1 Hz等明显倍频成分，可以判断此时轴承发生了内圈故障。因此仿真信号分析结果验证了本文方法在轴承振动信号特征提取中的可行性。

(a)最优高品质共振分量

(b)最优低品质共振分量

(c)残余分量图5 内圈故障仿真信号稀疏分解结果Fig.5 Results of inner race fault simulated signalssparse decomposition

为了展示本文方法的优势所在，将MED作为预处理方法，后处理以峭度作为度量指标，依然采用PSO对RSSD品质因子进行优化，与本文方法进行对比分析，结果如图6所示。在最佳低共振分量时域波形(图6c)中没有发现周期性故障瞬态冲击，且包络谱(图6d)中没有任何明显的故障特征频率成分，完全被高幅值干扰冲击所影响，无法判断滚动轴承是否存在故障,验证了本文所提方法的有效性。

3.2 外圈故障实验信号分析

(a)收敛曲线

(b)预处理后的信号

(c)低品质因子信号

(d)对比方法处理后的包络谱图6 对比方法的内圈故障仿真信号分析结果Fig.6 Results on inner race fault simulated signalfor comparison

图7 自制滚动轴承故障实验台Fig.7 Test rig for bearing fault detection

为了模拟实际铁路、机械等大型装备中轴承局部故障，首先采用图7所示的自制转子-轴承故障模拟实验台所产生的外圈故障实验信号进行分析。该实验台可以模拟不同的滚动轴承和转子故障。实验台包括伺服电机及控制器、支撑轴承、圆盘、轴承座、加速度传感器、计算机以及数据采集卡等，振动信号由加速度传感器采集并保存在计算机中。实验所用轴承型号为N205，为了模拟轴承局部故障，采用线切割技术在轴承外圈加工出宽度为0.5 mm的凹槽。实验过程中电机转速为1000 r/min，加速度传感器安装在实验轴承座正上方，图7已标出。采样频率为12 000 Hz，根据轴承各元件故障频率计算公式计算得到此时实验轴承外圈故障特征频率为87.51 Hz。

传感器采集到的轴承外圈故障信号时域波形如图8a所示，是故障导致周期性瞬态冲击与轴承系统传递函数卷积结果。图8a中故障冲击成分较为明显，这是因为人为加工的凹槽较为标准，导致故障冲击幅值较大。为了使所采集的振动信号更接近现场实际情况，添加幅值为4 m/s2的高斯随机噪声后信号如图8b所示。同时为了模拟外界偶然性冲击干扰，在时域信号图8b中采样点2281～2360范围内人为添加一段幅值为60 m/s2的随机振动，结果如图8c所示。偶然性冲击幅值远大于轴承故障冲击幅值，轴承外圈故障冲击特征在此干扰冲击和高斯噪声淹没下已经无法清晰识别。

(a)采集的原始加速度信号

(b)加噪后的加速度信号

(c)加入噪声和干扰冲击后的加速度信号图8 外圈故障仿真信号Fig.8 Simulated signals of outer race fault

(a)MOMEDA预处理后信号

(b)PSO-RSSD收敛曲线

(c)本文方法处理后的包络谱图9 外圈故障仿真信号分析结果Fig.9 Results on outer race fault using theproposed method

采用本文方法所得分析结果如图9所示。首先对原始信号进行MOMEDAD预处理，结果如图9a所示，与图8c相比，故障冲击得到初步增强。进一步采用PSO-RSSD方法进行瞬态冲击二次增强，设置PSO算法参数以及高低品质因子的取值范围后，寻优过程收敛曲线如图9b所示，可见适应度函数在第31次迭代后达到最大值。得到的最优高低品质因子分别为Q1=4.76，Q2=2.94。根据最优品质因子对预处理信号进行RSSD分解，得到的高低共振分量、残余分量如图10所示。图10a所示为信号的最优高品质共振分量，主要包含谐波成分；图10b所示为信号的最优低品质共振分量，瞬态冲击成分得到了明显增强。其包络谱如图9c所示，从包络谱中可以看到88.01 Hz的频率成分与外圈故障特征频率88 Hz非常接近，幅值突出且存在176 Hz、263 Hz等明显倍频成分，可以判断此时轴承发生了外圈故障。因此外圈故障实验数据分析结果验证了本文方法在高幅值偶然性冲击等强干扰下提取周期性故障特征的可行性。

(a)最优高品质共振分量

(b)最优低品质共振分量

(c)残余分量图10 外圈故障仿真信号稀疏分解结果Fig.10 Results of outer race fault simulated signalssparse decomposition

作为对比，以MED作为预处理方法，后处理以峭度作为优化指标，采用PSO对RSSD品质因子进行优化，结果如图11所示。最优低品质共振分量时域波形如图11c所示，其包络谱(图11d)中存在大量噪声，与本文方法相比，其特征提取效果不佳。

3.3 内圈故障实验信号分析

滚动轴承内圈故障信号同样来源于图7所示的自制转子-轴承故障模拟实验台，内圈故障如图12所示，通过公式计算可知轴承内圈故障特征频率fc=129.15 Hz。

(a)收敛曲线

(b)预处理后的信号

(c)最优低品质共振分量

(d)对比方法处理后的包络谱图11 对比方法的外圈故障信号分析结果Fig.11 Results on outer race fault for comparison

图12 内圈故障轴承实物图Fig.12 The physical drawing of inner ringfault bearing

图13a是由原始数据中随机截取的一段含有12 000个采样点的信号，可以看到原始信号杂乱无章，冲击成分完全被噪声覆盖,图13b为其包络谱。由图13b可见，由于噪声干扰严重，无法清晰地找到内圈故障特征频率。采用本文方法所得分析结果如图14所示。首先采用MOMEDA对原始信号进行预处理，结果如图14a所示，可见故障冲击成分得到初步增强，但不足以判断轴承故障发生情况。进一步采用PSO-RSSD方法进行瞬态冲击特征二次增强，以低共振分量包络谱的多点峭度(ESMK)作为目标函数，寻优过程收敛曲线如图14b所示，可见适应度函数在第3次迭代后达到最大值，得到的最优高低品质因子分别为Q1=8.04，Q2=2.34。最优低共振分量时域波形如图14c所示，冲击成分得到了明显的增强，其包络谱如图14d所示。从包络谱中可以看到130.4 Hz的频率成分与内圈故障特征频率fc=129.15 Hz吻合，幅值明显且存在260.9 Hz、391.3 Hz等显著倍频成分，可以判断此时轴承发生了内圈故障。内圈故障实验数据分析结果验证了本文方法在轴承振动信号特征提取中的有效性。

(a)原始信号

(b)包络谱图13 内圈故障信号及其包络谱Fig.13 The signal and its envelope spectrum ofinner race fault

(a)MOMEDA预处理后信号

(b)PSO-RSSD收敛曲线

(c)最优低共振分量

(d)本文方法处理后的包络谱图14 本文所提方法分析结果Fig.14 Results on inner race fault using theproposed method