
% *************************************************************************
%       GLAVNI PROGRAM ZA DEMONSTRACIJU RADA EZW KODERA I DEKODERA
%        Davor Zivko i Dario Zrno, NMDOS projekt, sijecanj 2005.
% *************************************************************************

% ucitavanje putanje slike
string = input('path -> ','s');

% ucitavanje slike    
[img map]=imread(string);
 

% konverzija RGB slike u sivu
img2 = rgb2gray(img);

% prikaz originalne slike
figure,imshow(img2),title('Originalna slika');

% SKALIRANJE -> kako bi se mogli postici i P i N simboli
X = double(img2) - 128;
X0=X;

% maksimalni elementi po stupcima
Y0=max(X);

% maksimalni element unutar cijele matrice
Y1=max(Y0); 


for i=0:20;                          % maksimalno 20 pragova (ako zadovoljavaju uvjete)
    if 2^i<=Y1 & 2^i>0.5*Y1;
        threshold=2^i;   
        initialthreshold=threshold; 
        %zapamti broj razina
        laststeplevel=i+1;
        break;
    end;
end;

sublist=[];
sub_list=[];

% dimenzije slike
[xx,yy]=size(X);

% mapping kreira Morton scale order matrix -> redoslijed skeniranja
A=mapping(xx);

% velicina Mortonove matrice
[m,n]=size(A);

% postavljanje vektora sa redoslijedom kodiranja kao globalne varijable
global N;

N=zeros(m*n,2);

% upisivanje indeksa skeniranja
for i=1:m,
    for j=1:n,
        N(A(i,j),1)=i;
        N(A(i,j),2)=j;
    end
end

% inicijalizacija reda skeniranja
order=1;  

% dok je prag razlicit od 0.5, obavlja se dominant-pass i subordinant-pass
while threshold ~= 0.5, 

    threshold

   % Dominant Pass
   [D,X,sublist,sub_list] = dominantpass(X,threshold,sublist,sub_list); 
   DD{order}=D

   significantlist{order}=sub_list;
   
   % Subordinate pass
   % prepolavlja se prag
   threshold=threshold/2;
   if threshold ==0.5,
       break;
   end
   
   % zapisuje preko subordinate-pass, u SS, preko S, je li neki koeficijent
   % jednak pragu ili nije
   S = subordinatepass(sublist,threshold);
   SS{order}=S
   order=order+1;
end


%**************************************************************************
%                          EZW dekoder
%**************************************************************************

global N;

% velicina inicijalne slike
[m,n]=size(N);

% inicijalizacija prostora i pocetnog praga
XX=zeros(sqrt(m)); 
threshold=initialthreshold; 

% kreiranje matrice indeksa pozicija za sve p i n koeficijente
sublist=[]; 
for level=1:laststeplevel, 
    RR=zeros(size(XX)); 
    [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
                % postavlja sve nasljednike nul-stabla na nulu
                RR=checkchildren(i,RR);
            end
            RR(N(i,1),N(i,2))=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

initialimage=X0;

reconstructedimage=XX;


razlika = XX - X0;

% prikaz rekonstruirane slike i razlike izmedju originala i rekonstrukcije
figure,imshow(uint8(reconstructedimage+128)),title('Rekonstruirana slika');
figure,imshow(uint8(razlika)),title('Razlika originalne i rekonstruirane slike');