CN111507221A - 基于vmd和最大重叠离散小波包变换的齿轮信号去噪方法 - Google Patents

Abstract

本发明公开了一种基于VMD和最大重叠离散小波包变换的齿轮信号去噪方法，将变分模态分解VMD和最大重叠离散小波包变换MODWPT相结合来对齿轮故障信号进行去噪，通过实测数据的实验结果表明，本发明的方法可有效滤除信号中的噪声干扰和迭代误差成分，实现信号特征的精确提取。

Description

基于VMD和最大重叠离散小波包变换的齿轮信号去噪方法

技术领域

本发明涉及齿轮故障信号去噪技术领域，特别是一种基于VMD和最大重叠离散小波包变换的齿轮信号去噪方法。

背景技术

齿轮作为最常用的传动装置被广泛应用于各类机械设备中。当齿轮出现故障时，会导致设备的整体性能下降，严重情况下会造成设备事故并带来巨大经济损失。在实际工况中，齿轮工作环境十分恶劣，环境噪声和信号衰减等因素对于早期故障特征有较大影响，而齿轮早期故障特征比较微弱，往往淹没在强背景噪声中难以提取。因此，对于齿轮故障信号去噪方法的研究具有重要意义。

当齿轮出现故障时，故障齿在啮合过程中，不仅会引起振动信号幅值和相位的变化，同时也会使振动信号的能量发生变化。相较于正常状态，齿轮的故障信号大多数是多分量的调幅-调频信号，且表现出非线性和非平稳性的特点。众多学者在该类信号降噪方面进行了大量研究。Huang等提出的经验模态分解(Empirical mode decomposition，EMD)方法是一种自适应信号处理方法，它能够根据信号的自身特性，将其分解成一些列表征信号特征的由高频到低频的固有模态函数(Intrinsic Mode Function，IMF)分量。

然而，由于自身算法的限制，EMD方法在分解过程中易产生模态混叠问题，影响信号的去噪效果。WU等针对该问题，提出一种噪声辅助分解的集合经验模态分解(EnsembleEmpirical Mode Decomposition，EEMD)。

实际计算中由于受限所添加白噪声的次数，EEMD分解过程中向信号内添加的噪声难以彻底消除，这将影响信号特征提取的准确性；同时，该方法的模态选取属于递归“筛选”，受采样频率影响较大，实际信号处理效果并不理想。近年来，Dragomiretskiy等提出一种可变尺度的非平稳信号分析方法——变分模态分解(Variational ModeDecomposition，VMD)，此方法可将复杂信号分解为多个单分量调幅-调频信号和的形式，分解过程中的模态个数可自适应地确定，有效避免了模态混叠问题的产生，并表现出更好的噪声鲁棒性。因此，针对齿轮故障信号的特性，VMD方法相较于传统递归式分解方法，具有更强的适用性。

基于上述分析，为有效地从含噪齿轮信号中提取相关故障特征信息，本专利提出一种基于VMD和最大重叠离散小波包变换(Maximal Overlap Discrete Wavelet PacketTransform，MODWPT)相结合的齿轮信号去噪方法。通过实测数据的实验结果表明，本文所提方法可有效滤除信号中的噪声干扰和迭代误差成分，实现信号特征的精确提取，为该类问题的解决提供了一种可行的方法。

发明内容

本发明的目的是要解决现有技术齿轮实际工况复杂、常伴随噪声干扰致使故障特征难以提取的问题，提供一种可有效滤除信号中的噪声干扰和迭代误差成分、实现信号特征的精确提取的基于VMD和最大重叠离散小波包变换的齿轮信号去噪方法。

为达到上述目的，本发明是按照以下技术方案实施的：

一种基于VMD和最大重叠离散小波包变换的齿轮信号去噪方法，包括以下步骤：

S1、采用变分模态分解VMD方法对采集到的齿轮故障信号进行分解，初始化预设尺度K值设为2，分解得到每个IMF分量的中心频率ω_k；

S2、判别预设尺度K值下分解所得各IMF分量中心频率ω_k是否出现相近值，若出现，则预设尺度K值设为K-1，否则，预设尺度K值设为K+1继续步骤S2；

S3、在选定最优预设尺度K值后，计算惩罚因子α在不同取值范围下，重构信号的模糊熵值FE，取min{FE}所对应的α作为最优惩罚因子用于齿轮故障信号的VMD分解；

S4、计算各IMF分量的能量熵增量Δq_i，根据Δq_i数值初定虚假模态分量；

S5、各IMF分量的频域互相关系数ρ_f，通过ρ_f数值综合判定虚假模态和噪声干扰成分，得到初步去噪信号x'(t)；

S6、针对含有高频噪声干扰的IMF分量采用最大重叠离散小波包变换MODWPT方法对其进行去噪处理，进一步提升去噪效果和性能指标；

S7、将MODWPT去噪的IMF分量和表征信号特征的各主模态分量进行重构，形成去噪信号x”(t)。

与现有技术相比，本发明将变分模态分解VMD和最大重叠离散小波包变换MODWPT相结合来对齿轮故障信号进行去噪，通过实测数据的实验结果表明，本发明的方法可有效滤除信号中的噪声干扰和迭代误差成分，实现信号特征的精确提取。

附图说明

图1为本发明的齿轮信号去噪流程图。

图2为本发明应用实施例中采集的齿轮故障信号时域波形及其包络谱：(a)为齿轮故障信号时域波形；(b)为齿轮故障信号包络谱。

图3为本发明应用实施例中齿轮故障信号经VMD分解后各IMF分量及其频谱。

图4为本发明应用实施例中采用EEMD方法对齿轮故障信号分解后的前5阶IMF分量的时域波形和频谱。

图5为本发明应用实施例中VMD分解得到各IMF分量的能量熵增量曲线。

图6为本发明应用实施例中VMD分解得到的各IMF分量的频域互相关系数。

图7为本发明应用实施例中高频噪声分量IMF1的MODWPT分解结果。

图8为本发明应用实施例中MODWPT去噪的IMF分量和表征信号特征的各主模态分量进行重构后去噪信号的时域波形及包络谱：(a)为重构后去噪信号时域波形；(b)为重构后去噪信号包络谱。

图9为本发明对比实施例中采用EEMD的小波阈值信号去噪方法对采集的齿轮故障信号去噪后齿轮故障信号时域波形及其包络谱：(a)为经EEMD小波阈值去噪后齿轮故障信号时域波形；(b)为经EEMD小波阈值去噪后齿轮故障信号包络谱。

具体实施方式

为使本发明的目的、技术方案及优点更加清楚明白，以下结合实施例，对本发明进行进一步的详细说明。此处所描述的具体实施例仅用于解释本发明，并不用于限定发明。

如图1所示，本实施例的一种基于VMD和最大重叠离散小波包变换的齿轮信号去噪方法，包括以下步骤：

S1、按照下述VMD方法对采集到的齿轮故障信号x(t)进行分解，由VMD算法可知，若模态分解个数K取值合理，则各IMF分量间的中心频率分布较为合理，数值产生出现相近或重叠的现象，即不会出现过分解或欠分解的模态混叠问题。因此，以各IMF分量中心频率是否出现相近项法确定模态分解个数K的最佳值。初始化预设尺度K值设为2，计算分解得到每个IMF分量的中心频率ω_k；

对采集到的齿轮振动信号x(t)其建立变分约束模型：

式中，K为模态分解个数；u_k(t)为第k个模态分量；ω_k为u_k(t)的中心频率。

引入增广拉格朗日函数ζ，将约束问题转化为非约束问题：

式中，α为二次惩罚因子；<>表示向量内积。

VMD采用乘法算子交替的方法求取上式，得到最优解将信号分解成为K个IMF，其实现流程如下:

1)初始化

λ¹、n，令其初始值为0；

2)执行循环n＝n+1；

3)令k＝0，k＝k+1，当k<K时，执行内层第1个循环，将u_k(t)更新为

4)令k＝0，k＝k+1，当k<K时，执行内层第2个循环，将ω_k(t)更新为

5)更新λ；

6)重复步骤2)～5)，直至满足迭代约束条件

则停止迭代，输出分解所得各IMF分量中心频率ω_k。

S2、根据公式(6)判别预设尺度K值下分解所得各IMF分量中心频率ω_k是否出现相近值，若出现，则预设尺度K值设为K-1，否则，预设尺度K值设为K+1继续步骤S2；

ω_k＝f_s×nor(ω_k) (6)

式中，f_s为信号的采样频率，nor(ω_k)为VMD分解所得IMF分量的归一化中心频率。

S3、当采用VMD方法对轴承振动信号进行分解重构后，当信号内的环境噪声和背景干扰成分较弱时，此时信号内包含可有效表征故障特征的主要信息，且应具有较明显的冲击特性，其模糊熵值较小。因此，在选定最优预设尺度K值后，计算惩罚因子α在不同取值范围下，根据公式(7)～(10)重构信号的模糊熵值FE，取min{FE}所对应的α作为最优惩罚因子用于齿轮故障信号的VMD分解；

1)设X＝[x(1),x(2),…,x(N)]为选定最优预设尺度K值下的重构信号，对X设定模式维数m，由X构造m维向量为

式中，

为(n(i),n(i+1),…,n(i+m-1))的均值，i∈[1,N-m]。

2)定义d_ij为任意两A^m间的最大欧式距离，并用指数函数u定义相似度D_ij

D_ij＝u(d_ij,n,r)＝exp[-(d_ij/r)ⁿ] (8)

式中，u为关于

的模糊隶属度函数；n和r分别为其边界的梯度和宽度。

3)定义

并重复步骤1)～3)，算出m+1时候的Φ^m+1(n,r)

4)在以上步骤的基础上，可定义模糊熵

式中，m为模式维数；n为相似容限；N为数据长度。

S4、信号经VMD分解所得IMF分量的中心频率由低频至高频分布，具有不同的能量。表征信号自身特征信息的主模态分量应占有主要能量，而特征信息不敏感的虚假模态分量的能量所占比例较小，故计算各IMF分量的能量熵增量Δq_i，根据Δq_i数值可初定虚假模态分量；

设u₁(t),u₂(t),…,u_n(t)为信号x(t)经VMD分解得到的IMF分量。基于IMF能量熵增量的敏感模态函数判别算法的具体过程如下：

计算各IMF分量的能量：

式中，u_i(t)为第i个IMF分量。

求解各IMF分量的归一化能量：

式中，N为IMF分量的总个数。

计算各IMF分量的能量熵增量值：

Δq_i＝-p(i)log₂(p(i)) (13)

S5、频域内噪声成分的功率谱密度较小，且功率谱的互相关系数受噪声干扰小，为此，再根据公式(14)计算各IMF分量的频域互相关系数ρ_f，通过ρ_f数值综合判定虚假模态和噪声干扰成分，得到初步去噪信号x'(t)；

根据时域互相关系数计算方法，设G_x、G_y分别为信号x_i和y_i功率谱，f_a为分析频率，则频域内信号x_i和y_i的互相关系数可表示为：

式中，∣ρ_f∣越大，表明两信号在频域上相关性越好；反之，∣ρ_f∣越小，说明两信号在频域上相关性越差。

S6、针对含有高频噪声干扰的IMF分量采用最大重叠离散小波包变换MODWPT方法对其进行去噪处理，进一步提升去噪效果和性能指标；

MODWPT的分解系数可由W_j,n＝{W_j,n,t,t＝0,…,N-1}表示，其中j为分解层数，n可看作一个随j变化的频率索引，则可计算出MODWPT的分解系数：

式中，若n mod 4＝0或3，则

若n mod4＝1或2，则

与传统DWT方法相比，MODWT具有小波系数和尺度系数的平移不变性、所有分解层数都保持相同时间分辨率、无相位扭曲等突出优点，非常适合处理非线性、非平稳信号。

S7、将MODWPT去噪的IMF分量和表征信号特征的各主模态分量进行重构，形成去噪信号x”(t)。

式中，x”(t)为最终去噪信号；x'(t)为经能量熵增量和频域互相关系数判定的对故障信息敏感的IMF分量所组成的重构信号；

为经MODWPT分解去噪后的IMF分量所构成的重构信号。

应用实施例

为验证上述实施例的方法的有效性，对实际采集到的齿轮断齿故障信号进行去噪处理与分析。所采用试验平台为PQZZ-II型机械故障模拟综合试验台。断齿故障模拟时，为驱动齿轮的一齿折断，折断线位于分度圆附近。

试验过程中，电动机转频为50Hz，数据采集系统以ADA16-8/2(LPCI)型高速多功能采集卡为基础进行搭建，传感器采用KD1001L型压电式加速度传感器，为全面监测齿轮箱运行状态，加速度传感器安装在安装在输出轴轴承X方向。信号采集时，设置采样频率f_s＝5120Hz，采样时长15s，采样精度为16bt，分析时长1s。齿轮箱中各齿轮均为标准直齿轮，其中，驱动齿轮和从动齿轮的齿数分别为55和75，齿轮材料为S45C，模数为2mm，输入轴的平均转速为880r/min。经计算，可得驱动齿轮的转频f₁＝14.67Hz，从动齿轮转频f₂＝10.76Hz，啮合频率f_z＝806.67Hz。

试验采集到的齿轮断齿故障振动信号的时域波形如图2(a)所示。对于齿轮故障的诊断，包络谱分析是最为常用的方法，因此，为诊断出齿轮故障类型，求解齿轮故障信号的包络谱，并示于图2(b)。

由图2可知，由于信号采集系统中未采用相应消噪装置，致使采集到的时域信号中含有较多噪声成分；由于噪声干扰成分的影响，该信号的包络谱中只在与驱动齿轮2倍转频相对应的29.4Hz处出现较为明显的峰值，根据齿轮振动机理，仅可判定与输入轴相连的驱动齿轮存在故障，但齿轮的故障特征以及所属故障类型从中难以读取和判别。

按照本文所提VMD分解参数选择方法对预设尺度K和惩罚参数α进行选取，经计算K＝4、α＝2100时，分解效果最佳。齿轮故障信号经VMD分解后各IMF分量及其频谱如图3所示。

由图3可知，VMD方法在一定程度上克服了EMD中存在的模态混叠问题，各IMF分量集中在各自中心频率附近，有效减少了各模态分量间的信息泄露。

为比较分解效果，采用EEMD方法对齿轮故障信号进行分解，得到12阶IMF分量，由于EEMD分解后的信号特征主要集中在前几个模态分量，在此仅列出前5阶模态分量的时域波形和频谱，并示于图4。其中，EEMD分解时，所添加白噪声的标准差ε₀＝0.2，总体实验次数I＝100。

由图4可知，EEMD方法在一定程度上抑制了模态混叠问题，但仍可看出分解效果并不理想，尤其对于IMF1分量，在其整个频率范围内机会都存在频率成分，且该方法都分解出了多阶迭代误差分量。这势必将影响后续信号故障特征提取的准确性与故障诊断的可靠性。

计算VMD分解得到各IMF分量的能量熵增量Δqi，结果如图5所示。

由图5可知，IMF1和IMF4的能量熵增量数值明显小于IMF2和IMF3分量，按本文算法将其暂定为虚假干扰的模态分量。进一步计算各IMF分量同原信号间的频域互相关系数，结果如图6所示。

由图6可知，IMF2和IMF3的频域互相关系数明显大于IMF1和IMF4，表明在降低噪声干扰的情况下，IMF2和IMF3包含原信号的信息较多，同图5的能量熵增量分析结果相一致。

根据上述分析并结合所分析信号的特性，判定IMF2和IMF3为表征信号自身特征的主模态分量；IMF1为包含分解迭代误差和背景信号的干扰分量；IMF4为高频噪声干扰成分。

为有效提升去噪效果和性能指标，采用MODWPT对包含高频噪声干扰成分的IMF1进行分解，分解过程中选择长度L＝22的Fejer-Korovkin小波滤波器，分解层数J＝4，低频分解结果如图7所示。

由图7可知，高频噪声分量IMF1经MODWPT分解后，其频率最低分量C1中的噪声干扰成分得到有效滤除，信号自身特征较为明显。将其作为高频去噪分量与VMD分解所得表征信号自身特征的主模态IMF2和IMF3分量进行重构，重构后去噪信号的时域波形和包络谱如图8所示。

由图8(a)可知，本文提出的基于VMD-MODWPT的去噪方法有效滤除了原信号内绝大多数无用的低频干扰和高频噪声成分，齿轮故障信号得到有效提纯；同时分析图8(b)的包络谱图可以发现在14.71Hz、20.63Hz、29.38Hz、44.15Hz、53.84Hz、58.75Hz、66.31Hz和73.13Hz出现明显峰值，这些峰值频率分别对应着驱动齿轮的1～5倍转频f1以及从动齿轮的2倍频、5倍频和6倍频，而其它频率区域峰值较低。

根据齿轮故障诊断原理以及齿轮断齿故障包络谱特征可知：当齿轮发生断齿故障时，在啮合过程中会激发瞬时冲击缺陷，致使断齿故障信号具有丰富的频域特征。断齿故障的包络谱中除包含故障齿轮的转频fr及其整数倍频外，同时还会有与其相啮合齿轮的高倍转频的出现。根据上述原理并结合图8(b)可诊断出与输入轴相连的驱动齿轮存在断齿故障。

对比实施例

为说明上述实施例的方法的有效性和优越性，采用EEMD的小波阈值信号去噪方法对齿轮故障信号进行去噪处理，去噪结果如图9所示。其中，小波阈值处理信号时采用db5小波基和“Heursure”规则进行软阈值降噪。

由图9可知，基于EEMD小波阈值信号去噪方法的去噪效果并不理想，去噪后信号中仍有大量噪声成分残留。而其包络谱中，虽能识别出驱动齿轮的转频及其部分整数倍频，但与其相啮合齿轮的高倍转频成分的谱线并不明显，从中无法准确判别齿轮的故障类型；同时，由于低频背景信号和迭代误差等虚假成分并未得到有效滤除，去噪信号包络谱的低频区域存在较多虚假频率成分，这也进一步增加了信号故障类型诊断的难度。

综述，通过实测数据的实验结果表明，本发明的方法可有效滤除信号中的噪声干扰和迭代误差成分，实现信号特征的精确提取。

本发明的技术方案不限于上述具体实施例的限制，凡是根据本发明的技术方案做出的技术变形，均落入本发明的保护范围之内。

Claims (8)

1.一种基于VMD和最大重叠离散小波包变换的齿轮信号去噪方法，其特征在于，包括以下步骤：

S1、采用变分模态分解VMD方法对采集到的齿轮故障信号进行分解，初始化预设尺度K值设为2，分解得到每个IMF分量的中心频率ω_k；

S2、判别预设尺度K值下分解所得各IMF分量中心频率ω_k是否出现相近值，若出现，则预设尺度K值设为K-1，否则，预设尺度K值设为K+1继续步骤S2；

S3、在选定最优预设尺度K值后，计算惩罚因子α在不同取值范围下，重构信号的模糊熵值FE，取min{FE}所对应的α作为最优惩罚因子用于齿轮故障信号的VMD分解；

S4、计算各IMF分量的能量熵增量Δq_i，根据Δq_i数值初定虚假模态分量；

S5、各IMF分量的频域互相关系数ρ_f，通过ρ_f数值综合判定虚假模态和噪声干扰成分，得到初步去噪信号x'(t)；

S6、针对含有高频噪声干扰的IMF分量采用最大重叠离散小波包变换MODWPT方法对其进行去噪处理，进一步提升去噪效果和性能指标；

S7、将MODWPT去噪的IMF分量和表征信号特征的各主模态分量进行重构，形成去噪信号x”(t)。

2.根据权利要求1所述的基于VMD和最大重叠离散小波包变换的齿轮信号去噪方法，其特征在于，所述步骤S1具体包括：

S101、对采集到的齿轮振动信号x(t)其建立变分约束模型：

式中，K为模态分解个数，K值设为2；u_k(t)为第k个模态分量；ω_k为u_k(t)的中心频率；

S102、引入增广拉格朗日函数ζ，将步骤S101的变分约束模型转化为非约束模型：

式中，α为二次惩罚因子；<>表示向量内积；

S103、采用乘法算子交替的方法求取步骤S102的非约束模型，得到最优解将信号分解成为K个IMF分量，乘法算子交替的方法的求解流程如下：

1)初始化

λ¹、n，令其初始值为0；

2)执行循环n＝n+1；

3)令k＝0，k＝k+1，当k<K时，执行内层第1个循环，将u_k(t)更新为

4)令k＝0，k＝k+1，当k<K时，执行内层第2个循环，将ω_k(t)更新为

5)更新λ：

6)重复步骤2)～5)，直至满足迭代约束条件

则停止迭代，输出分解所得各IMF分量中心频率ω_k。

3.根据权利要求2所述的基于VMD和最大重叠离散小波包变换的齿轮信号去噪方法，其特征在于：所述步骤S2中判别预设尺度K值下分解所得各IMF分量中心频率ω_k是否出现相近值是根据以下公式来判别的：

ω_k＝f_s×nor(ω_k) (6)

式中，f_s为信号的采样频率，nor(ω_k)为VMD分解所得IMF分量的归一化中心频率。

4.根据权利要求3所述的基于VMD和最大重叠离散小波包变换的齿轮信号去噪方法，其特征在于，所述步骤S3中重构信号的模糊熵值FE的具体步骤如下：

S301、设X＝[x(1),x(2),…,x(N)]为选定最优预设尺度K值下的重构信号，对X设定模式维数m，由X构造m维向量为：

式中，

为(n(i),n(i+1),···,n(i+m-1))的均值，i∈[1,N-m]；

S302、定义d_ij为任意两A^m间的最大欧式距离，并用指数函数u定义相似度D_ij：

D_ij＝u(d_ij,n,r)＝exp[-(d_ij/r)ⁿ]， (8)

式中，u为关于

的模糊隶属度函数；n和r分别为其边界的梯度和宽度；

S303、定义：

并重复步骤S301～S303)，算出m+1时候的Φ^m+1(n,r)；

S304、在以上步骤的基础上，定义模糊熵：

式中，m为模式维数；n为相似容限；N为数据长度。

5.根据权利要求4所述的基于VMD和最大重叠离散小波包变换的齿轮信号去噪方法，其特征在于，所述步骤S4具体包括：

S401、设u₁(t),u₂(t),…,u_n(t)为信号x(t)经VMD分解得到的IMF分量，基于IMF能量熵增量的敏感模态函数判别算法计算各IMF分量的能量：

式中，u_i(t)为第i个IMF分量；

S402、求解各IMF分量的归一化能量：

式中，N为IMF分量的总个数；

S403、计算各IMF分量的能量熵增量值：

Δq_i＝-p(i)log₂(p(i)) (13)。

6.根据权利要求5所述的基于VMD和最大重叠离散小波包变换的齿轮信号去噪方法，其特征在于，所述步骤S5具体包括：

根据时域互相关系数计算方法，设G_x、G_y分别为信号x_i和y_i功率谱，f_a为分析频率，则频域内信号x_i和y_i的互相关系数可表示为：

式中，∣ρ_f∣越大，表明两信号在频域上相关性越好；反之，∣ρ_f∣越小，说明两信号在频域上相关性越差。

7.根据权利要求6所述的基于VMD和最大重叠离散小波包变换的齿轮信号去噪方法，其特征在于，所述步骤S6中针对含有高频噪声干扰的IMF分量采用最大重叠离散小波包变换MODWPT方法对其进行去噪处理的具体步骤为：

MODWPT的分解系数由W_j,n＝{W_j,n,t,t＝0,…,N-1}表示，其中j为分解层数，n看作一个随j变化的频率索引，则计算出MODWPT的分解系数：

式中，若n mod 4＝0或3，则

若n mod4＝1或2，则

8.根据权利要求7所述的基于VMD和最大重叠离散小波包变换的齿轮信号去噪方法，其特征在于，所述步骤S7中去噪信号表示如下：

式中，x”(t)为最终去噪信号；x'(t)为经能量熵增量和频域互相关系数判定的对故障信息敏感的IMF分量所组成的重构信号；

为经MODWPT分解去噪后的IMF分量所构成的重构信号。

Images