现代信号处理教程（第三版）课件 第11章 正交与双正交小波构造、正交小波包

11.1正交小波概述11.2由递推求解的方法。11.3消失矩、规则性及支撑范围11.4Daubechies正交小波构造11.5接近于对称的正交小波及Coiflet小波11.6双正交滤波器组11.7双正交小波11.8双正交小波的构造11.9双正交样条小波11.10正交小波包第11章正交、双正交小波构造及正交小波包11.1正交小波概述

1.Haar小波Haar小波在时域的正交性是明显的正交、双正交；时域有限支撑，频域无限扩展2.Shannon小波Shannon小波的尺度函数，

中的正交基。如何求出Shannon小波？Shannon小波Haar和Shannon小波是正交小波中两个极端的例子。

Shannon小波及其尺度函数的频谱如图所示，显然，它们在频域有着明显的正交性。?3.Daubechies小波设计出更通用、时域和频域都具有紧支撑的正交小波设计思路：对要设计的小波提出要求；将这些要求转化为对滤波器的要求；设计出对应该小波的共轭滤波器组；由滤波器得到尺度函数和小波函数。关键是找到小波和滤波器之间的关系小波和滤波器之间的关系：二尺度差分方程时域关系CQMFB先讨论如何由滤波器递推求出尺度函数和小波函数；再讨论如何对小波提出各种要求。二尺度差分方程频域关系如何实现?目的:

11.2由递推求解的方法近似因为未知，所以该式无法利用考虑频率关系每两点插入一个零每两点插入个零依次类推点序列点序列点序列卷积的结果，序列的长度急剧增加。最后要压缩到一个有限长的区间中。令

例当趋近于无穷时，逼近，

“逼近”连续函数

这时需要将连成一个阶梯形的曲线，然后压缩成一条连续曲线。压缩办法，见“wavefun.m”在中，由于不知道，所以无法递推，因而采用的上述的递推方法。假定取某一简单的形式，可用另外的方法递推出所需要的。令

的起始值由递推求解的又一方法：假定的长度为为Haar函数，则：现研究一下支撑范围:即的支撑范围是0～N，也是紧支撑的。继续递推，有现举例说明上述递推过程。例令

其他求递推过程：令

例令

例

但是，并不是上述的递推都可以收敛到一个尺度函数，例如：为什么？11.3消失矩、规则性及支撑范围

定义：1．消失矩(Vanishingmoments)

的阶矩如果在处有阶重零点，即由傅里叶变换的性质，有？则：显然即是容许条件且则说该小波函数具有阶消失矩。上述说法具有阶消失矩。关于消失矩的解释：令并假定小波变换恒为零

由有若可展成一高阶的多项式(如用台劳级数)，如 阶，。那么其中阶次小于的多项式部分(对应低频)在小波变换中的贡献恒为零，反映在小波变换中,有的只是阶次大于

的多项式部分，它们对应高频端，这就有利于突出信号中的高频成分及信号中的突变点。从这个角度讲，我们希望能具有尽量高的消失矩。消失矩越高，在处越平滑地为零，越具有好的带通性质。构造小波要考虑的因数之一若的值大部分为零,有利数据压缩；去除噪声；突出信号中的奇异点。取决于信号自身的特点；小波是否紧支撑；小波是否具有高的消失矩。因为有关于例

以上关系为我们设计具有高阶消失矩的小波提供了一个切实可行的方法。下面的定理进一步明确了有关消失矩的几个相关概念。

令保证了在处有

阶重零点。所以，又保证了在处应有

阶重零点。都具有阶消失矩因此因为定理11.1(1)小波有阶消失矩如果在处阶连续可微，则下面三个说法是等效的：(2)和它的前阶导数在处等于0(3)和它的前阶导数在处等于0证明：两边微分，即有（3）和（2）等效所以（1）和（2）等效所以，如果在处有阶重零点，则上述结论也适用于离散信号：

具有阶消失矩

数学上用于描述函数局部特征；在信号处理中用于描述信号在某点，或某个区域内的平滑性和奇异性。规则性常用Lipschitz指数来度量。有关Lipschitz指数的定义及在信号检测中的应用将在第13章讨论。2．规则性(regularity)

粗略的说，函数在某一点，如的导数的阶次越高，函数在该点越平滑，其Lipschitz指数也越高。在区间

上都有高阶的导数，该函数在整个区间上是平滑的。用这样的函数去分析信号，更有利于突出信号的变化。几个特殊的情况(记李氏指数为）：若在处一次可微，但一阶导数不连续（分段线性函数），则在拐点处的阶跃函数在阶跃点是连续的，但不可微，所以单位冲激函数在处，既不连续，更不可导，所以因此，若想构造出高规则性的尺度函数和小波函数，还要从滤波器入手。讨论规则性的目的是构造出连续且具有高阶导数的尺度函数和小波函数。由上例可知：由递推求出的收敛，即平滑递推求出的不收敛，不平滑Daubechies给出了正交小波的构造方法，还给出了保证具有高的消失矩和规则性的方法。假定：

阶

阶保证具有 阶消失矩

函数，随衰减，大衰减快可导出1.如果Daubechies证明了：可保证递推出的

是收敛的2.再如果

不但是收敛的，且是次可微的是设计小波滤波器的基本结构。合适的选择和可得到不同性能的滤波器及对应的小波函数。例所以，递推出的收敛。例所以，递推出的不收敛。

如果是有限支撑的，则尺度函数和小波函数均是有限支撑的。若的支撑范围是则的支撑范围也是，而的支撑范围是：3.支撑范围定理11.2证明见教材11.4Daubechies正交小波构造思路：由要求：：正交、紧支撑：高阶消失矩

已给定：：好的规则性，即尽可能平滑对小波函数和尺度函数的要求转化为对滤波器的要求：1.低通特性的要求

在要有零点，且零点的阶次（）越高越好;2.

3.要保证规则性要求，越大越好若能确定，则 也就确定目标：如何求出具有最小阶次且满足上述 要求的多项式，并使得

满足CQMFB:阶次为:阶次为Daubechies证明了：满足要求的的最小阶次为：这样，的阶次为；即有个系数。Daubechies算法的导出：

的偶函数，可表为的函数

也是的偶函数，故也可表为的函数因为所以则令则再令

由有任务：求Bezout方程Bezout方程

BezoutEtienne（贝祖）（1730－1783）

法国数学家Bezout定理指出:若，为阶次，的多项式，且二者之间没有共同的零点，那么，唯一地存在两个阶次分别为,的多项式,使得则令无共同零点由

是一奇对称多项式，即的不同选择可构造出不同类型的小波，在构造正交小波时，Daubechies选择因为可求出：因此有：给定,右边多项式可求。做谱分解，可将单位圆内的零点赋给,则 是最小相位的。这样可以求出令，求db小波及尺度函数例由由求出：压缩到的范围内，得：由滤波器系数

由FB的关系，求出又可求出Haar小波令，求db小波及尺度函数

例求出：db小波在时及的波形

11.5接近于对称的正交小波 及Coiflet小波

问题

db小波是紧支撑的正交小波，但它们不是对称的,也即和不具有线性相位。

SymN小波和db小波构造的方法基本相同

使序列较为对称的办法是令为混合相位系统，即其零点有的在单位圆内，有的在单位圆外。

令，求db小波及尺度函数

例该多项式六个零点：db小波：sym4小波：

利用MATLAB的有关文件可给出N＝2～10时的滤波器的系数及相应的尺度函数、小波函数。关于coifN”小波的构造在“dbN”中，小波函数具有最大的消失矩()，但尺度函数的消失矩没有考虑。Coifman感到这一类小波在数值分析中的应用中尚不能满足要求，因此要求Daubechies构造一类使也具有和一样消失矩的正交紧支撑小波。Daubechies构造了这样的小波，并称之为“coiflet小波”。特点：对应频域等效Daubechies给出了时的求解方法

保证满足功率互补关系coifN小波和dbN小波相比，支撑范围变宽，但具有更好的对称性。正交小波不足：

和都不是对称的;symN,coifN接近对称，但不是真正的对称。

正交小波的优点：

尺度函数和小波函数都是紧支撑的，有着高的消失矩

和的不对称来自于和的不对称。因此，要设计出对称的尺度函数和小波函数就要首先设计具有对称的滤波器组。我们已在7.8节讨论了具有线性相位的双正交滤波器组的基本概念，给出了可准确重建的双正交滤波器组的设计方法。

将上述内容引入到小波分析，即可给出适合小波变换的双正交滤波器组准确重建的条件，同时也给出双正交条件下的多分辨率分析及双正交小波的构造方法。

双正交滤波器组一对滤波器用于信号的分解，一对滤波器用于信号的重建（对偶滤波器）两个滤波器的长度不再相同。11.6 双正交滤波器组准确重建的充要条件:普通两通道FBPR条件

双正交FBPR条件区别翻转来自不同的滤波器双正交滤波器组中的正交关系：同一支路滤波器系数偶数移位正交上下支路，滤波器系数偶数移位交叉正交注意两组正交的不同上下支路交叉正交:

正交于;正交于上下支路各自是正交的:

和其对偶正交;和其对偶正交11.7 双正交小波

离散小波变换:

分解重建对偶关系正交小波：二者相等正交小波二尺度差分方程频域关系在双正交滤波器组的情况下，分解滤波器(

)和重建滤波器(

)将产生两个尺度函数()和两个小波函数()。其中和对应信号的分解，而和对应信号的重建。

引出如下关系双正交小波：两尺度差分方程频域关系要求低通高通低通带通满足要求两套函数分解重建双正交时：不要求互相正交

但要求和其对偶正交

放宽正交性是使具有线性相位，从而使 更具有对称性，从而减小了相位失真。双正交情况下的多分辨率分析：

正交情况下：

双正交情况下：不再要求！两套空间互为对偶要求双正交！正交性在频谱上的反映双正交小波变换的快速算法和正交小波变换的快速算法基本相同，区别是在重建时使用的是对偶滤波器。两对滤波器的频率特性尺度、小波函数和其对偶函数的频率特性

分解出的概貌分解出的细节双正交情况下的多分辨率分解：双正交情况下的多分辨率重建：二插值在双正交小波中，分解和重建对应两套滤波器、尺度和小波函数，MATLAB中也是这样给出的。例如[RF,DF]=BIORWAVF(W)returnstwoscalingfiltersassociatedwithbiorthogonalwaveletspecifiedbythestringW.W='biorNr.Nd'wherepossiblevaluesforNrandNdare:Nr=1Nd=1,3or5；Nr=2Nd=2,4,6or8Nr=3Nd=1,3,5,7or9；Nr=4Nd=4Nr=5Nd=5；Nr=6Nd=8RFisthereconstructionfilterDFisthedecompositionfilter[LO_D,HI_D,LO_R,HI_R]=BIORFILT(DF,RF)computesfourfiltersassociatedwithbiorthogonalwaveletspecifiedbydecompositionfilterDFandreconstructionfilterRF.LO_DDecompositionlow-passfilterHI_DDecompositionhigh-passfilterLO_RReconstructionlow-passfilterHI_RReconstructionhigh-passfilter[LO_D1,HI_D1,LO_R1,HI_R1,LO_D2,HI_D2,LO_R2,HI_R2]=BIORFILT(DF,RF,'8')returnseightfilters,thefirstfourassociatedwiththedecompositionwaveletandthelastfourassociatedwiththereconstructionwavelet.Forabiorthogonalwavelet:[PHI1,PSI1,PHI2,PSI2,XVAL]=WAVEFUN('wname',ITER)returnsthescalingandwaveletfunctionsbothfordecomposition(PHI1,PSI1)andforreconstruction(PHI2,PSI2).11.8双正交小波的构造

双正交小波构造关键是的构造：1.支撑范围：在构造时，应尽量让它们在处有高阶的重零点，以保证2.消失矩：有阶消失矩；技术要求：

有着相同的规则性；的规则性随着在处零点数的增加而增加；的规则性也随着在处零点数的增加而增加；如果在处有不同的零点数，则的规则性也不相同。

3.规则性：4.对称性：要求：具有线性相位；结果：具有对称性；

为对称中心为对称；以位移中心反对称；

为对称中心对称方式为对称；以位移中心为对称；

及的构造线性相位要求滤波器系数是实的如果以为对称中心，则滤波器的系数仍是实的，增加一个移位因子为保证、具有阶消失矩，令类似正交小波构造现在的任务是寻找，使其满足：Cohen，Daubechies双正交小波的结构方法：

因为和都是的实系数多项式，所以可将它们分别写成的函数：（A）

情况1：为对称；当以为对称中心，即

情况2：（B）将（A）、（B）二式分别代入：:对应(A)式:对应(B)式有：可表为的函数，即同理再令：得到：令：Bezout方程

解得：以对称以对称问题1:如何得到？问题2：如何和

联系起来？在双正交小波中，有多少函数需要构造注意下面的符号？11.9 双正交样条小波

样条函数是分段光滑且在连结点处具有一定光滑性的一类函数，它在数值逼近方面有广泛的应用。其中基数B样条(CardinalB-Spline)函数具有最小的支撑范围且又容易在计算机上实现，因此被认为是构造小波函数的最佳候选。令：次卷积

：不连续：连续，一阶导数不连续：连续，一阶导数也连续随着的增大，越来越光滑用样条函数构造双正交小波的思路：有多少函数需要构造？注意下面的符号Step1

FT这样，在双正交小波中，一个尺度函数基本上不需求解，只需把相应阶次的样条函数作为该尺度函数即可。Step2令有（A）（A）Step3由求出再实现下面的分解：现在的问题：如何得到？综合Step1,Step2,对应（A）式，有于是可以求出（B）同理，对应（B）式,有：

：仅和有关，和无关，取决于

：和，都有关，取决于和(A)(B)给定不同的和，即可求出：例令由必有选中（B）式。由注意：双正交的情况下，分解滤波器和综合滤波器的系数一般不等长。现：实际尺度为2；可取不同的值，因此有不同的长度令由必有已有由其他Haar函数是正交的，同时也是双正交的MATLAB中记法令由必有已有再令必有已有可以求出：注意：和系数的长度不同，随着的增大，二者的差距急剧变大。例令由必有选中（A）式由有由例令必有求出令必有求出令必有求出：分解和重建滤波器是对称的，都具有线性相位；4、阶是阶（）样条函数；可求出：5、按照频域关系，有3、滤波器系数，在分母是2的整次幂；2、滤波器和其对偶滤波器的长度有较大的差别；双正交小波及滤波器组的特点：用于分解的用于分解的用于分解的上述双正交滤波器的改进：因为所以，的长度即是的长度。信号分解时：

两个滤波器的长度有较大的差别，一方面影响滤波器的性能，另一方面也使输出数据的长度有较大的差别，这是不希望的。如何解决长度差别太大的问题？长度差别太大的原因是由于分解：改变的方法是令：Daubechies的方法是：

的一阶实根

的共轭复根在保证、系数始终为实数的情况下，考虑对和的分配。例按原来的方法，即令六阶多项式比高六阶现在它们各由两部分组成令共同部分：长度等于5长度等于7长度等于9不同部分将零点“db9/7小波”长度接近相等！在已求出：的情况下，“db5/3小波”用于分解的

及用于重建的

例有两种分解方法，其中之一是：用于分解的

及用于重建的

11.10 正交小波包

Mallat多分辨率分析：特点：只对低频空间进行分解，对高频空间不分解。空间中给出离散系数，它反映了信号的概貌，空间中给出离散系数，它反映了信号的细节，也即小波系数。显然比有着好的时间分辨率，但差的频率分辨率。若在二者之间取折中，则可将也作分解。

空间的逐级分解概念：

在上面的分解中，任取一组空间进行组合，如果这组空间：①能将空间

覆盖；②相互之间不重合;则称这一组空间中的正交归一基的集合构造了一个小波包(waveletpacket)。都符合要求上述分解的实现：问题：1.在每一子空间，正交基是否存在？如何构造？2.如何选择最佳小波包？为了探讨小波包中正交基的存在及构造方法，再一次考察Mallat算法中正交基的行为：(1)

信号本身的性质；(2)

信号分解的目的；(3)“最佳”原则的选择。一个最佳小波包的选择取决于三个因素：二尺度差分方程CQMFB√?定理12.5

令是空间中的正交归一基， 是一对共轭正交滤波器，令则是中的正交归一基。该定理说明：这样，在对和的向下逐级分解中，每一个子空间中都存在着按二尺度差分方程生成的正交基第级，有个子空间统一命名：对应低频,正交基：对应高频,正交基小波包中的基本关系子空间总数的一半关于空间和基向量名称的说明：如此类推，可得到各个子空间的正交基(A)(B)(C)例Haar小波由（A）式：由（B）,(C)式：例令中的为db5对应的尺度函数，当时，求出中的正交基上面两例可以表示为下面的树状结构：对小波包的节点，记令则作为多分辨率分析的输入为在该结点(或子空间)处的小波包系数，下面的定理给出了小波包分解和重建算法：在小波包的分解中，在结点处的小波包系数由下式给出定理12.6在结点处的小波包系数由下式重建插入一个零后的新序列信号的小波包分解与重建信号是MATLAB中所给的信号noisdopp.mat，下图是其小波包分解：

例关于最佳分解“树”的选择：一个最佳的小波包应使信号在其各个子空间中的投影(即)尽可能地大，这取决于信号能量随频率的分布。

在决定“最佳小波包”的过程中，需要一个“代价函数”，从而使在各种小波包选择的可能中，选择一种具有最小代价的小波包。去噪：选择将信号和噪声最大程度分离的空间；数据压缩：选择信号能量最大集中的空间。

求出信号在每一个子空间的“熵”，或“范数”，或指定“阈值”，并把它们作为代价函数来决定小波包的选择。令为信号在某一子空间正交基上的投影。定义：waveletToolbox中有着丰富的有关小波包的命令：为的Shannon熵

为的范数为的对数能量熵MATLAB中有关小波包的文件：WaveletsPacketsAlgorithms.bestlevt-Bestleveltree(waveletpacket).besttree-Besttree(waveletpacket).entrupd-Entropyupdate(waveletpacket).wentropy-Entropy(waveletpacket).wp2wtree-Extractwavelettreefromwavelet packettree.wpcoef-Waveletpacketcoefficients.wpcutree-Cutwaveletpackettree.wpdec-Waveletpacketdecomposition1-D.wpdec2-Waveletpacketdecomposition2-D.wpfun-Waveletpacketfunctions.wpjoin-Recomposewaveletpacket.wprcoef-Reconstructwaveletpacketcoefficients.wprec-Waveletpacketreconstruction1-D.wprec2-Waveletpacketreconstruction2-D.wpsplt-Split(decompose)waveletpacket.function[wpt,ento,n2m]=bestlevt(wpt)%BESTLEVT：Bestlevelwaveletpackettree.%BESTLEVTcomputestheoptimalcompletesub-treeofan%initialtreewithrespecttoanentropytypecriterion.%Theresultingcompletetreemaybeofsmallerdepth%thantheinitialone.%T=BESTLEVT(T)computesthemodifiedtreeT%correspondingtothebestleveltreedecomposition.[T,E]=BESTLEVT(T)returnsbesttreeT%andinaddition,thebestentropyvalueE.%Theoptimalentropyofthenodewhoseindexisj-1isE(j).function[wpt,ento,n2m]=besttree(wpt)%BESTTREEBestwaveletpackettree.%BESTTREEcomputestheoptimalsub-treeofaninitialtree%withrespecttoanentropytypecriterion.%Theresultingtreemaybemuchsmallerthantheinitialone.%T=BESTTREE(T)computesthemodifiedtreeT%correspondingtothebestentropyvalue.%[T,E]=BESTTREE(T)returnsthebesttreeT%andinaddition,thebestentropyvalueE.%Theoptimalentropyofthenodewhoseindexisj-1functiont=entrupd(t,ent,in3)%ENTRUPDEntropyupdate(waveletpackettree).%T=ENTRUPD(T,ENT)or％T=ENTRUPD(T,ENT,PAR)%updatestheentropyofwaveletpackettreeT%usingtheentropyfunctionENTwithoptional%parameterPAR(seeWENTROPYformoreinformation).functionent=wentropy(x,t_ent,in3)%E=WENTROPY(X,T,P)returnstheentropyEoftheX.%Tisastringcontainingthetypeofentropy:%T='shannon','threshold','norm',%'logenergy'(or'logenergy'),'sure','user'.%PisanoptionalparameterdependingonthevalueofT:%IfT='shannon'or'logenergy',Pisnotused.%IfT='threshold'or'sure',Pisthethreshold，apositivenumber.%IfT='norm',Pisthepowerandmustbesuchthat1<=P.%IfT='user',PisastringcontainingtheM-filenameofyour％ ownentropyfunction,withasingleinputX.function[wpt,n2m]=wp2wtree(wpt)%WP2WTREEExtractwavelettreefromwaveletpackettree.%T=WP2WTREE(T)computesthemodifiedtreeT%correspondingtothewaveletdecompositiontree.functionx=wpcoef(wpt,node)%WPCOEFWaveletpacketcoefficients.%X=WPCOEF(T,N)returnsthecoefficientsassociated%withthenodeNofthewaveletpackettreeT.%IfNdoesn'texist,X=[];function[wpt,n2m]=wpcutree(wpt,level)%WPCUTREECutwaveletpackettree.%T=WPCUTREE(T,L)cutsthetreeTatlevelL.%Inaddition,[T,RN]=WPCUTREE(T,L)returnsthevectorRN%whichcontainstheindicesofthereconstructednodes.functionx=wprcoef(t,node)%WPRCOEFReconstructwaveletpacketcoefficients.%X=WPRCOEF(T,N)computesreconstructedcoefficients%ofthenodeNofthewaveletpackettreeT.functionvarargout=wpdec(x,depth,wname,type_ent,parameter)%WPDECWaveletpacketdecomposition1-D.%T=WPDEC(X,N,'wname',E,P)returnsawptreeobjectT%correspondingtoawaveletpacketdecompositionofthevect