% This is a test program that tests the EZW encoder and decoder

% matrix X is an example matrix same as Shapiro's paper given





% Written by xiangang li

%

% This function calls:

%    dominantpass.m

%    subordinatepass.m

%    checkdescents1.m

%    checkchildren.m

%    mapping.m







X=[63 -34 49 10 7 13 -12 7;

   -31 23 14 -13 3 4 6 -1;

   15 14 3 -12 5 -7 3 9;

   -9 -7 -14 8 4 -2 3 2;

   -5 9 -1 47 4 6 -2 2;

   3 0 -3 2 3 -2 0 4;

   2 -3 6 -4 3 6 3 6;

   5 11 5 6 0 3 -4 4];



X0=X;



Y0=max(X);

     Y1=max(Y0); 

     for i=0:20;

        if 2^i<=Y1 & 2^i>0.5*Y1;

           threshold=2^i;   % get initial threshold T0;

           initialthreshold=threshold; % get initial threshold T0;

           laststeplevel=i+1;% last step level

           break;

        end;

     end;

     

     

     sublist=[];

     sub_list=[];

     [xx,yy]=size(X);

     

  

A=mapping(xx);

[m,n]=size(A);

global N;      % Let Morton scanorder vector as a global variable

 N=zeros(m*n,2);



     for i=1:m,

        for j=1:n,

            

                 N(A(i,j),1)=i;

                 N(A(i,j),2)=j;

                

         end

     end



       order=1;  

while threshold ~= 0.5,  % if threshold~=1, do dominantpass and subordinatepass.

         

   threshold

   

   %Dominant Pass

   

   [D,X,sublist,sub_list] = dominantpass(X,threshold,sublist,sub_list); 

   

   DD{order}=D

   significantlist{order}=sub_list;

     

    %Subordinate pass

    

    threshold=threshold/2;

    if threshold ==0.5,

        break;

    end

    

    

   S = subordinatepass(sublist,threshold);

   SS{order}=S

   order=order+1;

    

end



  

%*********************************%   

%*********EZW decoder*************%

%*********************************%



 global N;

 [m,n]=size(N);% the size of initial image

               % m is the pixels of initial image

 XX=zeros(sqrt(m)); % initialize the reconstructured image to zero;

 threshold=initialthreshold; % initial theshold ;

 

 sublist=[]; % sublist is the new position matrix 

             % for all significant coefficients 'p' and 'n';

             

             

for level=1:laststeplevel, 

                

            RR=zeros(size(XX)); % reference matrix RR;

  

     [a,b]=size(DD{level}); % ?    

 

 % dominant pass

 

     

        i=1; j=1;

     while j<=b,

     

         if RR(N(i,1),N(i,2))==0

                      

              if DD{level}(j)=='p'

                  if threshold==1

                     XX(N(i,1),N(i,2))=threshold; 

                 else

                     XX(N(i,1),N(i,2))=1.5*threshold;

                 end

                

              end

     

              if DD{level}(j)=='n'

                 if threshold==1

                     XX(N(i,1),N(i,2))=-threshold; 

                 else

                     XX(N(i,1),N(i,2))=-1.5*threshold;

                 end

              end

         

                

              if DD{level}(j)=='t'& A(N(i,1),N(i,2))<=m/4

                 RR=checkchildren(i,RR);% all zerotree's descendants are set to 1.

              end

             

              RR(N(i,1),N(i,2))=1; %reference matrix =1;

              i=i+1;

              j=j+1;

               

        else i=i+1;

               

             

        end

         

    end

 

 

 % subordinate pass

 

 

 [xx,yy]=size(significantlist{level});

 threshold=threshold/2;

 

   for i=1:xx,

       

       if level==laststeplevel|threshold==0.5

           break

       end

       

      if SS{level}(i)==1

         if XX(sub_list(i,1),sub_list(i,2))>0;

           

            XX(sub_list(i,1),sub_list(i,2))= fix(XX(sub_list(i,1),sub_list(i,2))+ threshold/2);

            

         else 

             

            XX(sub_list(i,1),sub_list(i,2))= fix(XX(sub_list(i,1),sub_list(i,2))-threshold/2);

           

         end

      end

      

      if SS{level}(i)==0

         if XX(sub_list(i,1),sub_list(i,2))>0;

           

            XX(sub_list(i,1),sub_list(i,2))= fix(XX(sub_list(i,1),sub_list(i,2))-threshold/2);

            

         else 

             

            XX(sub_list(i,1),sub_list(i,2))= fix(XX(sub_list(i,1),sub_list(i,2))+threshold/2);

         

         end

      end

   end

   

    



      threshold

      level

      XX

      

end

initialiamge=X0

reconstructedimage=XX

difference=XX-X0