首页 理论教育 离散小波变换的多分辨率分析

离散小波变换的多分辨率分析

时间:2023-02-13 理论教育 版权反馈
【摘要】:小波变换作为能随频率的变化自动调整分析窗大小的分析工具,自20世纪80年代中期以来得到了迅猛的发展,并在信号处理、故障诊断、计算机视觉、图像处理、语音分析与合成等众多的领域得到应用。ψ又称为基本小波,或母小波。这两个式子中出现的a2是由于定义小波变换时在分母中出现了,而式中又要对a作积分所引入的。但小波变换的Parseval定理稍为复杂,它不但要有常数加权,而且以cψ的存在为条件。

第一章 小波变换基础

第一节 小波变换概述

把函数分解成一系列简单基函数的表示,无论是理论上,还是实际应用中都有重要意义。

1807年法国数学家傅里叶(J. Fourier 1768—1830)提出并证明了将周期函数展开为正弦级数的原理,奠定了傅里叶级数理论的基础。傅里叶级数用于分析周期性的函数或分布,理论分析时经常假定周期是2π,定义如式[1,2]

其中

然而,被分析函数的性质并不能完整地由傅里叶系数来刻划,即不可能仅根据傅里叶系数大小就预知函数的性质(如大小、正则性)。

傅里叶变换(Fourier Transform,FT)的定义如式(1.1.1)、(1.1.2):

通过引入广义函数或分布的概念,可获得奇异函数(如冲击函数)的傅里叶变换的存在。对于时域的常量函数,在频域将表现为冲击函数,表明具有很好的频域局部化性质。由式(1.1.1)可知,为了得到F(ω),必须有关于f(x)的过去和未来的所有知识,而且f(x)在时域局部值的变化会扩散到整个频域,也就是F(ω)的任意有限区域的信息都不足以确定任意小区域的f(x)。在时域,哈尔(Haar)基是一组具有最好的时域分辨能力的正交基,它在时域上是完全局部化的,但在频域的局部化却很不好,这是由于哈尔基的两个缺点:缺乏正则性与缺乏振动性。研究者们希望寻找关于空间变量(或时间变量)与频域变量都同时好的希尔伯特(Hilbert)基。

为此,人们提出了短时傅里叶变换(Short Time Fourier Transform,STFT)的概念:

定义1.1 若W∈L2(R)选择使得W与它的傅里叶变换满足:

那么使用W作为窗函数,在式(1.1.3)中引入的窗口傅里叶变换称为STFT:

当窗函数选择为高斯(Gaussian)函数时,则为Gabor变换。

STFT的缺点是分析窗的大小和形状是恒定的。因为频率与周期成反比,所以反映信号的高频成分需要窄的时间窗,而反映信号的低频成分需要宽的时间窗,STFT无法满足要求。此外,STFT的冗余很大,增加了不必要的计算量。

小波变换作为能随频率的变化自动调整分析窗大小的分析工具,自20世纪80年代中期以来得到了迅猛的发展,并在信号处理、故障诊断、计算机视觉、图像处理、语音分析与合成等众多的领域得到应用。

一、小波变换的定义[3,4,5]

给定一个基本函数ψ(t),令

式中a,b均为常数,且a>0。显然,ψa,b(t)是基本函数ψ(t)先作移位再作伸缩以后得到的。若a,b不断地变化,可得到一族函数ψa,b(t)。给定平方可积的信号x(t),即x(t)∈L2(R),则x(t)的小波变换(Wavelet Transform,WT)定义为

式中a,b和t均是连续变量,因此该式又称为连续小波变换(CWT)。如无特别说明,式中及以后各式中的积分都是从−∞到+∞。信号x(t)的小波变换WTx(a,b)是a和b的函数,b是时移,a是尺度因子。ψ(t)又称为基本小波,或母小波。ψa,b(t)是母小波经移位和伸缩所产生的一族函数,称之为小波基函数,或简称小波基。这样,(1.1.5)式的WT又可解释为信号x(t)和一族小波基的内积。

母小波可以是实函数,也可以是复函数。若x(t)是实信号,ψ(t)也是实的,则WTx(a,b)也是实的,反之,WTx(a,b)为复函数。

在(1.1.4)式中,b的作用是确定对x(t)分析的时间位置,也即时间中心。尺度因子a的作用是把基本小波ψ(t)作伸缩。由ψ(t)变成,当a>1时,若a越大,则的时域支撑范围(即时域宽度)较之ψ(t)变得越大,反之,当a<1时,a越小,则的宽度越窄。这样,a和b联合起来确定了对x(t)分析的中心位置及分析的时间宽度。

这样,(1.1.5)式的WT可理解为用一族分析宽度不断变化的基函数对x(t)作分析,这一变化正好适应了对信号分析时在不同频率范围所需要不同的分辨率这一基本要求。

(1.1.4)式中的因子是为了保证在不同的尺度a时,ψa,b(t)始终能和母函数ψ(t)有着相同的能量,即

,则dt=adt′,这样,上式的积分即等于

令x(t)的傅里叶变换为X(Ω),ψ(t)的傅里叶变换为Ψ(Ω),由傅里叶变换的性质,ψa,b(t)的傅里叶变换为:

由Parsevals定理,(1.1.5)式可重新表示为:

此式(1.1.7)即为小波变换的频域表达式。

二、小波反变换及小波容许条件

下述定理给出了连续小波反变换的公式及反变换存在的条件。

定理 1.1 设x1(t),x2(t)和ψ(t)∈L2(R),x1(t),x2(t)的小波变换分别是WTx1(a,b)和WTx2(a,b),则

式中

其中Ψ(Ω)为ψ(t)的傅里叶变换。

(1.1.8)式实际上可看作是小波变换的Parseval定理。该式又可写成更简单的形式,即

进一步,如果令x1(t)=x2(t)=x(t),有

该式更清楚地说明,小波变换的幅平方在尺度——位移平面上的加权积分等于信号在时域的总能量,因此,小波变换的幅平方可看作是信号能量时-频分布的一种表示形式。

(1.1.8)和(1.1.10)式中对a的积分是从0~∞,这是因为假定a总为正值。这两个式子中出现的a−2是由于定义小波变换时在分母中出现了,而式中又要对a作积分所引入的。

傅里叶变换中的Parseval定理,即时域中的能量等于频域中的能量。但小波变换的Parseval定理稍为复杂,它不但要有常数加权,而且以cψ的存在为条件。

定理1.2 设x(t),ψ(t)∈L2(R),记Ψ(Ω)为ψ(t)的傅里叶变换,若

则x(t)可由其小波变换WTx(a,b)来恢复,即

在定理1.1和定理1.2中,结论的成立都是以cψ<∞为前提条件的。(1.1.9)式又称为“容许条件(admissibility condition)。该容许条件含有多层的意思:

第一,并不是时域的任一函数ψ(t)∈L2(R)都可以充当小波。其可以作为小波的必要条件是其傅里叶变换满足该容许条件;

第二,由(1.1.9)式可知,若cψ<∞,则必有Ψ(0)=0,否则cψ必趋于无穷。小波函数ψ(t)必然是带通函数;

第三,由于,因此必有

这一结论指出,ψ(t)的取值必然是有正有负,也即它是振荡的。

以上三条勾画出了作为小波的函数所应具有的大致特征,即ψ(t)是一带通函数,它的时域波形应是振荡的。此外,从时-频定位的角度,总希望ψ(t)是有限支撑的,因此它应是快速衰减的。时域有限长且振荡的这一类函数即是被称作小波(wavelet)的原因。

由上述讨论,ψ(t)自然应和一般的窗函数一样满足:

尺度a常按a=2j来离散化,j∈Z。由(1.1.6)式,对应的傅里叶变换2j/2Ψ(2jΩ)e−jΩb,由于需要在不同的尺度下对信号进行分析,同时也需要在该尺度下由WTx(a,b)来重建x(t),因此要求是有界的,当j由−∞~+∞时,应有

(1.1.6)式中0<A≤B<∞。该式称为小波变换的稳定性条件,它是在频域对小波函数提出的又一要求。满足(1.1.12)式的小波称作“二进(dyadic)”小波。

三、连续小波变换的计算

在(1.1.5)式关于小波变换的定义中,变量t, a和b都是连续的,在计算机上实现一个信号的小波变换时,t,a和b均应离散化。对a离散化最常用的方法是取,j∈Z,并取a0=2,这样2ja=。对于a按2的整次幂取值所得到的小波习惯上称之为“二进(dyadic)”小波。对这一类小波的小波变换,可用有关离散小波变换的方法来实现。然而取a=2j,j∈Z,在实际工作中有时显得尺度跳跃太大。当希望a任意取值(a>0),也即在a>0的范围内任意取值时,这时的小波变换即是连续小波变换。

计算(1.1.5)式的最简单的方法是用数值积分的方法,即,令

由于在t=k~k+1的区间内,x(t)=x(k),所以上式又可写为:

由该式可以看出,小波变换WTx(a,b)可看作是x(k)和的卷积后的累加所得到的结果,卷积的中间变量是t,卷积后的变量为b及a。Matlab中的cwt.m即是按此思路来实现的。其具体过程大致如下:

第一,先由指定的小波名称得到母小波ψ(t)及其时间轴上的刻度,假定刻度长为0~N1−;

第二,从时间轴坐标的起点开始求积分,k=1,…,N1−;

第三,由尺度a确定对上述积分值选择的步长,a越大,上述积分值被选中的越多;

第四,求x(k)和所选中的积分值序列的卷积,然后再作差分,即完成(1.1.13)式。

本方法的不足之处是在a变化时,(1.1.13)式中括号内的积分、差分后的点数不同,也即和x(k)卷积后的点数不同。解决的方法是在不同的尺度下对ψ(t)作插值,使其在不同的尺度下,在其有效支撑范围内的点数始终相同。

例1.1 令x(t)为一正弦加噪声信号,它取自Matlab中的noissin.mat。对该信号作CWT,a分别等于2和128,a=2时,小波变换的结果对应信号中的高频成分,a=128时,小波变换对应信号中的低频成分。其源始信号及变换结果见图1.1(a),(b)和(c)。

图1.1 信号“noissin”的小波变换

图1.2 多尺度下小波变换的灰度表示

例1.2 仍然使用例1.1的信号“noissin”,对其作CWT时a分别取10,30,60,90,120及150。所得到的图1.2是在各个尺度下的小波系数的灰度图。颜色越深,说明在该尺度及该位移(水平轴)处的小波系数越大。此例旨在说明对小波变换的结果有不同表示方式。

四、尺度离散化的小波变换

在(1.1.7)式定义了信号x(t)的连续小波变换,式中a,b和t都是连续变量。为了在计算机上有效地实现小波变换,t自然应取离散值,a和b也应取离散值。从减少信息冗余的角度,a和b也没有必要连续取值。

a和b形成了一个二维的“尺度—位移”平面。前已述及,a越大,Ψ(aΩ)对应的频率越低,反之,对应的频率越高。因此,ab−平面也可视为“时-频平面”。对同一个信号x(t),已给出过不同的表示形式,如STFT,Gabor变换,WVD及本章的WT。

现重写几个有关的公式,即

其中,(1.1.15)式是用时-频平面离散栅格(m,n)上的点来表示x(t),即Gabor展开,(1.1.16)式是具有双线性变换的表示形式,它和其他三种表示形式有较大的区别。(1.1.14)和(1.1.17)式说明同一信号x(t)在时-频平面上具有不同的表示形式。在第二章会指出,(1.1.14)式的反变换是有信息冗余的,即不需要STFT(t,Ω)的所有的值即可恢复x(t)。同理,(1.1.17)式的小波变换也存在着信息冗余。在这两个式子中,只需取时-频平面上的离散栅格处的点即可。问题的关键是如何决定a和b抽样的步长以保证对x(t)的准确重建。

下面,首先考虑尺度a的离散化,然后再考虑a和b的同时离散化。

1. 尺度离散化的小波变换

目前通用的对a离散化的方法是按幂级数的形式逐步加大a,即令,j∈Z。若取a0=2,则

称为“半离散化二进小波”,而

称为二进小波变换。

设母小波ψ(t)的中心频率为Ω0,带宽为ΔΩ,当a=2j时,ψj,b(t)的中心频率变为(Ω,带宽ΔΩ若a=2j+1时,ψj+1,b(t)的中心频率和带宽分别是:−1。从对信号作频域分析的角度,希望当a由2j变成2j+1时,ψj,b(t)和ψj+1,b(t)在频域对应的分析窗[(Ωj)0−ΔΩj,(Ωj)0Ωj]和[(Ωj+1)0−ΔΩj+1,(Ωj+1)0Ωj+1]能够相连接。这样,当

j由0变至无穷时,ψj,b(t)的傅里叶变换可以覆盖整个Ω轴。显然,若令母小波ψ(t)的(Ω)0=3ΔΩ,则上面两个频域窗首尾相连,即

首尾相连。通过对母小波作合适的调制,可以方便地做到(Ω)0=3ΔΩ

现在,来讨论如何由(1.1.18)式的WTx(j,b)来恢复x(t),设是ψ(t)的对偶小波,并令取类似的形式,即

这样,通过对偶小波,希望能重建x(t):

为了寻找和ψ(t)应满足的关系,现对上式作如下改变:

式中ℑ代表求傅里叶变换。有

则(1.1.20)式的右边变成X(Ω)的傅里叶反变换,自然就是x(t)。

对于满足容许条件的小波ψ(t),当a=2j,j∈Z时,其二进制小波ψj,b(t)对应的傅里叶变换应满足(1.1.12)式的稳定性条件。这样,结合(1.1.12)和(1.1.21)式,可由下式得到对偶小波

由于(1.1.22)式的分母满足(1.1.12)式,因此有

这样,对偶小波也满足稳定性条件,也即,总可以找到一个“稳定的”对偶小波由(1.1.19)式重建出x(t)。下面的定理更完整地回答了在半离散二进小波变换情况下x(t)的重建问题。

定理1.3 如果存在常数A>0,B>0,使得

如果

满足

该定理指出,若ψ(t)的傅里叶变换满足稳定性条件,则x(t)在a=2j,j∈Z上的小波变换的幅平方的和是有界的。进而,ψ(t)和的傅里叶变换若满足(1.1.24)式(也即(1.1.21)式),则x(t)可由(1.1.25)式重建。

总之,若ψ(t)满足容许条件,且再满足稳定性条件,由二进小波变换WTx(j,b)总可以重建x(t),也即一个满足稳定性条件的对偶小波总是存在的。但是,满足稳定性条件的对偶小波不一定是唯一的。如何构造“好”的小波ψ(t)及得到唯一的对偶小波是小波理论中的重要内容。

若(1.1.23)式的稳定性条件满足,则(1.1.9)式的容许条件必定满足,且

从而,由连续小波变换WTx(a,b)总可以恢复x(t)。

以上讨论的是仅对a作二进制离散化情况,现考虑a和b同时离散化的相应理论问题。

2. 离散栅格上的小波变换

,j∈Z,可实现对a的离散化。若j=0,则ψj,b(t)=ψ(t−b)。欲对b离散化,最简单的方法是将b均匀抽样,如令b=kb0,b0的选择应保证能由WTx(j,k)来恢复出x(t)。当j≠0时,将a由变成时,即是将a扩大了a0倍,这时小波ψj,k(t)的中心频率比ψj−1,k(t)的中心频率下降了a0倍,带宽也下降了a0倍。因此,这时对b抽样的间隔也可相应地扩大a0倍。由此可以看出,当尺度a分别取时,对b的抽样间隔可以取,这样,对a和b离散化后的结果是:

对给定的信号x(t),(1.1.7)式的连续小波变换可变成如下离散栅格上的小波变换,即

此式称为“离散小波变换(Discrete Wavelet Transform,DWT)”。注意式中t仍是连续变量。

,可以仿照傅里叶级数和Gabor展开那样来重建x(t),即

该式称为小波级数,称为小波系数,的对偶函数,或对偶小波。

对任一周期信号x(t),若周期为T,且2x(t)∈L(0,T),则x(t)可展成傅里叶级数,即

式中X(kΩ0)是x(t)的傅里叶系数,它由下式求出:

小波级数和傅里叶级数形式上类似,但其物理概念却有着明显的不同:

(1)傅里叶级数的基函数,是一组正交基,即。而小波级数所用的一族函数不一定是正交基,甚至不一定是一组“基”;

(2)对傅里叶级数来说,基函数是固定的,且分析和重建的基函数是一样的,即都是(差一负号);对小波级数来说,分析所用的函数是可变的,且分析和重建所用的函数是不相同的,即分析时是,而重建时是

(3)在傅里叶级数中,时域和频域的分辨率是固定不变的,而小波级数在a,b轴上的离散化是不等距的,这正体现了小波变换“变焦”和“恒Q”性的特点。

第二节 离散小波变换的多分辨率分析

在实际应用中,特别是在计算机上实现小波变换时,信号总要取成离散的,因此,研究a,b及t都是离散值情况下的小波变换,进一步发展一套快速小波变换算法将更有意义。由Mallat和Meyer自20世纪80年代末期所创立的“多分辨率分析”技术在这方面起到了关键的作用。该算法和多抽样率信号处理中的滤波器组及图像处理中的金字塔编码等算法结合起来,构成了小波分析的重要工具。

本节将详细讨论多分辨率分析的定义、算法及应用[3,4,5]

一、多分辨率分析的定义

Mallat给出了多分辨率分析的定义:

设{Vj}j∈Z是L2(R)空间中的一系列闭合子空间,如果它们满足如下六个性质,则说{Vj},j∈Z是一个多分辨率近似。这六个性质是:

1. ∀(j,k)∈Z2,若x(t)∈Vj,则x(t−2jk)∈Vj (1.2.1)

2. ∀j∈Z,Vj⊃Vj+1,即…V0⊃V1⊃V2…Vj⊃Vj+1… (1.2.2)

3. j∀∈Z,若x(t)∈Vj,则

4.

5.

6. 存在一个基本函数θ(t),使得{θ(t−k)},k∈Z是V0中的Riesz基。

现对以上性质作一些直观的解释:

性质1说明,空间Vj对于正比于尺度2j的位移具有不变性,也即函数的时移不改变其所属的空间。在上一章对(a,b)作二进制离散化时曾说明,若令,则b0应取b=2jkb,将b0归一化为1,则

所以,(1.2.1)式实际上应等效为:

j

∀∈Z,若x(t)∈Vj,则x(t−k)∈Vj

这是因为j∀∈Z,必有2j∈Z;

性质2说明,在尺度2j(或j)时,对x(t)作的是分辨率为2j−的近似,其结果将包含在较低一级分辨率2−j−1时对x(t)近似的所有信息,此即空间的包含,也即(1.2.2)式;

性质3是性质2的直接结果。在Vj+1中,函数作了二倍的扩展,分辨率降为2−j−1,所以应属于Vj+1

性质 4 说明当j→∞时,分辨率2−j→0,这时将会丢失x(t)的所有信息,也即

从空间上讲,所有Vj(j=−∞~+∞)的交集为零空间;

性质5是性质4的另一面,即当j→−∞时,分辨率2j−→∞,那么信号x(t)在该尺度下的近似将收敛于它自身,即

从空间上讲,即是所有Vj(j=−∞~+∞)的并集收敛于整个L2(R)空间;

性质6说明了V0中Riesz基的存在性问题,并将要由此引出V0,V1,…,Vj中正交基的存在性问题,因此,需要着重加以解释。

设V0是一Hilbert空间(注:能量有限的空间L2(R)即是Hilbert空间),{θk=θ(t−k)},k∈Z是V0中的一组向量,其个数与V0的维数一致。自然,V0中的任一元素x都可表示为θk的线性组合,即

(1){θk=θ(t−k)},k∈Z之间是线性无关的,且

(2)存在常数0<A≤B<∞使得

则θ(t−k),k∈Z是V0中的Riesz基。

Riesz基本身是一个标架,但它比标架的要求要高,即kθ之间是线性无关的,但它又比正交基要求低,即并不要求kθ之间一定要两两正交。(1.2.7)式的能量约束关系保证了(1.2.6)式对x(t)表示的数值稳定性。

下述定理给出了在V0中存在Riesz基的充要条件。

定理 1.4 θ(t−k),k∈Z是V0中的 Riesz 基的充要条件是存在常数A>0, B>0

使得

式中是θ(t)的傅里叶变换。

在实际工作中,人们总偏爱正交基。下述定理给出了如何由Riesz基构造正交基的方法。

定理1.5 令{Vj},j∈Z是一多分辨分析,φ(t)是一尺度函数,若其傅里叶变换可由下式给出:

并令

则φj,k(t)是Vj中的正交归一基,对所有的j∈Z,式中是产生Riesz基的基本函数θ(t)的傅里叶变换。

定理1.4给出了V0空间中Riesz基的存在性,定理1.5给出了由Riesz基过渡到正交基的方法。在实际工作中,找到一个正交归一的基函数{φ(t−k)},k∈Z并不太容易,但找到一组Riesz基{θ(t−k)},k∈Z却比较容易。具体步骤是:第一,由θ(t)作FT得;第二,由(1.2.9)式求Φ(Ω);第三,由Φ(Ω)作逆傅里叶变换得φ(t),则{φ(t−k)},k∈Z即是一组正交基。

利用此方法构造Battle−Lemarie小波,其思路是令θ(t)为一个三角波,其傅里叶变换为

可以证明,{θ(t−k)},k∈Z构成一组Riesz基,但是,θ(t−k)之间并不正交,可以求出:

显然,是有界的,满足(1.2.8)式所要求的Riesz基的频域条件。按(1.2.9)式,可以求出

由Φ(Ω)作反变换即可得到φ(t)。{φ(t−k)},k∈Z即可构成一组正交基。

二、空间Vj、Wj中信号的分解

由上两节关于频率轴剖分的思想,φ(t)应是V0中的低通函数,ψ(t)应是W0中的带通函数。将φ(t)归一化,有

定理1.5已指出,{φ(t−k)},k∈Z是V0中的正交归一基,{φj,k(t),j,k∈Z2}是 Vj中的正交归一基。这样,可将x(t)∈L2(R)按此基函数逐级进行分解。

1. 子空间V0

令P0x(t)是在V0中的投影,则

式中a0(k)是加权系数,它应是一个离散序列。由φ(t−k)的正交性质,有

φ0,k(t)和P0x(t)作内积实质上是φ0,k(t)和x(t)作内积,即

这样

由于P0x(t)是时间t的函数,而φ(t−k)又具有低通性质,因此称P0x(t)是x(t)在V0中的“分段平滑”逼近,而称a0(k)为x(t)在V0中的“离散”逼近。它们都是x(t)在分辨率j=0时的“概貌”。

2. 子空间V1

由多分辨分析的定义,若φ(t)∈V0,则,由定理1.5,是V1中的正交归一基。仿照(1.2.10)~(1.2.12)式,有

3. 子空间W1

若在子空间W0中能找到一个带通函数ψ(t),使{ψ(t−k),k∈Z}是W0中的正交归一基,类似尺度函数φ(t),因ψ(t)∈W0,则也可构成W1中的正交归一基,即

依次类推,ψj,k(t)将是Wj中的正交归一基。称ψ(t)为小波函数,满足上述正交归一性质的正交小波的构造问题将在下一章详细讨论。这样,可依次将x(t)在Wj中作类似在Vj各空间中的分解。

D1x(t)是x(t)在子空间W1上的投影,它是时间t的函数。因为ψ(t)是带通函数,所以D1x(t)是x(t)的连续细节逼近。同理,d1(k)是x(t)在W1中的离散细节。由于V0=V1⊕W1,所以必有

P0x(t)=P1x(t)+D1x(t)

D1x(t)=P0x(t)−P1x(t)

这两个式子指出,x(t)在W1中的投影等于x(t)分别在V0和V1中的投影的差,它也是在j=0和j=1这两个分辨率水平上的逼近之差,因此,D1x(t)和d1(k)均被称为x(t)的“细节”。实际上,它们反映的也是x(t)的高频成份,且d1(k)就是尺度2ja=时的离散栅格上的小波变换。

4. 对子空间Vj,Wj

将上述的讨论加以推广,自然有如下结论:

一般,令j=1~∞,可依次实现对x(t)的多分辨率分析。

三、二尺度差分方程

前已指出,φj,k(t)是Vj中的正交归一基,ψj,k(t)是Wj中的正交归一基,并且Vj⊥Wj,Vj−1=Vj⊕Wj。在相邻尺度(如j和j1−)下的尺度函数和尺度函数之间、尺度函数和小波函数之间必然存在着一定的联系。

由于φj,0(t)=2−j/2φ(2−jt)∈Vj,而Vj−1包含在Vj中,这样,把φj,0(t)设想成是Vj−1中的一个元素,因此它当然可以表示为Vj−1中正交基的线性组合,即

式中h0(k)是加权系数,它是一个离散序列。将上式进一步展开,有

同理,由于Wj也包含在Vj−1中,因此,Wj中的ψj,0(t)也可表示为Vj−1中正交基φj−1,k(t)的线性组合,即

h1(k)也是加权系数。(1.2.15)和(1.2.16)两式被称为“二尺度差分方程”,它们揭示了在多分辨率分析中尺度函数和小波函数的相互关系。这一关系存在于任意相邻的两级之间,如j=1,有

该式又等效于

因此,二尺度差分方程是多分辨率分析中小波函数和尺度函数的一个重要性质。

由φj,k和ψj,k各自的正交性,h0(k),h1(k)可由下式求得:

,则

h0(k)=〈φ1,0(t),φ0,k(t)〉

同理

h1(k)=〈ψ1,0(t),φ0,k(t)〉

这两个式子揭示了一个重要的关系,即h0(k)和h1(k)与j无关,它对任意两个相邻级中的φ和ψ的关系都适用。这就是说,由j=0和j=1的二尺度差分方程求出的h0(k)和h1(k)适用于j取任何整数时的二尺度差分方程。h0(k)和h1(k)类似于两通道滤波器组,h0(k)对应低通滤波器H0(z),h1(k)对应高通滤波器H1(z),且在每一级,H0(z)和H1(z)保证不变。有:

对比(1.2.17)和(1.2.18)式,有

这是Haar小波及其尺度函数所对应的滤波器的系数。

对(1.2.17)式两边取傅里叶变换,有

该式是FT和DTFT的混合表示形式,式中φ(t−k)是φ(t)在b轴上离散取值所得到的,假定对b轴的抽样间隔为Ts,则上式积分的

式中ω=ΩTs,Ω是相对连续信号的角频率,Ω=2πf,而ω是相对离散信号的圆频率。后面的讨论以离散信号和离散系统为主,将Ω,ω都记为ω,并将e简记为ω。得

同理,有

Φ(ω)、Φ(2ω)和Ψ(2ω)都是连续信号的傅里叶变换(FT),而H0(ω),H1(ω)是离散信号的傅里叶变换(Discrete Time Fourier Transform,DTFT)。

将(1.2.17)和(1.2.18)两式的两边分别对t积分,由于∫φ(t)dt=1,∫ψ(t)dt=0,所以:

对应于频域,有

因此,H0(z)应是低通滤波器,H1(z)应是高通滤波器。

由(1.2.19)式,有

由于当J→∞时,,因此

式中。同理可由(1.2.20)式求出:

式中。这样,(1.2.21)和(1.2.22)式建立了H0(ω)、 H1(ω)分别和Φ(ω)、Ψ(ω)的直接关系。若H0(z),H1(z)已知,可由它们求出相应的Φ(ω)和Ψ(ω),进一步求出相应的φ(t)和ψ(t)。

四、二尺度差分方程与共轭正交滤波器组

(1.2.19)和(1.2.20)式给出了二尺度差分方程的频域关系,结合正交基的频域性质。在此基础上,可导出在二尺度差分方程中h0(k)和h1(k)的频域关系,从而把多分辨率分析和滤波器组结合起来。

定理 1.6 设φ∈L2(R),ψ∈L2(R)分别是多分辨率分析中的尺度函数和小波函数,h0(k),h1(k)分别是满足二尺度差分方程(1.2.17)和(1.2.18)式的滤波器系数,则

将(1.2.23)式写成Z变换的形式,有

发现:满足小波变换多分辨率分析中二尺度差分方程的H0(z)、H1(z)正是一对共轭正交滤波器组(CQMFB)。这样,就把小波变换和滤波器组联系了起来,从而为离散信号的小波变换的快速实现提供了有效的途径。

注意,满足(1.2.23c)式的H0(ω)和H1(ω)并不唯一,其中一个解是

,都可满足(1.2.23c)式。

(1.2.25b)式所对应的时域关系是:

h1(k)=(−1)kh0(1−k)

至此,给出了一系列重要的概念,它们分别是:

(1)在V0中总存在θ(t),使{θ(t−k),k∈Z}构成V0中的Riesz基;

(2)定理1.5说明如何由Riesz基θ(t−k)得到V0中的正交归一基,进而φj,k 是Vj中的正交归一基,即φ(t)是尺度函数。

(3)把Wj视为Vj的正交补,并假定在W0中存在小波函数ψ(t),使{ψ(t−k),k∈Z}是W0中的正交归一基,进而ψj,k是Wj中的正交归一基;

(4)由于假定Wj⊥Vj,所以假定ψj,k和φj,k是正交的;

(5)按ψj,k,φj,k分别对x(t)作分解,得到(1.2.13)~(1.2.14)式的分解(或投影)关系;

(6)由Vj⊃Vj+1,Vj+1⊕Wj+1=Vj这一包含关系,得到了(1.2.18)式的二尺度差分方程;其频域关系如(1.2.19)和(1.2.20)式所示;

(7)由定理1.6,满足二尺度差分方程的H0(z)和H1(z)恰是一对共轭正交滤波器组,即它们满足(1.2.23)式。

定理 1.7 令{Vj}是一多分辨率分析序列,φj,k(t)是Vj中的正交归一基,再令H0(z)和H1(z)是一对共轭正交滤波器组,记ψ(t)的傅里叶变换为Ψ(ω),

则存在基本小波函数ψ(t),使{ψ(t−k),k∈Z}是W0中的正交归一基,进而,{ψj,k(t),j,k∈Z}是Wj中的正交归一基。

五、Mallat算法

在上述多分辨率分析的基础上,下述两个定理给出了如何通过滤波器组实现信号的小波变换及反变换。

定理1.8 令aj(k),dj(k)是多分辨率分析中的离散逼近系数,h0(k),h1(k)是满足(1.2.17)和(1.2.18)式的二尺度差分方程的滤波器,则aj(k),dj(k)存在如下递推关系:

式中h(k)=h(−k)。

设j=0,a0(k)是x(t)在V0中由正交基φ(t−k)作分解的系数,它是在V0中对x(t)所作的离散平滑逼近。将a0(k)通过一滤波器后得到x(t)在V1中的离散平滑逼近a1(k)。该滤波器是将h0(k)先作一次翻转,得,然后a0(k)再和作卷积运算。

如果令j由0逐级增大,即得到多分辨率的逐级实现,如图1.3所示。该图所反映的过程(也即(1.2.26)式)即是Mallat算法,也即小波变换的快速实现。

图1.3 多分辨率分解的滤波器组实现

由(1.2.26)式及图1.3,可以看出Mallat多分辨率分析的思路:

(1)从滤波器组的角度看,若h0(−k)的频带在0~π/2,h1(−k)的频带在π/2~π,那么,a0(k)所处的频带是0~π,a1(k)在0~π/2,d1(k)在π/2~π;对a1(k)再分解后,a2(k)在0~π/4,而d2(k)在π/4~π/2。这就实现了对频带的逐级剖分。一方面保证了各子频带的恒 Q 性,另一方面又保证了H0(z)和H1(z)在各级的不变性;

(2)若记a0(k)所处的频带为空间V0,a1(k)处于V1,d1(k)处于W1,由它们频带的性质,显然,V0=V1⊕W1,V1⊥W1,Vj−1=Vj⊕Wj,Vj⊥Wj,j=1~∞,同时,有

V0=W1⊕W2⊕…⊕Wj⊕Vj,j∈Z

当j→∞时,aj(k),dj(k)占据的空间(也即频带)趋于无穷小,因此必有,这时的分辨率最差,因此。这就是在多分辨率分析中所讨论的主要思想;

(3)将多分辨率分解归结到图1.3来实现,这就把对离散信号的小波变换归结到逐级的线性卷积来实现。若h0(n),h1(n)的系数不是太长,卷积可在时域直接实现,否则,可用DFT及FFT来实现。

下面定理讨论了信号的重建问题,也即小波反变换。

定理1.9 若aj+1(k),dj+1(k)按(1.2.26)式得到,则aj(k)可由下式重建:

第三节 小波标架

将(1.1.2)式的连续小波变换改变成(1.1.26)式的离散小波变换,人们自然会问:一族小波函数ψj,k(t),j,k∈Z,在空间L2(R)上是否是完备的?所谓完备,是指对任一x(t)∈L2(R),它都可以由这一组函数(即ψj,k(t))来表示;如果ψj,k(t)是完备的,那么ψj,k(t)对x(t)的表示是否有信息的冗余?对a和b的抽样间隔如何选取才能保证对x(t)的表示不存在信息的冗余?

Daubechies对上述问题进行了深入的研究,给出了“小波标架”的理论[3,4,5]

一、小波标架的定义

当母小波为ψ(t)时,子小波为。对其尺度因子和平移因子进行离散,得到如下的小波标架:

如果小波标架的下上界分别为A和B,且时间采样步长为u0,a=2,每个倍频程有N个尺度因子,则可由下面的公式来估计上下界:

其中

根据上面估计的小波标架上下限,可以用迭代算法求取对偶标架:

其中要求

通过对偶标架可以重构信号:

二、计算对偶函数

当N=1时(N为倍频程的尺度因子),u0(时间采样步长)分别等于0.5、1.0、1.5时,计算墨西哥帽(Mexican_hat)小波对偶函数如图1.4所示。

图1.4 墨西哥帽(Mexican_hat)小波对偶函数

墨西哥帽的时域:,频域:

可以看出,u0越大,B/A越大,对偶函数的局域性越来越差,与原信号相差越大。

三、用小波标架重构信号

信号重构步骤如下:

(1)对信号做小波变换,得到WF(s,t);

(2)对得到的WF(s,t)进行采样,保持海森伯格盒(Heisenberg box)大小不变,得到WF(m,n);

(3)计算Morlet小波标架的上下界;

(4)利用迭代法计算每一个尺度因子s下的对偶小波;

(5)利用公式重构信号。

用小波标架重构正弦信号结果如图1.5所示,其中左上图为原始信号,右上图为墨西哥帽小波,左下图为重构出来的信号,右下图为对偶标架。

图1.5 用小波标架重构信号

四、多分辨率分析

1. 尺度函数与尺度空间

尺度函数定义:若函数φ(t)∈L2(R),由其整数平移而得到的函数序列φk(t)=φ(t−k)满足

则称φ(t)为尺度函数。

零尺度空间V0

为L2(R)的闭子空间,∀f(t)∈V0,则有

定义两个算子

其中φ(t)m,n是一个bump function,宽度为2mW,中心在2mt=n。

尺度为j的尺度空间Vj

对于j∈Z,有一系列的尺度空间{Vj}j∈Z

关于平移步长的说明:

−T≤2−jt−k≤T

−2jT≤t−2jk≤2jT

故尺度为j时平移步长为2j,即在时频平面上,时间采样间隔为2j

2. 多分辨率分析概念的引入

定义:多分辨率分析是满足下述性质的一系列闭子空间{Vj},j∈Z

(1)一致单调性:

…⊂V2⊂V1⊂V0⊂V−1⊂V−2⊂…

(2)渐近完备性:

(3)伸缩规则性:

(4)平移不变性:

(5)正交基存在性:∃φ∈V0,使得{φ(t−n)}是V0的正交基,即

3. 小波函数与小波空间

,不同尺度之间不具有正交性,故不可能成为L2(R)空间的正交基。为了寻找L2(R)空间中的一组正交基,从Vm+1在Vm中的正交补开始。

Vm称为Vm+1与Wm+1的直交和,Wm+1称为Vm+1在Vm中的正交补,记作:

Vm=Vm+1⊕Wm+1

∀fm∈Vm,有唯一的分解

fm=fm+1+dm+1,fm+1∈Vm+1,dm+1∈Wm+1

对Vm+1继续进行分解,可得

Vm=Vm+1⊕Wm+1=Vm+2⊕Wm+2⊕Wm+1=Vm+3⊕Wm+3⊕Wm+2⊕Wm+1=...

任给一个f∈L2(R),可进行如下分解:

函数f(t)在空间Vm中的投影可写为

cm,k称为低频分量,或称为光滑部分。

函数f(t)在空间Wk中的投影可写为

dk,j称为高频分量,或称为细节部分。

综合以上讨论,可得

令m→−∞,有

如果再令M→−∞,则有

4. 二尺度方程组及多分辨率滤波器组

其中:

(1.3.1),(1.3.2)式描述的是相邻两个尺度空间基函数之间的关系,称之为二尺度方程。

设H0(ω)与H1(ω)分别是h0与h1的傅里叶变换,其定义如下:

则二尺度方程(1.3.1)与(1.3.2)频域表示为:

其中,Φ(ω)与Ψ(ω)分别为φ(t)与ψ(t)的Fourier变换。

5. 滤波器组系数及其性质

(1)h0(n)和h1(n)的和式满足

(2)频域初值

H0(ω=0)=1

H1(ω=0)=0

(3)递推关系

(4)h0(n)满足偶次平移正交性

h1(n)满足偶次平移正交性

(5)滤波器H0(ω),H1(ω)的特性

以4阶Daubechies小波为例,对含有高频噪声信号进行分解和重构。使用正交小波变换的快速算法(塔式算法)。其中图1.6为源信号、重构信号和重构误差;图1.7为各阶平滑分量;图1.8为各阶细节分量。

图1.6 源信号,重构信号与重构误差

从图1.6可以看出,重构出来的信号与源信号的误差非常小(信号突变的地方误差明显大于平坦的地方),这是因为分解与重构过程中没有截断与补零的问题,重构信号与源信号只存在计算误差,误差可以忽略不计。

图1.7 平滑分量

图1.8 细节分量

五、实验

实验1:Gabor标架

Matlab程序:

实验2:小波标架

实验3:多分辨率分析

第四节 小波的提升构造

利用提升方案(Lifting Scheme)构造小波[4,5,6]

一、提升方案的基本原理

小波函数ψj,k(t)通常定义为一个属于L2(R)空间的母小波的二进伸缩(Dilates)和平移(Translate):

ψj,k(t)=2j/2ψ(2jt−k)

这样的小波称为第一代小波。在更一般的情况下,小波并不必须是彼此的伸缩与平移,但仍然具有第一代小波的特点,这样的小波称为第二代小波,利用提升方案可以构造它们。

第一代小波具有如下性质:

P1:是L2(R)空间的Riesz基,还是Lebesgue、Lipschitz、Sobolev和Besov空间的无条件基。

P2:小波及其对偶在空间和频域是局域化的,有些小波还是紧支的。

P3:小波分析可纳入多分辨分析的框架,这导致了快速小波变换算法。

在研究中常有如下需要:

G1:第一代小波提供了定义在nR上函数的基,但在如数据分割、在一般定义域上的微分和积分方程的求解,需要定义在任意的、可能不光滑的域上的小波。

G2:第一代小波典型地只提供具有不变测度的空间的基,而微分方程的对角化、在曲线或表面上的分析等需要可适应加权测度的基。

G3:第一代小波隐含对数据进行规则采样,而实际问题经常要处理不规则采样的数据。

具有性质P1~P3而又满足G1~G3性质的第一代小波的推广称为第二代小波。这里的关键问题是平移与伸缩并不是属性P1~P3所必需的,放弃平移和伸缩,隐含着傅里叶变换不能再用作构造工具。下面介绍利用提升方案构造第二代小波的方法。

考虑信号X={xk|xk∈R}k∈Z,把 X 分成二个不相交的集合:偶下标采样Xe={x2k}k∈Z和奇下标采样Xo={x2k+1}k∈Z,通常情况下这两个集合是紧密相关的,因而从一个集合能很好地建立另一个集合的预测P:

d=xe−P(xo)

知道了d和奇采样值,可立即恢复信号

xe=d+P(xo)

若P性能好,则d将是一个稀疏集。

利用相邻两偶采样对奇采样进行预测,记下差值

若信号是相关的,则大多数小波系数γ−1,k将很小。在理论上,可以继续通过对{λ−1,k}k∈Z施加以上操作,然而,上述简单的操作性能并不好,为此引入另一个条件,即希望λj,k系数的平均值在每一次分解时保持一致,或者说使,此前所进行的下采样很显然不具有这种特点,可通过借助于γ−1,k对λ−1,k进行提升来实现这点:

现在,每一级小波变换由两步构成:首先计算小波系数,其次提升下采样系数。逆变换可立即得到:只需把式(1.4.2)中的加号换成减号,再把式(1.4.1)的等式中的项作一下移动即可。在进行小波变换时,可进行同址运算,即不需要辅助存储器,这对硬件实现十分有利。

下面的定理给出提升方案的一般方法。

定理1.10 给定双正交滤波器算子的初始集合,那么可通过如下方法获得一个新的双正交滤波器算子集

式中Sj是一个从l2(M(j))到l2(K(j))的算子。

二、把小波变换分解成基本的提升步骤

已经证明所有 FIR 小波滤波器都有能分解成基本的提升步骤。用矩阵表示时,一个提升步骤对应一个单元(elementary)矩阵。分解的基本理论依据是矩阵代数,根据矩阵代数,任何具有多项式元素项且行列式为1的矩阵都可以分解成一系列的单元矩阵。

把求自然数的最大公约数的Euclidean算法推广到求两个多项式的最大公因子。两个多项式的公因子取决于因子pz,在多项式的情形下,解并不是唯一的。

定理1.11 (多项式的Euclidean 算法)。设有两个多项式a(z)和b(z)≠0,而且。令a0(z)=a(z),b0(z)=b(z),从 i = 0开始循环执行以下步骤

那么an(z)=gcd(a(z),b(z)),如果n是使得bn(z)=0的最小数。

定理中定义为:若,则

用矩阵形式表示为

相应地

式中n≤b(z)+1。这样,an(z)能整除a(z)、b(z),如果an(z)是一个单项式的话,那么a(z),b(z)是互素的。

为了把FIR小波滤波器(h, g)分解成基本的提升步骤,首先注意到he(z),ho(z)必须是互素的,而且,利用公约数的不唯一性,总是可以使公约数为常量K,即

对于给定的滤波器h,通过如下操作,总可以找到一个互补滤波器0g,

即令

在式(1.4.3)中

当i为奇时使用式(1.4.4)的第一个等式,当i为偶时使用第二个等式,有

通过一个提升步骤可获得滤波器g,

定理1.12 给定互补滤波器对(h, g),那么总是存在多项式si(z)和ti(z),1≤i≤n,以及一个常量K,使得

与对偶滤波器对相关的多相(polyphase)矩阵为

在正交小波滤波器时有P(z)=P~(z),这就对应着两种不同的分解,也就是说,把FIR小波滤波器分解成基本的提升步骤时,分解是不唯一的。

作为例子,下面给出对具有两阶消失矩的D4正交小波的分解:

其中

多相矩阵是

因式分解是

使用式(1.4.5)作为P(z) 的分解,则分析用的多相矩阵为

由此可得小波分解算法

由分解算法,通过反向进行操作,并改变相应的符号可得重构算法

利用提升方案进行小波变换具有可进行同址运算优点,这样在具体实现时可省去大量的存储器开销,在进行图像处理时,这个优点更为明显。它的另一个优点是可提高小波变换的速度。所以把现存的有限长小波滤波器分解成基本的提升步骤,可加快小波变换的进行,根据Daubechies的分析,随滤波器长度的增加,运算速度趋于常规小波变换的2倍,换言之,在同等的硬件条件下,对一维小波变换而言,运算时间降低一半,对二维小波变换则降为原来的四分之一。这个优点在实时性要求比较高的场合有很大的实用价值。

三、二代小波漫谈

现在举例,对一个8点序列,怎样实现第二代小波变换。

1. 奇偶分开

非常简单,就是[2,4,6,8]组成一列向量,[1,3,5,7]组成一列向量。

2. 预测

用[2,4,6,8]来预测[1,3,5,7]。比如说1, 3估计2;3, 5估计4;5, 7估计6;7,1估计8。(边缘处理,采用循环方法)。估计公式可以用别人的,也可以自己做。举一个线性的例子:2=1*a+3*b,4=3*a+5*b,…,其他的都一样。这样就可找到最优的a,b,使得(2−(1*a+3*b)).^2+(4−(3*a+5*b)).^2+…最小化。就是最小均方准则。若正好为零,说明偶可以完全预测奇,也就是只要存储偶数列向量a和b就可以了,压缩也就是实现了。对于信号很长序列,就等于压缩了一半。当然,可以采用更复杂的立方差值预测,多项式预测,或其他的准则来使其最小,这样的压缩也就得到了最优。

3. 提升

经过预测后,要存储的是偶数序列[2,4,6,8],新的奇数序列[n1,n3,n5,n7]= [2−(1*a+3*b),4−(3*a+5*b),…]和线性变换系数(a,b)。这里新的奇数序列就是高频分量。但偶数序列是不能完全代表信号的性质的,有所差距。所以要对偶数序列进行修正。即所谓的提升。这次用个简单的提升吧。[n2,n4,n6,n8]= [2,4,6,8]+k*[n1,n3,n5,n7]。[n2,n4,n6,n8]就是要分解的低频分量。那k怎么求呢?因为要保持n2,n4,n6,n8和原始信号 [1,2,3,4,5,6,7,8]一样的性质。一般就是均值和高阶矩。这里就一个未知数k,所以用均值相等就行了。1/8*(1+2+3+…+8)= 1/4*(n2+n4+n6+n8)。k 很容易就求出来了。最终存储的就是[n1,n3,n5,n7]和[n2,n4,n6,n8]以及a,b,k。

最后,对二代小波总结如下:

第一,反变换,就是3−>2−>1。

第二,处理二维图像:先行提升,再列提升。

第三,整数阶变换:就是加一个取整。

第四,多层或小波包提升,就是在对序列[n1,n3,n5,n7]或[n2,n4,n6,n8]再做1−>2−>3。

第五,灵活。不一定是a, b,也可能就一个a或a,b,c;不一定是一个k,也可能是k1, k2。但越多计算量太大。最好是用大师们做好的CDF,5/3,7/9等。

第六,最重要的是,任何一代小波,总可以通过一次或多次提升实现,它和一代小波没有本质区别。

四、整数小波变换

提升方案为扩展小波变换的应用领域提供了更多的灵活性。常规的小波变换都是采用浮点运算的,但利用提升方案所带来的便利,可十分方便地构造整数到整数的小波变换。将整数小波变换用于图像压缩就可以用小波变换进行无失真的图像压缩。

最早的整数到整数小波变换是S变换,是哈尔(Haar)变换的整数形式:

Said和Pearlman提出了S+P(tranSform + Prediction),就是在S变换之后,利

用低通滤波器的系数来产生一个新的高通滤波器系数,它的一般形式是:

式中α−1=0,α0=2/8,α1=3/8,β1=−2/8。

S和S+P变换的逆变换由式(1.4.6)、(1.4.7)把执行顺序变动一下,再改变一下相关项的符号就可以获得。

通常小波滤波器的系数都是浮点数,只能把整数映射成浮点数,要进行无失真变换,必须构造把整数映射成整数的小波变换,所有正交或双正交小波滤波器,用提升方案进行分解后,都可用与S+P类似的方式来构造变换的整数版本。

用提升方案构造小波变换有如下优点:

(1)同址计算:即不需要辅助存储器,源信号(图像)可被小波变换的结果覆盖。

(2)更快的小波变换:传统上,快速小波变换首先把信号分解成高通和低通成分,并进行下抽样,然后对低通成分重复进行该过程直到所需要的变换级数。提升方案可把变换速度提高1倍。

(3)不需借助傅氏分析便可获得逆变换。实际上,只要简单地调整一下正变换中的正负号即可。此优点使得不需要很强的傅氏分析的背景便可理解小波的特性和小波变换。

S−变换的逆变换非常容易得到:

免责声明:以上内容源自网络,版权归原作者所有,如有侵犯您的原创版权请告知,我们将尽快删除相关内容。

我要反馈