中国矿业大学2015-2016学年第一学期《数字视频技术》课程小设计考核设计题目:图像的算术编码研究专业班级:学生姓名:学生学号:指导教师:成绩:本人郑重声明:本人认真、独立完成了查找资料、完成作业、编写程序等考核任务,无抄袭行为。签字:日期:成绩评阅人一、设计任务、目的和要求:1.1设计任务:图像的算术编码1.2设计目的:1.了解图像压缩的意义方法,对比不同的压缩方法;2.熟悉算术编码的基本原理和特点;3.掌握改进的算术编码的方法与具体实例应用;4.掌握利用MATLAB编程实现数字图像的算术编码。1.3设计要求:要求实现灰度图像的算术编码和解码恢复图像;处理结果要求最终图像显示,且计算图像的信息熵,平均码字长度,编码效率,压缩比。二、总体方案设计2.1算术编码简介算术编码,是图像压缩的主要算法之一。是一种无损数据压缩方法,也是一种熵编码的方法。和其它熵编码方法不同的地方在于,其他的熵编码方法通常是把输入的消息分割为符号,然后对每个符号进行编码,而算术编码是直接把整个输入的消息编码为一个数,一个满足(0.0≤n1.0)的小数n。算术编码是一种到目前为止编码效率最高的统计熵编码方法,它比著名的Huffman编码效率提高10%左右,但由于其编码复杂性和实现技术的限制以及一些专利权的限制,所以并不象Huffman编码那样应用广泛。算术编码有两点优于Huffman码:①它的符号表示更紧凑;②它的编码和符号的统计模型是分离的,可以和任何一种概率模型协同工作。后者非常重要,因为只要提高模型的性能就可以提高编码效率。2.2软件运行环境系统运行环境:windows操作系统。软件编程平台:Matlab2014a。2.3编解码算法原理2.3.1编码算术编码将整个要编码的数据映射到一个位于[0,1)的实数区间中。并且输出一个小于1同时大于0的小数来表示全部数据。利用这种方法算术编码可以让压缩率无限的接近数据的熵值,从而获得理论上的最高压缩率。算术编码进行编码时,从实数区间[0,1)开始。按照符号的频度将当前的区间分割成多个子区间。根据当前输入的符号选择对应的子区间,然后从选择的子区间中继续进行下一轮的分割。不断的进行这个过程,直到所有符号编码完毕。对于最后选择的一个子区间,输出属于该区间的一个小数。这个小数就是所有数据的编码。子区间计算的迭代递推公式:StartN=StartB+LeftC*LEndN=StartB+RightC*L其中StartN表示新子区间的起始位置,EndN表示新子区间的结束位置,StartB表示前子区间的起始位置,LeftC表示当前符号的区间左端,RightC表示当前符号的区间右端,L表示前子区间宽度。2.3.2解码算术编码进行解码时仅输入一个小数。解码前首先需要对区间[0,1)按照初始时的符号频度进行分割。然后观察输入的小数位于那个子区间。输出对应的符号,选择对应的子区间,然后从选择的子区间中继续进行下一轮的分割。不断的进行这个过程,直到所有的符号都解码出来。整个过程相当于编码时的逆运算。2.3.3算术编码的改进以上描述的算法,在当前的计算机系统上是很难实现的。尤其是无限精度的实数运算。所以在实现的时候,需要对算法做一些改进。使得它可以在当前的计算机系统上较快的运行。当然,这种改进是以降低运算精度为代价的。也就是说,这种改进实际上会降低算法的压缩率。但是,它会使算法的实现成为可能。观察前面描述的算法过程可以发现,运算时区间的上下沿都是小于1的小数。那么我们可以省略0和小数点,仅仅使用小数的尾数来表示小数。省略0和小数点后的尾数,实际上就是一个无限大的整数。使用无限整数的部分高位来表示整数,并在这些整数上进行整数运算就可以模拟出实数运算。另外,分割区间、选择子区间的过程,相当于将一个区间映射到另一个更小的区间中(以下简称“映射区间”)。如果我们知道一个符号的频度。以及符号值小于该符号的其它符号的频度总计(以下简称“累积频度(CumulativeFrequency)”)。还有到目前为止所有符号频度的总计(以下简称“总计频度(TotalFrequency)”)。那么就可以根据这些频度信息,从当前区间中计算出映射区间。计算的公式如下:Range=High-Low+1High=Low+Range*(CumFreq+Freq)/Total–1Low=Low+Range*CumFreq/Total其中Low表示区间的下沿;High表示区间的上沿;Range表示区间的范围;Freq表示符号频度;CumFreq表示累积频度;Total表示总计频度。这些变量中保存的都是整数,并进行整数运算。其中”/”表示整除。另外需要注意一点,这里使用闭区间[Low,High],而不是使用右开区间[Low,High)。在解码的时候也可以进行整数运算。根据输入的整数数值、当前区间的下沿和总计频度,可以计算出一个估算出来的累积频度(以下简称“估算频度(EstimateFrequency)”)。其计算公式如下。Range=High-Low+1EstFreq=((Value-Low+1)*Total-1)divRange其中,Value表示输入的整数数值;EstFreq表示估算频度。利用估算频度在当前的累积频度表中查找,当满足CumFreq≤EstFreq<CumFreq+Freq的条件时,就可以解码出一个符号。利用解码出的符号可以得到对应的累积频度和频度。根据这些频度信息,可以从当前区间中计算出映射区间。这一点同编码时是一样的。计算出映射区间后,更新对应符号的频度,又可以进行新的一轮解码。2.3.4流程图(1)主程序(2)概率统计(3)编码函数三、设计与实现Matlab所有程序如下:(1)Matlab主程序closeall;clearall;clc;%关闭所有图形窗口,清除工作空间所有变量,清空命令行image=imread('lena.bmp');%读取图像,得到图像的数据矩阵X=imread('lena.bmp');subplot(1,3,1);imshow(X);title('原始图像')%对图像缩小N=16;seq=imresize(image,[N,N]);%B=imresize(A,[rowscols]),[rowscols]为图像调整后的尺寸%当使用公式seq=imresize(image,M),M指缩小倍数时,利用代码[m,n]=size(seq)[r,c]=size(A),当有两个输出参数时,size函数将矩阵的行数返回到第一个输出变量r,将矩阵的列数返回到第二个输出变量c。%二维数据转化成一维k=1;forj=1:Nfori=1:Nx(k)=seq(j,i);k=k+1;endend%计算信源每个灰度值的概率[alphacnt]=probmodel(x);%cnt(i)为信源符号出现的次数,alpha(i)信源符号%开始编码btag=arithintcod(alpha,cnt,x);%计算编码长度b=length(btag);%开始解码outseq=arithintdecod(btag,alpha,cnt,length(x));%对解码出的一维数据转换为二维图像k=1;forj=1:Nfori=1:Noutimg(j,i)=outseq(k);k=k+1;endend%计算信息熵sum=N*N;H=0;%初始化信息熵fori=0:255;[r,c]=find(seq==i);%统计每个灰度值的像素点总数num(i+1)=length(r);p(i+1)=num(i+1)/sum;%统计每个灰度值的概率ifp(i+1)~=0H=H-p(i+1)*log2(p(i+1));%计算信息熵endend%计算平均码字长度pjmc=b/sum;%计算编码效率bmxl=H/pjmc;%计算压缩比ysb=sum*8/b;%输出最终结果%输出编码前图像disp('原图像');disp(x);%输出编码disp(strcat('编码=',btag));%输出解码disp('解码图像');disp(outseq);disp('信息熵');disp(H);disp('平均码字长度');disp(pjmc);disp('编码效率');disp(bmxl);disp('压缩比');disp(ysb);subplot(1,3,2),imshow(seq),title('压缩图像');subplot(1,3,3),imshow(outimg),title('解码图像');d=length(x);(2)概率统计probmodel函数%概率统计程序function[alphacnt]=probmodel(seq)%输入seq原始图像%cnt(i)信源符号出现的次数%alpha(i)信源符号if~isempty(seq)alpha(1)=seq(1);%初始化cnt(1)=1;%初始化l=length(seq);%计算循环次数k=2;fori=2:l%遍历数据idx=find(seq(i)==alpha);%寻找alpha中是否有该符号ifisempty(idx)%若没有,则添加入alphaalpha(k)=seq(i);cnt(k)=1;k=k+1;else%若有,则该该符号出现次数加一cnt(idx)=cnt(idx)+1;endendelsealpha=0;cnt=0;endend(3)编码函数%实现编码的函数functiontag=arithintcod(alpha,cnt,seq)%算术编码%alpha表示信源符号%cnt表示信源符号出现次数,概率%seq表示输入图像数据ls=length(seq);%计算输入数据的总个数CC(1)=0;%计算累计概率,便于计算上下限fori=1:length(cnt)CC(i+1)=CC(i)+cnt(i);endtotcount=CC(i+1);%总符号数m=ceil(log2(totcount*4));%需要编码的长度l=0;%下限u=2^m-1;%上限tag='';%编码标签scale=0;%开始迭代过程fori=1:ls%开始迭代编码p=find(seq(i)==alpha);%该符号在信源符号的位置l1=l+floor(((u-l+1)*CC(p))/totcount);%计算上下限u=l+floor(((u-l+1)*CC(p+1))/totcount)-1;l=l1;lb=dec2bin(l,m);%转换为字符串二进制ub=dec2bin(u,m);E2=1;E3=1;%区间扩展的两个条件while(E2|E3)fbl=lb(1);%上下限第一位fbu=ub(1);iffbl==fbu%若首位相等,则可移出,并扩展区间E2=1;tag=strcat(tag,fbl);%移出lb(1)='';%左移,扩展区间lb(end+1)='0';ub(1)='';ub(end+1)='1';iffbl=='0'sc='1';elsesc='0';endwhilescale0tag=strcat(tag,sc);scale=scale-1;endelseE2=0;endsbl=lb(2);sbu=ub(2);fbl=lb(1);%上下限第一位fbu=ub(1);iffbl~=fbuifsbl=='1'&sbu=='0'%忽略次高位,并记录忽略次数lb(1)='';lb(end+1)='0';ub(1)='';ub(end+1)='1';lb(1)='0';ub(1)='1';scale=scale+1;E3=1;elseE3=0;endendendl=bin2dec(lb);%转换为10进制u=bin2dec(ub);endf=lb(1);%对最终区间编码lb(1)='';iff=='0'sc='1';elsesc='0';endtag=strcat(tag,f);whilescale0tag=strcat(tag,sc);scale=scale-1;endtag=strcat(tag,lb);end(4)解码函