clear all
clc
close all

fname='C:\skola\pokus\data\Svobodova\'; %cesta k datum

index=1;

[maska0,maska15,maska30,maska45,maska60,maska75,maska90,maska105,maska120,maska135,maska150,maska165]= cevy(); %zisk masek pro tlusté cévy
[maska_st0,maska_st15,maska_st30,maska_st45,maska_st60,maska_st75,maska_st90,maska_st105,maska_st120,maska_st135,maska_st150,maska_st165] = cevy2(); %zisk masek pro středně tlusté cévy
[maska_t0,maska_t15,maska_t30,maska_t45,maska_t60,maska_t75,maska_t90,maska_t105,maska_t120,maska_t135,maska_t150,maska_t165] = cevy3(); %zisk masek pro tenke cévy
num=0; %celkove pixelu v obrazku 

for (i=1:10)
    %%obrazky
    if(i<10)
    fname2=[fname '0' num2str(i)];
    else
    fname2=[fname num2str(i)]; 
    end
    A=dir(fname2);
    
    %%labely
         if(i==1)
    fname3=['C:\skola\pokus\PixelLabelData'];
    else
    fname3=['C:\skola\pokus\PixelLabelData_' num2str(i-1)]; 
         end
    
          
    B=dir(fname3);  
    delka =length(B)-2;
    
    for(j=3:length(A))
        images{index}= imresize(double(imread([fname2 '\' A(j).name]))/65536,1/3);
        labels{index}= imresize(imread([fname3 '\Label_' num2str(j-2) '.png']),1/3);
        labels{index}= labels{index}>1;
        
        max_odezva{index}=filtrace(images{index},maska0,maska15,maska30,maska45,maska60,maska75,maska90,maska105,maska120,maska135,maska150,maska165); %konvoluce obrazu s jednotlivými maskami a zisk maximální odezvy
        max_odezva2{index}=filtrace(images{index},maska_st0,maska_st15,maska_st30,maska_st45,maska_st60,maska_st75,maska_st90,maska_st105,maska_st120,maska_st135,maska_st150,maska_st165);
        max_odezva3{index}=filtrace(images{index},maska_t0,maska_t15,maska_t30,maska_t45,maska_t60,maska_t75,maska_t90,maska_t105,maska_t120,maska_t135,maska_t150,maska_t165);

    E{index} = entropyfilt(images{index});
    S{index} = stdfilt(images{index},ones(9));
    R{index} = rangefilt(images{index},ones(9));
    
    num=num+numel(E{index}/3);%počet prvku v poli numel
    
    
        index=index+1;
    end
    
   
end
% % 
% pom=0;
% for(k=[20:30])
% pom=pom+sum(sum(labels{k}));
% end

%numsamples=num;
numsamples=27250000;
priznaky=zeros(numsamples,7); %%podle příznaků -> 7 mam E,S,R,image,priznaky
anotace=zeros(numsamples,2); %%anotace bud ceva nebo pozadi 
ind=1;

numbg = 0;
numfw = 0;

for(k=131:length(E)) %test
    
    for(i=1:size(E{k},1))
       
        for(j=1:size(E{k},2))
            if(ind>numsamples)
                break
            end
            
            priznaky(ind,7)=(E{k}(i,j));
            priznaky(ind,6)=(S{k}(i,j));
            priznaky(ind,3)=(R{k}(i,j));
            priznaky(ind,2)=(images{k}(i,j));
% %             
            priznaky(ind,1)=(max_odezva{k}(i,j));
            priznaky(ind,5)=(max_odezva2{k}(i,j));
            priznaky(ind,4)=(max_odezva3{k}(i,j));
            
            anotace(ind,2)=labels{k}(i,j);
            anotace(ind,1)=~labels{k}(i,j);
            ind=ind+1;
        end
        
        
    end
end

%     load('sit')
    FF=sit(priznaky'); %vložit do site
    FF=vec2ind(FF);
    
   
  BW_out=cell(size(E));
ind=1;
for(i=131:length(E)) %prochazim obrazky
    BW_out{i}=zeros(size(images{i},[1,2]));
    for(j=1:size(images{i},1))
        for(k=1:size(images{i},2))
            BW_out{i}(j,k)=FF(ind)-1; %% na vystupu segmentovaného obrazku 0 nebo 1
            ind=ind+1;
        end
    end
      BW_out{i}=medfilt2(BW_out{i},[8,8]);
      BW_out{i}= bwareaopen(BW_out{i},1250);
end

% figure
% % % subplot(131)
% imshow(images{131},[]);
% % title('původní obrazek')
% % subplot(132)
% figure
% imshow(labels{131},[]);
% % title('Ground Truth')
% % subplot(133)
% figure
% imshow(BW_out{131},[]);
% % % % title('Segmentace')
% figure
% imshowpair(BW_out{133},labels{133})

%% vysledky
Podobnost=[]; %dice skore

for i = 131:length(E) %for i = 1:length(E)
    Podobnost(i) = dice(labels{i},BW_out{i});
end

FN=0; %pro vypocet se, sp 
TN=0;
TP=0;
FP=0;

for(i=131:length(E)) 
  for r=1:(size(labels{i},1))
    for s=1:(size(images{i},2))
        if labels{i}(r,s) == 1 
            if BW_out{i}(r,s) == 0  
                FN = FN + 1;
            else    
                TP = TP + 1;
            end
        else    
            if BW_out{i}(r,s) == 1 
                FP = FP + 1;
            else    
                TN = TN + 1;
            end
        end
    end
  end
end
SE = TP/(TP+FN);
SP = TN/(FP+TN);

% A=zeros(2,2);
% 
% for i=(131:length(E)) %for i=(1:length(E))
%     mat=confusionmat(reshape(labels{i},1,[]),reshape(BW_out{i},1,[]));
%     % confusionmat produkuje tabulku predikci - [TN, FP; FN, TP]
%     A=A+mat;
% end
