EZW算法的过程详解和Matlab代码（1）构建扫描次序表（修正小波树结构）

前段时间，我们讨论了嵌入式小波零树算法的基本原理。（http://blog.csdn.net/chenyusiyuan/archive/2007/11/16/1888968.aspx）一个多星期过去了，我根据自己对算法的理解，编写出EZW算法的Matlab程序，可以实现图像的任意级别的小波分解和重构、以及任意精度的EZW编解码过程。下面，我们以一幅16*16的Lena局部图像经过3级小波分解后的小波数据为例详细的说明EZW算法的编解码过程，并给出相应的Matlab代码。原始图像如下

分解后的小波图像数据为：

一、在开始编码之前，首先要求出初始阈值T1：

MaxDecIm=max(max(abs(DecIm)));

T=zeros(1,codeDim);

T(1)=2^flor(log2(MaxDecIm));

二、然后是建立小波树结构，构建扫描次序表。这个扫描次序表非常重要，后面的编码、解码过程都要按照扫描次序表逐个处理数据矩阵的各个元素。构建过程如下：用（r，c）表示数据矩阵上各元素的位置。row、col作为全局变量，表示数据矩阵的行、列数。

1、小波树结构的特点：

（1）对于LL-N低频子带的点(r,c)，有3个孩子：(r,c+W)、(r+H,c)和(r+H,c+W)，其中W、H分别是LL-N子带的宽和高；

（2）第N~2高频子带(LH、HL、HH)的点都有4个孩子，即： tp=[2*r-1,2*c-1;2*r-1,2*c;2*r,2*c-1;2*r,2*c];

（3）第1高频子带的点没有孩子。根据小波树的这个特点，可编写如下“小波树”函数treeMat()，输入矩阵内任一点的位置（r，c），给出该点的子孙列表cp 。

function cp=treeMat(r,c) %这个函数是一个递归函数

global row col dim % dim是小波分解级数

HLL=row/2^dim;

WLL=col/2^dim;

if (r<=HLL)&&(c<=WLL)

tp1=[r,c+WLL;r+HLL,c;r+HLL,c+WLL];

cP=[tp1;treeMat(r,c+WLL);treeMat(r+HLL,c);treeMat(r+HLL,c+WLL)];

elseif (r>row/2)||(c>col/2)

cP=[];

else

tp=[2*r-1,2*c-1;2*r-1,2*c;2*r,2*c-1;2*r,2*c];

tm1=[];tm2=[];tm3=[];tm4=[];

if (tp(4,1)<=row/2)&&(tp(4,2)<=col/2)

t1=treeMat(tp(1,1),tp(1,2));

tm1=[tm1;t1];

t2=treeMat(tp(2,1),tp(2,2));

tm2=[tm2;t2];

t3=treeMat(tp(3,1),tp(3,2));

tm3=[tm3;t3];

t4=treeMat(tp(4,1),tp(4,2));

tm4=[tm4;t4];

end

cP=[tp;tm1;tm2;tm3;tm4];

end

示例，当row=8，col=8，dim=2时，LL-N低频子带的点(1,1) (1,2) (2,1) (2,2) 的子孙分布如下：

由这个小波树列表tree_p，我们可以进一步构建函数childMat()，给出矩阵数据Mat和矩阵任一点的位置（r，c），返回该点的子孙数据列表chMat。

function chMat=childMat(Mat,chRows,chCols)

global row col dim

chPoint=treeMat(chRows,chCols);

chMat=[];

[mRows,mCols]=size(chPoint);

for iRows=1:mRows

chMat=[chMat;chPoint(iRows,1),chPoint(iRows,2),Mat(chPoint(iRows,1),chPoint(iRows,2))];

end

2、构建扫描次序表

本文EZW算法的扫描次序为Morton式，其特征是从（1，1）开始，每4个点组成一个“Z”型扫描单元，从微观到宏观上都是严格的“Z”型结构，可以用递归方法来构建扫描次序表。扫描次序表scanlist由两部分组成，一个是将数据矩阵Mat按照morton扫描次序转换成数据列表matlist，一个是按照扫描次序组成的矩阵各点位置的（r，c）列表lsorder。

function scanlist=morton(Mat)

global row col

matlist=mat2list(Mat);

scanorder=listorder(row,col,1,1);

scanlist=[];

for i=1:row*col

scanlist=[scanlist;i scanorder(i,:) matlist(i)];

end

function mls=mat2list(Mat) % 该函数为递归函数

[r,c]=size(Mat);

if (r==2)&&(c==2)

mls=[Mat(1,1);Mat(1,2);Mat(2,1);Mat(2,2)];

else

M1=Mat(1:r/2,1:c/2);

M2=Mat(1:r/2,c/2+1:c);

M3=Mat(r/2+1:r,1:c/2);

M4=Mat(r/2+1:r,c/2+1:c);

lt1=mat2list(M1);

lt2=mat2list(M2);

lt3=mat2list(M3);

lt4=mat2list(M4);

mls=[lt1;lt2;lt3;lt4];

end

function lsorder=listorder(mr,mc,pr,pc) % 该函数为递归函数

lso=[pr,pc;pr,pc+mc/2;pr+mr/2,pc;pr+mr/2,pc+mc/2];

mr=mr/2;mc=mc/2;

lm1=[];lm2=[];lm3=[];lm4=[];

if (mr>1)&&(mc>1)

ls1=listorder(mr,mc,lso(1,1),lso(1,2));

lm1=[lm1;ls1];

ls2=listorder(mr,mc,lso(2,1),lso(2,2));

lm2=[lm2;ls2];

ls3=listorder(mr,mc,lso(3,1),lso(3,2));

lm3=[lm3;ls3];

ls4=listorder(mr,mc,lso(4,1),lso(4,2));

lm4=[lm4;ls4];

end

lsorder=[lso;lm1;lm2;lm3;lm4];

len=length(lsorder);

lsorder=lsorder(len-mr*mc*4+1:len,1:2);