%--------------------------------------------------------------------------
%------------------------- DIPLOMOVÁ PRÁCE --------------------------------
%--------------------------------------------------------------------------
% Tvarová klasifikace pro detekci chybně segmentovaných kostí v CT datech
%--------------------------------------------------------------------------
%                                                     Autor: Tomáš Janovič
%--------------------------------------------------------------------------

%% Inicializace
close all;
clear all;
clc;

%% Načtení objemu dat a předzpracování
try
    data = loaddicom(); %načetní dat
catch ex
    error(ex.message) %v případě, že se nepodaří, výpis vyjímky
end
data = mat2gray(im2double(data)); %přetypování na double a normalizace v intervalu <0,1>
otsu_prah = graythresh(data); %práh dle Otsu algoritmu
masky = zeros(size(data)); %alokace paměti pro masky hlavy
for i=1:size(data,3)
    masky(:,:,i) = headmask(data(:,:,i),otsu_prah); %zisk masky hlavy
end
masky = logical(masky); %přetypování masek na logický typ
data = data.*masky; %násobení masek s původními obrazy

%% Histogram celého objemu dat a algoritmus fitování prav. rozdělení
histogram = imhist(data(masky),1024); %výpočet histogramu
[prah_kost,funkce,funkce_fit] = cortical_bones_thresh(histogram,'pearson','exhaustive_search'); %optimalizace prav. rozdělení a zjištění prahu pro segmentaci

figure(1) %zobrazení histogramu a optimalizované funkce
bar(0:1/1023:1,histogram,'b')
hold on
plot(0:1/1023:1,funkce,'r','LineWidth',2)
plot(0:1/1023:1,funkce_fit,'green','LineWidth',2)
plot(prah_kost*ones(1,5e5),1:5e5,'black','LineWidth',2)
hold off
axis tight
title('Histogram celého objemu dat');
xlabel('Normalizovaná intenzita');
ylabel('Pocet voxelu');
legend('Histogram','Výchozí funkce','Optimalizovaná funkce','Práh');

%% Segmentace
kosti_kort = (data>=prah_kost); %segmentace kortikálních kostí

% figure(2) %možnost zobrazení segmentace
% for i=1:size(data,3)
%     imshow(kosti_kort(:,:,i))
%     pause(0.01)
% end

%% Detekce potenciálně chybně segmentovaných úseků
[pozice_preruseni,img_pozice] = detection_missegmented_points(kosti_kort); %samotná detekce potenciálně přerušených kostí

% figure(3) %možnost zobrazení detekce
% for i=1:size(data,3)
%     imshow(img_pozice(:,:,:,i))
%     pause(0.1)
% end








