抑制EMD端点效应的方法与流程

本发明涉及信号处理方法技术领域，尤其涉及一种抑制EMD端点效应的方法。

背景技术：

EMD(经验模态分解)是由美国国家航空航天局的Huang等人于1998年提出的一种自适应的信号处理方法，其分解得到的IMF分量具有真实的物理意义，因此被广泛地应用于处理非线性、非平稳信号。但是，其在理论上还存在着端点效应问题。为此，许多专家学者提出几种有效的方法：邓拥军等人提出利用神经网络对原始信号进行延拓的算法，该算法对短信号具有较好的端点效应抑制效果，但神经网络学习时间长，影响EMD的分解速度，不适合实时信号处理；赵近平等人提出了镜像延拓法，该算法只对原始信号进行一次延拓，在处理较长周期信号时可以得到理想的抑制效果，但当信号的端点不确定是否是极值点时，要截去一部分数据，在处理短信号时效果欠佳；此外，还有B样条经验模式分解、ARMA模型、窗函数、多项式拟合等方法，这些方法大多数都可以改善端点效应问题，但由于各种方法延拓结果因人而异，得到的结果也不尽相同，加之非平稳信号的随机性，最终结果都不可避免地存在一些问题。

技术实现要素：

本发明所要解决的技术问题是提供一种抑制EMD端点效应的方法，所述方法采用SVM(支持向量机模型)延拓和窗函数相结合的方法，弥补了SVM延拓依然找不到端点以及窗函数会改变原始信号的缺点，解决了EMD分解过程中端点效应引起的失真现象。

为解决上述技术问题，本发明所采取的技术方案是：一种抑制EMD端点效应的方法，其特征在于所述方法包括如下步骤：

利用SVM对原始信号两端进行延拓；

对未延拓的数据加矩形窗处理，延拓数据加海明窗处理；

对加窗处理后的数据进行EMD分解，把分解得到的IMF分量的两端延拓部分去掉，得到原始信号真实的IMF分量。

进一步的技术方案在于，所述方法还包括：利用传感器采集齿轮箱齿轮裂纹故障信号的步骤。

进一步的技术方案在于，所述方法还包括：通过得到的真实的IMF分量进行齿轮故障特征提取，进而判断齿轮故障的步骤。

优选的，所述的步骤利用SVM对原始信号两端进行延拓中支持向量机模型SVM的构造方法包括如下步骤：

参数设置：选择精度参数ε、惩罚参数C、损失函数e和适当的核函数k(x_i,x_j)；

最小化：

即在和0≤a_i,的约束条件下对a_i和a_i^*求解式(2)所示目标函数的最小值，即

得到最优解

构造支持向量回归模型：

式中，

优选的，所述的步骤利用SVM对原始信号两端进行延拓的方法包括如下步骤：

先考虑向后延拓的情况，对于一个给定的数据序列s(1),s(2),…,s(N)，其中N为采样点数，首先确定训练样本数l，按一定的规则产生一个训练集L＝{(x₁,y₁),…,(x_l,y_l)}；

按照上述方法构造如式(3)所示的回归模型，利用该回归模型得到边界外的第一个预测值s(N+1)，再将s(N+1)作为原始数据新的边界点，得到第2个数据序列的延拓值s(N+2)，以此类推，根据所需要延拓的数据个数M，得到全部延拓序列s(N),s(N+1),…,s(N+M)；

最后，按照向后延拓的方法对原始数据序列向前延拓，最终得到原始信号延拓数据。

优选的，所述的对未延拓的数据加矩形窗处理，延拓数据加海明窗处理的方法如下：

将得到的延拓数据分为延拓的和未延拓的两个部分；

向前、向后延拓的数据分别加海明窗处理，即延拓数据与海明窗函数相乘；

未延拓的数据加矩形窗处理，即该数据与矩形窗函数相乘；

将上述延拓加窗处理的数据整合到一起形成新的信号序列。

优选的，所述的对加窗处理后的数据进行EMD分解，把分解得到的IMF分量的两端延拓部分去掉，得到原始信号真实的IMF分量的方法如下：

将加窗处理后的信号序列记为x(t)，然后确定信号x(t)的所有极值点，接着用三次样条曲线将所有的极值点连接起来形成上、下包络线；

将上、下包络线的平均值记为m，求出h＝x(t)-m；

把h当作新的x(t)，重复上述步骤直到满足IMF条件时，记c₁＝h；

把c₁视为一个IMF，r＝x(t)-c₁，r当作新的x(t)，重复上述步骤得到其他的IMF分量，由此可得

式中，r_n称为残余函数，代表信号的平均趋势；将得到的IMF分量去除延拓的部分，得到可以反映原始信号的IMF分量。

采用上述技术方案所产生的有益效果在于：所述方法采用SVM延拓和窗函数相结合的方法，弥补了SVM延拓依然找不到端点以及窗函数会改变原始信号的缺点，解决了EMD分解过程中端点效应引起的失真现象。得到的IMF分量更加真实的反映了原始信号的特征，为后续齿轮故障诊断中精确地提取故障特征频率提供了帮助。

附图说明

图1是本发明第一个实施例所述方法的流程图；

图2是本发明第二个实施例所述方法的流程图；

图3是本发明第三个实施例所述方法的流程图；

图4是验证本发明所述方法时仿真信号时域图；

图5是验证本发明所述方法时仿真信号加窗时域图；

图6a是仿真信号未加窗EMD分解结果图；

图6b是仿真信号加窗EMD分解结果图；

图7是原始信号向前、向后延拓波形图；

图8a是SVM延拓加窗处理后的波形图；

图8b是SVM延拓加窗处理后的EMD分解图；

图9a是使用实验测试平台及传感器采集的实测信号图；

图9b是实测信号EMD分解图；

图10a是实测信号SVM延拓加窗处理后的波形图；

图10b是实测信号SVM延拓加窗处理后的EMD分解图。

具体实施方式

下面结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例仅仅是本发明的一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。

在下面的描述中阐述了很多具体细节以便于充分理解本发明，但是本发明还可以采用其他不同于在此描述的其它方式来实施，本领域技术人员可以在不违背本发明内涵的情况下做类似推广，因此本发明不受下面公开的具体实施例的限制。

实施例一

如图1所示，本发明公开了一种抑制EMD端点效应的方法，所述方法包括如下步骤：

S101：利用SVM对原始信号两端进行延拓；

利用SVM对原始信号两端进行延拓的具体方法如下：

S1011：选择精度参数ε、惩罚参数C、损失函数e和适当的核函数k(x_i,x_j)。

S1012：最小化

即在和0≤a_i,的约束条件下对a_i和a_i^*求解式(2)所示目标函数的最小值，即

得到最优解

S1013：构造回归模型：

式中，

S1014：构造支持向量回归模型后，就可以利用它对数据点分别向前和向后进行延拓。先考虑向后延拓的情况，对于一个给定的数据序列s(1),s(2),…,s(N)，其中N为采样点数，首先确定训练样本数l，按一定的规则产生一个训练集L＝{(x₁,y₁),…,(x_l,y_l)}。

S1015：按照上述方法构造如式(3)所示的回归模型，利用该回归模型就可以得到边界外的第一个预测值s(N+1)，再将s(N+1)作为原始数据新的边界点，就可以得到第2个数据序列的延拓值s(N+2)，以此类推，根据所需要延拓的数据个数M，就可以得到全部延拓序列s(N),s(N+1),…,s(N+M)。

S1016：按照上面所述方法对原始数据序列向前延拓，最终得到原始信号延拓数据。

S102：对未延拓的数据加矩形窗处理，延拓数据加海明窗处理，具体方法如下：

S1021：将得到的延拓数据分为延拓的和未延拓的两个部分。

S1022：向前、向后延拓的数据分别加海明窗处理，即延拓数据与海明窗函数相乘。

S1023：未延拓的数据加矩形窗处理，即该数据与矩形窗函数相乘。由于原始信号为未延拓的那一部分，因此不会改变原始信号的特征。

S1024：将上述延拓加窗处理的数据整合到一起形成新的信号序列。

S103：对加窗处理后的数据进行EMD分解，把分解得到的IMF分量的两端延拓部分去掉，得到原始信号真实的IMF分量，具体方法如下：

S1031：将上述处理后的数据序列记为x(t)，然后确定信号x(t)的所有极值点，接着用三次样条曲线将所有的极值点连接起来形成上、下包络线(包络线要包络所有极值点)。

S1032：将上、下包络线的平均值记为m，求出h＝x(t)-m。

S1033：把h当作新的x(t)，重复上述步骤直到满足IMF条件时，记c₁＝h。

S1034：把c₁视为一个IMF，r＝x(t)-c₁，r当作新的x(t)，重复上述步骤得到其他的IMF分量，由此可得

式中，r_n称为残余函数，代表信号的平均趋势。

S1035：将得到的IMF分量去除延拓的部分，得到可以反映原始信号的IMF分量。

实施例二

如图2所示，本发明实施例公开了一种抑制EMD端点效应的方法，所述方法包括如下步骤：

S201：利用传感器采集齿轮箱齿轮裂纹故障信号；

S202：利用SVM对原始信号两端进行延拓，具体方法如下：

利用SVM对原始信号两端进行延拓的具体方法如下：

S2021：选择精度参数ε、惩罚参数C、损失函数e和适当的核函数k(x_i,x_j)。

S2022：最小化

即在和0≤a_i,的约束条件下对a_i和a_i^*求解式(2)所示目标函数的最小值，即

得到最优解

S2023：构造回归模型：

式中，

S2024：构造支持向量回归模型后，就可以利用它对数据点分别向前和向后进行延拓。先考虑向后延拓的情况，对于一个给定的数据序列s(1),s(2),…,s(N)，其中N为采样点数，首先确定训练样本数l，按一定的规则产生一个训练集L＝{(x₁,y₁),…,(x_l,y_l)}。

S2025：按照上述方法构造如式(3)所示的回归模型，利用该回归模型就可以得到边界外的第一个预测值s(N+1)，再将s(N+1)作为原始数据新的边界点，就可以得到第2个数据序列的延拓值s(N+2)，以此类推，根据所需要延拓的数据个数M，就可以得到全部延拓序列s(N),s(N+1),…,s(N+M)。

S2026：按照上面所述方法对原始数据序列向前延拓，最终得到原始信号延拓数据。

S203：对未延拓的数据加矩形窗处理，延拓数据加海明窗处理，具体方法如下：

S2031：将得到的延拓数据分为延拓的和未延拓的两个部分。

S2032：向前、向后延拓的数据分别加海明窗处理，即延拓数据与海明窗函数相乘。

S2033：未延拓的数据加矩形窗处理，即该数据与矩形窗函数相乘。由于原始信号为未延拓的那一部分，因此不会改变原始信号的特征。

S2034：将上述延拓加窗处理的数据整合到一起形成新的信号序列。

S204：对加窗处理后的数据进行EMD分解，把分解得到的IMF分量的两端延拓部分去掉，得到原始信号真实的IMF分量，具体方法如下：

S2041：将上述处理后的数据序列记为x(t)，然后确定信号x(t)的所有极值点，接着用三次样条曲线将所有的极值点连接起来形成上、下包络线(包络线要包络所有极值点)。

S2042：将上、下包络线的平均值记为m，求出h＝x(t)-m。

S2043：把h当作新的x(t)，重复上述步骤直到满足IMF条件时，记c₁＝h。

S2044：把c₁视为一个IMF，r＝x(t)-c₁，r当作新的x(t)，重复上述步骤得到其他的IMF分量，由此可得

式中，r_n称为残余函数，代表信号的平均趋势。

S2045：将得到的IMF分量去除延拓的部分，得到可以反映原始信号的IMF分量。

实施例三

如图3所示，本发明实施例公开了一种抑制EMD端点效应的方法，所述方法包括如下步骤：

S301：利用传感器采集齿轮箱齿轮裂纹故障信号。

S302：利用SVM对原始信号两端进行延拓，具体方法如下：

S3021：选择精度参数ε、惩罚参数C、损失函数e和适当的核函数k(x_i,x_j)。

S3022：最小化

即在和0≤a_i,的约束条件下对a_i和a_i^*求解式(2)所示目标函数的最小值，即

得到最优解

S3023：构造回归模型：

式中，

S3024：构造支持向量回归模型后，就可以利用它对数据点分别向前和向后进行延拓。先考虑向后延拓的情况，对于一个给定的数据序列s(1),s(2),…,s(N)，其中N为采样点数，首先确定训练样本数l，按一定的规则产生一个训练集L＝{(x₁,y₁),…,(x_l,y_l)}。

S3025：按照上述方法构造如式(3)所示的回归模型，利用该回归模型就可以得到边界外的第一个预测值s(N+1)，再将s(N+1)作为原始数据新的边界点，就可以得到第2个数据序列的延拓值s(N+2)，以此类推，根据所需要延拓的数据个数M，就可以得到全部延拓序列s(N),s(N+1),…,s(N+M)。

S3026：按照上面所述方法对原始数据序列向前延拓，最终得到原始信号延拓数据。

S303：对延拓后的数据中间加矩形窗，延拓数据加海明窗处理，具体方法如下：

S3031：、将得到的延拓数据分为延拓的和未延拓的两个部分。

S3032：向前、向后延拓的数据分别加海明窗处理，即延拓数据与海明窗函数相乘。

S3033：未延拓的数据加矩形窗处理，即该数据与矩形窗函数相乘。由于原始信号为未延拓的那一部分，因此不会改变原始信号的特征。

S3034：将上述延拓加窗处理的数据整合到一起形成新的信号序列。

S304：对加窗处理后的数据进行EMD分解，把分解得到的IMF分量的两端延拓部分去掉，得到原始信号真实的IMF分量，具体方法如下：

S3041：将上述处理后的数据序列记为x(t)，然后确定信号x(t)的所有极值点，接着用三次样条曲线将所有的极值点连接起来形成上、下包络线(包络线要包络所有极值点)。

S3042：将上、下包络线的平均值记为m，求出h＝x(t)-m。

S3043：把h当作新的x(t)，重复上述步骤直到满足IMF条件时，记c₁＝h。

S3044：把c₁视为一个IMF，r＝x(t)-c₁，r当作新的x(t)，重复上述步骤得到其他的IMF分量，由此可得

式中，r_n称为残余函数，代表信号的平均趋势。

S3045：将得到的IMF分量去除延拓的部分，得到可以反映原始信号的IMF分量。

S305：通过得到的真实的IMF分量进行齿轮故障特征提取。

验证方法：

为了说明EMD的端点效应问题和窗函数法虽然可以解决端点效应问题，但是它却改变了原始信号特征的问题，本发明以仿真信号x＝sin(2*pi*10*t)+sin(2*pi*50*t)+sin(2*pi*100*t)+0.1*randn(1,length(t))为例来说明。该信号是由一段正弦和随机噪声叠加而成的信号，其中正弦由三个频率成分组成，分别为10Hz、50Hz和100Hz。

图4所示为仿真信号波形图，图5所示为仿真信号加海明窗之后的波形图，由图可以看出原始信号已经发生了改变。图6a所示为仿真信号EMD分解波形图，可以看出信号存在严重的端点效应，并且向内扩散“污染”了整个信号序列，对后续信号的处理会产生不利的影响，图6b为加窗之后的EMD分解波形图，可以看出端点效应得到了有效抑制，但是却改变了原始信号，这对后续信号的处理也是不利的。

为了验证本发明的有效性，仍然以上面的仿真信号为例，其时域波形如图1，未进行任何处理的EMD分解结果如图6a所示，可以看出各个IMF存在严重的端点效应问题。为了解决端点效应，采用基于SVM延拓和窗函数相结合的方法进行处理。实验中将原始信号序列分成100个训练组，并对信号分别向前、后各延拓50个点，原始信号延拓后的结果见图7，延拓之后加窗(中间加矩形窗，两端加海明窗)的结果见图8a，对应的EMD分解结果见图8b。

与未加任何处理的EMD分解结果相比较，可以看出图8b中各个IMF分量在端点处出现的发散现象均得到了改善，表明本文提出的方法可以有效的抑制EMD分解中出现的端点效应问题，并且各个IMF分量也都可以反映原始信号的特征，基本符合实际情况。

为了进一步地确认该发明的有效性，本发明以实测齿轮箱齿轮裂纹故障信号进行分析。实验采用江苏千鹏诊断系统—QPZZ-II，其中输入轴转速为990r/min，采样频率为20000Hz，采样点数为262144，故障齿轮为3.2mm的裂纹程度。为了方便研究截取前800个点，并采用SVM向前、向后延拓50个点，利用本文介绍的方法处理研究。

图9a所示为测试信号未进行任何处理的时域波形，直接对其EMD分解的(本文取了前4个IMF分量)结果见图9b，由图可以看出信号两端存在着严重的端点发散现象，并向内“污染”了整个数据序列，导致分解得到的IMF分量不能真实的反映原始信号特征。

为了解决端点效应问题，将该数据序列按本发明所述方法处理，其SVM延拓加窗处理结果见图10a所示。由图可以看出处理后的信号与原始信号存在着显著地区别，这是在信号两端加窗处理的原因，因此最后EMD分解得到的各个IMF分量要剔除两端延拓部分，这样就可以真实的反映原始信号特征。处理后的数据EMD分解结果见图10b，由图可以看出各个IMF分量的端点效应都得到了有效的抑制，进一步的验证了本发明的有效性。

本发明所述方法利用SVM延拓和窗函数相结合的方法抑制了齿轮箱齿轮裂纹故障信号在进行EMD分解时出现的端点效应问题。采用SVM延拓和窗函数相结合的方法，弥补了SVM延拓依然找不到端点以及窗函数会改变原始信号的缺点。解决了由EMD端点效应引起的齿轮箱齿轮裂纹故障信号分解后的失真现象。本发明抑制端点效应之后得到的IMF分量可以更加真实的反映原始信号特征，可以准确地分离出故障频率及其倍频，进一步准确地判断出齿轮裂纹故障。