ipcount=ipcount+1;
else
codetable（icodecount）=ipcount;
codetable（icodecount+1）=p1;
icodecount=icodecount+2;
p1=p2;
ipcount=1;
end;
end;
codetable（icodecount）=ipcount;
codetable（icodecount+1）=p1;
icodecount=icodecount+2;
codetable（icodecount）=nn;%行结束符号
icodecount=icodecount+1;
end;
codetable（icodecount）=65535;%码表结束符号
（2）连续灰度图像的霍夫曼编码程序代码
functions=reduce（p）
s=cell（length（p）,1）;
fori=1:
length（p）
s{i}=i;
end;
n=size（s,1）;
whilen>2
[p,i]=sort（p）;
p
（2）=p
（1）+p
（2）;
p
（1）=[];
s=s（i）;
s{2}={s{1},s{2}};
s
（1）=[];
n=size（s,1）;
end;
functionmakecode（sc,codeword）
globalCODE
ifisa（sc,'cell'）
makecode（sc{1},[codeword0]）;
makecode（sc{2},[codeword1]）;
else
CODE{sc}=char（'0'+codeword）;
end;
functionCODE=huffman（p）
globalCODE
CODE=cell（length（p）,1）;
iflength（p）>1
p=p/sum（p）;
s=reduce（p）;
makecode（s,[]）;
else
CODE={'1'};
end;
functiony=mat2huff（x）
y.size=uint32（size（x））;
x=round（double（x））;
xmin=min（x（:
））;
xmax=max（x（:
））;
pmin=double（int16（xmin））;
pmin=uint16（pmin+32768）;
y.min=pmin;
x=x（:
）';
h=histc（x,xmin:
xmax）;
maxh=max（h）;
ifmaxh>65535
h=65535*h/maxh;
end;
h=uint16（h）;
y.hist=h;
map=huffman（double（h））;
hx=map（x（:
）-xmin+1）;
hx=char（hx）';
hx=hx（:
）';
hx（hx==''）=[];
ysize=ceil（length（hx）/16）;
hx16=repmat（'0',1,ysize*16）;
hx16（1:
length（hx））=hx;
hx16=reshape（hx16,16,ysize）;
hx16=hx16'-'0';
twos=pow2（15:
-1:
0）;
%y.code=uint16（sum（t3,2））';
t1=ones（ysize,1）;
t2=twos（t1,:
）;
t3=hx16.*t2;
t4=sum（t3,2）;
y.code=uint16（t4）';
实验四图像DCT变换编码与压缩
一、实验题目：
图像DCT变换编码与压缩
二、实验目的：
（1）掌握离散余弦变换DCT的实现方法，了解DCT的幅度分布特性，从而加深对DCT变换的认识。
（2）掌握图像DCT变换编码的实现方法，从而加深对变换编码压缩图像原理的理解。
三、实验内容：
编程实现图像DCT变换编码。
四、预备知识：
（1）熟悉DCT原理。
（2）熟悉变换编码原理。
（3）熟悉在MATLAB环境下对图像文件的I/O操作。
五、实验原理：
变换编码是在变换域进行图像压缩的一种技术。
图2.2.1显示了一个典型的变换编码系统。
图4.1变换编码系统
在变换编码系统中，如果正变换采用DCT变换就称为DCT变换编码系统。
DCT用于把一幅图像映射为一组变换系数，然后对系数进行量化和编码。
对于大多数的正常图像来说，多数系数具有较小的数值且可以被粗略地量化（或者完全抛弃），而产生的图像失真较小。
DCT是仅次于K-L变换的次最佳正交变换，且以获得广泛应用，并成为许多图像编码国际标准的核心。
离散余弦变换的变换核为余弦函数，计算速度快，有利于图像压缩和其他处理。
对于N×N的数字图像，二维DCT变换的正反变换可表示为：
（4.1）
其中，
MATLAB图像处理工具箱实现离散余弦变换有两种方法：
（1）使用函数dct2，该函数用一个基于FFT的算法来提高当输入较大的方阵时的计算速度。
（2）使用由dctmtx函数返回的DCT变换矩阵，这种方法较适合于较小的输入方阵（例如8×8或16×16）。
①函数：
dct2
实现图像的二维离散余弦变换。
调用格式为：
B=dct2（A）
B=dct2（A,[MN]）
B=dct2（A,M,N）
式中A表示要变换的图像，M和N是可选参数，表示填充后的图像矩阵大小，B表示变换后得到的图像矩阵。
②函数：
dctmtx
除了用dct2函数实现二维离散余弦变换，还可用dctmtx函数来计算变换矩阵，调用格式为：
D=dctmtx（N）
式中D是返回N×N的DCT变换矩阵，如果矩阵A是N×N方阵，则A的DCT变换可用D×A×D’来计算。
这在有时比dct2计算快，特别是对于A很大的情况。
③函数：
idct2
实现图像的二维离散余弦反变换。
调用格式为：
B=idct2（A）
B=idct2（A,[MN]）
B=idct2（A,M,N）
式中参数同dct2。
此外，为了实现8×8子块的DCT图像变换还要用到MATLAB中的blkproc函数。
将这个函数和函数dctmtx一起用于块处理可以大大简化运算。
调用函数blkproc的格式为：
B=blkpro（A,[M,N],FUN,P1,P2,…）
其中，A表示原图像，[M,N]指定了大小为M×N的滑动邻域，FUN是对M×N的矩阵进行计算的函数，Pi是传递给FUN的附加参数。
该函数自动实现图像块处理的整个过程。
Blkproc把A分成M×N个块，对每个块调用参数为P1，P2，…的函数FUN，并重新将结果组合到输出图像B。
blkproc函数实现n×n矩阵的DCT变换和反变换。
编程中可写成：
Y=blkproc（F，[88]，’P1*x*P2’，H，H’）。
同样的道理，blkproc函数还用于量化和反量化。
显示误差直方图可能用到的MATLAB函数有：
Max%找图像差最大值
[]=hist%用于生成直方图数据
Bar%显示图像差值直方图
以上函数用MATLAB的help查看具体使用方法。
图4.2显示了采用JPEG标准化矩阵进行DCT变换编码的结果。
图4.2DCT变换编码
六、实验步骤：
DCT变换编码流程如下:
步骤1：
设置JPEG标准化数组;
步骤2：
求8×8快的DCT变换矩阵;
步骤3:
计算8×8快的DCT变换；
步骤4：
对DCT系数量化和反量化；
步骤5：
求反量化系数的逆DCT变换；
步骤6：
重新显示重建图像、误差图像和误差图像的直方图。
量化时可采用JPEG标准推荐的归一化数组，如表4.1所示。
表4.1JPEG标准化数组
16
11
10
16
24
40
51
61
12
12
14
19
26
58
60
55
14
13
16
24
40
57
69
56
14
17
22
29
51
87
80
62
18
22
37
56
68
109
103
77
24
35
55
64
81
104
113
92
49
64
78
87
103
121
120
101
72
92
95
98
112
100
103
99
七、思考题目：
（1）观察图像8×8子块的DCT系数的分布，并分析其特点。
（2）将量化步长分别增大为初始值的2倍、4倍、8倍后再进行DCT变换编码，显示不同量化步长条件下的重建图像、误差图像以及误差图像的直方图。
分析重建图像质量和量化步长的关系。
八、实验程序代码：
functiony=dctcoder（x,quality）
m=[1611101624405161
1212141926586055
1413162440576956
1417222951878062
182237566810910377
243555648110411392
49647887103121120101
7292959811210010399]*quality;
[xm,xn]=size（x）;
xx=x;
x=double（x）;
t=dctmtx（8）;
y1=blkproc（x,[88],'P1*x*P2',t,t'）;%DCT
y2=blkproc（y1,[88],'round（x./P1）',m）;%Q
y3=blkproc（y2,[88],'x.*P1',m）;%IQ
y=blkproc（y3,[88],'P1*x*P2',t',t）;%IDCT
subplot（2,2,1）;imshow（xx）;title（'原始图像'）;
subplot（2,2,2）;imshow（mat2gray（y））;title（'重建图像'）;%reconstrutedimage
d=x-y;%original-reconstruted
subplot（2,2,3）;imshow（mat2gray（d））;title（'误差图像'）;
[h,k]=hist（d（:
）,512）;
subplot（2,2,4）;bar（k,h,'k'）;title（'误差图像直方图'）;
实验五图像小波变换编码与压缩
一、实验题目：
图像小波变换编码与压缩
二、实验目的：
掌握图像压缩编码的原理，熟悉小波变换的基本知识，了解嵌入零数小波编码原理，并应用MATLAB编程实现嵌入零数小波编、解码程序。
三、实验内容：
（1）对图像进行小波分解求取小波系数；
（2）编程实现嵌入零数小波的编码程序；
（3）编程实现嵌入零数小波的解码程序。
四、预备知识：
（1）熟悉图像读写和显示；
（2）理解图像压缩编码的原理；
（3）理解小波变换的原理；
（4）理解嵌入零数小波编码原理。
五、实验原理：
随着多媒体和网络技术的迅猛发展,人们对基于网络的图像传输、浏览和检索等应用提出了愈来愈广泛的需求,研究具有高效的压缩能力、支持用户个性化浏览并适合网络渐进传输的图像编码方法已成为目前静态图像编码领域的研究热点。
在大量的图像编码方法中,Shapiro提出的嵌入式零树小波编码算法（EmbeddedZerotreeWaveletencoder，简称EZW编码）是公认的效率最高的图像渐进式编码（progressiveencoding）方法之一,这种算法得到的比特流中的比特按其重要性排序。
使用这种算法，编码者能够在任一点结束编码，允许精确到任一个目标比特率或目标失真率。
（1）嵌入零数小波编码算法原理
①内嵌编码（EmbeddedCoding）
内嵌编码就是编码器将待编码的比特流按重要性的不同进行排序，根据目标码率或失真度大小要求随时结束编码；同样，对于给定码流解码器也能够随时结束解码，并可以得到相应码流截断处的目标码率的恢复图像。
内嵌编码的次序是从最重要的位（最高位）到最不重要的位（最低位）逐个发送，直到达到所需码率时停止。
内嵌编码的输出信息主要包括两部分：
排序信息和重要像素的位信息。
位信息是编码必不可少的有效信息；排序信息是辅助信息，按其重要性从左到右排列，反映了重要像素在原图上的空间位置，用于恢复原始的数据结构。
一幅图像经过三级小波分解后形成了10个子带，如图5-1所示。
小波系数的分布特点是越往低频子带系数值越大，包含的图像信息越多，对视觉比较重要，如图7-1中的LL3子带。
越往高频子带系数值越小，包含的图像信息越少，对视觉来说不太重要。
这样对相同数值的系数选择先传较低频的系数的重要比特，后传较高频系数的重要比特。
正是由于小波系数具有这些特点，它非常适合于嵌入式图像的编码算法。
在JPEG2000标准中以小波变换作为图像编码的变换方法。
图5.1小波系数间的父子关系
②几个基本概念
一幅原始图像经过变换后得到的是空间分布的小波系数，引入字符来表示每一个小波系数状态，这样更有利于将来的量化和熵编码。
零数根（ZTR）：
在本次量化中，该小波系数的幅值小于给定的阈值T，并且它的子孙后代中任何一个的幅值都比T小。
与此同时，还要求该系数的父节点系数的幅值大于T；
负大系数（NEG）：
在本次量化中，该小波系数的幅值大于给定的阈值T。
但是，它的最高位为1，即它实际的值是一个负数；
正大系数（POS）：
在本次量化中，该小波系数的幅值大于给定的阈值T。
并且，它的最高位为0，即它实际的值是一个正数；
孤独零（IZ）：
在本次量化中，该小波系数的幅值小于给定阈值T。
但是，它的子孙后代中至少存在一个系数的幅值大于T，即该小波系数的子孙后代中存在有大系数（无论正负）。
与此同时，该系数的父节点的幅值可以为任意值。
零数，在某一量化级中它的任何一个节点都不是大系数，其中，最上层
升级会员 