CN104931197A - 基于eemd的自动平衡机振动信号处理方法 - Google Patents

Abstract

本发明公开了一种基于EEMD的自动平衡机振动信号处理方法，首先采用数学形态学滤波器对提取到的振动信号进行预处理，消除信号中存在的脉冲干扰，然后针对经验模态分解中存在的模态混叠问题，采用总体经验模态分解的方法，提取出不平衡量的频率和相位。包括以下步骤：一、利用数学形态学方法对传感器测得的振动信号进行预处理，滤除脉冲干扰；二、利用集合经验模态分解法对经过步骤一数学形态学方法处理后的信号进行分解，得到信号的本征模态函数分量；三、通过希尔伯特变换得到振动信号幅值相位。

Description

基于EEMD的自动平衡机振动信号处理方法

技术领域

本发明属于自动平衡机振动信号处理相关领域，具体涉及一种基于EEMD(集合经验模态分解)的自动平衡机振动信号处理方法。

背景技术

旋转机械遍及生活生产各个领域，在众多机械设备中占比例巨大。但由于设计和制造的原因导致转子本身的质量分布不均匀，机器在生产组装的过程中产生误差，机器在运行中对转子产生不同程度的磨损，不同的环境影响等原因，转子在实际旋转中的主惯性轴和旋转轴线有一定的偏差，这种现象称为转子质量不平衡。不平衡的转子在旋转的过程中，不平衡质量就产生了不平衡力，不平衡力会使旋转机械产生振动和噪声，长时间下去机械就会损坏严重。因此及时发现和处理不平衡问题至关重要。

对转子进行不平衡测试的技术就是针对于转子在实际的运转中容易产生不平衡振动的现象产生的，而能够通过测试技术确定不平衡量的大小和位置并对其进行校正的仪器就是所谓的平衡机。在生产中，工作人员可以通过平衡机时刻监测转子的运转情况，当发生不平衡振动的情况时，可以根据平衡机测量出的不平衡量及时对转子进行平衡校正。衡量平衡机测试技术的性能好坏和精确度的标准是：对故障转子校正后不平衡量的减少的比例，或者校正后整个运转的机器可以达到的最小的剩余不平衡量的多少。本发明主要针对动平衡机测量***中振动信号处理方法进行研究，也就是精确地对不平衡量的大小和位置进行确定。

关于自动平衡机振动信号处理方法的研究，现有常用的方法主要为频域分析方法和时频分析方法。频域分析方法主要是利用傅里叶变换的方法。傅里叶变换及其快速算法(fast Fourier transform,FFT)是信号处理的理论基础，在振动信号包含成分比较简单的时候应用比较广泛，算法实现简单，运算速度快，也能比较准确地得到信号中的各个频率分量，但是在傅里叶变换中，采用的是时域截断的方法，其抗干扰能力有所降低，会引起能量泄露和栅栏现象，影响了提取精度。研究人员对于傅里叶变换方法的改进主要手段是给信号加窗，增加窗序列的长度，并合理地选择采样点数和采样频率。时频分析方法能够同时获得信号在时域和频域上的分布特征。时频分析方法方面，主要有短时傅里叶变换、分数阶傅里叶变换、Wigner-Ville分布和小波变换等方法。小波变换方法也是在傅里叶变换的基础上产生的，它能够克服傅里叶变换对复杂信号分析方面存在的不足和缺陷，对信号进行分解和重构，采用预先设定的小波基函数，对原始信号进行局部特征信息分析，通过不同尺度对应的带宽将信号进行分解，能够很好的去除干扰噪声。但是在小波变换中，小波基函数需要根据信号进行调整，如何满 足不同信号的分析要求是个关键点，在整个信号分析过程中选定的小波基是无法改变的，若小波基仅为局部最优，在整体范围内效果并不一定十分理想。小波基的设定不同，对同一个信号的分析结果差别很大，限制了小波变换的应用。希尔伯特-黄变换(Hilbert-Huang transform,HHT)方法是在傅立叶变换的基础上发展起来的线性和稳态谱分析，其处理的过程是先对待处理的信号进行经验模态分解(Empirical Mode Decomposition,EMD)，然后将分解后的每个本征模态函数(Intrinsic Mode Function,IMF)进行希尔伯特变换，就可以对信号进行时频分析。通过经验模态分解可以得到想要频率成分的信号，滤除异频噪声的干扰，此方法在机械故障诊断方面运用的非常多，但是对于振动信号中含有大量脉冲干扰的情况下，经验模态分解过程中容易出现模态混叠问题。

发明内容

为解决现有方法存在的以上不足，本发明提供一种基于EEMD的自动平衡机振动信号处理方法，首先采用数学形态学滤波器对提取到的振动信号进行预处理，消除信号中存在的脉冲干扰，然后针对经验模态分解中存在的模态混叠问题，采用总体经验模态分解的方法，提取出不平衡量的频率和相位。

本发明的目的是通过以下方案实现的：一种基于EEMD的自动平衡机振动信号处理方法，包括采用数学形态学滤波器消除不平衡振动信号中的脉冲干扰，并利用经验模态分解方法提取不平衡量的频率、幅值和相位；其步骤如下：

步骤一、利用数学形态学方法对传感器测得的振动信号进行预处理，滤除脉冲干扰。数学形态学以集合来描述目标信号，通过结构元素在信号内部移动收集信号的信息，对信号进行匹配，从而达到提取信号、保持细节和抑制噪声的目的。利用数学形态学开闭组合滤波器对传感器测得的振动信号进行预处理，能够滤除信号中的脉冲干扰。

步骤二、利用集合经验模态分解法对经过步骤一数学形态学方法处理后的信号进行分解，得到信号的本征模态函数(IMF)分量：利用集合经验模态分解法结合白噪声均匀的分布性，通过向经过所述步骤一数学形态学方法处理后的信号添加高斯白噪声，然后对加噪信号进行多次经验模态分解(EMD)，从而在不同尺度上，消除原信号由于跳变产生的间断，然后对获得的本征模态函数分量求取平均值，可以去除添加的白噪声，获得输入信号的本征模态函数组。

步骤三、通过希尔伯特变换得到振动信号幅值相位：将经过所述步骤二获得的本征模态函数组通过希尔伯特变换得到每个本征模态函数(IMF)分量的瞬时频率和瞬时幅值。

其中，步骤三通过希尔伯特变换提取振动信号幅值相位的具体步骤如下：

对本征模态函数h_i(t)做希尔伯特变换，可以得到：

${\hat{h}}_i\left(t\right)=\frac{1}{\mathrm{\π}}{\∫}_{-\mathrm{\∞}}^{\mathrm{\∞}}\frac{h_i\left(\mathrm{\τ}\right)}{t-\mathrm{\τ}}d\mathrm{\τ}$

从而得到信号x(t)的希尔伯特谱可以表示为：

$H(\mathrm{\ω},t)=\mathrm{Re}\underset{i=1}{\overset{n}{\Σ}}{a}_i\left(t\right)e^{j\∫{\mathrm{\ω}}_i\left(t\right)dt}$

再定义希尔伯特边际谱：

$h\left(\mathrm{\ω}\right)={\∫}_0^{T}H(\mathrm{\ω},t)dt$

H(ω,t)表示信号的幅值在整个频率段上与时间和频率的分布关系；h(ω)反映了信号的幅值随频率的变化情况。

本发明具有以下优点：

傅里叶变换是信号处理领域经典的方法之一，在不平衡振动信号提取过程中，主要利用傅里叶变换的快速算法。但是变换中采取的时域截断方法会产生能量泄漏使得精度下降，而且，傅里叶变换在整体上将信号分解为不同的频率分量，而没有局域性信息，无法得知某种频率分量出现、消失在哪些时间内，只适用于分析信号组成分量的频率不随时间变化的平稳信号。

本发明采用数学形态学滤波器对提取到的振动信号进行预处理，消除信号中存在的脉冲干扰，然后针对经验模态分解中存在的模态混叠问题，采用总体经验模态分解的方法，提取出不平衡量的频率和相位，该方法测量的精确较高，计算速度快。

附图说明

图1本发明提出的方法流程图

图2经验模态分解(EMD)流程图

图3集合经验模态分解(EEMD)流程图

图4电机转速为3600r/min时，原始含噪声振动信号幅值谱

图5电机转速为3600r/min时，数学形态学去噪后的振动信号幅值谱

图6电机转速为3600r/min时，EEMD得到信号的各IMF分量及余项的幅值谱，其中，图6(a)为IMF1分量幅值谱，图6(b)为IMF2分量幅值谱，图6(c)为IMF3分量幅值谱，图6(d)为IMF4分量幅值谱，图6(e)为IMF5分量幅值谱，图6(f)为IMF6分量幅值谱，图6(g)为IMF7分量幅值谱，图6(h)为余项幅值谱

图7电机转速为3600r/min时，各IMF分量瞬时相位谱，其中，图7(a)为IMF1分量的瞬时相位，图7(b)为IMF2分量的瞬时相位，图7(c)为IMF3分量的瞬时相位，图7(d)为IMF4分量的瞬时相位，图7(e)为IMF5分量的瞬时相位，图7(f)为IMF6分量的瞬时相位，图7(g)为IMF7分量的瞬时相位，图7(h)为IMF8分量的瞬时相位，图7(i)为IMF9分量的瞬时相位

图8电机转速为3600r/min时，处理后的振动信号的希尔伯特三维谱

图9电机转速为3600r/min时，处理后的振动信号的希尔伯特时频谱

图10电机转速为7200r/min时，处理后的振动信号的希尔伯特三维谱

图11电机转速为12000r/min时，处理后的振动信号的希尔伯特三维谱

具体实施方式

步骤一、利用数学形态学方法对传感器测得的振动信号进行预处理，滤除脉冲干扰。

设原始信号f(n)是定义在F(n)＝(0,1,...,N-1)上的离散函数，定义结构元素g(n)为G＝(0,1,...,M-1)上的离散函数，且N≥M，则：

f(n)对g(n)的开运算定义为：

f(n)对g(n)的闭运算定义为：

$f\left(n\right)\·g\left(n\right)=(f\left(n\right)\⊕[-g\left(n\right)\left]\right)\mathrm{\Θ}[-g\left(n\right)]---\left(2\right)$

从式(1)和式(2)可以看出，闭运算相当于开运算的对偶运算。开运算和闭运算具有幂等性，即运算的结果是不受运算次数的影响，都是相同的。

开、闭运算可以理解为由结构元素g(n)构成的“球”在信号内部“滚动”。开运算可以使含有毛刺的信号变得光滑，抑制信号中存在的正向脉冲噪声；闭运算则可以填平信号轮廓上的“沟壑”，可以抑制信号中存在的负向脉冲噪声。因此，利用形态学开、闭运算构成的数学形态滤波器对于信号中存在的正负脉冲信号干扰可以进行有效地抑制。

根据开、闭运算的特质，为了同时去除不平衡振动信号中存在的正、负脉冲噪声，滤波器的设计采用级联形式，包括滤除负脉冲干扰的形态开滤波器和滤除正脉冲干扰的形态闭滤波器。Maragos采用同型结构元素，调整开、闭运算级联顺序，定义了如式(3)和(4)所示的开-闭(open—closing)滤波器和闭-开(close—opening)滤波器滤除信号的正负脉冲噪声。

F_OC(f(n))＝(fοg·g)(n)  (3)

F_CO(f(n))＝(f·gοg)(n)  (4)

经滤波处理后的信号存在一定的误差，即偏倚现象。原因是开运算的收缩性导致通过闭-开滤波器之后的信号的输出偏小，即统计左偏移；闭运算的扩张性导致通过闭-开滤波器之后的信号输出偏大，即统计右偏移，单独一种滤波器使用会造成信号处理的不对称性。为了有效地抑制噪声，同时避免出现偏倚现象，本发明将开-闭和闭-开两种滤波器级联在一起，形成组合级联形式的形态滤波器，对通过两个滤波器的信号求取其平均值，用来消除误差。具体表达式如式(5)所示:

$y\left(n\right)=\frac{1}{2}\[F_{OC}\left(f\left(n\right)\right)+F_{CO}\left(f\left(n\right)\right)\]---\left(5\right)$

图4是原始含噪声振动信号，可以看出信号中含有大量正负脉冲噪声，这些噪声的存在会干扰频率，影响经验模态分解结果。利用数学形态学滤波器对不平衡振动信号进行去早处理，得到如图3所示的信号模型，可以看到脉冲干扰得到有效地抑制，表明了数学形态学在抑制脉冲噪声干扰上的有效性。

步骤二、利用集合经验模态分解法(Ensemble Empirical Mode Decomposition,EEMD)对经过步骤一数学形态学方法处理后的信号进行分解，得到信号的本征模态函数(IMF)分量：

经验模态分解(EMD)法的提出是基于这样的假设：任何信号都包含不同的单一的本征振动模态分量，无论线性还是非线性的分量都要含有相同数量的极值点和过零点，并且过零点与相邻过零点之间必须只有一个极值点，每个模态都互相独立。这些假设如果成立，那么每个信号都可以被分解成为几个本征模态函数(IMF)，且每个本征模态函数必须满足以下的定义：(1).整个序列内的极值点和过零点数量相差不超过一个；(2).在任意时间点上，由局部最大值定义的(上)包络线和由局部最小值定义的(下)包络线的均值为零。

集合经验模态分解(EEMD)法利用白噪声均匀的分布性，一方面顺滑了由于跳变产生的异常扰动；另一方面添入了均匀分布的随机尺度，从而有效的抑制模态混叠现象。因此，EEMD方法的实质为：通过向信号添加高斯白噪声，然后进行多次EMD分解，从而在不同尺度上，消除原信号由于跳变产生的间断，然后对获得的本征模态函数分量求取平均值，以此来去除添加的白噪声，达到更好的分解效果。该方法作为EMD的改进型方法，不仅可以成功抑制模态混叠，还可以使得求得的本征模态函数更具实际意义。EEMD算法的具体实现过程如下：

(1).初始化。设置总体平均次数M和噪声幅值系数k，并且使m＝1；

(2).执行第m次EMD，详细过程如下：

a.将高斯白噪声n_m(t)添加到输入信号x(t)中：

x_m(t)＝x(t)+kn_m(t)  (6) 

b.用EMD分解加噪信号x_m(t)，得到一组IMF分量，记作c_n,m，表示第m次EMD分解得到的第n个本征模态函数分量。

c.若m＜M，则返回步骤(1)，令m＝m+1，直至m＝M。

(3).计算M次分解得到的IMF c_n,m的总体均值：

$y_n=\frac{1}{M}\underset{m=1}{\overset{M}{\Σ}}{c}_{n,m}---\left(7\right)$

式(7)中，n＝1,2,…,N，m＝1,2,…,M。

(4).把M次分解得到的n个IMF分量的总体均值y_n作为最终的结果，从而获得输入信号 的本征模态函数组Y_IMF＝[y₁,y₂,…y_n]

EEMD分解后的输入信号x(t)可表示为Y_IMF＝[y₁,y₂,…y_n]之后，即

$x\left(t\right)=\underset{n=1}{\overset{N}{\Σ}}{y}_n\left(t\right)+r_n---\left(8\right)$

在EMD分解的过程中，只要对分解参数设置完成，其每次分解得到的IMF都是恒定不变的，但是由于EEMD在分解之前要添加白噪声，而每次添加的白噪声是随机的，幅度并不完全相同，所以分解出来的IMF波形每次都会有些许的区别。因此幅度的改变会对每次分解出的IMF结果造成影响。此外，即使设置的参数一样，可是添加的白噪声是随机产生的，因此计算出来的结果仍然存在差异。但是，随着总体平均次数的越来越多，达到一定数量时，最终会相互抵消每次添加的白噪声，直至完全消除。因此，当总体平均的次数达到一定数量时，不仅能有效消除白噪声造成的干扰，而且能分解出更加精确的IMF，从而真实的反映出信号内部的特征。

Z.Wu给出了白噪声引入的数据误差经验公式:

${\mathrm{\ϵ}}_n=\frac{\mathrm{\ϵ}}{\sqrt{M}}---\left(9\right)$

式中：ε为白噪声的幅度，ε_n为最终误差的标准差，M为总体平均次数。从式中可以知道，当总体平均次数M一定时，误差ε_n会随着噪声幅度ε的升高而变大。

因此，在实际应用中可根据精度的要求以及计算时间的长短来对变量进行控制，使最终分解得到的结果更加保真。

步骤三、通过希尔伯特变换得到每个IMF分量的瞬时频率和瞬时幅值

通过步骤二后，经验模态分解利用良好的自适应性，将信号分解为若干个具有物理意义的IMF分量之和，这样使得瞬时频率具有了实际的物理意义，从而可以通过希尔伯特变换得到每个IMF分量的瞬时频率和瞬时幅值。

对本征模态函数(IMF)h_i(t)做希尔伯特变换，可以得到

${\hat{h}}_i\left(t\right)=\frac{1}{\mathrm{\π}}{\∫}_{-\mathrm{\∞}}^{\mathrm{\∞}}\frac{h_i\left(\mathrm{\τ}\right)}{t-\mathrm{\τ}}d\mathrm{\τ}---\left(10\right)$

构造解析信号 

于是得到幅值函数为

$a_i\left(t\right)=\sqrt{h_i^2\left(t\right)+\widehat{{h_i}^2}\left(t\right)}---\left(12\right)$

相位函数为 

进一步可以求出瞬时频率为

可以得到

式中省略了残量r_n；Re表示取实部；信号x(t)的希尔伯特谱可以表示为

$H(\mathrm{\ω},t)=\mathrm{Re}\underset{i=1}{\overset{n}{\Σ}}{a}_i\left(t\right)e^{j\∫{\mathrm{\ω}}_i\left(t\right)dt}---\left(16\right)$

再定义希尔伯特边际谱

$h\left(\mathrm{\ω}\right)={\∫}_0^{T}H(\mathrm{\ω},t)dt---\left(17\right)$

式(17)中，T表示信号的总长度；H(ω,t)表示信号的幅值在整个频率段上与时间和频率的分布关系；h(ω)反映了信号的幅值随频率的变化情况。

利用希尔伯特变换求出各IMF分量的瞬时相位谱，对瞬时相位求导得到信号的瞬时频率，继而得到信号的希尔伯特三维谱和希尔伯特时频谱。根据式(10)和式(13)，求得各IMF瞬时相位谱如图7所示，其中θ1～θ9表示IMF1～IMF9的瞬时相位，相位单位为rad。通过零点的相位值，可以得到每个IMF分量的初相位。

步骤四、算法验证

本发明在对于算法性能验证上，用Matlab软件仿真对信号进行滤除脉冲干扰，利用EEMD对信号进行分解成IMF分量，再利用希尔伯特变换得出信号的瞬时幅值和瞬时相位。

(1)为验证步骤一提到的数学形态学滤波器的有效性，用Matlab软件对振动信号进行仿真实验。

图4为原始含噪声信号，利用数学形态滤波器对信号进行去噪处理。滤波器采用形态开闭组合滤波器，结构体元素采用正弦型，幅值为0.01V，长度为5个采样点。处理结果如图5所示。从图5和图4的对比可以观察出，原始含噪声的信号中，经过数学形态学去噪后，正负脉冲干扰得到有效地抑制，图像中“毛刺”减少，曲线较之图4变得平滑。因此可以证明数学形态学在抑制脉冲噪声干扰上的有效性。

经过形态学开闭组合滤波器处理过的信号，没有了频率成分比较复杂的正负脉冲信号，同时也不会因为脉冲干扰的存在造成的时间尺度上的跳变，为以后的经验模态分解创造了有利条件。

(2)验证集合经验模态分解方法在抑制模态混叠上的有效性

对图5中经过数学形态学滤波后的信号进行所述步骤二中的EEMD分解，可以得到各个分量如图6中的各IMF分量，EEMD方法很好地将高频扰动信号和低频信号分开，克服了模 态混叠现象，得到的IMF分量也具有了明显的物理意义。同时保证了经验模态分解的自适应特性，显示了EEMD的优越性。

(3)验证集对各IMF分量进行希尔伯特变换得到瞬时幅值和瞬时相位的有效性

利用希尔伯特变换求出各IMF分量的瞬时相位谱，对瞬时相位求导得到信号的瞬时频率，继而得到信号的希尔伯特三维谱和希尔伯特时频谱。根据式(10)和式(13)，求得各IMF瞬时相位谱如图7所示，其中θ1～θ9表示IMF1～IMF9的瞬时相位，相位单位为rad。通过零点的相位值，可以得到每个IMF分量的初相位。

已知电机转速为3600r/min，即与转子同频信号频率为60Hz，采样点为1400点。不平衡量为与转子同频的IMF分量，振幅为0.35V，初相位为π/6。由图6可以看出IMF2为与转子同频的不平衡分量，测得的幅值为0.3504V。根据图7的各个IMF分量的相位曲线可以得到，IMF2的相位为0.5226rad。其余的频率分量是由于机械本身和传感器内部噪声产生，与本发明中所需的不平衡分量在统计上相互独立，并无影响。

实验测得的结果证明了本发明提供的基于EEMD的自动平衡机振动信号处理方法的有效性。

图8为采用本发明的方法处理后的振动信号的希尔伯特三维谱图。通过三维谱图，可以清晰得出各个频率对应的幅值和时间。

图9表示的是采用本发明的方法处理后的振动信号的希尔伯特时频谱，从图中可以看出各个频率在时间上的分布情况。

改变电机转速，分别取7200r/min和12000r/min，分别对其进行所述步骤一、二、三的处理，当电机转数为7200r/min时，转子频率为120Hz，可以从图10中看出频率为120Hz处不平衡信号分量的相位和幅值；当电机转数为12000r/min时，转子频率为200Hz，可以从图11中看出频率为200Hz处不平衡信号分量的相位和幅值。

可以从以上3组实测振动信号所得结果来看，EEMD方法能够比较准确地提取出与转子同频的不平衡量幅值和相位，分析结果和误差如表1所示：

表1 不同信噪比下信号幅值相位测量值

由表1可以看出，对于不同转速的转子，即不平衡量位于不同频率的特征量提取方面，EEMD方法得到的相位和幅值在数值上很接近真实值，相对误差较小，证明了方法的有效性。具有实际应用意义。

Claims (2)

1.一种基于EEMD的自动平衡机振动信号处理方法，其特征在于，采用数学形态学滤波器消除不平衡振动信号中的脉冲干扰，并采用经验模态分解方法提取不平衡量的频率、幅值和相位；包括以下步骤：

步骤一、利用数学形态学方法对传感器测得的振动信号进行预处理，滤除脉冲干扰：利用数学形态学开闭组合滤波器对传感器测得的振动信号进行预处理，滤除信号中的脉冲干扰；

步骤二、利用集合经验模态分解法对经过所述步骤一数学形态学方法处理后的信号进行分解，得到信号的本征模态函数分量：利用集合经验模态分解法，通过向经过所述步骤一数学形态学方法处理后的信号添加高斯白噪声，然后对加噪信号进行多次经验模态分解，从而在不同尺度上，消除原信号由于跳变产生的间断，然后对获得的本征模态函数分量求取平均值，去除添加的白噪声，获得输入信号的本征模态函数组；

步骤三、通过希尔伯特变换得到振动信号幅值相位：将经过所述步骤二获得的本征模态函数组通过希尔伯特变换得到每个本征模态函数分量的瞬时频率和瞬时幅值。

2.如权利要求1所述的一种基于EEMD的自动平衡机振动信号处理方法，其特征在于，所述步骤三通过希尔伯特变换提取振动信号幅值相位的具体步骤如下：

对本征模态函数h_i(t)做希尔伯特变换，可以得到：

${\hat{h}}_i\left(t\right)=\frac{1}{\mathrm{\π}}{\∫}_{-\mathrm{\∞}}^{\mathrm{\∞}}\frac{h_i\left(\mathrm{\τ}\right)}{t-\mathrm{\τ}}d\mathrm{\τ}$

从而得到信号x(t)的希尔伯特谱可以表示为：

$H(\mathrm{\ω},t)=\mathrm{Re}\underset{i=1}{\overset{n}{\Σ}}{a}_i\left(t\right)e^{{}^j\∫{\mathrm{\ω}}_i\left(t\right)dt}$

再定义希尔伯特边际谱：

$h\left(\mathrm{\ω}\right)={\∫}_0^{T}H(\mathrm{\ω},t)dt$

H(ω,t)表示信号的幅值在整个频率段上与时间和频率的分布关系；h(ω)反映了信号的幅值随频率的变化情况。