第六章基于小波变换的故障诊断方法ppt课件.ppt
第六章 基于小波变换的故障诊断方法,小波变换的基本原理 奇异性的检测 基于小波变换的原油管道泄漏检测,一、小波变换的基本原理,小波变换是由法国理论物理学家Grossmann与法国数学家Morlet共同提出的。小波分析是近20多年来发展起来的新兴学科,其基础是平移和伸缩下的不变性,这使得能将一个信号分解成对空间和尺度的独立贡献,同时又不丢失原有信号的信息。,小波的由来,小波变换是一种能够在时间频率两域对信号进行分析的方法,具有可以对信号在不同范围、不同的时间区域内进行分析,对噪声不敏感,能够分析到信号的任意细节等优点,在信号处理领域获得越来越广泛的应用,被誉为“数学显微镜”。,小波分析和Fourier分析,傅立叶变换是一个十分重要的工具,无论是在一般的科学研究中,还是在工程技术的应用中,它都发挥着基本工具的作用。,从历史发展的角度来看,自从法国科学家J.Fourier在1807年为了得到热传导方程简便解法而首次提出著名的傅立叶分析技术以来,傅立叶变换首先在电气工程领域得到成功应用,之后,傅立叶变换迅速得到越来越广泛的应用,而且理论上也得到了深入研究。,傅立叶变换最重要的意义是它引进了频率的概念,他把一个函数展开成各种频率的谐波的线性叠加,由此引出了一系列频谱分析的理论。很多在时域中看不清的问题,在频域中却能一目了然 。,因此,长期以来,Fourier分析理论不论在数学中还是工程科学中一直占领着极其重要的地位。,傅立叶分析的实质在于将一个任意的函数f(t)表示为具有不同频率的谐波函数的线性叠加。即一族标准函数 的加权求和,从而将对原来函数的研究转化为对这个叠加的权系数的研究:,其中,权函数:,就是原来函数f(t)的傅里叶变换。,经过以上的变换,就将对,的研究,转化为对权系数,即其傅氏变换,的研究。,从以上分析可知,经典的傅氏分析是一种纯频域分析。,上式中,各符号的含义: 表示频域函数; 表示对原函数f(t)的傅里叶变换; 表示对频域函数 的傅里叶反变换。,傅里叶变换是时域到频域互相转化的工具,从物理意义上讲,傅里叶变换的实质是把f(t)这个波形分解成许多不同频率的正弦波的叠加和。,从傅里叶变换中可以看出,这些标准基是由正弦波及其高次谐波组成的,因此它在频域内是局部化的。,例:假设一信号的主要频率成分是100Hz和400Hz,如下图所示,通过傅里叶变换对其频率成分进行频域分析。,上图为原始信号,从图中看不出100Hz和400Hz的任何频域信息。但从下图的信号频谱分析中,可以明显看出信号的频率特性。,从上例中可知,虽然傅里叶变换能够将信号的时域特征和频域特征联系起来,能分别从信号的时域和频域进行观察,但却不能把两者有机地结合起来。,信号的时域波形中不包含任何频域信息;而其傅里叶谱是信号的统计特性,从其表达式中也可以看出,它是整个时间域内的积分,没有局部化分析信号的功能,完全不具备时域信息。,也就是说,对于傅里叶谱中的某一频率,不知道这个频率是在什么时侯产生的。这样,在信号分析中就面临一对最基本的矛盾:时域和频域的局部化矛盾。,在实际的信号处理过程中,尤其是对非平稳信号的处理中,信号在任一时刻附近的频域特征都很重要。,如在故障诊断中,故障点(机械故障、控制系统故障、电力系统故障等)一般都对应于测试信号的突变点。,对于这些时变信号进行分析,通常需要提取某一时间段(或瞬间)的频率信息或某一频率段所对应的时间信息。因此,需要寻求一种具有一定的时间和频率分辨率的基函数来分析时变信号。,为了研究信号的局部特征,科学家们提出了一些对傅里叶变换进行改进的算法,其中短时傅里叶变换(Short Time Fourier TransformSTFT)就是比较有代表性的一种。,短时傅里叶变换是一种折衷的信号时、频信息分析方法,它是Dennis Gabor于1946年提出的。,短时傅里叶变换的基本思想是:通过给信号加一个小窗,将信号划分为许多小的时间间隔,用傅里叶变换来对每一个时间间隔内的信号进行分析,以便确定该时间间隔内的频率信息。,它假定非平稳信号在分析窗函数g(t)的这个短时间间隔内是平稳的(伪平稳),并移动分析窗函数,使f(t)g(t- )在不同的有限时间宽度内是平稳信号,从而计算出各个不同时刻的功率谱。,短时傅里叶变换定义如下:,其中,f(t)是待分析的信号; 函数 是 的复共轭函数; g(t)是固定的紧支集函数,称为窗口函数。,随着时间的变化,g(t)所确定的“时间窗”在t轴上移动,使f(t)“逐渐”进行分析。,短时傅里叶变换 大致反映了f(t)在时刻时,频率为的“信号成分”的相对含量。,这样,信号在窗函数上的展开就可以表示为在,这一区域内的状态,,并把这一区域称为窗口, 和分别称为窗口的时宽和频宽,表示了时频分析中的分辨率,窗宽越小则分辨率越高。,为了得到更好的时频分析效果,希望和都非常小,但是由海森堡测不准定理(Heisenberg Uncertainty Principle)可知, 和是互相制约的,两者不可能同时都任意小。,(事实上, 0.5,且仅当g(t)为高斯函数时,等号成立。),由此可见,短时傅里叶变换虽然在一定程度上克服了标准傅里叶变换不具有局部分析能力的缺陷,但它也存在着自身不可克服的缺陷,即当窗函数g(t)确定后,矩形窗口的形状就确定了, 和只能改变窗口在相平面上的位置,而不能改变窗口的形状。,可以说,短时傅里叶变换是具有单一分辨率的分析,这对分析信号来说是很不利的。因为,一般来说高频信号持续的时间比较短,低频信号持续的时间比较长。为了更好地分析信号,信号的高频成分需要窄的时间窗,而信号的低频成分需要宽的时间窗。而单一分辨率无法满足这种要求。,正是由于傅立叶分析理论存在上述缺陷,人们一直在寻找更好的基来展开和描绘任意函数,经过多年的探索和总结,逐渐发展成为小波分析理论。,小波变换继承和发展了短时傅里叶变换的局部化思想,并且克服了其窗口大小和形状固定不变的缺点。它不但可以同时从时域和频域观测信号的局部特征,而且时间分辨率和频率分辨率都是可以变化的,是一种比较理想的信号处理方法。,1984年,法国地球物理学家Morlet在分析地震波的局部性质时,发现传统的Fourier变换难以达到要求,因而引入小波概念用于对信号进行分解。,小波变换理论发展过程中的重要阶段,1985年,Meyer构造了具有一定衰减性质的光滑函数,它的二进制伸缩与平移构成了L2(R)的规范正交基,这一发展标志着小波热的开始。,1986年,Lemarie和Battle分别提出了具有指数衰减的小波函数。,1987年,法国马赛召开第一次有关小波的国际会议。,1990年,崔锦泰和王建忠构造了基于样条函数的单正交小波函数。,1988年,Mallat与Meyer合作提出了多分辨分析的框架。,1988年,Daubechies构造了具有有限支集的正交小波基。在美国PureAppl.Math.发表一篇长达87页的论文,被公认是小波分析的经典文献。,1989年,Mallat在多分辨率分析基础上,构造了Mallat算法。为此,Mallat于1989年荣获IEEE论文奖。,1990年,Meyer等出版第一部小波系统性专著小波与算子,共三卷。尤众、王耀东、邓东皋等译校成中文本(共两册)。这套书详细研究了各种小波基的构造,小波基与函数空间的关系,小波分析在复分析、算子论、偏微分方程与分线性分析等方面的应用。,1991年,邓东皋等在数学进展上发表“小波分析”国内第一篇小波论文。对国内小波的研究和应用起了很大的推动作用。,1992年,Daubechies的小波10讲系统论述了正交小波的紧支性、正则性、对称性及时频特性,介绍了离散小波变换和连续小波变换等。,到此,经典小波理论已基本成熟,1992年以后,在国际上,重点转向小波的推广和应用。,在国内,由于对小波的研究起步较晚,20世纪90年代以来,可以说小波的理论研究和应用研究几乎同时开始。1994年,形成国内的小波高潮。,近十年来,小波理论一直在各个不同研究领域扮演着重要的角色。主要集中在数学物理(如分形、混沌、求解方程等)、图像与数据压缩、信号处理、神经网络、故障诊断与检测、石油地质勘探等方面。,定义1:,称满足 的函数f(x)为平方可积函数,并把这类函数的集合记为L2(R)。其中,R表示实数集合。,若f(x),g(x) L2(R),为常数,则f(x)g(x) L2(R)。因此,L2(R)构成了一个线性空间。我们称其为平方可积函数空间。,预备知识,定义2:,在L2(R)空间中的内积定义为:,其中, 表示g(x)的共扼。,定义3:,在L2(R)空间,函数f(x)的范数f(x)定义为:,定义4:,在L2(R)空间,若:内积0,则称函数f与函数g正交。,定义5:,在L2(R)空间,两个函数f(x)与g(x)的卷积定义为:,定义6:,函数f(x)的傅里叶变换 定义为:,定义7:,对任意函数f(x),其扩张函数fs(x)定义为:,其中,s为尺度因子(scale factor),或简称为尺度。,定义8:,把希尔伯特空间(Hilbert space)中的可测的、平方可积的两维函数构成的子空间记作:L2(R2)。,函数f(x,y) L2(R2)的经典范数f(x,y)定义为:,定义9:,f(x,y) L2(R2)的傅里叶变换f(x,y)定义为:,定义10:,定义11:,设f(t)为在R上定义的函数,我们称集合为函数f(t)的支集(即f(t) 0的点所构成的集合的闭包)。,具有紧支集的函数就是在有限区间外恒等于零的函数。,小波与小波变换,我们称满足条件,定义12:,的平方可积函数(x)(即(x) L2(R))为基本小波,或小波母函数。,函数f(x) L2(R)的连续小波变换定义为:,定义13:,其中,*表示卷积。,因此,Wf(s,x)关于x的傅里叶变换可以表示为:,由定义13可知,小波变换Wf(s,x)是尺度s与空间位置x的函数。小波变换通过(x)在尺度上的伸缩和空间域(时域)上的平移来分析信号。,尺度s增大时,s在空间域(时域)上伸展,小波变换的空间域分辨率降低; s()在频域上收缩,其中心频率降低,变换的频域分辨率升高。反之,尺度s减小时, s在空间域(时域)上收缩,小波变换的空间域分辨率升高; s()在频域上伸展,其中心频率升高,变换的频域分辨率降低。,连续小波变换的定义,也即:当检测低频信号时(即对于大的s0),时间窗会自动变宽,以便在低频域用低频对信号进行轮廓分析。反之,当检测高频信息时,(即对于小的s0),时间窗会自动变窄,以便在频率域用较高的频率对信号进行细节分析。因而,小波分析具有“数学显微镜”的美誉。,图 小波变换的时频窗口,例:,图 联合时频分析,小波变换可以对信号做联合时-频域分析得到其特征。,最下面的图是信号在时域的波形,右上图为该信号的频谱,左上的大图为联合时频分析一种算法的结果,前后两个400Hz的频率成分通过联合时频分析可以清楚地看到,而传统傅立叶变换则只能分辨出含有400Hz的信号,不能从时域上分辨出包括两个400Hz频率信号。,通常使用的小波母函数有:Daubechies小波、Harr小波以及Morlet小波。,Morlet小波函数由下式描述:,小波变换具有多分辨即多尺度特点,可以由粗及精的观察信号。,可以将小波变换看成基本频率特性为 的带通滤波器在不同尺度a下对信号作滤波。,小波变换带通滤波器的带宽 与中心频率f成正比即 ,亦即滤波器有一个恒定的相对带宽,即品质因数恒定(称之为等Q结构,Q为滤波器的品质因数)。,适当的选择基本小波,使 在时域上为有限支撑, 在频域上也比较集中,便可以使小波变换在时、频两域都具有表征信号局部特征的能力,因此非常适合于检测信号的瞬态或奇异点。,二、奇异性的检测,通常用李普西兹指数(Lipschitz)来描述函数的局部奇异性。,设n是一非负整数, ,我们说f(x)在点x0为李普西兹,如果存在两个常数A和h00,及n次多项式Pn(h),使得对任意的 ,均有:,如果上式对所有 均成立,且 ,,称f(x)在(a,b)上是一致李普西兹。,信号的奇异性表征与小波变换的模极大值,李普西兹指数越大,函数越光滑。,函数在一点连续、可微,则在该点的李普西兹指数为1;,在一点可导,而导数有界但不连续时,李普西兹指数仍为1。,如果f(x)在x0李普西兹指数小于1,则称函数f(x)在x0点是奇异的。,阶跃信号的李普西兹指数为0,脉冲函数的李普西兹指数为-1。,设实函数 满足 且 ,,如果我们选择小波函数为它的一阶导数,即,,同时记 ,这时,小波变换:,即小波变换 可表示为信号f(x)在尺度s被,平滑后的一阶导数。,图 信号突变点与其小波变换模极大值的关系,例:,x0,x2是信号f(x)的突变点;x1是f(x)慢变区间的转折点。,x0,x2是 的快变化点;x1对应 的慢变化点。,这两种拐点可以通过观察 的极值点是极大点还是极小点分辨出来。x0,x2对应 的极大点;x1对应 的极小点。,函数的奇异点可以从其小波变换的模极大值检测出来。小波变换的模极大值都是出现在信号有突变的地方。信号突变越大,其小波变换的模极大值就越大。,对于边沿检测或奇异点检测来说,我们只是对的极大点感兴趣。,事实上, 的局部极大值通常刻画了信号非正规性的Lipschitz指数。,结论:,奇异点检测的小波的选择,下面以阶跃式边沿和 函数式尖峰这两类突变为例,介绍小波变换的过零点和极值点来检测信号的局部突变的特性。,图 用 , 作小波对阶跃输入和脉冲输入的处理结果,由以上分析可得,突变点的位置有时是由小波变换的过零点反映的,有时是由小波变换的极值点反映的。,一般地说,根据过零点作检测不如根据极值点。因为过零点易受噪声干扰,而且有时过零点反映的不是突变点,而是信号在慢变区间的转折点。,检测边沿宜采用如 的反对称小波;检测尖峰脉冲宜采用如 的对称小波。,结论:,要使奇异检测有效,必须满足适当条件:,, 应是某一平滑函数的一、二阶导数;,尺度a必须适当,以便使y(t)的突变点基本上能反映待分析信号x(t)的突变点;且只有在适当尺度下各突变点引起的小波变换才能避免交叠干扰。,三、基于小波变换的原油管道泄漏检测,1、泄漏检测原理,当流体输送管道因为机械、人为破坏、材料失效等原因发生泄漏时,由于管道内流体压力很高而管道外一般为大气压力,管内输送的流体在内外压差的作用下迅速流失,泄漏部位产生物质损失,这会引起发生泄漏场所的流体的密度减小,进而引起管道内此处流体的压力降低。,瞬态负压波法泄漏检测原理及定位公式,由于流体的连续性,管道中的流体速度不会立即发生改变,流体在泄漏点和与其相邻的两边的区域之间的压力产生差异,这种差异导致泄漏点上下游区域内的高压流体流向泄漏点处的低压区域,从而又引起与泄漏点相邻区域流体的密度减小和压力降低。,这种现象从泄漏点处沿管道依次向上、下游方向扩散,在水力学上称为负压波(又称为减压波)。,泄漏在管道中的总体反映就是从泄漏点处产生了同时向上、下游端传播的瞬态负压波,它的传播过程类似于声波在介质中的传播,它的传播速度是声波在管道输送流体中的传播速度,原油管道中负压力波的传播速度约在10001200米/秒之间。,在管道两端安装压力传感器能够捕捉到包含泄漏信息的瞬态负压波,就可以检测泄漏的发生,并根据泄漏产生的瞬态负压波传播到管道两端的时间差进行漏点定位。,沿管道传播的瞬态负压波中包含有泄漏的信息,由于管道的波导作用,它能够传播数十公里以上的远端。,该方法即为瞬态负压波法,它具有快速的反应速度和很高的定位精度,能够及时检测出泄漏,防止泄漏事故扩大,减少流体损失赢得宝贵的时间,是一种受到广泛重视的泄漏检测方法。,瞬态负压波泄漏定位示意图,其中: x 泄漏点距上游站测压点的距离,单位:m; L 上下游站间距,单位:m; a 负压波的传播速度,单位:m/s; t 上游站压力突变时间与下游站压力突变时间 差,单位:s。,泄漏点的计算公式为:,2、瞬态负压波泄漏定位准确的关键,由泄漏点的定位公式:,可以看出,负压波传播到上、下游传感器的时间差的精确确定,和管内负压波速度的确定是瞬态负压波定位方法的两项关键所在。,在分析泄漏引发的负压波信号序列,确定负压波信号传到管道首、末端的时刻时,一个显然的要求是首、末端压力信号序列起始时刻应该一致,这就要求统一担任首、末端数据采集系统的工控机的系统时间,可以采用全球定位系统(GPS)来定时统一各站工控机的系统时钟。 这个方案即满足了泄漏监测系统对统一时标的要求,实施也很方便,造价低廉,有很广泛的应用场所。,GPS是英文Global Positioning System的缩写,意即全球定位系统。全球定位系统利用导航卫星进行测时和测距,使在地球上任何地方的用户,都能计算出他们所处的方位。,GPS系统包括以下三大部分:(1)GPS卫星(空间部分);(2)地面支撑系统(地面监控部分);(3)GPS接收机(用户部分)。,在泄漏监测系统中使用到的是GPS的精确授时功能。,GPS开始时只用于军事目的,后来由于GPS接收机技术的发展,超大规模芯片的应用,使接收机成本不断下降,现在也已广泛应用于航海、航空、科学研究、交通运输、石油勘探、地形测量以及商业、旅游业等一切行业,甚至要渗透到个人生活的各个方面。,它可以准确测定用户的三维位置、三维速度,并给出精确的时间基准。由于GPS具有定位精度高、使用范围广、可全天候应用、用户设备简单等优点,GPS问世以来,已充分显示了其在导航、定位、授时领域的霸主地位。,管内压力波的传播速度决定于液体的弹性、液体的密度和管材的弹性,液体的体积弹性系数随品种、温度、压力的不同而不同。,式中:a 管内压力波的传播速度,m/s;K 液体的体积弹性系数,Pa; 液体的密度,kg/m3;E 管材的弹性,Pa; D 管道直径,m;e 管壁厚度,m;C1 与管道约束条件有关的修正系数。,负压波传播速度公式如下:,在实际的泄漏监测系统中,总是采集压力传感器送来的数据,再分析采集到的数据序列,从中寻找泄漏信息。精确确定泄漏引发的负压波传播到上、下游传感器的时间差,就必须先确定瞬态负压波传到管道首、末端的时刻,即需要准确地捕捉到泄漏负压波传到首、末端信号序列的对应特征点。,而由于不可避免的工业现场的电磁干扰、输油泵的振动等因素的存在,采集到的压力波形序列附加着大量的噪声(如下图所示) 。,图 原始压力信号,上图为在中石化管道储运公司沧州输油公司沧州临邑线长60公里的沧州东光输油管线的30公里处做的一次泄漏放油实验时在管道一端采集到的压力信号序列,噪声信号很强,由此信号根本无法确定负压波的边沿,因而也就不能对泄漏点定位。,如何在强噪声干扰中提取信号的特征拐点是泄漏检测与定位中必须解决的问题,下面提出采用离散小波变换确定负压波信号的特征拐点,滤波器组计算小波变换的方法。,此方法在管道实际运行中做到实时监测,取得了良好的效果。,离散小波变换及滤波器组,1、多分辨分析,1988年S.Mallat在构造正交小波基时,将计算机视觉领域内的多分辨率思想率先引入小波变换,提出了多分辨分析(Multi-Resolution Analysis)的概念,在空间的概念上形象地说明了小波变换的多分辨特性,并使用多分辨分析将此之前Meyer等提出的各种具体小波基的构造法统一起来。在Burt and Adelson图像分解和重构的塔式算法的启发下,基于多分辨率框架,提出了塔式多分辨率分解和重构算法,给出了正交小波的构造方法以及正交小波变换的快速算法Mallat快速小波分解和重构算法。,对于多分辨分析的理解,可以用一个三层的分解进行说明,其小波分解树如下图所示。,图 三层多分辨分析小波分解树结构图,多分辨分析只是对低频部分进行进一步分解,而高频部分则不予考虑。分解关系式为:S=A3+D3D2D1。在图中,只是以一个层分解进行说明,如果要进行进一步的分解,则可以把低频部分A3分解成低频部分A4和高频部分D4,以下再分解依此类推。,多分辨分析分解的最终目的是力求构造一个在频率上高度逼近 空间的正交小波基,这些频率分辨率不同的正交小波基相当于带宽各异的带通滤波器。,从上图的多分辨分析树型结构可以看出,多分辨率分析只对低频空间进行进一步的分解,使频率的分辨率变得越来越高。,下面分析多分辨分析是如何构造正交小波基的。,空间L2(R)中的多分辨分析是指 中满足如下条件的一个空间序列 :, 单调性: ,对任意jZ;, 逼近性: ;, 伸缩性: ,伸缩性体现了尺度的变化、逼近正交小波函数的变化和空间的变化具有一致性;, 平移不变性:对任意kZ,有, Riesz基存在性:,可以称(t)为尺度函数(Scaling Function),并可以定义如下函数:,该函数系 是规范正交的。,设以Vj表示小波分解树结构图中的低频部分Aj,Wj表示分解中的高频部分Dj,则多分辨分析的子空间V0是一个有限个子空间的逼近,即:,上式中的空间列 具有如下性质:,设法找出一个确定的函数 ,使得对每个,,函数系 构成空间Wj的规范正交基,,其中,,(1),若令 代表分辨率为2-j的函数 的逼近(即函数f的低频部分或“粗糙像”),而,代表逼近的误差(即函数f的高频部分或,“细节”部分),则下式,意味着:,注意到f=f0,所以上式可简写为:,(2),(3),上式(3)表明,任何函数 都可以根据分辨率为2-N时f(t)的低频部分(“粗糙像”)和分辨率2-j(1jN)下f(t)的高频部分(“细节”部分)完全重构,这也从另一侧面表示Mallat塔式重构算法。,(3),从包容关系 ,可很容易得到尺度函数(t)的一个极为有用的性质。,注意到 ,所以 可以用V-1子空间的基函数 展开,令展开系数为hk,则,这就是尺度函数的双尺度方程。,(4),另一方面,由于 ,故,这就意味着小波基函数 可以用V-1子空间的正交基 展开,令展开系数为gk,即有,这就是小波函数的双尺度方程。,(5),由双尺度方程式(4)和式(5)可知,尺度函数与小波基函数的构造归结为系数 的设计。,小波基函数 可由尺度函数 的平移和伸缩的线性组合获得。,若令,则把尺度函数和小波基函数的设计可以归结为滤波器H(),G()的设计。,(4),(5),尺度函数的双尺度方程:,小波函数的双尺度方程:,构造正交小波时,滤波器H(),G()必须满足以下三个条件:,(6),(7),(8),联合求解式(6)、(7)、(8)可得,(9),所以,要设计正交小波,只需要设计滤波器H()。,综上分析,为了使 构成Vj子空间的正交基, 应该具有下列基本性质:, 尺度函数的容许条件, 。, 能量归一化条件: 。, 尺度函数 具有正交性,即, 跨尺度的尺度函数 和 相关。, 尺度函数 与基小波函数 正交,即, 基小波函数 和 相关。,将尺度函数的容许条件与小波的容许条件,作比较可知,尺度函数的傅里叶变换 具有低通滤波特性(相当于一个低通滤波器),而小波基函数的傅里叶变换 则具有高通滤波特性(相当于一个带通滤波器)。,多分辨分析可以对信号进行有效的时频分解,它不同于短时傅里叶变换对信号频带的等间隔划分,它的尺度是按二进制划分的,所以在高频频段其频率分辨率较差,而在低频频段其时间分辨率较差,即对信号的频带进行指数等间隔划分(具有等Q结构)。,2、一维Mallat算法,S.Mallat在Burt and Adelson图像分解和重构的塔式算法的启发下,基于多分辨率框架,提出了一种具有完美数学描述的塔式多分辨率离散小波分解与重构算法,即Mallat算法。,Mallat算法的本质是:不需要知道小波函数和尺度函数的具体结构,仅根据系数就可实现信号的分解和重构。而且这种算法可以使每次小波分解后信号长度减半,大大减少了小波变换的复杂度,因而它是一种快速算法。,Mallat快速算法在小波变换中的地位,与快速傅里叶变换在傅里叶变换中的地位相当。,设(t)为尺度函数,(t)为基本小波,设信号f(t) L2(R),将信号进行分解,用A表示低频,D表示高频,下标数字表示小波分解的尺度层数,已得到f(t)在2-j分辨率下的粗糙像Ajf Vj,Vjj Z构成L2(R)的多分辨分析,从而有:,式中,,(10),式(10)可写成如下形式:,由尺度函数的双尺度方程可得:,利用尺度函数的正交性,有:,(12),(13),同理,由小波函数的双尺度方程可得:,(14),(11),由式(11)(14),可得:,(15),(16),(17),引入无穷矩阵 ,,其中 ,,则式(15)、(16)、(17),可分别表示为更简洁的形式:,其中 和 分别是H和G对应的对偶算子,或分别理解为H和G的共扼转置矩阵,满足正交性条件,,。H和G分别称为低通滤波器和带通滤波器。,(18),(19),Mallat一维分解算法,Mallat一维重构算法,式(18)即为Mallat一维分解算法;式(19)即为Mallat一维重构算法,它们都是递推的快速算法,可由下图来表示:,(b)重构算法,图 Mallat小波分解和重构算法示意图,2:表示每个样本点中入 1个零值,(a)分解算法,2:表示2个样本点中取1个,将实时采集获得的管道压力信号通过两通道滤波器分为低频概貌和高频细节输出,定位负压波的位置实际就是确定高频细节当中最小值的位置。,在实际应用当中,基于虚拟仪器LabVIEW开发了原油管道泄漏实时监测系统,其中包括泄漏定位模块,其检测定位原理即基于离散小波变换和Mallat算法。如下图所示。,图 泄漏点定位模块,泄漏点定位部分精确地定位出泄漏点的位置,同时将泄漏报警时间、泄漏点位置等参数写入泄漏信息数据库以供日后总结分析。定位部分包括自动定位和手动定位两部分,手动定位需要操作人员根据经验,通过移动光标捕捉两端波形对应拐点。,补充:MATLAB中的小波工具包简介,MATLAB中常用小波函数介绍,1、Haar小波,Haar小波是小波分析发展过程中用的最早的小波,也是最简单的小波。Haar小波本身是一个阶跃函数如下图所示。,Haar小波,Haar小波的支集长度为1,滤波器长度为2。,2、Daubechies小波,Daubechies小波是由著名小波学者Ingrid Daubechies所创造的。,Daubechies系列的小波简写为dbN,其中N表示阶数,db是小波名字的前缀,除db1(等同于Haar小波)外,其余的db系列小波函数都没有解析的表达式。,db系列小波函数,db4,db8,3、SymletsA(SymN)小波族,sym小波的构造类似于db小波族,两者的差别在于sym小波有更好的对称性。,sym系列小波函数,sym4,sym8,同前面db小波族相比,sym小波族有更好的对称性。其它的性质如连续性、支集长度、滤波器长度都同db小波族一致。,4、Biorthogonal(BiorNr.Nd)小波族,这是一族双正交小波(也叫半正交小波),与正交小波的区别在于:正交小波的伸缩和平移构成的基函数完全正交;双正交小波对不同尺度伸缩下的小波函数之间也有正交性,而同尺度之间通过平移得到的小波函数系之间没有正交性。,双正交小波用于分解和重构的小波不是同一个函数,相应的滤波器也不能由同一个小波生成。,bior系列小波函数,bior2.4,bior4.4,在双正交小波的命名中,Nr和Nd分别是和重构和分解滤波器长度有关的参数。下图为bior2.4和bior4.4小波函数。,5、Coiflet(coifN)小波族,Coiflet小波族是Daubechies提出的另一个小波系,它有更长的支集长度和更大的消失矩,对称性比较好。,coif系列小波函数,coif3,coif5,下图为coif3和coif5小波函数。从图中可以看出,它们在对称性上的特点。,6、Morlet小波,Morlet小波是一个具有解析表达式的小波,但它不具备正交性。Morlet小波的小波函数如下图所示。,morlet小波函数,7、Mexican Hat小波,类似于morlet小波,Mexican Hat小波同样是有解析的小波函数,但同样不具有正交性。 Mexican Hat小波的小波函数如下图所示。,Mexican Hat小波函数,8、Meyer小波,Meyer小波是在频域定义的具有解析形式的正交小波。Meyer小波的尺度函数和小波函数如下图所示。,Meyer小波,尺度函数,小波函数,一维离散小波变换,一维离散小波变换实现的算法一般是mallat算法。,即先对较大尺度的信号进行小波变换,再选取其中的低频部分在原尺度的1/2尺度上再进行小波变换。,对于一个长度为N的信号s,第一步从原始信号开始,产生两组参数,一组是作用低通滤波器得到的近似信号,另一组是作用高通滤波器得到的细节信号。,1、一维离散小波分解算法,在物理信号中,低频部分是表征信号本身特征的,而高频部分则是表征信号细微差别的。,比如,声音信号,如果只保留低频信号,仍可以辨别出说话的内容,但是可能不太容易辨别出说话的人。但如果去除了低频部分,就只能听到一些噪声。,在MATLAB小波工具箱中,实现多尺度分解的函数是:wavedec,该函数的调用格式为:,C,L = wavedec(X,N,wname)C,L = wavedec(X,N,Lo_D,Hi_D),式中,X:待分解的信号; N:分解层数; wname:使用的小波函数; Lo_D:给定的低通滤波器; Hi_D:给定的高通滤波器。,即:在第二种调用方式中,信号X在给定的低通滤波器和高通滤波器下进行多尺度分解。,在MATLAB中,离散小波变换分解算法主要使用如下几个常用指令:,dwt; 用于信号的单层分解wavedec; 用于信号的多层分解wmaxlev; 在多层分解前求最大的分解层数,例:对白噪声信号noissin进行3层小波分解,所用小波函数为db4。,load noissin; 读入白噪声信号,s=noissin(1:1000); 取信号的前1000个采样点c,l=wavedec(s,3,db4); 对信号做3层多尺度分解cd1,cd2,cd3=detcoef(c,l,1,2,3);得到三个尺度的细节系数ca3=appcoef(c,l,db4,3); 得到尺度3的近似系数figure;subplot(511);plot(1:1000,s);title(s); 绘制原始信号subplot(512);plot(ca3);title(ca3); 绘制尺度3的近似系数subplot(513);plot(cd3);title(cd3); 绘制尺度3的细节系数,subplot(514);plot(cd2);title(cd2); 绘制尺度2的细节系数subplot(515);plot(cd1);title(cd1); 绘制尺度1的细节系数,Noissin信号在db4小波下三层分解的近似系数和细节系数,2、一维离散小波重建算法,重建运算是小波变换的逆变换,也就是把分解得到的近似系数和细节系数叠加得到原始信号。,与小波分解过程类似,重构过程首先从尺度最低的近似系数cAj和细节系数cDj开始,通过作用低频和高频重构滤波器恢复出上一尺度的近似信号cAj-1,把这个过程继续下去,直到得到原始信号s。,单单做这个过程没有什么意义,只是把分解的信号又重建一下。实际中,往往会在重建之前都要对分解系数做各种处理,以达到预期的目的。,在MATLAB中,用于离散小波重建算法的命令主要有如下几个:,idwt; 用于单层小波重建waverec; 用于多层小波重建原始信号,要求输入 参数同小波分解得到结果的格式一致upcoef; 用于重建小波系数至上一层次,要求输入 参数同小波分解得到结果的格式一致,用于得到某一层次的小波系数的命令主要有如下几个:,detcoef; 求得某一层次的细节系数appcoef; 求得某一层次的近似系数upwlev; 重建组织小波系数的排列形式,load sumsin;读入信号sumsin,该信号为不同频率正弦波信号的叠加s=sumsin;c,l=wavedec(s,4,sym1);把信号s用sym1小波分解到第四层,分解的系数存到数组c中,各层分解后的长度存到数组l中;nc,nl=upwlev(c,l,sym1);通过第四层小波系数重建第三层小波近似系数,把三层的系数存放在数组nc中,三层分解的长度存放到数组nl中figure;subplot(311);plot(s);title(原始信号); 绘制原始信号,例:利用upwlev对正弦波信号进行重构。,subplot(312);plot(c);title(做4层wavedec得到的结果);subplot(313);plot(nc);title(做3层wavedec得到的结果);,用upwlev重建小波系数示意图,信号的扩展,在用计算机实现小波变换的过程中,一个实际的问题就是信号长度和滤波器长度的匹配问题。,对有限支集长度的滤波器,在离散小波变换中,对每个点的变换都需要用到其它点的信息,所以如果点与边界的距离不足滤波器的长度,就无法得到做变换所需要的完全信息,从而导致滤波器不可逆。,为了保证滤波器构成的矩阵是可逆的,就需要对滤波器做调整,而这种调整也等价于对信号的调整,调整的方法就是在每层的变换中都把信号边界做扩展,使得最靠近边界的点也能获得小波变换所需要的完全信息。,以求细节系数的高通滤波为例,过程如下图所示。下图中,高通滤波器假设支集长度为4,那么在求第一个细节系数d0的时侯,就需要用到信号的扩展s-1。在信号的另一侧,如为了求第n个细节系数dn,同样需要用到信号的扩展sn+1。,图 小波变换做高通滤波时对信号的扩展,MATLAB中,对信号的扩展方式主要有以下几种。,零填充:是最简单的一种方式。把边界以外到滤波器支集以内所需要的信息都置为0。这种方式最大的问题是在边界处不连续。 对称扩展:令s-1=s1,s-2=s2,以此类推,在另一侧的边界做同样处理,它是MATLAB小波工具箱中缺省的扩展方式,很好地解决了边界点的不连续问题。但是它的一阶导数仍然是不连续的。 一阶光滑扩展:保证边界以外的信息满足一阶导数的连续,那么所有边界外的点应该落在通过s0,s1的直线上,另一侧边界也做同样处理。 零阶光滑扩展:保证边界以外的信号与原信号连续,那么只需要令:si=s0(in)。,在对信号的扩展中,还有一点需要注意,为了保持小波变换的可逆性,分解阶段和重建阶段必须指定相同的扩展方式。,在MATLAB中,指定信号扩展方式的命令是dwtmode。其调用格式如下:,ST = dwtmode(mode),ST = dwtmode,主要的扩展方式如下:,sym:对称扩展;zpd:零填充;spd或sp1:一阶光滑扩展;sp0:零阶光滑扩展。,小波变换用于信号降噪,前面介绍了如何利用小波工具箱把信号分解以及重建的方法,而小波分析之所以强大,就在于它能将信号中我们关心的成分尽可能详细地展示出来,从这个意义上讲,单单分解和重建是没有意义的,小波分析最重要的过程就是对小波分解系数的处理。,小波分析的一般过程就是先把信号分解为小波系数,然后对分解出来的系数根据问题的要求做一些处理,再用小波重建方法恢复信号。,下面通过小波变换用于信号降噪来描述小波分析的主要方法。,光滑性:在大部分情况下,降噪后的信号应该至少和原信号具有同等的光滑性; 相似性:降噪后的信号和原信号的方差估计应该是最坏情况下的方差最小。,1、信号降噪的准则,2、小波分析用于降噪的过程,分解过程:选定一种小波,对信号进行N层小波分解; 作用阈值过程:对分解得到的各层系数选择一个阈值,并对细节系数作用软阈值处理; 重建过程:将处理后的系数通过小波重建恢复原始信号。,携带信息的原始信号在频域或小波域的能量相对集中,表现为能量密集区域的信号分解系数的绝对值比较大,而噪声信号的能量谱相对分散,所以其系数的绝对值小。这样,可以通过作用阈值的方法过滤掉绝对值小于一定阈值的小波系数,从而达到降噪的效果。,3、从原始信号确定各级阈值,在小波分析用于降噪的过程中,核心的步骤就是在系数上作用阈值。因为阈值的选取直接影响降噪的质量,所以人们提出了各种理论的和经验的