G3dPatchFromPolyProcs.mesa
Copyright © 1986 by Xerox Corporation. All rights reserved.
Last Edited by: Crow, September 30, 1989 6:24:26 pm PDT
Glassner, July 5, 1989 6:33:58 pm PDT
Bloomenthal, June 19, 1989 12:13:08 pm PDT
DIRECTORY Atom, Basics, G2dVector, G3dBasic, G3dMappedAndSolidTexture, G3dMatrix, G3dRender, G3dScanConvert, G3dClipXfmShade, G3dShape, G3dSortandDisplay, G3dSpline, G3dVector, RealFns;
G3dPatchFromPolyProcs: CEDAR MONITOR
IMPORTS Atom, Basics, G2dVector, G3dMappedAndSolidTexture, G3dMatrix, G3dRender, G3dClipXfmShade, G3dSortandDisplay, G3dSpline, G3dVector, RealFns
= BEGIN
Types
Ray:      TYPE ~ G3dBasic.Ray;
NatSequence:   TYPE ~ G3dRender.NatSequence;
NatSequenceRep:  TYPE ~ G3dBasic.NatSequenceRep;
Pair:     TYPE ~ G3dRender.Pair;
Triple:    TYPE ~ G3dRender.Triple;
TripleSequence:  TYPE ~ G3dRender.TripleSequence;
TripleSequenceRep: TYPE ~ G3dBasic.TripleSequenceRep;
RGB:     TYPE ~ G3dRender.RGB;
RealSequence:  TYPE ~ G3dRender.RealSequence;
RealSequenceRep: TYPE ~ G3dBasic.RealSequenceRep;
Context:    TYPE ~ G3dRender.Context;
SixSides:    TYPE ~ G3dRender.SixSides;
Matrix:    TYPE ~ G3dRender.Matrix;
Shape:     TYPE ~ G3dRender.Shape;
Vertex:     TYPE ~ G3dShape.Vertex;
Patch:     TYPE ~ G3dRender.Patch;
PatchSequence:  TYPE ~ G3dRender.PatchSequence;
PatchSequenceRep: TYPE ~ G3dRender.PatchSequenceRep;
PatchProc:   TYPE ~ G3dRender.PatchProc;
CtlPoint:    TYPE ~ G3dRender.CtlPoint;
CtlPtInfo:    TYPE ~ G3dRender.CtlPtInfo;
RECORD[coord: CtlPoint, shade: Shading, vtxPtr: NAT, data: REF];
CtlPtInfoSequence: TYPE ~ G3dRender.CtlPtInfoSequence;
CtlPtInfoProc:  TYPE ~ G3dRender.CtlPtInfoProc;
Shading:    TYPE ~ G3dRender.Shading;
RenderStyle:   TYPE ~ G3dRender.RenderStyle;
RenderData:   TYPE ~ G3dRender.RenderData;
ShapeClass:   TYPE ~ G3dRender.ShapeClass;
ShapeProc:   TYPE ~ G3dRender.ShapeProc;
ClipState:    TYPE ~ G3dRender.ClipState;
OutCode:    TYPE ~ G3dRender.OutCode;
AllOut:    OutCode ~ G3dRender.AllOut;
NoneOut:    OutCode ~ G3dRender.NoneOut;
FacingDir:   TYPE ~ G3dRender.FacingDir;
LORA:     TYPE = LIST OF REF ANY;
MidPtProc: TYPE ~ PROC[v0, v1: REF CtlPtInfo] RETURNS[REF CtlPtInfo];
nullTriple: Triple ~ [0.0, 0.0, 0.0];
Corner: TYPE ~ RECORD [ inVtx, outVtx: NAT ← 0,
        inDir, outDir, normal, interiorKnot: Triple ← nullTriple,
        concave: BOOLEANFALSE ];
Represents a corner of a polygon,
- inVtx is the vertex at the other end of the incoming edge (previous vertex in order),
- outVtx is other vertex on outgoing edge,
- inDir, outDir are the corresponding outward direction vectors,
- normal is the carefully determined normal vector
- concave = TRUE indicates corner is concave vertex needing reversed cross product.
CornerSeq: TYPE ~ RECORD [ length: NAT ← 0, s: SEQUENCE maxLength: NAT OF Corner ];
CornerSeqSeq: TYPE ~ RECORD [ length: NAT ← 0,
          s: SEQUENCE maxLength: NAT OF REF CornerSeq ];
TangentSet: TYPE ~ RECORD [t0, et0, t1, et1: Triple ← nullTriple];
Represents the tangents at the endpoints of an edge (in object space and eyespace), the edge is ordered in the vertex order of the polygon (t0 at leading endpoint, t1 at trailing endpoint)
TangentSeq: TYPE ~ RECORD [length: NAT ← 0, s: SEQUENCE maxLength: NAT OF TangentSet];
TangentTriple: TYPE ~ ARRAY [0..3) OF TangentSet;
TangentQuad: TYPE ~ ARRAY [0..4) OF TangentSet;
TangentSeqSeq: TYPE ~ RECORD [
       length: NAT ← 0, s: SEQUENCE maxLength: NAT OF REF TangentSeq];
       
Triangle: TYPE ~ RECORD [ v: ARRAY[0..3) OF CtlPtInfo, t: ARRAY[0..3) OF TangentSet ];
NatSequenceSequence: TYPE ~ RECORD [
       length: NAT ← 0, s: SEQUENCE maxLength: NAT OF NatSequence ];
BoolSequence: TYPE ~ RECORD [length: NAT ← 0, s: SEQUENCE maxLength: NAT OF BOOLEAN];
Renamed Procedures
Add: PROC[v1, v2: Triple] RETURNS[Triple] ~ G3dVector.Add;
Mul: PROC[v: Triple, s: REAL] RETURNS[Triple] ~ G3dVector.Mul;
Cross: PROC[v1, v2: Triple] RETURNS[Triple] ~ G3dVector.Cross;
Div: PROC[v: Triple, s: REAL] RETURNS[Triple] ~ G3dVector.Div;
Dot: PROC[v1, v2: Triple] RETURNS[REAL] ~ G3dVector.Dot;
Length: PROC[v: Triple] RETURNS[REAL] ~ G3dVector.Length;
Negate: PROC[v: Triple] RETURNS[Triple] ~ G3dVector.Negate;
Nmlize: PROC[v: Triple] RETURNS[Triple] ~ G3dVector.Unit;
Sub3: PROC[v1, v2: Triple] RETURNS[Triple] ~ G3dVector.Sub;
ReleasePatch: PROC [p: REF Patch] ~ G3dClipXfmShade.ReleasePatch;
GetPatch: PROC [size: NAT] RETURNS [REF Patch] ~ G3dClipXfmShade.GetPatch;
GetProp: PROC [propList: Atom.PropList, prop: REF ANY] RETURNS [REF ANY] ~
                     Atom.GetPropFromList;
PutProp: PROC [propList: Atom.PropList, prop: REF ANY, val: REF ANY]
   RETURNS
[Atom.PropList] ~ Atom.PutPropOnList;
Global Variables
maxDeviation: REAL ← .5;        -- subdivision tolerance with antialiasing
maxJaggyDeviation: REAL ← 2.0;      -- jaggy subdivision tolerance (pixels)
closenessFactor: REAL ← 10.0;   -- times maxDeviation gets clipping cutoff for straightness      
recurseLimit: NAT ← 14;       -- safety valve on recursion
minCosToAlign: REAL ← -0.866;   -- Align edges if they meet at > 150 degrees
minCosToAlignOpen: REAL ← 0.707;  -- Align open edges if they meet at > 45 degrees
stopIfStraight: BOOLEANTRUE;     -- set false to defeat termination algorithm
unitNormals: BOOLEANFALSE; -- treat all polygons equally when making vertex normal
showLines: BOOLEANFALSE;        -- debug and pedagogical aid
stopAtPoly: NAT ← 0;      -- debug aid, allows breakpoint on given polygon #
useManhattan: BOOLEANTRUE;  -- switches termination procedures in EdgeStraight
Utility Procedures
PutPropSafely: PROC[propList: Atom.PropList, prop, val: REF ANY] RETURNS[Atom.PropList] ~{
Put property on new property list, to avoid clobbering other lists inherited from same place
newProps: Atom.PropList ← NIL;
FOR list: Atom.PropList ← propList, list.rest UNTIL list = NIL DO-- new proplist
element: Atom.DottedPair ← NEW[Atom.DottedPairNode ← list.first^];
newProps ← CONS[element, newProps];
ENDLOOP;
RETURN[ PutProp[ newProps, prop, val ] ];
};
Sqr: PROCEDURE [number: REAL] RETURNS [REAL] ~ INLINE { RETURN[number * number]; };
Sub: PROC[f, g: REAL] RETURNS[REAL] ~ {
 Difference f - g, returns zero if less than 3 decimal places
IF RealFns.AlmostEqual[f, g, -10] THEN RETURN [0.0] ELSE RETURN[f-g];
};
Sub1: PROC[f, g: REAL] RETURNS[REAL] ~ {  
Screen space f - g, returns zero if less than noticeable
result: REAL ← f-g;
IF ABS[result] < G3dScanConvert.justNoticeable THEN RETURN [0.0] ELSE RETURN[result];
};
DiffPosnsVtx: PROC[vtx1, vtx2: Vertex, space: ATOMNIL] RETURNS[Triple] ~ {
Subtracts vtx2 from vtx1
RETURN[ [
Sub[vtx1.point.x, vtx2.point.x],
Sub[vtx1.point.y, vtx2.point.y],
Sub[vtx1.point.z, vtx2.point.z]
] ];
};
DiffPosnsCtlPt: PROC[vtx1, vtx2: CtlPoint, space: ATOMNIL] RETURNS[Triple] ~ {
Subtracts vtx2 from vtx1 in selected space
SELECT space FROM
$Eye  => RETURN[[Sub[vtx1.ex, vtx2.ex], Sub[vtx1.ey, vtx2.ey], Sub[vtx1.ez, vtx2.ez]]];
$Screen=> RETURN[[Sub1[vtx1.sx, vtx2.sx], Sub1[vtx1.sy, vtx2.sy], Sub1[vtx1.sz, vtx2.sz]]];
ENDCASE =>      -- object space
RETURN[ [ Sub[vtx1.x, vtx2.x], Sub[vtx1.y, vtx2.y], Sub[vtx1.z, vtx2.z] ] ];
};
GetSlopeVec: PROC[normal, edge: Triple, hermite: BOOLFALSE] RETURNS[slope: Triple] ~ {
Returns a vector normal to the 1st vector and in the plane defined by both vectors
Magnitude is scaled to equal that of second vector ("edge")
Magnitude is scaled to to proper length for Bezier inner control pt approximating circular arc
slope ← Cross[ Cross[normal, edge], normal ];
IF hermite
THEN slope ← G3dVector.Mul[ Nmlize[slope], Length[edge] ]
THEN slope ← G3dVector.Mul[ ScaleTangent[ slope, edge ], 3.0 ]
ELSE slope ← ScaleTangent[ slope, edge ];
};
GetNmlVec: PROC[vec: Triple, p: CtlPtInfo, reverse: BOOLFALSE]
    RETURNS
[nmlVec: Triple] ~ {
Returns a unitd vector normal to the 1st vector and in the plane defined by 2nd
normal: Triple ← [p.shade.exn, p.shade.eyn, p.shade.ezn];
nmlVec ← IF reverse
THEN Nmlize[ Cross[normal, vec] ]
ELSE Nmlize[ Cross[vec, normal] ];
};
ScaleTangent: PROC[tangent, edgeDir: Triple] RETURNS[Triple] ~ {
Scale tangent vector to proper length for Bezier inner control pt approximating circular arc
adjSide, oppSide, scale: REAL;
edgeLength: REAL ← Length[edgeDir];
hypotenuse: REAL ← Length[tangent];
IF edgeLength = 0.0 OR hypotenuse = 0.0 THEN RETURN [[0.0, 0.0, 0.0]];
adjSide ← ABS[Dot[edgeDir, tangent] / edgeLength];
oppSide ← IF adjSide >= hypotenuse
THEN 0.0
ELSE RealFns.SqRt[hypotenuse*hypotenuse - adjSide*adjSide];
scale ← IF oppSide/hypotenuse > .01
THEN edgeLength * 2.0 * (hypotenuse - adjSide) / (3.0 * oppSide * oppSide)
ELSE .333333 * edgeLength / hypotenuse;
RETURN[ G3dVector.Mul[ tangent, scale ] ];
};
TooClose: PROC[ context: Context, v: CtlPoint, outCode: OutCode, tol: REAL ]
   RETURNS[BOOLEAN] ~ {
Find max distance from screen edge on inside
maxDist: REAL ← 0.0;
IF v.sz = 0.0 THEN RETURN [FALSE];     -- invalid screen coord, ignore
IF outCode.left
THEN maxDist ← MAX[maxDist, v.sx - context.viewPort.x];
IF outCode.right
THEN maxDist ← MAX[maxDist, context.viewPort.w + context.viewPort.x - v.sx];
IF outCode.bottom
THEN maxDist ← MAX[maxDist, v.sy - context.viewPort.y];
IF outCode.top
THEN maxDist ← MAX[maxDist, context.viewPort.h + context.viewPort.y - v.sy];
IF maxDist < tol * closenessFactor THEN RETURN[TRUE] ELSE RETURN[FALSE];
};
EdgeStraight: PROC[context: Context, v1, v2: CtlPtInfo, slope1, slope2: Triple, tol: REAL]
     RETURNS[BOOL] ~ {
distance: REAL;
IF NOT useManhattan THEN {
Tests for perpendicular distance from slope vectors to edge vector (has bug somewhere!!)
pos1: Triple ← G3dClipXfmShade.XfmPtToDisplay[ context, -- near slope + near end, screen
Add[[v1.coord.ex, v1.coord.ey, v1.coord.ez], slope1 ]
];
pos2: Triple ← G3dClipXfmShade.XfmPtToDisplay[ context, -- far slope + near end, screen
Add[[v1.coord.ex, v1.coord.ey, v1.coord.ez], slope2 ]
];
Get unit length vector in edge direction
edge: Pair ← G2dVector.Unit[ [v2.coord.sx - v1.coord.sx, v2.coord.sy - v1.coord.sy] ];
Use cross product to get area which = height when base is unit length
distance1: REAL ← G2dVector.Cross[
[pos1.x - v1.coord.sx, pos1.y - v1.coord.sy], [edge.x, edge.y]
];
distance2: REAL ← G2dVector.Cross[
[pos2.x - v1.coord.sx, pos2.y - v1.coord.sy], [edge.x, edge.y]
];
distance ← (ABS[distance1] + ABS[distance2]) / 2;
}
ELSE {
Tests for slopes within n pixels of edge on screen
pos1: Triple ← G3dClipXfmShade.XfmPtToDisplay[ context,  -- near slope plus near end
Add[[v1.coord.ex, v1.coord.ey, v1.coord.ez], slope1 ]
];
pos2: Triple ← G3dClipXfmShade.XfmPtToDisplay[ context,  -- far slope plus near end
Add[[v1.coord.ex, v1.coord.ey, v1.coord.ez], slope2 ]
];
Slopes roughly same length as edge, so distance nears zero when edge straight
distance ← ABS[pos1.x - v2.coord.sx] + ABS[pos1.y - v2.coord.sy]
   + ABS[pos2.x - v2.coord.sx] + ABS[pos2.y - v2.coord.sy];
distance ← distance / 4.0;   -- summed manhattan distances / 4 estimates deviation
};
distance ← distance / 2.0;   -- heuristic (fudge factor)
IF distance <= tol THEN RETURN[TRUE] ELSE RETURN[FALSE];
};
PatchDepthSort: PROC[context: Context, p: PatchSequence]
     RETURNS[PatchSequence] ~ {
z: RealSequence ← NEW[RealSequenceRep[p.length]];
pOut: PatchSequence ← NEW [PatchSequenceRep[p.length]];
FOR i: NAT IN [0..p.length) DO      -- Get average of depths at vertices
z[i] ← 0.0;
FOR j: NAT IN [0..p[i].nVtces) DO z[i] ← z[i] + p[i][j].coord.ez; ENDLOOP;
z[i] ← z[i] / p[i].nVtces;
pOut[i] ← p[i];
ENDLOOP;
FOR i: NAT IN [1..p.length) DO
FOR j: NAT DECREASING IN [0..i) DO    -- bubble sort to increasing depth order
IF z[j+1] < z[j] THEN {
t: REAL ← z[j+1]; p: REF Patch ← pOut[j+1];
z[j+1] ← z[j];  pOut[j+1] ← pOut[j];
z[j] ← t;    pOut[j] ← p;
};
ENDLOOP;
ENDLOOP;
IF NOT context.antiAliasing THEN FOR i: NAT IN [0..p.length/2) DO-- re-order back-to-front
j: NAT ← p.length-1 - i;
tmpP: REF Patch ← pOut[i]; pOut[i] ← pOut[j]; pOut[j] ← tmpP;
ENDLOOP;
pOut.length ← p.length;
RETURN[pOut];
};
Display Procedures
DisplayNothing: PatchProc ~ {      -- dummy routine for timing tests
RETURN[ patch ];
};
DisplayPatchEdges: PatchProc ~ {
PROC[ context: Context, patch: REF Patch, data: REF ANYNIL ] RETURNS[REF Patch]
shape: Shape ← NARROW[ GetProp[patch.props, $Shape] ];
patchNo: NATNARROW[ GetProp[patch.props, $PatchNo], REF NAT ]^;
renderData: REF RenderData ← patch.renderData;
tangents: REF TangentSeqSeq ← NARROW[ GetProp[renderData.props, $PatchTangents] ];
corners: REF CornerSeqSeq ← NARROW[ GetProp[renderData.props, $PatchCorners] ];
xfm: Matrix ← G3dMatrix.Mul[shape.matrix, context.eyeSpaceXfm];
clr: RGB ← renderData.shadingClass.color;
npts: NATIF data # NIL THEN NARROW[data, REF NAT]^ ELSE 8;
FOR i: NAT IN [0..patch.nVtces) DO
j: NAT ← (i+1) MOD patch.nVtces;
outP: REF Patch ← GetPatch[npts];
outState: OutCode ← AllOut; inState: OutCode ← NoneOut;
spline: G3dSpline.Spline;
SELECT patch.type FROM
$PolygonWithNormals => spline ← GetEdgeCurveNmls[ patch[i], patch[j] ];
$PolygonWithTangents, $PolygonToTnsrPatch => spline ← GetEdgeCurveTngs[
patch[i], patch[j], tangents[patchNo][i]
];
ENDCASE => SIGNAL G3dRender.Error[$Unimplemented, "Unknown surface type"];
FOR k: NAT IN [0..npts) DO
OPEN outP[k].coord;
t: REAL ← 1.0 * k / (npts-1);
[[x, y, z]] ← G3dSpline.Position[spline, t];
[[ex, ey, ez], clip] ← G3dClipXfmShade.XfmPtToEyeSpace[context, [x, y, z], xfm];
clip ← G3dClipXfmShade.GetClipCodeForPt[ context, [ex, ey, ez] ];
IF clip = NoneOut
THEN [[sx, sy, sz]] ← G3dClipXfmShade.XfmPtToDisplay[context, [ex, ey, ez]];
outState ← LOOPHOLE[ Basics.BITAND[LOOPHOLE[outState], LOOPHOLE[clip]] ];
inState ← LOOPHOLE[ Basics.BITOR[LOOPHOLE[inState], LOOPHOLE[clip]] ];
[outP[k].shade.er, outP[k].shade.eg, outP[k].shade.eb] ← clr;  -- use shape color
ENDLOOP;
outP.type ← $PolyLine;
outP.nVtces ← npts;
outP.clipState ← IF inState = NoneOut
THEN in
ELSE IF outState # NoneOut THEN out ELSE clipped;
outP.renderData ← patch.renderData;
outP.props ← patch.props;
[outP] ← G3dSortandDisplay.OutputPolygon[context, outP];
ENDLOOP;
RETURN[NIL];
};
Initialization of Classes
InitClasses: PROC[] ~ {    -- register procedures for basic surface types
standardClass: ShapeClass ← G3dRender.GetShapeClass[$ConvexPolygon];
standardClass.type ← $PolygonWithNormals;
standardClass.validate ← ValidatePolyWNmls;
standardClass.displayPatch ← DisplayPatchNmls;
G3dRender.RegisterShapeClass[standardClass, $PolygonWithNormals];
standardClass.type ← $PolygonWithTangents;
standardClass.validate ← ValidatePolyWTangents;
standardClass.displayPatch ← DisplayPatchTngs;
G3dRender.RegisterShapeClass[standardClass, $PolygonWithTangents];
standardClass.type ← $PolygonToTnsrPatch;
standardClass.displayPatch ← DisplayPatchTnsr;
G3dRender.RegisterShapeClass[standardClass, $PolygonToTnsrPatch];
standardClass.type ← $PolygonNoImage;
standardClass.validate ← G3dSortandDisplay.ValidatePolyhedron;
standardClass.displayPatch ← DisplayNothing;
G3dRender.RegisterShapeClass[standardClass, $PolygonNoImage];
};
Utilities for expansion of non-planar polygons
TriangulateAndDisplay: PatchProc ~ {
Subdivide polygon to triangular patches, depth sort, then call TriangleDisplay on each one
patchNo: NAT ← 0;
tol: REALIF context.antiAliasing THEN maxDeviation ELSE maxJaggyDeviation; 
tangents: REF TangentSeq ← NARROW[ data ];
IF patch.nVtces < 3 THEN SIGNAL G3dRender.Error[$MisMatch, "Not enough vertices"];
IF NARROW[ GetProp[patch.props, $PatchNo], REF NAT ]^ = stopAtPoly
THEN patchNo ← stopAtPoly;   -- debug breakpoint here to stop at a given polygon
FOR i: NAT IN [0..patch.nVtces) DO    -- normalize normal vectors just to be sure
OPEN patch[i].shade;
[[exn, eyn, ezn]] ← Nmlize[[exn, eyn, ezn]];
ENDLOOP;
IF patch.nVtces = 3
THEN {
IF patch.type = $PolygonWithTangents THEN patch.props ← PutPropSafely[
patch.props, $Tangents,
NEW[TangentTriple ← [
NmlizeTangentSet[ patch[0], patch[1], tangents[0] ],
NmlizeTangentSet[ patch[1], patch[2], tangents[1] ],
NmlizeTangentSet[ patch[2], patch[0], tangents[2] ]
]]
];
TriangleDisplay[context, patch, 0, tol];
}
ELSE {
outPatch: PatchSequence ← NEW [PatchSequenceRep[patch.nVtces]];
midPt: CtlPtInfo ← GetCenterPt[context, patch, tangents];
midTangents: REF TangentSeq ← NARROW[ midPt.data ];
FOR i: NAT IN [0..patch.nVtces) DO
j: NAT ← (i + 1) MOD patch.nVtces;  -- next vertex
outPatch[i] ← GetPatch[3];  -- released by display action
outPatch[i].type ← patch.type;
outPatch[i].oneSided ← patch.oneSided;
outPatch[i].nVtces ← 3;
outPatch[i].clipState ← patch.clipState;
outPatch[i].dir ← unknown;
outPatch[i].renderData ← patch.renderData;
outPatch[i].props ← patch.props;
outPatch[i][0] ← patch[i];
outPatch[i][1] ← patch[(i+1) MOD patch.nVtces];
outPatch[i][2] ← midPt;
IF patch.type = $PolygonWithTangents THEN outPatch[i].props ← PutPropSafely[
outPatch[i].props, $Tangents,
NEW[TangentTriple ← [
NmlizeTangentSet[ outPatch[i][0], outPatch[i][1], tangents[i] ], -- orig. edge
midTangents[j],         -- in to middle, from next vertex
[ t0: midTangents[i].t1, et0: midTangents[i].et1,  -- back out to current vertex
t1: midTangents[i].t0, et1: midTangents[i].et0 ]
]]
];
IF outPatch[i].clipState # in THEN G3dClipXfmShade.GetPatchClipState[ outPatch[i] ];
ENDLOOP;
outPatch.length ← patch.nVtces;
outPatch ← PatchDepthSort[ context, outPatch ];    -- sort to depth order
FOR i: NAT IN [0..patch.nVtces) DO
TriangleDisplay[context, outPatch[i], 0, tol];   ReleasePatch[outPatch[i]];
ENDLOOP;
};
RETURN[patch];
};
GetCenterPt: PROC[context: Context, patch: REF Patch, tangents: REF TangentSeq ← NIL]
    RETURNS[midPt: CtlPtInfo] ~ {
Find central point in polygon for use in triangulation
tol: REALIF context.antiAliasing THEN maxDeviation ELSE maxJaggyDeviation; 
IF patch.nVtces <= 3 THEN SIGNAL G3dRender.Error[$MisMatch, "Not enough vertices"];
{
midTangents: REF TangentSeq ← IF patch.type = $PolygonWithTangents
         OR
patch.type = $PolygonToTnsrPatch
THEN NEW[ TangentSeq[patch.nVtces] ] ELSE NIL;
halfVtces: NAT ← patch.nVtces / 2;
loopEnd: NATIF halfVtces * 2 # patch.nVtces THEN patch.nVtces
              ELSE
patch.nVtces/2;
midPt.shade.r ← midPt.shade.g ← midPt.shade.b ← midPt.shade.t ← 0.0;
FOR i: NAT IN [0..loopEnd) DO -- guess at middle point, (odd # of Vtces will be too flat)
pt: CtlPtInfo; flat: BOOLEAN;  -- sum vertices given by each opposing vertex pair
t, t0, t1: TangentSet;
j: NAT ← (i + halfVtces) MOD patch.nVtces;  -- vertex across poly from this one
IF patch.type = $PolygonWithTangents OR patch.type = $PolygonToTnsrPatch
THEN { -- get sum of outgoing and incoming tangents at vertex for midpt direction
offst0, offst1: Triple;
k: NAT ← (i + patch.nVtces-1) MOD patch.nVtces; -- previous vertex in polygon
offst0 ← tangents[i].et0;
t.et0 ← Add[ tangents[i].et0, tangents[k].et1 ];
k ← (j + patch.nVtces-1) MOD patch.nVtces; -- previous vertex across polygon
offst1 ← tangents[k].et1;
t.et1 ← Add[ tangents[j].et0, tangents[k].et1 ];
t ← NmlizeTangentSet[ patch[i], patch[j], t ];
[pt, t0, t1, flat] ← CurveDivideTan[
context, patch[i], patch[j], t.et0, t.et1, offst0, offst1, 0.5, tol];
In following, t0 toward center at vertex, t1 away from center at center
midTangents[i] ← t0;
midTangents[j] ← [t0: t1.t1, et0: t1.et1, t1: t1.t0, et1: t1.et0 ]; -- reverse order
}
ELSE {
[pt, flat] ← TriangleCurveDivideNml[
context, patch[i], patch[j], patch.renderData, 0, tol
];
midPt.shade.exn ← midPt.shade.exn + pt.shade.exn;
midPt.shade.eyn ← midPt.shade.eyn + pt.shade.eyn;
midPt.shade.ezn ← midPt.shade.ezn + pt.shade.ezn;
};
midPt.coord.ex ← midPt.coord.ex + pt.coord.ex;
midPt.coord.ey ← midPt.coord.ey + pt.coord.ey;
midPt.coord.ez ← midPt.coord.ez + pt.coord.ez;
midPt.shade.r ← midPt.shade.r + pt.shade.r;
midPt.shade.g ← midPt.shade.g + pt.shade.g;
midPt.shade.b ← midPt.shade.b + pt.shade.b;
midPt.shade.t ← midPt.shade.t + pt.shade.t;
IF patch.renderData.shadingClass.texture # NIL THEN {    -- for textures
midPt.coord.x ← midPt.coord.x + pt.coord.x;
midPt.coord.y ← midPt.coord.y + pt.coord.y;
midPt.coord.z ← midPt.coord.z + pt.coord.z;
midPt.shade.txtrX ← midPt.shade.txtrX + pt.shade.txtrX;
midPt.shade.txtrY ← midPt.shade.txtrY + pt.shade.txtrY;
midPt.shade.xn ← midPt.shade.xn + pt.shade.xn;
midPt.shade.yn ← midPt.shade.yn + pt.shade.yn;
midPt.shade.zn ← midPt.shade.zn + pt.shade.zn;
};
ENDLOOP;
midPt.coord.ex ← midPt.coord.ex / loopEnd;  -- get average by dividing by no. edges
midPt.coord.ey ← midPt.coord.ey / loopEnd;
midPt.coord.ez ← midPt.coord.ez / loopEnd;
IF patch.type = $PolygonWithTangents OR patch.type = $PolygonToTnsrPatch
THEN {    -- get mid-normal by cross-products on mid-tangents
FOR i: NAT IN [0..loopEnd-1) DO
[[midPt.shade.exn, midPt.shade.eyn, midPt.shade.ezn]] ← Add[
[midPt.shade.exn, midPt.shade.eyn, midPt.shade.ezn],
Cross[midTangents[i].et1, midTangents[i+1].et1]
];
ENDLOOP;
[[midPt.shade.exn, midPt.shade.eyn, midPt.shade.ezn]] ← Nmlize[
[midPt.shade.exn, midPt.shade.eyn, midPt.shade.ezn]
];
FOR i: NAT IN [0..patch.nVtces) DO -- Rebuild midtangents using middle point
midTangents[i].et0 ← GetSlopeVec[
[patch[i].shade.exn, patch[i].shade.eyn, patch[i].shade.ezn],
DiffPosnsCtlPt[midPt.coord, patch[i].coord, $Eye]
];
midTangents[i].et1 ← GetSlopeVec[
[midPt.shade.exn, midPt.shade.eyn, midPt.shade.ezn],
DiffPosnsCtlPt[patch[i].coord, midPt.coord, $Eye]
];
ENDLOOP;
}
ELSE {
midPt.shade.exn ← midPt.shade.exn / loopEnd;
midPt.shade.eyn ← midPt.shade.eyn / loopEnd;
midPt.shade.ezn ← midPt.shade.ezn / loopEnd;
};
midPt.shade.r ← midPt.shade.r / loopEnd;
midPt.shade.g ← midPt.shade.g / loopEnd;
midPt.shade.b ← midPt.shade.b / loopEnd;
midPt.shade.t ← midPt.shade.t / loopEnd;
IF patch.renderData.shadingClass.texture # NIL THEN {    -- for solid textures
midPt.coord.x ← midPt.coord.x / loopEnd;
midPt.coord.y ← midPt.coord.y / loopEnd;
midPt.coord.z ← midPt.coord.z / loopEnd;
midPt.shade.txtrX ← midPt.shade.txtrX / loopEnd;
midPt.shade.txtrY ← midPt.shade.txtrY / loopEnd;
IF patch.type = $PolygonWithTangents OR patch.type = $PolygonToTnsrPatch
THEN {    -- get mid-normal by cross-products on mid-tangents
FOR i: NAT IN [0..loopEnd-1) DO
[[midPt.shade.xn, midPt.shade.yn, midPt.shade.zn]] ← Add[
[midPt.shade.xn, midPt.shade.yn, midPt.shade.zn],
Cross[midTangents[i].t1, midTangents[i+1].t1]
];
ENDLOOP;
[[midPt.shade.xn, midPt.shade.yn, midPt.shade.zn]] ← Nmlize[
[midPt.shade.xn, midPt.shade.yn, midPt.shade.zn]
];
FOR i: NAT IN [0..patch.nVtces) DO-- Rebuild midtangents using middle point
midTangents[i].t0 ← GetSlopeVec[
[patch[i].shade.xn, patch[i].shade.yn, patch[i].shade.zn],
DiffPosnsCtlPt[midPt.coord, patch[i].coord]
];
midTangents[i].t1 ← GetSlopeVec[
[midPt.shade.xn, midPt.shade.yn, midPt.shade.zn],
DiffPosnsCtlPt[patch[i].coord, midPt.coord]
];
ENDLOOP;
}
ELSE {
midPt.shade.xn ← midPt.shade.xn / loopEnd;
midPt.shade.yn ← midPt.shade.yn / loopEnd;
midPt.shade.zn ← midPt.shade.zn / loopEnd;
};
};
{ OPEN midPt.coord;
shape: Shape ← NARROW[ GetProp[patch.props, $Shape] ];
clip ← G3dClipXfmShade.GetClipCodeForPt[context, [ ex, ey, ez] ];
IF clip = NoneOut
THEN [ [sx, sy, sz] ] ← G3dClipXfmShade.XfmPtToDisplay[ context, [ex, ey, ez], shape ];
};
midPt.data ← midTangents;
};
};
TriangleDisplay: PROC[ context: Context, p: REF Patch, level: NAT, tol: REAL] ~ {
Recursively divide triangular patches to straight edges, then display
subP: PatchSequence ← NIL;
IF context.stopMe^ THEN RETURN[];  -- shut down if stop signal received
IF p.dir = unknown AND p.type = $PolygonWithTangents
THEN [] ← TriangleTngsBackFacing[p];  -- backface test
IF ( p.clipState # out ) AND ( NOT p.oneSided OR NOT p.dir = back ) THEN {
IF level < recurseLimit THEN subP ← SubdivideTriangle[context, p, level, tol];
IF subP = NIL
THEN {          -- recursion limit hit or edges all straight
IF showLines THEN p.type ← $PolyLine ELSE p.type ← $ConvexPolygon;
G3dClipXfmShade.ShadePoly[context, p];
IF (p[0].coord.sz = 0.0 OR p[1].coord.sz = 0.0 OR p[2].coord.sz = 0.0 )
THEN SIGNAL G3dRender.Error[$MisMatch, "Zeroed vertex"];
[p] ← G3dSortandDisplay.OutputPolygon[context, p];      -- display
}
ELSE {
subP ← PatchDepthSort[ context, subP ]; -- sort to display order
TriangleDisplay[context, subP[0], level+1, tol];  ReleasePatch[subP[0]];
TriangleDisplay[context, subP[1], level+1, tol];  ReleasePatch[subP[1]];
TriangleDisplay[context, subP[2], level+1, tol];  ReleasePatch[subP[2]];
TriangleDisplay[context, subP[3], level+1, tol];  ReleasePatch[subP[3]];
};
};
};
SubdivideTriangle: PROC[context: Context, p: REF Patch, level: NAT, tol: REAL]
    RETURNS[PatchSequence] ~ {
Divide triangular patch into four subtriangles
shape: Shape ← NARROW[ GetProp[p.props, $Shape] ];
v0, v1, v2: CtlPtInfo;  flat0, flat1, flat2: BOOLEAN;
t: REF TangentTriple ← NARROW[GetProp[p.props, $Tangents] ];
t00, t01, t10, t11, t20, t21, tm0, tm1, tm2: TangentSet;
outPatch: PatchSequence ← NEW[PatchSequenceRep[4]];
IF p.type = $PolygonWithTangents AND t # NIL
THEN {  -- get midpoints and edge tangents
[v0, t00, t01, flat0] ← CurveDivideTan[
context, p[0], p[1], t[0].et0, t[0].et1,
GetNmlVec[t[0].et0, p[0]], GetNmlVec[t[0].et1, p[1], TRUE], 0.5, tol ];
[v1, t10, t11, flat1] ← CurveDivideTan[
context, p[1], p[2], t[1].et0, t[1].et1,
GetNmlVec[t[1].et0, p[1]], GetNmlVec[t[1].et1, p[2], TRUE], 0.5, tol];
[v2, t20, t21, flat2] ← CurveDivideTan[
context, p[2], p[0], t[2].et0, t[2].et1,
GetNmlVec[t[2].et0, p[2]], GetNmlVec[t[2].et1, p[0], TRUE], 0.5, tol];
get inner edge tangents based on midpoints
tm0.et0 ← GetSlopeVec[ [v0.shade.exn, v0.shade.eyn, v0.shade.ezn],
       DiffPosnsCtlPt[v1.coord, v0.coord, $Eye] ];
tm0.et1 ← GetSlopeVec[ [v1.shade.exn, v1.shade.eyn, v1.shade.ezn],
       DiffPosnsCtlPt[v0.coord, v1.coord, $Eye] ];
tm1.et0 ← GetSlopeVec[ [v1.shade.exn, v1.shade.eyn, v1.shade.ezn],
       DiffPosnsCtlPt[v2.coord, v1.coord, $Eye] ];
tm1.et1 ← GetSlopeVec[ [v2.shade.exn, v2.shade.eyn, v2.shade.ezn],
       DiffPosnsCtlPt[v1.coord, v2.coord, $Eye] ];
tm2.et0 ← GetSlopeVec[ [v2.shade.exn, v2.shade.eyn, v2.shade.ezn],
       DiffPosnsCtlPt[v0.coord, v2.coord, $Eye] ];
tm2.et1 ← GetSlopeVec[ [v0.shade.exn, v0.shade.eyn, v0.shade.ezn],
       DiffPosnsCtlPt[v2.coord, v0.coord, $Eye] ];
}
ELSE {
[v0, flat0] ← TriangleCurveDivideNml[context, p[0], p[1], p.renderData, level, tol];
[v1, flat1] ← TriangleCurveDivideNml[context, p[1], p[2], p.renderData, level, tol];
[v2, flat2] ← TriangleCurveDivideNml[context, p[2], p[0], p.renderData, level, tol];
};
IF flat0 AND flat1 AND flat2 AND stopIfStraight THEN RETURN[NIL]; -- if all straight, done
{ -- Reverse mid-normals if on wrong side of triangle plane by > 30 deg.
normal: Triple ← Nmlize[ Cross[
DiffPosnsCtlPt[p[1].coord, p[0].coord, $Eye], DiffPosnsCtlPt[p[2].coord, p[0].coord, $Eye]
] ];
check: REAL;
check ← Dot[[v0.shade.exn, v0.shade.eyn, v0.shade.ezn], normal ];
IF check < -0.5 THEN [[v0.shade.exn, v0.shade.eyn, v0.shade.ezn]] ←
      Negate[[v0.shade.exn, v0.shade.eyn, v0.shade.ezn]];
check ← Dot[[v1.shade.exn, v1.shade.eyn, v1.shade.ezn], normal ];
IF check < -0.5 THEN [[v1.shade.exn, v1.shade.eyn, v1.shade.ezn]] ←
      Negate[[v1.shade.exn, v1.shade.eyn, v1.shade.ezn]];
check ← Dot[[v2.shade.exn, v2.shade.eyn, v2.shade.ezn], normal ];
IF check < -0.5 THEN [[v2.shade.exn, v2.shade.eyn, v2.shade.ezn]] ←
      Negate[[v2.shade.exn, v2.shade.eyn, v2.shade.ezn]];
};
FOR i: NAT IN [0..4) DO
outPatch[i] ← GetPatch[3];    -- allocate 3 point patch
outPatch[i].type ← p.type;
outPatch[i].oneSided ← p.oneSided;
outPatch[i].nVtces ← 3;
outPatch[i].clipState ← p.clipState;
outPatch[i].dir ← p.dir;
outPatch[i].renderData ← p.renderData;
outPatch[i].props ← p.props;
ENDLOOP;
{ OPEN v0.coord; clip ← G3dClipXfmShade.GetClipCodeForPt[context, [ ex, ey, ez] ];
IF clip = NoneOut
THEN [ [sx, sy, sz] ] ← G3dClipXfmShade.XfmPtToDisplay[ context, [ex, ey, ez], shape ]; };
{ OPEN v1.coord; clip ← G3dClipXfmShade.GetClipCodeForPt[context, [ ex, ey, ez] ];
IF clip = NoneOut
THEN [ [sx, sy, sz] ] ← G3dClipXfmShade.XfmPtToDisplay[ context, [ex, ey, ez], shape ]; };
{ OPEN v2.coord; clip ← G3dClipXfmShade.GetClipCodeForPt[context, [ ex, ey, ez] ];
IF clip = NoneOut
THEN [ [sx, sy, sz] ] ← G3dClipXfmShade.XfmPtToDisplay[ context, [ex, ey, ez], shape ]; };
outPatch[0][0] ← p[0]; outPatch[0][1] ← v0; outPatch[0][2] ← v2;
outPatch[1][0] ← p[1]; outPatch[1][1] ← v1; outPatch[1][2] ← v0;
outPatch[2][0] ← p[2]; outPatch[2][1] ← v2; outPatch[2][2] ← v1;
outPatch[3][0] ← v0; outPatch[3][1] ← v1; outPatch[3][2] ← v2;
IF p.type = $PolygonWithTangents AND t # NIL THEN {
outPatch[0].props ← PutPropSafely[ outPatch[0].props, $Tangents,
       NEW[TangentTriple ← [t00, [tm2.t1, tm2.et1, tm2.t0, tm2.et0], t21] ] ];
outPatch[1].props ← PutPropSafely[ outPatch[1].props, $Tangents,
       NEW[TangentTriple ← [t10, [tm0.t1, tm0.et1, tm0.t0, tm0.et0], t01] ] ];
outPatch[2].props ← PutPropSafely[ outPatch[2].props, $Tangents,
      NEW[TangentTriple ← [t20, [tm1.t1, tm1.et1, tm1.t0, tm1.et0], t11] ] ];
outPatch[3].props ← PutPropSafely[ outPatch[3].props, $Tangents,
          NEW[TangentTriple ← [tm0, tm1, tm2] ] ];
};
FOR i: NAT IN [0..4) DO G3dClipXfmShade.GetPatchClipState[ outPatch[i] ]; ENDLOOP; --bad!!
outPatch.length ← 4;
RETURN[ outPatch ]; -- return four sub-patches
};
Procedures for expansion of non-planar polygons using normals at vertices
DisplayPatchNmls: PatchProc ~ {
renderStyle: RenderStyle;
patch.type ← $PolygonWithNormals;
WITH patch.renderData.shadingClass.renderMethod SELECT FROM
style: REF RenderStyle => {
renderStyle ← style^;
SELECT renderStyle FROM
lines => { [] ← DisplayPatchEdges[context, patch]; RETURN[patch]; };
ENDCASE;
};
ENDCASE;
Do we have texture?
IF patch.renderData.shadingClass.texture # NIL THEN {
txtrRange: REF Pair ← NARROW[
GetProp[patch.renderData.shadingProps, $TxtrCoordRange]
];
IF txtrRange # NIL THEN G3dMappedAndSolidTexture.AdjustTexture[ -- fix texture seams
patch, patch.renderData.shadingClass.texture, txtrRange^
]
};
patch ← TriangulateAndDisplay[context, patch]; -- make triangles, recursively subdivide
RETURN[ patch ];
};
ValidatePolyWNmls: ShapeProc ~ {
PROC[context: Context, shape: Shape, data: REF] RETURNS[Shape];
render: REF RenderData;
shape ← G3dSortandDisplay.ValidatePolyhedron[context, shape]; -- Update shading, xfm
render ← G3dRender.RenderDataFrom[shape];
IF shape.clipState # out THEN {
FOR i: NAT IN [0..render.patch.length) DO    -- check backfacing for curved surface
neg, pos: BOOLEANFALSE;
FOR j: NAT IN [0..render.patch[i].nVtces) DO
OPEN render.patch[i][j];
awayness: REAL ← G3dVector.Dot[   -- normal front or back facing?
[coord.ex, coord.ey, coord.ez], [shade.exn, shade.eyn, shade.ezn]
];
IF awayness > 0.0 THEN pos ← TRUE ELSE IF awayness < 0.0 THEN neg ← TRUE;
ENDLOOP;
IF pos AND NOT neg THEN render.patch[i].dir ← back
ELSE IF neg AND NOT pos THEN render.patch[i].dir ← front
ELSE render.patch[i].dir ← unknown;
ENDLOOP;
};
RETURN[ shape ];
};
GetEdgeCurveNmls: PROC[p1, p2: CtlPtInfo] RETURNS[spline: G3dSpline.Spline] ~ {
edge: Triple ← DiffPosnsCtlPt[p2.coord, p1.coord, $Object];
pos1: Triple ← [p1.coord.x, p1.coord.y, p1.coord.z];
pos2: Triple ← [p2.coord.x, p2.coord.y, p2.coord.z];
slope1: Triple ← GetSlopeVec[[p1.shade.xn, p1.shade.yn, p1.shade.zn], edge, TRUE];
slope2: Triple ← GetSlopeVec[[p2.shade.xn, p2.shade.yn, p2.shade.zn], edge, TRUE];
spline ← G3dSpline.SplineFromHermite[[pos1, pos2, slope1, slope2]];
};
TriangleNmlsBackFacing: PROC[p: REF Patch] RETURNS [BOOLEAN] ~ {
Use slopes to make corner triangles and internal hexagon, evaluate at corners of hexagon
FacingFromTrio: PROC[p0, p1, p2: Triple] ~ {
direction: Triple ← Cross[
[p1.x - p0.x, p1.y - p0.y, p1.z - p0.z],
[p2.x - p0.x, p2.y - p0.y, p2.z - p0.z]
];
dotDir: REAL ← Dot[p0 , direction];
IF dotDir > 0.0
THEN back ← TRUE
ELSE IF dotDir < 0.0 THEN front ← TRUE ELSE unknown ← TRUE;
};
back, front, unknown: BOOLEANFALSE;
t: REF TangentTriple ← NARROW[ GetProp[p.props, $Tangents] ];
pts: ARRAY [0..9) OF Triple;
FOR i: NAT IN [0..3) DO
pts[i*3] ← [ p[i].coord.ex, p[i].coord.ey, p[i].coord.ez ];
ENDLOOP;
FOR i: NAT IN [0..3) DO
i1: NAT ← i*3+1; i2: NAT ← i*3+2;
pts[i*3+1] ← Add[ pts[i*3], t[i].et0];
pts[i*3+2] ← Add[ pts[((i+1) MOD 3) * 3], t[i].et1];
ENDLOOP;
FOR i: NAT IN [0..3) DO
j: NAT ← i*3; j1: NAT ← j+1; j2: NAT ← (j+8) MOD 9;
FacingFromTrio[ pts[j], pts[j1], pts[j2] ];
FacingFromTrio[ pts[j1], pts[j1+1], pts[j2] ];
FacingFromTrio[ pts[j2], pts[j1], pts[j2-1] ];
ENDLOOP;
IF (back AND front) OR unknown
THEN p.dir ← unknown
ELSE IF back
THEN p.dir ← back
ELSE IF front
THEN p.dir ← front
ELSE p.dir ← unknown;
IF p.dir = back THEN RETURN[TRUE] ELSE RETURN[FALSE];
};
TriangleCurveDivideNml: PROC[ context: Context, vtx0, vtx1: CtlPtInfo,
           renderData: REF RenderData, level: NAT, tol: REAL ]
         RETURNS[pt: CtlPtInfo, flat: BOOLEAN] ~ {
Returns parametric midpoint of curve given by v0-v1
t: REAL ~ 0.5;
t2: REAL ~ 0.25;
t3: REAL ~ 0.125;
b0: REAL ~ 1.0 + 0.0*t - 3.*t2 + 2.*t3; -- value of basis functions at t=1/2
b1: REAL ~ 0.0 + 0.0*t + 3.*t2 - 2.*t3;
b2: REAL ~ 0.0 + 1.0*t - 2.*t2 + 1.*t3;
b3: REAL ~ 0.0 + 0.0*t - 1.*t2 + 1.*t3;
s0: REAL ~    0.0 - 3.*2.*t + 2.*3.*t2; -- slope of basis functions at t=1/2
s1: REAL ~   0.0 + 3.*2.*t - 2.*3.*t2;
s2: REAL ~ 1.0 - 2.*2.*t + 1.*3.*t2;
s3: REAL ~   0.0 - 1.*2.*t + 1.*3.*t2;
pos1, pos2, slope1, slope2, slopeMid: Triple;
Ensure consistent evaluation for arithmetic stability
v0: CtlPtInfo ← IF vtx0.coord.x < vtx1.coord.x THEN vtx0 ELSE vtx1;
v1: CtlPtInfo ← IF vtx0.coord.x < vtx1.coord.x THEN vtx1 ELSE vtx0;
clipStraight, tooShort: BOOLEANFALSE;
IF v1.coord.clip # NoneOut AND v0.coord.clip # NoneOut
THEN clipStraight ← TRUE
ELSE IF v0.coord.clip # NoneOut
THEN clipStraight ← TooClose[context, v1.coord, v0.coord.clip, tol]
ELSE IF v1.coord.clip # NoneOut
THEN clipStraight ← TooClose[context, v0.coord, v1.coord.clip, tol];
{ edge: Triple ← DiffPosnsCtlPt[v1.coord, v0.coord, $Eye];
IF Length[edge] <= G3dScanConvert.justNoticeable
THEN tooShort ← TRUE
ELSE {
slope1 ← GetSlopeVec[[v0.shade.exn, v0.shade.eyn, v0.shade.ezn], edge, TRUE];
slope2 ← GetSlopeVec[[v1.shade.exn, v1.shade.eyn, v1.shade.ezn], edge, TRUE];
};
};
pos1 ← [v0.coord.ex, v0.coord.ey, v0.coord.ez];
pos2 ← [v1.coord.ex, v1.coord.ey, v1.coord.ez];
IF clipStraight OR tooShort OR EdgeStraight[context, v0, v1, slope1, slope2, tol]
THEN {         -- straight edge, average endpoints for midpoint
pt.coord.ex ← (pos1.x + pos2.x) / 2; pt.shade.exn ← (v0.shade.exn + v1.shade.exn) / 2;
pt.coord.ey ← (pos1.y + pos2.y) / 2; pt.shade.eyn ← (v0.shade.eyn + v1.shade.eyn) / 2;
pt.coord.ez ← (pos1.z + pos2.z) / 2; pt.shade.ezn ← (v0.shade.ezn + v1.shade.ezn) / 2;
[[pt.shade.exn, pt.shade.eyn, pt.shade.ezn]] ← Nmlize[
[pt.shade.exn, pt.shade.eyn, pt.shade.ezn]
];
flat ← TRUE;
}
ELSE {         -- not straight, evaluate curve
slopeMid.x ← pos1.x*s0 + pos2.x*s1 + slope1.x*s2 + slope2.x*s3;
slopeMid.y ← pos1.y*s0 + pos2.y*s1 + slope1.y*s2 + slope2.y*s3;
slopeMid.z ← pos1.z*s0 + pos2.z*s1 + slope1.z*s2 + slope2.z*s3;
pt.coord.ex ← pos1.x*b0 + pos2.x*b1 + slope1.x*b2 + slope2.x*b3;
pt.coord.ey ← pos1.y*b0 + pos2.y*b1 + slope1.y*b2 + slope2.y*b3;
pt.coord.ez ← pos1.z*b0 + pos2.z*b1 + slope1.z*b2 + slope2.z*b3;
[[pt.shade.exn, pt.shade.eyn, pt.shade.ezn]] ← Nmlize[Add[
Do with both ends to get consistent mid normal, curve is likely to be non planar
Nmlize[GetSlopeVec[slopeMid, [v0.shade.exn, v0.shade.eyn, v0.shade.ezn], TRUE]],
Nmlize[GetSlopeVec[slopeMid, [v1.shade.exn, v1.shade.eyn, v1.shade.ezn], TRUE]]
]];
flat ← FALSE;
};
pt.shade.r ← (v0.shade.r + v1.shade.r) / 2;
pt.shade.g ← (v0.shade.g + v1.shade.g) / 2;
pt.shade.b ← (v0.shade.b + v1.shade.b) / 2;
pt.shade.t ← (v0.shade.t + v1.shade.t) / 2;
IF renderData.shadingClass.texture # NIL THEN {    -- for solid textures
pos1 ← [v0.coord.x, v0.coord.y, v0.coord.z];
pos2 ← [v1.coord.x, v1.coord.y, v1.coord.z];
SELECT renderData.class.type FROM
$PolygonWithNormals => {
edge: Triple ← DiffPosnsCtlPt[v1.coord, v0.coord, $Object];
slope1 ← GetSlopeVec[[v0.shade.xn, v0.shade.yn, v0.shade.zn], edge, TRUE];
slope2 ← GetSlopeVec[[v1.shade.xn, v1.shade.yn, v1.shade.zn], edge, TRUE];
};
ENDCASE => SIGNAL G3dRender.Error[$Unimplemented, "Unknown surface type"];
IF clipStraight OR EdgeStraight[context, v0, v1, slope1, slope2, tol]
THEN {         -- straight edge, average endpoints for midpoint
pt.coord.x ← (pos1.x + pos2.x) / 2; pt.shade.xn ← (v0.shade.xn + v1.shade.xn) / 2;
pt.coord.y ← (pos1.y + pos2.y) / 2; pt.shade.yn ← (v0.shade.yn + v1.shade.yn) / 2;
pt.coord.z ← (pos1.z + pos2.z) / 2; pt.shade.zn ← (v0.shade.zn + v1.shade.zn) / 2;
}
ELSE {         -- not straight, evaluate curve
slopeMid.x ← pos1.x*s0 + pos2.x*s1 + slope1.x*s2 + slope2.x*s3;
slopeMid.y ← pos1.y*s0 + pos2.y*s1 + slope1.y*s2 + slope2.y*s3;
slopeMid.z ← pos1.z*s0 + pos2.z*s1 + slope1.z*s2 + slope2.z*s3;
pt.coord.x ← pos1.x*b0 + pos2.x*b1 + slope1.x*b2 + slope2.x*b3;
pt.coord.y ← pos1.y*b0 + pos2.y*b1 + slope1.y*b2 + slope2.y*b3;
pt.coord.z ← pos1.z*b0 + pos2.z*b1 + slope1.z*b2 + slope2.z*b3;
[[pt.shade.xn, pt.shade.yn, pt.shade.zn]] ← Add[
Do with both ends to get stable mid normal, must be a better way
GetSlopeVec[slopeMid, [v0.shade.xn, v0.shade.yn, v0.shade.zn], TRUE],
GetSlopeVec[slopeMid, [v1.shade.xn, v1.shade.yn, v1.shade.zn], TRUE]
];
};
pt.shade.txtrX ← (v0.shade.txtrX + v1.shade.txtrX) / 2; -- average texture coords
pt.shade.txtrY ← (v0.shade.txtrY + v1.shade.txtrY) / 2;
};
pt.data ← vtx0.data;
};
Procedures for expansion of non-planar polygons using tangent vectors at vertices
DisplayPatchTngs: PatchProc ~ {
renderStyle: RenderStyle;
WITH patch.renderData.shadingClass.renderMethod SELECT FROM
style: REF RenderStyle => {
renderStyle ← style^;
SELECT renderStyle FROM
lines => { [] ← DisplayPatchEdges[context, patch]; RETURN[patch]; };
ENDCASE;
};
ENDCASE;
patch.type ← $PolygonWithTangents;
IF data = NIL
THEN {
patchNo: NATNARROW[ GetProp[patch.props, $PatchNo], REF NAT ]^;
tangents: REF TangentSeqSeq ← NARROW[
GetProp[patch.renderData.props, $PatchTangents]
];
patch ← TriangulateAndDisplay[ context, patch, tangents[patchNo] ]; -- subdivide proc.
}
ELSE {
tangents: REF TangentSeq ← NARROW[ data ];
patch ← TriangulateAndDisplay[ context, patch, tangents ]; -- subdivide proc.
};
RETURN[ patch ];
};
ValidatePolyWTangents: ShapeProc ~ {
PROC[context: Context, shape: Shape, data: REF] RETURNS[Shape];
shape ← G3dSortandDisplay.ValidatePolyhedron[context, shape]; -- Update shading, xfm
IF GetProp[G3dRender.RenderDataFrom[shape].props, $PatchTangents] = NIL
THEN shape ← GetTangents[context, shape];  -- calculate tangent vectors if not read in
IF shape.clipState # out THEN {
xfm: Matrix ← G3dMatrix.Mul[shape.matrix, context.eyeSpaceXfm];
tangent: REF TangentSeqSeq ← NARROW[
GetProp[G3dRender.RenderDataFrom[shape].props, $PatchTangents]
];
FOR i: NAT IN [0..shape.surfaces.length) DO      -- transform tangent vectors
FOR j: NAT IN [0..tangent[i].length) DO OPEN tangent[i][j];
[ et0 ] ← G3dMatrix.TransformVec[ t0, xfm];
[ et1 ] ← G3dMatrix.TransformVec[ t1, xfm];
ENDLOOP;
ENDLOOP;
};
RETURN[ shape ];
};
GetEdgeCurveTngs: PROC[p1, p2: CtlPtInfo, tangent: TangentSet]
       RETURNS
[spline: G3dSpline.Spline] ~ {
b0, b1, b2, b3: Triple;
[b0, b1, b2, b3] ← TngsToBezKnots[ p1, p2, tangent ];
spline ← G3dSpline.SplineFromBezier[ [b0, b1, b2, b3] ];
};
TngsToBezKnots: PROC[p1, p2: CtlPtInfo, tangent: TangentSet]
          RETURNS[b0, b1, b2, b3: Triple] ~ {
Convert endpoints with tangents to Bezier knots
edgeDir: Triple;
b0 ← [ p1.coord.x, p1.coord.y, p1.coord.z ];
b3 ← [ p2.coord.x, p2.coord.y, p2.coord.z ];
edgeDir ← Sub3[b3, b0];
b1 ← Add[ b0, ScaleTangent[tangent.t0, edgeDir] ];
b2 ← Add[ b3, ScaleTangent[tangent.t1, Negate[edgeDir]] ];
};
TriangleTngsBackFacing: PROC[p: REF Patch] RETURNS [BOOLEAN] ~ {
Use slopes to make corner triangles and internal hexagon, evaluate at corners of hexagon
FacingFromTrio: PROC[p0, p1, p2: Triple] ~ {
direction: Triple ← Cross[
[p1.x - p0.x, p1.y - p0.y, p1.z - p0.z],
[p2.x - p0.x, p2.y - p0.y, p2.z - p0.z]
];
dotDir: REAL ← Dot[p0 , direction];
IF dotDir > 0.0
THEN back ← TRUE
ELSE IF dotDir < 0.0 THEN front ← TRUE ELSE unknown ← TRUE;
};
back, front, unknown: BOOLEANFALSE;
t: REF TangentTriple ← NARROW[ GetProp[p.props, $Tangents] ];
pts: ARRAY [0..9) OF Triple;
FOR i: NAT IN [0..3) DO
pts[i*3] ← [ p[i].coord.ex, p[i].coord.ey, p[i].coord.ez ];
ENDLOOP;
FOR i: NAT IN [0..3) DO
i1: NAT ← i*3+1; i2: NAT ← i*3+2;
pts[i*3+1] ← Add[ pts[i*3], t[i].et0];
pts[i*3+2] ← Add[ pts[((i+1) MOD 3) * 3], t[i].et1];
ENDLOOP;
FOR i: NAT IN [0..3) DO
j: NAT ← i*3; j1: NAT ← j+1; j2: NAT ← (j+8) MOD 9;
FacingFromTrio[ pts[j], pts[j1], pts[j2] ];
FacingFromTrio[ pts[j1], pts[j1+1], pts[j2] ];
FacingFromTrio[ pts[j2], pts[j1], pts[j2-1] ];
ENDLOOP;
IF (back AND front) OR unknown
THEN p.dir ← unknown
ELSE IF back
THEN p.dir ← back
ELSE IF front
THEN p.dir ← front
ELSE p.dir ← unknown;
IF p.dir = back THEN RETURN[TRUE] ELSE RETURN[FALSE];
};
NmlizeTangentSet: PROC[v0, v1: CtlPtInfo, tangent: TangentSet]
         RETURNS[TangentSet] ~ {
Convert endpoints with tangents to Properly scaled tangents
ep0: Triple ← [ v0.coord.ex, v0.coord.ey, v0.coord.ez ];
ep1: Triple ← [ v1.coord.ex, v1.coord.ey, v1.coord.ez ];
p0: Triple ← [ v0.coord.x, v0.coord.y, v0.coord.z ];
p1: Triple ← [ v1.coord.x, v1.coord.y, v1.coord.z ];
edgeDir: Triple ← Sub3[p1, p0];
t0: Triple ← ScaleTangent[ tangent.t0, edgeDir ];
t1: Triple ← ScaleTangent[ tangent.t1, Negate[edgeDir] ];
et0, et1: Triple;
edgeDir ← Sub3[ep1, ep0];
et0 ← ScaleTangent[ tangent.et0, edgeDir ];
et1 ← ScaleTangent[ tangent.et1, Negate[edgeDir] ];
RETURN[ [t0, et0, t1, et1] ];
};
IsStraight: PROC[context: Context, p0, p1, s0, s1: Triple, tol: REAL] RETURNS[BOOLEAN] ~{
DistOffEdge: PROC[v: Triple, line: Ray] RETURNS[REAL] ~ {
ptOnEdge: Triple ← G3dClipXfmShade.XfmPtToDisplay[
context, G3dVector.NearestToLine[line, v]
];
vs: Triple ← G3dClipXfmShade.XfmPtToDisplay[ context, v ];
RETURN[ RealFns.SqRt[ Sqr[ptOnEdge.x - vs.x] + Sqr[ptOnEdge.y - vs.y] ] ];
};
line: Ray ← [ p0, Sub3[p1, p0] ];
Get screen distance of inner knots from edge formed by outer knots
dist0: REAL ← DistOffEdge[Add[p0, s0], line];
dist1: REAL ← DistOffEdge[Add[p1, s1], line];
IF dist0 > tol OR dist1 > tol THEN RETURN[FALSE] ELSE RETURN[TRUE];
};
CurveDivideTan: PROC[ context: Context, vtx0, vtx1: CtlPtInfo,
           slope1, slope2, offst1, offst2: Triple, parameter, tol: REAL ]
         RETURNS[ pt: CtlPtInfo, t0, t1: TangentSet, flat: BOOLEAN ] ~ {
Returns point on curve given by v0-v1 and tangents and parameter
p1: REAL ← 1.0 - parameter; p2: REAL ← parameter;
pos1, pos2, edge, oPt: Triple;
v0: CtlPtInfo ← vtx0;
v1: CtlPtInfo ← vtx1;
clipStraight, tooShort: BOOLEANFALSE;
IF v1.coord.clip # NoneOut AND v0.coord.clip # NoneOut    -- not a good test!!
THEN clipStraight ← TRUE
ELSE IF v0.coord.clip # NoneOut
THEN clipStraight ← TooClose[context, v1.coord, v0.coord.clip, tol]
ELSE IF v1.coord.clip # NoneOut
THEN clipStraight ← TooClose[context, v0.coord, v1.coord.clip, tol];
edge ← DiffPosnsCtlPt[v1.coord, v0.coord, $Eye];
IF Length[edge] <= G3dScanConvert.justNoticeable THEN tooShort ← TRUE;
pos1 ← [v0.coord.ex, v0.coord.ey, v0.coord.ez];
pos2 ← [v1.coord.ex, v1.coord.ey, v1.coord.ez];
IF clipStraight OR tooShort OR IsStraight[ context, pos1, pos2, slope1, slope2, tol ]
THEN {         -- straight edge, average endpoints for midpoint
[[pt.coord.ex, pt.coord.ey, pt.coord.ez]] ← Add[ Mul[pos1, p1], Mul[pos2, p2] ];
t0.et0 ← Mul[ Sub3[pos2, pos1], p2/3.0]; t1.et0 ← Mul[ Sub3[pos2, pos1], p1/3.0];
t0.et1 ← Negate[ t0.et0 ]; t1.et1 ← Negate[ t1.et0 ];
oPt ← Add[ Mul[Add[pos1, offst1], p1], Mul[Add[pos2, offst2], p2] ];
flat ← TRUE;
}
ELSE {         -- not straight, evaluate midpoint by subdivision
b0: Triple ← pos1; b3: Triple ← pos2;
b1: Triple ← Add[pos1, slope1]; b2: Triple ← Add[pos2, slope2];
m01, m12, m23, mm0, mm1: Triple;
m01 ← Add[Mul[b0, p1], Mul[b1, p2]];
m12 ← Add[Mul[b1, p1], Mul[b2, p2]];
m23 ← Add[Mul[b2, p1], Mul[b3, p2]];
mm0 ← Add[Mul[m01, p1], Mul[m12, p2]]; mm1 ← Add[Mul[m23, p1], Mul[m12, p2]];
[[pt.coord.ex, pt.coord.ey, pt.coord.ez]] ← Add[ Mul[mm1, p1], Mul[mm0, p2] ];
t0.et0 ← Mul[ slope1, p2 ]; t0.et1 ← Sub3[mm0, [pt.coord.ex, pt.coord.ey, pt.coord.ez] ];
t1.et0 ← Sub3[mm1, [pt.coord.ex, pt.coord.ey, pt.coord.ez] ]; t1.et1 ← Mul[ slope2, p1 ];
b0 ← Add[pos1, offst1]; b3 ← Add[pos2, offst2]; -- get offset midpoint
b1 ← Add[b0, slope1]; b2 ← Add[b3, slope2];
m01 ← Add[Mul[b0, p1], Mul[b1, p2]];
m12 ← Add[Mul[b1, p1], Mul[b2, p2]];
m23 ← Add[Mul[b2, p1], Mul[b3, p2]];
mm0 ← Add[Mul[m01, p1], Mul[m12, p2]]; mm1 ← Add[Mul[m23, p1], Mul[m12, p2]];
oPt ← Add[ Mul[mm1, p1], Mul[mm0, p2] ];
flat ← FALSE;
};
[[pt.shade.exn, pt.shade.eyn, pt.shade.ezn]] ← Nmlize[
Cross[ Sub3[ oPt, [pt.coord.ex, pt.coord.ey, pt.coord.ez] ], t1.et0 ]
];
pt.shade.r ← p1 * v0.shade.r + p2 * v1.shade.r;
pt.shade.g ← p1 * v0.shade.g + p2 * v1.shade.g;
pt.shade.b ← p1 * v0.shade.b + p2 * v1.shade.b;
pt.shade.t ← p1 * v0.shade.t + p2 * v1.shade.t;
IF shape.shadingClass.texture # NIL THEN {    -- for solid textures
lerpProc: CtlPtInfoProc ← shape.shadingClass.lerpVtxAux;
pos1 ← [v0.coord.x, v0.coord.y, v0.coord.z];
pos2 ← [v1.coord.x, v1.coord.y, v1.coord.z];
SELECT shape.class.type FROM
$PolygonWithTangents => {
edge: Triple ← DiffPosnsCtlPt[v1.coord, v0.coord, $Object];
slope1 ← IF tan.t0 # nullTriple
THEN tan.t0
ELSE GetSlopeVec[ [v0.shade.xn, v0.shade.yn, v0.shade.zn], edge ];
slope2 ← IF tan.t1 # nullTriple
THEN tan.t1
ELSE GetSlopeVec[ [v1.shade.xn, v1.shade.yn, v1.shade.zn], Negate[edge] ];
};
ENDCASE => SIGNAL G3dRender.Error[$Unimplemented, "Unknown surface type"];
pt.shade.xn ← (v0.shade.xn + v1.shade.xn) / 2;
pt.shade.yn ← (v0.shade.yn + v1.shade.yn) / 2;
pt.shade.zn ← (v0.shade.zn + v1.shade.zn) / 2;
IF flat
THEN {         -- straight edge, average endpoints for midpoint
pt.coord.x ← (pos1.x + pos2.x) / 2;
pt.coord.y ← (pos1.y + pos2.y) / 2;
pt.coord.z ← (pos1.z + pos2.z) / 2;
}
ELSE {         -- not straight, evaluate curve
b0: Triple ← pos1; b3: Triple ← pos2;
b1: Triple ← Add[pos1, slope1];
b2: Triple ← Add[pos2, slope2];
m01, m12, m23, mm0, mm1: Triple;
m01.x ← (b0.x + b1.x)/2; m01.y ← (b0.y + b1.y)/2; m01.z ← (b0.z + b1.z)/2;
m12.x ← (b1.x + b2.x)/2; m12.y ← (b1.y + b2.y)/2; m12.z ← (b1.z + b2.z)/2;
m23.x ← (b2.x + b3.x)/2; m23.y ← (b2.y + b3.y)/2; m23.z ← (b2.z + b3.z)/2;
mm0.x ← (m01.x+m12.x)/2; mm0.y ← (m01.y+m12.y)/2; mm0.z ← (m01.z+m12.z)/2;
mm1.x ← (m23.x+m12.x)/2; mm1.y ← (m23.y+m12.y)/2; mm1.z ← (m23.z+m12.z)/2;
pt.coord.x ← (mm1.x + mm0.x) / 2;
pt.coord.y ← (mm1.y + mm0.y) / 2;
pt.coord.z ← (mm1.z + mm0.z) / 2;
t0.t0 ← [ slope1.x / 2, slope1.y / 2, slope1.z / 2 ];
t0.t1 ← [ mm0.x - pt.coord.x, mm0.y - pt.coord.y, mm0.z - pt.coord.z ];
t1.t0 ← [ mm1.x - pt.coord.x, mm1.y - pt.coord.y, mm1.z - pt.coord.z ];
t1.t1 ← [ slope2.x / 2, slope2.y / 2, slope2.z / 2 ];
};
IF lerpProc # NIL AND v0.aux # NIL THEN { -- get auxiliary info from supplied proc
data: LORALIST[ v0.aux, v1.aux, NEW[REAL ← 0.5], NEW[REAL ← 0.5] ];
pt ← lerpProc[ NIL, pt, data];     -- average, for lack of better strategy
};
};
pt.data ← v0.data
};
Procedures for tensor product quadrilaterals using tangent vectors at vertices
DisplayPatchTnsr: PatchProc ~ {
renderStyle: RenderStyle;
WITH patch.renderData.shadingClass.renderMethod SELECT FROM
style: REF RenderStyle => {
renderStyle ← style^;
SELECT renderStyle FROM
lines => { [] ← DisplayPatchEdges[context, patch]; RETURN[patch]; };
ENDCASE;
};
ENDCASE;
patch.type ← $PolygonToTnsrPatch;
IF data = NIL
THEN {
patchNo: NATNARROW[ GetProp[patch.props, $PatchNo], REF NAT ]^;
tangents: REF TangentSeqSeq ← NARROW[
GetProp[patch.renderData.props, $PatchTangents]
];
patch ← SubDivideTnsr[ context, patch, tangents[patchNo] ]; -- subdivide proc.
}
ELSE {
tangents: REF TangentSeq ← NARROW[ data ];
patch ← SubDivideTnsr[ context, patch, tangents ]; -- subdivide proc.
};
IF patch # NIL THEN G3dClipXfmShade.ReleasePatch[patch];    -- end of life for patch
RETURN[ NIL ];
};
SubDivideTnsr: PatchProc ~ {
Divide polygon into quadrilateral and triangular patches, depth sort, divide triangles into quadrilaterals then call QuadrilateralDisplay on each one
tol: REALIF context.antiAliasing THEN maxDeviation ELSE maxJaggyDeviation; 
tangents: REF TangentSeq ← NARROW[ data ];
IF patch.nVtces < 3 THEN SIGNAL G3dRender.Error[$MisMatch, "Not enough vertices"];
IF patch.nVtces = 3
THEN {
patch.props ← PutPropSafely[
patch.props, $Tangents,
NEW[TangentTriple ← [
NmlizeTangentSet[ patch[0], patch[1], tangents[0] ],
NmlizeTangentSet[ patch[1], patch[2], tangents[1] ],
NmlizeTangentSet[ patch[2], patch[0], tangents[2] ]
]]
];
TnsrTriangleDivide[context, patch, tol];
}
ELSE IF patch.nVtces > 4 THEN {
outPatch: PatchSequence ← NEW [PatchSequenceRep[patch.nVtces]];
midPt: CtlPtInfo ← GetCenterPt[context, patch, tangents];
midTangents: REF TangentSeq ← NARROW[ midPt.data ];
FOR i: NAT IN [0..patch.nVtces) DO
j: NAT ← (i + 1) MOD patch.nVtces;  -- next vertex
outPatch[i] ← GetPatch[3];  -- released by display action
outPatch[i].type ← patch.type;
outPatch[i].oneSided ← patch.oneSided;
outPatch[i].nVtces ← 3;
outPatch[i].clipState ← patch.clipState;
outPatch[i].dir ← unknown;
outPatch[i].renderData ← patch.renderData;
outPatch[i].props ← patch.props;
outPatch[i][0] ← patch[i];
outPatch[i][1] ← patch[(i+1) MOD patch.nVtces];
outPatch[i][2] ← midPt;
outPatch[i].props ← PutPropSafely[
outPatch[i].props, $Tangents,
NEW[TangentTriple ← [
NmlizeTangentSet[ outPatch[i][0], outPatch[i][1], tangents[i] ], -- orig. edge
midTangents[j],         -- in to middle, from next vertex
[ t0: midTangents[i].t1, et0: midTangents[i].et1,  -- back out to current vertex
t1: midTangents[i].t0, et1: midTangents[i].et0 ]
]]
];
IF outPatch[i].clipState # in THEN G3dClipXfmShade.GetPatchClipState[ outPatch[i] ];
ENDLOOP;
outPatch.length ← patch.nVtces;
outPatch ← PatchDepthSort[ context, outPatch ];    -- sort to depth order
FOR i: NAT IN [0..patch.nVtces) DO      -- subdivide and display
TnsrTriangleDivide[context, outPatch[i], tol];  ReleasePatch[outPatch[i]];
ENDLOOP;
}
ELSE {
patch.props ← PutPropSafely[
patch.props, $Tangents,
NEW[TangentQuad ← [
NmlizeTangentSet[ patch[0], patch[1], tangents[0] ],
NmlizeTangentSet[ patch[1], patch[2], tangents[1] ],
NmlizeTangentSet[ patch[2], patch[3], tangents[2] ],
NmlizeTangentSet[ patch[3], patch[0], tangents[3] ]
]]
];
TnsrQuadDisplay[context, patch, 0, tol];
};
RETURN[patch];
};
TnsrQuadDisplay: PROC[ context: Context, p: REF Patch, level: NAT, tol: REAL] ~ {
Recursively divide quadrilateral patches to straight edges, then display
subP: PatchSequence ← NIL;
IF context.stopMe^ THEN RETURN[];  -- shut down if stop signal received
IF p.dir = unknown THEN [] ← TnsrQuadBackFacing[p];  -- backface test
IF ( p.clipState # out ) AND ( NOT p.oneSided OR NOT p.dir = back ) THEN {
IF level < recurseLimit THEN subP ← SubdivideTnsrQuad[context, p, level, tol];
IF subP = NIL
THEN {          -- recursion limit hit or edges all straight
IF showLines THEN p.type ← $PolyLine ELSE p.type ← $ConvexPolygon;
G3dClipXfmShade.ShadePoly[context, p];
[p] ← G3dSortandDisplay.OutputPolygon[context, p];      -- display
}
ELSE {
subP ← PatchDepthSort[ context, subP ]; -- sort to display order
TnsrQuadDisplay[context, subP[0], level+1, tol];  ReleasePatch[subP[0]];
TnsrQuadDisplay[context, subP[1], level+1, tol];  ReleasePatch[subP[1]];
TnsrQuadDisplay[context, subP[2], level+1, tol];  ReleasePatch[subP[2]];
TnsrQuadDisplay[context, subP[3], level+1, tol];  ReleasePatch[subP[3]];
};
};
};
TnsrTriangleDivide: PROC[ context: Context, p: REF Patch, tol: REAL] ~ {
Divide triangular patch into three subquadrilaterals and send to display
shape: Shape ← NARROW[ GetProp[p.props, $Shape] ];
v0, v1, v2, vCtr, vCtr1, vCtr2, vCtr3: CtlPtInfo;  flat0, flat1, flat2: BOOLEAN;
t: REF TangentTriple ← NARROW[GetProp[p.props, $Tangents] ];
t00, t01, t10, t11, t20, t21, tm0, tm1, tm2: TangentSet;
outPatch: PatchSequence ← NEW[PatchSequenceRep[3]];
Find midpoints and midslopes for each side
[v0, t00, t01, flat0] ← CurveDivideTan[
context, p[0], p[1], t[0].et0, t[0].et1,
GetNmlVec[t[0].et0, p[0]], GetNmlVec[t[0].et1, p[1], TRUE], 0.5, tol ];
[v1, t10, t11, flat1] ← CurveDivideTan[
context, p[1], p[2], t[1].et0, t[1].et1,
GetNmlVec[t[1].et0, p[1]], GetNmlVec[t[1].et1, p[2], TRUE], 0.5, tol];
[v2, t20, t21, flat2] ← CurveDivideTan[
context, p[2], p[0], t[2].et0, t[2].et1,
GetNmlVec[t[2].et0, p[2]], GetNmlVec[t[2].et1, p[0], TRUE], 0.5, tol];
IF flat0 AND flat1 AND flat2 AND stopIfStraight     -- all straight? then done
THEN TnsrQuadDisplay[context, p, recurseLimit, tol];     -- display as polygon
Get inner edges from midpoint to opposite vertex curves
tm0.et0 ← GetSlopeVec[ [v0.shade.exn, v0.shade.eyn, v0.shade.ezn], Add[t10.et0, t21.et1] ];
tm0.et0 ← ScaleTangent[ tm0.et0, DiffPosnsCtlPt[p[2].coord, v0.coord, $Eye] ];
tm0.et1 ← GetSlopeVec[ [p[2].shade.exn, p[2].shade.eyn, p[2].shade.ezn], Add[t20.et0, t11.et1]];
tm0.et1 ← ScaleTangent[ tm0.et1, DiffPosnsCtlPt[v0.coord, p[2].coord, $Eye] ];
tm1.et0 ← GetSlopeVec[ [v1.shade.exn, v1.shade.eyn, v1.shade.ezn], Add[t20.et0, t01.et1]];
tm1.et0 ← ScaleTangent[ tm1.et0, DiffPosnsCtlPt[p[0].coord, v1.coord, $Eye] ];
tm1.et1 ← GetSlopeVec[ [p[0].shade.exn, p[0].shade.eyn, p[0].shade.ezn], Add[t00.et0, t21.et1]];
tm1.et1 ← ScaleTangent[ tm1.et1, DiffPosnsCtlPt[v1.coord, p[0].coord, $Eye] ];
tm2.et0 ← GetSlopeVec[ [v2.shade.exn, v2.shade.eyn, v2.shade.ezn], Add[t00.et0, t11.et1]];
tm2.et0 ← ScaleTangent[ tm2.et0, DiffPosnsCtlPt[p[1].coord, v2.coord, $Eye] ];
tm2.et1 ← GetSlopeVec[ [p[1].shade.exn, p[1].shade.eyn, p[1].shade.ezn], Add[t10.et0, t01.et1]];
tm2.et1 ← ScaleTangent[ tm2.et1, DiffPosnsCtlPt[v2.coord, p[1].coord, $Eye] ];
Split inner edges at 1/3 point and average for ctr. point
[vCtr1, tm0, , ] ← CurveDivideTan[
context, v0, p[2], tm0.et0, tm0.et1,
GetNmlVec[tm0.et0, v0], GetNmlVec[tm0.et1, p[2], TRUE], 1.0/3.0, tol ];
[vCtr2, tm1, , ] ← CurveDivideTan[
context, v1, p[0], tm1.et0, tm1.et1,
GetNmlVec[tm1.et0, v1], GetNmlVec[tm1.et1, p[0], TRUE], 1.0/3.0, tol ];
[vCtr3, tm2, , ] ← CurveDivideTan[
context, v2, p[1], tm2.et0, tm2.et1,
GetNmlVec[tm2.et0, v2], GetNmlVec[tm2.et1, p[1], TRUE], 1.0/3.0, tol ];
[[vCtr.coord.ex, vCtr.coord.ey, vCtr.coord.ez]] ← Div[     -- average positions
Add[
[vCtr1.coord.ex, vCtr1.coord.ey, vCtr1.coord.ez],
Add[
[vCtr2.coord.ex, vCtr2.coord.ey, vCtr2.coord.ez],
[vCtr3.coord.ex, vCtr3.coord.ey, vCtr3.coord.ez]
]
],
3.0
];
[[vCtr.shade.exn, vCtr.shade.eyn, vCtr.shade.ezn]] ← Div[    -- average normals
Add[
[vCtr1.shade.exn, vCtr1.shade.eyn, vCtr1.shade.ezn],
Add[
[vCtr2.shade.exn, vCtr2.shade.eyn, vCtr2.shade.ezn],
[vCtr3.shade.exn, vCtr3.shade.eyn, vCtr3.shade.ezn]
]
],
3.0
];
FOR i: NAT IN [0..3) DO
outPatch[i] ← GetPatch[4];     -- allocate 4 point patch
outPatch[i].type ← p.type;
outPatch[i].oneSided ← p.oneSided;
outPatch[i].nVtces ← 4;
outPatch[i].clipState ← p.clipState;
outPatch[i].dir ← p.dir;
outPatch[i].renderData ← p.renderData;
outPatch[i].props ← p.props;
ENDLOOP;
{ OPEN v0.coord; clip ← G3dClipXfmShade.GetClipCodeForPt[context, [ ex, ey, ez] ];
IF clip = NoneOut
THEN [ [sx, sy, sz] ] ← G3dClipXfmShade.XfmPtToDisplay[ context, [ex, ey, ez], shape ]; };
{ OPEN v1.coord; clip ← G3dClipXfmShade.GetClipCodeForPt[context, [ ex, ey, ez] ];
IF clip = NoneOut
THEN [ [sx, sy, sz] ] ← G3dClipXfmShade.XfmPtToDisplay[ context, [ex, ey, ez], shape ]; };
{ OPEN v2.coord; clip ← G3dClipXfmShade.GetClipCodeForPt[context, [ ex, ey, ez] ];
IF clip = NoneOut
THEN [ [sx, sy, sz] ] ← G3dClipXfmShade.XfmPtToDisplay[ context, [ex, ey, ez], shape ]; };
{ OPEN vCtr.coord; clip ← G3dClipXfmShade.GetClipCodeForPt[context, [ ex, ey, ez] ];
IF clip = NoneOut
THEN [ [sx, sy, sz] ] ← G3dClipXfmShade.XfmPtToDisplay[ context, [ex, ey, ez], shape ]; };
outPatch[0][0] ← p[0]; outPatch[0][1] ← v0; outPatch[0][2] ← vCtr; outPatch[0][3] ← v2;
outPatch[1][0] ← p[1]; outPatch[1][1] ← v1; outPatch[1][2] ← vCtr; outPatch[1][3] ← v0;
outPatch[2][0] ← p[2]; outPatch[2][1] ← v2; outPatch[2][2] ← vCtr; outPatch[2][3] ← v1;
outPatch[0].props ← PutPropSafely[ outPatch[0].props, $Tangents, NEW[TangentQuad ←
[t00, tm0, [tm2.t1, tm2.et1, tm2.t0, tm2.et0], t21] ] ];
outPatch[1].props ← PutPropSafely[ outPatch[1].props, $Tangents, NEW[TangentQuad ←
[t10, tm1, [tm0.t1, tm0.et1, tm0.t0, tm0.et0], t01] ] ];
outPatch[2].props ← PutPropSafely[ outPatch[2].props, $Tangents, NEW[TangentQuad ←
[t20, tm2, [tm1.t1, tm1.et1, tm1.t0, tm1.et0], t11] ] ];
FOR i: NAT IN [0..3) DO G3dClipXfmShade.GetPatchClipState[ outPatch[i] ]; ENDLOOP; --bad!!
outPatch.length ← 3;
FOR i: NAT IN [0..3) DO
TnsrQuadDisplay[ context, outPatch[i], 0, tol ];  ReleasePatch[outPatch[i]];
ENDLOOP;
};
SubdivideTnsrQuad: PROC[context: Context, p: REF Patch, level: NAT, tol: REAL]
       RETURNS[PatchSequence] ~ {
Divide quadrangular patch into four subquadrangles
shape: Shape ← NARROW[ GetProp[p.props, $Shape] ];
v0, v1, v2, v3, vCtr, vCtr2: CtlPtInfo;  flat0, flat1, flat2, flat3: BOOLEAN;
t: REF TangentQuad ← NARROW[GetProp[p.props, $Tangents] ];
t00, t01, t10, t11, t20, t21, t30, t31, tm0, tm1, tm2, tm3: TangentSet;
outPatch: PatchSequence ← NEW[PatchSequenceRep[4]];
Find midpoints and midslopes for each side
[v0, t00, t01, flat0] ← CurveDivideTan[
context, p[0], p[1], t[0].et0, t[0].et1,
GetNmlVec[t[0].et0, p[0]], GetNmlVec[t[0].et1, p[1], TRUE], 0.5, tol ];
[v1, t10, t11, flat1] ← CurveDivideTan[
context, p[1], p[2], t[1].et0, t[1].et1,
GetNmlVec[t[1].et0, p[1]], GetNmlVec[t[1].et1, p[2], TRUE], 0.5, tol];
[v2, t20, t21, flat2] ← CurveDivideTan[
context, p[2], p[3], t[2].et0, t[2].et1,
GetNmlVec[t[2].et0, p[2]], GetNmlVec[t[2].et1, p[3], TRUE], 0.5, tol];
[v3, t30, t31, flat3] ← CurveDivideTan[
context, p[3], p[0], t[3].et0, t[3].et1,
GetNmlVec[t[3].et0, p[3]], GetNmlVec[t[3].et1, p[0], TRUE], 0.5, tol];
IF flat0 AND flat1 AND flat2 AND flat3 AND stopIfStraight  -- all straight? then done
THEN RETURN[NIL];
Get inner edge endpoint tangents from opposite midpoints
- first project direction given by endpoint cross-tangents on plane given by normal
- Then scale direction by distance to opposite midpoint
tm0.et0 ← GetSlopeVec[ [v0.shade.exn, v0.shade.eyn, v0.shade.ezn], Add[t10.et0, t31.et1] ];
tm0.et0 ← ScaleTangent[ tm0.et0, DiffPosnsCtlPt[v2.coord, v0.coord, $Eye] ];
tm0.et1 ← GetSlopeVec[ [v2.shade.exn, v2.shade.eyn, v2.shade.ezn], Add[t11.et1, t30.et0] ];
tm0.et1 ← ScaleTangent[ tm0.et1, DiffPosnsCtlPt[v0.coord, v2.coord, $Eye] ];
tm1.et0 ← GetSlopeVec[ [v1.shade.exn, v1.shade.eyn, v1.shade.ezn], Add[t20.et0, t01.et1] ];
tm1.et0 ← ScaleTangent[ tm1.et0, DiffPosnsCtlPt[v3.coord, v1.coord, $Eye] ];
tm1.et1 ← GetSlopeVec[ [v3.shade.exn, v3.shade.eyn, v3.shade.ezn], Add[t21.et1, t00.et0] ];
tm1.et1 ← ScaleTangent[ tm1.et1, DiffPosnsCtlPt[v1.coord, v3.coord, $Eye] ];
tm0.et0 ← GetSlopeVec[ [v0.shade.exn, v0.shade.eyn, v0.shade.ezn],
       DiffPosnsCtlPt[v2.coord, v0.coord, $Eye] ];
tm0.et1 ← GetSlopeVec[ [v2.shade.exn, v2.shade.eyn, v2.shade.ezn],
       DiffPosnsCtlPt[v0.coord, v2.coord, $Eye] ];
tm1.et0 ← GetSlopeVec[ [v1.shade.exn, v1.shade.eyn, v1.shade.ezn],
       DiffPosnsCtlPt[v3.coord, v1.coord, $Eye] ];
tm1.et1 ← GetSlopeVec[ [v3.shade.exn, v3.shade.eyn, v3.shade.ezn],
       DiffPosnsCtlPt[v1.coord, v3.coord, $Eye] ];
Get center point, normal and cross tangents
[vCtr, tm0, tm2, ] ← CurveDivideTan[
context, v0, v2, tm0.et0, tm0.et1,
GetNmlVec[tm0.et0, v0], GetNmlVec[tm0.et1, v2, TRUE], 0.5, tol ];
[vCtr2, tm1, tm3, ] ← CurveDivideTan[
context, v1, v3, tm1.et0, tm1.et1,
GetNmlVec[tm1.et0, v1], GetNmlVec[tm1.et1, v3, TRUE], 0.5, tol ];
[[vCtr.coord.ex, vCtr.coord.ey, vCtr.coord.ez]] ← Div[       -- average
Add[
[vCtr.coord.ex, vCtr.coord.ey, vCtr.coord.ez],
[vCtr2.coord.ex, vCtr2.coord.ey, vCtr2.coord.ez]
],
2
];
[[vCtr.shade.exn, vCtr.shade.eyn, vCtr.shade.ezn]] ← Nmlize[ Cross[ tm2.et0 , tm3.et0 ] ];
FOR i: NAT IN [0..4) DO
outPatch[i] ← GetPatch[4]; -- 4 point patch, released by display action
outPatch[i].type ← p.type;
outPatch[i].oneSided ← p.oneSided;
outPatch[i].nVtces ← 4;
outPatch[i].clipState ← p.clipState;
outPatch[i].dir ← p.dir;
outPatch[i].renderData ← p.renderData;
outPatch[i].props ← p.props;
ENDLOOP;
{ OPEN v0.coord; clip ← G3dClipXfmShade.GetClipCodeForPt[context, [ ex, ey, ez] ];
IF clip = NoneOut
THEN [ [sx, sy, sz] ] ← G3dClipXfmShade.XfmPtToDisplay[ context, [ex, ey, ez], shape ]; };
{ OPEN v1.coord; clip ← G3dClipXfmShade.GetClipCodeForPt[context, [ ex, ey, ez] ];
IF clip = NoneOut
THEN [ [sx, sy, sz] ] ← G3dClipXfmShade.XfmPtToDisplay[ context, [ex, ey, ez], shape ]; };
{ OPEN v2.coord; clip ← G3dClipXfmShade.GetClipCodeForPt[context, [ ex, ey, ez] ];
IF clip = NoneOut
THEN [ [sx, sy, sz] ] ← G3dClipXfmShade.XfmPtToDisplay[ context, [ex, ey, ez], shape ]; };
{ OPEN v3.coord; clip ← G3dClipXfmShade.GetClipCodeForPt[context, [ ex, ey, ez] ];
IF clip = NoneOut
THEN [ [sx, sy, sz] ] ← G3dClipXfmShade.XfmPtToDisplay[ context, [ex, ey, ez], shape ]; };
{ OPEN vCtr.coord; clip ← G3dClipXfmShade.GetClipCodeForPt[context, [ ex, ey, ez] ];
IF clip = NoneOut
THEN [ [sx, sy, sz] ] ← G3dClipXfmShade.XfmPtToDisplay[ context, [ex, ey, ez], shape ]; };
outPatch[0][0] ← p[0]; outPatch[0][1] ← v0; outPatch[0][2] ← vCtr; outPatch[0][3] ← v3;
outPatch[1][0] ← p[1]; outPatch[1][1] ← v1; outPatch[1][2] ← vCtr; outPatch[1][3] ← v0;
outPatch[2][0] ← p[2]; outPatch[2][1] ← v2; outPatch[2][2] ← vCtr; outPatch[2][3] ← v1;
outPatch[3][0] ← p[3]; outPatch[3][1] ← v3; outPatch[3][2] ← vCtr; outPatch[3][3] ← v2;
IF t # NIL THEN {
outPatch[0].props ← PutPropSafely[ outPatch[0].props, $Tangents, NEW[TangentQuad ←
[t00, tm0, tm3, t31] ] ];
outPatch[1].props ← PutPropSafely[ outPatch[1].props, $Tangents, NEW[TangentQuad ←
[t10, tm1, [tm0.t1, tm0.et1, tm0.t0, tm0.et0], t01] ] ];
outPatch[2].props ← PutPropSafely[ outPatch[2].props, $Tangents, NEW[TangentQuad ←
[t20, [tm2.t1, tm2.et1, tm2.t0, tm2.et0], [tm1.t1, tm1.et1, tm1.t0, tm1.et0], t11] ] ];
outPatch[3].props ← PutPropSafely[ outPatch[3].props, $Tangents, NEW[TangentQuad ←
[t30, [tm3.t1, tm3.et1, tm3.t0, tm3.et0], tm2, t21] ] ];
};
FOR i: NAT IN [0..4) DO G3dClipXfmShade.GetPatchClipState[ outPatch[i] ]; ENDLOOP; --bad!!
outPatch.length ← 4;
RETURN[ outPatch ]; -- return four sub-patches
};
TnsrQuadBackFacing: PROC[p: REF Patch] RETURNS [BOOLEAN] ~ {
p.dir ← unknown;        -- don't know how to do this yet
IF p.dir = back THEN RETURN[TRUE] ELSE RETURN[FALSE];
};
Procedures for computing tangent vectors at vertices
GetTangents: ShapeProc ~ {
Build endpoint tangents for each edge of a polygon
GetSlope: PROC[normal, edge1, edge2: Triple, reverse: BOOL] RETURNS[slope: Triple] ~ {
Returns a vector normal to the 1st vector and in the plane defined by both vectors
cornerNorm: Triple ← Cross[edge1, edge2];
IF reverse THEN cornerNorm ← Negate[ cornerNorm ];
slope ← Cross[ Cross[normal, edge2], normal ];
IF Length[slope] < .0001 THEN IF Dot[normal, edge2] < 0.0  -- catch degenerate cases
THEN slope ← cornerNorm ELSE slope ← Negate[ cornerNorm ];
IF Dot[ normal, cornerNorm ] < 0.0 THEN
IF Dot[ Nmlize[edge2], Nmlize[normal] ] > Dot[ Nmlize[edge1], Nmlize[normal] ]
THEN slope ← Negate[ slope ]; -- wrong side of plane, high-curvature on one side
SIGNAL G3dRender.Error[$Unimplemented, "surface ambiguously curved"];
slope ← ScaleTangent[ slope, edge2 ];
};
nullTriple: Triple ← [0.0, 0.0, 0.0];
corners: REF CornerSeqSeq;
Pass through polygon structure, getting vertex normals robustly while building corner structure.
corners ← GetNormalsAndDirections[shape]; -- Get robust vertex normals and edge directions
Order contiguous polygons to clockwise around vertex, seen from outside
FOR i: NAT IN [0..corners.length) DO-- sort corners to clockwise order by matching edges
IF corners[i] # NIL THEN corners[i] ← SortCtlPoint[ corners[i] ];
ENDLOOP;
Build tangents at corners around each vertex. By projecting edge directions on plane given by normal
FOR i: NAT IN [0..corners.length) DO
nCorners: NATIF corners[i] # NIL THEN corners[i].length ELSE 0;
outDirs: TripleSequence ← NEW[TripleSequenceRep[nCorners]]; -- temp tangents
inDirs: TripleSequence ← NEW[TripleSequenceRep[nCorners]];
FOR this: NAT IN [0..nCorners) DO
last: NAT ← (this + nCorners-1) MOD nCorners;
outDirs[this] ← GetSlope[
corners[i][this].normal, corners[i][this].inDir, corners[i][this].outDir,
corners[i][this].concave
];
IF corners[i][last].inVtx = corners[i][this].outVtx   -- matches previous corner
THEN inDirs[last] ← outDirs[this]      -- copy vector
ELSE inDirs[last] ← GetSlope[       -- open edge, compute vector
corners[i][last].normal, corners[i][last].outDir, corners[i][last].inDir,
NOT corners[i][this].concave
];
ENDLOOP;
FOR this: NAT IN [0..nCorners) DO     -- store new tangents
corners[i][this].outDiroutDirs[this];
corners[i][this].inDirinDirs[this];
ENDLOOP;
ENDLOOP;
Make tangents continuous if within 30 degrees of being so. Always make continuous if open edge
FOR i: NAT IN [0..corners.length) DO
IF corners[i] # NIL THEN corners[i] ← FindContinuousEdges[ shape, i, corners[i] ];
ENDLOOP;
{ -- build tangent structure for polygons, using the vertex tangent structure just constructed
renderData: REF RenderData ← G3dRender.RenderDataFrom[shape];
tangents: REF TangentSeqSeq ← NEW[ TangentSeqSeq[shape.surfaces.length] ];
FOR i: NAT IN [0..shape.surfaces.length) DO
nVtces: NAT ← shape.surfaces[i].vertices.length;
tangents[i] ← NEW[ TangentSeq[nVtces] ];
FOR cVtx: NAT IN [0..nVtces) DO-- get direction vectors for each edge at vtx.
vtx: NAT ← shape.surfaces[i].vertices[cVtx];       -- current vertex
nVtx: NAT ← shape.surfaces[i].vertices[(cVtx + 1) MOD nVtces]; -- other end of edge
Find outgoing direction vector by pointer to opposing vertex
IF corners[vtx][corners[vtx].length-1].inVtx # corners[vtx][0].outVtx -- open edge?
 AND corners[vtx][corners[vtx].length-1].inVtx = nVtx -- check last inDir
THEN tangents[i][cVtx].t0 ← corners[vtx][corners[vtx].length-1].inDir
ELSE FOR j: NAT IN [0..corners[vtx].length) DO-- otherwise check all outDirs
IF corners[vtx][j].outVtx = nVtx THEN
{ tangents[i][cVtx].t0 ← corners[vtx][j].outDir; EXIT; };
ENDLOOP;
Find direction vector at other end of edge by looking back from opposing vertex
IF corners[nVtx][corners[nVtx].length-1].inVtx # corners[nVtx][0].outVtx -- open?
 AND corners[nVtx][corners[nVtx].length-1].inVtx = vtx -- last inDir
THEN tangents[i][cVtx].t1 ← corners[nVtx][corners[nVtx].length-1].inDir
ELSE FOR j: NAT IN [0..corners[nVtx].length) DO -- otherwise check all outDirs
IF corners[nVtx][j].outVtx = vtx THEN
{ tangents[i][cVtx].t1 ← corners[nVtx][j].outDir; EXIT; };
ENDLOOP;
IF tangents[i][cVtx].t0 = nullTriple OR tangents[i][cVtx].t1 = nullTriple THEN
SIGNAL G3dRender.Error[$Fatal, "Null tangent vector, bad topology?"];
ENDLOOP;
tangents[i].length ← nVtces;
ENDLOOP;
tangents.length ← shape.surfaces.length;
renderData.props ← PutProp[renderData.props, $PatchTangents, tangents];
};
RETURN[shape];
};
GetNormalsAndDirections: PROC[shape: Shape] RETURNS[REF CornerSeqSeq] ~ {
Get vertex normals robustly ( using robust polygon normal )
Find incoming and outgoing directions for each polygon at a vertex
corners: REF CornerSeqSeq ← NEW[ CornerSeqSeq[shape.vertices.length] ];
Clear vertex normals
IF NOT shape.vertices.valid.normal THEN FOR i: NAT IN [0..shape.vertices.length) DO
shape.vertices[i].normal ← nullTriple;
ENDLOOP;
FOR i: NAT IN [0..shape.surfaces.length) DO
nVtces: NAT ← shape.surfaces[i].vertices.length;
normal: Triple ← nullTriple;
concave: REF BoolSequence ← NEW[ BoolSequence[nVtces] ];
cNmls: TripleSequence ← NEW[ TripleSequenceRep[nVtces] ];
Get polygon normal robustly ( allowing slightly concave polygons ). Sum of normals given at corners is assumed to be the dominant direction. Normals over 90 degrees off that are assumed to come from concave vertices and are therefore reversed.
FOR cVtx: NAT IN [0..nVtces) DO -- get normal for each corner
vtx: NAT ← shape.surfaces[i].vertices[cVtx];
nVtx: NAT ← shape.surfaces[i].vertices[(cVtx + 1) MOD nVtces];
lVtx: NAT ← shape.surfaces[i].vertices[(cVtx + nVtces - 1) MOD nVtces];
cNmls[cVtx] ← Cross[  -- in object space so do right-handed
DiffPosnsVtx[ shape.vertices[lVtx], shape.vertices[vtx] ],
DiffPosnsVtx[ shape.vertices[nVtx], shape.vertices[vtx] ]
];
IF unitNormals THEN cNmls[cVtx] ← Nmlize[cNmls[cVtx]];   
normal ← Add[ normal, cNmls[cVtx] ];  -- sum normals
ENDLOOP;
FOR cVtx: NAT IN [0..nVtces) DO -- Get normals robustly
IF Dot[ cNmls[cVtx], normal ] < 0.0
THEN {         -- more than 90 degrees off poly norm.
cNmls[cVtx] ← Negate[ cNmls[cVtx] ];
normal ← Add[ normal, cNmls[cVtx] ];  -- cancel 1st entry
normal ← Add[ normal, cNmls[cVtx] ]; -- contribute to sum
concave[cVtx] ← TRUE;    -- put concave tag in here
}
ELSE concave[cVtx] ← FALSE;
ENDLOOP;
FOR cVtx: NAT IN [0..nVtces) DO
vtx: NAT ← shape.surfaces[i].vertices[cVtx];
IF shape.vertices[vtx] # NIL THEN {
OPEN shape.vertices[vtx];
IF shape.vertices.valid.normal
THEN cNmls[cVtx] ← normal      -- Get predefined normal
ELSE normal ← Add[normal, cNmls[cVtx]];  -- sum normals
};
ENDLOOP;
Build corner structure
FOR cVtx: NAT IN [0..nVtces) DO-- get direction vectors for each edge at vtx.
vtx: NAT ← shape.surfaces[i].vertices[cVtx];
nVtx: NAT ← shape.surfaces[i].vertices[(cVtx + 1) MOD nVtces];
lVtx: NAT ← shape.surfaces[i].vertices[(cVtx + nVtces - 1) MOD nVtces];
corners[vtx] ← UpdateCornerVecSeq[
corners[vtx],
[ lVtx, nVtx,     -- store number of incoming vertex and outgoing vertex
DiffPosnsVtx[shape.vertices[lVtx], shape.vertices[vtx]], -- dir to incoming vtx
DiffPosnsVtx[shape.vertices[nVtx], shape.vertices[vtx]],  -- outgoing direction
cNmls[cVtx],              -- normal at corner
nullTriple,              -- interior knot
concave[cVtx]             -- concave tag
]
];
ENDLOOP;
ENDLOOP;
corners.length ← shape.vertices.length;
FOR i: NAT IN [0..shape.vertices.length) DO -- unit new vertex normals, store in corners
OPEN shape.vertices[i];
IF Length[ normal ] > shape.sphereExtent.radius * .0001
THEN normal ← Nmlize[ normal ]
ELSE normal ← nullTriple;  -- likely unused, set default to stop trouble
IF corners[i] # NIL THEN FOR j: NAT IN [0..corners[i].length) DO
corners[i][j].normal ← normal;
ENDLOOP;
ENDLOOP;
RETURN[ corners ];
};
UpdateCornerVecSeq: PROC[corner: REF CornerSeq, entry: Corner]
       RETURNS
[REF CornerSeq] ~ {
newSeq: REF CornerSeq;
IF corner = NIL
THEN newSeq ← NEW[ CornerSeq[6] ]    -- get brand new sequence
ELSE IF corner.length = corner.maxLength
THEN {           -- expand filled sequence
newSeq ← NEW[ CornerSeq[ corner.maxLength * 2 ] ];
FOR i: NAT IN [0..corner.maxLength) DO
newSeq[i] ← corner[i];    -- copy into new sequence
ENDLOOP;
}
ELSE newSeq ← corner;      -- no changes required
corner ← newSeq;
corner[corner.length] ← entry;
corner.length ← corner.length + 1;
RETURN[corner];
};
FindContinuousEdges: PROC[ shape: Shape, vtx: NAT, corner: REF CornerSeq ]
        RETURNS[ REF CornerSeq ] ~ {
Make tangents continuous if within 30 degrees of being so. Always make continuous if open edge, adjacent edges not allowed to be
nCorners: NAT ← corner.length;
acrossFrom: REF NatSequenceSequence ← NEW[NatSequenceSequence[nCorners]];
FOR this: NAT IN [0..nCorners) DO
acrossFrom[this] ← NEW[NatSequenceRep[nCorners]];
ENDLOOP;
Check for open edge, if there, align tangents unless included angle less than 45 degrees
IF corner[nCorners-1].inVtx # corner[0].outVtx
THEN { -- open edge (should only be one)
cosAng: REAL ← Dot[Nmlize[corner[nCorners-1].inDir], Nmlize[ corner[0].outDir]];
IF cosAng < minCosToAlignOpen THEN { -- included angle over limit (aligned => -1.0)
IF unitNormals THEN {   -- unit to give equal weight to direction
corner[nCorners-1].inDir ← Nmlize[corner[nCorners-1].inDir];
corner[0].outDir    ← Nmlize[corner[0].outDir];
};
corner[nCorners-1].inDir ← Add[ -- sum to get aligned direction
corner[nCorners-1].inDir, Negate[corner[0].outDir]
];
corner[0].outDir ← Negate[corner[nCorners-1].inDir]; -- opposite dir
corner[nCorners-1].inDir ← ScaleTangent[     -- scale inDir to Bezier length
corner[nCorners-1].inDir,            -- outdir done below
DiffPosnsVtx[shape.vertices[corner[nCorners-1].inVtx], shape.vertices[vtx]]
];
};
}
ELSE SELECT nCorners FROM  -- no open edge
1, 2 => SIGNAL G3dRender.Error[$MisMatch, "too few edges at vertex"];
3 => {};        -- leave things alone if 3
4 => FOR i: NAT IN [0..2) DO        -- align opposing vertices if 4
o: NAT ← i+2; o1: NAT ← i+1; i1: NAT ← (i+3) MOD 4; -- opposite, next, and prev.
corner[o].outDir ← corner[o1].inDir ← Add[ -- sum for aligned dir
corner[o].outDir, Negate[corner[i].outDir]
];
corner[i].outDir ← corner[i1].inDir ← Negate[corner[o].outDir];
ENDLOOP;
ENDCASE => {
};
This will attempt to align edges in the more complex cases
IF nCorners > 50 THEN {    -- for corners over 4
FOR this: NAT IN [0..nCorners) DO  -- tag nearly continuous pairs
FOR other: NAT IN (this+1..nCorners) DO -- don't test adjacent edges, musn't be aligned
cosAng: REAL ← Dot[Nmlize[corner[other].outDir], Nmlize[corner[this].outDir]];
IF this = 0 AND other = nCorners-1 AND corner[other].inVtx = corner[this].outVtx
THEN EXIT;       -- skip if adjacent edge wrapping around zero
IF cosAng < minCosToAlign THEN { -- included angle over limit, mark both edges
acrossFrom[this][acrossFrom[this].length] ← other;
acrossFrom[this].length ← acrossFrom[this].length + 1;
acrossFrom[other][acrossFrom[other].length] ← this;
acrossFrom[other].length ← acrossFrom[other].length + 1;
};
ENDLOOP;
ENDLOOP;
FOR i: NAT IN [0..nCorners) DO          -- align paired tangents
Sum: PROC[ k: NAT ] ~ {
FOR j: NAT IN [1..acrossFrom[k].length) DO-- sum tangents across from this one
corner[acrossFrom[k][0]].outDir ← Add[
corner[acrossFrom[k][0]].outDir, corner[acrossFrom[k][j]].outDir
];
ENDLOOP;
};
Copy: PROC[ k: NAT ] ~ {
FOR j: NAT IN [1..acrossFrom[k].length) DO    -- copy sum to other tangents
corner[acrossFrom[k][j]].outDir ← corner[acrossFrom[k][0]].outDir;
ENDLOOP;
};
Kill: PROC[ k, other: NAT ] ~ {
FOR j: NAT IN [1..acrossFrom[k].length) DO-- kill pointers to prevent further use
IF j # other THEN acrossFrom[acrossFrom[k][j]].length ← 0;
ENDLOOP;
acrossFrom[k].length ← 0;
};
this, last, bigOne: NAT ← i;
biggest: BOOLEANFALSE;
WHILE NOT biggest DO     -- chain through to largest set of paired tangents
biggestFound: NAT ← acrossFrom[this].length;
biggest ← TRUE; bigOne ← this;
FOR j: NAT IN [0..acrossFrom[this].length) DO -- check all edges paired with this one
IF acrossFrom[acrossFrom[this][j]].length >= biggestFound THEN {
bigOne ← acrossFrom[this][j];     -- found one the same size at least
IF acrossFrom[acrossFrom[this][j]].length > biggestFound THEN {-- a bigger one
biggestFound ← acrossFrom[acrossFrom[this][j]].length; biggest ← FALSE;
};
};
ENDLOOP;
last ← this; this ← bigOne;
ENDLOOP;
IF acrossFrom[bigOne].length > 1 AND acrossFrom[last].length > 1 THEN {
FOR j: NAT IN [0..acrossFrom[bigOne].length) DO  -- cull any neighboring edges
IF (acrossFrom[bigOne][j] = (last + 1) MOD nCorners)
OR (acrossFrom[bigOne][j] = (last + nCorners-1) MOD nCorners)
THEN {
FOR k: NAT IN [j..acrossFrom[bigOne].length) DO
acrossFrom[bigOne][k] ← acrossFrom[bigOne][k+1];
ENDLOOP;
acrossFrom[bigOne].length ← acrossFrom[bigOne].length - 1;
};
ENDLOOP;
};
Sum[bigOne]; Sum[last];   -- sum all directions in set
IF bigOne # last THEN {
corner[acrossFrom[last][0]].outDir ← Add[ 
corner[acrossFrom[last][0]].outDir,         -- get aligned direction
Negate[corner[acrossFrom[bigOne][0]].outDir]
];
corner[acrossFrom[bigOne][0]].outDir ← Negate[ corner[acrossFrom[last][0]].outDir ];
};
Copy[bigOne]; Copy[last];       -- copy directions to others in sets
Kill[bigOne, last]; Kill[last, bigOne];        -- kill sets
ENDLOOP;
};
FOR this: NAT IN [0..nCorners) DO    -- go back over tangents and scale properly
last: NAT ← (this + nCorners -1) MOD nCorners;
corner[this].outDir ← ScaleTangent[
corner[this].outDir,
DiffPosnsVtx[shape.vertices[corner[this].outVtx], shape.vertices[vtx]]
];
IF corner[last].inVtx = corner[this].outVtx THEN corner[last].inDir ← corner[this].outDir;
ENDLOOP;
RETURN[ corner ];
};
SortCtlPoint: PROC[corner: REF CornerSeq] RETURNS[REF CornerSeq] ~ {
Sort polygon entries to clockwise order at each vertex, ensure discontinuities are at extremes
Build chain of corners by checking both ends of chain against all remaining unattached corners, only one chain allowed, since only one open edge allowed at vertex
nCrnrs: NAT ← corner.length;
someLeft, someFound: BOOLEANTRUE;
newSeq: REF CornerSeq ← NEW[CornerSeq[nCrnrs]];
newSeq[0] ← corner[0]; newSeq.length ← 1;
WHILE someLeft DO
someLeft ← someFound ← FALSE;
FOR j: NAT IN [1..nCrnrs) DO
IF newSeq[newSeq.length-1].inVtx = corner[j].outVtx -- match last incoming edge
THEN {             -- next corner in clockwiser order
newSeq[newSeq.length] ← corner[j]; newSeq.length ← newSeq.length+1;
corner[j].outVtx ← corner[j].inVtx ← LAST[NAT]; -- no further use
someFound ← TRUE;
}
ELSE IF newSeq[0].outVtx = corner[j].inVtx    -- match 1st outgoing edge
THEN {           -- previous corner in clockwiser order
FOR k: NAT DECREASING IN [0..newSeq.length) DO
newSeq[k+1] ← newSeq[k];
ENDLOOP;
newSeq[0] ← corner[j]; newSeq.length ← newSeq.length+1;
corner[j].outVtx ← corner[j].inVtx ← LAST[NAT]; -- no further use
someFound ← TRUE;
}
ELSE IF corner[j].outVtx # LAST[NAT] THEN someLeft ← TRUE; -- not done
ENDLOOP;
IF NOT someFound AND someLeft THEN  -- more than one chain, apparently
SIGNAL G3dRender.Error[$Fatal, "Multiple open edges at a vertex not allowed"];
ENDLOOP;
IF newSeq.length > 0 THEN corner ← newSeq;
RETURN[ corner ];
};
InitClasses[];
END.