利用粒子滤波与谱峭度的滚动轴承故障诊断方法技术

本发明专利技术公开了一种利用粒子滤波与谱峭度的滚动轴承故障诊断方法，涉及粒子滤波降噪处理以及谱峭度计算。该方法提出了在快速谱峭度的基础上，利用粒子滤波降噪处理提高信噪比，解决了快速谱峭度在低信噪比情况下可行性差的问题。首先建立振动信号的状态方程；然后提取背景噪声，将其与状态方程之和作为观测方程；最后联立状态方程与观测方程建立状态空间模型；采用粒子滤波算法对信号重新估计，得到新的序列即是降噪之后的信号；最后用快速谱峭度的方法，获取最佳分析频带，得出故障频率。本发明专利技术的滚动轴承故障诊断方法，降低故障信号中的噪声干扰，提高信噪比，实现了滚动轴承早期微弱故障的诊断。

【技术实现步骤摘要】
利用粒子滤波与谱峭度的滚动轴承故障诊断方法
本专利技术涉及一种利用粒子滤波与谱峭度的滚动轴承故障诊断方法，涉及粒子滤波降噪处理以及谱峭度计算。
技术介绍
滚动轴承是机械设备尤其是旋转机械中重要的零部件之一，但是其寿命的随机性较大，容易损坏，目前还无法准确预测其寿命的长短。近年来，科技和工业生产都在快速的发展，机械设备逐渐向高速化、大型化和自动化发展，这在提高生产力的同时对设备的安全维护也提出了更高的要求，在某一环节出现一个微小的问题都将造成不可挽回的巨额损失。鉴于滚动轴承在机械设备中的地位，可以得出滚动轴承能否正常运行工作关系着机械设备乃至生产线的正常运转，因而掌握滚动轴承运行的工作状态以及故障的形成和发展，是目前机械故障诊断领域中所研究的重要课题之一。根据故障诊断技术机理的不同，诊断滚动轴承故障的常用方法有振动分析法、噪声分析法、油样分析法、温度分析法、油膜电阻法、声发射诊断等。这些分析方法中，振动分析法是运用最广泛的一种方法。振动信号包含丰富的设备异常或故障的信息，通过对振动信号的实时在线监测和分析，可以判断机械系统及其部件是否正常运行，并以此确定机械故障产生的原因，判断故障产生的具体部位和损伤程度。较常用的振动诊断方法有倒频谱分析、特征参数分析法、冲击脉冲法、包络分析法、小波分析等。虽然诊断方法很多，但是对早期的微弱故障诊断仍然是一个难题。因为在故障早期，所需频段信号能量非常微弱、特征信号成分极易被噪声和其他频段信号俺没，未经任何处理的原始频段微弱信号中很难提取出有用的特征信号成分，如何提取被淹没在噪声无用频段内的微弱故障特征信号是目前国内外研究焦点。
技术实现思路
振动信号是典型的非平稳、非线性的时序信号，粒子滤波在解决非线性滤波问题时，随机量不会受到高斯分布的制约，比高斯模型具有更广泛的分布，对非线性特性有很强的建模能力。本专利技术提供了一种利用粒子滤波与谱峭度的滚动轴承故障诊断方法，解决了普通的快速谱峭度方法难以实现在低信噪比的情况下诊断出轴承故障的问题。一种利用粒子滤波与谱峭度的滚动轴承故障诊断方法，具体过程如下：粒子滤波是以状态空间方程为主对信号滤波，所以首先要建立轴承振动信号的状态空间方程。对于一个非线性、非高斯系统，状态空间模型的一般表达式为：其中，vk和nk分别是过程噪声和观测噪声，并且都是相互独立的，协方差分别是Qk和Rk的零均值噪声序列，f和h为状态转移方程和观测方程。众所周知，轴承振动信号是典型的非线性、非高斯时间序列。所以用AR时序模型作为状态方程，模型的残差方差即为过程噪声。提取原始信号的背景噪声作为观测噪声，将其与状态方程加载一起作为观测方程。综上所述，振动信号的状态空间表达式可以写为：式中，p为AR模型阶数，at(t)是模型参数序列，为模型残差方差，即过程噪声，为信号背景噪声，具体获得方法为：用小波分解函数对采集的振动信号进行多尺度分解,得到小波变换系数。再设定合适的阈值,将大于阈值的小波系数置为零,小于阈值的小波系数保留下来。对经过处理的小波系数进行小波重构,就得到所需的噪声信号，即观测噪声。et和nt都是服从N～(0,1)的标准正太分布。模型参数估计方法为Yule-Walker法，模型定界准则是AIC准则。粒子滤波的本质是利用蒙特卡罗积分法求解贝叶斯积分运算，它的关键是核心粒子的权值大小要能不断被修正和调整，最后获得精确的估计值，下文将标准粒子滤波的具体流程。首先用一组随机粒子和这些粒子的权值大小来近似表达后验概率密度函数p(xk|z1:k)：式中，δ(·)是狄拉克函数，z1:k是观测量。如果这些粒子是通过真实概率密度函数p(xk|z1:k)获得，那么式(3)可写为：然而，在很多情况下获得真实概率密度函数粒子是非常困难的，所以就利用重要性概率函数q(xk|z1:k)获取粒子。假设p(xk|z1:k)与一个分析函数π(xk|z1:k)是线性关系，则是(3)中的权值大小可以通过下式计算：假设重要性函数可以分解为：将式(5)和时(6)代入到式(4)中更新粒子权值大小当满足条件:的时候，式(7)可以写为：又因为先验概率密度函数通常选择为：所以式(8)可以写为：粒子滤波结束后结合快速谱峭度法提取特征故障频率。谱峭度的定义式为：式中，为信号的n阶谱矩，其中，E<·>——均值运算，E|·|取模运算。X(f,t)是信号x(t)在频率f处的复包络。快速谱峭度算法采用的是在全频带范围内1/3和二叉树结构设计的带通滤波器分析信号。算法如下：(1)构造一个截止频率fc＝1/8+ε的低通滤波器h(n)，其中fs＝1,ε＞0。低通滤波器的截止频率增大ε是为了消除过度带的影响。并基于h(n)构造准解析低通滤波器h0(n)和准解析高通滤波器h1(n)，其分析频带分别为[0,1/4]、[1/4,1/2]，即：(2)分别以h0(n)和h1(n)低通和高通滤波，对滤波结果做2倍降采样，如图1所示。如此迭代进行，可以获得如图2所示的滤波器树及相应的滤波结果。其中表示第k层第i个滤波器产生的滤波结果，其中i取值为0到2k-1；同时也是信号在中心频带fi＝(i+2-1)2k-1和带宽(Δf)k＝2k-1的复包络，2倍降采样的目的是保证滤波器每一层中的数据长度与原始数据相同。(3)对于每一个滤波器结果，按照式(10)计算峭度：(4)将所有的谱峭度汇总，从而得到“快速谱峭度图”。根据峭度图选择峭度最大值频带区间解调分析，计算功率谱，对功率谱做FFT变换得到频谱图，根据频谱信息做出故障类型的诊断。本专利技术提供了一种利用粒子滤波与谱峭度的滚动轴承故障诊断方法，降低故障信号中的噪声干扰，提高信噪比，实现了滚动轴承早期微弱故障的诊断。附图说明图1为低通/高通滤波器分解图。图2为二叉树滤波器组及相应的滤波结果。图3内环故障时域振动信号。图4外环故障时域振动信号。图5内环振动信号的AIC曲线图。图6外环振动信号的AIC曲线图。图7滤波之后的内环时域信号。图8滤波之后的外环时域信号。图9为内环信号降噪后频谱图。图10为外环信号降噪后频谱图。图11内环信号滤波之前的频谱图。图12外环信号滤波之前的频谱图。具体实施方式以下结合技术方案和附图详细叙述本专利技术的实施例。实施例的数据来自美国西储大学，数据为电机驱动端轴承故障数据,轴承型号为6205-2RSJEMSKF，深沟球轴承。轴承参数如下表所示：表1轴承参数表内圈直径外圈直径厚度滚动体直径Dd节经dm接触角α滚珠数Z25mm52mm15mm7.94mm139.04mm60°9实验电机转速1797rpm/m，故障尺寸为0.1778mm，内外环的故障特征频率计算公式如下：内环外环观察图3和图4外所示的内环、外环时域振动信号，无法看出冲击信号，所以不能判断是否存在故障。要得到故障特征频率只有作进一步的分析处理。图5和图6分别是内环和外环的AIC曲线图。AIC准则认为：AIC值越小，模型精度越高。但是这样模型的阶数就会变大，计算量就会上升。在上升到一定阶段时，模型的精确度度提高非常缓慢。综合考虑两者因素，我们选定当AIC值下降缓慢时候的阶数为模型阶数。内环的AIC曲线在阶数为20时，值已经趋于平缓，下降幅度很小；而外环在第8阶时下降幅度最小。所以模型阶数分别定本文档来自技高网...

【技术保护点】
利用粒子滤波与谱峭度的滚动轴承故障诊断方法，其特征在于如下步骤：(一)滤波降噪，过程如下：(1)初始化；由先验概率p(x0)产生粒子群所有粒子权值为1/N；(2)更新；在k时刻，更新粒子权值，wki=wk-1ip(zk|xki)=wk-1ip(zk-h(xki)),i=1,2,···N,]]>并且归一化权值则得k时刻未知参数x最小均方估计为：(3)重采样；得到新的粒子集合(4)预测；利用状态方程f预测未知参数(5)时刻k＝k+1，转到步骤(2)；(二)滤波降噪之后，采用快速谱峭度的算法计算得到峭度图，过程如下：(6)构造一个截止频率fc＝1/8+ε的低通滤波器h(n)，其中fs＝1,ε＞0；低通滤波器的截止频率增大ε是为了消除过度带的影响；基于h(n)构造准解析低通滤波器h0(n)和准解析高通滤波器h1(n)，其分析频带分别为[0,0.25]、[0.25,0.5]，即：h0(n)=h(n)ejπn/4h1(n)=h(n)ejπn/4]]>(7)分别以h0(n)和h1(n)进行低通和高通滤波，对滤波结果做2倍降采样；如此迭代进行，获得滤波器树及相应的滤波结果；其中表示第k层第i个滤波器产生的滤波结果，其中i取值为0到2k‑1；同时也是信号在中心频带fi＝(i+2‑1)2k‑1和带宽(Δf)k＝2k‑1的复包络，2倍降采样是保证滤波器每一层中的数据长度与原始数据相同；(8)对于每一个滤波器结果，按照式f≠0计算峭度：K=E<|cki(n)|4>{E<|cki(n)|2>}2-2]]>(9)将所有的谱峭度汇总，从而得到“快速谱峭度图”；(10)选择角度最大的频带做解调分析，计算功率谱得到频谱图。...

【技术特征摘要】
1.利用粒子滤波与谱峭度的滚动轴承故障诊断方法，其特征在于如下步骤：(一)滤波降噪，过程如下：(1)初始化；由先验概率p(x0)产生粒子群所有粒子权值为1/N；(2)更新；在k时刻，更新粒子权值，并且归一化权值则得k时刻未知参数x最小均方估计为：(3)重采样；得到新的粒子集合(4)预测；利用状态方程f预测未知参数(5)时刻k＝k+1，转到步骤(2)；(二)滤波降噪之后，采用快速谱峭度的算法计算得到峭度图，过程如下：(6)构造一个截止频率fc＝1/8+ε的低通滤波器h(n)，其中fs＝1,ε＞0；低通滤波器的截止频率增大ε是为了消除过度带的影响；基于h(n)构造准解析低通滤波器h0(n)和准解析...

【专利技术属性】
技术研发人员：李宏坤，任远杰，
申请(专利权)人：大连理工大学，
类型：发明
国别省市：辽宁;21