function segmentace_final = ivd_segmentation(magnitude, image, spine_curve, pozice_ivd_boxu)
% funkce provede segmentaci meziobratlove plotenky podle Zhu et al.:
% A method of localization and segmentation of intervertebral discs in
% spine MRI based on Gabor filter bank

obr_velikost = size(magnitude{1, 1});
meritka = size(magnitude, 1);
smery = size(magnitude, 2);

%% Segmentace meziobralových plotének
% Sumace obrazů pro každý z pěti scales
suma_meritek = zeros(obr_velikost(1), obr_velikost(2), meritka);
for i = 1:smery
    temp = sum(reshape(cell2mat(magnitude(:, i)'), [obr_velikost(1), obr_velikost(2), meritka]), 3);
    temp = func.trim_spine(temp, spine_curve);
    suma_meritek(:, :, i) = temp ./ smery;
end
% Vybrání maximální hodnoty ze všech pěti obrazů
maximum_meritek = max(suma_meritek, [], 3);

maximum_meritek_binary = maximum_meritek > graythresh(maximum_meritek(maximum_meritek>0));

% Inicializace prahovych hodnot a segmentaci
ivd_boxy = cell(1, size(pozice_ivd_boxu, 1));
ivd_segmentace = cell(1, length(ivd_boxy));
maximum_meritek_binary_boxy = cell(1, length(ivd_boxy));
for i = 1:size(ivd_segmentace, 2)
    right_bound = pozice_ivd_boxu(i, 2) + pozice_ivd_boxu(i, 4);
    bottom_bound = pozice_ivd_boxu(i, 1) + pozice_ivd_boxu(i, 3);
    indices_vert = pozice_ivd_boxu(i, 2):right_bound;
    indices_hor = pozice_ivd_boxu(i, 1):bottom_bound;
    
    ivd_boxy{i} = image(indices_vert, indices_hor);
    maximum_meritek_binary_boxy{i} = maximum_meritek_binary(indices_vert, indices_hor);
end

% Nastavení vstupů pro iterativní postup hledání prahové hodnoty
max_iterace = 30;
element = strel('rectangle', [8 12]); % pozor na prepocet
prah_step = 0.01;

vysledek = cell(1, length(ivd_boxy));
% Cyklus itarativního hledání lokálního prahu pro každý z IVD
for i = 1:length(vysledek)
    % Nastavení začátku iterace
    iterace = 0;
    
    prah = graythresh(ivd_boxy{i}); % prvotni prah
    while iterace < max_iterace
        % segmentace
        ivd_segmentace{i} = ivd_boxy{i} < prah;
        
%         imshow(ivd_segmentace{i}, []); pause
        
        A1 = abs(sum(ivd_segmentace{i}(:)) - sum(maximum_meritek_binary_boxy{i}(:)));
        A2 = sum(maximum_meritek_binary_boxy{i}(:)) / 2;
        T1 = min(A1, A2);
        T2 = max(A1, A2);
        
        ivd_segmentace{i} = imclose(ivd_segmentace{i}, element);
        
        area_before_filling = sum(ivd_segmentace{i}(:));
        
        ivd_segmentace{i} = imfill(ivd_segmentace{i},'holes');
        
        area_after_filling = sum(ivd_segmentace{i}(:));
        
        area_diff = abs(area_before_filling - area_after_filling);
        
        if area_diff >= T1 && area_diff <= T2
            break
        elseif abs(area_after_filling) > T2
            prah = prah - prah_step;
        else
            prah = prah + prah_step;
        end
        
        iterace = iterace + 1;
    end
    
    ivd_segmentace{i} = bwareafilt(ivd_segmentace{i}, 1, 'largest');
    
%     ivd_segmentace{i} = imclearborder(ivd_segmentace{i});
end

segmentace_final = func.superpozice_ivd_segmentace(obr_velikost, ivd_segmentace, pozice_ivd_boxu);

end

