%% Trénování SVM klasifikátoru a predikce tříd - RITE databáze
% Skript slouží pro natrénování klasifikátoru SVM a následnou predikci
% tříd cév.
% Skript je po natrénování klasifikátoru rozdělen do dvou částí, první část
% vyhodnocuje výsledek klasifikace na "raw" testovací množině. Druhá část
% skriptu slouží pro zpracování výsledku pro určitý testovací obraz, kde je
% zahrnuta i vizualizace výsledku.

clc; clear;

%% Načtení trénovacích dat
load('rite_train_label_mat_skeleton.mat');
load('rite_train_priznak_mat_skeleton.mat');
load('rite_train_label_mat_vess.mat');
load('rite_train_priznak_mat_vess.mat');

%% Načtení již natrénovaného klasifikátoru SVM
load('SVMModel_rbf.mat');
 
%% Trénování SVM
t = templateSVM('Standardize', 1, 'KernelFunction', 'rbf');
SVMModel_rbf = cell(1, 2);
SVMModel_rbf{1} = fitcecoc(rite_train_priznak_mat_skeleton, rite_train_label_mat_skeleton, 'Learners', t, 'FitPosterior', 1, 'ClassNames', [2, 3]);
SVMModel_rbf{2} = fitcecoc(rite_train_priznak_mat_vess, rite_train_label_mat_vess, 'Learners', t, 'FitPosterior', 1, 'ClassNames', [2, 3]);

%% Pro celou množinu testovacích dat
% Načtení testovacích dat
load('rite_test_label_mat_skeleton.mat');
load('rite_test_priznak_mat_skeleton.mat');
load('rite_test_label_mat_vess.mat');
load('rite_test_priznak_mat_vess.mat');

% Predikce tříd cév pomocí natrénovného SVM modelu
for i = 1:2
    [label{i}, score{i}] = predict(SVMModel_rbf{i}, rite_test_priznak_mat_vess);
end
    
% load('label.mat');
% load('score.mat');

% Pro RITE databázi - v labelech zahruny i cévy, kt. nebyly klasifikovány
% expertem, tyto pixely zavádí zkreslení výsledku a snižuje úspěšnost. Tyto
% neklas. pixely jsou odstraněny.
label{1}(rite_test_label_mat_skeleton == 1) = [];
label{2}(rite_test_label_mat_vess == 1) = [];
rite_test_label_mat_skeleton(rite_test_label_mat_skeleton == 1) = [];
rite_test_label_mat_vess(rite_test_label_mat_vess == 1) = [];

% Výpočet statistických metrik klasifikace
acc_all = double([]);
error_rate_all = double([]);
spravne = 0;
nespravne = 0;
for i = 1:length(label)
    if (label{2}(i) == rite_test_label_mat_vess(i)) == 1
        spravne = spravne + 1;
    else
        nespravne = nespravne + 1;
    end
end

acc_pom_all = spravne / (length(label{2})) * 100;            % procento přesnosti klasifikace
acc_all = [acc_all; acc_pom_all];

error_rate_pom_all = nespravne / (length(label{2})) * 100;   % procento chybovosti klasifikace
error_rate_all = [error_rate_all; error_rate_pom_all];
   
% vytvoření matice záměn
C_all = confusionmat(rite_test_label_mat_vess, label{2}, 'Order', [2, 3]);
CCH_all = confusionchart(C_all);

%% Pro jednotlivé testovací obrazy
% Načtení dat daného obrazu
load('rite_test_priznak_mat_skeleton1.mat');
load('rite_test_priznak_mat_vess1.mat');

% segmentace cév
seg_location = 'C:\Users\jana-\Documents\FEKT\DP\data\RITE\test\vessels';   % cesta k maskám segmentovaných cév
addpath(seg_location);
seg_data = imageDatastore(seg_location);

% labely cév
class_location = 'C:\Users\jana-\Documents\FEKT\DP\data\RITE\test\classes'; % cesta k labelům cév označneých experty
addpath(class_location);
class_data = imageDatastore(class_location);

% body OD
mask_location = 'C:\Users\jana-\Documents\FEKT\DP\data\RITE\test\masks';    % cesta k ROI maskám
addpath(mask_location);
mask_data = fileDatastore(mask_location, 'ReadFcn', @load, 'FileExtensions', '.mat');

% načtení snímků a masek 
segmentace = im2double(imread('rite08_test_seg.png'));
klas = imread('rite08_test_class.png');
load('rite08_test_mask.mat');

%% Predikce tříd naučeného SVM modelu
[label_skel, ~] = predict(SVMModel_rbf{1}, rite_test_priznak_mat_skeleton8);
[label_vess, ~] = predict(SVMModel_rbf{2}, rite_test_priznak_mat_vess8);

%% Vytvoření skeletonu, vymezení oblasti zájmu jen na cévy
BW = imbinarize(segmentace); % binarizace segmentovaného obrazu
skeleton = bwskel(BW);       % vytvoření skeletonu cévního stromu
branch = bwmorph(skeleton, 'branchpoints');
skeleton(branch) = 0;        % odstranění bifurkací a bodů křížení

ROI_cevy = immultiply(segmentace, maska_rozdil);  % násobení segmentační masky a ROI okolo OD
ROI_skeleton = immultiply(ROI_cevy, skeleton);    % násobení ROI a skeletonu

%% VESS predikce - převod vektoru do matice
ROI_cevy_pom = reshape(ROI_cevy, [], 1);
pom_lab_vess = zeros(size(ROI_cevy));
pom_lab_vess = reshape(pom_lab_vess, [], 1);
pom_ind1 = 1;
for m = 1:length(pom_lab_vess)
    if ROI_cevy_pom(m) == 1
        pom_lab_vess(m) = label_vess(pom_ind1);
        pom_ind1 = pom_ind1 + 1;
    end
end
pom_lab_vess = reshape(pom_lab_vess, size(ROI_cevy, 1), size(ROI_cevy, 2));

%% SKELETON predikce - převod vektoru do matice
ROI_skeleton_pom = reshape(ROI_skeleton, [], 1);
pom_lab_skel = zeros(size(ROI_skeleton));
pom_lab_skel = reshape(pom_lab_skel, [], 1);
pom_ind2 = 1;
for n = 1:length(pom_lab_skel)
    if ROI_skeleton_pom(n) == 1
        pom_lab_skel(n) = label_skel(pom_ind2);
        pom_ind2 = pom_ind2 + 1;
    end
end
pom_lab_skel = reshape(pom_lab_skel, size(ROI_skeleton, 1), size(ROI_skeleton, 2));
        
%% SPOJENÍ VESS A SKELETONU
pom_final = pom_lab_vess;
for k = 1:size(pom_final, 1)
    for l = size(pom_final, 2)
        if pom_lab_skel(k, l) ~= 0
            pom_final(k, l) = pom_lab_skel(k, l);
        end
    end
end

%% Vizualizace výsledku klasifikace
vysledek_klas_R = size(pom_final);
vysledek_klas_G = size(pom_final);
vysledek_klas_B = size(pom_final);
for s = 1:size(pom_final, 1)
    for t = 1:size(pom_final, 2)
        if pom_final(s, t) == 0             % pixely pozadí -> černá
            vysledek_klas_R(s, t) = 0;
            vysledek_klas_G(s, t) = 0;
            vysledek_klas_B(s, t) = 0;
        elseif pom_final(s, t) == 2         % pixely tepen -> červená
            vysledek_klas_R(s, t) = 255;
            vysledek_klas_G(s, t) = 0;
            vysledek_klas_B(s, t) = 0;
        elseif pom_final(s, t) == 3         % pixely žil -> modrá
            vysledek_klas_R(s, t) = 0;
            vysledek_klas_G(s, t) = 0;
            vysledek_klas_B(s, t) = 255;
        end
    end
end      
   
vysledek_klas = cat(3, vysledek_klas_R, vysledek_klas_G, vysledek_klas_B); % vytvoření RGB obrazu  

% úprava klasifikovaných labelů expertů - odstranění neklasifikovaných, zavádějí chybu klasifikace
klas_pom = im2double(klas);
[r_pom, g_pom, b_pom] = imsplit(klas_pom);
r_pom = immultiply(r_pom, maska_rozdil);
r_pom = immultiply(r_pom, segmentace);
g_pom = immultiply(g_pom, maska_rozdil);
g_pom = immultiply(g_pom, segmentace);
b_pom = immultiply(b_pom, maska_rozdil);
b_pom = immultiply(b_pom, segmentace);
klas_viz = cat(3, r_pom, g_pom, b_pom);

figure
subplot 121; imshow(klas_viz, []); title('Ground truth'); subplot 122; imshow(vysledek_klas, []); title('Klasifikace pomocí SVM');

%% Výpočet statistických metrik
% výpočet správně a špatně klasifikovaných pixelů
pom_final = labels_rite(vysledek_klas);
label_test = labels_rite(klas_viz);

label_test(label_test == 1) = 0;
pom_final(label_test == 0) = 0;

pom_final = reshape(pom_final, [], 1);
pom_final(pom_final == 0) = [];

label_test = reshape(label_test, [], 1);
label_test(label_test == 0) = [];

acc = double([]);
error_rate = double([]);
spravne = 0;
nespravne = 0;
for i = 1:length(pom_final)
    if (pom_final(i) == label_test(i)) == 1
        spravne = spravne + 1;
    else
        nespravne = nespravne + 1;
    end
end

acc_pom = spravne / (length(pom_final)) * 100;           % procento přesnosti klasifikace
acc = [acc; acc_pom];

error_rate_pom = nespravne / (length(pom_final)) * 100;  % procento chybovosti klasifikace
error_rate = [error_rate; error_rate_pom];
   
% vytvoření matice záměn
C = confusionmat(label_test, pom_final, 'Order', [2, 3]);
CCH = confusionchart(C);