ImplicitValueImpl.mesa
Copyright Ó 1985, 1990 by Xerox Corporation. All rights reserved.
Bloomenthal, November 21, 1992 7:00 pm PST
DIRECTORY G3dBasic, G3dPlane, G3dVector, ImplicitDefs, ImplicitValue, RealFns;
ImplicitValueImpl: CEDAR MONITOR
IMPORTS G3dPlane, G3dVector, RealFns
EXPORTS ImplicitValue
~ BEGIN
Pair:      TYPE ~ G3dBasic.Pair;
Triple:     TYPE ~ G3dBasic.Triple;
TripleSequence:   TYPE ~ G3dBasic.TripleSequence;
TripleSequenceRep:  TYPE ~ G3dBasic.TripleSequenceRep;
Plane:      TYPE ~ G3dPlane.Plane;
NearSegment:   TYPE ~ G3dVector.NearSegment;
DistanceMode:   TYPE ~ ImplicitDefs.DistanceMode;
Sample:     TYPE ~ ImplicitDefs.Sample;
Source2d:     TYPE ~ ImplicitValue.Source2d;
Source2dSequence:  TYPE ~ ImplicitValue.Source2dSequence;
Source2dSequenceRep: TYPE ~ ImplicitValue.Source2dSequenceRep;
Source3d:     TYPE ~ ImplicitValue.Source3d;
Source3dSequence:  TYPE ~ ImplicitValue.Source3dSequence;
Source3dSequenceRep: TYPE ~ ImplicitValue.Source3dSequenceRep;
ValueProc:    TYPE ~ ImplicitValue.ValueProc;
Support
SqDist3d: PROC [i, j: Triple] RETURNS [REAL] ~ {
d: Triple ~ [j.x-i.x, j.y-i.y, j.z-i.z];
RETURN[d.x*d.x+d.y*d.y+d.z*d.z];
};
Dist3d: PROC [i, j: Triple] RETURNS [REAL] ~ {
d: Triple ~ [j.x-i.x, j.y-i.y, j.z-i.z];
RETURN[RealFns.SqRt[d.x*d.x+d.y*d.y+d.z*d.z]];
};
SqDist2d: PROC [i, j: Pair] RETURNS [REAL] ~ {
d: Pair ~ [j.x-i.x, j.y-i.y];
RETURN[d.x*d.x+d.y*d.y];
};
Dist2d: PROC [i, j: Pair] RETURNS [REAL] ~ {
d: Pair ~ [j.x-i.x, j.y-i.y];
RETURN[RealFns.SqRt[d.x*d.x+d.y*d.y]];
};
Point and Line Segment Sources
OfPoint: PUBLIC PROC [
q, source: Triple,
strength: REAL,
distanceMode: DistanceMode]
RETURNS [sample: Sample]
~ {
sample.near ¬ source;
sample.vector ¬ [source.x-q.x, source.y-q.y, source.z-q.z];
sample.squareDistance ¬
sample.vector.x*sample.vector.x+sample.vector.y*sample.vector.y+
sample.vector.z*sample.vector.z;
SELECT distanceMode FROM
inverse, wtInverse => {
sample.distance ¬ RealFns.SqRt[sample.squareDistance];
sample.distanceSet ¬ TRUE;
sample.value ¬ strength/MAX[0.0001, sample.distance];
};
inverseSqrd, wtInverseSqrd =>
sample.value ¬ strength*strength/MAX[0.0001, sample.squareDistance];
ENDCASE;
};
Miscellany
CombineSamples: PUBLIC PROC [samples: LIST OF Sample] RETURNS [sample: Sample] ~ {
value: REAL ¬ 0.0;
FOR s: LIST OF Sample ¬ samples, s.rest WHILE s # NIL DO
IF s.first.value > sample.value THEN sample ¬ s.first;
value ¬ value+s.first.value;
ENDLOOP;
sample.value ¬ value;
};
SourceSequences
nScratchSourcer2d: NAT ~ 6;
nScratchSources3d: NAT ~ 6;
scratchSourcer2d: ARRAY [0..nScratchSourcer2d) OF Source2dSequence ¬ ALL[NIL];
scratchSources3d: ARRAY [0..nScratchSources3d) OF Source3dSequence ¬ ALL[NIL];
ObtainSources2d: PUBLIC ENTRY PROC [nSources: NAT] RETURNS [Source2dSequence] ~ {
FOR i: NAT IN [0..nScratchSourcer2d) DO
sources: Source2dSequence ¬ scratchSourcer2d[i];
IF sources # NIL THEN {
scratchSourcer2d[i] ¬ NIL;
IF sources.maxLength < nSources THEN sources ¬ NEW[Source2dSequenceRep[nSources]];
RETURN [sources];
};
ENDLOOP;
RETURN[NEW[Source2dSequenceRep[nSources]]];
};
ObtainSources3d: PUBLIC ENTRY PROC [nSources: NAT] RETURNS [Source3dSequence] ~ {
FOR i: NAT IN [0..nScratchSources3d) DO
sources: Source3dSequence ¬ scratchSources3d[i];
IF sources # NIL THEN {
scratchSources3d[i] ¬ NIL;
IF sources.maxLength < nSources THEN sources ¬ NEW[Source3dSequenceRep[nSources]];
RETURN [sources];
};
ENDLOOP;
RETURN[NEW[Source3dSequenceRep[nSources]]];
};
ReleaseSources2d: PUBLIC ENTRY PROC [sources: Source2dSequence] ~ {
FOR i: NAT IN [0..nScratchSourcer2d) DO
IF scratchSourcer2d[i] = NIL THEN {
scratchSourcer2d[i] ¬ sources;
RETURN;
};
ENDLOOP;
};
ReleaseSources3d: PUBLIC ENTRY PROC [sources: Source3dSequence] ~ {
FOR i: NAT IN [0..nScratchSources3d) DO
IF scratchSources3d[i] = NIL THEN {
scratchSources3d[i] ¬ sources;
RETURN;
};
ENDLOOP;
};
WtFunction: TYPE ~ PROC [strength, sqDist, dist: REAL] RETURNS [val: REAL];
OfSources2d: PUBLIC PROC [
q: Pair,
sources: Source2dSequence,
distanceMode: DistanceMode]
RETURNS [sum: REAL ¬ 0.0]
~ {
Inverse: PROC ~ {
FOR n: NAT IN [0..sources.length) DO
s: Source2d ~ sources[n];
sum ¬ sum+s.strength/MAX[0.0001, Dist2d[s.p, q]];
ENDLOOP;
};
InverseSquared: PROC ~ {
FOR n: NAT IN [0..sources.length) DO
s: Source2d ~ sources[n];
sum ¬ sum+s.strength*s.strength/MAX[0.0001, SqDist2d[s.p, q]];
ENDLOOP;
};
Wt: PROC [f: WtFunction] ~ {
info: TripleSequence ¬ NEW[TripleSequenceRep[sources.length]]; -- [sqDist, dist, val]
FOR i: NAT IN [0..sources.length) DO
info[i].y ¬ RealFns.SqRt[info[i].x ¬ SqDist2d[sources[i].p, q]];
info[i].z ¬ f[sources[i].strength, info[i].x, info[i].y];
ENDLOOP;
FOR n: NAT IN [0..sources.length) DO
FOR nn: NAT IN [n+1..sources.length) DO
d: REAL ~ Dist2d[sources[n].p, sources[nn].p];
sum ¬ sum+(d/(info[n].y+info[nn].y))*(info[n].z+info[nn].z);
ENDLOOP;
ENDLOOP;
};
WtInverse: PROC ~ {
F: WtFunction ~ {val ¬ strength/MAX[0.0001, sqDist]};
Wt[F];
};
WtInvSq: PROC ~ {
F: WtFunction ~ {val ¬ strength*strength/MAX[0.0001, dist]};
Wt[F];
};
SELECT distanceMode FROM
inverse   => Inverse[];
inverseSqrd  => InverseSquared[];
wtInverse  => IF sources.length = 1 THEN Inverse[] ELSE WtInverse[];
wtInverseSqrd => IF sources.length = 1 THEN InverseSquared[] ELSE WtInvSq[];
ENDCASE;
};
OfSources3d: PUBLIC PROC [
q: Triple,
sources: Source3dSequence,
distanceMode: DistanceMode]
RETURNS [sum: REAL ¬ 0.0]
~ {
Inverse: PROC ~ {
FOR n: NAT IN [0..sources.length) DO
s: Source3d ~ sources[n];
sum ¬ sum+s.strength/MAX[0.0001, Dist3d[s.p, q]];
ENDLOOP;
};
InverseSquared: PROC ~ {
FOR n: NAT IN [0..sources.length) DO
s: Source3d ~ sources[n];
sum ¬ sum+s.strength*s.strength/MAX[0.0001, SqDist3d[s.p, q]];
ENDLOOP;
};
Wt: PROC [f: WtFunction] ~ {
info: TripleSequence ¬ NEW[TripleSequenceRep[sources.length]]; -- [sqDist, dist, val]
FOR i: NAT IN [0..sources.length) DO
info[i].y ¬ RealFns.SqRt[info[i].x ¬ SqDist3d[sources[i].p, q]];
info[i].z ¬ f[sources[i].strength, info[i].x, info[i].y];
ENDLOOP;
FOR n: NAT IN [0..sources.length) DO
FOR nn: NAT IN [n+1..sources.length) DO
d: REAL ~ Dist3d[sources[n].p, sources[nn].p];
sum ¬ sum+(d/(info[n].y+info[nn].y))*(info[n].z+info[nn].z);
ENDLOOP;
ENDLOOP;
};
WtInverse: PROC ~ {
F: WtFunction ~ {val ¬ strength/MAX[0.0001, sqDist]};
Wt[F];
};
WtInvSq: PROC ~ {
F: WtFunction ~ {val ¬ strength*strength/MAX[0.0001, dist]};
Wt[F];
};
SELECT distanceMode FROM
inverse   => Inverse[];
inverseSqrd  => InverseSquared[];
wtInverse  => IF sources.length = 1 THEN Inverse[] ELSE WtInverse[];
wtInverseSqrd => IF sources.length = 1 THEN InverseSquared[] ELSE WtInvSq[];
ENDCASE;
};
Simple Geometric Shapes
OfPlane: PUBLIC PROC [q: Triple, plane: Plane] RETURNS [sample: Sample] ~ {
sample.distance ¬ G3dPlane.DistanceToPoint[q, plane];
sample.distanceSet ¬ TRUE;
sample.squareDistance ¬ sample.distance*sample.distance;
};
OfSphere: PUBLIC PROC [q, center: Triple, radiusSqd: REAL] RETURNS [sample: Sample] ~ {
sample.near ¬ center;
sample.squareDistance ¬ G3dVector.SquareDistance[center, q];
sample.value ¬ radiusSqd*sample.squareDistance;
};
OfCylinder: PUBLIC PROC [q: Triple, end0, end1: Triple, radiusSqd: REAL]
RETURNS [sample: Sample]
~ {
n: NearSegment ~ G3dVector.NearestToSegment[end0, end1, q];
sample.near ¬ n.point;
sample.squareDistance ¬ G3dVector.SquareDistance[q, sample.near];
sample.value ¬ radiusSqd*sample.squareDistance;
};
OfTorus: PUBLIC PROC [q: Triple, minorRadius, majorRadius: REAL]
RETURNS [sample: Sample]
~ {
r2:  REAL ~ majorRadius*majorRadius;
rr2: REAL ~ minorRadius*minorRadius;
r4:  REAL ~ r2*r2;
rr4: REAL ~ rr2*rr2;
a:  REAL ~ r2+rr2;
b:  REAL ~ r2-rr2;
c:  REAL ~ r2*rr2;
t:  Triple ~ G3dVector.MulVectors[q, q];
sample.value ¬ t.x*t.x+t.y*t.y+t.z*t.z+r4+rr4+2.0*(t.x*t.y+t.x*t.z+t.y*t.z-a*t.x+b*t.y-a*t.z-c);
};
END.
..
OfSources2d: PUBLIC PROC [
q: Pair,
sources: Source2dSequence,
distanceMode: DistanceMode,
threshold: REAL ¬ 1.0]
RETURNS [sum: REAL ¬ 0.0]
~ {
Inverse: PROC ~ {
FOR n: NAT IN [0..sources.length) DO
s: Source2d ~ sources[n];
sum ¬ sum+s.strength*threshold/MAX[0.0001, Dist2d[s.p, q]];
ENDLOOP;
};
InvSq: PROC ~ {
FOR n: NAT IN [0..sources.length) DO
s: Source2d ~ sources[n];
sum ¬ sum+s.strength*s.strength*threshold/MAX[0.0001, SqDist2d[s.p, q]];
ENDLOOP;
};
Wt: PROC [f: WtFunction] ~ {
info: TripleSequence ¬ NEW[TripleSequenceRep[sources.length]]; -- [sqDist, dist, val]
FOR i: NAT IN [0..sources.length) DO
info[i].y ¬ RealFns.SqRt[info[i].x ¬ SqDist2d[sources[i].p, q]];
info[i].z ¬ f[sources[i].strength, info[i].x, info[i].y];
ENDLOOP;
FOR n: NAT IN [0..sources.length) DO
FOR nn: NAT IN [n+1..sources.length) DO
d: REAL ~ Dist2d[sources[n].p, sources[nn].p];
sum ¬ sum+(d/(info[n].y+info[nn].y))*(info[n].z+info[nn].z);
ENDLOOP;
ENDLOOP;
};
WtInverse: PROC ~ {
F: WtFunction ~ {val ¬ strength*threshold/MAX[0.0001, sqDist]};
Wt[F];
};
WtInvSq: PROC ~ {
F: WtFunction ~ {val ¬ strength*strength*threshold/MAX[0.0001, dist]};
Wt[F];
};
SELECT distanceMode FROM
inverse  => Inverse[];
invSq   => InvSq[];
wtInverse => IF sources.length = 1 THEN Inverse[] ELSE WtInverse[];
wtInvSq  => IF sources.length = 1 THEN InvSq[] ELSE WtInvSq[];
ENDCASE;
};
OfMesh: PUBLIC PROC [
q: Triple,
mesh: Mesh,
distanceMode: DistanceMode,
threshold: REAL ¬ 1.0]
RETURNS [sample: Sample]
~ {
Notes: an inside value is derived from the length of a perpendicular dropped to a primitive;
in the case of a segment, this is a perpendicular to the line containing the segment, in the
case of triangles and polygons, it is a perpendicular dropped to the plane containing the
triangle or polygon. If the intersection of perpendicular with line (plane) does not lie within
the segment (triangle or polygon) distance to the closest segment endpoint (triangle or
polygon edge) is used to compute a value.
s: SegmentSequence ~ mesh.segments;
t: TriangleSequence ~ mesh.triangles;
p: PolygonSequence ~ mesh.polygons;
temp: Sample;
sInside: Sample ¬ OfInsideSegments[q, s, distanceMode, threshold];
tInside: Sample ¬ OfInsideTriangles[q, t, distanceMode, threshold];
pInside: Sample ¬ OfInsidePolygons[q, p, distanceMode, threshold];
SELECT TRUE FROM
sInside.value = 0.0 AND tInside.value = 0.0 AND pInside.value = 0.0 => {
Nothing registers inside, so find closest edge:
temp ¬ OfEdgeSegments[q, s, distanceMode, threshold];
IF temp.value > sample.value THEN sample ¬ temp;
temp ¬ OfEdgeTriangles[q, t, distanceMode, threshold];
IF temp.value > sample.value THEN sample ¬ temp;
temp ¬ OfEdgePolygons[q, p, distanceMode, threshold];
IF temp.value > sample.value THEN sample ¬ temp;
};
tInside.value = 0.0 AND pInside.value = 0.0 => {
A segment registers inside, so add in closest edge of triangles and polygons:
sample ¬ sInside;
temp ¬ OfEdgeTriangles[q, t, distanceMode, threshold];
IF temp.value > sample.value THEN sample ¬ temp;
temp ¬ OfEdgePolygons[q, p, distanceMode, threshold];
IF temp.value > sample.value THEN sample ¬ temp;
};
ENDCASE => {
sample ¬ sInside;
IF tInside.value > sample.value THEN sample ¬ tInside;
IF pInside.value > sample.value THEN sample ¬ pInside;
sample.value ¬ sInside.value+tInside.value+pInside.value;
};
};
Skirt
nScratchItemSequences:  NAT ~ 10;
scratchItemSequences:  ARRAY [0..nScratchItemSequences) OF ItemSequence ¬ ALL[NIL];
ObtainItems: ENTRY PROC [nItems: NAT] RETURNS [ItemSequence] ~ {
IF nItems = 0 THEN RETURN[NIL];
FOR i: NAT IN [0..nScratchItemSequences) DO
items: ItemSequence ¬ scratchItemSequences[i];
IF items = NIL THEN LOOP;
scratchItemSequences[i] ¬ NIL;
IF items.maxLength < nItems THEN items ¬ NEW[ItemSequenceRep[nItems]];
items.length ¬ 0;
RETURN [items];
ENDLOOP;
RETURN[NEW[ItemSequenceRep[nItems]]];
};
ReleaseItems: ENTRY PROC [items: ItemSequence] ~ {
FOR i: NAT IN [0..nScratchItemSequences) DO
IF scratchItemSequences[i] # NIL THEN LOOP;
scratchItemSequences[i] ¬ items;
RETURN;
ENDLOOP;
};
NItems: PROC [mesh: Mesh] RETURNS [n: NAT ¬ 0] ~ {
IF mesh = NIL THEN RETURN;
IF mesh.segments # NIL THEN n ¬ n+mesh.segments.length;
IF mesh.triangles # NIL THEN n ¬ n+mesh.triangles.length;
IF mesh.polygons # NIL THEN n ¬ n+mesh.polygons.length;
};
TriangleStrength: PROC [triangle: Triangle, near: NearTriangle] RETURNS [REAL] ~ {
RETURN[
(triangle.s0*near.w0+triangle.s1*near.w1+triangle.s2*near.w2)/(near.w0+near.w1+near.w2)];
};
DoItems: PROC [
q: Triple,
mesh: Mesh,
distanceMode: DistanceMode,
segmentScale: REAL ¬ 1.0,
triangleScale: REAL ¬ 1.0,
threshold: REAL ¬ 1.0]
RETURNS [items: ItemSequence]
~ {
AddToItemSequence: PROC [item: Item, scale: REAL] ~ {
item.sample.value ¬ item.sample.value*scale;
items[items.length] ¬ item;
items.length ¬ items.length+1;
};
IF (items ¬ ObtainItems[NItems[mesh]]) = NIL THEN RETURN;
IF mesh.triangles # NIL THEN FOR n: NAT IN [0..mesh.triangles.length) DO
triangle: Triangle ~ mesh.triangles[n];
near: NearTriangle ~ G3dPolygon.NearestToTriangle[triangle, q];
sample: Sample ¬ IF near.inside
THEN ImplicitValue.OfPoint[
q, near.point, TriangleStrength[triangle, near], distanceMode, threshold]
ELSE ImplicitValue.OfTriangleEdges[q, triangle, distanceMode, threshold];
AddToItemSequence[[sample, near.inside, tripoly], triangleScale];
ENDLOOP;
IF mesh.polygons # NIL THEN FOR n: NAT IN [0..mesh.polygons.length) DO
polygon: Polygon ~ mesh.polygons[n];
near: NearPolygon ~ G3dPolygon.NearestToPolygon[polygon, q];
sample: Sample ¬ IF near.inside
THEN ImplicitValue.OfPoint[
q, near.point, ImplicitValue.PolygonStrength[polygon, q], distanceMode, threshold]
ELSE ImplicitValue.OfPolygonEdges[q, polygon, distanceMode, threshold];
AddToItemSequence[[sample, near.inside, tripoly], triangleScale];
ENDLOOP;
IF mesh.segments # NIL THEN FOR n: NAT IN [0..mesh.segments.length) DO
seg: Segment3d ~ mesh.segments[n];
near: NearSegment ~ G3dVector.NearestTo3dSegment[seg.p0, seg.p1, q, seg.acc];
sample: Sample ¬ ImplicitValue.OfPoint[
q, near.point, seg.s0*near.w0+seg.s1*near.w1, distanceMode, threshold];
AddToItemSequence[[sample, near.inside, segment], segmentScale];
ENDLOOP;
};
MaxItem: PROC [items: ItemSequence] RETURNS [index: INTEGER ¬ -1] ~ {
Inside or edge, whichever max
max: REAL ¬ -100000.0;
IF items # NIL THEN FOR n: NAT IN [0..items.length) DO
IF items[n].sample.value > max THEN {max ¬ items[n].sample.value; index ¬ n};
ENDLOOP;
};
MaxSegmentItem: PROC [items: ItemSequence] RETURNS [index: INTEGER ¬ -1] ~ {
max: REAL ¬ -100000.0;
IF items # NIL THEN FOR n: NAT IN [0..items.length) DO
i: Item ¬ items[n];
IF i.type = segment AND i.sample.value > max THEN {max ¬ i.sample.value; index ¬ n};
ENDLOOP;
};
MaxTripolyItem: PROC [items: ItemSequence] RETURNS [index: INTEGER ¬ -1] ~ {
max: REAL ¬ -100000.0;
IF items # NIL THEN FOR n: NAT IN [0..items.length) DO
i: Item ¬ items[n];
IF i.type = tripoly AND i.sample.value > max THEN {max ¬ i.sample.value; index ¬ n};
ENDLOOP;
};
SimpleValue: PROC [
q: Triple,
mesh: Mesh,
distanceMode: DistanceMode,
segmentScale: REAL ¬ 1.0,
triangleScale: REAL ¬ 1.0,
threshold: REAL ¬ 1.0]
RETURNS [sample: Sample]
~ {
items: ItemSequence ¬ DoItems[q, mesh, distanceMode, segmentScale, triangleScale, threshold];
IF items # NIL THEN {
Sum: PROC [n: INTEGER, v: Triple] RETURNS [REAL ¬ 0.0] ~ {
s: Sample ¬ items[n].sample;
l: REAL ¬ IF s.distanceSet THEN s.distance ELSE RealFns.SqRt[s.squareDistance];
IF l # 0.0 THEN {
acos: REAL ¬ (v.x*s.vector.x+v.y*s.vector.y+v.z*s.vector.z)/l;
RETURN[0.5*(1.0-acos)*s.value];
};
};
SumSegments: PROC RETURNS [sample: Sample] ~ {
nMax: NAT ¬ MaxSegmentItem[items];
max: Item ¬ items[nMax];
sample ¬ max.sample;
IF mesh.segments.length = 1
THEN RETURN[max.sample]
ELSE {
v: Triple ¬ G3dVector.Normalize[max.sample.vector];
FOR n: NAT IN [0..items.length) DO
IF n # nMax AND items[n].type = segment AND items[n].inside
THEN sample.value ¬ sample.value+Sum[n, v];
ENDLOOP;
};
};
SumTripolys: PROC RETURNS [sample: Sample] ~ {
nTripolys: INTEGER ¬
(IF mesh.triangles = NIL THEN 0 ELSE mesh.triangles.length)+
(IF mesh.polygons = NIL THEN 0 ELSE mesh.polygons.length);
nMax: NAT ¬ MaxTripolyItem[items];
max: Item ¬ items[nMax];
sample ¬ max.sample;
IF nTripolys = 1
THEN RETURN[max.sample]
ELSE {
v: Triple ¬ G3dVector.Normalize[max.sample.vector];
FOR n: NAT IN [0..items.length) DO
IF n # nMax AND items[n].type = tripoly AND items[n].inside
THEN sample.value ¬ sample.value+Sum[n, v];
ENDLOOP;
};
};
SELECT TRUE FROM
mesh.segments = NIL => RETURN[SumTripolys[]];
mesh.triangles = NIL AND mesh.polygons = NIL => RETURN[SumSegments[]];
ENDCASE => {
sSample: Sample ¬ SumSegments[];
tSample: Sample ¬ SumTripolys[];
sample ¬ IF sSample.value > tSample.value THEN sSample ELSE tSample;
sample.value ¬ sSample.value+tSample.value;
};
};
ReleaseItems[items];
};
CleverValue: PROC [
q: Triple,
mesh: Mesh,
distanceMode: DistanceMode,
segmentScale: REAL ¬ 1.0,
triangleScale: REAL ¬ 1.0,
threshold: REAL ¬ 1.0]
RETURNS [sample: Sample]
~ {
items: ItemSequence ¬ DoItems[q, mesh, distanceMode, segmentScale, triangleScale, threshold];
sample ¬ DoClever[q, items, distanceMode, threshold];
ReleaseItems[items];
};
CleverValueDebug: PROC [
q: Triple,
mesh: Mesh,
distanceMode: DistanceMode,
segmentScale: REAL ¬ 1.0,
triangleScale: REAL ¬ 1.0,
threshold: REAL ¬ 1.0]
RETURNS [sample: Sample, items: ItemSequence]
~ {
items ¬ DoItems[q, mesh, distanceMode, segmentScale, triangleScale, threshold];
sample ¬ DoClever[q, items, distanceMode, threshold];
};
DoClever: PROC [
q: Triple,
items: ItemSequence,
distanceMode: DistanceMode,
threshold: REAL]
RETURNS [sample: Sample]
~ {
IF items # NIL THEN {
Sum: PROC [n: INTEGER] ~ {
IF n # -1 THEN {
s: Sample ¬ items[n].sample;
l: REAL ¬ IF s.distanceSet THEN s.distance ELSE RealFns.SqRt[s.squareDistance];
IF l # 0.0 THEN {
acos: REAL ¬ G3dVector.Dot[v, s.vector]/l;
sample.value ¬ sample.value+0.5*(1.0-acos)*s.value;
};
};
};
nMax: NAT ¬ MaxItem[items];
max: Item ¬ items[nMax];
v: Triple ¬ G3dVector.Normalize[max.sample.vector];
sample ¬ max.sample;
SELECT max.type FROM
segment => {
FOR n: NAT IN [0..items.length) DO
IF n # nMax AND items[n].inside AND items[n].type = segment THEN Sum[n];
ENDLOOP;
Sum[MaxTripolyItem[items]];
};
tripoly => IF max.inside
THEN {
FOR n: NAT IN [0..items.length) DO
IF n # nMax AND items[n].inside AND items[n].type = tripoly THEN Sum[n];
ENDLOOP;
Sum[MaxSegmentItem[items]];
}
ELSE {
someInsideSegments: BOOL ¬ FALSE;
FOR n: NAT IN [0..items.length) DO
IF items[n].type # segment OR NOT items[n].inside THEN LOOP;
someInsideSegments ¬ TRUE;
Sum[n];
ENDLOOP;
IF NOT someInsideSegments THEN Sum[MaxSegmentItem[items]];
};
ENDCASE;
};
};