小波分析(Wavelet Analysis)是数字信号处理中非常有力的一种工具。它是在 20世纪 80 年代初,由Morlet在分析研究地球物理信号时提出来的。传统的傅里叶分析中,信号完全是在频域展开的,不包含任何时域信息。但在故障诊断时这些丢弃的时域信息同样非常重要,因此许多学者对傅里叶分析进行了推广,提出了很多能表征时域和频域信息的信号分析方法。小波分析便属于时频分析方法的一种。它是一种信号的时间-尺度(时间-频率)分析法,具有多分辨率分析特点,而且在时域和频域都具有表征信号局部特征的能力,是一种窗口大小固定不变,但其形状可变,时间窗和频率窗都可改变的时频局部化分析方法。小波分析的基本思想是用一族函数去表示或逼近某一信号或函数,这一族函数成为小波函数系,它是通过基本小波函数的不同尺度的平移和伸缩构成的 [137] 。
“小波”即小的波形,“小”是指具有衰减性,“波”即具有波动性。小波分析的原始思想形成于 20 世纪初,Haar提出 L 2 ( R )函数空间的一组正交基,后来被认为是最早的小波基。“小波”这一概念是由法国地质学家Morlet在 20 世纪 80 年代研究地下岩石油层分布时正式提出的,并成功应用于地质数据处理中。
小波变换(wavelet transform,WT)是一种新的变换分析方法,它继承和发展了短时傅里叶变换局部化的思想,同时又克服了窗口大小不随频率变化等缺点,能够提供一个随频率改变的“时间-频率”窗口,是进行信号时频分析和处理的理想工具。小波变换在低频部分具有较高的频率分辨率和较低的时间分辨率,在高频部分具有较高的时间分辨率和较低的频率风变绿,很适合于探测正常信号中夹带的瞬态反常现象并展示其成分。它克服了传统傅里叶变换的缺陷,具有良好的时、频局部化性能,使得小波理论在信号去噪、信号处理、图像处理、计算机分类与识别、语言合成、数值分析、医学成像与诊断、故障诊断、地震勘探、分形理论、天体力学等领域得到了广泛应用。它的主要特点是通过变换能够充分突出问题某些方面的特征,能对时间(空间)频率局部化分析,通过伸缩平移运算对信号(函数)逐步进行多尺度细化,最终达到高频处时间细分,低频处频率细分,能自动适应时频信号分析的要求,从而可聚焦到信号的任意细节,解决了傅里叶变换/Fourier变换的困难问题。小波变换的实质是原信号与小波基函数具有相似性,小波系数就是小波基函数与原信号相似的系数。
小波变换主要包括连续小波变换和离散小波变换 [138] 。
设 L 2 ( R )为能量有限的信号空间,即: ,若 ψ ( t )∈ L 2 ( R ),其傅里叶变换满足容许性条件:
即 C ψ 有界,则称 ψ 为基小波或母小波。将 ψ 经伸缩和平移后,可得小波序列:
其中, a , b ∈ R ,且 a ≠ 0。称 a 是伸缩因子, b 是平移因子。
对于任意的函数 f ( t )∈ L 2 ( R ),可以定义其关于基小波 ψ 的连续小波变换式:
其中, 为 ψ 的共轭运算。
原来的一维信号经小波变换转换成二维信号,便于更好地分析其时频特性。小波逆变换把二维信号重构回一维,变换公式如下:
f ( t )的小波变换是一个二元函数,其实质是将 L 2 ( R )空间的任意函数 f ( t )表示为该函数在具有不同平移因子和伸缩因子的 ψ a , b ( t )上的投影叠加。从形式上看,是 f ( t )在 t = b 附近按 ψ a , b ( t )进行加权平均,体现了以 ψ a , b ( t )为标准的 f ( t )的变化快慢。由此,可以说小波变换是一个“变焦镜头”,平移因子 b 是一个时间中心参数,相当于镜头相对目标的平行移动,而伸缩因子 a 相当于调焦旋钮,实现焦距的调节。小波变换将一维的时域函数映射至二维的“时间-尺度”域,使信号在小波基上的展开具有多尺度性,调整平移因子和伸缩因子获得不同时频宽度的小波,匹配原始信号的任意位置,实现对信号的时频局特性分析。
连续小波变换是一种线性变换,它具有以下几个方面的性质 [139] :
已知 x ( t )和 y ( t )∈ L 2 ( R ), k 1 , k 2 为任意常数,且 x ( t )的连续小波变换为 W fx ( a , τ ), y ( t )的连续小波变换为 W fy ( a , τ ),则 z ( t )= k 1 W fx ( a , τ )+ k 2 W fy ( a , τ )的连续小波变换为:
设 x ( t )的连续小波变换为 W fx ( a , τ ), x ( t- t 0 )的连续小波变换为 W fx ( a , τ -t 0 ),即延时后的信号 x ( t- t 0 )的小波系数是将原信号 x ( t )的小波系数在τ轴上进行同样时移。
设 x ( t )的连续小波变换为 W fx ( a , τ ),则 的连续小波变换为 , λ > 0。
若 x 1 ( t )和 x 2 ( t )∈ L 2 ( R ),它们的连续小波变换分别为 W fx 1 ( a , τ )和 W fx 2 ( a , τ ),也即:
此外,连续小波变换还有能量不变性等其他一些性质。
在连续小波变换中,平移因子和伸缩因子均为连续变化的实数,具体应用时需要计算连续积分,不适于实际应用中数字信号的处理。因此,实际应用中常采用离散小波变换,由离散化的平移因子和伸缩因子来完成。一般取:
将其代入式(2-2),可得离散小波:
相应的离散小波变换为:
在实际应用中,为使小波变换的计算更有效,构造的小波函数均具有正交性,即:
从理论上可以证明,将连续小波变换离散成离散小波变换,信号的基本信息不会丢失。相反,由于小波基函数的正交性,使得小波空间中两点之间因冗余度造成的关联得以消除;同时,计算误差更小,变换结果使时频函数更能反映信号本身的特性。
对于 f ( t )∈ L 2 ( R n )( n >1),公式 f ( t )= 存在几种扩展的可能性。一种可能性是选择小波 f ( t )∈ L 2 ( R n )使其为球对称,其傅里叶变换也为球对称,即:
其相容性条件变为:
对所有的 f , g ∈ L 2 ( g n )有:
其中, W f ( a , b )=< ψ a , b >, ψ a , b ( t )= ,( a ∈ R + , a ≠ 0, b ∈ R n )。
公式(2-1)可改写成:
若选择的小波不是球对称的,但可以用旋转进行同样的扩展与平移。在二维时,可定义:
其中, a > 0, b ∈ R 2 , R θ = 。
相容性条件变为:
对应的重构公式为:
消噪问题是信号处理的经典问题,在实际的工程中所采集的信号都存在背景噪声,这对信号的分析是极为不利的,严重的甚至会影响到检测仪器对设备性能的检测,使之做出错误的判断,所以对信号进行消噪处理是信号处理必备的一环。传统的消噪方法多采用平均或线性方法,不适于非平稳信号,去噪效果差,难以识别混杂于信号中的噪声,从而影响故障检测过程及故障诊断结果。小波变换以其良好的时频特性,实现非线性去噪处理,在消噪处理中受到越来越多的青睐。
小波消噪的实质即从含噪信号中准确恢复出原始信号,其本质是函数逼近的问题,即在小波基函数伸缩和平移构成的函数空间中,依据某种准则,实现对原始信号的最佳逼近。小波消噪法的基本思想是:由小波变换将含噪的信号进行多尺度分解,取得一组小波系数,在每个尺度下对小波系数进行阈值量化处理,最大可能地去除噪声的小波系数,而保留并增强信号的小波系数,对处理后的小波系数重构,从而得到消噪后的信号。小波消噪的一般过程如图 2-1 所示。
图 2-1 小波消噪法原理框图
常用的小波消噪方法有:基于小波变换的模极大值消噪、基于相邻尺度小波系数相关性消噪、基于小波变换的阈值消噪 [140] 。
模极大值去噪算法是根据信号和噪声在多尺度空间上小波变换系数的模极值传播规律的不同而发展起来的一种去噪算法。理论上只要信号与噪声的奇异性有差异,就能产生很好的去噪效果。一般信号小波系数的模极大值将随着小波分解层数的增大而增大;而对于白噪声信号,其模值随着分解层数的增大而减小。因此,观察不同尺度间小波变换模极大值变化的规律,去除幅度随尺度的增加而减小的点(对应噪声的极值点),保留幅度随尺度增加而增入的点(对应于有用信号的极值点)。再由保留的模极大值点用交替投影法进行重建,即可以达到去噪的目的。但是,交替投影法算法复杂,容易造成投影信号的偏差,难以在实际应用中对信号进行实时处理。
相关性去噪算法是根据信号经小波变换后,其小波系数在各尺度上有较强的相关性,尤其是在信号的边缘附近,其相关性更加明显,而噪声对应的小波系数在各尺度间却是没有这种明显的相关性来去噪的。在尺度空间上的相关运算能使噪声的幅值大为减小,从而抑制了噪声和小的边缘,增强了信号的主要边缘,更好地刻画了原始信号。并且在小尺度上,这种作用明显大于在大尺度上的作用。由于噪声能量主要是分布在小尺度上,因而这种随尺度增大而作用强度递减的性质,恰好滤除了噪声,很好地保留了有用信号。
依据小波变换的去相关性,变换使得信号的能量集中在幅值较大的少量小波系数中,而噪声的能量分布在幅值较小的大量小波系数中,并且将覆盖整个小波域。由此使得幅值较大的小波系数很大程度上对应于有用信号,幅值较小的系数可认为是噪声,设定一个合适的阈值,把信号的小波系数保留或增强,而消除噪声系数,即可得到消噪后的信号。基于这一思想,Donoho和Johnstone提出硬阈值和软阈值去噪方法 [138] ,即将小波系数中幅值较小的系数置零,保留或收缩幅值较大的系数,得到估计小波系数并重构,从含噪信号中估计真实信号。
设阈值为 λ ,小波变换多尺度分解后的系数为 w ,阈值处理后的小波系数为 w λ 。
当小波系数的幅值小于给定阈值时,令其为 0;系数幅值大于阈值时,保持不变,其阈值函数:
当小波系数的幅值小于给定阈值时,令其为 0;系数幅值大于阈值时,将其减去阈值,其阈值函数:
在小波阈值消噪过程中,最为关键的是阈值函数及阈值的确定,它们从某种程度上将会直接影响消噪的质量。硬阈值和软阈值的阈值处理函数如图 2-2 所示。
图 2-2 两种阈值处理函数示意图
由图 2-2 可知,硬阈值函数在 w = λ 处不连续,容易使去噪信号在奇异点附近出现伪吉布斯(Pseudo-Gibbs)现象。
阈值 λ 的确定是消噪的另一个难题。阈值选择过小,小波变换系数中将会包含较多的噪声成分,而阈值太大,则可能将有用信号一同消除,造成信号失真。常用的选取阈值的规则有以下几种:无偏似然估计、固定阈值估计、极值阈值估计和启发式阈值估计。其中,极值阈值估计法和无偏似然估计法较保守,当噪声在高频区分布得较少时,二者具有良好的去噪效果,能够在较低信噪比的条件下提取信号。启发式阈值估计方法和固定阈值估计方法的去噪结果较有效、彻底,但容易把有用高频信号误判成噪声而舍弃,造成信息丢失。目前,各种阈值公式不断涌现,综合考虑算法的复杂度及其消噪效果,在实际应用中通常采用Johnstone和Donoho给出的统一的阈值计算公式:
其中, N 是信号的长度, σ 是噪声的标准方差。
式(2-9)是在正态高斯噪声模型下,针对多维独立正态变量的联合分布在维数趋于无穷时的情况得出的结论。当系数大于该阈值时,含有噪声的概率近似于零。由于阈值 λ 正比于信号长度 N 的对数的平方根,当 N 的取值很大时, λ 具有将所有小波系数置为零的趋势,此时的小波滤波器相当于一个低通滤波器。在实际应用中,噪声的标准方差通常是未知的,常用下式的估算方法来确定:
其中, Med ()表示求中值运算, d j 为第 j 层分解的小波系数。
硬、软阈值消噪法,计算方法较复杂,消噪效果也不十分令人满意。因此,本文提出一种基于风电轴承振动信号的自适应阈值小波消噪方法。
自适应阈值消噪法是在阈值法基础上的改进。采用软阈值消噪时,总体效果较好,但当含噪信号不规则时仍会失去一部分信号的细节信息;当采用硬阈值消噪时,消噪效果并不理想,仍然含有明显的噪声。这说明当噪声为时变时,传统的消噪方法效果很有限。采用自适应阈值消噪法,可以克服上述缺陷。
依据模式识别问题中的分类方法,将风力发电机组轴承的振动信号和干扰信号看作两类事物,它们具有类内密集、类间分离的特性,以最小错误率为分类准则进行模式分类,这种轴承信号的自适应阈值小波消噪法具体流程如图 2-3 所示,包括多尺度分解、阈值量化处理和小波重构三个重要步骤。
图 2-3 风电轴承振动信号的自适应阈值小波消噪法流程图
假设将某一叠加了噪声的含噪信号 hz ( i )表示为:
其中, s ( i )为实际信号, z ( i )为噪声, i 表示采样序号。
消噪处理的目的即准确地从含噪信号 hz ( i )中恢复实际信号 s ( i )。采用Mallat算法将含噪信号 hz ( i )进行正交小波变换,得到各层高、低频分解系数,如下:
式中, j 为分解层数, N 为采样点数; a j , k , d j , k 为尺度系数和小波系数; h , g 分别为低通和高通滤波器,二者相互正交。运用该算法将使得信号在每次分解时,其长度减半, n 表示每层的信号长度。
传统的硬阈值、软阈值等小波变换阈值处理方法,主要适用于含有高斯白噪声的信号,采用单一的阈值对小波函数进行处理,不能在每一尺度上将信号和噪声有效分离。针对传统软硬阈值方法的缺点,本文采用自适应阈值法,对不同尺度的小波系数选用不同的阈值,具体步骤如下:
①设定初始阈值:
其中, d max 和 d min 为小波系数中的最大值和最小值。
②由计算的阈值将小波系数分为两类:
计算每个小波系数 d j , k 出现的概率 P i ,若两个分类 W 1 和 W 2 出现的概率分别为 P W 1 和 P W 2 ,则两个分类的均值 、 和方差 、 :
③计算类内方差 、类间方差 和总体方差 ,考查“类内密集、类间分离”特性:
④重新计算阈值 λ j , n + 1 以及类间方差和总体方差之比 K n :
重复执行①—④步,若 K n + 1 < K n ,循环结束,此时得到的即为最佳阈值 λ m 。
得到最佳阈值后,另一个关键是阈值函数的确定。对软、硬阈值方法扬长避短,按照式(2-17)进行阈值的量化处理:
其中, d 为小波系数, α 为调整参数, w 为处理后新的小波系数。
信号各尺度分解系数与能量分布相关,由式(2-17)可知,能量分布因素参与了阈值处理过程,能够根据轴承振动信号的能量分布特性进行自适应降噪。
将各层小波系数进行阈值处理后,按照式(2-18)进行重构:
其中, j 为分解层数, N 为采样点数, n 为每层分解的信号长度; a j , k 尺度系数, w j , k 为处理后的新小波系数; h , g 为式(2-12)中的正交滤波器。
至此,完成了信号多尺度分解、阈值量化处理和小波重构整个去噪过程。总结其思路为:采用Mallat算法对信号 hz ( i )进行 N 层小波分解,得到各层的尺度系数(即低频系数)和小波系数(即高频系数);对 1~ N 层,每层选取一个最佳阈值 λ m ,按照式(2-17)对该层的小波系数进行阈值量化处理。最后,由小波分解的第 N 层尺度系数和经过阈值处理的 1 至 N 层的新小波系数,按照式(2-18)重构,从而得到去噪信号。