function [G, Y]= delE00(gren, yelow, laf) %D6

green = rgb2lab(gren); %C1
yellow = rgb2lab(yelow); 
leaf = rgb2lab(laf); %C2

% green= meanCircle(green);
% yellow= meanCircle(yellow);
% leaf = mean_leaf(leaf);

deltaLbarg = leaf(1) - green(1); 
lbar = (green(1)+leaf(1))/2;

c1barg = sqrt(green(2)^2 + green(3)^2);
c2barg = sqrt(leaf(2)^2 + leaf(3)^2);


cinversebarg= (c1barg +c2barg)/2;
delCbar = c2barg - c1barg;


a1barg = green(2) + (green(2)/2)* (1- sqrt(cinversebarg^(7)/(cinversebarg^(7)+25^(7))));
a2barg = leaf(2) + (leaf(2)/2)* (1- sqrt(cinversebarg^(7)/(cinversebarg^(7)+25^(7))));

h1= atan2d(green(3),a1barg);
h1bar = mod(h1, 360);
%h1bar = mod(atan2d(green(3), a1barg) , 360);
h2= atan2d(leaf(3),a2barg);
h2bar = mod(h2, 360);
%h2bar = mod(atan2d(leaf(3), a2barg) , 360);
diff = abs(h1bar- h2bar);
if(diff <= 180)
    deltahbarg = h2bar-h1bar;
elseif (diff >180 && (h2bar <= h1bar))
    deltahbarg = (h2bar-h1bar) + 360;
   
elseif (diff >180 && (h2bar > h1bar))
     deltahbarg = (h2bar-h1bar) - 360;
     
end

delHbarg = 2* sqrt(c1barg*c2barg)*sind(deltahbarg/2);

if (diff <= 180)
    Hinversebarg = (h1bar+h2bar)/2;
elseif (diff >180 && (h1bar+h2bar)<360)
    Hinversebarg = (h1bar+h2bar +360)/2;
elseif (diff >180 && (h1bar+h2bar)>=360)
    Hinversebarg = (h1bar+h2bar -360)/2;
end

T = 1- 0.17*cosd(Hinversebarg - 30) + 0.24 *cosd(2*Hinversebarg) + 0.32*cosd(3*Hinversebarg +6) - 0.20* cosd(4*Hinversebarg - 63);

sl = 1+ (0.015* (lbar-50)^2/ sqrt(20+ (lbar - 50)^2));
sc = 1+ 0.045*cinversebarg;
sh = 1+ 0.015*cinversebarg*T;
rt = -2 * sqrt(cinversebarg^(7)/(cinversebarg^(7)+ 25^7)) * sind(60 * exp(-1* ((Hinversebarg - 275)/25)^2));


d1 = sqrt((deltaLbarg/sl)^2 + (delCbar/sc)^2 + (delHbarg / sh)^2 + rt*(delCbar/sc)*(delHbarg / sh)); %d1


deltaLbarg = leaf(1) - yellow(1);
lbar = (yellow(1)+leaf(1))/2;

c1bary = sqrt(yellow(2)^2 + yellow(3)^2);
c2bary = sqrt(leaf(2)^2 + leaf(3)^2);


cinversebary= (c1bary +c2bary)/2;
delCbar = c2bary - c1bary;


a1bary = yellow(2) + (yellow(2)/2)* (1- sqrt(cinversebary^(7)/(cinversebary^(7)+25^(7))));
a2bary = leaf(2) + (leaf(2)/2)* (1- sqrt(cinversebary^(7)/(cinversebary^(7)+25^(7))));


h1= atan2d(yellow(3),a1bary);
h1bar = mod(h1, 360);
%h1bar = mod(atan2d(yellow(3), a1barg) , 360);
h2= atan2d(leaf(3),a2bary);
h2bar = mod(h2, 360);
%h2bar = mod(atan2d(leaf(3), a2barg) , 360);
diff = abs(h1bar- h2bar);
if(diff <= 180)
    deltahbary = h2bar-h1bar;
elseif (diff >180 && (h2bar <= h1bar))
    deltahbary = (h2bar-h1bar) + 360;
   
elseif (diff >180 && (h2bar > h1bar))
     deltahbary = (h2bar-h1bar) - 360;
     
end

delHbary = 2* sqrt(c1bary*c2bary)*sind(deltahbary/2);

if (diff <= 180)
    Hinversebary = (h1bar+h2bar)/2;
elseif (diff >180 && (h1bar+h2bar)<360)
    Hinversebary = (h1bar+h2bar +360)/2;
elseif (diff >180 && (h1bar+h2bar)>=360)
    Hinversebary = (h1bar+h2bar -360)/2;
end

T = 1- 0.17*cosd(Hinversebary - 30) + 0.24 *cosd(2*Hinversebary) + 0.32*cosd(3*Hinversebary +6) - 0.20* cosd(4*Hinversebary - 63);

sl = 1+ (0.015* (lbar-50)^2/ sqrt(20+ (lbar - 50)^2));
sc = 1+ 0.045*cinversebary;
sh = 1+ 0.015*cinversebary*T;
rt = -2 * sqrt(cinversebary^(7)/(cinversebary^(7)+ 25^7)) * sind(60 * exp(-1* ((Hinversebary - 275)/25)^2));


d2 = sqrt((deltaLbarg/sl)^2 + (delCbar/sc)^2 + (delHbary / sh)^2 + rt*(delCbar/sc)*(delHbary / sh)); %d2

G = 100 *(1- (d1/ (d1+d2)));       %Green percentage
Y = 100 *(1- (d2/(d1+d2)));        %Yellow percentage
end
