function detection(SelImage)

% sluzi na detekciu castic. Uzivatel prvym kruhom ktory sa zobrazi vo novom
% okne definuje pozadie a druhym kruhom definuje nejaku referencnu
% casticu (polohu nastavi ako okolie co najkvalitnejsej castice v obraze).
% Rozmer a tvar kruhu je volitelny LT (myši) kde sa  dá meniť tvar aj
% rozmery a pomocou klávesnice sa potvrdzuje vyber. Potom sa v danom vyseku v
% ktorom bola definovana castica postupne zvysuje prah. Ako hodnota prahu
% je zvolena hodnota, kde je hmotnost v zavislosti na prahu priblizne
% linearna = vyberam hodnotu minima z gradientu vahy. Druhy pristup je
% vyber minima z gradientu histogramu obrazu. V skripte je ale pre vypocet
% prahovej hodnoty pouzity prvy pristup.
global T myGLOBAL INFO;
myGLOBAL = [];

X=SelImage;
%% filtracia a skalovanie obrazu
% pouzity priemerovaci gaussovsky filter, skalovanie prebieha podelenim
% obrazu jeho najvyssou hodnotou
% J = medfilt2(X,[3 3]);
h1 = [1 2 1;2 4 2;1 2 1]./16;
% h1 = [1 2 4 2 1; 2 4 8 4 2; 4 8 16 8 4;2 4 8 4 2; 1 2 4 2 1]./16;
J = imfilter(X,h1);
b = max(J);
b=max(b);
a = min(J);
a=min(a);
J= J./b;
f=figure; imshow(J); title('Select background')
%% vytvara elipsu nahodne velku s nahodnou poziciou a vklada ju do obrazu
% na zaklade zvolenej oblasti a po kliknuti na ENTER sa vytvori maska, teda
% obraz kde vsetko okrem plochy daneho kruhu je 0 a plocha je 1. v tejto
% poloche sa potom z filtrovaneho obrazu pocita priemerna hodnota pozadia.
h = imellipse(gca, [100 100 100 100]);
fcn = makeConstrainToRectFcn('imellipse',get(gca,'XLim'),get(gca,'YLim'));
setPositionConstraintFcn(h,fcn);
pause
Bc = createMask(h);
background = mean(J(Bc==1));
%% vytvara elipsu a robi to iste co predchadzajuca cast ale pre plochu kde
% je castica
close(f)
fa=figure; imshow(J); title('Select Particle')
g = imellipse(gca, [100 100 100 100]);
fun = makeConstrainToRectFcn('imellipse',get(gca,'XLim'),get(gca,'YLim'));
setPositionConstraintFcn(g,fun);
pause
Part = createMask(g);
close(fa)
Centr = bwlabel(Part);  % tvorba znaceneho obrazu Centr
% vypocet vlastnosti masky ktora obsahuje casticu:
maska = regionprops(Centr,'PixelIdxList','Area','Centroid');
[INDa INDb]=find(Part==1); % hlada idnexy pozicie vybraneho kruhu
IND = [INDa INDb];
j=1;
a=1;
SE = strel('disk',2);
Vaha = [];
% histogram(J)
%% for cyklus, postupne zvysovanie prahu od pozadia vyssie
% for cyklus prebieha od hodnoty pozadia smerom ku vyssim hodnotam.
% Postupne zvysuje prah a prahuje obraz masky. V maske sa potom pocita
% suma jasovych urovni a odpocitava sa priemerna jasova uroven pozadia 
% (pocitana z prvej masky) 

ab = 1;
for i = background:0.001:0.5
    prah = i;
    BIN = im2bw(J,prah);
    %     BIN = imopen(BIN,SE);
    BIN = BIN.*Part;
    %     figure; imshow(BIN)
    VahaB(ab) = sum(J(BIN==1))-background;
    %     hmot2(ab) = sum(J(BIN==1))-background;
    prahy(ab) = prah;
    ab = ab+1;
end
prahy(VahaB<=1)=[];
VahaB(VahaB<=1)=[];
dif = abs(diff(VahaB));
grad = abs(gradient(VahaB));% vypocet gradientu Vahy
grad = abs(gradient(grad));
%% prelozenie krivky gradientom
% aproximácia krivky gradientu (byva casto ovplyvnena sumom a tym
% skokovita). Dochadza tym ku vyhladeniu krivky a lahsiemu vyberu minima,
% ktore reprezentuje dany prah
for i = 1:length(grad)
    MAT(i,1) = prahy(i);
    MAT(i,2) = grad(i);
end
MAT = sortrows(MAT,1);
coeff = polyfit(MAT(:,1),MAT(:,2),9);
val = polyval(coeff,MAT(:,1));
%% pracne hladanie lokalnych minim aproximovaneho gradientu
% porovnava 5 predchadzajucih a 5 nasledovnych hodnot v krivke - testovanim
% sa ukazalo ako ovela spolahlivejsie nez hladanie prikazom ,,min''
as=1;
for i = 5:length(val)-5
    
    if (val(i-4)>val(i-3))==1 && (val(i-3)>val(i-2))==1 && (val(i-2)>val(i-1))==1 && (val(i-1)>val(i))==1 && (val(i+4)>val(i+3))==1 && (val(i+3)>val(i+2))==1 && (val(i+2)>val(i+1))==1 && (val(i+1)>val(i))==1
        pom(as) = i;
        as=as+1;
    end
    
end
trAPROX = prahy(pom(1)-2); % prah (najvyhovujucejsie lokalne minimum. 
% Hodnota je posunuta o dve hodnoty spať vzhladom na klesajucu tendenciu krivky a najdeniu pociatku linearity krivky I/p)


% figure
% subplot(2,1,1); plot(prahy,VahaB); title('Vaha')
% subplot(2,1,2); plot(MAT(:,1),val,'g',MAT(:,1),MAT(:,2),'r')

% figure
% subplot(2,2,1); plot(prahyB,VahaB); title('Vaha')
% subplot(2,2,3); plot(prahyB(1:end-1), dif); title('dif')
% subplot(2,2,4); plot(prahyB(2:end-1),grad(2:end-1)); title('grad')

BINAR2 = im2bw(J,trAPROX); %binarny obraz vzniknuty prahovanim

% figure ; imshow(BINAR2) ; title('aprox-grad')

%% Uprava vzniknuteho binarneho obrazu
% imclearborder maže castice dotykajuce sa hrany, imopen maze male sumove
% castice. Pomocou bwboundaries su najdene indexy hranic objektov a su
% ulozene do strukturovanej promennej HA. Nasledne su pocita rozne
% vlastnosti castic pomcoou prikazu regionprops a vsetky tieto info su
% ulozene do globalnej promennej myGLOBAL.
FINAL = imclearborder(BINAR2); %odstranenie hranicnych castic
Se = strel('disk',3);
FINAL = imopen(FINAL,Se);
Se = strel('disk',2);
FINAL = imdilate(FINAL,Se);
% FINAL = imfill(FINAL,'holes');
HA = bwboundaries(FINAL);
[LAB,number] = bwlabel(FINAL);
properties = regionprops(LAB,'Area','Centroid','EquivDiameter', ...
    'EulerNumber','BoundingBox','Extent','Perimeter',...
    'PixelIdxList', 'FilledArea','PixelList','FilledImage',...
    'Solidity','Eccentricity','MajorAxisLength','MinorAxisLength');
for i = 1:length(properties)
    vect(i) = properties(i).Perimeter;
end

%% odtranenie pozostatkov sumu
% berie sa hodnota 0,25 priemernej velkosti obvodu ako hranica pre
% odstranenie prilis malych castic - pozostatkov impulzneho rusenia
avPer = mean(vect);
limit = avPer*0.25;
for j = 1:length(properties)
    if properties(j).Perimeter<limit
        HA(j)=[];
    end
end
% myGLOBAL.Images.LAB = LAB;
% myGLOBAL.Images.FINAL = FINAL;
% myGLOBAL.Properties = properties;
% myGLOBAL.Boundaries = HA;
% masscalc
% klasif
% ab = figure
% imshow(FINAL)

%% vykreslenie vysledku detekcie

aa=figure
imshow(J)
hold on
for ii = 1:length(HA)
    hranice1 = HA{ii};
    plot(hranice1(:,2), hranice1(:,1), 'g', 'LineWidth', 1) %Vykreslenie hranic kazdej castice
end
hold off

%% ulozenie prahu do globalnej premennej
threshold = trAPROX;
INFO.tresholdONE = threshold;
trAPROX = trAPROX*b;
INFO.Tresh = trAPROX;

