function [ newMask ] = BayesClassification( mask, img )
%UNTITLED Summary of this function goes here
%   Detailed explanation goes here
[m n] = size(mask);
CountF = sum(sum(mask));

CntNF = (m*n) - CountF;

pF = CountF / (m*n);
pNF = 1 - pF;


% countR_F = double(0);
% countG_F = double(0);
% countB_F = double(0);
% 
% countR_NF = double(0);
% countR_NF = double(0);
% countR_NF = double(0);

Rf(1,CountF) = double(0);
Gf(1,CountF) = double(0);
Bf(1,CountF) = double(0);

% R_nf = zeros;
% G_nf = zeros;
% B_nf = zeros;

counterF = double(1);
% countNF = double(1);


for i= 1:m
    for j= 1:n

        if(mask(i,j) == 1)

            Rf(1,counterF) =   img(i,j,1);
            Gf(1,counterF) =  img(i,j,2);
            Bf(1,counterF) =  img(i,j,3);
            counterF = counterF + 1;
        end
    end
end


mewR_F = sum(Rf) / (counterF-1);
mewG_F = sum(Gf) / (counterF-1);
mewB_F = sum(Bf) / (counterF-1);


% Coveriance of fence pixels

L = 1/counterF;
Xrr = Rf-mewR_F;
Xgg = Gf-mewG_F;
Xbb = Bf-mewB_F;

Srr = L*sum(Xrr.*Xrr);
Srg = L*sum(Xrr.*Xgg);
Srb = L*sum(Xrr.*Xbb);

Sgr = L*sum(Xgg.*Xrr);
Sgg = L*sum(Xgg.*Xgg);
Sgb = L*sum(Xgg.*Xbb);

Sbr = L*sum(Xbb.*Xrr);
Sbg = L*sum(Xbb.*Xgg);
Sbb = L*sum(Xbb.*Xbb);

Sigma_F =   [Srr Srg Srb;
            Sgr Sgg Sgb;
            Sbr Sbg Sbb];
        
iSigma_F = inv(Sigma_F);



newMask(m,n) = double(0);
prob(m,n) = double(0);

k1 = ( 1/( (2*22/7)^(3/2) ) * sqrt(det(Sigma_F) ));


for i=1:m
    for j=1:n
        cR = double(img(i,j,1));   %Current color
        cG = double(img(i,j,2));
        cB = double(img(i,j,3));
        A = [cR - mewR_F;
             cG - mewG_F;
             cB - mewB_F];
         
       d = -0.5 * (A'*iSigma_F*A);
       pf =  k1 * exp(d);
       prob(i,j) = pf;
       
       if(pf > 0.5)
           newMask(i,j) = 1;
       else
           newMask(i,j) = 0;
       end
       
        
    end
end

end

