unit OvUnit;
{
calculate areas A, B and AB overlap
}
interface
uses dialogs,sysutils,graphics;
type TrealPoint = record
x,y : double;
end;
Tedge = record
x1,y1,x2,y2 : double; //start and finish of edge
end;
TPoly12 = record
angles : byte;
coord :array[1..12] of TrealPoint;
end;
TEdge12 = record
edges : byte;
coord : array[1..12] of TEdge;
end;
TIntersection = record
fa,fb : double; //mult factor of vector to intersection
valid : boolean; //false if parallel
end;
TISlist = record //intersection list of A,B
Acount : byte;
Bcount : byte;
isdata : array[1..12,1..12] of TIntersection;
end;
TInnerPoint = record
x,y : double; //intersection coordinates
ia,ib : byte; //edges with (x,y)
end;
TIPlist = record
top : byte;
ipdata : array[1..50] of TInnerpoint;
end;
TTrack = record
x,y : double; //ref to point in TInnerpoint list
Ematch : LongInt; //bit set = point on this edge
visited : boolean; //used in innerpath
end; //B points index + 12
TTrackList = record
top : byte;
trackdata : array[1..50] of TTrack;
end;
TPolyList = record
top : byte;
coord : array[1..50] of TRealpoint;
end;
function overlap(var pA, pB : TPoly12) : double;
function PolyArea(var poly : TPoly12) : double;
function TestConvex(var poly : TPoly12) : boolean;
var Apoints : TPoly12; //angle coordinates of A
Bpoints : TPoly12; //angle coordinates of B
APosOK : boolean = false; //valid shape of A
BPosOK : boolean = false; //...............B
polyList : TPolyList; //final angle coordinates of inner polygon
implementation
uses unit1,paramunit;
{--------forward declarations------------}
procedure makeEdgeList(var Elist:Tedge12 ; crd: TPoly12);forward;
procedure makeIntersections;forward;
procedure makeIPlist;forward;
procedure makeTrackList;forward;
procedure makePolyList;forward;
function calcArea : double;forward;
function distance(i,j : byte) : double;forward;
procedure paintPolylist;forward;
procedure addtrackpoint(n : byte);forward;
function Polydistance(var po : TPoly12; i,j : byte) : double;forward;
{------------ tables used during calculation of areas -------}
var Aedges : TEdge12; //list of edges of polygon A
Bedges : TEdge12; //.........................B
ISList : TISlist; //table with intersection data
IPlist : TIPlist; //list of innerpoints, point may appear > once
TrackList : TTracklist; //bit encoded edges of inner polygon
{------------ main calls -------------}
function PolyArea(var poly : TPoly12) : double;
//calculate area of polygon
//this is done by breaking down the polygon into triangles
//with edges a,b,c and s = (a+b+c)/2
//apply herons lemma: area = sqrt(s*(s-a)*(s-b)*(s-c))
var a,b,c,s : double;
i : byte;
begin
result := 0;
with poly do
for i := 1 to angles-2 do
begin
a := Polydistance(poly,1,i+1);
b := Polydistance(poly,i+1,i+2);
c := Polydistance(poly,1,i+2);
s := (a+b+c)2;
result := result + sqrt(s*(s-a)*(s-b)*(s-c));
end;
end;
function Testconvex(var poly : TPoly12) : boolean;
//return true is polygon is convex
var i,j,iscount : byte;
b1,b2 : boolean;
begin
makeEdgeList(Aedges,poly);
makeEdgeList(Bedges,poly);
makeIntersections;
for j := 1 to ISList.Acount do
begin
iscount := 0;
for i := 1 to ISList.Acount do
with ISList.isdata[i,j] do
if valid then
begin
b1 := (fa >= 0) and (fa <= 1);
b2 := (fb >= 0) and (fb <= 1);
if b1 or b2 then inc(iscount);
end;
if iscount > 2 then //may intersect with only 2 other edges
begin
result := false;
exit;
end;
end;//for
result := true;
end;
function Overlap(var pA, pB : TPoly12) : double;
//calculate A,B areas and overlap
//A,B,shape must be OK
begin
result := 0;
makeEdgeList(Aedges , pA); //list of projected edges
makeEdgeList(Bedges , pB); //list of reference rectangle edges
makeIntersections;
makeIPlist;
makeTrackList; //list of inner points
makePolyList; //final result, global list
if PolyList.Top >= 3 then
begin
paintpolyList;
result := calcArea; //area of polylist
end;
end;
{----------- sub procs of proc. "overlap" ----------}
procedure makeEdgeList(var Elist:Tedge12 ; crd: TPoly12);
//use Tpoly12 array to build EdgeList
var i,n : byte;
begin
n := crd.angles;
Elist.edges := n;
for i := 1 to n do
begin
Elist.coord[i].x1 := crd.coord[i].x;
Elist.coord[i].y1 := crd.coord[i].y;
if i < n then
begin
Elist.coord[i].x2 := crd.coord[i+1].x;
Elist.coord[i].y2 := crd.coord[i+1].y;
end
else begin
Elist.coord[i].x2 := crd.coord[1].x;
Elist.coord[i].y2 := crd.coord[1].y;
end;
end;//for
end;
procedure makeIntersections;
//use edge lists of poly12
//to build the intersections list
var i,j : byte;
a,b,c,d,e,f : double;
x,y,discr : double;
begin
with ISlist do
begin
Acount := Aedges.edges;
Bcount := Bedges.edges;
end;
for i := 1 to ISlist.Acount do //find intersections of edges
for j := 1 to ISlist.Bcount do
begin
a := Aedges.coord[i].x2 - Aedges.coord[i].x1;
b := Bedges.coord[j].x1 - Bedges.coord[j].x2; {ax + by = e}
c := Aedges.coord[i].y2 - Aedges.coord[i].y1; {cx + dy = f}
d := Bedges.coord[j].y1 - Bedges.coord[j].y2; {fa=x,fb=y}
e := Bedges.coord[j].x1 - Aedges.coord[i].x1;
f := Bedges.coord[j].y1 - Aedges.coord[i].y1;
discr := a*d-b*c;
if discr <> 0 then
begin
ISList.isdata[i,j].fa := (e*d-b*f)discr; //vector mult factor
ISList.isdata[i,j].fb := (a*f-c*e)discr; //....
ISList.isdata[i,j].valid := true;
end
else
begin
ISList.isdata[i,j].valid := false; //edges parallel or coincident
end;
end;//for i,j
end;
procedure makeIPlist;
//use intersections to calculate IPlist (inner points list)
var i,j,n : byte;
Poscode : word;//1xx:minus out; x1x:inner; xx1: positive out
IPhit : boolean;
begin
for i := 1 to 50 do //clear list
with IPlist.ipdata[i] do
begin
ia := 0;
ib := 0;
end;
//
n := 0;
for i := 1 to ISList.Acount do //find intersections
for j := 1 to ISList.Bcount do
with ISList.isdata[i,j] do
if valid and (fb >= 0) and (fb <= 1) and (fa >= 0) and (fa <=1) then
begin
inc(n);
with Aedges.coord[i] do
with IPList.ipdata[n] do
begin
x := x1 + fa*(x2 - x1);
y := y1 + fa*(y2 - y1);
ia := i;
ib := j;
end;
end;
//
for i := 1 to ISList.Acount do //find inner A points
begin
poscode := 0;
IPhit := false;
for j := 1 to ISList.Bcount do
with ISList.isdata[i,j] do
if valid and (fb >=0) and (fb <=1) and (IPhit = false) then
begin
if (fa > 1) then poscode := poscode or $0001
else if (fa < 0) then poscode := poscode or $0100
else poscode := poscode or $0010;
if (poscode=$0110) or (poscode=$0101) then //first point
begin
inc(n);
IPhit := true;
with IPlist.ipdata[n] do
begin
x := Aedges.coord[i].x1;
y := Aedges.coord[i].y1;
ia := i;
end;
end;//if
if (poscode=$0011) or (poscode=$0101) then
begin //second point
inc(n);
IPhit := true;
with IPlist.ipdata[n] do
begin
x := Aedges.coord[i].x2;
y := Aedges.coord[i].y2;
ia := i;
end;
end;//if
end;//if valid
end;//for j
//
for j := 1 to ISList.Bcount do //find inner B points
begin
poscode := 0;
IPhit := false;
for i := 1 to ISList.Acount do
with ISList.isdata[i,j] do
if valid and (fa >= 0) and (fa <= 1) and (IPhit = false) then
begin
if fb > 1 then poscode := poscode or $0001
else if fb < 0 then poscode := poscode or $0100
else poscode := poscode or $0010;
if (poscode=$0110) or (poscode=$0101) then //first point
begin
inc(n);
IPhit := true;
with IPlist.ipdata[n] do
begin
x := Bedges.coord[j].x1;
y := Bedges.coord[j].y1;
ib := j;
end;
end;//if
if (poscode=$0011) or (poscode=$0101) then //second point
begin
inc(n);
IPhit := true;
with IPlist.ipdata[n] do
begin
x := Bedges.coord[j].x2;
y := Bedges.coord[j].y2;
ib := j;
end;
end;//if
end;//if valid
end;//for i
IPlist.top := n;
end;
procedure makeTrackList;
//use IP list to built the trackList
var i : byte;
begin
TrackList.top := 0;
//
for i := 1 to 50 do //clear the TrackList
with TrackList.trackdata[i] do
begin
x := 0; y := 0;
Ematch := 0;
visited := false;
end;
//
for i := 1 to IPlist.top do addTrackpoint(i);
end;
procedure makePolyList;
var oldEdges : longInt; //bit set = on edge, PEdges bias + 4
trackNr,i : byte;
hit : boolean;
begin
polylist.top := 0;
for i := 1 to 50 do //clear polylist
with polylist.coord[i] do
begin
x := 0;
y := 0;
end;
if Tracklist.Top < 3 then exit;
//
Polylist.Top := 1; //set 1st point
with PolyList.coord[1] do
begin
x := TrackList.trackdata[1].x;
y := Tracklist.trackdata[1].y;
end;
trackList.trackdata[1].visited := true;
oldEdges := TrackList.trackdata[1].Ematch;
// --- find connection in TrackList---
repeat
trackNr := 0; //entry of trackList
hit := false;
repeat
inc(trackNr);
with TrackList.trackdata[trackNr] do
if (visited = false) and (OldEdges and EMatch <> 0) then hit := true;
until hit or (trackNr = Tracklist.Top);
//
if hit then
with trackList.trackdata[tracknr] do
begin
visited := true;
oldEdges := Ematch;
inc(Polylist.Top);
PolyList.coord[Polylist.top].x := x;
PolyList.coord[Polylist.Top].y := y;
end;
until hit = false;
end;
procedure paintPolylist;
//paint in map2, copy to paintbox
//called in test mode to show IP list
var i : byte;
xx,yy : integer;
begin
with map2 do with canvas do
begin
pen.color := $000000;
brush.color := $000000;
brush.Style := bsSolid;
for i := 1 to Polylist.Top do
with PolyList.coord[i] do
begin
xx := XToPixel(x);
yy := YToPixel(y);
ellipse(xx-4,yy-4,xx+4,yy+4);
end;
end;//with
map2ToBox;
end;
function calcArea : double;
//use polygon list to calculate overlapping area
//this is done by breaking down the polygon into triangles
//with edges a,b,c and s = (a+b+c)/2
//apply herons lemma: area = sqrt(s*(s-a)*(s-b)*(s-c))
var a,b,c,s : double;
i : byte;
begin
result := 0;
if polylist.top < 3 then exit;
//
for i := 1 to Polylist.Top-2 do
begin
a := distance(1,i+1);
b := distance(i+1,i+2);
c := distance(1,i+2);
s := (a+b+c)2;
result := result + sqrt(s*(s-a)*(s-b)*(s-c));
end;
end;
//----- help for polylist generation , area calculation -----------
procedure addtrackpoint(n : byte);
//adds entry n of IP list to the tracklist
var i,j,aa,bb : byte;
x,y : double;
hit : boolean;
begin
i := 0;
x := IPlist.ipdata[n].x;
y := IPlist.ipdata[n].y;
hit := false;
j := 0;
if tracklist.top > 0 then
repeat
inc(j);
if (Tracklist.trackdata[j].x = x) and
(TrackList.trackdata[j].y = y) then hit := true;
until hit or (j = TrackList.Top);
//
if hit = false then
begin
inc(TrackList.Top); //create new entry
j := Tracklist.Top;
Tracklist.trackdata[j].x := x;
Tracklist.trackdata[j].y := y
end;
//
aa := IPlist.ipdata[n].ia; //point on edge ia
bb := IPlist.ipdata[n].ib; //..............ib
with TrackList.trackdata[j] do
begin
if aa <> 0 then Ematch := Ematch or (1 shl aa);
if bb <> 0 then Ematch := Ematch or (1 shl (bb+12));
end;
end;
function Polydistance(var po : TPoly12; i,j : byte) : double;
var dx, dy : double;
begin
with po do
begin
dx := (coord[i].x - coord[j].x);
dy := (coord[i].y - coord[j].y);
end;
result := sqrt(dx*dx + dy*dy);
end;
function distance(i,j : byte) : double;
//pythagoras, distance between (x1,y1)...(x2,y2) of polylist
var dx, dy : double;
begin
dx := (PolyList.coord[i].x - PolyList.coord[j].x);
dy := (PolyList.coord[i].y - PolyList.coord[j].y);
result := sqrt(dx*dx + dy*dy);
end;
end.