实验四-DCT变换HUFFman编码图像压缩

实验四图像压缩姓名：学号：邮箱：一、实验目的1.掌握DCT变换的原理2.了解DCT变化在图像压缩中的应用3.掌握图像压缩的基本原理及方法4.了解霍夫曼编码原理5.熟悉图像压缩的MATLAB编程二、实验原理DCT是目前比较好的图像变换，它有很多优点。DCT是正交变换，它可以将8x8图像空间表达式转换为频率域，只需要用少量的数据点表示图像；DCT产生的系数很容易被量化，因此能获得好的块压缩；DCT算法的性能很好，它有快速算法，如采用快速傅立叶变换可以进行高效的运算，因此它在硬件和软件中都容易实现；而且DCT算法是对称的，所以利用逆DCT算法可以用来解压缩图像。由于DCT主要应用在数据和图像的压缩，因此希望原信号的能量在变换后能尽量集中在少数系数上，且这些大能量的系数能处在相对集中的位置，这将有利于进一步的量化和编码。但是如果对整段的数据或整幅图像来做DCT，那就很难保证大能量的系数能处在相对集中的位置。因此，在实际应用中，一般都是将数据分成一段一段来做，一般分成8x8或16x16的方块来做。二维DCT正交变换的公式为：二维DCT逆变换公式：其中三、实验要求利用DCT变换对图像进行压缩，对比不同压缩比下的结果，对比不同压缩比下图像大小的变化。压缩过程如下图所示：四、实验过程与结果实验程序如下：（先给出主程序，然后给出各功能子函数的程序）主程序：clearload('lena.mat')%调入170*170大小的一幅彩色lena图像l=imresize(lena,[256256]);%将图像变换为8的整数倍大小X=rgb2gray(l);Y1=double(X);%读入图像数据lianghua=[1611101624405161;%61;%量化矩阵，量化的程度序决定压缩比1212141926586055;1413162440576956;1417222951878062;182237566810910377;243555648110411392;49647887103121120101;7292959811210010399];ilianghua=lianghua;%----------------------------------------------------------%图像压缩%----------------------------------------------------------t=dctmtx(8);J=blkproc(Y1,[88],'P1*x*P2',t,t');%分成8*8块进行DCT变换M=blkproc(J,[88],'round(x./P1)',lianghua);%量化u=abs(min(min(M1));M=(M1./u)+1;data=uint8(M);%Huffman编码要求为无符号整形数组读入图像DCT变换量化huffman编码M2=M-double(data);[zipped,info]=huffencode(data);%调用Huffman编码程序进行压缩unzipped=huffdecode(zipped,info,data);%调用Huffman解码程序进行解压缩k=1;fori=1:256forj=1:256unzippedray(i,j)=unzipped(k);k=k+1;endendunzippedray=unzippedray';%对解压缩后得到的一维数组进行变换，得到无损的量化后%二维数组，其值与data数组值是一致的，体现了Huffman编码是一种无损编码unzippedray=(double(unzippedray)-1+M2).*u;T=blkproc(unzippedray,[88],'x.*P1',ilianghua);%反量化I=blkproc(T,[88],'P1*x*P2',t',t);%8*8DCT反变换%----------------------------------------------------------%调用Huffman编码程序进行解码%显示原始图像和经编码后的图像，显示压缩比，并计算均方根误差得erms=0，表示是Huffman是无失真编码figuresubplot(221);imshow(Y1,[]);axissquare;xlabel('原256*256灰度图像');subplot(222);imshow(I,[]);axissquare;xlabel('Huffman解压缩后图像');subplot(223);imshow((Y1-I),[]);axissquare;xlabel('量化后损失的图像部分');[h,k]=hist((Y1-I),256);%生成直方图数据subplot(224);bar(k,h,'k');title('误差图像直方图');%subplot(224);imshow(I,[]);axissquare;xlabel('压缩图像');%erms=compare(data(:),unzipped(:))cr=info.ratiowhosdataunzippedzipped%huffencode函数对输入矩阵vector进行Huffman编码，返回%编码后的向量（压缩后数据）及相关信息Huffman编码子程序：function[zipped,info]=huffencode(vector)%输入和输出都是unit8格式%info返回解码需要的机构信息%info.pad是添加的比特数%info.huffcodes是Huffman码字%info.rows是原始图像行数%info.cols是原始图像行数%info.length是原始图像数据长度%info.maxcodelen是最长码长if~isa(vector,'uint8')error('inputargumentmustbeauint8vector');end[m,n]=size(vector);vector=vector(:)';f=frequency(vector);%计算各符号出现的概率symbols=find(f~=0);f=f(symbols);[f,sortindex]=sort(f);%将符号按照出现的概率大小排序symbols=symbols(sortindex);len=length(symbols);symbols_index=num2cell(1:len);codeword_tmp=cell(len,1);whilelength(f)1%生产Huffman树，得到码字编码表index1=symbols_index{1};index2=symbols_index{2};codeword_tmp(index1)=addnode(codeword_tmp(index1),uint8(0));codeword_tmp(index2)=addnode(codeword_tmp(index2),uint8(1));f=[sum(f(1:2))f(3:end)];symbols_index=[{[index1,index2]}symbols_index(3:end)];[f,sortindex]=sort(f);symbols_index=symbols_index(sortindex);endcodeword=cell(256,1);codeword(symbols)=codeword_tmp;len=0;forindex=1:length(vector)%得到整个图像所有比特数len=len+length(codeword{double(vector(index))+1});endstring=repmat(uint8(0),1,len);pointer=1;forindex=1:length(vector)%对输入图像进行编码code=codeword{double(vector(index))+1};len=length(code);string(pointer+(0:len-1))=code;pointer=pointer+len;endlen=length(string);pad=8-mod(len,8);%非8整数倍时，最后补pad个0ifpad0string=[stringuint8(zeros(1,pad))];endcodeword=codeword(symbols);codelen=zeros(size(codeword));weights=2.^(0:23);maxcodelen=0;forindex=1:length(codeword)len=length(codeword{index});iflenmaxcodelenmaxcodelen=len;endiflen0code=sum(weights(codeword{index}==1));code=bitset(code,len+1);codeword{index}=code;codelen(index)=len;endendcodeword=[codeword{:}];%计算压缩后的向量cols=length(string)/8;string=reshape(string,8,cols);weights=2.^(0:7);zipped=uint8(weights*double(string));%码表存储到一个稀疏矩阵huffcodes=sparse(1,1);forindex=1:nnz(codeword)huffcodes(codeword(index),1)=symbols(index);end%填写解码时所需的结构信息info.pad=pad;info.huffcodes=huffcodes;info.ratio=cols./length(vector);info.length=length(vector);info.maxcodelen=maxcodelen;info.rows=m;info.cols=n;%huffdecode函数对输入矩阵vector进行Huffman编码，%返回解压后的图像数据endHuffman解码子程序：functionvector=huffdecode(zipped,info,image)if~isa(zipped,'uint8')error('inputargumentmustbeauint8vector');end%产生0，1序列，每位占一个字节len=length(zipped);string=repmat(uint8(0),1,len.*8);bitindex=1:8;forindex=1:lenstring(bitindex+8.*(index-1))=uint8(bitget(zipped(index),bitindex));endstring=logical(string(:)');len=length(string);%开始解码weights=2.^(0:51);vector=repmat(uint8(0),1,info.length);vectorindex=1;codeindex=1;code=0;forindex=1:lencode=bitset(code,codeindex,string(index));codeindex=codeindex+1;byte=decode(bitset(code,codeindex),info);ifbyte0vector(vectorindex)=byte-1;codeindex=1;code=0;vectorindex=vectorindex+1;endend%vector=reshape(vector,info.rows,info.cols);%函数addnode添加节点endaddnode子程序：functioncodeword_new=addnode(codeword_old,item)codeword_new=cell(size(