首页 | 新闻 | 新品 | 文库 | 方案 | 视频 | 下载 | 商城 | 开发板 | 数据中心 | 座谈新版 | 培训 | 工具 | 博客 | 论坛 | 百科 | GEC | 活动 | 主题月 | 电子展
返回列表 回复 发帖

在Matlab小波工具箱中测试滤波器组和数据延拓模式

在Matlab小波工具箱中测试滤波器组和数据延拓模式

在Matlab的小波工具箱中,可以设置9种数据边界的控制模式。除per(periodization)之外的模式,其做256点数据的8级分解的结果的数据量,可能比256大几倍。居士早年的处理,结果始终是256点,如果要扔掉旧程序、转入Matlab小波工具箱,那么应对模式per最感兴趣,因为在Matlab中其变换结果数据量最小,所以可能最易兼容。此外,Matlab也声明:其1维、2维平稳小波变换(swt),也是用per模式定义的。
       对滤波器序列,以偶数为周期做周期化,然后用其偶数平移系分解单位矩阵的测试,可以给类似于per模式中的各级分解与合成,提供条件。而且,在小波包树上,不同块上使用的滤波器序列、周期长度可以完全独立,可以有多种混合。
       如果,不强求相间隔的尺度上的二进栅格采样的关系,在每个尺度上都只求周期为偶数,那么至多只需要补充一个点就可以把数据长度变为偶数。当然,在逆变换时,要相应地除去延拓部分,这是技巧,而不是分解单位矩阵的问题。延拓量被最小化,就使变换结果的数据最少了。
       如果,原始数据的长度,能被2的L(正整数)次幂整除,那么做L级周期化分解的小波包树上的任一组基上的系数数目都是相同的、最小的。这应该是,离散小波变换、小波包的某种正典;其形式工整美观;背后的思想观念清晰;操作利索;在有限长数据的处理中,也简便实用,因为,地球上最著名的、应用最广的变换技术,FFT,常用的数据长度就是2的正整数次幂。其方便,还例如,可以用列举、归纳来论述实际的应用问题。
       居士在讨论小波(包)的离散变换时,坚持使用了单位矩阵分解、有限维矢量空间的多分辨分析、有限维基向量、有限维小波包等概念,以强调一种理论观念和角度,而不仅限于实践技巧。Matlab使用的是普通滤波器的观点,per模式作为一种实用技巧,对付的是被处理数据,但是,这不必造成问题。不过,居士好奇的是,其per模式未被列入帮助文档(dwtmode)的表格中,而且,紧随表格的几行文字里“略有冗余的(slightly redundant,指某些非per模式)”、两个“But”和“whatever”以及“must be”对于“完美的重建(perfectreconstruction)”的说明,使per模式与其它模式形成对比,也很容易让人担心per模式不良。

       有趣的是,《维基百科》(en.wikipedia.org)中“Daubechieswavelet”的“Implementation”部分,也把某种Periodization与其它模式(被称为高级(尖端,sophisticated)的方法)“相对”而论,而且强调只影响“very ends”。但是,居士对整个滤波器系统的周期化的测试,并不管支集影响有多远;在256点数据的8级分解中,最后一个近似系数,实际上,涉及到了所有的256个数据点。

       在Matlab中,做些测试,是应该的!

       Matlab系统默认的模式是sym,它在表格中排在第一位。此外,系统有函数计算“合理的最大”分解深度。在命令窗中,运行

wmaxlev(256,'db33'), wmaxlev(64,'db20'),

得到的最大分解深度为1、0,即是说:只该分解一级、或干脆就不能做小波变换。分解深度越小,由滤波器系数本身的近似性和计算中的数字截断所引入的误差(不是理论上的重建问题!),一般应该更小。所以,测试时,用sym模式、系统计算的合理分解深度,作为比较标准。

       本文后的图片所示(右边)的结果,用单位能量的正态分布的随机序列(左边),比较了Matlab的小波变换中的3种模式、它的65组滤波器(依次为Daubechies1至Daubechies 45、 Coiflets 1至Coiflets5、15组双正交)。模式sym、zpd、per的误差(未做任何变换域数据修正)结果,几乎完全重合,而且,也与居士在Scilab中周期化滤波器系统以分解单位矩阵的结果一致。用模式zpd、per时,分解深度都被强行设置为8;然而,用sym模式时的平均分解深度,仅约为3,不到8的一半,可以说sym模式很占便宜。

       下面是居士写的计算图片中结果的Matlab程序temp8x.m:

function [Test8x]=temp8x

%test Matlab-Wavelet-Toolbox: dwtmode and itsFilters

%in Matlab-R2011a. 2013-August.

%      tic;[Test8x]=temp8x;toc,

%------prepare

close all force;

biof1=['1.1'; '1.3'; '1.5'; '2.2'; '2.4'; '2.6';'2.8';...

      '3.1'; '3.3'; '3.5'; '3.7'; '3.9'; '4.4'; '5.5'; '6.8'];

dwtmod1=['sym';'zpd';'per']; %3 modes

rng(1.0e9); x=randn(1,256);  x=x/norm(x(); %test data

ModeSymLevel=0;

Test8x=zeros(65,3);

%------compute

fori1=1:65;    % 65filters

    if(i1>0)&&(i1<46)

       FilterNam1=['db'int2str(i1)];     %Daubechies     

    elseif(i1>45)&&(i1<51)

       FilterNam1=['coif' int2str(i1-45)]; %Coiflets

    elseif(i1>50)&&(i1<66)

       FilterNam1=['bior' biof1(i1-50,];%Biorthogonal                 

    else

       error('Undefined Filter.');  %unknown

    end

    forj1=1:3;    %3mode error, from dwt and idwt

       if j1==1

           dl0=wmaxlev(length(x),FilterNam1);%level,only for Mode sym

           ModeSymLevel=ModeSymLevel+dl0; %memorize the levels

       else

           dl0=8;

       end

       dwtmode(dwtmod1(j1,); %set mode

       [xc1,xl1]=wavedec(x,dl0,FilterNam1);%transform      

       xr1=waverec(xc1,xl1,FilterNam1);  %inverse

       Test8x(i1,j1)=norm(xr1-x);    %norm of error

   end

end

%------display

ModeSymLevel=ModeSymLevel/length(Test8x(:,1)), %'sym' average level

subplot(1,2,1);  plot(x);  axis tight; %display raw data

title('The Normalized randn(1,256) forwavedec');

subplot(1,2,2); %display the result

semilogy(Test8x(:,1),'-g');  holdon;

semilogy(Test8x(:,2),'ob'); semilogy(Test8x(:,3),'+r');  hold off;

grid on; axis tight; legend(dwtmod1,'Location','NorthWest');

title('temp8x.m: Matlab dwt Modes and ItsFilters');

xlabel('db1-db45,coif1-coif5,15 BiorthogonalFilters');

ylabel('Reconstruction-Error Norm');

end

继承事业,薪火相传
返回列表