DIRECTORY
RealFns, SV2d, SVAngle, SVLines2d, SVPolygon2d, SVVector3d;

SVPolygon2dImpl: CEDAR PROGRAM
IMPORTS RealFns, SVAngle, SVLines2d, SVVector3d
EXPORTS SVPolygon2d =
BEGIN

Point2d: TYPE = SV2d.Point2d;
Path: TYPE = REF PathObj;
PathObj: TYPE = SV2d.PathObj;
Polygon: TYPE = REF PolygonObj;
PolygonObj: TYPE = SV2d.PolygonObj;
TrigLine: TYPE = REF TrigLineObj;
TrigLineObj: TYPE = SV2d.TrigLineObj;
TrigLineSeg: TYPE = SV2d.TrigLineSeg;
TrigPolygon: TYPE = REF TrigPolygonObj;
TrigPolygonObj: TYPE = SV2d.TrigPolygonObj;
Vector3d: TYPE = SVVector3d.Vector3d;


unitSquare: Polygon;
root3, root3Over3, twoRoot3Over3: REAL;-- constants computed at load time in Init (see below).
globalVLine, globalCLine: TrigLine;


CreatePoly: PUBLIC PROC [len: NAT] RETURNS [poly: Polygon] = {
poly _ NEW[PolygonObj[len]];
poly.len _ 0;
};
CopyPoly: PUBLIC PROC [poly: Polygon] RETURNS [copy: Polygon] = {
copy _ CreatePoly[poly.maxVerts];
FOR i: NAT IN [0..poly.len) DO
copy[i] _ poly[i];
ENDLOOP;
copy.len _ poly.len;
};

CircumHexagon: PUBLIC PROC [r: REAL] RETURNS [hex: Polygon] = {
hex _ CreatePoly[6];
hex _ AddPolyPoint[hex, [r*root3Over3, r]];-- front right
hex _ AddPolyPoint[hex, [r*twoRoot3Over3, 0]];-- right
hex _ AddPolyPoint[hex, [r*root3Over3, -r]];-- back right
hex _ AddPolyPoint[hex, [-r*root3Over3, -r]];-- back left
hex _ AddPolyPoint[hex, [-r*twoRoot3Over3, 0]];-- left
hex _ AddPolyPoint[hex, [-r*root3Over3, r]];-- front left
};
ClearPoly: PUBLIC PROC [poly: Polygon] = {
poly.len _ 0;
};
AddPolyPoint: PUBLIC PROC [poly: Polygon, point: Point2d] RETURNS [polyPlusPoint: Polygon] = {
IF poly.len = poly.maxVerts THEN {
polyPlusPoint _ CreatePoly[poly.maxVerts+1];
polyPlusPoint _ PartPolyGetsPartPoly[poly, 0, polyPlusPoint, 0, poly.maxVerts];
polyPlusPoint[poly.maxVerts] _ point;
polyPlusPoint.len _ poly.maxVerts+1;
}
ELSE {
poly[poly.len] _ point;
poly.len _ poly.len + 1;
polyPlusPoint _ poly;
};
}; -- end of AddPolyPoint
PutPolyPoint: PUBLIC PROC [poly: Polygon, index: NAT, point: Point2d] RETURNS [newPoly: Polygon] = {
IF index+1 > poly.maxVerts THEN {
newPoly _ CreatePoly[index+1];
newPoly _ PartPolyGetsPartPoly[poly, 0, newPoly, 0, poly.maxVerts];
newPoly[index] _ point;
newPoly.len _ index+1;
}
ELSE {
IF index+1 > poly.len THEN poly.len _ index+1;
poly[index] _ point;
newPoly _ poly};
};
IsClockwisePoly: PUBLIC PROC [poly: Polygon] RETURNS [BOOL] = {
RETURN [SignedArea[poly]>0];
};
InvertPoly: PUBLIC PROC [poly: Polygon] RETURNS [ylop: Polygon] = {
ylop _ NEW[PolygonObj[poly.len]];
FOR i: NAT IN[0..poly.len) DO
ylop[poly.len - i -1] _ poly[i];
ENDLOOP;
};
InvertPolyInPlace: PUBLIC PROC [poly: Polygon] = {
tempPoint: Point2d;
FOR i: NAT IN[0..poly.len/2) DO
tempPoint _ poly[poly.len - i - 1];
poly[poly.len - i - 1] _ poly[i];
poly[i] _ tempPoint;
ENDLOOP;
};
PartPolyGetsPartPath: PUBLIC PROC [fromPath: Path, fromStart: NAT, toPoly: Polygon, toStart: NAT, duration: NAT] RETURNS [newPoly: Polygon] = {
j: NAT _ fromStart;
IF toStart+duration > toPoly.maxVerts THEN {
newPoly _ CreatePoly[toStart+duration];
newPoly _ PartPolyGetsPartPoly[toPoly, 0, newPoly, 0, toPoly.len];
newPoly.len _ toStart+duration}
ELSE {
newPoly _ toPoly;
IF toStart+duration > newPoly.len THEN newPoly.len _ toStart+duration};
IF fromStart+duration > fromPath.len THEN ERROR; -- don't use bad data
FOR i: NAT IN[toStart..toStart+duration-1] DO
newPoly[i] _ fromPath[j];
j _ j + 1;
ENDLOOP;
};
PartPolyGetsPartPoly: PUBLIC PROC [fromPoly: Polygon, fromStart: NAT, toPoly: Polygon, toStart: NAT, duration: NAT] RETURNS [newPoly: Polygon] = {
j: NAT _ fromStart;
IF toStart+duration > toPoly.maxVerts THEN {
newPoly _ CreatePoly[toStart+duration];
newPoly _ PartPolyGetsPartPoly[toPoly, 0, newPoly, 0, toPoly.len];
newPoly.len _ toStart+duration}
ELSE {
newPoly _ toPoly;
IF toStart+duration > newPoly.len THEN newPoly.len _ toStart+duration};
IF fromStart+duration > fromPoly.len THEN ERROR; -- don't use bad data
FOR i: NAT IN[toStart..toStart+duration-1] DO
newPoly[i] _ fromPoly[j];
j _ j + 1;
ENDLOOP;
};
PathToPolygon: PUBLIC PROC [path: Path] RETURNS [poly: Polygon] = {
poly _ NEW[PolygonObj[path.len]];
FOR i: NAT IN[0..path.len) DO
poly[i] _ path[i];
ENDLOOP;
poly.len _ path.len;
}; -- end of PathToPolygon


PolygonToTrigPolygon: PUBLIC PROC [poly: Polygon] RETURNS [trigPoly: TrigPolygon] = {
lastPoint, thisPoint: Point2d;
newSeg: TrigLineSeg;
trigPoly _ NEW[TrigPolygonObj[poly.len]];
lastPoint _ poly[0];
FOR i: NAT IN[1..poly.len) DO
thisPoint _ poly[i];
newSeg _ SVLines2d.CreateTrigLineSeg[lastPoint, thisPoint];
trigPoly[i-1] _ newSeg;
lastPoint _ thisPoint;
ENDLOOP;
thisPoint _ poly[0];
newSeg _ SVLines2d.CreateTrigLineSeg[lastPoint, thisPoint];
trigPoly[poly.len-1] _ newSeg;
};
CreatePath: PUBLIC PROC [len: NAT] RETURNS [path: Path] = {
path _ NEW[PathObj[len]];
path.len _ 0;
};
CopyPath: PUBLIC PROC [path: Path] RETURNS [copy: Path] = {
copy _ CreatePath[path.maxVerts];
FOR i: NAT IN [0..path.len) DO
copy[i] _ path[i];
ENDLOOP;
copy.len _ path.len;
};
ClearPath: PUBLIC PROC [path: Path] = {
path.len _ 0;
};
AddPathPoint: PUBLIC PROC [path: Path, point: Point2d] RETURNS [pathPlusPoint: Path] = {
IF path.len = path.maxVerts THEN {
pathPlusPoint _ CreatePath[path.maxVerts+1];
pathPlusPoint _ PartPathGetsPartPath[path, 0, pathPlusPoint, 0, path.maxVerts];
pathPlusPoint[path.maxVerts] _ point;
pathPlusPoint.len _ path.maxVerts+1;
}
ELSE {
path[path.len] _ point;
path.len _ path.len + 1;
pathPlusPoint _ path;
};
}; -- end of AddPathPoint
InsertPathPoint: PUBLIC PROC [path: Path, point: Point2d] RETURNS [newPath: Path] = {
IF path.len = path.maxVerts THEN {
newPath _ CreatePath[path.maxVerts+1];
newPath[0] _ point;
newPath _ PartPathGetsPartPath[path, 0, newPath, 1, path.maxVerts];
newPath.len _ path.maxVerts+1;
}
ELSE {
ShiftUpPath[path, 0, 1];
path[0] _ point;
path.len _ path.len + 1;
newPath _ path;
};
}; -- end of InsertPathPoint
SplicePathPoint: PUBLIC PROC [path: Path, point: Point2d, index: INTEGER] RETURNS [newPath: Path]= {
IF index > path.len-1 THEN ERROR AttemptToSplicePastEndOfPath;
IF index = -1 THEN {newPath _ InsertPathPoint[path, point]; RETURN};
IF index = path.len-1 THEN {
newPath _ AddPathPoint[path, point]; RETURN};
IF path.len = path.maxVerts THEN {
newPath _ CreatePath[path.maxVerts+1];
newPath _ PartPathGetsPartPath[path, 0, newPath, 0, path.maxVerts];
ShiftUpPath[newPath, index+1, 1];
newPath[index+1] _ point;
}
ELSE {
ShiftUpPath[path, index+1, 1];
path[index+1] _ point;
newPath _ path;
};
}; -- end of SplicePathPoint

AttemptToSplicePastEndOfPath: PUBLIC ERROR = CODE;

DeletePathPoint: PUBLIC PROC [path: Path, index: NAT] RETURNS [newPath: Path] = {
IF index> path.len-1 THEN ERROR AttemptToDeleteNonExistentPoint;
IF path.len = 0 THEN ERROR PathEmpty;
IF index = path.len-1 THEN {path.len _ path.len -1; newPath _ path; RETURN};
ShiftDownPath[path, index+1, 1];
newPath _ path;
}; -- end of DeletePathPoint

AttemptToDeleteNonExistentPoint: PUBLIC ERROR = CODE;
PathEmpty: PUBLIC ERROR = CODE;

PutPathPoint: PUBLIC PROC [path: Path, index: NAT, point: Point2d] = {
IF index+1 > path.maxVerts THEN ERROR;
IF index+1 > path.len THEN path.len _ index+1;
path[index] _ point;
};
ConcatPath: PUBLIC PROC [path1, path2: Path] RETURNS [cat: Path] = {
cat _ CreatePath[path1.len+path2.len];
FOR i: NAT IN[0..path1.len) DO
cat[i] _ path1[i];
ENDLOOP;
FOR i: NAT IN[path1.len..path1.len+path2.len) DO
cat[i] _ path2[i-path1.len];
ENDLOOP;
};
SubPath: PUBLIC PROC [path: Path, lo, hi: NAT] RETURNS [subpath: Path] = {
IF hi < lo THEN ERROR;
subpath _ CreatePath[hi-lo+1];
FOR i: NAT IN[lo..hi] DO
subpath[i - lo] _ path[i];
ENDLOOP;
subpath.len _ hi-lo+1;
};
SubPathOfPoly: PUBLIC PROC [poly: Polygon, lo, hi: NAT] RETURNS [subpath: Path] = {
IF hi < lo THEN ERROR;
subpath _ CreatePath[hi-lo+1];
FOR i: NAT IN[lo..hi] DO
subpath[i - lo] _ poly[i];
ENDLOOP;
subpath.len _ hi-lo+1;
};
ShiftUpPath: PUBLIC PROC [path: Path, startAt: NAT, by: NAT] = {
IF path.len + by > path.maxVerts THEN ERROR;
FOR i: NAT DECREASING IN [startAt..path.len) DO
path[i+by]  _ path[i];
ENDLOOP;
path.len _ path.len + by;
}; -- end of ShiftUpPath
ShiftDownPath: PUBLIC PROC [path: Path, startAt: NAT, by: NAT] = {
IF startAt - by < 0 THEN ERROR ShiftingDataOffLeftEnd;
IF startAt > path.len-1 THEN ERROR ShiftingDownNonExistentElements;
FOR i: NAT IN [startAt..path.len) DO
path[i-by]  _ path[i];
ENDLOOP;
path.len _ path.len - by;
}; -- end of ShiftDownPath

ShiftingDataOffLeftEnd: PUBLIC ERROR = CODE;
ShiftingDownNonExistentElements: PUBLIC ERROR = CODE;

PolygonToPath: PUBLIC PROC [poly: Polygon] RETURNS [path: Path] = {
path _ NEW[PathObj[poly.len]];
FOR i: NAT IN[0..poly.len) DO
path[i] _ poly[i];
ENDLOOP;
path.len _ poly.len;
}; -- end of PolygonToPath
PartPathGetsPartPath: PUBLIC PROC [fromPath: Path, fromStart: NAT, toPath: Path, toStart: NAT, duration: NAT] RETURNS [newPath: Path] = {
j: NAT _ fromStart;
IF toStart+duration > toPath.maxVerts THEN {
newPath _ CreatePath[toStart+duration];
newPath _ PartPathGetsPartPath[toPath, 0, newPath, 0, toPath.len];
newPath.len _ toStart+duration}
ELSE {
newPath _ toPath;
IF toStart+duration > newPath.len THEN newPath.len _ toStart+duration};
IF fromStart+duration > fromPath.len THEN ERROR; -- don't use bad data
FOR i: NAT IN[toStart..toStart+duration-1] DO
newPath[i] _ fromPath[j];
j _ j + 1;
ENDLOOP;
};

PointPolyClass: TYPE = SVPolygon2d.PointPolyClass;-- {in, on, out};

CirclePointInPoly: PUBLIC PROC [point: Point2d, poly: Polygon] RETURNS [class: PointPolyClass] = {
v, c: Point2d;
thisTheta: REAL;
iPrime: NAT;
totalTheta: REAL _ 0;
v _ poly[0];
SVLines2d.FillTrigLineFromPoints[point, v, globalVLine];
FOR i: NAT IN[1..poly.len] DO
iPrime _ IF i = poly.len THEN 0 ELSE i;
c _ poly[iPrime];
SVLines2d.FillTrigLineFromPoints[point, c, globalCLine];
thisTheta _ SVAngle.ShortestDifference[globalCLine.theta, globalVLine.theta];
IF thisTheta = 180 OR thisTheta = -180 THEN RETURN[on];
totalTheta _ totalTheta + thisTheta;
v _ c;
SVLines2d.CopyTrigLine[from: globalCLine, to: globalVLine];
ENDLOOP;
IF -3.6 < totalTheta AND totalTheta < 3.6 THEN RETURN [out];
IF -363.6 < totalTheta AND totalTheta < -356.4 THEN RETURN [in]
ELSE RETURN [in]; -- **** should be an error
};
SquarePointInPoly: PUBLIC PROC [point: Point2d, poly: Polygon] RETURNS [class: PointPolyClass] = {
RETURN[in];
};

UpDown: TYPE = {up, on, down};
LeftRight: TYPE = {left, on, right};

YComparedToHorizLine: PRIVATE PROC [testY: REAL, yLine: REAL] RETURNS [UpDown] = {
IF testY > yLine THEN RETURN[up];
IF testY = yLine THEN RETURN[on];
RETURN[down];
};
XComparedToVertLine: PRIVATE PROC [testX: REAL, xLine: REAL] RETURNS [LeftRight] = {
IF testX > xLine THEN RETURN[right];
IF testX = xLine THEN RETURN[on];
RETURN[left];
};
LineSegHitsHorizLineLeftOfT: PRIVATE PROC [loPoint: Point2d, hiPoint: Point2d, t, lineY: REAL] RETURNS [LeftRight] = {
x: REAL;
IF loPoint[1] = hiPoint[1] THEN RETURN[XComparedToVertLine[loPoint[1], t]];
x _ loPoint[1] + (lineY - loPoint[2])*(hiPoint[1]-loPoint[1])/(hiPoint[2]-loPoint[2]);
RETURN[XComparedToVertLine[x, t]];
};

IsEven: PRIVATE PROC [n: NAT] RETURNS [BOOL] = {
RETURN [(n/2)*2 = n];
};
BoysePointInPoly: PUBLIC PROC [point: Point2d, poly: Polygon] RETURNS [class: PointPolyClass] = {
currentX, nextX: LeftRight;
currentY, nextY, resolveY: UpDown;
i: NAT;
crossings: NAT _ 0;
currentPoint: Point2d _ poly[0];
intersect: LeftRight;
currentX _ XComparedToVertLine[currentPoint[1], point[1]];
currentY _ YComparedToHorizLine[currentPoint[2], point[2]];
resolveY _ currentY;
FOR j: NAT IN [1..poly.len] DO
i _ IF j = poly.len THEN 0 ELSE j;
nextX _ XComparedToVertLine[poly[i][1], point[1]];
nextY _ YComparedToHorizLine[poly[i][2], point[2]];
SELECT TRUE FROM
currentX = left AND currentY = up =>
SELECT TRUE FROM
nextX = left AND nextY = on => {currentY _ on; resolveY _ up};  -- we stay left so no right crossings.
nextX = left => {resolveY _ currentY _ nextY};  -- we stay left so no right crossings
nextX = on AND nextY = up => {currentX _ on; resolveY _ up};
nextX = on AND nextY = on => {RETURN[on]};-- p is on a vertex of poly
nextX = on AND nextY = down => {currentX _ on; resolveY _ currentY _ down}; -- cross to the left
nextX = right AND nextY = up => {currentX _ right};-- cross vert line only
nextX = right AND nextY = on => {currentX _ right;  currentY _ on;  resolveY _ up}; -- we still consider ourselves up
nextX = right AND nextY = down => {-- one of the hard cases
currentX _ right; resolveY _ currentY _ down;
intersect _ LineSegHitsHorizLineLeftOfT[poly[i], currentPoint, point[1], point[2]];
IF intersect = right THEN crossings _ crossings + 1
ELSE IF intersect = on THEN RETURN[on];
};
ENDCASE => ERROR;
currentX = left AND currentY = on =>
SELECT TRUE FROM
nextX = left AND nextY = on => {currentY _ on}; -- resolve y unchanged
nextX = left => resolveY _ currentY _ nextY;
nextX = on AND nextY = up => {currentX _ on; resolveY _ currentY _ up};
nextX = on AND nextY = on => {RETURN[on]};-- p is on a vertex of poly
nextX = on AND nextY = down => {currentX _ on; resolveY _ currentY _ down};
nextX = right AND nextY = up => {currentX _ right; resolveY _ currentY _ up}; -- cross vert line only
nextX = right AND nextY = on => {RETURN[on]}; -- we crossed thru p.
nextX = right AND nextY = down => {currentX _ right;  resolveY _ currentY _ down};
ENDCASE => ERROR;
currentX = left AND currentY = down =>
SELECT TRUE FROM
nextX = left AND nextY = on => {currentY _ on};  -- we stay left. resolveY unchanged
nextX = left => {resolveY _ currentY _ nextY};  -- we stay left.  Take upDown of next
nextX = on AND nextY = up => {currentX _ on; resolveY _ currentY _ up};
nextX = on AND nextY = on => {RETURN[on]};-- p is on a vertex of poly
nextX = on AND nextY = down => {currentX _ on};
nextX = right AND nextY = up => {-- one of the hard cases
currentX _ right; resolveY _ currentY _ up;
intersect _ LineSegHitsHorizLineLeftOfT[currentPoint, poly[i], point[1], point[2]];
IF intersect = right THEN crossings _ crossings + 1
ELSE IF intersect = on THEN RETURN[on];
};
nextX = right AND nextY = on => {currentX _ right;  currentY _ on;  resolveY _ down}; -- we still consider ourselves down
nextX = right AND nextY = down => {currentX _ right};
ENDCASE => ERROR;
currentX = on AND currentY = up =>
SELECT TRUE FROM
nextX = left AND nextY = on => {currentX _ left; currentY _ on; resolveY _ up};  -- we go left.
nextX = left => {currentX _ left; resolveY _ currentY _ nextY};
nextX = on AND nextY = up => {};-- no change
nextX = on AND nextY = on => RETURN[on];-- p is on a vertex of poly
nextX = on AND nextY = down => RETURN[on]; -- we just crossed thru p
nextX = right AND nextY = up => {currentX _ right};
nextX = right AND nextY = on => {currentX _ right; currentY _ on;  resolveY _ up}; -- we still consider ourselves up
nextX = right AND nextY = down => {-- an automatic crossing!
currentX _ right; resolveY _ currentY _ down;
crossings _ crossings + 1};
ENDCASE => ERROR;
currentX = on AND currentY = on => RETURN[on]; -- should never occur.  Should be caught sooner.
currentX = on AND currentY = down =>
SELECT TRUE FROM
nextX = left AND nextY = on => {currentX _ left; currentY _ on; resolveY _ down};  -- we go left.
nextX = left => {currentX _ left; resolveY _ currentY _ nextY};
nextX = on AND nextY = up => RETURN[on]; -- we just crossed thru p
nextX = on AND nextY = on => RETURN[on];-- p is on a vertex of poly
nextX = on AND nextY = down => {}; -- no change
nextX = right AND nextY = up => {-- automatic crossing!
currentX _ right; resolveY _ currentY _ up;
crossings _ crossings + 1};
nextX = right AND nextY = on => {currentX _ right; currentY _ on;  resolveY _ down}; -- we still consider ourselves down
nextX = right AND nextY = down => {currentX _ right};
ENDCASE => ERROR;
currentX = right AND currentY = up =>
SELECT TRUE FROM
nextX = left AND nextY = up => {currentX _ left};
nextX = left AND nextY = on => {currentX _ left; currentY _ on; resolveY _ up};
nextX = left AND nextY = down => {-- one of the hard cases
currentX _ left; resolveY _ currentY _ down;
intersect _ LineSegHitsHorizLineLeftOfT[poly[i], currentPoint, point[1], point[2]];
IF intersect = right THEN crossings _ crossings + 1
ELSE IF intersect = on THEN RETURN[on];
};
nextX = on AND nextY = up => {currentX _ on};
nextX = on AND nextY = on => {RETURN[on]};-- p is on a vertex of poly
nextX = on AND nextY = down => {-- automatic crossing!
currentX _ on; resolveY _ currentY _  down;
crossings _ crossings + 1;
};
nextX = right AND nextY = up => {}; -- no change
nextX = right AND nextY = on => {currentY _ on;  resolveY _ up}; -- we still consider ourselves up
nextX = right AND nextY = down => { -- automatic crossing!
resolveY _ currentY _  down;
crossings _ crossings + 1;
};
ENDCASE => ERROR;
currentX = right AND currentY = on =>
SELECT TRUE FROM
nextX = left AND nextY = up => {-- this is a crossing if resolveY is down
IF resolveY = down THEN crossings _ crossings + 1;
currentX _ left; resolveY _ currentY _ up};
nextX = left AND nextY = on => {RETURN[on]};-- we pass thru p
nextX = left AND nextY = down => {-- this is a crossing if resolveY is up
IF resolveY = up THEN crossings _ crossings + 1;
currentX _ left; resolveY _ currentY _ down};
nextX = on AND nextY = up => {-- this is a crossing if resolveY is down
IF resolveY = down THEN crossings _ crossings + 1;
currentX _ on; resolveY _ currentY _ up};
nextX = on AND nextY = on => {RETURN[on]};-- p is on a vertex of poly
nextX = on AND nextY = down => {-- this is a crossing if resolveY is up
IF resolveY = up THEN crossings _ crossings + 1;
currentX _ on; resolveY _ currentY _ down};
nextX = right AND nextY = up => {-- this is a crossing if resolveY is down
IF resolveY = down THEN crossings _ crossings + 1;
resolveY _ currentY _ up};
nextX = right AND nextY = on => {}; -- no change
nextX = right AND nextY = down => {-- this is a crossing if resolveY is up
IF resolveY = up THEN crossings _ crossings + 1;
resolveY _ currentY _ down};
ENDCASE => ERROR;
currentX = right AND currentY = down =>
SELECT TRUE FROM
nextX = left AND nextY = up => {-- one of the hard cases
currentX _ left; resolveY _ currentY _ up;
intersect _ LineSegHitsHorizLineLeftOfT[currentPoint, poly[i], point[1], point[2]];
IF intersect = right THEN crossings _ crossings + 1
ELSE IF intersect = on THEN RETURN[on];
};
nextX = left AND nextY = on => {currentX _ left; currentY _ on; resolveY _ down};
nextX = left AND nextY = down => {currentX _ left};
nextX = on AND nextY = up => {-- automatic crossing!
currentX _ on; resolveY _ currentY _  up;
crossings _ crossings + 1};
nextX = on AND nextY = on => {RETURN[on]};-- p is on a vertex of poly
nextX = on AND nextY = down => {currentX _ on};
nextX = right AND nextY = up => {-- automatic crossing!
resolveY _ currentY _  up;
crossings _ crossings + 1};
nextX = right AND nextY = on => {currentY _ on;  resolveY _ down}; -- we still consider ourselves down
nextX = right AND nextY = down => {};-- no change
ENDCASE => ERROR;
ENDCASE => ERROR;
currentPoint _ poly[i];
ENDLOOP;
IF IsEven[crossings] THEN RETURN[out] ELSE RETURN[in];
};

PointToVector: PRIVATE PROC [point: Point2d] RETURNS [vector: Vector3d] = {
vector[1] _ point[1];
vector[2] _ point[2];
vector[3] _ 0;
};

SignedArea: PUBLIC PROC [poly: Polygon] RETURNS [area: REAL] = {
D1, D2: Vector3d;
iPlusOne: NAT;
areaVector: Vector3d;
partialArea: REAL;
area _ 0;
FOR i: NAT IN [0..poly.len-1] DO
iPlusOne _ IF i = poly.len-1 THEN 0 ELSE i + 1;
D1 _ PointToVector[poly[i]];
D2 _ PointToVector[poly[iPlusOne]];
areaVector _ SVVector3d.CrossProduct[D2, D1];-- will only have a z component
partialArea _ areaVector[3];
area _ area + partialArea;
ENDLOOP;
}; -- end of signed area
ClockwisePerimeterAroundUnitSquare: PUBLIC PROC [from, to: Point2d] RETURNS [perim: REAL] = {
Side: TYPE = NAT;
side: Side;
fromSide, toSide: Side;
perimPart: REAL;
SELECT TRUE FROM
from[1] = -0.5 => fromSide _ 0;
from[2] = 0.5  => fromSide _ 1;
from[1] = 0.5 => fromSide _ 2;
from[2] = -0.5 => fromSide _ 3;
ENDCASE => ERROR;
SELECT TRUE FROM
to[1] = -0.5 => toSide _ 0;
to[2] = 0.5  => toSide _ 1;
to[1] = 0.5 => toSide _ 2;
to[2] = -0.5 => toSide _ 3;
ENDCASE => ERROR;
side _ fromSide;
IF toSide = fromSide THEN {
IF IsClockwisePtoQAlongSquareSide[from, to, toSide] THEN 
perim _ DistanceClockwiseAlongSquareSide[from, to, toSide]
ELSE perim _ -DistanceClockwiseAlongSquareSide[to, from, toSide];
RETURN};
perim _ DistanceClockwiseAlongSquareSide[from, unitSquare[fromSide+1], fromSide];
FOR i: NAT IN[1..3] DO
side _ fromSide + i;
IF side > 3 THEN side _ side -4;
IF side = toSide THEN GOTO ToFound;
perim _ perim + 1;
REPEAT
ToFound => {-- find distance from last vertex to this point
perimPart _ DistanceClockwiseAlongSquareSide[unitSquare[side], to, side];
perim _ perim + perimPart;
};
ENDLOOP;
IF perim > 2 THEN perim _ perim - 4;
}; -- end of ClockwisePerimeterAroundUnitSquare
IsClockwisePtoQAlongSquareSide: PROC [p, q: Point2d, side: NAT] RETURNS [BOOL] = {
SELECT side FROM
0 => -- left side.  p to q is clockwise if q above p;
RETURN[q[2] > p[2]];
1 => -- top side.  p to q is clockwise if q right of p;
RETURN[q[1] > p[1]];
2 => -- right side.  p to q is clockwise if q below p;
RETURN[q[2] < p[2]];
3 => -- bottom side.  p to q is clockwise if q left of p;
RETURN[q[1] < p[1]];
ENDCASE => ERROR;
}; -- end of IsClockwisePtoQAlongSquareSide

DistanceClockwiseAlongSquareSide: PROC [p, q: Point2d, side: NAT] RETURNS [REAL] = {
SELECT side FROM
0 => -- left side.  Distance p to q is difference in y's;
RETURN[q[2] - p[2]];
1 => -- top side.  Distance p to q is difference in x's;
RETURN[q[1] - p[1]];
2 => -- right side.  Distance p to q is difference in y's;
RETURN[p[2] - q[2]];
3 => -- bottom side.  Distance p to q is difference in x's;
RETURN[p[1] - q[1]];
ENDCASE => ERROR;
}; -- end of DistanceClockwiseAlongSquareSide

Init: PROC = {
unitSquare _ CreatePoly[4];
unitSquare _ AddPolyPoint[unitSquare, [-.5, -.5]];
unitSquare _ AddPolyPoint[unitSquare, [-.5, .5]];
unitSquare _ AddPolyPoint[unitSquare, [.5, .5]];
unitSquare _ AddPolyPoint[unitSquare, [.5, -.5]];
root3 _ RealFns.SqRt[3];
root3Over3 _ root3/3.0;
twoRoot3Over3 _ 2*root3Over3;
globalVLine _ SVLines2d.CreateEmptyTrigLine[];
globalCLine _ SVLines2d.CreateEmptyTrigLine[];
}; -- end of Init


Init[];

END.


File: SVPolygon2dImpl.mesa
Last edited by Bier on June 1, 1984 4:51:36 pm PDT
Author: Eric Bier on February 18, 1987 6:21:46 pm PST
Contents: Routines for creating, manipulating, and playing with polygons.  Part of the SolidViews package.
GLOBAL VARIABLES
POLYGONS
Makes a hexagon with two of its sides parallel to the x axis.
toPoly in range [toStart..toStart+duration-1] gets fromPath in range [fromStart..fromStart+duration-1];
toPoly in range [toStart..toStart+duration-1] gets fromPoly in range [fromStart..fromStart+duration-1];
TRIG POLYGONS (a richer [and slower] structure)
IF index is -1, splice onto beginning, if index is path.len-1, splice onto end.
startAt is the index of the lowest element which should be shifted up
startAt  the index of the highest element which should be shifted down.  Must be greater than zero.
toPath in range [toStart..toStart+duration-1] gets fromPath in range [fromStart..fromStart+duration-1];
Take point.  Take a vertex V and its clockwise neighbor C from poly.  Find the angle from V to C.  Add up these angles for all such pairs V and C in poly.  If the total is -360 degress, the point is in the poly.  If the total is zero, the point is outside.  Say +or- 3.6 degrees is close enough to zero and -360 +or- 3.6 is good enough for 360.  Otherwise, error.
to find the angles, I could take the cross product of the vectors from point to V and to C, normalize, and the take ArcSin, but Cedar doesn't have ArcSin so for now, I will use the ArcTan machinery available from SV2d.CreateTrigLineSeg to find the angle of each vector with respect to the x-axis and then take the difference of these two angles.
We wish to know if the (infinite) line passing through loPoint and hiPoint hits the horizontal line y = lineY to the left, right, or smack dab on the point (t, lineY).
Let (x0,y0) = loPoint and (x1,y1) = hiPoint then the line equation can be expressed as:
y - y0 = (y1-y0)/(x1-x0)(x - x0) unless x1 = x0.
IF x1 = x0 then we compare x0 to t and return result.
Otherwise, we solve for x where y = lineY:
lineY - y0 = [(y1-y0)/(x1-x0)](x - x0)
(lineY - y0)(x1-x0) = (y1-y0)(x - x0)
(lineY - y0)(x1-x0)/(y1-y0) = x - x0
x0 + (lineY - y0)(x1-x0)/(y1-y0) = x.
So, x > t => x0 + (lineY - y0)(x1-x0)/(y1-y0) > t.
taking 1 mult and 1 div which =approx.  3 mults.
Consider a horiz line and a vert line passing through p.  How many edges of poly cross the horizontal line to the right of the vertical line?  Call this number "hits".  If hits is even then p is outside of poly.  If hits is odd, the p is inside poly.  If we find an edge of poly which passes through p, then p is on poly;
Call the vertical line vLine, (x = p[1]) and the horizontal line hLine (y = p[2]).
Method:  There are nine starting cases.  The first vert of poly may be up, on, or down of hLine.  It may be left, on or right of vLine.  Similarly, there are nine cases for the other endpoint.  Most cases make it easy to determine whether this edge of the poly cross hLine to the right of vLine or not.  The hard cases are:
1) up-left to down-right
2) down-left to up-right
For these cases, we use LineSegHitsHorizLineLeftOfT where t is p[1]
Hopefully, there will be some way of coding this without mentioning all 81 cases explicitly.
Chose an arbitrary point P in the plane.  Chose two consecutive vertices (v1, v2) on the polygon.  Call the vector P to v1, D1.  Call the distance P to v2, D2.  The area of the triangle made by these three points is |D1|*|D2|*sin(angle between D1 and D2)/2.
This is just D2 x D1 where "x" is the cross product operator. 
This will be positive if D2 is clockwise of D1.
If P is in the interior of poly for poly convex, and if we add up all of the areas we get by taking vertices pairwise in clockwise order, area will be positive.  If we take them in counter-clockwise order, area will be negative.  This turns out to be true whatever the position of P and even if poly is concave part of the time.
let P be the origin so the vector is the same as its endpoint
Find out which sides from and to are on.  Number the sides as follows:  0 = left, 1 = top, 2 = right, 3 = bottom.
unitsquare is a global polygon. unitsquare[0] = lower left; [1] = upper left; [2] = upper right; [3] = lower right;
toSide # fromSide so add on the clockwise distance to the next vertex and keep going.
Assumes that q is clockwise from p;
ΚΗ–

"cedar" style˜Iheadš
œ™Iprocš
œ2™2Lš
œ5™5Lš
œj™jL˜š
Οk	˜	Lš
œ;˜;—L˜Lšœ
˜Lšœ(˜/Lšœ˜Lš
˜˜Lšœ	œ˜Lšœœœ	˜Lšœ	œ˜Lšœ	œœ˜Lšœœ˜#Lšœ
œœ
˜!Lšœ
œ˜%Lšœ
œ˜%Lšœ
œœ˜'Lšœœ˜+Lšœ
œ˜%L˜Lš
œ™L˜Lš
œ˜Lšœ"œ8˜^Lš
œ#˜#L˜Lš
œ™L˜—š
Οn
œœ
œœœ˜>Lšœœ˜Lš
œ
˜
Lš
œ˜—šžœœ
œœ˜ALš
œ!˜!šœœ
œ˜Lš
œ˜—Lšœ
˜Lš
œ˜Lš
œ˜—L˜š
ž
œœ
œœœ˜?Lš
œ=™=Lš
œ˜Lš
œ9˜9Lš
œ6˜6Lš
œ9˜9Lš
œ9˜9Lš
œ6˜6Lš
œ9˜9Lš
œ˜—šž	œœ
œ˜*Lš
œ
˜
Lš
œ˜—šžœœ
œ!œ˜^šœœ˜"Lš
œ,˜,Lš
œO˜OLš
œ%˜%Lš
œ$˜$Lš
œ
˜
—šœ˜Lš
œ˜Lš
œ˜Lš
œ˜Lš
œ˜—Lš
œ˜—š
žœœ
œœœ˜dšœœ˜!Lš
œ˜Lš
œC˜CLš
œ˜Lš
œ˜Lš
œ
˜
—šœ˜Lšœœ˜.Lš
œ˜Lš
œ˜—Lš
œ˜—š
žœœ
œœœ˜?Lšœ˜Lš
œ˜—šž
œœ
œœ˜CLšœœ˜!šœœ
œ˜Lš
œ ˜ —Lšœ
˜Lš
œ˜—šžœœ
œ˜2Lš
œ˜šœœ
œ˜Lš
œ#˜#Lš
œ!˜!Lš
œ˜—Lšœ
˜Lš
œ˜—šžœœ
œœœœœ˜
Lš
œg™gLšœœ
˜šœ$œ˜,Lš
œ'˜'Lš
œB˜BLš
œ˜—šœ˜Lš
œ˜Lšœ œ!˜G—Lšœ#œ
œ˜Fšœœ
œ˜-Lš
œ˜Lš
œ
˜
—Lšœ
˜Lš
œ˜—šžœœ
œ œœœœ˜’
Lš
œg™gLšœœ
˜šœ$œ˜,Lš
œ'˜'Lš
œB˜BLš
œ˜—šœ˜Lš
œ˜Lšœ œ!˜G—Lšœ#œ
œ˜Fšœœ
œ˜-Lš
œ˜Lš
œ
˜
—Lšœ
˜Lš
œ˜—šž
œœ
œœ˜CLšœœ˜!šœœ
œ˜Lš
œ˜—Lšœ
˜Lš
œ˜Lš
œ˜—L˜Lš
œ/™/L˜šžœœ
œœ˜ULš
œ˜Lš
œ˜Lšœœ˜)Lš
œ˜šœœ
œ˜Lš
œ˜Lš
œ;˜;Lš
œ˜Lš
œ˜—Lšœ
˜Lš
œ˜Lš
œ;˜;Lš
œ˜Lš
œ˜—š
ž
œœ
œœœ˜;Lšœœ˜Lš
œ
˜
Lš
œ˜—šžœœ
œœ˜;Lš
œ!˜!šœœ
œ˜Lš
œ˜—Lšœ
˜Lš
œ˜Lš
œ˜—šž	œœ
œ˜'Lš
œ
˜
Lš
œ˜—šžœœ
œœ˜Xšœœ˜"Lš
œ,˜,Lš
œO˜OLš
œ%˜%Lš
œ$˜$Lš
œ
˜
—šœ˜Lš
œ˜Lš
œ˜Lš
œ˜Lš
œ˜—Lš
œ˜—šžœœ
œœ˜Ušœœ˜"Lš
œ&˜&Lš
œ˜Lš
œC˜CLš
œ˜Lš
œ
˜
—šœ˜Lš
œ˜Lš
œ˜Lš
œ˜Lš
œ˜Lš
œ˜—Lš
œ˜—š
žœœ
œ%œœ˜dLš
œO™OLšœœ
œ˜>Lšœœ*œ˜Dšœœ˜Lšœ%œ˜-—šœœ˜"Lš
œ&˜&Lš
œC˜CLš
œ!˜!Lš
œ˜Lš
œ
˜
—šœ˜Lš
œ˜Lš
œ˜Lš
œ˜Lš
œ˜—Lš
œ˜—L˜Lšœœ
œœ
˜2L˜š
žœœ
œœœ˜QLšœœ
œ!˜@Lšœœ
œ˜%Lšœœ*œ˜LLš
œ ˜ Lš
œ˜Lš
œ˜—L˜Lšœ!œ
œœ
˜5Lšœœ
œœ
˜L˜šžœœ
œœ˜FLšœœ
œ
˜&Lšœœ˜.Lš
œ˜Lš
œ˜—šž
œœ
œœ˜DLš
œ&˜&šœœ
œ˜Lš
œ˜—Lšœ
˜šœœ
œ!˜0Lš
œ˜—Lšœ
˜Lš
œ˜—š
žœœ
œœœ˜JLšœ	œ
œ
˜Lš
œ˜šœœ
œ	˜Lš
œ˜—Lšœ
˜Lš
œ˜Lš
œ˜—š
ž
œœ
œœœ˜SLšœ	œ
œ
˜Lš
œ˜šœœ
œ	˜Lš
œ˜—Lšœ
˜Lš
œ˜Lš
œ˜—š
žœœ
œœœ˜@Lš
œE™ELšœœ
œ
˜,š	œœ

œ
œ˜/Lš
œ˜—Lšœ
˜Lš
œ˜Lš
œ˜—š
ž
œœ
œœœ˜BLš
œc™cLšœœ
œ˜6Lšœœ
œ!˜Cšœœ
œ˜$Lš
œ˜—Lšœ
˜Lš
œ˜Lš
œ˜—L˜Lšœœ
œœ
˜,Lšœ!œ
œœ
˜5L˜šž
œœ
œœ˜CLšœœ˜šœœ
œ˜Lš
œ˜—Lšœ
˜Lš
œ˜Lš
œ˜—šžœœ
œœœœœ˜‰
Lš
œg™gLšœœ
˜šœ$œ˜,Lš
œ'˜'Lš
œB˜BLš
œ˜—šœ˜Lš
œ˜Lšœ œ!˜G—Lšœ#œ
œ˜Fšœœ
œ˜-Lš
œ˜Lš
œ
˜
—Lšœ
˜Lš
œ˜—L˜Lšœœ/˜CL˜šžœœ
œ!œ˜bLš
œλ™λLš
œΩ™ΩLš
œ˜Lšœœ
˜Lšœœ
˜Lšœœ˜Lš
œ˜Lš
œ8˜8šœœ
œ˜Lšœ	œœœ˜'Lš
œ˜Lš
œ8˜8Lš
œM˜MLšœœœ
œ˜7Lš
œ$˜$Lš
œ˜Lš
œ;˜;—Lšœ
˜Lšœœœ
œ˜<Lšœœœ
œ˜?Lšœ
œ!˜,Lš
œ˜—šžœœ
œ!œ˜bLšœ˜Lš
œ˜—L˜Lšœœ˜Lšœœ˜$L˜šžœœ
œ	œ	œœ
˜RLšœœ
œ˜!Lšœœ
œ˜!Lšœ˜
Lš
œ˜—šžœœ
œ	œ	œœ˜TLšœœ
œ˜$Lšœœ
œ˜!Lšœ˜
Lš
œ˜—š
žœœ
œ0œœ˜vLš
Ϥ
™§
Lš
œW™WLš
œ0™0Lš
œ5™5Lš
œ*™*Lš
œ&™&Lš
œ%™%Lš
œ$™$Lš
œ%™%Lš
œ2™2Lš
œ0™0Lšœœ
˜Lšœœ
œ%˜KLš
œV˜VLšœ˜"Lš
œ˜—L˜šžœœ
œœœœ˜0Lšœ˜Lš
œ˜—šžœœ
œ!œ˜aLš
œΑ™ΑLš
œR™RLš
œΓ™ΓLš
œ™Lš
œ™Lš
œC™CLš
œ\™\Lš
œ˜Lš
œ"˜"Lšœœ
˜Lšœœ˜Lš
œ ˜ Lš
œ˜Lš
œ:˜:Lš
œ;˜;Lš
œ˜šœœ
œ˜Lšœœœœ˜"Lš
œ2˜2Lš
œ3˜3Lšœ
œ
˜šœœ˜$Lšœ
œ
˜Lšœ
œV˜fLš
œU˜ULšœœ.˜<Lšœœœ!˜ELšœœR˜`Lšœœ9˜JLšœœd˜uLšœœ*˜;Lš
œ-˜-Lš
œS˜SLšœœ˜3Lšœ
œœ
œ˜'Lš
œ˜Lšœœ
˜—šœœ˜$Lšœ
œ
˜Lšœ
œ6˜FLš
œ,˜,Lšœœ9˜GLšœœœ!˜ELšœœ=˜KLšœœT˜eLšœœœ˜CLšœœA˜RLšœœ
˜—šœœ˜&Lšœ
œ
˜Lšœ
œD˜TLš
œU˜ULšœœ9˜GLšœœœ!˜ELšœœ!˜/šœœ(˜9Lš
œ+˜+Lš
œS˜SLšœœ˜3Lšœ
œœ
œ˜'Lš
œ˜—Lšœœh˜yLšœœ$˜5Lšœœ
˜—šœœ˜"Lšœ
œ
˜Lšœ
œO˜_Lš
œ?˜?Lšœœ˜,Lšœœœ ˜CLšœœœ˜DLšœœ"˜3Lšœœc˜tLšœœ+˜<Lš
œ-˜-Lš
œ˜Lšœœ
˜—Lšœœœ6˜_šœœ˜$Lšœ
œ
˜Lšœ
œQ˜aLš
œ?˜?Lšœœœ˜BLšœœœ ˜CLšœœ!˜/Lšœœ&˜7Lš
œ+˜+Lš
œ˜Lšœœg˜xLšœœ$˜5Lšœœ
˜—šœœ˜%Lšœ
œ
˜Lšœ
œ!˜1Lšœ
œ?˜Ošœ
œ*˜:Lš
œ,˜,Lš
œS˜SLšœœ˜3Lšœ
œœ
œ˜'Lš
œ˜—Lšœœ˜-Lšœœœ!˜Ešœœ(˜6Lš
œ+˜+Lš
œ˜Lš
œ˜—Lšœœ˜0LšœœQ˜bLšœœ)˜:Lš
œ˜Lš
œ˜Lš
œ˜Lšœœ
˜—šœœ˜%Lšœ
œ
˜šœ
œ9˜ILšœœ˜2Lš
œ+˜+—Lšœ
œœ˜=šœ
œ9˜ILšœœ˜0Lš
œ-˜-—šœœ9˜GLšœœ˜2Lš
œ)˜)—Lšœœœ!˜Ešœœ9˜GLšœœ˜0Lš
œ+˜+—šœœ9˜JLšœœ˜2Lš
œ˜—Lšœœ˜0šœœ9˜JLšœœ˜0Lš
œ˜—Lšœœ
˜—šœœ˜'Lšœ
œ
˜šœ
œ(˜8Lš
œ*˜*Lš
œS˜SLšœœ˜3Lšœ
œœ
œ˜'Lš
œ˜—Lšœ
œA˜QLšœ
œ#˜3šœœ&˜4Lš
œ)˜)Lš
œ˜—Lšœœœ!˜ELšœœ!˜/šœœ&˜7Lš
œ˜Lš
œ˜—LšœœU˜fLšœœ ˜1Lšœœ
˜—Lšœœ
˜Lš
œ˜—Lšœ
˜Lš
œœ
œœ
œ˜6Lš
œ˜—L˜šž
œœ
œœ˜KLš
œ˜Lš
œ˜Lš
œ˜Lš
œ˜—L˜š
ž
œœ
œœœ˜@Lš
œ™Lš
œ>™>Lš
œ/™/Lš
œΘ™ΘLš
œ˜Lšœ
œ
˜Lš
œ˜Lšœ
œ
˜Lš
œ	˜	Lš
œ=™=Lšœœ
œ˜ Lšœœœœ˜/Lš
œ˜Lš
œ#˜#Lš
œL˜LLš
œ˜Lš
œ˜Lšœ
˜Lš
œ˜—š
ž"œœ
œœ	œ˜]Lš
œq™qLš
œs™sLšœœœ
˜Lš
œ˜Lš
œ˜Lšœœ
˜Lšœ
œ
˜Lš
œ˜Lš
œ˜Lš
œ˜Lš
œ˜Lšœœ
˜Lšœ
œ
˜Lš
œ˜Lš
œ˜Lš
œ˜Lš
œ˜Lšœœ
˜Lš
œ˜šœœ˜Lšœ2œ
˜9Lš
œ:˜:Lšœ=˜ALšœ˜—Lš
œU™ULš
œQ˜Qšœœ
œ˜Lš
œ˜Lšœ
œ˜ Lšœœ
œ	˜#Lš
œ˜—š
˜Lš
œ;˜;Lš
œI˜ILš
œ˜Lš
œ˜—Lšœ
˜Lšœœ˜$Lš
œ/˜/—š
žœœœœœ˜RLšœ˜š
œ5˜5Lšœ˜—š
œ7˜7Lšœ˜—š
œ6˜6Lšœ˜Lš
œ9˜9—Lšœ˜Lšœœ
˜Lš
œ+˜+—L˜š
ž œœœœœ˜TLš
œ#™#Lšœ˜š
œ9˜9Lšœ˜—š
œ8˜8Lšœ˜—š
œ:˜:Lšœ˜—š
œ;˜;Lšœ˜—Lšœœ
˜Lš
œ-˜-—L˜šžœœ˜š
œ˜Lš
œ2˜2Lš
œ1˜1Lš
œ0˜0Lš
œ1˜1Lš
œ˜Lš
œ˜Lš
œ˜Lš
œ.˜.Lš
œ.˜.Lš
œ˜—L˜—L˜Lš
œ˜L˜Lšœ
˜L˜L˜—…—W‚Ϋ