TubeGeometryImpl.mesa
Copyright © 1985 by Xerox Corporation. All rights reserved.
Bloomenthal, February 24, 1987 6:14:43 pm PST
DIRECTORY Contours, Cubic2, Matrix3d, Real, RealFns, Spline3d, TubeContour, Vector3d, TubeDefs, TubeGeometry, TubeMisc;
TubeGeometryImpl: CEDAR PROGRAM
IMPORTS Contours, Cubic2, Matrix3d, Real, RealFns, Spline3d, TubeContour, TubeMisc, Vector3d
EXPORTS TubeGeometry
~ BEGIN
OPEN TubeDefs;
Setting Spline Coefficients
SetCoeffs: PUBLIC PROC [tube: Tube] ~ {
tubeProc: TubeProc ~ {SetSectionCoeffs[tube]};
TubeMisc.ApplyToTube[tube, tubeProc];
};
SetSectionCoeffs: PUBLIC PROC [tube: Tube] ~ {
tube.coeffs ← Spline3d.CoeffsFromHermite[
[tube.p0, tube.p1, Vector3d.Mul[tube.v0, tube.tens0], Vector3d.Mul[tube.v1, tube.tens1]],
tube.coeffs];
};
Setting Tube Lengths
SetLengths: PUBLIC PROC [tube: Tube] ~ {
tubeProc: TubeProc ~ {SetSectionLengths[tube]};
TubeMisc.ApplyToTube[tube, tubeProc];
};
SetSectionLengths: PUBLIC PROC [tube: Tube] ~ {
IF tube = NIL THEN RETURN;
IF tube.prev # NIL THEN tube.length0 ← tube.prev.length1;
tube.length1 ← tube.length0+Spline3d.Length[tube.coeffs];
};
Creating Tube Frames
Make: PUBLIC PROC [
tube: Tube,
epsilon: REAL ← 0.03,
scale: REAL ← 1.0,
taper: REAL ← 0.0,
skin: BOOL ← FALSE,
roundTip: BOOL ← FALSE,
view: Matrix ← NIL,
frameProc: FrameProc ← NIL]
~ {
SetCoeffs[tube];
SetLengths[tube];
MakeFrames[tube, epsilon, scale, taper, skin, roundTip, view, frameProc];
};
MakeFrames: PUBLIC PROC [
tube: Tube,
epsilon: REAL ← 0.03,
scale: REAL ← 1.0,
taper: REAL ← 0.0,
skin: BOOL ← FALSE,
roundTip: BOOL ← FALSE,
view: Matrix ← NIL,
frameProc: FrameProc ← NIL]
~ {
tubeProc: TubeProc ~ {
scale0, scale1: REAL ← scale;
IF tube = NIL THEN RETURN;
IF taper # 0.0 THEN {
scale0 ← scale*(1.0-taper*tube.length0);
scale1 ← scale*(1.0-taper*tube.length1);
};
MakeSectionFrames[tube, epsilon, scale0, scale1, skin, roundTip, view, frameProc];
};
TubeMisc.ApplyToTube[tube, tubeProc];
};
MakeSectionFrames: PUBLIC PROC [
tube: Tube,
epsilon: REAL ← 0.03,
scale0, scale1: REAL ← 1.0,
skin: BOOL ← FALSE,
roundTip: BOOL ← FALSE,
view: Matrix ← NIL,
frameProc: FrameProc ← NIL]
~ {
bez: Bezier ← Spline3d.BezierFromCoeffs[tube.coeffs];
maxLen: REAL ← IF view = NIL THEN 0.01*Spline3d.ConvexHullLength[bez] ELSE 0.0;
vv: Triple ← IF Vector3d.Equal[bez.b1, bez.b0, 0.0001]
THEN Vector3d.Sub[bez.b2, bez.b1]
ELSE Vector3d.Sub[bez.b1, bez.b0];
circle: Contour ← Contours.Circle[tube.circleRes];
c0: Contour ← TubeContour.TContour[tube.contours, circle, 0.0];
c1: Contour ← TubeContour.TContour[tube.contours, circle, 1.0];
MakeFrame: PROC [position, velocity: Triple, t: REAL, c: Contour] ~ {
tube.frames ← TestFrameSequence[tube.frames];
IF frameProc # NIL
THEN tube.frames[tube.frames.length] ← frameProc[position, velocity, t]
ELSE {
f: Frame ← tube.frames[tube.frames.length];
f.t ← t;
f.position ← position;
f.triad ← Basis[velocity, vv, tube.refVec];
Keep scale non-zero (to avoid subsequent surface normal computation problems):
f.scale ← MAX[.000001, (scale0+t*(scale1-scale0))*TubeMisc.Radius[tube, t, roundTip]];
f.twist ← tube.tw0+t*(tube.tw1-tube.tw0);
f.matrix ← RefMatrix[f.position, f.triad, f.scale, f.twist, f.matrix];
IF skin THEN {
f.contour ← c;
f.normals ← SELECT TRUE FROM
c.circle => c.pairs,
c.normals # NIL => c.normals,
ENDCASE => Contours.Normals[c];
};
};
tube.frames.length ← tube.frames.length+1;
vv ← velocity;
tube.refVec ← tube.frames[tube.frames.length-1].triad.n;
};
Walker: PROC [bez: Bezier, t0, t1: REAL, c0, c1: Contour] ~ {
IF DividedEnough[bez, t0, t1, c0, c1]
THEN MakeFrame[
bez.b0,
IF Vector3d.Equal[bez.b1, bez.b0, 0.0001]  -- degenerate end tangent check
THEN Vector3d.Sub[bez.b2, bez.b1]
ELSE Vector3d.Sub[bez.b1, bez.b0],
t0,
c0]
ELSE {
left, rite: Bezier;
t: REAL ← 0.5*(t0+t1);
c: Contour ← IF skin THEN TubeContour.TContour[tube.contours, circle, t] ELSE NIL;
[left, rite] ← Spline3d.SplitBezier[bez];
Walker[left, t0, t, c0, c];
Walker[rite, t, t1, c, c1];
};
};
DividedEnough: PROC [bez: Bezier, t0, t1: REAL, c0, c1: Contour] RETURNS [BOOL] ~ {
BezSmallEnough: PROC [bez: Bezier] RETURNS[BOOL] ~ {
RETURN[IF view = NIL
THEN Spline3d.ConvexHullLength[bez] < maxLen OR Spline3d.FlatBezier[bez, epsilon]
ELSE Cubic2.Flat[
[Matrix3d.TransformD[bez.b0, view], Matrix3d.TransformD[bez.b1, view],
Matrix3d.TransformD[bez.b2, view], Matrix3d.TransformD[bez.b3, view]],
100.0*epsilon]];
};
TwistSmallEnough: PROC [t0, t1: REAL] RETURNS[BOOL] ~ {
RETURN[ABS[(t1-t0)*(tube.tw1-tube.tw0)] < 20.0];
};
ContoursCloseEnough: PROC [c0, c1: Contour] RETURNS [ret: BOOL] ~ {
AllCircles: PROC [contours: ContourSequence] RETURNS [BOOL] ~ {
IF contours # NIL THEN FOR n: NAT IN [0..contours.length) DO
IF NOT contours[n].circle THEN RETURN [FALSE];
ENDLOOP;
RETURN[TRUE];
};
IF AllCircles[tube.contours]
THEN {
IF roundTip
THEN {
r0: REAL ~ TubeMisc.Radius[tube, t0, TRUE];
r1: REAL ~ TubeMisc.Radius[tube, t1, TRUE];
IF r0 = 0.0 OR r1 = 0.0 THEN RETURN[FALSE];
RETURN[1.0-MIN[r0, r1]/MAX[r0, r1] < 100.0*epsilon];
}
ELSE RETURN[TRUE];
};
RETURN[Contours.Similar[c0, c1] > 1.0-10.0*epsilon];
};
IF t1-t0 < 0.005 THEN RETURN[TRUE];
RETURN[
BezSmallEnough[bez] AND
TwistSmallEnough[t0, t1] AND
ContoursCloseEnough[c0, c1]
];
};
IF tube # NIL AND tube.frames # NIL THEN tube.frames.length ← 0;
tube.refVec ← IF tube.prev # NIL
THEN tube.prev.refVec
ELSE Spline3d.RefVec[tube.coeffs, 0.0];
Walker[bez, 0.0, 1.0, c0, c1];
MakeFrame[
bez.b3,
IF Vector3d.Equal[bez.b3, bez.b2, 0.0001] -- test for degenerate end tangent
THEN Vector3d.Sub[bez.b2, bez.b1]
ELSE Vector3d.Sub[bez.b3, bez.b2],
1.0,
c1];
IF view # NIL THEN tube.circleRes ←
Real.RoundI[TubeContour.ScreenCircleRes[tube, view]/MAX[0.01, 100*epsilon]];
};
Altering Tube Frames
ReScale: PUBLIC PROC [tube: Tube, scale: REAL ← 1.0, taper: REAL ← 0.0] ~ {
tubeProc: TubeProc ~ {
scale0: REAL ← scale*tube.r0;
scale1: REAL ← scale*tube.r1;
IF tube = NIL THEN RETURN;
IF taper # 0.0 THEN {
scale0 ← scale0*(1.0-taper*tube.length0);
scale1 ← scale1*(1.0-taper*tube.length1);
};
FOR n: NAT IN [0..NFrames[tube]) DO
f: Frame ← tube.frames[n];
f.scale ← MAX[0.000001, scale0+f.t*(scale1-scale0)];
f.matrix ← RefMatrix[f.position, f.triad, f.scale, f.twist, f.matrix];
ENDLOOP;
};
TubeMisc.ApplyToTube[tube, tubeProc];
};
SetCircles: PUBLIC PROC [tube: Tube, circleRes: NAT] ~ {
tubeProc: TubeProc ~ {
tube.circleRes ← circleRes;
FOR n: NAT IN [0..NFrames[tube]) DO
f: Frame ← tube.frames[n];
f.contour ← circle;
f.normals ← circle.normals;
ENDLOOP;
};
circle: Contour ← Contours.Circle[circleRes];
TubeMisc.ApplyToTube[tube, tubeProc];
};
Miscellaneous Frame Operations
NFrames: PUBLIC PROC [tube: Tube] RETURNS [INTEGER] ~ {
RETURN[IF tube.frames # NIL THEN tube.frames.length ELSE 0];
};
TestFrameSequence: PROC [frames: FrameSequence] RETURNS [FrameSequence] ~ {
IF frames = NIL OR frames.length >= frames.maxLength
THEN frames ← LengthenFrameSequence[frames];
RETURN[frames];
};
LengthenFrameSequence: PROC [frames: FrameSequence] RETURNS [FrameSequence] ~ {
nFrames: INTEGER ← IF frames = NIL THEN 0 ELSE frames.length;
temp: FrameSequence ← NEW[FrameSequenceRep[MAX[5, Real.RoundI[1.3*nFrames]]]];
FOR i: NAT IN [0..nFrames) DO temp[i] ← frames[i]; ENDLOOP;
FOR i: NAT IN [nFrames..temp.maxLength) DO temp[i] ← NEW[FrameRep]; ENDLOOP;
temp.length ← nFrames;
RETURN[temp];
};
GetFrame: PUBLIC PROC [tube: Tube, t: REAL] RETURNS [Frame] ~ {
nFrames: INTEGER ← NFrames[tube];
frames: FrameSequence ← tube.frames;
SELECT TRUE FROM
nFrames = 0 =>
RETURN[NIL];
t < frames[0].t =>
RETURN[frames[0]];
t > tube.frames[nFrames-1].t =>
RETURN[frames[nFrames-1]];
ENDCASE => {
n: NAT ← 0;
WHILE frames[n].t < t DO n ← n+1; ENDLOOP;
IF frames[n].t = t
THEN RETURN[frames[n]]
ELSE {
frame: Frame ← NEW[FrameRep];
frame0: Frame ← frames[n-1];
frame1: Frame ← frames[n];
alpha: REAL ← (frame.t-frame0.t)/(frame1.t-frame0.t);
frame.t ← t;
frame.position ← Spline3d.Position[tube.coeffs, t];
frame.triad ← Basis[Spline3d.Velocity[tube.coeffs, t], frame0.triad.v, frame0.triad.n];
frame.scale ← frame0.scale+alpha*(frame1.scale-frame0.scale);
frame.twist ← frame0.twist+alpha*(frame1.twist-frame0.twist);
frame.matrix ← RefMatrix[frame.position, frame.triad, frame.scale, frame.twist];
frame.contour ← Contours.Interpolate[frame0.contour, frame1.contour, alpha];
frame.normals ← Contours.Normals[frame.contour];
RETURN[frame];
};
};
};
PutFrame: PUBLIC PROC [tube: Tube, frame: Frame] ~ {
InsertFrame: PROC [frame: Frame, frames: FrameSequence, n: NAT] ~ {
FOR i: NAT DECREASING IN [n..frames.length) DO frames[n+1] ← frames[n]; ENDLOOP;
frames[n] ← frame;
};
tube.frames ← TestFrameSequence[tube.frames];
SELECT TRUE FROM
tube.frames.length = 0 =>
tube.frames[0] ← frame;
frame.t < tube.frames[0].t =>
InsertFrame[frame, tube.frames, 0];
frame.t > tube.frames[tube.frames.length-1].t =>
tube.frames[tube.frames.length] ← frame;
ENDCASE => {
n: NAT ← 0;
WHILE tube.frames[n].t < frame.t DO n ← n+1; ENDLOOP;
IF tube.frames[n].t = frame.t THEN {
tube.frames[n] ← frame;
RETURN;
};
InsertFrame[frame, tube.frames, n];
};
tube.frames.length ← tube.frames.length+1;
};
CopyFrame: PUBLIC PROC [frame: Frame] RETURNS [Frame] ~ {
copy: Frame ← NIL;
IF frame # NIL THEN {
copy ← NEW[FrameRep];
copy.t ← frame.t;
copy.position ← frame.position;
copy.triad ← frame.triad;
copy.scale ← frame.scale;
copy.matrix ← Matrix3d.CopyMatrix[frame.matrix];
copy.contour ← Contours.Copy[frame.contour];
};
RETURN[copy];
};
CopyFrames: PUBLIC PROC [frames: FrameSequence] RETURNS [FrameSequence] ~ {
copy: FrameSequence ← NIL;
IF frames # NIL THEN {
copy ← NEW[FrameSequenceRep[frames.length]];
copy.length ← frames.length;
FOR n: NAT IN [0..copy.length) DO copy[n] ← CopyFrame[frames[n]]; ENDLOOP;
};
RETURN[copy];
};
Frame Geometric Operations
Basis: PUBLIC PROC [v, vv, rv: Triple] RETURNS [Triad] ~ {
dot: REAL;
b, n: Triple;
vv ← Vector3d.Normalize[vv];
v ← Vector3d.Normalize[v];
rv ← Vector3d.Normalize[rv];
dot ← Vector3d.Dot[v, vv];
IF ABS[dot] > 0.9999 THEN {        -- tangents v, vv colinear
n ← IF dot > 0.0 THEN rv ELSE [-rv.x, -rv.y, -rv.z];
n ← Vector3d.V90[v, n];       -- ensure orthogonality
}
ELSE {
a: REAL ← RealFns.ArcTanDeg[RealFns.SqRt[1.0-dot*dot], dot]; -- == Acos[]
axis: Triple ← Vector3d.Cross[v, vv];
n ← Vector3d.Normalize[Vector3d.RotateAbout[rv, axis, a]]; -- incremental kludge
};
b ← Vector3d.Normalize[Vector3d.Cross[v, n]];
RETURN[[v, n, b]];
};
RefMatrix: PUBLIC PROC [position: Triple, triad: Triad, scale, twist: REAL, out: Matrix ← NIL]
RETURNS [Matrix] ~ {
IF ABS[twist] > 0.001 THEN {
out ← Matrix3d.MakePureRotate[triad.v, twist, , out];
triad.n ← Matrix3d.TransformVec[triad.n, out];
triad.b ← Matrix3d.TransformVec[triad.b, out];
};
RETURN[
Matrix3d.LocalScale[Matrix3d.MakeFromTriad[triad.n, triad.b, triad.v, position], scale, out]];
};
XYFromFrame: PUBLIC PROC [point: Triple, tube: Tube, t: REAL] RETURNS [Pair] ~ {
ThreeTriples: TYPE ~ RECORD [o, x, y: Triple]; -- origin, x and y vectors
GetOXY: PROC [tube: Tube, t: REAL] RETURNS [tt: ThreeTriples] ~ {
OXYFromFrame: PROC [f: Frame] RETURNS [tt: ThreeTriples] ~ {
m: Matrix ~ f.matrix;
tt.o ← f.position;
tt.x ← Vector3d.Normalize[[m[0][0], m[0][1], m[0][2]]];
tt.y ← Vector3d.Normalize[[m[1][0], m[1][1], m[1][2]]];
};
nFrames: INTEGER ~ NFrames[tube];
SELECT TRUE FROM
nFrames = 0 =>      RETURN[[origin, origin, origin]];
t < tube.frames[0].t =>    RETURN[OXYFromFrame[tube.frames[0]]];
t > tube.frames[nFrames-1].t => RETURN[OXYFromFrame[tube.frames[nFrames-1]]];
ENDCASE => {
n: NAT ← 0;
WHILE tube.frames[n].t < t DO n ← n+1; ENDLOOP;
IF tube.frames[n].t = t THEN RETURN[OXYFromFrame[tube.frames[n]]]
ELSE {
f0: Frame ← tube.frames[n-1];
f1: Frame ← tube.frames[n];
alpha: REAL ~ (t-f0.t)/(f1.t-f0.t);
x0: Triple ~ [f0.matrix[0][0], f0.matrix[0][1], f0.matrix[0][2]];
x1: Triple ~ [f1.matrix[0][0], f1.matrix[0][1], f1.matrix[0][2]];
y0: Triple ~ [f0.matrix[1][0], f0.matrix[1][1], f0.matrix[1][2]];
y1: Triple ~ [f1.matrix[1][0], f1.matrix[1][1], f1.matrix[1][2]];
tt.o ← Spline3d.Position[tube.coeffs, t];
tt.x ← Vector3d.Normalize[Vector3d.Combine[x0, 1.0-alpha, x1, alpha]];
tt.y ← Vector3d.Normalize[Vector3d.Combine[y0, 1.0-alpha, y1, alpha]];
};
};
};
tt: ThreeTriples ← GetOXY[tube, t];
q: Triple ← Vector3d.Sub[point, tt.o];   -- q now in tt space
RETURN[[Vector3d.Dot[q, tt.x], Vector3d.Dot[q, tt.y]]];
};
END.
..
MakeMatrices: PUBLIC PROC [s: SplineData, scale: REAL] ~ {
sc0: REAL ← s.rad0*scale;
sc1: REAL ← s.rad1*scale;
tw0: REAL ← s.twist0;
tw1: REAL ← s.twist1;
vv: Triple ← s.v0;
mats: MatrixSeq ← s.mats;
s.refVec ← IF s.prev # NIL THEN s.prev.refVec ELSE Spline3d.RefVec[s.c, 0.0];
mats[0] ←Misc3d.RefMatrix[s.p0, s.refVec, Vector3d.Cross[vv, s.refVec], vv, sc0, tw0, mats[0]];
FOR i: NAT IN[1..s.axialRes) DO
n, b: Triple;
t: REAL ← i/REAL[s.axialRes-1.0];
p: Triple ← Spline3d.Position[s.c, t];
v: Triple ← Spline3d.Tangent[s.c, t];
[n, b] ← Misc3d.Basis[v, vv, s.refVec];
mats[i] ← Misc3d.RefMatrix[p, n, b, v, (1.0-t)*sc0+t*sc1, (1.0-t)*tw0+t*tw1, mats[i]];
vv ← v;
s.refVec ← n;
ENDLOOP;
};
GetFrame: PUBLIC PROC [tube: Tube, t: REAL] RETURNS [Frame] ~ {
frames: FrameSequence ← tube.frames;
numFrames: INT ← IF frames # NIL THEN frames.maxLength ELSE 0;
SELECT TRUE FROM
numFrames = 0 =>
RETURN[NIL];
frame.t < frames[0].t =>
RETURN[tube.frames[0]];
frame.t > tube.frames[numFrames-1].t =>
RETURN[tube.frames[numFrames-1];
ENDCASE => {
n: NAT ← 0;
WHILE frames[n].t < frame.t DO tube.frames[n] ← frames[n]; n ← n+1; ENDLOOP;
IF tube.frames[n].t = t
THEN RETURN[tube.frames[n]]
ELSE {    -- the following is sloppy:
frame0: Frame ← tube.frames[n-1];
frame1: Frame ← tube.frames[n];
m0: Matrix ← frame0.m;
m1: Matrix ← frame1.m;
alpha: REAL ← (frame.t-frame0.t)/(frame1.t-frame0.t);
beta: REAL ← 1.0-alpha;
x: Triple ← Vector3d.Combine[
[m0[0][0], m0[0][1], m0[0][2]], alpha, [m1[0][0], m1[0][1], m1[0][2]], beta];
y: Triple ← Vector3d.Combine[
[m0[1][0], m0[1][1], m0[1][2]], alpha, [m1[1][0], m1[1][1], m1[1][2]], beta];
z: Triple ← Vector3d.Combine[
[m0[2][0], m0[2][1], m0[2][2]], alpha, [m1[2][0], m1[2][1], m1[2][2]], beta];
v: Triple ← Vector3d.SameMag[z, Spline3d.Velocity[tube.coeffs, t]];
n: Triple ← Vector3d.SameMag[x, Vector3d.Cross[v, y]];
b: Triple ← Vector3d.SameMag[y, Vector3d.Cross[v, n]];
p: Triple ← Spline3d.Position[tube.coeffs, t];
frame: Frame ← NEW[FrameRep];
frame.t ← t;
frame.c ← Contour.Interpolate[frame0.c, frame1.c, alpha];
frame.m ← RefMatrix[p, n, b, v, 1.0, 0.0];
RETURN[frame];
};
};
};