多分辨率分析:离散小波变换为什么需要离散小波变换虽然离散化的连续小波变换(即小波级数)使得连续小波变换的运算可以用计算机来实现,但这还不是真正的离散变换。事实上,小波级数仅仅是CWT的采样形式。即便是考虑到信号的重构,小波级数所包含的信息也是高度冗余的。这些冗余的信息同样会占用巨大的计算时间和资源。而离散小波变换(DWT)则不仅提供了信号分析和重构所需的足够信息,其运算量也大为减少。 相比CWT,DWT的实现要容易得多。本小节将介绍DWT的基本概念及其性质,以及用来实现其计算的算法。如前面的内容一样,会举一些应用实例来帮助理解DWT。 离散小波变换(DWT)历史DWT的建立要追溯到1976年。当时,Croiser, Esteban, 和 Galand发明了一种分解离散时间信号的新技术。几乎在同时,Crochiere, Weber, 和 Flanagan在语音信号编码上也做了类似的工作。他们将其命名为子带编码。1983年,Burt定义了一种与子带编码非常类似的新方法,并取名为金字塔编码。现在,这两种编码方法都又称为多分辨分析。到1989年,Vetterli 和 Le Gall对子带编码方法进行了一些改进,并且去除了金字塔编码中的冗余。下面将会简要介绍子带编码。离散小波变换及多分辨分析理论的详细讨论可在很多相关的论文及专著中找到,这里不详细展开。 子带编码和多分辨分析主要思想与CWT完全一样。利用数字滤波技术,可以得到数字信号的时间-尺度表示。回忆一下前面关于CWT的讨论中,CWT可以看做是不同尺度的小波和信号的相关运算,其结果表征了信号与不同尺度小波的相似程度。改变分析窗的尺度,在时间上移动窗函数,并将其与信号相乘来计算连续小波变换。在离散的情况下,用不同截止频率的滤波器在不同的尺度上对信号进行分析。信号通过一组高通滤波器来分析其高频分量,通过一组低通滤波器分析其低通分量。信号的分辨率这个信号详细信息的表征方法,会随着滤波方式而变化。尺度通过上采样(upsampling)和下采样(downsampling)而变化。下采样一个信号指的是降低信号采样率,或者说去除一些信号中的采样点。例如,2倍的下采样指的是信号每隔一个采样点保留一个值。n倍的下采样则是将采样点降为1/n。 上采样一个信号指的则是通过对原始信号增加采样点数来提高信号采样率。例如,2倍的上采样指的是在每两个相邻的采样点之间增加一个新的采样点,其值通常为0。n倍的上采样则是在通过增加n倍的采样点数来将采样率提高n倍。 虽然不是唯一的选择,但DWT的系数通常是以二进栅格的方式从CWT采样得到的。也即是说,s0 = 2 ,tau0 = 1,从而有s=2^j ,tau = k*2^j,这点在前面已经有所讨论。由于信号是离散时间函数,在下面的讨论中将不再区分函数与序列这两个概念。序列通常用x[n]表示,其中n为整数。 处理过程从信号通过一个冲激响应为h[n]的半带数字低通滤波器开始。对信号滤波的过程在数学上等效为信号与滤波器冲激响应的卷积。离散时间的卷积定义如下:
半带低通滤波器滤除了信号中超过最高频率一半的所有频率分量。例如,如果信号的最高频率分量为1000Hz,则半带低通滤波器滤除所有超过500Hz的频率分量。【译者注:这里可能有点小问题,应该是若滤波器采样频率为1000Hz,则信号中超过500Hz的频率分量全部被滤除。】 这里要特别注意频率单位的问题。对离散信号来说,频率的单位通常用弧度(radian)来表示。相应地,根据圆周频率,信号的采样频率为2pi弧度。于是,如果信号以奈奎斯特频率(信号最高频率分量的2倍)采样,则经半带低通滤波器滤波后得到的信号其最高频率为pi弧度。这也就是说,在数字频率域,奈奎斯特频率对应为pi rad/s。因此,对离散信号用Hz为单位便不再合适。但是,人们在讨论中已经习惯用Hz来表示频率的单位。但要牢记,离散时间信号的频率单位是弧度。【译者注:满足采样定理的情况下,采样频率对应为2pi弧度,奈奎斯特频率对应pi弧度】 信号通过半带低通滤波器之后,根据奈奎斯特采样规则,由于信号现在的最高频率为 pi/2弧度而不是 pi弧度,因此可以扔掉一半的采样点。只要每隔一个点扔掉一个点就可实现对原信号2倍的下采样。这样信号长度变为原来的一半。信号的尺度变为原来的2倍。需要说明的是,低通滤波器仅是滤除了高频信息,但并没有改变尺度。只有下采样过程才改变了尺度分辨率。但从另外一方面来说,下采样又与信号所包含的信息有关,受滤波运算的影响。半带低通滤波器滤除了一半的频率信息,这可解释为丢失了一半的信息。因此,滤波处理以后,分辨率降为一半。需要说明的是,滤波之后的下采样处理并不影响分辨率,因为从信号中滤除一半的频率分量使得信号中一半的采样点为冗余数据。丢掉一半的采样点不会丢失任何信息。小结一下,低通滤波将分辨率降低一半,但是。不能改变尺度。然后再进行2倍的下采样处理,扔掉一半的冗余数据 上述过程用数学公式可以表示如下:
说到这里,我们来看一下DWT到底是如何计算的:DWT通过在不同的频段利用不同的分辨率来将信号分解为粗略的近似和详细信息。DWT用到了两组函数,即尺度函数和小波函数,分别对应了低通滤波器和高通滤波器。将信号分解到不同的频段可简单地在时域将信号不断通过高通和低通滤波来实现。原始信号x[n]首先通过一个半带高通滤波器g[n]和一个半带低通滤波器h[n]。滤波之后,由于信号现在的最高频率是pi/2而不是pi,根据奈奎斯特规则,有一半的采样点可以扔掉。于是接下来进行2倍的下采样处理,即每隔一个点扔掉一个点。这就是一层的分解,并可用数学公式表示如下:
其中yhigh[k] 和 ylow[k]分别表示经过2倍下采样处理后的高通和低通滤波器输出。 上述分解过程将时间分辨率变为一半,因为仅需用一半的采样点即可表示原始信号的全部特征。但是,同时频率分辨率变为原来的两倍,因为信号的频带只有原始信号频带的一半,从而有效降低了一半的频率不确定性。上述分解过程,也即是通常所说的子带编码,还可以重复下去不断分解。在每一层的分解中,滤波和下采样会导致采样点数不断减半(因此使时间分辨率减半)和频带减半(因此使频率分辨率加倍)。整个的分解过程,即子带编码过程如图4.1所示。其中x[n]为待分解的原始信号,h[n]和g[n]分别为低通和高通滤波器。每一层的信号带宽在图中用"f"来标识。
来看一个例子,假定原始信号x[n]的长度为512,频带范围为0到pi弧度。在第一层分解中,信号通过高通和低通滤波器,接下来再进行2倍的下采样处理。高通滤波器的输出为256个采样点(因此时间分辨率降低一半),但是频带范围为pi/2到pi(因此提高了一采样点倍的频率分辨率)。这256个采样点即为第一层DWT的系数。低通滤波器的输出也是256个采样点,其频带范围为0到pi/2。在下一层的分解中,低通滤波之后的信号再经过同样的低通和高通滤波器。经第二个低通滤波及下采样之后,信号只有128个采样点,频带范围从0到pi/4。经第二个高通滤波及下采样之后,信号也只有128个采样点,频带范围从pi/4到pi/2。第二个高通滤波器的输出信号就是第二层DWT的系数。相比第一层,此时信号的时间分辨率降低一半,频率分辨率提高一倍。相比原始信号,时间分辨率降为原来的1/4,频率分辨率提高为原来的4倍。低通滤波器的输出继续作为下一层分解的输入信号,再进行滤波及下采样的处理。这个过程一直继续下去,直到剩下2个采样点才结束分解。对于刚才的这个512点的例子,总共需要进行8层分解,每一层的输出点数均为前面一层点数的一半。最后,原始信号的DWT包括从最后一层开始的每层的高通输出。信号的DWT的数据点数与原始信号的点数是一样多的。 原始信号中主要的频率分量,在DWT中会出现在包含这些特定频率的信号中。与傅立叶变换不同的是,这种变换保留了信号的时间信息。当然,时间的分辨率与频率所处的层数有关。如果信息的频段在高频部分,那么时间分辨率就很高,因为此时信号的采样点数还比较多。如果信息的频段在低频部分,那么时间分辨率就很低,因为此时信号只有很少的几个采样点。DWT在高频部分能获得很好的时间分辨率,在低频部分能获得很好的频率分辨率。许多实际的信号正是需要这种处理。 原始信号中不重要的频段幅度很小,于是扔掉这部分信号几乎不会丢失信息,这样就可以减少数据量。图4.2给出了一个DWT的例子,并演示了如何减小数据量。图4.2a是一个512点的信号,其幅度进行了归一化处理。横轴为采样点,纵轴为归一化幅度。图4.2b为4.2a所示信号的DWT,总共需要8层分解。最后256个点对应着原始信号最高频带,即pi/2到pi的频率。往前的128个点对应着次高频带,即pi/4到pi/2的频率。以此类推。由图可以看出,仅有前面64个点的值比较大,这对应着信息主要在低频部分。其余部分的值都很小,几乎不携带任何信息。因此,只要保留最前面的64个点,其余的点可以全部扔掉,仍然不会损失信息。这正是DWT可以有效减少数据量的原理所在。
图4.2所示的DWT的例子后面还会再讨论,因为它为解释小波变换提供了一个很重要的视角。但现在,我们还是先来简单小结一下DWT的数学形式。 DWT的一个重要性质是高通滤波器和低通滤波器单位冲激响应之间存在某种对应关系。高通滤波器和低通滤波器的关系如下:
其中g[n]是高通滤波器,h[n]为低通滤波器,L为滤波器长度(点数)。这也即是说,低通和高通滤波器是互为翻转奇数点反转的对称关系。【译者补充:这种对称关系也就是先在时间轴上反转,然后再对奇数点的值求反】。低通到高通的反转通过(-1)n项来实现。信号处理中,通常将满足上述条件的滤波器称为正交镜像滤波器(QMF)。滤波及下采样过程可以用如下公式表示:
这种情况下信号的重构非常容易,因为半带滤波器构成一组正交基。重构的时候,按照上述过程完全相反的次序即可。每层的信号先进行2倍的上采样,然后再通过合成滤波器 g’[n]和h’[n](分别对应高通和低通滤波器),然后再相加。非常有意思的是,分析和合成滤波器系数是完全一样的,只是在时间上是翻转的。因此,每层的重构公式均可表示为:
要注意的是,如果滤波器不是理想的半带滤波器,那么就不可能实现信号的完全重构。虽然在实际中不可能实现理想的滤波器,但在某种条件下设计出实现完全重构的滤波器还是有可能的。其中最著名的是由数学家Ingrid Daubechies设计的一组滤波器,被称为Daubechies(多贝西)小波。 由于要不断地进行2倍的下采样,因此为了有效进行处理,要保证信号长度为2的幂次方,或者至少是2的幂次方的倍数。信号的长度觉得了需要分解的层数。例如,信号长度为1024,那么分解层数为10。 解释DWT系数并不是一件容易的事,因为DWT系数的表示方式很特别。长话短说,每层的DWT系数组合起来,从最后一层的开始到第一层为止,所得结果即为DWT。下面用一个例子来更清晰地说明这个概念: 假定有一个256长度的信号,采样率为10MHz,需要计算其DWT。由于信号的采样频率为10MHz,信号的最高频率分量为5MHz。在第一层,信号通过低通滤波器h[n]和高通滤波器g[n],滤波器的输出均再进行2倍的下采样。经下采样处理后的高通滤波器的输出即为第一层的DWT系数。这部分数据量为128点,包括了信号中频率从2.5MHz到5MHz的信息。这128个数据按顺序排在最终结果的最后128个位置。经下采样处理后的低通滤波器的输出,同样也是128个数据点,包括了信号中频率从0到2.5MHz的信息。低频的输出再通过同样的滤波器h[n]和g[n]来进行下一层的分解。第二层经下采样处理后的高通滤波器的输出为第二层的DWT系数,数据点数为64点。在最终结果中的位置是在第一层128个系数之前。也是按顺序排列,占据64个位置。第二层经下采样处理后的低通输出继续进行分解,同样是通过滤波器h[n]和g[n]。经下采样处理后的高通滤波器输出为第三层的DWT系数,数据点数为32点。在最终结果中的位置是在第二层64个系数之前。也是按顺序排列,占据32个位置。 这个过程一直继续下去,直到在第9层仅有1个DWT系数需要计算。这个系数排列在DWT最终结果的第一的位置。接下来是第8层分解得到的2个系数,接下来是第7层分解得到的4个系数,第6层分解得到的8个系数,第5层分解得到的16个系数,第4层分解得到的32个系数,第3层分解得到的64个系数,第2层分解得到的128个系数,第1层分解得到的256个系数。【译者注:256点长度的数据仅需8层分解即可,这里好像有点小问题,但原文如此,读者要注意】。由上面的讨论可知,在随着频段越来越低,数据点数也越来越少,也即是说,随着时间分辨率随着频率的降低而降低。但是,由于在低频段频率间隔越来越小,因而频率分辨率越来越高。显然,最前面很少的几个系数并不能携带全部的信息,因为此时时间分辨率太低。 为了说明DWT这种异乎寻常的能力,我们还是来看一个现实世界中的信号实例。原始信号为256点长的超声波信号,采样频率为25MHz。信号是利用2.25MHz的变频器来产生的,因此最主要的频谱分量在2.25MHz。最后128个采样点对应着[6.25 12.5]频率范围的信号。由图可见,这128个数据点的数值均很小,几乎不包含任何信息。因此扔掉这些点不会造成任何信息的损失。【译者注:这里图可能指图4.1,但原文中并没有明确说明】。前面64个数据点对应着[3.12 6.25] MHz频率范围的信号,同样幅度很小,没有携带什么有用的信息。很小的毛刺对应着信号中的高频噪声分量。再前面的32点对应着[1.5 3.1] MHz频率范围的信号。由图可以看到,信号最主要的能量就是集中在这32个数据点,这与期望的是一致的。再前面的16点对应着[0.75 1.5]MHz频率范围的信号,这层中的峰值可能表示信号的低频包络。再往前的采样点看起来也没有包含什么很重要的信息。于是,我们只要保留第3层和第4层的系数即可。这也就是说,我们可以仅用16+32=48个数据点来表示这个256点长的信号,这大大降低了信号的数据量。 从小波变换特性中受益最大的可能是图像处理。众所周知,图像,特别是高分辨率图像,占用了大量的磁盘空间。其实,要是下载本教程花费了较长的时间,那么最主要的原因可能就是来自图片。DWT可以用来降低图像的尺寸,但同时几乎不降低分辨率。原因如下: 对于一幅给定的图像,可以计算其DWT,比如说对每行计算DWT。然后设定一个门限,凡是DWT系数低于这个门限的数据点全部扔掉。于是对于每行来说,只保留哪些高于门限的DWT系数。当需要重构时,那些被扔掉的系数用0来填充即可。通过逆DWT,即可得到原始的图像。我们同样可以在不同的频段上对图像进行分析,并且仅用某一特定频段的系数来重构原始图像。希望能在不久的将来用一幅实际的例子来说明这点。 另外一个越来越受关注的问题是,当进行分解(子带编码)时,不仅只有低频分量进一步分解,高频分量同样也还要进一步分解。换句话说,对信号的低频和高频部分同时放大。也即是说图4.1中两边同时出现树形结构。这也就是所谓的小波包。这里并不准备讨论小波包,因为它超出了本教程的范围。读者若对小波包感兴趣,或者还想更详细地了解DWT,可以参考一些市面上小波方面的著名教材。 最后小结一下这个微型系列的小波教程。如果它能对试图理解小波的读者有点帮助的话,作者会感到为这个教程所花费的时间和努力都是值得的。需要提醒的是,这个教程并不覆盖全部的小波变换的内容。这仅仅是小波概念概述性的教程,希望作为现有的还显复杂的小波教材的一份补充阅读资料。对其中还存在的结构上或者技术上的错误,欢迎指正。 感谢阅读本教程。 |