CN105699082B - 一种稀疏化的最大谐噪比解卷积方法 - Google Patents

Abstract

一种稀疏化的最大谐噪比解卷积方法，首先对采集的信号进行截断和去均值处理，然后对未提前给定精确周期的情况进行估计周期操作，再对信号进行解卷积处理将谐噪比作为目标函数，对滤波器系数进行求导，进而得到迭代表达式，并在每次迭代滤波的过程中对滤波信号进行稀疏处理，稀疏处理的阈值和周期都会随滤波后的信号进行更新，最后对解卷积后的信号进行包络分析，从包络谱中能提取故障特征频率，本发明在特征频率的提取过程中不需要人为参与，有利于实现故障特征提取和诊断监测的自动化，节约时间，效率更高。

Description

一种稀疏化的最大谐噪比解卷积方法

技术领域

本发明涉及机械设备故障诊断技术领域，特别涉及一种稀疏化的最大谐噪比解卷积(Sparse Maximum Harmonics-to-noise ratio Deconvolution，SMHD)方法。

背景技术

振动分析为现阶段机械设备故障诊断最为有效的途径之一，机械设备的状态劣化往往表现为振动信息的变化或异常。目前基于振动信息的信号处理方法，例如时域法、频域法还有时频域法，这些都已成功应用于轴承故障诊断中，并产生了很好的效果。然而，在滚动轴承故障诊断领域仍然面临着很多挑战，轴承故障的提取仍然具有很多困难。1、测试传感器与故障源之间冗长而复杂的传递路径能严重影响传递函数，进而使冲击信号的幅值降低、时间拉长，从而导致故障引起的脉冲很容易被噪声覆盖。2、轴承中滚子的随机波动会导致原本呈准周期的故障冲击的频谱包络谱进一步模糊。3、来自机械系统中的非周期性噪声及周期性干扰的影响，为提取轴承故障冲击增加了更多挑战。

解卷积方法以其能消除传递路径影响并能增强故障冲击的优点被广泛应用。2007年最小熵解卷积(Minimum entropy deconvolution，MED)方法被Sawalhi和Randall等学者首次应用于滚动轴承的故障诊断领域，并取得了一定的效果。MED是一种不需要任何先验假设的信号时域盲解卷积技术，通过迭代选择一个(finite impulse response)FIR滤波器来最小化滤波信号的输出熵(即峭度最大化)，旨在最小化噪声的同时提取出故障冲击，因此其在高噪声下也能得到理想的诊断结果。虽然MED对冲击性的增强与提取具有明显的效果，但是其目标函数仅仅是寻求滤波信号峭度值最大化，所以易受随机孤立冲击的干扰。即当故障信号中同时存在周期性冲击序列和随机孤立冲击时，通过MED增强技术有可能只能增强孤立冲击，而对真正反映故障特征的周期性冲击序列没有效果。鉴于上述问题，2012年McDonald等学者提出了最大相关峭度解卷积(Maximum correlated kurtosisdeconvolution,MCKD)算法，提出了相关峭度的概念，兼顾了故障的冲击性和周期性特征，减小了随机冲击的干扰，并成功地将其应用到齿根裂纹故障诊断上。然而MCKD方法对提前准确估计故障周期具有极大的依赖性，而在工程实际中，由于测速设备的限制，导致测速装置难以安装或者成本很高，而且复杂的工况使得设备转速不可能保持恒定。这些都会影响测速的精度进而导致周期估计存在误差，同时MCKD对由随机滑动引起的周期波动容忍性差，这些不足给MCKD的应用推广带来了诸多不便。

发明内容

为了克服上述现有技术的缺点，本发明的目的在于提供一种稀疏化的最大谐噪比解卷积方法，在不提供精准预估故障周期的情况下，也能实现精确的故障诊断。

为实现上述目的，本发明采取的技术方案为：

一种稀疏化的最大谐噪比解卷积方法，包括以下步骤：

步骤一：将振动加速度传感器吸附于被测试滚动轴承的轴承座上，并对测试得到的信号进行高频采样、截断和去均值处理，将信号记为x；

步骤二：对提供的周期信息进行验证，判断是否提供了精准的故障周期，如果提供精准的故障周期，则利用这个周期；如果没有提供精准周期，则需要估计周期，利用谐噪比(Harmonics-to-Noise Ratio，HNR)的定义，确定信号x的包络信号的自相关函数中除0位置以外的局部极大值为周期；

步骤三：将谐噪比为目标函数对滤波系数进行求偏导：

其中t为时间，T为周期，对目标函数进行离散化后，使用其中f(l)为滤波器系数，l＝1,2,……,L,L为滤波器长度，求导后的结果为：

将上式写成矩阵的形式：

Af＝b (3)

其中：

b--逆滤波器的输入信号x、输出信号y的互相关，b为L维列向量；

A--输入信号x的自相关，A为L×L维矩阵；

f--逆滤波器的滤波器系数，f为L维列向量；

首先计算自相关矩阵A；再假设逆滤波器的初始值f⁽⁰⁾，设置滤波器长度L＝100，并给定初始滤波器系数为[0 0…… 1-1 ……0 0]，用y⁽⁰⁾和x⁽⁰⁾计算列向量b⁽¹⁾；然后求解新的滤波器系数f⁽¹⁾＝A^-1b⁽¹⁾；另外在每次更新滤波器系数后，对滤波信号进行稀疏处理即：

其中σ为阈值常数，初始阈值为原信号的均方根值或者绝对均值，为稀疏变换后的滤波信号；

然后求滤波后的信号的峭度，比较前后两次滤波信号的峭度，以此来设置新的阈值σ，当峭度增加时，增大σ；当峭度减小时，减小σ，并且对滤波后的信号求包络，计算包络信号的周期，以此更新周期，设定最大迭代次数为30次，得到经过SMHD处理后的信号y_k；

步骤四：对SMHD处理后的信号y_k进行包络分析并得到包络谱，对包络谱进行分析，进而提取故障特征频率。

本发明相比于现有技术，具有以下有益效果：

a)本发明源于传统MED方法，用谐噪比代替峭度作为目标函数，能兼顾故障的冲击性和周期性特征，极大地减小了随机冲击的干扰。

b)本发明不需要任何先验知识，也无需对系统故障特征频率进行预估，方法具有较好的鲁棒性。

c)本发明能提取故障周期具有一定波动的冲击信号，尤其是带有随机滑动的滚动轴承故障信号。

d)本发明对能提供预估周期的信号具有更好的效果，并且对给定的初始周期具有很大的容忍性。

附图说明

图1为实施例机车轮对轴承试验台。

图2为本发明方法流程图。

图3为实施例中对原信号进行截断处理后的信号x。

图4为实施例中信号x的包络信号。

图5为实施例中原信号的频谱图。

图6为实施例中原信号的包络谱图。

图7为实施例中经SMHD方法处理后信号y_k。

图8为实施例中经SMHD方法处理后信号y_k的包络谱图。

图9为实施例中经MCKD滤波后的信号。

图10为实施例中经MCKD滤波后的信号的包络谱图。

具体实施方式

下面结合附图和实施例对本发明做详细描述。

以某机车轴承试验台为例，该试验台由液压马达、驱动轮、轴承和机车轮对等组成，如图1所示，液压马达带动驱动轮运动进而驱动轴承外圈运动，轴承内圈固定在机车轮对的车轴上，加速度传感器固定在轴承一端，测量轴承的振动信号。

由于试验台速度不是恒定的，所以无法精确估计轴承的故障特征频率，显然MCKD方法不太适用该场合。因此可以先通过本发明提出的方法提取出精确的故障特征频率，然后带入MCKD方法，进行效果对比。

如图2所示，稀疏化的最大谐噪比解卷积方法，包括以下步骤：

步骤一：将振动加速度传感器吸附于被测试滚动轴承的轴承座上，并对测试得到的信号进行高频采样、截断和去均值处理，其中采样频率为76.8kHz，在使用数据时需要去掉起始噪声部分，截取整段信号中共2s的数据，如图3所示，将信号记为x；

步骤二：对提供的周期信息进行验证，判断是否提供了精准的故障周期，显然没有提供精准周期，则需要估计周期，利用谐噪比的定义，计算信号x的包络信号，如图4所示，并确定它的自相关函数中除0位置以外的局部极大值为周期，经过计算周期为501个采样点；

步骤三：将谐噪比为目标函数对滤波系数进行求偏导：

其中t为时间，T为周期，对目标函数进行离散化后，使用其中f(l)为滤波器系数，l＝1,2,……,L,L为滤波器长度，求导后的结果为：

将上式写成矩阵的形式：

Af＝b (3)

其中：

b--逆滤波器的输入信号x、输出信号y的互相关，b为L维列向量；

A--输入信号x的自相关，A为L×L维矩阵；

f--逆滤波器的滤波器系数，f为L维列向量；

首先计算自相关矩阵A；再假设逆滤波器的初始值f⁽⁰⁾，设置滤波器长度L＝100，并给定初始滤波器系数为[0 0…… 1-1 ……0 0]，用y⁽⁰⁾和x⁽⁰⁾计算列向量b⁽¹⁾；然后求解新的滤波器系数f⁽¹⁾＝A^-1b⁽¹⁾；另外在每次更新滤波器系数后，对滤波信号进行稀疏处理即：

其中σ为阈值常数，初始阈值为原信号的绝对均值0.0871，为稀疏变换后的滤波信号；

然后求滤波后信号的峭度，比较前后两次滤波信号的峭度，以此来设置新的阈值σ，当峭度增加时，增大σ；当峭度减小时，减小σ。并且对滤波后的信号求包络，计算包络信号的周期，以此更新周期，设定最大迭代次数为30次，得到经过SMHD处理后的信号y_k；

步骤四：对SMHD处理后的信号y_k进行包络分析并得到包络谱，对包络谱进行分析，进而提取故障特征频率。

参照图5和图6，图5和图6分别是原信号的频谱图和包络谱图，这两幅图无法对应上轴承故障特征频率，参照图7，图7为本发明提出的方法SMHD滤波后的信号，由于该试验台与一般的内圈旋转的试验台不同，该试验台是外圈旋转，所以从图7可以看出，轴承故障为外圈故障，该滤波信号效果十分突出，信噪比很高，直接从时域波形就能判断故障类型。参照图8，图8为对应的包络谱，包络谱也十分干净清楚，从图8中可以得到，该轴承故障类型为外圈故障，故障特征频率为26.07Hz，将其带入MCKD方法，按照McDonald等提出的MCKD方法对重采样后的信号进行处理，对该方法的参数按照McDonald提出的推荐设置，其中参数设置为精确计算故障特征频率f＝26.07Hz，移位距离T＝f_s*1/f＝76800/26.07，移位次数M＝3，滤波器长度设置为100，最大滤波次数设置为30。参照图9和图10，图9为MCKD滤波后的信号，图10为对应的包络谱，对比两种方法的处理结果，本发明提出的方法优势十分明显，时域波形比MCKD滤波后的信号具有更高的信噪比，而且包络解调后的频谱也突显了这一优势。而且本发明不需要提取轴承的精确预估故障周期，在工程实际中更加适用。

Claims (1)

1.一种稀疏化的最大谐噪比解卷积(Sparse Maximum Harmonics-to-noise ratioDeconvolution，SMHD)方法，其特征在于，包括以下步骤：

步骤一：将振动加速度传感器吸附于被测试滚动轴承的轴承座上，并对测试得到的信号进行高频采样、截断和去均值处理，将信号记为x；

步骤二：对提供的周期信息进行验证，判断是否提供了精准的故障周期，如果提供精准的故障周期，则利用这个周期；如果没有提供精准周期，则需要估计周期，利用谐噪比(Harmonics-to-Noise Ratio，HNR)的定义，确定信号x的包络信号的自相关函数中除0位置以外的局部极大值为周期；

步骤三：将谐噪比为目标函数对滤波系数进行求偏导：

<mrow> <msub> <mi>O</mi> <mi>k</mi> </msub> <mo>&amp;lsqb;</mo> <mi>f</mi> <mrow> <mo>(</mo> <mi>l</mi> <mo>)</mo> </mrow> <mo>&amp;rsqb;</mo> <mo>=</mo> <mfrac> <mrow> <mo>&amp;Integral;</mo> <mi>y</mi> <mrow> <mo>(</mo> <mi>t</mi> <mo>)</mo> </mrow> <mi>y</mi> <mrow> <mo>(</mo> <mi>t</mi> <mo>+</mo> <mi>T</mi> <mo>)</mo> </mrow> <mi>d</mi> <mi>t</mi> </mrow> <mrow> <mo>&amp;Integral;</mo> <msup> <mi>y</mi> <mn>2</mn> </msup> <mrow> <mo>(</mo> <mi>t</mi> <mo>)</mo> </mrow> <mi>d</mi> <mi>t</mi> <mo>-</mo> <mo>&amp;Integral;</mo> <mi>y</mi> <mrow> <mo>(</mo> <mi>t</mi> <mo>)</mo> </mrow> <mi>y</mi> <mrow> <mo>(</mo> <mi>t</mi> <mo>+</mo> <mi>T</mi> <mo>)</mo> </mrow> <mi>d</mi> <mi>t</mi> </mrow> </mfrac> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>1</mn> <mo>)</mo> </mrow> </mrow>

其中t为时间，T为周期，对目标函数进行离散化后，使用其中f(l)为滤波器系数，l＝1,2,……,L,L为滤波器长度，求导后的结果为：

将上式写成矩阵的形式：

Af＝b (3)

其中：

b--逆滤波器的输入信号x、输出信号y的互相关，b为L维列向量；

A--输入信号x的自相关，A为L×L维矩阵；

f--逆滤波器的滤波器系数，f为L维列向量；

首先计算自相关矩阵A；再假设逆滤波器的初始值f⁽⁰⁾，设置滤波器长度L＝100，并给定初始滤波器系数为[00……1-1……00]，用y⁽⁰⁾和x⁽⁰⁾计算列向量b⁽¹⁾；然后求解新的滤波器系数f⁽¹⁾＝A^-1b⁽¹⁾；另外在每次更新滤波器系数后，对滤波信号进行稀疏处理即：

<mrow> <mover> <mi>y</mi> <mo>&amp;OverBar;</mo> </mover> <mo>=</mo> <mi>y</mi> <mo>&amp;times;</mo> <mrow> <mo>(</mo> <mn>1</mn> <mo>-</mo> <mi>exp</mi> <mo>(</mo> <mfrac> <mrow> <mo>-</mo> <msup> <mi>y</mi> <mn>2</mn> </msup> </mrow> <mrow> <mn>2</mn> <msup> <mi>&amp;sigma;</mi> <mn>2</mn> </msup> </mrow> </mfrac> <mo>)</mo> <mo>)</mo> </mrow> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>4</mn> <mo>)</mo> </mrow> </mrow>

其中σ为阈值常数，初始阈值为原信号的均方根值或者绝对均值，为稀疏变换后的滤波信号；

然后求滤波后的信号的峭度，比较前后两次滤波信号的峭度，以此来设置新的阈值σ，当峭度增加时，增大σ；当峭度减小时，减小σ，并且对滤波后的信号求包络，计算包络信号的周期，以此更新周期，设定最大迭代次数为30次，得到经过SMHD处理后的信号y_k；

步骤四：对SMHD处理后的信号y_k进行包络分析并得到包络谱，对包络谱进行分析，进而提取故障特征频率。