function integral = project(x1,y1,x2,y2,shift,grow_hrt,grow_lng)
%
%  integral = project(x1,y1,x2,y2,shift,grow_hrt,grow_lng)
%
%  Project calculates the integral on a line from point
%  (X1,Y1) to (X2,Y2) of a synthetic 2-D slice model
%  of a human chest.  SHIFT allows a horizontal translation
%  while GROW_HRT and GROW_LNG allows a dilation of the "heart"
%  and "lungs & chest wall", respectively.  The model is contained
%  within a square box with vertices at (0,0) and (1,1).  Any points
%  can be given for the line, including inside the model.
%  SHIFT is limited to +/-0.04, GROW_HRT is limited
%  to 0 to 0.04, and GROW_LNG is limited from 0 to 0.08.
%
%  USE  modgraph(shift,grow_hrt,grow_lng,res) to
%  display the model.  RES = 0.004 is suggested.
%
%  SEE modgraph.m, projgraf.m
%
%  Developed by Timothy D. Dorney
%               Rice University
%               February, 1999
%               tdorney@ieee.org
%
%  Coded using MATLAB 5.X.X with NO additional toolboxes.
%
%	REVISION HISTORY
%
%	VERSION 1.0.0		FEB. 15, 1999	TIM DORNEY
%	VERSION 1.0.1		FEB. 16, 1999	TIM DORNEY
%	  Added feature so that if x1==x2 & y1==y2
%         the value of the density at the point is returned.
%	  This allows a method to graph the 2-D slice without
%         using MODGRAPH.
%	VERSION 1.0.2		FEB. 28, 1999	TIM DORNEY
%	  MATLAB rounding errors are causing problems with
%	  hroizontal and vertical projections.  Projection 
%         end points will be limited to 7 significant digits
%         after the decimal point.
%

if (nargin ~= 7)
        disp('PROJECT requires 7 input arguments!')
	return;
end
 
if ((shift < -0.04) | (shift > 0.04))
	shift = 0;
	disp('SHIFT limited to +/- 0.04.  The limit was exceeded.  SHIFT set to 0!');
end
if ((grow_hrt < 0) | (grow_hrt > 0.04))
	grow_hrt = 0;
	disp('GROW_HRT limited from 0 to 0.04.  The limit was exceeded.  GROW_HRT set to 0!');
end
if ((grow_lng < 0) | (grow_lng > 0.08))
	grow_lng = 0;
	disp('GROW_LNG limited from 0 to 0.08.  The limit was exceeded.  GROW_LNG set to 0!');
end

% DENSITY		REGION		PART			GROWTH
DA = 1.00;		% A		center spine circle     N
DB = 0.99;		% B             left spine ellipse	N
DC = 0.98;		% C		right spine ellipse	N
DD = 0.5;		% D		center heart circle	Symmetric
DE = 0.6;		% E		left lung ellipse	Asymmetric left
DF = 0.61;		% F		right lung ellipse	Asymmetric right
DG = 0.4;		% G		outer skin ellipse	Asymmetric vertical

% LIMIT PROJECTION END POINT ACCURACY TO 7 SIG. FIGS. AFTER DECIMAL POINT

x1 = (round(x1 * 10^7))./(10^7);
x2 = (round(x2 * 10^7))./(10^7);
y1 = (round(y1 * 10^7))./(10^7);
y2 = (round(y2 * 10^7))./(10^7);

% ========== GGGGGGGGGG ==========

if (((x1-x2) ~= 0) & ((y1-y2) ~= 0))
	M = (y1-y2)/(x1-x2);
	YC = y1-M*x1;
	xx = x1;	yy = y1;
% REGION G
	A = 2.25;	B = 1.95;	C = 0.03;	D = (sqrt(C+grow_lng/10) - sqrt(C))*B;
	[e,f] = find( ( ((xx-0.5-shift)./A).^2 + ((yy-0.35-D)./B).^2 ) <= C+grow_lng/10);

	if (~isempty(f))	% IF NOT NULL THEN (X1,Y1) ARE INSIDE THE BODY.
		GX1 = x1;
		GY1 = y1;
	else
		AA = 1/(A^2) + (M^2)/(B^2);
		BB = (2*(-0.5-shift))/(A^2) + 2*M*(YC-0.35-D)/(B^2);
		CC = ((-0.5-shift)^2)/(A^2) + ((YC-0.35-D)^2)/(B^2) - C - grow_lng/10;
		XANS1 = (-BB + sqrt(BB^2 - 4*AA*CC))/(2*AA);
		XANS2 = (-BB - sqrt(BB^2 - 4*AA*CC))/(2*AA);
		if ((imag(XANS1) ~= 0) & (imag(XANS2) ~= 0))	% NEVER INTERSECT.
			integral = 0;
			return;
		else
			YANS1 = M*XANS1 + YC;
			YANS2 = M*XANS2 + YC;
			DIST1 = sqrt((XANS1-x1)^2 + (YANS1-y1)^2);
			DIST2 = sqrt((XANS2-x1)^2 + (YANS2-y1)^2);
			[Y,I] = min([DIST1 DIST2]);

			if (Y == 1)
				GX1 = XANS1; GY1 = YANS1;
			else
				GX1 = XANS2; GY1 = YANS2;
			end
		end
	end
	xx = x2;	yy = y2;
% REGION G
	A = 2.25;	B = 1.95;	C = 0.03;	D = (sqrt(C+grow_lng/10) - sqrt(C))*B;
	[e,f] = find( ( ((xx-0.5-shift)./A).^2 + ((yy-0.35-D)./B).^2 ) <= C+grow_lng/10);

% IF NOT NULL THEN (X2,Y2) ARE INSIDE THE BODY.
	if (~isempty(f))
		GX2 = x2;
		GY2 = y2;
	else
		AA = 1/(A^2) + (M^2)/(B^2);
		BB = (2*(-0.5-shift))/(A^2) + 2*M*(YC-0.35-D)/(B^2);
		CC = ((-0.5-shift)^2)/(A^2) + ((YC-0.35-D)^2)/(B^2) - C - grow_lng/10;
		XANS1 = (-BB + sqrt(BB^2 - 4*AA*CC))/(2*AA);
		XANS2 = (-BB - sqrt(BB^2 - 4*AA*CC))/(2*AA);
		if ((imag(XANS1) ~= 0) & (imag(XANS2) ~= 0))	% NEVER INTERSECT.
			integral = 0;
			return;
		else
			YANS1 = M*XANS1 + YC;
			YANS2 = M*XANS2 + YC;
			DIST1 = sqrt((XANS1-x1)^2 + (YANS1-y1)^2);
			DIST2 = sqrt((XANS2-x1)^2 + (YANS2-y1)^2);
			[Y,I] = min([DIST1 DIST2]);

			if (Y == 1)
				GX2 = XANS2; GY2 = YANS2;
			else
				GX2 = XANS1; GY2 = YANS1;
			end
		end
	end
	GDIST = sqrt((GX1-GX2)^2 + (GY1-GY2)^2);
elseif (((x1-x2) == 0) & ((y1-y2) ~= 0))
	xx = x1;	yy = y1;
% REGION G
	A = 2.25;	B = 1.95;	C = 0.03;	D = (sqrt(C+grow_lng/10) - sqrt(C))*B;
	[e,f] = find( ( ((xx-0.5-shift)./A).^2 + ((yy-0.35-D)./B).^2 ) <= C+grow_lng/10);

	if (~isempty(f))	% IF NOT NULL THEN (X1,Y1) ARE INSIDE THE BODY.
		GX1 = x1;
		GY1 = y1;
	else
		T1 = sqrt((C+grow_lng/10) - ((xx-0.5-shift)./A).^2);
		YANS1 = (B*T1)+ 0.35 + D;
		YANS2 = (B*-T1)+ 0.35 + D;
		if ((imag(YANS1) ~= 0) & (imag(YANS2) ~= 0))	% NEVER INTERSECT.
			integral = 0;
			return;
		else
			XANS1 = x1;
			XANS2 = x2;
			DIST1 = sqrt((XANS1-x1)^2 + (YANS1-y1)^2);
			DIST2 = sqrt((XANS2-x1)^2 + (YANS2-y1)^2);
			[Y,I] = min([DIST1 DIST2]);

			if (Y == 1)
				GX1 = XANS1; GY1 = YANS1;
			else
				GX1 = XANS2; GY1 = YANS2;
			end
		end
	end
	xx = x2;	yy = y2;
% REGION G
	A = 2.25;	B = 1.95;	C = 0.03;	D = (sqrt(C+grow_lng/10) - sqrt(C))*B;
	[e,f] = find( ( ((xx-0.5-shift)./A).^2 + ((yy-0.35-D)./B).^2 ) <= C+grow_lng/10);

% IF NOT NULL THEN (X2,Y2) ARE INSIDE THE BODY.
	if (~isempty(f))
		GX2 = x2;
		GY2 = y2;
	else
		T1 = sqrt((C+grow_lng/10) - ((xx-0.5-shift)./A).^2);
		YANS1 = (B*T1)+ 0.35 + D;
		YANS2 = (B*-T1)+ 0.35 + D;
		if ((imag(YANS1) ~= 0) & (imag(YANS2) ~= 0))	% NEVER INTERSECT.
			integral = 0;
			return;
		else
			XANS1 = x1;
			XANS2 = x2;
			DIST1 = sqrt((XANS1-x1)^2 + (YANS1-y1)^2);
			DIST2 = sqrt((XANS2-x1)^2 + (YANS2-y1)^2);
			[Y,I] = min([DIST1 DIST2]);

			if (Y == 1)
				GX2 = XANS2; GY2 = YANS2;
			else
				GX2 = XANS1; GY2 = YANS1;
			end
		end
	end
	GDIST = sqrt((GX1-GX2)^2 + (GY1-GY2)^2);
elseif (((x1-x2) ~= 0) & ((y1-y2) == 0))
	xx = x1;	yy = y1;
% REGION G
	A = 2.25;	B = 1.95;	C = 0.03;	D = (sqrt(C+grow_lng/10) - sqrt(C))*B;
	[e,f] = find( ( ((xx-0.5-shift)./A).^2 + ((yy-0.35-D)./B).^2 ) <= C+grow_lng/10);

	if (~isempty(f))	% IF NOT NULL THEN (X1,Y1) ARE INSIDE THE BODY.
		GX1 = x1;
		GY1 = y1;
	else
		T1 = sqrt((C+grow_lng/10) - ((yy-0.35-D)./B).^2);
		XANS1 = (A*T1)+ 0.5 + shift;
		XANS2 = (A*-T1)+ 0.5 + shift;
		if ((imag(XANS1) ~= 0) & (imag(XANS2) ~= 0))	% NEVER INTERSECT.
			integral = 0;
			return;
		else
			YANS1 = y1;
			YANS2 = y2;
			DIST1 = sqrt((XANS1-x1)^2 + (YANS1-y1)^2);
			DIST2 = sqrt((XANS2-x1)^2 + (YANS2-y1)^2);
			[Y,I] = min([DIST1 DIST2]);

			if (Y == 1)
				GX1 = XANS1; GY1 = YANS1;
			else
				GX1 = XANS2; GY1 = YANS2;
			end
		end
	end
	xx = x2;	yy = y2;
% REGION G
	A = 2.25;	B = 1.95;	C = 0.03;	D = (sqrt(C+grow_lng/10) - sqrt(C))*B;
	[e,f] = find( ( ((xx-0.5-shift)./A).^2 + ((yy-0.35-D)./B).^2 ) <= C+grow_lng/10);

% IF NOT NULL THEN (X2,Y2) ARE INSIDE THE BODY.
	if (~isempty(f))
		GX2 = x2;
		GY2 = y2;
	else
		T1 = sqrt((C+grow_lng/10) - ((yy-0.35-D)./B).^2);
		XANS1 = (A*T1)+ 0.5 + shift;
		XANS2 = (A*-T1)+ 0.5 + shift;
		if ((imag(XANS1) ~= 0) & (imag(XANS2) ~= 0))	% NEVER INTERSECT.
			integral = 0;
			return;
		else
			YANS1 = y1;
			YANS2 = y2;
			DIST1 = sqrt((XANS1-x1)^2 + (YANS1-y1)^2);
			DIST2 = sqrt((XANS2-x1)^2 + (YANS2-y1)^2);
			[Y,I] = min([DIST1 DIST2]);

			if (Y == 1)
				GX2 = XANS2; GY2 = YANS2;
			else
				GX2 = XANS1; GY2 = YANS1;
			end
		end
	end
	GDIST = sqrt((GX1-GX2)^2 + (GY1-GY2)^2);
else	% x1 == x2 & y1 == y2
	xx = x1;	yy = y1;
% REGION G
	A = 2.25;	B = 1.95;	C = 0.03;	D = (sqrt(C+grow_lng/10) - sqrt(C))*B;
	[e,f] = find( ( ((xx-0.5-shift)./A).^2 + ((yy-0.35-D)./B).^2 ) <= C+grow_lng/10);

% IF NOT NULL THEN (X2,Y2) ARE INSIDE THE BODY.
	if (~isempty(f))
		GDIST = 999;
	else
		GDIST = 0;
	end
end
% GDIST

% ========== AAAAAAAAAA ==========

AFLAG = 0;
if (((x1-x2) ~= 0) & ((y1-y2) ~= 0))
	M = (y1-y2)/(x1-x2);
	YC = y1-M*x1;
	xx = x1;	yy = y1;
% REGION A
	[e,f] = find (((xx-0.5-shift).^2 + (yy-0.1).^2) <= 0.0009);

	if (~isempty(f))	% IF NOT NULL THEN (X1,Y1) ARE INSIDE REGION A.
		GX1 = x1;
		GY1 = y1;
	else
		AA = 1 + (M^2);
		BB = (2*(-0.5-shift)) + 2*M*(YC-0.1);
		CC = ((-0.5-shift)^2) + ((YC-0.1)^2) - 0.0009;
		XANS1 = (-BB + sqrt(BB^2 - 4*AA*CC))/(2*AA);
		XANS2 = (-BB - sqrt(BB^2 - 4*AA*CC))/(2*AA);
		if ((imag(XANS1) ~= 0) & (imag(XANS2) ~= 0))	% NEVER INTERSECT.
			ADIST = 0;
			AFLAG = 1;
		else
			YANS1 = M*XANS1 + YC;
			YANS2 = M*XANS2 + YC;
			DIST1 = sqrt((XANS1-x1)^2 + (YANS1-y1)^2);
			DIST2 = sqrt((XANS2-x1)^2 + (YANS2-y1)^2);
			[Y,I] = min([DIST1 DIST2]);

			if (Y == 1)
				GX1 = XANS1; GY1 = YANS1;
			else
				GX1 = XANS2; GY1 = YANS2;
			end
		end
	end
	if (AFLAG == 0)
		xx = x2;	yy = y2;
% REGION A
		[e,f] = find (((xx-0.5-shift).^2 + (yy-0.1).^2) <= 0.0009);

% IF NOT NULL THEN (X2,Y2) ARE INSIDE REGION A.
		if (~isempty(f))
			GX2 = x2;
			GY2 = y2;
		else
			AA = 1 + (M^2);
			BB = (2*(-0.5-shift)) + 2*M*(YC-0.1);
			CC = ((-0.5-shift)^2) + ((YC-0.1)^2) - 0.0009;
			XANS1 = (-BB + sqrt(BB^2 - 4*AA*CC))/(2*AA);
			XANS2 = (-BB - sqrt(BB^2 - 4*AA*CC))/(2*AA);
			YANS1 = M*XANS1 + YC;
			YANS2 = M*XANS2 + YC;
			DIST1 = sqrt((XANS1-x1)^2 + (YANS1-y1)^2);
			DIST2 = sqrt((XANS2-x1)^2 + (YANS2-y1)^2);
			[Y,I] = min([DIST1 DIST2]);

			if (Y == 1)
				GX2 = XANS2; GY2 = YANS2;
			else
				GX2 = XANS1; GY2 = YANS1;
			end
		end
		ADIST = sqrt((GX1-GX2)^2 + (GY1-GY2)^2);
	end
elseif (((x1-x2) == 0) & ((y1-y2) ~= 0))
	xx = x1;	yy = y1;
% REGION A
	[e,f] = find (((xx-0.5-shift).^2 + (yy-0.1).^2) <= 0.0009);

	if (~isempty(f))	% IF NOT NULL THEN (X1,Y1) ARE INSIDE REGION A.
		GX1 = x1;
		GY1 = y1;
	else
		T1 = sqrt(0.0009 - ((xx-0.5-shift).^2));
		YANS1 = T1 + 0.1;
		YANS2 = -T1 + 0.1;
		if ((imag(YANS1) ~= 0) & (imag(YANS2) ~= 0))	% NEVER INTERSECT.
			ADIST = 0;
			AFLAG = 1;
		else
			XANS1 = x1;
			XANS2 = x2;
			DIST1 = sqrt((XANS1-x1)^2 + (YANS1-y1)^2);
			DIST2 = sqrt((XANS2-x1)^2 + (YANS2-y1)^2);
			[Y,I] = min([DIST1 DIST2]);

			if (Y == 1)
				GX1 = XANS1; GY1 = YANS1;
			else
				GX1 = XANS2; GY1 = YANS2;
			end
		end
	end
	if (AFLAG == 0)
		xx = x2;	yy = y2;
% REGION A
		[e,f] = find (((xx-0.5-shift).^2 + (yy-0.1).^2) <= 0.0009);

% IF NOT NULL THEN (X2,Y2) ARE INSIDE REGION A.
		if (~isempty(f))
			GX2 = x2;
			GY2 = y2;
		else
			T1 = sqrt(0.0009 - ((xx-0.5-shift).^2));
			YANS1 = T1 + 0.1;
			YANS2 = -T1 + 0.1;
			XANS1 = x1;
			XANS2 = x2;
			DIST1 = sqrt((XANS1-x1)^2 + (YANS1-y1)^2);
			DIST2 = sqrt((XANS2-x1)^2 + (YANS2-y1)^2);
			[Y,I] = min([DIST1 DIST2]);

			if (Y == 1)
				GX2 = XANS2; GY2 = YANS2;
			else
				GX2 = XANS1; GY2 = YANS1;
			end
		end
		ADIST = sqrt((GX1-GX2)^2 + (GY1-GY2)^2);
	end
elseif (((x1-x2) ~= 0) & ((y1-y2) == 0))
	xx = x1;	yy = y1;
% REGION A
	[e,f] = find (((xx-0.5-shift).^2 + (yy-0.1).^2) <= 0.0009);

	if (~isempty(f))	% IF NOT NULL THEN (X1,Y1) ARE INSIDE REGION A.
		GX1 = x1;
		GY1 = y1;
	else
		T1 = sqrt(0.0009 - ((yy-0.1).^2));
		XANS1 = T1 + 0.5 + shift;
		XANS2 = -T1 + 0.5 + shift;
		if ((imag(XANS1) ~= 0) & (imag(XANS2) ~= 0))	% NEVER INTERSECT.
			ADIST = 0;
			AFLAG = 1;
		else
			YANS1 = y1;
			YANS2 = y2;
			DIST1 = sqrt((XANS1-x1)^2 + (YANS1-y1)^2);
			DIST2 = sqrt((XANS2-x1)^2 + (YANS2-y1)^2);
			[Y,I] = min([DIST1 DIST2]);

			if (Y == 1)
				GX1 = XANS1; GY1 = YANS1;
			else
				GX1 = XANS2; GY1 = YANS2;
			end
		end
	end
	if (AFLAG == 0)
		xx = x2;	yy = y2;
% REGION A
		[e,f] = find (((xx-0.5-shift).^2 + (yy-0.1).^2) <= 0.0009);

% IF NOT NULL THEN (X2,Y2) ARE INSIDE REGION A.
		if (~isempty(f))
			GX2 = x2;
			GY2 = y2;
		else
			T1 = sqrt(0.0009 - ((yy-0.1).^2));
			XANS1 = T1 + 0.5 + shift;
			XANS2 = -T1 + 0.5 + shift;
			YANS1 = y1;
			YANS2 = y2;
			DIST1 = sqrt((XANS1-x1)^2 + (YANS1-y1)^2);
			DIST2 = sqrt((XANS2-x1)^2 + (YANS2-y1)^2);
			[Y,I] = min([DIST1 DIST2]);

			if (Y == 1)
				GX2 = XANS2; GY2 = YANS2;
			else
				GX2 = XANS1; GY2 = YANS1;
			end
		end
		ADIST = sqrt((GX1-GX2)^2 + (GY1-GY2)^2);
	end
else	% x1 == x2 & y1 == y2
	xx = x1;	yy = y1;
% REGION A
	[e,f] = find (((xx-0.5-shift).^2 + (yy-0.1).^2) <= 0.0009);

% IF NOT NULL THEN (X2,Y2) ARE INSIDE REGION A.
	if (~isempty(f))
		ADIST = 999;
	else
		ADIST = 0;
	end
end
% ADIST

% ========== BBBBBBBBBB ==========

BFLAG = 0;
if (((x1-x2) ~= 0) & ((y1-y2) ~= 0))
	M = (y1-y2)/(x1-x2);
	YC = y1-M*x1;
	xx = x1;	yy = y1;
% REGION B
	A = 2;	B = 1;
	[e,f] = find (( ((xx-0.415-shift)./A).^2 + ((yy-0.1)./B).^2) <= 0.000625);

	if (~isempty(f))	% IF NOT NULL THEN (X1,Y1) ARE INSIDE REGION B.
		GX1 = x1;
		GY1 = y1;
	else
		AA = 1/(A^2) + (M^2)/(B^2);
		BB = (2*(-0.415-shift))/(A^2) + 2*M*(YC-0.1)/(B^2);
		CC = ((-0.415-shift)^2)/(A^2) + ((YC-0.1)^2)/(B^2) - 0.000625;
		XANS1 = (-BB + sqrt(BB^2 - 4*AA*CC))/(2*AA);
		XANS2 = (-BB - sqrt(BB^2 - 4*AA*CC))/(2*AA);
		if ((imag(XANS1) ~= 0) & (imag(XANS2) ~= 0))	% NEVER INTERSECT.
			BDIST = 0;
			BFLAG = 1;
		else
			YANS1 = M*XANS1 + YC;
			YANS2 = M*XANS2 + YC;
			DIST1 = sqrt((XANS1-x1)^2 + (YANS1-y1)^2);
			DIST2 = sqrt((XANS2-x1)^2 + (YANS2-y1)^2);
			[Y,I] = min([DIST1 DIST2]);

			if (Y == 1)
				GX1 = XANS1; GY1 = YANS1;
			else
				GX1 = XANS2; GY1 = YANS2;
			end
		end
	end
	if (BFLAG == 0)
		xx = x2;	yy = y2;
% REGION B
		A = 2;	B = 1;
		[e,f] = find (( ((xx-0.415-shift)./A).^2 + ((yy-0.1)./B).^2) <= 0.000625);

% IF NOT NULL THEN (X2,Y2) ARE INSIDE REGION B.
		if (~isempty(f))
			GX2 = x2;
			GY2 = y2;
		else
			AA = 1/(A^2) + (M^2)/(B^2);
			BB = (2*(-0.415-shift))/(A^2) + 2*M*(YC-0.1)/(B^2);
			CC = ((-0.415-shift)^2)/(A^2) + ((YC-0.1)^2)/(B^2) - 0.000625;
			XANS1 = (-BB + sqrt(BB^2 - 4*AA*CC))/(2*AA);
			XANS2 = (-BB - sqrt(BB^2 - 4*AA*CC))/(2*AA);
			YANS1 = M*XANS1 + YC;
			YANS2 = M*XANS2 + YC;
			DIST1 = sqrt((XANS1-x1)^2 + (YANS1-y1)^2);
			DIST2 = sqrt((XANS2-x1)^2 + (YANS2-y1)^2);
			[Y,I] = min([DIST1 DIST2]);

			if (Y == 1)
				GX2 = XANS2; GY2 = YANS2;
			else
				GX2 = XANS1; GY2 = YANS1;
			end
		end
		BDIST = sqrt((GX1-GX2)^2 + (GY1-GY2)^2);
	end
elseif (((x1-x2) == 0) & ((y1-y2) ~= 0))
	xx = x1;	yy = y1;
% REGION B
	A = 2;	B = 1;
	[e,f] = find (( ((xx-0.415-shift)./A).^2 + ((yy-0.1)./B).^2) <= 0.000625);

	if (~isempty(f))	% IF NOT NULL THEN (X1,Y1) ARE INSIDE REGION B.
		GX1 = x1;
		GY1 = y1;
	else
		T1 = sqrt(0.000625 - ((xx-0.415-shift)./A).^2);
		YANS1 = (B*T1) + 0.1;
		YANS2 = (B*-T1) + 0.1;
		if ((imag(YANS1) ~= 0) & (imag(YANS2) ~= 0))	% NEVER INTERSECT.
			BDIST = 0;
			BFLAG = 1;
		else
			XANS1 = x1;
			XANS2 = x2;
			DIST1 = sqrt((XANS1-x1)^2 + (YANS1-y1)^2);
			DIST2 = sqrt((XANS2-x1)^2 + (YANS2-y1)^2);
			[Y,I] = min([DIST1 DIST2]);

			if (Y == 1)
				GX1 = XANS1; GY1 = YANS1;
			else
				GX1 = XANS2; GY1 = YANS2;
			end
		end
	end
	if (BFLAG == 0)
		xx = x2;	yy = y2;
% REGION B
		A = 2;	B = 1;
		[e,f] = find (( ((xx-0.415-shift)./A).^2 + ((yy-0.1)./B).^2) <= 0.000625);

% IF NOT NULL THEN (X2,Y2) ARE INSIDE REGION B.
		if (~isempty(f))
			GX2 = x2;
			GY2 = y2;
		else
			T1 = sqrt(0.000625 - ((xx-0.415-shift)./A).^2);
			YANS1 = (B*T1) + 0.1;
			YANS2 = (B*-T1) + 0.1;
			XANS1 = x1;
			XANS2 = x2;
			DIST1 = sqrt((XANS1-x1)^2 + (YANS1-y1)^2);
			DIST2 = sqrt((XANS2-x1)^2 + (YANS2-y1)^2);
			[Y,I] = min([DIST1 DIST2]);

			if (Y == 1)
				GX2 = XANS2; GY2 = YANS2;
			else
				GX2 = XANS1; GY2 = YANS1;
			end
		end
	BDIST = sqrt((GX1-GX2)^2 + (GY1-GY2)^2);
	end
elseif (((x1-x2) ~= 0) & ((y1-y2) == 0))
	xx = x1;	yy = y1;
% REGION B
	A = 2;	B = 1;
	[e,f] = find (( ((xx-0.415-shift)./A).^2 + ((yy-0.1)./B).^2) <= 0.000625);

	if (~isempty(f))	% IF NOT NULL THEN (X1,Y1) ARE INSIDE REGION B.
		GX1 = x1;
		GY1 = y1;
	else
		T1 = sqrt(0.000625 - ((yy-0.1)./B).^2);
		XANS1 = (A*T1) + 0.415 + shift;
		XANS2 = (A*-T1) + 0.415 + shift;
		if ((imag(XANS1) ~= 0) & (imag(XANS2) ~= 0))	% NEVER INTERSECT.
			BDIST = 0;
			BFLAG = 1;
		else
			YANS1 = y1;
			YANS2 = y2;
			DIST1 = sqrt((XANS1-x1)^2 + (YANS1-y1)^2);
			DIST2 = sqrt((XANS2-x1)^2 + (YANS2-y1)^2);
			[Y,I] = min([DIST1 DIST2]);

			if (Y == 1)
				GX1 = XANS1; GY1 = YANS1;
			else
				GX1 = XANS2; GY1 = YANS2;
			end
		end
	end
	if (BFLAG == 0)
		xx = x2;	yy = y2;
% REGION B
		A = 2;	B = 1;
		[e,f] = find (( ((xx-0.415-shift)./A).^2 + ((yy-0.1)./B).^2) <= 0.000625);

% IF NOT NULL THEN (X2,Y2) ARE INSIDE REGION B.
		if (~isempty(f))
			GX2 = x2;
			GY2 = y2;
		else
			T1 = sqrt(0.000625 - ((yy-0.1)./B).^2);
			XANS1 = (A*T1) + 0.415 + shift;
			XANS2 = (A*-T1) + 0.415 + shift;
			YANS1 = y1;
			YANS2 = y2;
			DIST1 = sqrt((XANS1-x1)^2 + (YANS1-y1)^2);
			DIST2 = sqrt((XANS2-x1)^2 + (YANS2-y1)^2);
			[Y,I] = min([DIST1 DIST2]);

			if (Y == 1)
				GX2 = XANS2; GY2 = YANS2;
			else
				GX2 = XANS1; GY2 = YANS1;
			end
		end
		BDIST = sqrt((GX1-GX2)^2 + (GY1-GY2)^2);
	end
else	% x1 == x2 & y1 == y2
	xx = x1;	yy = y1;
% REGION B
	A = 2;	B = 1;
	[e,f] = find (( ((xx-0.415-shift)./A).^2 + ((yy-0.1)./B).^2) <= 0.000625);

% IF NOT NULL THEN (X2,Y2) ARE INSIDE REGION B.
	if (~isempty(f))
		BDIST = 999;
	else
		BDIST = 0;
	end
end
% BDIST

% ========== CCCCCCCCCC ==========

CFLAG = 0;
if (((x1-x2) ~= 0) & ((y1-y2) ~= 0))
	M = (y1-y2)/(x1-x2);
	YC = y1-M*x1;
	xx = x1;	yy = y1;
% REGION C
	A = 1;	B = 2;
	[e,f] = find (( ((xx-0.56-shift)./A).^2 + ((yy-0.1)./B).^2) <= 0.000625);

	if (~isempty(f))	% IF NOT NULL THEN (X1,Y1) ARE INSIDE REGION C.
		GX1 = x1;
		GY1 = y1;
	else
		AA = 1/(A^2) + (M^2)/(B^2);
		BB = (2*(-0.56-shift))/(A^2) + 2*M*(YC-0.1)/(B^2);
		CC = ((-0.56-shift)^2)/(A^2) + ((YC-0.1)^2)/(B^2) - 0.000625;
		XANS1 = (-BB + sqrt(BB^2 - 4*AA*CC))/(2*AA);
		XANS2 = (-BB - sqrt(BB^2 - 4*AA*CC))/(2*AA);
		if ((imag(XANS1) ~= 0) & (imag(XANS2) ~= 0))	% NEVER INTERSECT.
			CDIST = 0;
			CFLAG = 1;
		else
			YANS1 = M*XANS1 + YC;
			YANS2 = M*XANS2 + YC;
			DIST1 = sqrt((XANS1-x1)^2 + (YANS1-y1)^2);
			DIST2 = sqrt((XANS2-x1)^2 + (YANS2-y1)^2);
			[Y,I] = min([DIST1 DIST2]);

			if (Y == 1)
				GX1 = XANS1; GY1 = YANS1;
			else
				GX1 = XANS2; GY1 = YANS2;
			end
		end
	end
	if (CFLAG == 0)
		xx = x2;	yy = y2;
% REGION C
		A = 1;	B = 2;
		[e,f] = find (( ((xx-0.56-shift)./A).^2 + ((yy-0.1)./B).^2) <= 0.000625);

% IF NOT NULL THEN (X2,Y2) ARE INSIDE REGION C.
		if (~isempty(f))
			GX2 = x2;
			GY2 = y2;
		else
			AA = 1/(A^2) + (M^2)/(B^2);
			BB = (2*(-0.56-shift))/(A^2) + 2*M*(YC-0.1)/(B^2);
			CC = ((-0.56-shift)^2)/(A^2) + ((YC-0.1)^2)/(B^2) - 0.000625;
			XANS1 = (-BB + sqrt(BB^2 - 4*AA*CC))/(2*AA);
			XANS2 = (-BB - sqrt(BB^2 - 4*AA*CC))/(2*AA);
			YANS1 = M*XANS1 + YC;
			YANS2 = M*XANS2 + YC;
			DIST1 = sqrt((XANS1-x1)^2 + (YANS1-y1)^2);
			DIST2 = sqrt((XANS2-x1)^2 + (YANS2-y1)^2);
			[Y,I] = min([DIST1 DIST2]);

			if (Y == 1)
				GX2 = XANS2; GY2 = YANS2;
			else
				GX2 = XANS1; GY2 = YANS1;
			end
		end
		CDIST = sqrt((GX1-GX2)^2 + (GY1-GY2)^2);
	end
elseif (((x1-x2) == 0) & ((y1-y2) ~= 0))
	xx = x1;	yy = y1;
% REGION C
	A = 1;	B = 2;
	[e,f] = find (( ((xx-0.56-shift)./A).^2 + ((yy-0.1)./B).^2) <= 0.000625);

	if (~isempty(f))	% IF NOT NULL THEN (X1,Y1) ARE INSIDE REGION C.
		GX1 = x1;
		GY1 = y1;
	else
		T1 = sqrt(0.000625 - ((xx-0.56-shift)./A).^2);
		YANS1 = (B*T1) + 0.1;
		YANS2 = (B*-T1) + 0.1;
		if ((imag(YANS1) ~= 0) & (imag(YANS2) ~= 0))	% NEVER INTERSECT.
			CDIST = 0;
			CFLAG = 1;
		else
			XANS1 = x1;
			XANS2 = x2;
			DIST1 = sqrt((XANS1-x1)^2 + (YANS1-y1)^2);
			DIST2 = sqrt((XANS2-x1)^2 + (YANS2-y1)^2);
			[Y,I] = min([DIST1 DIST2]);

			if (Y == 1)
				GX1 = XANS1; GY1 = YANS1;
			else
				GX1 = XANS2; GY1 = YANS2;
			end
		end
	end
	if (CFLAG == 0)
		xx = x2;	yy = y2;
% REGION C
		A = 1;	B = 2;
		[e,f] = find (( ((xx-0.56-shift)./A).^2 + ((yy-0.1)./B).^2) <= 0.000625);

% IF NOT NULL THEN (X2,Y2) ARE INSIDE REGION C.
		if (~isempty(f))
			GX2 = x2;
			GY2 = y2;
		else
			T1 = sqrt(0.000625 - ((xx-0.56-shift)./A).^2);
			YANS1 = (B*T1) + 0.1;
			YANS2 = (B*-T1) + 0.1;
			XANS1 = x1;
			XANS2 = x2;
			DIST1 = sqrt((XANS1-x1)^2 + (YANS1-y1)^2);
			DIST2 = sqrt((XANS2-x1)^2 + (YANS2-y1)^2);
			[Y,I] = min([DIST1 DIST2]);

			if (Y == 1)
				GX2 = XANS2; GY2 = YANS2;
			else
				GX2 = XANS1; GY2 = YANS1;
			end
		end
		CDIST = sqrt((GX1-GX2)^2 + (GY1-GY2)^2);
	end
elseif (((x1-x2) ~= 0) & ((y1-y2) == 0))
	xx = x1;	yy = y1;
% REGION C
	A = 1;	B = 2;
	[e,f] = find (( ((xx-0.56-shift)./A).^2 + ((yy-0.1)./B).^2) <= 0.000625);

	if (~isempty(f))	% IF NOT NULL THEN (X1,Y1) ARE INSIDE REGION C.
		GX1 = x1;
		GY1 = y1;
	else
		T1 = sqrt(0.000625 - ((yy-0.1)./B).^2);
		XANS1 = (A*T1) + 0.56 + shift;
		XANS2 = (A*-T1) + 0.56 + shift;
		if ((imag(XANS1) ~= 0) & (imag(XANS2) ~= 0))	% NEVER INTERSECT.
			CDIST = 0;
			CFLAG = 1;
		else
			YANS1 = y1;
			YANS2 = y2;
			DIST1 = sqrt((XANS1-x1)^2 + (YANS1-y1)^2);
			DIST2 = sqrt((XANS2-x1)^2 + (YANS2-y1)^2);
			[Y,I] = min([DIST1 DIST2]);

			if (Y == 1)
				GX1 = XANS1; GY1 = YANS1;
			else
				GX1 = XANS2; GY1 = YANS2;
			end
		end
	end
	if (CFLAG == 0)
		xx = x2;	yy = y2;
% REGION C
		A = 1;	B = 2;
		[e,f] = find (( ((xx-0.56-shift)./A).^2 + ((yy-0.1)./B).^2) <= 0.000625);

% IF NOT NULL THEN (X2,Y2) ARE INSIDE REGION C.
		if (~isempty(f))
			GX2 = x2;
			GY2 = y2;
		else
			T1 = sqrt(0.000625 - ((yy-0.1)./B).^2);
			XANS1 = (A*T1) + 0.56 + shift;
			XANS2 = (A*-T1) + 0.56 + shift;
			YANS1 = y1;
			YANS2 = y2;
			DIST1 = sqrt((XANS1-x1)^2 + (YANS1-y1)^2);
			DIST2 = sqrt((XANS2-x1)^2 + (YANS2-y1)^2);
			[Y,I] = min([DIST1 DIST2]);

			if (Y == 1)
				GX2 = XANS2; GY2 = YANS2;
			else
				GX2 = XANS1; GY2 = YANS1;
			end
		end
		CDIST = sqrt((GX1-GX2)^2 + (GY1-GY2)^2);
	end
else	% x1 == x2 & y1 == y2
	xx = x1;	yy = y1;
% REGION C
	A = 1;	B = 2;
	[e,f] = find (( ((xx-0.56-shift)./A).^2 + ((yy-0.1)./B).^2) <= 0.000625);

% IF NOT NULL THEN (X2,Y2) ARE INSIDE REGION C.
	if (~isempty(f))
		CDIST = 999;
	else
		CDIST = 0;
	end
end
% CDIST

% ========== DDDDDDDDDD ==========

DFLAG = 0;
if (((x1-x2) ~= 0) & ((y1-y2) ~= 0))
	M = (y1-y2)/(x1-x2);
	YC = y1-M*x1;
	xx = x1;	yy = y1;
% REGION D
	[e,f] = find (((xx-0.5-shift).^2 + (yy-0.25).^2) <= 0.003+grow_hrt/8);

	if (~isempty(f))	% IF NOT NULL THEN (X1,Y1) ARE INSIDE REGION D.
		GX1 = x1;
		GY1 = y1;
	else
		AA = 1 + (M^2);
		BB = (2*(-0.5-shift)) + 2*M*(YC-0.25);
		CC = ((-0.5-shift)^2) + ((YC-0.25)^2) - 0.003 - grow_hrt/8;
		XANS1 = (-BB + sqrt(BB^2 - 4*AA*CC))/(2*AA);
		XANS2 = (-BB - sqrt(BB^2 - 4*AA*CC))/(2*AA);
		if ((imag(XANS1) ~= 0) & (imag(XANS2) ~= 0))	% NEVER INTERSECT.
			DDIST = 0;
			DFLAG = 1;
		else
			YANS1 = M*XANS1 + YC;
			YANS2 = M*XANS2 + YC;
			DIST1 = sqrt((XANS1-x1)^2 + (YANS1-y1)^2);
			DIST2 = sqrt((XANS2-x1)^2 + (YANS2-y1)^2);
			[Y,I] = min([DIST1 DIST2]);

			if (Y == 1)
				GX1 = XANS1; GY1 = YANS1;
			else
				GX1 = XANS2; GY1 = YANS2;
			end
		end
	end
	if (DFLAG == 0)
		xx = x2;	yy = y2;
% REGION D
		[e,f] = find (((xx-0.5-shift).^2 + (yy-0.25).^2) <= 0.003+grow_hrt/8);

% IF NOT NULL THEN (X2,Y2) ARE INSIDE REGION D.
		if (~isempty(f))
			GX2 = x2;
			GY2 = y2;
		else
			AA = 1 + (M^2);
			BB = (2*(-0.5-shift)) + 2*M*(YC-0.25);
			CC = ((-0.5-shift)^2) + ((YC-0.25)^2) - 0.003 - grow_hrt/8;
			XANS1 = (-BB + sqrt(BB^2 - 4*AA*CC))/(2*AA);
			XANS2 = (-BB - sqrt(BB^2 - 4*AA*CC))/(2*AA);
			YANS1 = M*XANS1 + YC;
			YANS2 = M*XANS2 + YC;
			DIST1 = sqrt((XANS1-x1)^2 + (YANS1-y1)^2);
			DIST2 = sqrt((XANS2-x1)^2 + (YANS2-y1)^2);
			[Y,I] = min([DIST1 DIST2]);

			if (Y == 1)
				GX2 = XANS2; GY2 = YANS2;
			else
				GX2 = XANS1; GY2 = YANS1;
			end
		end
		DDIST = sqrt((GX1-GX2)^2 + (GY1-GY2)^2);
	end
elseif (((x1-x2) == 0) & ((y1-y2) ~= 0))
	xx = x1;	yy = y1;
% REGION D
	[e,f] = find (((xx-0.5-shift).^2 + (yy-0.25).^2) <= 0.003+grow_hrt/8);

	if (~isempty(f))	% IF NOT NULL THEN (X1,Y1) ARE INSIDE REGION D.
		GX1 = x1;
		GY1 = y1;
	else
		T1 = sqrt((0.003+grow_hrt/8) - (xx-0.5-shift).^2);
		YANS1 = T1 + 0.25;
		YANS2 = -T1 + 0.25;
		if ((imag(YANS1) ~= 0) & (imag(YANS2) ~= 0))	% NEVER INTERSECT.
			DDIST = 0;
			DFLAG = 1;
		else
			XANS1 = x1;
			XANS2 = x2;
			DIST1 = sqrt((XANS1-x1)^2 + (YANS1-y1)^2);
			DIST2 = sqrt((XANS2-x1)^2 + (YANS2-y1)^2);
			[Y,I] = min([DIST1 DIST2]);

			if (Y == 1)
				GX1 = XANS1; GY1 = YANS1;
			else
				GX1 = XANS2; GY1 = YANS2;
			end
		end
	end
	if (DFLAG == 0)
		xx = x2;	yy = y2;
% REGION D
		[e,f] = find (((xx-0.5-shift).^2 + (yy-0.25).^2) <= 0.003+grow_hrt/8);

% IF NOT NULL THEN (X2,Y2) ARE INSIDE REGION D.
		if (~isempty(f))
			GX2 = x2;
			GY2 = y2;
		else
			T1 = sqrt((0.003+grow_hrt/8) - (xx-0.5-shift).^2);
			YANS1 = T1 + 0.25;
			YANS2 = -T1 + 0.25;
			XANS1 = x1;
			XANS2 = x2;
			DIST1 = sqrt((XANS1-x1)^2 + (YANS1-y1)^2);
			DIST2 = sqrt((XANS2-x1)^2 + (YANS2-y1)^2);
			[Y,I] = min([DIST1 DIST2]);

			if (Y == 1)
				GX2 = XANS2; GY2 = YANS2;
			else
				GX2 = XANS1; GY2 = YANS1;
			end
		end
		DDIST = sqrt((GX1-GX2)^2 + (GY1-GY2)^2);
	end
elseif (((x1-x2) ~= 0) & ((y1-y2) == 0))
	xx = x1;	yy = y1;
% REGION D
	[e,f] = find (((xx-0.5-shift).^2 + (yy-0.25).^2) <= 0.003+grow_hrt/8);

	if (~isempty(f))	% IF NOT NULL THEN (X1,Y1) ARE INSIDE REGION D.
		GX1 = x1;
		GY1 = y1;
	else
		T1 = sqrt((0.003+grow_hrt/8) - (yy-0.25).^2);
		XANS1 = T1 + 0.5 + shift;
		XANS2 = -T1 + 0.5 + shift;
		if ((imag(XANS1) ~= 0) & (imag(XANS2) ~= 0))	% NEVER INTERSECT.
			DDIST = 0;
			DFLAG = 1;
		else
			YANS1 = y1;
			YANS2 = y2;
			DIST1 = sqrt((XANS1-x1)^2 + (YANS1-y1)^2);
			DIST2 = sqrt((XANS2-x1)^2 + (YANS2-y1)^2);
			[Y,I] = min([DIST1 DIST2]);

			if (Y == 1)
				GX1 = XANS1; GY1 = YANS1;
			else
				GX1 = XANS2; GY1 = YANS2;
			end
		end
	end
	if (DFLAG == 0)
		xx = x2;	yy = y2;
% REGION D
		[e,f] = find (((xx-0.5-shift).^2 + (yy-0.25).^2) <= 0.003+grow_hrt/8);

% IF NOT NULL THEN (X2,Y2) ARE INSIDE REGION D.
		if (~isempty(f))
			GX2 = x2;
			GY2 = y2;
		else
			T1 = sqrt((0.003+grow_hrt/8) - (yy-0.25).^2);
			XANS1 = T1 + 0.5 + shift;
			XANS2 = -T1 + 0.5 + shift;
			YANS1 = y1;
			YANS2 = y2;
			DIST1 = sqrt((XANS1-x1)^2 + (YANS1-y1)^2);
			DIST2 = sqrt((XANS2-x1)^2 + (YANS2-y1)^2);
			[Y,I] = min([DIST1 DIST2]);

			if (Y == 1)
				GX2 = XANS2; GY2 = YANS2;
			else
				GX2 = XANS1; GY2 = YANS1;
			end
		end
		DDIST = sqrt((GX1-GX2)^2 + (GY1-GY2)^2);
	end
else	% x1 == x2 & y1 == y2
	xx = x1;	yy = y1;
% REGION D
	[e,f] = find (((xx-0.5-shift).^2 + (yy-0.25).^2) <= 0.003+grow_hrt/8);

% IF NOT NULL THEN (X2,Y2) ARE INSIDE REGION D.
	if (~isempty(f))
		DDIST = 999;
	else
		DDIST = 0;
	end
end
% DDIST

% ========== EEEEEEEEEE ==========

EFLAG = 0;
y1s = y1 - 0.420 - (sqrt(0.008+grow_lng/20) - sqrt(0.008))*1.3;
x1s = x1 - 0.308 - shift;
x1a = x1s*cos(pi/3)+y1s*sin(pi/3);
y1a = -x1s*sin(pi/3)+y1s*cos(pi/3);
y2s = y2 - 0.420 - (sqrt(0.008+grow_lng/20) - sqrt(0.008))*1.3;
x2s = x2 - 0.308 - shift;
x2a = x2s*cos(pi/3)+y2s*sin(pi/3);
y2a = -x2s*sin(pi/3)+y2s*cos(pi/3);

if (((x1a-x2a) ~= 0) & ((y1a-y2a) ~= 0))
	M = (y1a-y2a)/(x1a-x2a);
	YC = y1a-M*x1a;
	xx = x1a;	yy = y1a;
	xs = x1a;	ys = y1a;
% REGION E
	[e,f] = find( ( (xx./2.2).^2 + (yy./1.3).^2 )  <= 0.008+grow_lng/20);

	if (~isempty(f))	% IF NOT NULL THEN (X1,Y1) ARE INSIDE REGION E.
		GX1 = xs;
		GY1 = ys;
	else
		AA = 1/(2.2^2) + (M^2)/(1.3^2);
		BB = 2*M*YC/(1.3^2);
		CC = (YC^2)/(1.3^2) - 0.008 - grow_lng/20;
		XANS1 = (-BB + sqrt(BB^2 - 4*AA*CC))/(2*AA);
		XANS2 = (-BB - sqrt(BB^2 - 4*AA*CC))/(2*AA);
		if ((imag(XANS1) ~= 0) & (imag(XANS2) ~= 0))	% NEVER INTERSECT.
			EDIST = 0;
			EFLAG = 1;
		else
			YANS1 = M*XANS1 + YC;
			YANS2 = M*XANS2 + YC;
			DIST1 = sqrt((XANS1-xs)^2 + (YANS1-ys)^2);
			DIST2 = sqrt((XANS2-xs)^2 + (YANS2-ys)^2);
			[Y,I] = min([DIST1 DIST2]);

			if (Y == 1)
				GX1 = XANS1; GY1 = YANS1;
			else
				GX1 = XANS2; GY1 = YANS2;
			end
		end
	end
	if (EFLAG == 0)
		xs = x2a;	ys = y2a;
		xx = x2a;	yy = y2a;
% REGION E
		[e,f] = find( ( (xx./2.2).^2 + (yy./1.3).^2 )  <= 0.008+grow_lng/20);

% IF NOT NULL THEN (X2,Y2) ARE INSIDE REGION E.
		if (~isempty(f))
			GX2 = xs;
			GY2 = ys;
		else
			AA = 1/(2.2^2) + (M^2)/(1.3^2);
			BB = 2*M*YC/(1.3^2);
			CC = (YC^2)/(1.3^2) - 0.008 - grow_lng/20;
			XANS1 = (-BB + sqrt(BB^2 - 4*AA*CC))/(2*AA);
			XANS2 = (-BB - sqrt(BB^2 - 4*AA*CC))/(2*AA);
			YANS1 = M*XANS1 + YC;
			YANS2 = M*XANS2 + YC;
			DIST1 = sqrt((XANS1-xs)^2 + (YANS1-ys)^2);
			DIST2 = sqrt((XANS2-xs)^2 + (YANS2-ys)^2);
			[Y,I] = min([DIST1 DIST2]);

			if (Y == 1)
				GX2 = XANS2; GY2 = YANS2;
			else
				GX2 = XANS1; GY2 = YANS1;
			end
		end
		EDIST = sqrt((GX1-GX2)^2 + (GY1-GY2)^2);
	end
elseif (((x1a-x2a) == 0) & ((y1a-y2a) ~= 0))
	xs = x1a;	ys = y1a;
	xx = x1a;	yy = y1a;
% REGION E
	[e,f] = find( ( (xx./2.2).^2 + (yy./1.3).^2 )  <= 0.008+grow_lng/20);

	if (~isempty(f))	% IF NOT NULL THEN (X1,Y1) ARE INSIDE REGION E.
		GX1 = xs;
		GY1 = ys;
	else
		T1 = sqrt((0.008+grow_lng/20) - (xx/2.2).^2);
		YANS1 = 1.3*T1;
		YANS2 = 1.3*-T1;
		if ((imag(YANS1) ~= 0) & (imag(YANS2) ~= 0))	% NEVER INTERSECT.
			EDIST = 0;
			EFLAG = 1;
		else
			XANS1 = xs;
			XANS2 = xs;
			DIST1 = sqrt((XANS1-xs)^2 + (YANS1-ys)^2);
			DIST2 = sqrt((XANS2-xs)^2 + (YANS2-ys)^2);
			[Y,I] = min([DIST1 DIST2]);

			if (Y == 1)
				GX1 = XANS1; GY1 = YANS1;
			else
				GX1 = XANS2; GY1 = YANS2;
			end
		end
	end
	if (EFLAG == 0)
		xs = x2a;	ys = y2a;
		xx = x2a;	yy = y2a;
% REGION E
		[e,f] = find( ( (xx./2.2).^2 + (yy./1.3).^2 )  <= 0.008+grow_lng/20);

% IF NOT NULL THEN (X2,Y2) ARE INSIDE REGION E.
		if (~isempty(f))
			GX2 = xs;
			GY2 = ys;
		else
			T1 = sqrt((0.008+grow_lng/20) - (xx./2.2).^2);
			YANS1 = 1.3*T1;
			YANS2 = 1.3*-T1;
			XANS1 = xs;
			XANS2 = xs;
			DIST1 = sqrt((XANS1-xs)^2 + (YANS1-ys)^2);
			DIST2 = sqrt((XANS2-xs)^2 + (YANS2-ys)^2);
			[Y,I] = min([DIST1 DIST2]);

			if (Y == 1)
				GX2 = XANS2; GY2 = YANS2;
			else
				GX2 = XANS1; GY2 = YANS1;
			end
		end
		EDIST = sqrt((GX1-GX2)^2 + (GY1-GY2)^2);
	end
elseif (((x1a-x2a) ~= 0) & ((y1a-y2a) == 0))
	xx = x1a;	yy = y1a;
	xs = x1a;	ys = y1a;
% REGION E
	[e,f] = find( ( (xx./2.2).^2 + (yy./1.3).^2 )  <= 0.008+grow_lng/20);

	if (~isempty(f))	% IF NOT NULL THEN (X1,Y1) ARE INSIDE REGION E.
		GX1 = xs;
		GY1 = ys;
	else
		T1 = sqrt((0.008+grow_lng/20) - (yy./1.3).^2);
		XANS1 = 2.2*T1;
		XANS2 = 2.2*-T1;
		if ((imag(XANS1) ~= 0) & (imag(XANS2) ~= 0))	% NEVER INTERSECT.
			EDIST = 0;
			EFLAG = 1;
		else
			YANS1 = ys;
			YANS2 = ys;
			DIST1 = sqrt((XANS1-xs)^2 + (YANS1-ys)^2);
			DIST2 = sqrt((XANS2-xs)^2 + (YANS2-ys)^2);
			[Y,I] = min([DIST1 DIST2]);

			if (Y == 1)
				GX1 = XANS1; GY1 = YANS1;
			else
				GX1 = XANS2; GY1 = YANS2;
			end
		end
	end
	if (EFLAG == 0)
		xx = x2a;	yy = y2a;
		xs = x2a;	ys = y2a;
% REGION E
		[e,f] = find( ( (xx./2.2).^2 + (yy./1.3).^2 )  <= 0.008+grow_lng/20);

% IF NOT NULL THEN (X2,Y2) ARE INSIDE REGION E.
		if (~isempty(f))
			GX2 = xs;
			GY2 = ys;
		else
			T1 = sqrt((0.008+grow_lng/20) - (yy./1.3).^2);
			XANS1 = 2.2*T1;
			XANS2 = 2.2*-T1;
			YANS1 = ys;
			YANS2 = ys;
			DIST1 = sqrt((XANS1-xs)^2 + (YANS1-ys)^2);
			DIST2 = sqrt((XANS2-xs)^2 + (YANS2-ys)^2);
			[Y,I] = min([DIST1 DIST2]);

			if (Y == 1)
				GX2 = XANS2; GY2 = YANS2;
			else
				GX2 = XANS1; GY2 = YANS1;
			end
		end
		EDIST = sqrt((GX1-GX2)^2 + (GY1-GY2)^2);
	end
else	% x1a == x2a & y1a == y2a
	xx = x1a;	yy = y1a;
	xs = x1a;	ys = y1a;
% REGION E
	[e,f] = find( ( (xx./2.2).^2 + (yy./1.3).^2 )  <= 0.008+grow_lng/20);

% IF NOT NULL THEN (X2,Y2) ARE INSIDE REGION E.
	if (~isempty(f))
		EDIST = 999;
	else
		EDIST = 0;
	end
end
% EDIST


% ========== FFFFFFFFFF ==========

FFLAG = 0;
y1s = y1 - 0.420 - (sqrt(0.008+grow_lng/20) - sqrt(0.008))*1.3;
x1s = x1 - 0.692 - shift;
x1a = x1s*cos(2*pi/3)+y1s*sin(2*pi/3);
y1a = -x1s*sin(2*pi/3)+y1s*cos(2*pi/3);
y2s = y2 - 0.420 - (sqrt(0.008+grow_lng/20) - sqrt(0.008))*1.3;
x2s = x2 - 0.692 - shift;
x2a = x2s*cos(2*pi/3)+y2s*sin(2*pi/3);
y2a = -x2s*sin(2*pi/3)+y2s*cos(2*pi/3);

if (((x1a-x2a) ~= 0) & ((y1a-y2a) ~= 0))
	M = (y1a-y2a)/(x1a-x2a);
	YC = y1a-M*x1a;
	xs = x1a;	ys = y1a;
	xx = x1a;	yy = y1a;
% REGION F
	[e,f] = find( ( (xx./2.2).^2 + (yy./1.3).^2 )  <= 0.008+grow_lng/20);

	if (~isempty(f))	% IF NOT NULL THEN (X1,Y1) ARE INSIDE REGION F.
		GX1 = xs;
		GY1 = ys;
	else
		AA = 1/(2.2^2) + (M^2)/(1.3^2);
		BB = 2*M*YC/(1.3^2);
		CC = (YC^2)/(1.3^2) - 0.008 - grow_lng/20;
		XANS1 = (-BB + sqrt(BB^2 - 4*AA*CC))/(2*AA);
		XANS2 = (-BB - sqrt(BB^2 - 4*AA*CC))/(2*AA);
		if ((imag(XANS1) ~= 0) & (imag(XANS2) ~= 0))	% NEVER INTERSECT.
			FDIST = 0;
			FFLAG = 1;
		else
			YANS1 = M*XANS1 + YC;
			YANS2 = M*XANS2 + YC;
			DIST1 = sqrt((XANS1-xs)^2 + (YANS1-ys)^2);
			DIST2 = sqrt((XANS2-xs)^2 + (YANS2-ys)^2);
			[Y,I] = min([DIST1 DIST2]);

			if (Y == 1)
				GX1 = XANS1; GY1 = YANS1;
			else
				GX1 = XANS2; GY1 = YANS2;
			end
		end
	end
	if (FFLAG == 0)
		xs = x2a;	ys = y2a;
		xx = x2a;	yy = y2a;
% REGION F
		[e,f] = find( ( (xx./2.2).^2 + (yy./1.3).^2 )  <= 0.008+grow_lng/20);

% IF NOT NULL THEN (X2,Y2) ARE INSIDE REGION F.
		if (~isempty(f))
			GX2 = xs;
			GY2 = ys;
		else
			AA = 1/(2.2^2) + (M^2)/(1.3^2);
			BB = 2*M*YC/(1.3^2);
			CC = (YC^2)/(1.3^2) - 0.008 - grow_lng/20;
			XANS1 = (-BB + sqrt(BB^2 - 4*AA*CC))/(2*AA);
			XANS2 = (-BB - sqrt(BB^2 - 4*AA*CC))/(2*AA);
			YANS1 = M*XANS1 + YC;
			YANS2 = M*XANS2 + YC;
			DIST1 = sqrt((XANS1-xs)^2 + (YANS1-ys)^2);
			DIST2 = sqrt((XANS2-xs)^2 + (YANS2-ys)^2);
			[Y,I] = min([DIST1 DIST2]);

			if (Y == 1)
				GX2 = XANS2; GY2 = YANS2;
			else
				GX2 = XANS1; GY2 = YANS1;
			end
		end
		FDIST = sqrt((GX1-GX2)^2 + (GY1-GY2)^2);
	end
elseif (((x1a-x2a) == 0) & ((y1a-y2a) ~= 0))
	xs = x1a;	ys = y1a;
	xx = x1a;	yy = y1a;
% REGION F
	[e,f] = find( ( (xx./2.2).^2 + (yy./1.3).^2 )  <= 0.008+grow_lng/20);

	if (~isempty(f))	% IF NOT NULL THEN (X1,Y1) ARE INSIDE REGION F.
		GX1 = xs;
		GY1 = ys;
	else
		T1 = sqrt((0.008+grow_lng/20) - (xx/2.2).^2);
		YANS1 = 1.3*T1;
		YANS2 = 1.3*-T1;
		if ((imag(YANS1) ~= 0) & (imag(YANS2) ~= 0))	% NEVER INTERSECT.
			FDIST = 0;
			FFLAG = 1;
		else
			XANS1 = xs;
			XANS2 = xs;
			DIST1 = sqrt((XANS1-xs)^2 + (YANS1-ys)^2);
			DIST2 = sqrt((XANS2-xs)^2 + (YANS2-ys)^2);
			[Y,I] = min([DIST1 DIST2]);

			if (Y == 1)
				GX1 = XANS1; GY1 = YANS1;
			else
				GX1 = XANS2; GY1 = YANS2;
			end
		end
	end
	if (FFLAG == 0)
		xs = x2a;	ys = y2a;
		xx = x2a;	yy = y2a;
% REGION F
		[e,f] = find( ( (xx./2.2).^2 + (yy./1.3).^2 )  <= 0.008+grow_lng/20);

% IF NOT NULL THEN (X2,Y2) ARE INSIDE REGION F.
		if (~isempty(f))
			GX2 = xs;
			GY2 = ys;
		else
			T1 = sqrt((0.008+grow_lng/20) - (xx./2.2).^2);
			YANS1 = 1.3*T1;
			YANS2 = 1.3*-T1;
			XANS1 = xs;
			XANS2 = xs;
			DIST1 = sqrt((XANS1-xs)^2 + (YANS1-ys)^2);
			DIST2 = sqrt((XANS2-xs)^2 + (YANS2-ys)^2);
			[Y,I] = min([DIST1 DIST2]);

			if (Y == 1)
				GX2 = XANS2; GY2 = YANS2;
			else
				GX2 = XANS1; GY2 = YANS1;
			end
		end
		FDIST = sqrt((GX1-GX2)^2 + (GY1-GY2)^2);
	end
elseif (((x1a-x2a) ~= 0) & ((y1a-y2a) == 0))
	xx = x1a;	yy = y1a;
	xs = x1a;	ys = y1a;
% REGION F
	[e,f] = find( ( (xx./2.2).^2 + (yy./1.3).^2 )  <= 0.008+grow_lng/20);

	if (~isempty(f))	% IF NOT NULL THEN (X1,Y1) ARE INSIDE REGION F.
		GX1 = xs;
		GY1 = ys;
	else
		T1 = sqrt((0.008+grow_lng/20) - (yy./1.3).^2);
		XANS1 = 2.2*T1;
		XANS2 = 2.2*-T1;
		if ((imag(XANS1) ~= 0) & (imag(XANS2) ~= 0))	% NEVER INTERSECT.
			FDIST = 0;
			FFLAG = 1;
		else
			YANS1 = ys;
			YANS2 = ys;
			DIST1 = sqrt((XANS1-xs)^2 + (YANS1-ys)^2);
			DIST2 = sqrt((XANS2-xs)^2 + (YANS2-ys)^2);
			[Y,I] = min([DIST1 DIST2]);

			if (Y == 1)
				GX1 = XANS1; GY1 = YANS1;
			else
				GX1 = XANS2; GY1 = YANS2;
			end
		end
	end
	if (FFLAG == 0)
		xx = x2a;	yy = y2a;
		xs = x2a;	ys = y2a;
% REGION F
		[e,f] = find( ( (xx./2.2).^2 + (yy./1.3).^2 )  <= 0.008+grow_lng/20);

% IF NOT NULL THEN (X2,Y2) ARE INSIDE REGION F.
		if (~isempty(f))
			GX2 = xs;
			GY2 = ys;
		else
			T1 = sqrt((0.008+grow_lng/20) - (yy./1.3).^2);
			XANS1 = 2.2*T1;
			XANS2 = 2.2*-T1;
			YANS1 = ys;
			YANS2 = ys;
			DIST1 = sqrt((XANS1-xs)^2 + (YANS1-ys)^2);
			DIST2 = sqrt((XANS2-xs)^2 + (YANS2-ys)^2);
			[Y,I] = min([DIST1 DIST2]);

			if (Y == 1)
				GX2 = XANS2; GY2 = YANS2;
			else
				GX2 = XANS1; GY2 = YANS1;
			end
		end
		FDIST = sqrt((GX1-GX2)^2 + (GY1-GY2)^2);
	end
else	% x1a == x2a & y1a == y2a
	xx = x1a;	yy = y1a;
	xs = x1a;	ys = y1a;
% REGION F
	[e,f] = find( ( (xx./2.2).^2 + (yy./1.3).^2 )  <= 0.008+grow_lng/20);

% IF NOT NULL THEN (X2,Y2) ARE INSIDE REGION F.
	if (~isempty(f))
		FDIST = 999;
	else
		FDIST = 0;
	end
end
% FDIST

integral = 0;
if (GDIST == 999)
	integral = DG;
	if (ADIST == 999)
		integral = DA;
	elseif (BDIST == 999)
		integral = DB;
	elseif (CDIST == 999)
		integral = DC;
	elseif (DDIST == 999)
		integral = DD;
	elseif (EDIST == 999)
		integral = DE;
	elseif (FDIST == 999)
		integral = DF;
	end
elseif (GDIST ~= 0)
	GGDIST = GDIST - ADIST - BDIST - CDIST - DDIST - EDIST - FDIST;
	integral = GGDIST*DG + ADIST*DA + BDIST*DB + CDIST*DC + DDIST*DD + EDIST*DE + FDIST*DF;
end

break;

GDIST
ADIST
BDIST
CDIST
DDIST
EDIST
FDIST