DSPRelated.com
Forums

DCT & DWT

Started by seth5 July 17, 2006
This is a MATLAB code for DCT and DWT Compression. I tried this code with
512X512 TIF image. The DCT part worked fine, but I faced some problems
with DWT compression. The output image was just blank. Any help regarding
this will be highly appreciated!!



inp_image=imread('c:\zenacolor.tif');
inp_imaged=double(inp_image);
siz = size(inp_imaged) ;
no_row = siz(1) ;
no_col  = siz(2) ;
energy_image = sum(sum(inp_imaged.^2))        % Energy of the input image
%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% 2-D DCT Calculation and Experiment
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
dctcoeffd=zeros(siz(1),siz(2)) ; %%%%%%%%%%%%%%%%%%%555
zeros(no_row,no_col)
dctrearr=zeros(siz(1),siz(2)) ;
dctrearrs=zeros(siz(1),siz(2)) ;
rmse = zeros(8,8) ; 
nbr = no_row/8 ; 
nbc = no_col/8 ; 
for m=0:(no_row/8-1) %%use nbr 
   for n=0:(no_col/8-1)%%% use nbc
      x=m*8 ; y=n*8 ;
      block=inp_imaged(x+1:x+8,y+1:y+8); %%%%%%%%%%%%%%% change
      dctcoeffd(x+1:x+8,y+1:y+8)=dct2(block) ;%%%%%%%%%%%%%%% change
   end
end
%DCT calculation is over. Now rearranging the data
for m=0:7
   for n=0:7
      x=m*nbr ; y=n*nbc ;
      block =dctcoeffd(m+1:8:end,n+1:8:end) ;
      dctrearr(x+1:x+nbr,y+1:y+nbr) = block ;
      dctrearrs(x+1:x+nbr,y+1:y+nbr) = (1+25*(m+n))*block ;
      rmse(m+1,n+1)=sqrt(sum(sum(block.*block))/(nbr*nbc)) ;
   end
end
maxv = max(max(abs(dctrearr))) ;
%imshow(dctcoeffd) ;
dctrearrs = abs(dctrearrs) ;  % The DCTC > 0 now
imshow(dctrearrs/maxv) ;
imshow(dctrearrs(1:2*nbr,1:2*nbc)/maxv) ;
print -dtiff plot.tiff
%
fprintf('\nRMS Energy cof the DCT coefficients are as follows:\n') ;
round(rmse)
%print -dtiff plot.tiff
%
%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% 2-D DCT Calculation and Experiment Over
% Now 2-D DWT Calculation and Experiment
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
nstage   = 3 ;          % No. of Decomposition stage
wtype    = 1 ;          % Wavelet filter type, 0=Haar, 1=Daub-4

wtcoeffd  = dwt2(inp_imaged,no_row,no_col,nstage,wtype,1) ;

out_image = dwt2(wtcoeffd,no_row,no_col,nstage,wtype,-1) ;
wtcoeffdisp = zeros(no_row,no_col) ;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Calculating RMSE of the subbands.
for level=nstage:-1:1
   k=2^level ;
   res = nstage-level ;
   npix = no_row*no_col/(k^2) ;
   if (res==0)
      block = wtcoeffd(1:no_row/k, 1:no_col/k) ;
      rmse = sqrt(sum(sum((block).^2))/npix);   % LH subband
      fprintf('\nRMSE in subband-%d = %d', 3*(nstage-level),rmse) 
   end
   k=2*k;
   block1 = wtcoeffd(1:no_row/k, (no_col/k)+1:no_col/(k/2)) ;
   block2 = wtcoeffd(no_row/k+1:no_row/(k/2), 1:no_col/(k)) ;
   block3 = wtcoeffd(no_row/k+1:no_row/(k/2), no_col/k+1:no_col/(k/2)) ;
   
   
   %RMSE Calculation
   rmse1 = sqrt(sum(sum((block1).^2))/npix) ;  % LH subband
   rmse2 = sqrt(sum(sum((block2).^2))/npix) ;  % HL subband
   rmse3 = sqrt(sum(sum((block3).^2))/npix) ;  % HH subband
   fprintf('\nRMSE in subband-%d = %d', 3*res+1,rmse1) ;
   fprintf('\nRMSE in subband-%d = %d', 3*res+2,rmse2) ;
   fprintf('\nRMSE in subband-%d = %d', 3*res+3,rmse3) ;
   %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
   %%The following code scales up the coefficients for better viewing
   % The scaling is adhoc
   
   wtcoeffdisp(1:no_row/k, no_col/k+1:no_col/(k/2)) = (2^(res+1))*block1
;
   wtcoeffdisp(no_row/k+1:no_row/(k/2), 1:no_col/k) = (2^(res+1))*block2
;
   wtcoeffdisp(no_row/k+1:no_row/(k/2), no_col/k+1:no_col/(k/2)) =
(2.5^(res+1))*block3 ;

   k=k/2;
   if (res==0)
       wtcoeffdisp(1:no_row/k, 1:no_col/k) = block/5;
   end
   %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
end
energy_percent = (rmse^2)*(no_row*no_col/(4^nstage))*100/energy_image ; 
fprintf('\nEnergy compacted in the LL subband = %d\n', energy_percent) ;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%
% A large number of DWT coefficients are likely to exceed the 8-bit
% dynamic range i.e., 0-255. The following code will clip the
coefficients
% to [0,255] range
wtcoeffdisp = abs(wtcoeffdisp) ;  % The DCTC > 0 now
mask = wtcoeffdisp > 255 ;
wtcoeffdisp1 = wtcoeffdisp.*(~mask) + 255.* mask ; % all coeffs are now <=
255
% Converting the image to uint8 format
wtcoeffdisp8b = uint8(wtcoeffdisp1) ;
imshow(wtcoeffdisp8b) ;
print -dtiff plot.tiff
imshow(wtcoeffdisp8b(1:no_row,1:no_col)) ;
%print -dtiff plot.tiff
imwrite(wtcoeffdisp8b,'C:\wtcoeffdisp.tif','tif') ;