基于信号子空间的新型盲解卷积方法

2022-02-22欧曙东

振动与冲击 2022年3期

周涛，赵明，郭栋，欧曙东

(1．西安交通大学机械制造系统工程国家重点实验室，西安 710049；2．重庆理工大学汽车零部件先进制造技术教育部重点实验室，重庆 400054)

在各类回转机械中，例如涡轮发动机、离心机、风力发电机、汽车变速箱等，滚动轴承是极为关键的承载部件，轴承失效可能导致机械设备异常停机等严重后果，因此滚动轴承组件早期损伤诊断具有重要研究意义。在失效进程前期，故障冲击振动信号往往被背景噪声和其他干扰所淹没，加之复杂的传递路径的影响，失效特征难以准确提取，给滚动轴承故障诊断带来极大挑战。

解卷积方法可有效消除复杂传递路径的影响，增强故障冲击振动特征。Wiggins[1]在地震信号处理领域首次提出了最小熵解卷积(minimum entropy deconvolution,MED)方法，以滤波结果峭度最大化为目标，通过迭代算法从传递路径未知的测试信号中提取冲击信号。Cabrelli[2]使用D范数提出了一种类似MED的解卷积方法，称之为最优最小熵解卷积(optimal MED,OMED)方法。OMED与MED十分相似，但OMED方法不需要迭代过程。Endo等[3]将基于自回归模型的滤波器技术和MED相结合应用于齿轮故障诊断中，并证明了该方法的有效性。尽管MED结合其他方法在多种应用[4-8]中都取得了良好的效果，但MED仍存在一些主要缺点。MED在迭代过程中以峭度最大化为目标，这使迭代结果更加倾向于提取具有较大值的单脉冲信号，但在旋转机械中，更加希望提取连续的脉冲信号。针对MED的缺陷，Mcdonald等[9]引入了相关峭度，并提出了一种基于最大相关峭度的解卷积方法(maximum correlated kurtosis deconvolution,MCKD)。但MCKD在实际使用过程中仍存在以下不足，例如：移位系数难以优选；需要预知准确的故障特征周期；以及当故障特征周期为信号间隔的非整数倍时需要额外进行重采样。针对上述问题，Miao等[10]利用自相关函数实现了故障特征周期的自适应估计，提出了改进的MCKD方法。Mcdonald等[11]对Toeplitz矩阵进行了调整，消除了由数据不连续带来的虚假脉冲，并针对特征周期为信号间隔非整数倍时故障信号的提取问题，提出了多点最优最小熵解卷积调整(multi-point optimal minimum entropy deconvolution adjusted,MOMEDA)方法。Buzzoni等[12]提出了一种基于循环平稳性的盲解卷积方法，称之为最大二阶循环平稳盲解卷积(maximum second-order cyclostationarity blind deconvolution,CYCBD)方法，并给出了MCKD是基于循环平稳标准的论证。该方法可以使用测得的瞬时转速在转速非恒定时精确地提取故障脉冲。上述改进方法取得了广泛应用[13-15]，但在很多工业场景下，设备往往都会经历不同程度的转速变化，瞬时转速无法准确测量并且轴承在运行过程中存在一定程度的打滑现象，这使得故障特征先验知识不再准确有效，上述方法对故障冲击信号的提取效果不佳，为滚动轴承故障诊断技术的应用带来了新的挑战。因此急需一种无需预知准确故障特征周期，并对随机冲击具有良好鲁棒性的滚动轴承故障特征提取方法。同时目前现有的解卷积方法将测试信号重构相空间后并没有将相空间分解，这会使得滤波结果受到噪声影响以及其他信号干扰。

信号子空间技术可以将相空间分解为若干个子空间，通过筛选准则保留子空间中的若干个分量重构信号，抑制噪声对信号的影响。奇异值分解(singular value decomposition,SVD)方法作为一种信号子空间技术在旋转机械故障诊断中得到了广泛的应用[16-18]。赵学智等[19]提出一种小波-SVD差分谱组合模式，并成功应用于滚动轴承故障诊断。曾鸣等[20]提出了μ-SVD降噪算法用于强噪声情况下有效提取齿轮故障特征。Xu等[21]提出了奇异谱检测编码器信号中的微弱波动。Zhao等[22]引入周期调制强度评估分量信号的信息量，通过调整权重重构信号。现有可用的基于SVD的方法通常以奇异值阈值、重构信号能量或信息量进行筛选，很少以奇异值矩阵分量的信息量为指标建立筛选准则。奇异值矩阵分量不仅可以直接反应矩阵分量中故障信息量，无需重构分量信号进行筛选，而且矩阵分量是标准化的，消除了能量对信息量的影响。相空间重构和信号子空间技术对解卷积方法提供了新的方向，在迭代提取故障冲击信号之前筛选有效故障信号分量，加强故障信号的提取效果。

在此背景下，本文提出了一种基于信号子空间的新型盲解卷积方法，称之为最优稀疏矩阵分量解卷积(optimal sparse matrix component deconvolution,OSMCD)方法。该方法通过SVD分解原始信号构成的相空间，通过稀疏编码收缩抑制[23]子空间噪声，以脉冲稀疏指数为指标筛选有效子空间，最后迭代实现故障冲击的提取。通过变转速工况下的仿真信号和列车轴承试验结果表明，所提方法不仅可以有效消除随机冲击和噪声，而且在缺乏准确的故障特征周期情况下仍能有效实现故障脉冲的提取，为故障诊断技术提供了有价值的参考。

1 解卷积方法原理

轴承故障冲击为振荡衰减波形，由于振动传感器和损伤冲击位置的传递路径特性未知，振动采样信号的冲击特性并不明显。设多分量的测试信号x可表示为

x=hb*b+hr*r+hd*d+hn*n

(1)

式中:b,r,d,n分别为故障冲击信号，随机冲击信号，离散谐波信号以及高斯白噪声；hb,hr,hd,hn分别为上述信号各自的传递函数的FIR滤波器近似值。解卷积问题的目的在于采用一个逆滤波器从测试信号中消除传递路径的影响并提取故障冲击信号，基于不同的目标函数提出了不同的现有解卷积方法，如图1所示。通常，逆滤波器f可表示为长度为L的FIR滤波器，滤波信号y可以表示为

y=x*f=(hb*b+hr*r+hd*d+hn*n)*f≈b

(2)

现有的解卷积方法优化的目标函数有基于峭度、相关峭度、D范数、多D范数和循环平稳指标等，以MEDA(minimum entropy deconvolution adjusted)为例，假定滤波后信号为零均值的信号，则使用的目标函数为

峭度

(3)

基于最大化目标函数，求解获得逆滤波器系数，进而获得最终的滤波信号。逆滤波器系数的迭代公式可表示为

(4)

其中

(5)

滤波信号的迭代公式即为

(6)

2 基于相空间重构的盲解卷积方法

2.1 方法与原理

由原始信号构成的L×(N-L+1)相空间实矩阵X0可以分解为一系列奇异值对应的加噪信号子空间和噪声子空间。奇异值分解结果可表示为

(7)

其中

(8)

(9)

Σ=[diag(σ1,σ2,…,σl),0]∈RL×N-L+1

(10)

(11)

式中：U，V为单位正交矩阵;Σ为对角矩阵;并且σ1>σ2>…>σl;U1、Σ1、V1为加噪信号子空间对应的奇异值对角阵与奇异值矩阵;U2、V2和Σ2分别为噪声子空间对应的奇异值对角阵与奇异值矩阵。

式(7)给出了一种相空间分解方法，将等式代入式(6)可得

(12)

(13)

即

(14)

式中,VL为V矩阵的前L个子向量成的子矩阵。式(6)计算形式由整个相空间转换到了子空间计算，与MEDA等价。并且SVD得到的奇异值矩阵分量与信号频率成分相对应，转速变化不改变共振频率，因此该计算方法在变转速工况下仍然适用。一般来说，VL矩阵组成的相空间中不仅噪声子空间会对结果有负面影响，而且加噪信号子空间中的随机冲击信号、离散谐波信号以及分量噪声对解卷积结果的负面影响更为严重，因此需要对子空间降噪和筛选。

2.2 最优稀疏矩阵分量的筛选准则

稀疏编码收缩利用非高斯分量的统计特性获得软阈值收缩函数对矩阵分量降噪，并且可以滤除频率与目标信号相同或相似的噪声。假定VL矩阵分量可表示为

x=s+n

(15)

式中:s为纯信号;n为噪声。降噪结果可表示为

(16)

…

(17)

关联积分可表示为

(18)

其中,

(19)

(20)

d[X(i),X(j)]=

(21)

脉冲稀疏指数可以计算为

(22)

脉冲稀疏指数用来度量子空间时间序列中脉冲的稀疏程度，通常故障冲击信号成分多的子空间的脉冲稀疏指数大于其他子空间。在实际应用中，各子空间的嵌入维度取同一维度m，一般在3～10，ε取各子空间的标准差，取前n个脉冲稀疏指数大的子空间作为最优子空间，通常取滤波器长度的5%～10%。为避免式(14)中迭代过程发生幅值变化，因此对最优子空间的幅值进行修正

(23)

3 仿真结果分析

为了评价OSMCD方法的性能，本节中进行了与MED、MOMED、MCKD、CYCBD等方法的对比仿真研究。

如第1章所述，测试信号由故障冲击信号、随机冲击信号、离散谐波信号以及高斯白噪声组成

(24)

式中：第一部分为滚动轴承故障引起的冲击信号；Ai为第i次故障脉冲的幅值；Tb为两个故障冲击之间的时间间隔，与传统模式下平均周期T不同，该参数随着转速变化而改变；τi为由于滑动影响下第i次故障冲击相对于Tb的微小波动，通常是Tb的1%～2%。滚动轴承的脉冲响应函数表示为

si=e-atsin(2πfnt+φi)

(25)

式中：a为谐振阻尼系数；fn为共振频率；φi为相位。

式(24)第二部分表示在测量过程中产生的随机冲击干扰，由于随机冲击干扰的不确定性，因此随机冲击的时间Tr以及振幅Dj都是随机变量。式(24)第三部分表示来及与轴等其他部件的离散谐波信号，Pk、fk、βk分别表示第k次离散谐波信号的幅值、频率以及相位。式(24)第四部分表示高斯白噪声。

在此仿真中，采样频率设为10 000 Hz，采样时间为1 s，仿真工况转速变化曲线如图2所示。

图2 仿真工况转速变化曲线

(a) 故障冲击信号

各方法滤波器长度L均选择为50,迭代次数选择为200。对于MCKD以及CYCBD，需要输入给定故障特征周期，在仿真中输入的故障特征周期为平均转速对应的故障特征周期，即为0.028 1 s。OSMCD方法取脉冲稀疏指数最大的两个子空间为最优子空间。MED、MOMED、MCKD、CYCBD 4种方法对仿真信号的解卷积结果如图4所示，以及滤波结果的包络谱如图5所示。

从图4及图5中可以看出，MED、MCKD、CYCBD、MOMED 4种方法提取故障冲击信号效果都较差，在各方法的包络谱中不能明显地发现故障特征阶次。MED对随机冲击信号鲁棒性不足，滤波后随机冲击信号的幅值明显大于其他信号。MOMED虽然较MED对随机冲击信号的幅值明显较小，但解卷积结果仍然不能明显看出故障冲击信号。MCKD因为仿真信号是变转速工况下产生，故障冲击之间的时间间隔并不固定，因此依赖故障周期先验知识的MCKD并不能得到很好的提取效果。CYCBD在提供瞬时转速的情况下可以对变转速时的故障冲击有较好的提取效果，但仿真信号并没有提供瞬时转速，因此CYCBD同样不能得到很好的提取效果。

(a) MED

OSMCD方法得到的各子空间脉冲稀疏指数,如图6所示。

图6 各子空间脉冲稀疏指数

取脉冲稀疏指数最大的两个子空间作为最优子空间，则滤波结果如图7所示。

(a) 滤波信号

将图7中的解卷积结果和图4及图5中的结果相对比可以发现，OSMCD对故障冲击信号有更好的提取效果，对随机冲击的鲁棒性较好，噪声被有效消除，并且在解卷积过程中不需要预知故障特征周期。下面将以现场实测数据对OSMCD方法的有效性和先进性进一步验证。

4 现场测试验证

在机车试验台架对CRI-2692型机车轴箱轴承进行振动测试，整车滚动综合试验台可进行单节车辆的滚动试验，每次试验可搭载8个试验轴承。选用的试验轴承的结构参数,如表1所示。经计算可得轴承各特征阶次，如表2所示。OO,OI,Ob分别表示外圈故障特征阶次、内圈故障特征阶次和滚子故障特征阶次。

表1 轴承设计参数

表2 轴承特征频率及阶次

4.1 轴承外圈故障试验验证

首先分析实际故障为外圈故障的轴承试验。采样频率fs为5 000 Hz，共采集5 s，实际工况为从75 r/min匀加速到150 r/min，图8为外圈故障信号及其包络谱。

(a) 原始信号

各方法滤波器长度L均选择为50,迭代次数选择为200。MCKD以及CYCBD输入的故障特征周期为0.065 s。MED、MCKD、CYCBD、MOMED 4种方法的滤波结果如图9所示以及滤波结果的包络谱如图10所示。

从图9、图10可以看出，由于有很强的随机冲击的干扰，4种方法都不能有效的提取故障冲击信号，在包络谱中也难以发现故障特征阶次。OSMCD方法滤波器长度L选择为50,迭代次数选择为200。OSMCD方法得到各子空间脉冲稀疏指数,如图11所示。

(a) MED

图11 各子空间脉冲稀疏指数

选择脉冲稀疏指数最大的两个子空间进行信号的解卷积，得到滤波结果以及滤波结果的包络谱如图12所示。

从图12中可以看到，OSMCD方法相比较于其他4种方法，明显地提取了故障冲击信号,并有效消除了噪声，对随机冲击具有良好的鲁棒性，在包络谱中能明显发现外圈故障特征阶次。

(a) 滤波信号

4.2 轴承滚动体故障试验验证

首先分析实际故障为滚动体故障的轴承试验。采样频率fs为5 000 Hz，共采集5 s，实际工况为从60 r/min匀加速到105 r/min，图13为滚动体故障信号及其包络谱。

(a) 原始信号

各方法滤波器长度L均选择为50,迭代次数选择为200。MCKD以及CYCBD输入的故障特征周期为0.210 s。MED、MCKD、CYCBD、MOMED 4种方法的滤波结果如下图14所示以及滤波结果的包络谱如图15所示。

从图14以及图15中可以看出，MED、MCKD、MOMED 3种方法都不能有效的提取故障冲击信号，在包络谱中也难以发现故障特征阶次。CYCBD方法提取出了一部分故障冲击，在包络谱中可以发现一阶故障特征阶次，但仍有较大的噪声。OSMCD方法选择滤波器长度L选择为50,迭代次数选择为200。OSMCD方法得到各子空间脉冲稀疏指数,如图16所示。