谐波小波与近似熵相结合的噪声分析Word格式.docx
- 文档编号:20993398
- 上传时间:2023-01-26
- 格式:DOCX
- 页数:50
- 大小:816.03KB
谐波小波与近似熵相结合的噪声分析Word格式.docx
《谐波小波与近似熵相结合的噪声分析Word格式.docx》由会员分享,可在线阅读,更多相关《谐波小波与近似熵相结合的噪声分析Word格式.docx(50页珍藏版)》请在冰豆网上搜索。
谐波小波与近似熵相结合的噪声信号分析
孙伟杰0221057
0引言
近似熵(Approximateentropy,简称ApEn)是最近发展起来的一种度量序列的复杂性和统计量化的规则。
它是在20世纪90年代初由Pincus为了克服混沌现象中求解熵的困难提出的。
近似熵是对非线性时间序列复杂度的一种非负的定量描述,它对于相对较短的(大于100个数据点)、含噪声的时间序列显示出潜在的应用价值,这是因为产生近似熵的主要的技术思想是:
它并不是企图完全重构吸引子(吸引子是一个数学概念,用于描写运动的收敛类型),而是用一种有效的统计方式——边缘概率的分布来区分各种过程(边缘概率在数学概念中是指当实验所获取的事例按不同的标准进行分类时,忽略掉某些分类标准而只考虑在某一种分类标准下某事件出现的概率)。
在应用的过程中,近似熵表现出以下主要的特点[1]:
(1)只需要比较短的数据就能估计出比较稳定的统计值。
所需的数据点大致在100~5000点,一般在1000点左右。
(2)有较好的抗干扰和抗噪的能力。
在实际应用中,常把它作为一个诊断的判据,已经在生物系统,生理电信号、机械设备故障诊断等领域进行了尝试并获得了良好的效果。
(3)对于随机信号或是确定性信号都可以使用,也可以应用于由随机成分和确定性成分混合的信号。
若一个非线性的物理过程复杂程度越高,那么近似熵将越大。
(4)近似熵更主要的是从衡量时间序列复杂性的角度来度量信号中产生新模式的概率大小,产生新模式的概率越大,序列的复杂性越大,相应的近似熵也越大。
可用近似熵来描述振动信号的不规则性和复杂性,通过比较一定条件下振动信号在不同噪声干扰下的近似熵的相对变化,可以直接反映该振动信号在此期间的运行状况。
本论文着重介绍了近似熵的一般算法和快速算法,并应用实例来说明近似熵在检测振动信号的复杂性方面具有很强的能力。
最后通过编程提取两组不同状态下包含噪声的振动信号的近似熵值,并进行分析比较,以说明近似熵在分析振动信号是否包含振动噪声方面效果良好。
谐波小波是一种复小波,在频域紧支,有明确的函数表达式,其伸缩与平移构成了L2(R)空间的规范正交基。
谐波小波分解算法是通过信号的快速傅立叶变换(FFT)及其逆变换(IFFT)实现的,算法速度快,精度高,因而具有很好的工程应用价值。
对于非平稳振动信号的分析,利用传统的信号分析方法——Fourier变换难以完成。
Fourier变换得到的功率谱密度仅能给出振动信号的平均统计结果,且Fourier变换为纯频域的,它并不能描述信号局部的时频特征。
对于非平稳振动信号中奇异点(突变点)的确定,时变、变频信号的动态分析,单从时域或频域难以识别的非周期信号等,皆需要对信号进行精细的时频分解。
小波变换的出现,使得对于上述非周期信号的时频细节的分解成为可能。
在工程应用中更关心离散小波变换,Mallat[7]和Daubechies[8]给出了计算二进离散小波变换的多分辨分析和塔形算法,为小波变换的广泛应用提供了条件。
但此算法分解得到的不同尺度上信号频率中心为二进的,时域带宽也为二进的。
(所谓二进小波及其离散变换,是指在实际运用中,特别是在计算机实现上,将连续小波及其变换进行二进制离散化的小波和相应的小波变换。
)而Newland[9]提出的谐波小波(harmonicwavelet)除继承了通常意义下小波函数的优点外,另外其还存在以下优点:
小波函数具有明确的函数表达式,无需通过繁冗的尺度函数迭代;
谐波小波变换的时频分解更加灵活,没有上面提到的二进限制;
算法实现简单。
小波分解是基于小波函数的阶段性滤波特性,而谐波小波函数具有频域盒形紧支谱特性及良好的相位定位能力,因此,国内外部分科研人员用谐波小波变换用于振动信号的分析。
正是基于以上考虑,本文将用以谐波小波为主的复小波变换来分析振动信号。
1近似熵的定义及其算法
1.1近似熵的定义
近似熵是用一个非负数来表示某时间序列的复杂性,越复杂的时间序列对应的近似熵越大,换句话而言,近似熵是从衡量时间序列复杂性的角度来度量信号中产生新模式的概率大小,产生新模式的概率越大,序列的复杂性越大,相应的近似熵也越大。
即说明系统越趋近于随机状态,包含频率成分越丰富、系统越复杂、而近似熵越低则信号越趋于周期性、信号包含的频谱越窄。
1.2近似熵的算法与实用快速算法
以下为近似熵具体的算法:
计算近似熵时,需输入两个参数m、r(其中m称为模式维数,r称为相似容限)且这两个参数在整个计算过程中固定不变。
m可以认为是比较序列的长度,即窗口长度,r可以认为是一个有效的阈值。
给出N个点u
(1),u
(2),…u(N),对固定的m和r定义两个参数,一个是极限值ApEn(m,r),一个是这N个点的统计估计值ApEn(m,r,N)。
下面结合matlab程序来说明近似熵的算法步骤。
(1)设原始数据为u
(1),u
(2),…,u(N),共N个点。
(2)将序列{u(i)}按顺序组成m维矢量X(i),即按序号连续顺序组成一组m维矢量:
从X
(1)到X(N-m+1),其中:
X(i)=[u(i),u(i+1),…,u(i+m-1)]。
i=1~N-m+1。
①
这些矢量代表着从第i个点开始的连续m个u值。
对每一个i值计算矢量X(i)与其余矢量X(j)之间的距离:
d[X(i),x(j)]=maxβu(i+k)-u(j+k)β②
其在matlab中通过一个循环用aux1=repmat(X(j,:
),N-m+1,1)表达。
(3)定义矢量X(i)和X(j)间的距离d[X(i),X(j)]为两者对应元素中差值最大的一个,即:
d[x(i),X(j)]=max[|u(i+k)-u(j+k)|] k=0~m-1
(此时X(i),X(j)中其他对应元素间差值自然都小于d)。
并对每一个i值计算X(i)与其余矢量X(j)(j=1~N-m+1,但j≠i)间的距离为
d[X(i),X(j)]。
其在matlab中通过max(abs(dif_aux(k,:
))表达。
(4)给定阈值r,对每个i≤N-m+1的值,统计d[X(i),X(j)]小于r的数目及此数目与距离总数N-m的比值,记作Cmi(r)。
即:
Cmi(r)=1/(N-m){d[X(i),X(j)]<
r}的数目i=1~N-m+1。
)))<
r_factor*std(signal)判断表达。
(5)先将Cmi(r)取对数,再求其对所有i的平均值,记作φm(r)。
Фm(r)=1/(N-m+1)*lnCmi(r)i=1~N-m+1
(6)再把维数m加1,变成m+1重复步骤2~5得Cm+1i(r)和φm+1(r)。
(7)理论上此序列的近似熵为:
ApEn(m,r)=lim[φm(r)-φm+1(r)]N→∞
一般言之,此极限值以概率1存在。
实际工作时N不可能为∞。
当N为有限值时按上述步骤得出的是序列长度为N时ApEn的估计值。
记作:
ApEn(m,r,N)=φm(r)-φm+1(r)
ApEn的值显然与m,r的取值有关。
Pincus根据实践,建议取m=2,r=0.1~0.2STD[STD是原始数据u(i),i=1~N的标准差(Standarddeviation)]。
其在matlab中用[ApEn_value]=ApEn(signal,m,r),r=r_factor*std(signal)表示。
在实际应用中,当需要用近似熵提取振动信号的特征时,其计算算过程对模式维数m,相似容限r都选用m=2,r=0.25*STD(u)作为标准。
为了进一步说明该近似熵算法在度量对于机械信号复杂性方面的能力,下文将引用文献[1]中的一个例子:
某一电厂50MW汽轮发电机组在大修后的重新开机过程中出现了低负荷下振动峰峰值突然增大并导致跳闸停机的现象作为例子,图1.1中XG1时域波形为机组在故障发生时6MW负荷下低压缸4M瓦的振动信号,XG2的小波包分解是小波包分解2层后4个独立频带内的时域波形图,其中XG12,i的频带分别为0~250Hz、250~500Hz、500~750Hz和750~1000Hz,时间历程均为0~0.128s.
图1.1 故障发生时6MW负荷下的振动信号XG1时域波形及其小波包分解
可以看出原始信号及其小波包分解后各频带内的波形杂乱,上下明显不对称,信号的能量主要分布在0~250Hz、250~500Hz的频带内.结合设备运行工况,诊断认为是轴瓦紧力不足和支撑不善造成的松动故障。
按照设计标准,预紧位移应达到0.25mm,然而停机检查却发现预紧位移仅为0.11mm,且该瓦的右后垫铁间隙为0.4mm×
25mm,左后垫铁间隙为0.05mm×
30mm,确实存在着支撑不善和紧力不足的现象.在增加轴瓦紧力和改善左右垫铁支撑后重新开机,机组整体振动有了明显好转。
图1.2、图1.3分别为机组修理后6MW负荷和45MW负荷下的振动信号时域波形和小波包分解后各频带时域波形,频带范围均同上。
可以看出,维修后该瓦的振动峰峰值已经降至43Lm,且所有频带内振动幅值较故障时显著下降,非平稳杂波分量明显减少,尤其是0~250Hz和250~500Hz的频带内的变化最为突出,这表明由松动引起轴瓦振动的摩擦、撞击故障等非平稳因素减少,运行工况得到改善,而且这两个频带为松动故障的敏感频带。
当负荷增大到45MW时,振动幅值有不同程度增大,但0~250Hz和250~500Hz频带内的振动幅值增长幅度较小,表明增加紧力和改善垫铁支撑后较明显地控制这两个频带的振动
图1.2维修后6MW负荷下的振动信号XG2时域波形及其小波包分解
Fig1.2VibrationsignalXG2(6MW)anditswaveletpacketdecompositionafterrepairing
图1.3维修后增大负荷到45MW的振动信号XG3时域波形及小波包分解
Fig1.3VibrationsignalXG3(45MW)anditswaveletpacketdecompositionafterrepairing
观察机组维修前后该瓦的原始振动信号及其小波包分解后各频带内的时域波形可以发现,其复杂性有不同程度的变化,因此通过比较各频带内时域波形的近似熵在不同运行工况下的变化可以很好反映出松动故障的特征.表1.1修理前后及增大负荷时该瓦振动信号及其小波包分解4个频带的信号的近似熵.其中ApEn为原始信号近似熵,ApEn2,1为0~250Hz频带信号近似熵,ApEn2,2为250~500Hz频带信号近似熵,ApEn2,3为500~750Hz频带信号近似熵,ApEn2,4为750~1000Hz频带信号近似熵。
表1.1维修前后及增加负荷时振动信号及其小波包分解的近似熵
ApEn
ApEn2,1
ApEn2,2
ApEn2,3
ApEn2,4
6MW故障信号XG1
1.1519
0.3869
0.9839
1.1318
1.1447
修理后6MW信号XG2
0.8142
0.1765
0.8551
1.0732
1.1379
增大负荷45MW信号XG3
0.9921
0.3791
0.9214
1.1218
1.1725
从表1.1可以看出,松动故障发生时振动信号XG1的近似熵(1.1519)明显大于维修后的振动信号XG2的近似熵(0.8142),这表明松动故障导致振动信号不规则化或复杂化。
维修后各频带信号的近似熵都有所降低,尤其是0~250Hz和250~500Hz频带内近似熵变化最大,这说明基于小波包分解后的各频带时域波形的近似熵对轴瓦松动故障非常敏感。
增大负荷后通频及各频带的信号近似熵比维修后6MW相应的近似熵都有所增大,这说明振动信号中又出现非平稳的成分。
这主要是由于没有对轴瓦松动故障进行彻底维修造成的。
当时虽然用小波分析法准确地找到故障原因,但开机期限非常紧迫,来不及吊出转子,连夜抢修时只增加了轴瓦预紧力和改善左右垫铁支撑,没有调整下瓦垫铁间隙,轴瓦松动故障没有得到彻底根除,机组“带病”运行。
由于高频不稳定性是松动故障的典型特征,所以在负荷增大时各频带内的近似熵都有所增大。
因此,从这个实例中,通过比较机械设备不同运行期间频带内的时域波形的近似熵的变化,可以有效地监测故障的发生和发展;
通过对机械设备同一运行期间不同频带近似熵的比较,可以确定故障的特征频带。
高频不稳定性是松动故障的典型特征。
通过近似熵在敏感频带内的变化可以有效提取出松动故障特征。
由上述实例可以看出,近似熵的计算实际上是在确定一个时间序列在模式上的自相似程度有多大,从另外一个角度讲,就是在衡量当维数变化时该时间序列中产生新模式的概率的大小,产生新模式的概率越大,序列就越复杂。
因此从理论上讲,近似熵能够表征信号的不规则性(复杂性),越复杂的信号近似熵应该越大。
近似熵只是希望从统计的角度来区别时间过程的复杂性,而不企图描述或重建奇异吸引子的全貌,因此只用较短的数据就可以估计出合理的近似熵。
近似熵大致相当于维数变化时新模式出现的对数条件概率的均值,在衡量时间序列的复杂性方面具有一般意义,而不仅仅是一个非线性动力学参数,因此近似熵的估计对随机过程和确定性过程都适用。
同时,当噪声的幅度低于相似容限r时,该噪声将被抑制,若时间序列中存在较大的瞬时状态的干扰时,干扰产生的数据(即所谓的‘野点’)与相邻数据组成的矢量与X(i)的距离必定很大,因而在阈值检测中将被去除,因此,近似熵具有很好的抗噪、抗干扰能力。
以下为近似熵的实用快速算法:
近似熵的计算可以按照上面定义的步骤去进行,然而其中有很多的冗余计算,降低了计算效率,不利于实时运用。
Pincus根据实践,建议取m=2,r=0.1~0.2STD[STD是原始数据u(i),i=1~N的标准差(Standarddeviation)],并给出了一种实用快速算法,可将计算速度提高到定义算法的5倍左右,其算法如下所示:
第一步:
对N点序列,先计算N×
N的距离矩阵D,D的第i行第j列元素记为dij,其定义为节点i,j之间的距离。
dij=βu(i)-u(j)β<
r
或者
dij=βu(i)-u(j)β≧r
i=1~N,j=1~N,i不等于j
其在matlab中是通过D=abs(signal*ones(1,N)-ones(N,1)*signal'
);
与S=zeros(N,N);
S(find(D<
=r_factor*std(signal)))=1;
来实现的,即计算两个矩阵。
第二步:
利用矩阵D中的元素,可以方便地计算得到Ci2(r)和Ci3(r)(假设m=2)。
Ci2(r)=∑dij∩d(i-1)(j-1)j=1...N-2
Ci2(r)=∑dij∩d(i-1)(j-1)∩d(i-1)(j-2)j=1...N-2
其在matlab中是通过Cmr_ij=S(k,1:
N-1).*S(k+1,2:
N);
%(*)
与Nm_i=sum(Cmr_ij);
Cmr_i=Nm_i/(N-(m-1));
来实现的,即计算统计d[X(i),X(j)]小于r的数目及此数目与距离总数N-m的比值,即Cmi(r)
第三步:
由Ci2(r)和Ci3(r)分别计算φm(r)和φm+1(r)。
在matlab中可以通过以下四句话Cmr=[Cmr;
Cmr_i]
和Cmr_ij(end)=[];
Cmr_ij_1=Cmr_ij.*S(k+2,3:
N)
和Nm_i_1=sum(Cmr_ij_1)和Cmr_i_1=Nm_i_1/(N-m)
以及Cmr_1=[Cmr_1;
Cmr_i_1]来实现计算φm(r)和φm+1(r)
第四步:
ApEn(m,r)=φm(r)-φm+1(r)
在matlab中通过phi_m=mean(log(Cmr));
phi_m_1=mean(log(Cmr_1))来实现此功能。
对照一下近似熵算法的步骤和近似熵实用快速算法的步骤,可以看出近似熵实用快速算法与近似熵算法最主要的不同点在于近似熵快速算法将近似熵定义算法中的第一步骤构造矢量的过程给省略了,同时近似熵快速算法不再分别计算模式维数m=2或模式维数m=3时各矢量之间的距离而用以求解时间序列中各数据点的差值替代,即避免了相同维数的矢量之间距离的重复计算,也减少了维数变化时的计算距离过程中的不必要计算,从而提高了运算效果,对于工程而言具有极强的实用价值。
下面还是以实例来说明快速近似熵在度量信号复杂性方面的能力,在图1.2中所示为周期信号:
x=sin(2*pi*50*t)+sin(2*pi*300*t)
其中,取时间间隔t为0.001即说明采样频率为1000Hz;
该信号产生的是主要频率为50Hz和300Hz的信号。
又图1.2中为信号中加入白噪声r后的波形,直观上可以看出信号x+r比信号x要复杂的多,按上述快速算法求其相应的近似熵分别为0.8511和0.2079,前者几乎是后者的3倍多,即越复杂的信号近似熵越大,从而表明近似熵可以很好地用来显示信号的复杂性。
[程序见附录3]
图1.2周期信号的波形与周期信号叠加白噪声的波形
2谐波小波与复morlet小波的定义及性质
2.1谐波小波的定义
小波是满足允许条件的函数,如果一个小波具有完全“盒形”的频谱将是非常理想的。
从这一考虑出发,设有实偶数Wo(t),它们的傅利叶变换
- 配套讲稿:
如PPT文件的首页显示word图标,表示该PPT已包含配套word讲稿。双击word图标可打开word文档。
- 特殊限制:
部分文档作品中含有的国旗、国徽等图片,仅作为作品整体效果示例展示,禁止商用。设计者仅对作品中独创性部分享有著作权。
- 关 键 词:
- 谐波 近似 相结合 噪声 分析
![提示](https://static.bdocx.com/images/bang_tan.gif)