实验五 huffman编码

Huffman 编码
clear all
close all
clc
I = [0 1 3 2 1 3 2 1;0 5 7 6 2 5 6 7;1 6 0 6 1 6 3 4;2 6 7 5 3 5 6 5;3 2 2 7 2 6 1 6;
2 6 5 0 2 7 5 0;1 2
3 2 1 2 1 2;3 1 2 3 1 2 2 1];%读入一幅图像的灰度值
%load I;
[m,n]=size(I);
%将矩阵的不同灰度值统计在数组c 的第一列中
p1=1;
s=m*n;
for k=1:m
for L=1:n
f = 0;
for b = 1:p1-1;
if c(b,1)==I(k,L)
f = 1;
break;
end
end
if f==0
c(p1,1)=I(k,L);
p1 = p1+1;
end
end
end
%将相同的灰度值占总像素的比例统计在数组p 中
for g=1:p1-1
p(g)=0;
c(g,2)=0;
for k=1:m
for L=1:n
if c(g,1)==I(k,L)
p(g)=p(g)+1;
end
end
end
p(g)=p(g)/s;
end
p11=p;
%找到最小的概率，相加直到等于1，把最小概率的序号存在tree第一列中，次小的放在第二列，和放在p像素比例之后
pn=0;
po=1;
while(1)
if(pn>=1)
break;
else
[pm,p2]=min(p(1:p1-1));
p(p2)=1.1;
[pm2,p3]=min(p(1:p1-1));
p(p3)=1.1;
pn=pm+pm2;
p(p1)=pn;
tree(po,1)=p2;
tree(po,2)=p3;
po = po+1;
p1=p1+1;
end
end
%数组第一维表示值，第二维表示代码数值大小，第三维表示代码的位数
for k=1:po-1
tt=k;
m1=1;
if or(tree(k,1)<=g,tree(k,2)<=g)
if tree(k,1)<=g
c(tree(k,1),2)=c(tree(k,1),2)+m1;
m2=1;
while(tt<po-1)
m1=m1*2;
for L=tt:po-1
if tree(L,1)==tt+g
c(tree(k,1),2)=c(tree(k,1),2)+m1;
m2=m2+1;
tt=L;
break;
elseif tree(L,2)==tt+g
m2=m2+1;
tt=L;
break;
end
end
end
c(tree(k,1),3)=m2;
end
tt=k;
m1=1;
if tree(k,2)<g
m2=1;
while tt<po-1
m1=m1*2;
for L=tt:po-1
if tree(L,1)==tt+g
c(tree(k,2),2)=c(tree(k,2),2)+m1;
m2=m2+1;
tt=L;
break;
elseif tree(L,2)==tt+g
m2 =m2+1;
tt=L;
break;
end
end
end
c(tree(k,2),3)=m2;
end
end
end
%把概率小的值用1标记，概率大的值用0标记
[M,N]=size(c);
disp('编码')
A1=dec2bin(c(1,2),c(1,3))%说明：这里可以把编码存在高维数组或构造数组、细胞数组同时显示
A2=dec2bin(c(2,2),c(2,3))
A3=dec2bin(c(3,2),c(3,3))
A4=dec2bin(c(4,2),c(4,3))
A5=dec2bin(c(5,2),c(5,3))
A6=dec2bin(c(6,2),c(6,3))
A7=dec2bin(c(7,2),c(7,3))
A8=dec2bin(c(8,2),c(8,3))
for m=1:M
if p11(m)~=0
H(m)=-p11(m)*log2(p11(m));
end
end
disp('信源的熵');
H1=sum(H)
NN=0;
for i =1:M
NN=NN+p11(1,i)*c(i,3);%平均码长
end
disp('平均码长')
NN
disp('编码效率')
yita=H1/(NN*log2(2))
disp('冗余度')
Rd=1-yita
Dctde编码
function dctdemo(action,varargin)
%DCTDEMO 2-D DCT image compression demo.
%
% This demo lets you explore image compression using the 2-D
% discrete cosine transform (DCT). The original image is
% transformed in 8-by-8 blocks and then inverse transformed in
% 8-by-8 blocks to create the reconstructed image. The inverse
% DCT is performed using the subset of DCT coefficients that
% are shown in white. The error image (the difference between
% the original and reconstructed image) is also displayed.
%
% Use the slider to select a subset of the DCT coefficients and then
% press the "Apply" button to produce the corresponding reconstructed % image. The DCT coefficients with small variance are discarded
% first, as they have the least effect on the reconstruction.
% Try to find out how few coefficients are needed to create a
% reasonable approximation to the original.
%
% The following images are courtesy of J. C. Russ, The Image Processing % Handbook, Third Edition, 1998, CRC Press: Saturn and Flower.
%
% See also DCT2, BLKPROC.
% Copyright 1993-2002 The MathWorks, Inc.
% $Revision: 1.17 $ $Date: 2002/04/04 18:03:38 $
% Callback Routines:
%
% InitializeDCTDEMO() sets up the demo. It is run only once upon
% the demo comes up.
%
% ApplyDCT(DemoFig) will call MakeDCT if it hasn't been called yet
% and it will apply the DCT Matrix to the image and update the
% reconstructed image. When ApplyDCT is called from the initialization % routine, it is called with an argument which is the figure handle of
% the demo. From a callback routine, it gets the figure handle with
% gcbf.
%
% SliderUpdate() updates the mask which is applied to the DCT
% coefficients before reconstruction. Update the surface which is
% the DCT Coefficient display.
%
% MakeDCT() will grab the current image and calculate its
% 8x8 blocked DCT. It puts it in the userdata of the
% original image. This is not used as a callback.
%
% ShowInfo() will pull up helpwin with the info for the demo in it.
% Images we can use:
% pout.tif
% saturn.tif
% eight.tif
% circuit.tif
if nargin<1,
action='InitializeDCTDEMO';
end;
feval(action,varargin{:});
return;
%%%
%%% Sub-function - InitializeDCTDEMO
%%%
function InitializeDCTDEMO()
% If dctdemo is already running, bring it to the foreground.
h = findobj(allchild(0), 'tag', '2D DCT Image Compression Demo'); if ~isempty(h)
figure(h(1))
return
end
screenD = get(0, 'ScreenDepth');
if screenD>8
grayres=256;
else
grayres=128;
end
DctDemoFig=figure( ...
'Name','2-D DCT Image Compression Demo', ...
'NumberTitle','off', 'HandleV isibility', 'on', ...
'tag', '2D DCT Image Compression Demo', ...
'V isible','off', 'Resize', 'off',...
'BusyAction','Queue','Interruptible','off',...
'IntegerHandle', 'off', ...
'Doublebuffer', 'on', ...
'Colormap', gray(grayres));
%====================================
% Information for all buttons (and menus)
top=0.92;
left=0.785;
btnWid=0.175;
btnHt=0.06;
% Spacing between the button and the next command's label spacing=0.025;
%====================================
% The CONSOLE frame
frmBorder=0.02; frmBottom=0.05; frmHeight = 0.46;
yPos=frmBottom-frmBorder;
frmPos=[left-frmBorder yPos btnWid+2*frmBorder frmHeight+2*frmBorder]; h=uicontrol( ...
'Parent', DctDemoFig, ...
'BusyAction','Queue','Interruptible','off',...
'Style','frame', ...
'Units','normalized', ...
'Position',frmPos, ...
'BackgroundColor',[0.45 0.45 0.45]);
%====================================
% The Apply button
labelStr='Apply';
callbackStr='dctdemo(''ApplyDCT'')';
yPos=frmBottom+frmHeight-btnHt-spacing;
applyHndl=uicontrol( ...
'Parent', DctDemoFig, ...
'BusyAction','Queue','Interruptible','off',...
'Style','pushbutton', ...
'Units','normalized', ...
'Position',[left yPos btnWid btnHt], ...
'String',labelStr, ...
'Enable', 'off', ...
'Callback',callbackStr);
%====================================
% The Images Popup
yPos=frmBottom;
hImgPop=uicontrol( ...
'Parent', DctDemoFig, ...
'BusyAction','Queue','Interruptible','off',...
'Style','popupmenu', ...
'BackgroundColor',[.8 .8 .8], ...
'Units','normalized', ...
'Position',[left yPos+3*btnHt+2*spacing btnWid btnHt], ...
'String','Saturn|Flower|Pout|Trees|Quarter|Circuit', ...
'Tag','ImagesPop', ...
'Callback','dctdemo(''LoadNewImage'')');
%====================================
% The Images Label
h = uicontrol( ...
'Parent', DctDemoFig, ...
'BusyAction','Queue','Interruptible','off',...
'Style','text', ...
'Units','normalized', ...
'Position',[left yPos+4*btnHt+2*spacing btnWid .035 ], ...
'Horiz','center', ...
'Background',[0.45 0.45 0.45], ...
'Foreground','white', ...
'String','Select Image:');
%====================================
% The Info button
labelStr='Info';
callbackStr='helpwin dctdemo';
yPos=frmBottom;
helpHndl=uicontrol( ...
'Parent', DctDemoFig, ...
'BusyAction','Queue','Interruptible','off',...
'Style','pushbutton', ...
'Units','normalized', ...
'Position',[left yPos+btnHt+2*spacing btnWid btnHt], ...
'String',labelStr, ...
'Enable', 'off', ...
'Callback',callbackStr);
%====================================
% The CLOSE button
labelStr='Close';
callbackStr='close(gcf)';
yPos=frmBottom;
closeHndl=uicontrol( ...
'Parent', DctDemoFig, ...
'BusyAction','Queue','Interruptible','off',...
'Style','pushbutton', ...
'Units','normalized', ...
'Position',[left yPos btnWid btnHt], ...
'Enable', 'off', ...
'String',labelStr, ...
'Callback',callbackStr);
%===================================
% Set up scroll bar
barHeight = 0.04; barWidth = btnWid+2*frmBorder; sbBottom=top-barWidth-2*barHeight;
sbLeft = left-frmBorder;
textWidth = barWidth/2;
scrollPos = [sbLeft sbBottom barWidth barHeight];
callbackStr = 'dctdemo(''SliderUpdate'')';
hSlider = uicontrol( ...
'Parent', DctDemoFig, ...
'BusyAction','Queue','Interruptible','off',...
'Style','slider', ...
'Units','normalized', ...
'Position',scrollPos, ...
'V alue',64, ...
'min',1, ...
'max',64, ...
'Interruptible','off', ...
'Callback',callbackStr);
% Left and right range indicators
c = get(DctDemoFig,'Color');
if [.298936021 .58704307445 .114020904255]*c'<.5,
fgColor = [1 1 1];
else
fgColor = [0 0 0];
end
rangePos = [sbLeft sbBottom-barHeight textWidth barHeight];
uicontrol( ...
'Parent', DctDemoFig, ...
'BusyAction','Queue','Interruptible','off',...
'Style','text', ...
'Units','normalized', ...
'Position',rangePos, ...
'Horiz','left', ...
'Background',c, ...
'Foreground',fgColor, ...
'String','1');
rangePos = [sbLeft+barWidth/2 sbBottom-barHeight textWidth barHeight]; uicontrol( ...
'Parent', DctDemoFig, ...
'BusyAction','Queue','Interruptible','off',...
'Style','text', ...
'Units','normalized', ...
'Position',rangePos, ...
'Horiz','right', ...
'Background',c, ...
'Foreground',fgColor, ...
'String','64');
%==================================
% Set up the 8-by-8 array
aPos = [sbLeft sbBottom+2*barHeight barWidth barWidth];
h8x8 = axes( ...
'Parent', DctDemoFig, ...
'BusyAction','Queue','Interruptible','off',...
'Units','normalized', ...
'Position',aPos, ...
'XTick',[],'YTick',[], ...
'ydir','reverse', ...
'Box','on', ...
'xlim',[1 9], ...
'ylim',[1 9], ...
'clim',[0 1], ...
'DataAspectRatio',[1 1 1]);
title('DCT coefficients');
rangePos = [sbLeft sbBottom-(2*barHeight) barWidth barHeight]; hNumCof = uicontrol( ...
'Parent', DctDemoFig, ...
'BusyAction','Queue','Interruptible','off',...
'Style','text', ...
'Units','normalized', ...
'Position',rangePos, ...
'Horiz','center', ...
'Background',c, ...
'Foreground',fgColor, ...
'String','64 Coefficients Selected');
hSurf=surface(ones(9),'facecolor','flat','edgecolor',[1 1 1]/2);
%==================================
% Set up the image axes
figpos = get(DctDemoFig, 'position');
row = figpos(4); col = figpos(3); % dimensions of figure window
vertSpac = (row-256)/3;
leftSide = col*0.25;
axHeight = (top-(frmBottom-frmBorder)-3*spacing)/2;
hTrueAx = axes('Parent', DctDemoFig, ...
'units', 'pixels', ...
'BusyAction','Queue','Interruptible','off',...
'ydir', 'reverse', ...
'XLim', [.5 128.5], ...
'YLim', [.5 128.5],...
'CLim', [0 1], ...
'XTick',[],'YTick',[], ...
'Position', [leftSide row-118-vertSpac 128 128]);
hNewAx = axes('Parent', DctDemoFig, ...
'units', 'pixels', ...
'BusyAction','Queue', 'Interruptible','off',...
'ydir', 'reverse', ...
'XLim', [.5 128.5], ...
'YLim', [.5 128.5],...
'CLim', [0 1], ...
'XTick',[],'YTick',[], ...
'Position', [leftSide-64-vertSpac/2 vertSpac+10 128 128]); title('Reconstructed Image');
hErrAx = axes('Parent', DctDemoFig, ...
'units', 'pixels', ...
'BusyAction','Queue', 'Interruptible','off',...
'ydir', 'reverse', ...
'XLim', [.5 128.5], ...
'YLim', [.5 128.5],...
'CLim', [-0.8 0.8], ...
'XTick',[],'YTick',[], ...
'Position', [leftSide+64+vertSpac/2 vertSpac+10 128 128]); title('Error Image');
blank = repmat(uint8(0),128,128);
hTrueImage = image('Parent', hTrueAx,...
'CData', blank, ...
'BusyAction','Queue','Interruptible','off',...
'CDataMapping', 'scaled', ...
'Xdata', [1 128],...
'Ydata', [1 128],...
'EraseMode', 'none');
hNewImage = image('Parent', hNewAx,...
'CData', blank, ...
'BusyAction','Queue','Interruptible','off',...
'CDataMapping', 'scaled', ...
'Xdata', [1 128],...
'Ydata', [1 128],...
'EraseMode', 'none');
hErrImage = image('Parent', hErrAx,...
'CData', blank, ...
'BusyAction','Queue','Interruptible','off',...
'CDataMapping', 'scaled', ...
'Xdata', [1 128],...
'Ydata', [1 128],...
'EraseMode', 'none');
% Status bar
% rangePos = [64 3 280 15];
rangePos = [0 .01 .75 .05];
hStatus = uicontrol( ...
'Parent', DctDemoFig, ...
'BusyAction','Queue','Interruptible','off',...
'Style','text', ...
'Units','normalized', ...
'Position',rangePos, ...
'Horiz','center', ...
'Background',c, ...
'Foreground',[.8 0 0], ...
'Tag', 'Status', ...
'String','Status bar');
set(DctDemoFig,'defaultaxesposition',[0.10 0.1 0.60 0.85]) setstatus(DctDemoFig, 'Initializing DCT Demo...');
set(DctDemoFig, 'Pointer', 'watch');
drawnow
set(DctDemoFig, 'V isible','on');
% Put handles to graphics objects and controls in the figure's userdata handles.TrueImage = hTrueImage;
handles.NewImage = hNewImage;
handles.ErrImage = hErrImage;
handles.h8x8 = h8x8;
handles.Slider = hSlider;
handles.Surface = hSurf;
handles.ImgPop = hImgPop;
handles.TrueAx = hTrueAx;
handles.NewAx = hNewAx;
handles.ErrAx = hErrAx;
handles.NumCof = hNumCof;
set(DctDemoFig, 'UserData',handles);
LoadNewImage(DctDemoFig) % Load the defaultImage
set(DctDemoFig, 'HandleV isibility','Callback')
set([helpHndl closeHndl applyHndl], 'Enable', 'on'); return
%%%
%%% Sub-Function - LoadNewImage
%%%
function LoadNewImage(DemoFig)
if nargin<1
DemoFig = gcbf;
end
set(DemoFig,'Pointer','watch');
handles = get(DemoFig,'UserData');
hTrueImage=handles.TrueImage;
hTrueAx = handles.TrueAx;
hImgPop=handles.ImgPop;
v = get(hImgPop,{'value','String'});
name = deblank(v{2}(v{1},:));
setstatus(DemoFig, ['Loading the ' name ' image...']); drawnow
switch name
case 'Saturn',
saturn = []; % Make sure saturn is parsed as a variable load imdemos saturn
img = double(saturn)/255;
set(handles.ErrAx, 'CLim', [-.45 .45]);
case 'Pout',
pout = []; % Make sure pout is parsed as a variable
load imdemos pout
img = double(pout)/255;
set(handles.ErrAx, 'CLim', [-.33 .33]);
case 'Trees',
trees = []; % Make sure trees is parsed as a variable
load imdemos trees
img = double(trees)/255;
set(handles.ErrAx, 'CLim', [-.72 .72]);
case 'Quarter',
quarter = []; % Make sure quarter is parsed as a variable load imdemos quarter
img = double(quarter)/255;
set(handles.ErrAx, 'CLim', [-.52 .52]);
case 'Circuit',
circuit = []; % Make sure circuit is parsed as a variable load imdemos circuit
img = double(circuit)/255;
set(handles.ErrAx, 'CLim', [-.32 .32]);
case 'Flower',
flower = []; % Make sure flower is parsed as a variable load imdemos flower
img = double(flower)/255;
set(handles.ErrAx, 'CLim', [-.65 .65]);
otherwise
error('DCTDEMO: Unknown Image Option!');
end
set(hTrueImage, 'Cdata', img);
set(get(hTrueAx,'title'),'string',['Original ' name ' Image']); blank = repmat(uint8(0),128,128);
set(handles.NewImage, 'Cdata', blank);
set(handles.ErrImage, 'Cdata', blank);
drawnow
MakeDCT(DemoFig);
SliderUpdate(DemoFig);
ApplyDCT(DemoFig);
return;
%%%
%%% Sub-Function - ApplyDCT
%%%
function ApplyDCT(DemoFig)
if nargin<1
DemoFig = gcbf;
end
set(DemoFig,'Pointer','watch'); drawnow;
handles = get(DemoFig,'UserData');
hTrueImage=handles.TrueImage;
hNewImage=handles.NewImage;
hErrImage=handles.ErrImage;
h8x8=handles.h8x8;
hSlider=handles.Slider;
hSurf=handles.Surface;
hNewAx=handles.NewAx;
imageDCT = get(hTrueImage,'UserData');
if isempty(imageDCT), % Create imageDCT and DCTvar if necessary MakeDCT(DemoFig);
imageDCT = get(hTrueImage,'UserData');
end
setstatus(DemoFig, 'Reconstructing image...'); drawnow;
% Determine 8x8 image
mask = get(hSurf,'Cdata');
% Get 8-by-8 DCT matrix
dctm = dctmtx(8);
% Reconstruct and display image
% newImage = blkproc(imageDCT,[8 8],'idct2(x.*P1)',mask(1:8,1:8)); newImage = blkproc(imageDCT,[8 8],'P1*(x.*P2)*P3',dctm.',mask(1:8,1:8),dctm);
set(hNewImage, 'CData', newImage);
drawnow;
trueImage = getimage(hTrueImage);
set(hErrImage, 'CData', trueImage-newImage);
error = (trueImage.^2 - newImage.^2);
MSE = sum(error(:))/prod(size(trueImage));
if MSE 1e-10, MSE=0; end
msestr = sprintf('%0.3g',MSE);
setstatus(DemoFig, ['The MSE (with images normalized) is ' msestr ' .']);
set(DemoFig,'Pointer','arrow'); drawnow
return
%%%
%%% Sub-function - SliderUpdate
%%%
function SliderUpdate(DemoFig)
if nargin<1
DemoFig = gcbf;
end
handles = get(DemoFig,'UserData');
hSlider=handles.Slider;
hSurf=handles.Surface;
axHndl=gca;
v = 64-floor(get(hSlider,'value'));
var = get(hSlider,'UserData');
DCTvar = var(1,:); order = var(2,:);
if v==63
set(handles.NumCof, 'String', [int2str(64-v) ' Coefficient Selected']);
else
set(handles.NumCof, 'String', [int2str(64-v) ' Coefficients Selected']);
end
% Determine 8x8 image and display
mask = ones(8,8); mask(order(1:v)) = zeros(1,v);
im8x8 = zeros(9,9); im8x8(1:8,1:8) = mask;
set(hSurf,'Cdata',im8x8);
% If we were called from LoadNewImage (nargin==1), we don't want to update the % status bar since LoadNewImage will call ApplyDCT anyways.
if nargin==0
setstatus(DemoFig, 'Press ''Apply'' to reconstruct the image.');
end
return
%%%
%%% Sub-Function - MakeDCT
%%%
%%% This function is not used as a call back. It is called
%%% only from LoadNewImage.
function MakeDCT(DemoFig)
setstatus(DemoFig, 'Computing DCT...');
handles = get(DemoFig,'UserData');
hTrueImage=handles.TrueImage;
hNewImage=handles.NewImage;
h8x8=handles.h8x8;
hSlider=handles.Slider;
hSurf=handles.Surface;
I = getimage(hTrueImage);
% Get 8-by-8 DCT matrix
dctm = dctmtx(8);
% Old call: imageDCT = blkproc(I,[8 8],'dct2'); imageDCT = blkproc(I,[8 8],'P1*x*P2',dctm,dctm.'); DCTvar = im2col(imageDCT,[8 8],'distinct').';
n = size(DCTvar,1);
DCTvar = (sum(DCTvar.*DCTvar) - (sum(DCTvar)/n).^2)/n; [dum,order] = sort(DCTvar);
set(hTrueImage,'UserData',imageDCT)
set(hSlider,'UserData',[DCTvar;order])
return。