www.pudn.com > SPIHT_bandelet.rar > func_SPIHT_Main.m

%function func_SPIHT_bandelet_Main 
% Matlab implementation of SPIHT (without Arithmatic coding stage) 
% 
% Main function  
% 
% input:    Orig_I : the original image. 
%           rate : bits per pixel 
% output:   img_spiht  
% 
% Jing Tian 
% Contact me : scuteejtian@hotmail.com 
% This program is part of my undergraduate project in GuangZhou, P. R. China. 
% April - July 1999 
clear all; 
clc; 
warning off; 
%-----------   Input   ---------------- 
Orig_I = func_ReadRaw('lena512.raw', 512*512 , 512 , 512); 
%Orig_I =rgb2gray(imread('lena256.bmp')); 
%Orig_I=mat2gray(Orig_I1); 
figure;imshow(Orig_I,[]); 
 
M=Orig_I; 
rate =0.5; 
 
%-----------   Pre-processing   ---------------- 
 
OrigSize = size(Orig_I, 1); 
max_bits = floor(rate * OrigSize^2); 
 
OutSize = OrigSize; 
image_spiht = zeros(size(Orig_I)); 
% "image " is the input of codec 
[nRow, nColumn] = size(Orig_I); 
 
 
%-----------  bandelet Decomposition   ---------------- 
options.use_single_qt = 1;  
options.wavelet_vm=[4,4]; 
%set to: 
%   1 if you want the 3 wavelets bands H/V/D to share the same quadtree 
%   0 for 3 different quadtree (more adaptivity but more bits for 
%   geometry).% set to: 
[m,n] = size(Orig_I); 
s=Inf; 
Jmin =4; 
Jmax = log2(n)-1; 
j_min = 3; 
j_max = 4; 
level=Jmax-Jmin+1; 
T =0.9; 
%rates=[ 0.2, 0.25, 0.4, 0.5, 1 ]; 
 
block_size=2; 
 
%for  i=1:5 
%   rate=rates(i); 
    disp('computing quadtree'); 
    [QT,Theta] = compute_wavelet_quadtree_xhl(M,Jmin,T,j_min,j_max,s, options,'biorthogonal'); 
    [MB,Rbandg] = perform_wavelet_bandelet_transform_xhl(M,Jmin,QT,Theta,1, options,'biorthogonal'); 
 
    %-----------   Coding   ---------------- 
    img_enc = func_SPIHT_Enc(MB, max_bits, nRow*nColumn, level);   
    [p,q]=size(QT); 
    temp_QT1=zeros(1,p*q/8); 
    QT1=reshape(QT,p*q/8,8); 
    for i=1:p*q/8 
        for j=1:8 
          if QT1(i,j)==3 
              temp_QT1(1,i)=bitset(temp_QT1(1,i),j,0); 
          else 
              temp_QT1(1,i)=bitset(temp_QT1(1,i),j,1); 
          end 
        end 
    end 
            
%     temp_QT= rle(QT); 
%     %temp_Theta=btc_encode(Theta,2); 
%     [tmp1,Mu,Mi]=btc_encode(Theta,block_size); 
    %-----------   Decoding   ---------------- 
    img_dec = func_SPIHT_Dec(img_enc); 
    QT2=zeros(p*q/8,8); 
    
     for i=1:p*q/8 
        for j=1:8 
          if bitget(temp_QT1,j)==0 
           QT2(i,j)=3; 
          else 
           QT2(i,j)=4; 
          end   
        end 
    end 
   QT3=reshape(QT2,p,q); 
        
%   QT1= rle(temp_QT); 
    %Theta1=btc_image_decode(temp_Theta,2,m,n); 
%   Theta1=btc_decode(tmp1,Mu,Mi,block_size); 
    %-----------  bandelet Reconstruction   ---------------- 
    MB1=img_dec; 
    disp('Computing inverse bandelet transform.'); 
    [Mb,Rbandg] = perform_wavelet_bandelet_transform_xhl(MB1,Jmin,QT,Theta,-1, options,'biorthogonal'); 
    img_spiht=Mb; 
    %-----------   PSNR analysis   ---------------- 
    Q = 255; 
    OrigSize = nRow*nColumn; 
    MSE = sum(sum((img_spiht - Orig_I) .^ 2)) / OrigSize; 
    psnr = 10*log10(Q*Q/MSE); 
    % --------------compute the rate of the image compression ---------- 
     
    %total_bit=max_bits+size 
    %----------------------------------------------------- 
    figure; 
    imshow(img_spiht,[]); 
    title(['imge spiht psnr=',num2str(psnr),' rate=',num2str(rate),' T=',num2str(T)]); 
%end