双自寻优小波去噪方法及其在滚动轴承故障诊断中的应用
2013-01-29肖方煜
汤 伟, 肖方煜, 陈 曦
(陕西科技大学 自动化应用研究所, 陕西 西安 710021)
0 引言
滚动轴承在各类机械中应用很广泛,且其运行状态直接影响整台机器的性能、寿命、功能、功率及安全.在旋转机械中滚动轴承发生故障的概率很高,据统计大约有30% 的机械故障都是由滚动轴承引起的[1].但在实际的故障诊断过程中,所获取的滚动轴承故障信号往往夹杂着各种各样的噪声,给故障的判断造成了很大的影响.所以在诊断故障前,必须去除噪声.小波阈值去噪方法由于原理简练、可操作性强与去噪效果显著的特点,已在信号去噪领域中获得了成功应用.
小波阈值去噪方法中,决定去噪效果最为关键的两个因素是阈值和阈值函数的选取[2].但是如何选取阈值和阈值函数,并且它们之间如何搭配,使其去噪的效果最优成为小波阈值去噪方法最关键的问题[3].本文正是基于这个问题,提出一种阈值和阈值函数各自寻优又相互结合的选取策略予以解决.通过模拟信号消噪仿真与滚动轴承故障信号消噪仿真实例表明,该方法无论在信噪比和去噪图像的直观效果上都优于传统的去噪方法.
1 现有阈值和阈值函数选取方法简介
关于阈值和阈值函数的选取思想分为3类:第一,阈值和阈值函数的选取相互独立.第二,阈值的选取基于阈值函数.第三,阈值函数的选取基于阈值.第一类中阈值和阈值函数相互独立选取策略忽略了阈值和阈值函数之间耦合的关系.因为根据风险函数的定义[4]:
(1)
去噪的目标是使估计信号与真实信号的风险最小化[5].而估计信号是含噪信号通过由阈值和阈值函数映射处理后得到的,因此可以把阈值和阈值函数当做一个耦合的系统.第一类思想的选取方法忽略了两者之间的耦合关系,因此很难使去噪效果达到最优化.第二类思想有基于某种阈值函数的阈值选取策略,代表有基于史坦的无偏似然估计方法、文献[6]提出的基于硬阈值函数的阈值选取策略等.这种策略认同了阈值和阈值函数的耦合关系,会提高去噪的效果.但是,每种阈值函数都有它们的缺点,文献[7]中证明硬阈值函数存在较大的方差,软阈值函数存在较大的偏差,基于第二类思想选取策略势必也存在着缺点,会影响去噪的效果.基于阈值函数存在的这些缺点,很多学者创造出了含有变量参数的阈值函数,典型的有半软阈值函数、指数型阈值函数等.第三类方法是在这些含有参数的阈值函数中建立起来的,通过某一确定阈值选取策略,对阈值函数按照试凑或自寻优的方法得到能使去噪效果达到最优的参数[8].
上述的3种方法在实际应用中都展露出了较为出色的效果,但这3种方法都没能把阈值和阈值函数的选取完全联合在一起.不仅从阈值和阈值函数之间耦合的角度来看存在着片面性,而且从数与数的组合(阈值和阈值函数的参数都是具体的数值)角度来看也存在片面性,因为上述3类方法中皆都运用数理统计决策论的方法估算阈值和阈值函数,这样从数与数的组合角度来看估算到最优的搭配的概率极低,会对产生最优去噪效果埋下伏笔.
2 双自寻优方法
双自寻优方法的目标是所有阈值和所有阈值函数参数都能互相结合,从各种结合中寻找出最优的搭配.这样的好处不仅解决了耦合问题,而且也涵盖了所有数与数之间组合性的所有可能性.
2.1 阈值自寻优区间及高频死去值现象
阈值是一个不小于零的数[9],它的自寻优区间是从0到一个未知数,数学表示为[0,x].在确定自寻优区间这个未知数之前,先介绍在试验中得到的一个现象--高频死去值现象[10].
高频死去值现象的表述:小波阈值去噪方法中,含噪信号经过小波分解后,对任意一层阈值使其从零按某一较小步长递增,当阈值增大到某一数值之后,一切衡量去噪效果的数值(信噪比等)不会再随阈值的增大而改变,定义此阈值为这层的高频死去值.
证明:设vj表示小波分解后信号的低频部分,Wj表示对应的高频部分,{vj}j∈z是L2(R)(L2(R)是全体能量有限信号的集合)的一个OMRA,φ(t)是相应的生成元,相应的双尺度方程为:
(2)
(3)
φjk(t)=2j/2φ(2jt-k)
(4)
ψjk(t)=2j/2ψ(2jt-k)
(5)
而Vj与Wj分别是由{φjk(t),k∈z}和{ψjk(t),k∈z}生成,即
(6)
(7)
设fj(t)和wj(t)分别是f(t)在vj和wj中的逼近与细节,即
(8)
(9)
因为Vj+1=Vj⊕Wj,Vj⊥Wj,所以
cj,k=〈fj,φjk〉=〈fj+wj,φjk〉=〈fj+1,φjk〉
(10)
而
(11)
代入(10)式,得
(12)
同理,得
(13)
经过分解后信号有如下表达式:
(14)
假定对fJ(t)在第j空间上的细节小波分量wj(t) 进行映射处理,也就是对dj,k进行阈值处理.因为dj,k是一个具体的值,所以在第j空间上必定可以找到一个阈值λ,使λ≥dj,k.因此,定义λ=dj,k的λ值为高频死去值.
根据{φj+1,k(t),k∈z}的正交性,得
cj+1,k=〈fj+1,φj+1,k〉=〈fj,φj+1,k〉+〈wj,φj+1,k〉
(15)
因为选择高频死去值λ=dj,k作为阈值,使得dj,k被置零,即〈wj,φj+1,k〉 项为零.因此
cj+1,k=〈fj+1,φj+1,k〉
(16)
证毕.
2.2 高频死去值的获取及阈值自寻优区间的确定
因为不同阈值函数的映射性质不同,相同的含噪信号经过不同阈值函数产生的信噪比数值也就不同.根据高频死去值性质中大于高频死去值后的小波系数不含有信号成分这一特性,可以使用两个不同阈值函数进行自寻优搜索,使它们自寻优得到的信噪比数值相减.当相减的结果在某一值之后为零,则这个值就为高频死去值.图1为基于同一含噪信号两种不同阈值函数消噪的信噪比-阈值示意图.图2为经过matlab的仿真数据图,此图是基于Blocks[11]信号,小波基为sym8,信噪比为4 db,在第一层时模平方阈值函数和软阈值函数产生的信噪比数值矩阵相减结果图.此图表明当阈值小于54时,其差值曲线一直变化,当大于54之后,其差值曲线为零.又因为所选的步长为十分位,所以这个信号在第一层的高频死去值为5.4.
图1 基于同一含噪信号两种不同阈值函数消噪的信噪比-阈值图
图2 基于不同阈值函数的信噪比相减结果图
2.3 基于软硬阈值折衷法阈值函数的自寻优区间的确定
针对含有参数的阈值函数的寻优过程方法比较简单,只需对其参数在某一区间进行自寻优的过程.例如针对软硬阈值折衷法阈值函数[12].
其定义为:
(17)
参数α的自寻优区间是[0,1],假设以0.1为其步长进行寻优过程,因此经过11次去噪算法后将得到11个衡量去噪效果的数值(例如信噪比),在这十一个数值中经过大小的比较,可以得到一个最优的数值.
因为软阈值函数代表收缩力最强的阈值函数,而硬阈值函数代表保留性最强的阈值函数[13].因此,阈值函数中参数选取的区间范围,是使阈值函数能从软阈值函数变化到硬阈值函数为标准的.
2.4 双自寻优方法程序的实现
整个双自寻优的过程方法如下.
首先,根据确定的分解层数,运用2.2节的方法确定最高分解层的高频死去值.由此确定阈值的自寻优范围是0到最高分解层高频死去值.(因为分解层数越高,则对应的高频部分能量越大[14],则高频死去值也越大,所以这里是最高层的高频死去值)
其次,确定含有参数的软硬阈值折衷法阈值函数的自寻优区间.这个范围是[0,1].
最后,编写整个自寻优去噪程序.
程序如下(本文程序为matlab语言程序,且以分解层数为一层,含有参数的阈值函数是软硬阈值折衷法来编写):
for i=0∶1:最高层高频死去值(阈值自寻优过程)
for j=0∶1:已知值(阈值函数自寻优过程)
t=i/步长1;
u=j/步长2;
e=含噪信号;e=e';
[c,l]=wavedec(e,1,'sym8');
a1=appcoef(c,l,'sym8',1);
d1=detcoef(c,l,1);
softd1=wthresh(d1,'p',t);
c2=[a1 softd1];
s3=waverec(c2,l,'sym8');
p1=1/length(e)*norm(e)^2;
p2h=1/length(e)*norm(e-s3)^2;
snrh=10*log(p1/p2h);
rmseh=sqrt(p2h);
plot(s3),axis([0 1200 -3 3]);
snr= rmse=[rmseh]
end
end
其中wthresh调用的程序是
function y = wthresh(x,sorh,t)
switch sorh
case 'p'
tmp=(abs(x)-u*t);
tmp=(tmp+abs(tmp))/2;
y=sign(x)*tmp;
最后,由max(max(snr)就可得到衡量去噪效果最优的信噪比值,并可以对应的找出具体的u和t是多少.从而得到最优的去噪效果.
3 模拟信号去噪
通过构造一个数据长度为1 024且受到强噪声污染的冲击信号,以模拟滚动轴承故障信号,利用Symmlet 8小波对含噪冲击信号进行3尺度静态小波分解.本文的自寻优算法将分别基于软阈值函数、硬阈值函数与经典的Heursure、Sureshrink、Visushrink、Mnimaxi方法进行比较,比较结果在表1中,表1所得数据为模拟信号与所选阈值与阈值函数共同决定的信噪比,单位为db;产生的图形如图3所示.
表1 所得信噪比试验结果
图3 模拟信号去噪结果图形的比较
4 滚动轴承故障信号去噪
试验数据由美国凯斯西储大学电气工程实验室网站下载[15].试验装置由一个1.5 kW 电动机、一个扭矩传感器/译码器、一个功率测试计与相应电器控制装置组成.实验中使用加速传感器采集振动信号,加速度传感器安装在风扇端6点钟方向.振动信号通过16通道的DAT记录器采集的MAT格式,数字信号的采集信号为120 000 S/s.采用6205-2RS JEM SKF 型轴承,轴承的内径直径为25.001 2 mm ,外径直径为51.998 9 mm ,厚度为0.590 6 mm,节径为39.039 8 mm,滚动体直径为15.001 2 mm,其外圈损伤直径为0.177 8 mm,故障特征频率为104 Hz.选择时间长度1 s(对应的是1 200个数据点)的故障数据,如图4所示.
图4 原始信号
图4所显示的信号难以判断其表示故障特征的冲击信号.图5为本文方法去噪后的故障信号,其去噪的时域图形中故障冲击信号特征明显.可以清晰的分辨出由故障引起的冲击信号大约为115个点位周期出现,这与故障特征频率为104 Hz(104 Hz对应的点为115个点)比较一致.特别是在400~600处的点和1 100~1 200处的点信号恢复的情况尤为明显,可以明显辨认出故障冲击信号.图6为经minimaxi与模平方阈值函数去噪后的故障信号.500~600点处的信号恢复较好但600~700和1 100~1 200点处的信号仍然难以辨认.图7为经rigrsure和软阈值函数去噪后的故障信号,600~700和1 100~1 200点处的信号难以辨认.图8为经rigrsure和软硬阈值折衷法函数去噪后的故障信号.其较好的得到了故障特征信号,但在1 100~1 200点处的信号仍然难以辨认.
图5 本文方法
图6 minimaxi与模平方阈值函数方法
图7 rigrsure和软阈值函数方法
图8 rigrsure和软硬阈值折衷法函数
上述三种方法中,第一种minimaxi与模平方阈值函数的结合是本文第一章中提到的第一类方法,阈值与阈值函数相互独立的选取方法.因为阈值与阈值函数相互独立的选取方法忽略了两者之间的耦合关系因此去噪效果不佳.第二种rigrsure和软阈值函数的结合是第二类方法,因为软阈值函数自身存在的偏差大的特点使得去噪效果也不是很好.第三种rigrsure和软硬阈值折衷法函数结合的方法属于第三类方法,因为这种方法只是一个因素在变化,所以去噪效果仍然不是很好.
5 结束语
阈值和阈值函数是小波阈值去噪方法中最重要的两个因素.为改善小波阈值去噪方法的去噪性能,本文从传统方法中归纳其优缺点,并从阈值和阈值函数相互耦合的思想出发,运用双自寻优方法对含噪故障信号进行去噪处理.从模拟的故障冲击信号和来自西储大学的滚动轴承信号实例表明,基于软硬阈值折衷法的双寻优方法无论在信噪比还是得到的恢复图形中皆优于传统小波阈值去噪方法,具有一定的理论研究与工程实际应用价值.
[1] 钟秉林.黄 仁.机械故障诊断学[M].北京:机械工业出版社,2007.
[2] 张德丰.Matlab小波分析与工程应用[M].北京:国防工业出版社,2008:42-49.
[3] D.L.Donoho,I.M.Johnstone,Ideal spatial adaptation by wavelet shrinkage[J].Biometrika,1994,81(12):425-455.
[4] Ronggen Yang.Wavelet denoising using principal component analysis[J].Expert Systems with Applications,2011,11(38):1 073-1 076.
[5] 葛哲学.小波分析理论与MATLAB R2007实现[M].北京:电子工业出版社,2007:67-69.
[6] 赵治栋,潘 敏.小波收缩中统一阈值函数及其偏差、方差与风险分析[J].电子与信息学报,2005,27(4):536-540.
[7] 张 弦.进化小波消噪方法及其在滚动轴承故障诊断中的应用[J].机械工程学报,2010,46(15):76-81.
[8] Hilton de Oliveira Mota.Partial discharge signal denoising with spatially adaptive wavelet thresholding and support vector machines[J].Electric Power Systems Resarch,2011,23(81):644-659.
[9] D.L.Donoho,I.M.Johnstone,Wavelet Shrinkage:Asymptopia[J].Journal of the Royal Statistical Society,1995,46(57):301-369.
[10] 肖方煜,汤 伟.自寻优阈值小波去噪方法[J].信号处理,2012,28(4):577-586.
[11] 祝海龙,郭天佑,屈梁生.基于二进小波变化和软阈值改进的信号消噪[J].自动化学报,2004,3(2):199-206.
[12] L.-K.Shark,Gear fault detection using customized multiwavelet lifting schemes[J]. Mechanical Systems and Signal Processing,2010,24(5):1 509-1 528.
[13] Andrew G.Bruce.Understanding waveshrink:variance and bias estimation[J].Biometrika,1996,83(4):727-745.
[14] 孟志强,王鼎顺,周华安.小波图像消噪中阈值与信噪比的单峰规律及阈值试探方法[J].电子测量与仪器学报,2006,20(5):63-67.
[15] The Case Western Reserve University Bearing Data Center Website[EB/OL].http://www.eecs.cwru.edu/laboratory/bearing/.