G3dCurveImpl.mesa
Copyright Ó 1984, 1992 by Xerox Corporation. All rights reserved.
Bloomenthal, October 21, 1992 5:01 pm PDT
DIRECTORY Convert, G3dBasic, G3dCurve, G3dMatrix, G3dSpline, IO, Real, Rope;
G3dCurveImpl: CEDAR PROGRAM
IMPORTS Convert, G3dMatrix, G3dSpline, IO, Real, Rope
EXPORTS G3dCurve
~ BEGIN
Errors
Error:    PUBLIC SIGNAL [code: ATOM, reason: ROPE ¬ NIL] = CODE;
Definitions
Triple:    TYPE ~ G3dBasic.Triple;    -- RECORD [x, y, z: REAL]
RealSequence:  TYPE ~ G3dBasic.RealSequence;  -- sequence of reals
Matrix:    TYPE ~ G3dMatrix.Matrix;    -- a cubic spline
Spline:    TYPE ~ G3dSpline.Spline;    -- a cubic spline
SplineSequence:  TYPE ~ G3dSpline.SplineSequence; -- a sequence of cubic splines
ROPE:     TYPE ~ Rope.ROPE;
STREAM:    TYPE ~ IO.STREAM;
Curve:    TYPE ~ G3dCurve.Curve;
CurveRep:   TYPE ~ G3dCurve.CurveRep;
Section:     TYPE ~ G3dCurve.Section;
Fragment:   TYPE ~ G3dCurve.Fragment;
Fragments:   TYPE ~ G3dCurve.Fragments;
Spot:     TYPE ~ G3dCurve.Spot;
CurveSequence:  TYPE ~ G3dCurve.CurveSequence;
CurveSequenceRep: TYPE ~ G3dCurve.CurveSequenceRep;
FSpot:     TYPE ~ RECORD [f: Fragment, t: REAL];
Creation
CurveFromSplines: PUBLIC PROC [splines: SplineSequence, res: NAT ¬ 100] RETURNS [c: Curve]
~ {
IF splines # NIL THEN {
tmp: Spline ¬ G3dMatrix.ObtainMatrix[];
nSegs: NAT ¬ MAX[1, res];
lenSegs: NAT ¬ MAX[1, nSegs/2];
c ¬ NEW[CurveRep[splines.length]];
c.length ¬ splines.length;
FOR s: NAT IN [0..splines.length) DO
t0: REAL ¬ 0.0;
spline: Spline ¬ splines[s];
fragments: Fragments ¬ NEW[G3dCurve.FragmentsRep[nSegs]];
c[s] ¬ [spline: spline, accumLength: c.totalLength, fragments: fragments];
FOR n: NAT IN [1..nSegs] DO
t1: REAL ¬ REAL[n]/REAL[nSegs];
len: REAL ¬ G3dSpline.Length[G3dSpline.Reparameterize[spline, t0, t1, tmp], lenSegs];
fragments[n-1] ¬ [t0, t1, t1-t0, len, c.totalLength];
c.totalLength ¬ c.totalLength+len;
t0 ¬ t1;
ENDLOOP;
c[s].length ¬ c.totalLength-c[s].accumLength;
ENDLOOP;
G3dMatrix.ReleaseMatrix[tmp];
};
};
Tranformation
Transform: PUBLIC PROC [curve: Curve, matrix: Matrix] ~ {
IF curve # NIL THEN FOR n: NAT IN [0..curve.length) DO
curve[n].spline ¬ G3dMatrix.Mul[curve[n].spline, matrix];
ENDLOOP;
};
Support
FSpotFromT: PROC [c: Curve, t: REAL] RETURNS [FSpot] ~ {
n: NAT ¬ Real.Floor[t];          -- index to curve section
tt: REAL ¬ t-n;            -- tt = parameter into section
fragments: Fragments ¬ c[MIN[n, c.length-1]].fragments;
nf: NAT ¬ Real.Round[tt*fragments.length];    -- guess at fragment
f: Fragment ¬ fragments[MIN[nf, fragments.length-1]];
WHILE tt NOT IN [f.t0..f.t1] DO        -- ensure correct fragment 
f ¬ fragments[nf ¬ nf+(IF tt < f.t0 THEN -1 ELSE 1)];
ENDLOOP;
RETURN[[f, tt]];
};
Conversions Between Arc-length and t
SpotFromT: PUBLIC PROC [c: Curve, t: REAL] RETURNS [Spot] ~ {
n: NAT ¬ Real.Floor[t];
RETURN[[c[n].spline, t-n]];
};
SpotFromArcLength: PUBLIC PROC [c: Curve, arcLength: REAL] RETURNS [Spot ¬ []] ~ {
SectionSpot: PROC [section: Section] RETURNS [spot: Spot] ~ {
FragmentSpot: PROC [f: Fragment] RETURNS [Spot] ~ {
RETURN[[section.spline, f.t0+f.dt*((arcLength-f.accumLength)/f.length)]];
};
f: Fragments ¬ section.fragments;
spot.spline ¬ section.spline;
FOR nf: NAT IN [1..f.length) DO
IF f[nf].accumLength > arcLength THEN RETURN[FragmentSpot[f[nf-1]]];
ENDLOOP;
RETURN[FragmentSpot[f[f.length-1]]];
};
IF arcLength < 0.0 THEN Error[$BadArcLength];
IF arcLength > c.totalLength THEN {
s: Section ¬ c[c.length-1];
f: Fragment ¬ s.fragments[s.fragments.length-1];
IF arcLength-c.totalLength > f.length*0.1  -- arbitrary
THEN Error[$BadArcLength]
ELSE RETURN[[s.spline, f.t1]];
};
FOR n: NAT IN [1..c.length) DO
IF c[n].accumLength > arcLength THEN RETURN[SectionSpot[c[n-1]]];
ENDLOOP;
RETURN[SectionSpot[c[c.length-1]]];
};
TFromArcLength: PUBLIC PROC [c: Curve, arcLength: REAL] RETURNS [t: REAL] ~ {
RETURN[SpotFromArcLength[c, arcLength].t];
};
ArcLengthFromT: PUBLIC PROC [c: Curve, t: REAL] RETURNS [REAL] ~ {
fs: FSpot ¬ FSpotFromT[c, t];
RETURN[fs.f.accumLength+((fs.t-fs.f.t0)/fs.f.dt)*fs.f.length];
};
ArcLengthFromTs: PUBLIC PROC [c: Curve, t0, t1: REAL] RETURNS [REAL] ~ {
RETURN[ArcLengthFromT[c, t1]-ArcLengthFromT[c, t0]];
};
Evaluation Given Parameter t
PositionT: PUBLIC PROC [c: Curve, t: REAL] RETURNS [Triple] ~ {
s: Spot ¬ SpotFromT[c, t];
RETURN[G3dSpline.Position[s.spline, s.t]];
};
VelocityT: PUBLIC PROC [c: Curve, t: REAL] RETURNS [Triple] ~ {
s: Spot ¬ SpotFromT[c, t];
RETURN[G3dSpline.Position[s.spline, s.t]];
};
AccelerationT: PUBLIC PROC [c: Curve, t: REAL] RETURNS [Triple] ~ {
s: Spot ¬ SpotFromT[c, t];
RETURN[G3dSpline.Acceleration[s.spline, s.t]];
};
CurvatureT: PUBLIC PROC [c: Curve, t: REAL] RETURNS [Triple] ~ {
s: Spot ¬ SpotFromT[c, t];
RETURN[G3dSpline.Curvature[s.spline, s.t]];
};
CurvatureMagT: PUBLIC PROC [c: Curve, t: REAL] RETURNS [REAL] ~ {
s: Spot ¬ SpotFromT[c, t];
RETURN[G3dSpline.CurvatureMag[s.spline, s.t]];
};
Evaluation Given Arc Length
PositionA: PUBLIC PROC [c: Curve, arcLength: REAL] RETURNS [Triple] ~ {
s: Spot ¬ SpotFromArcLength[c, arcLength];
RETURN[G3dSpline.Position[s.spline, s.t]];
};
VelocityA: PUBLIC PROC [c: Curve, arcLength: REAL] RETURNS [Triple] ~ {
s: Spot ¬ SpotFromArcLength[c, arcLength];
RETURN[G3dSpline.Velocity[s.spline, s.t]];
};
AccelerationA: PUBLIC PROC [c: Curve, arcLength: REAL] RETURNS [Triple] ~ {
s: Spot ¬ SpotFromArcLength[c, arcLength];
RETURN[G3dSpline.Acceleration[s.spline, s.t]];
};
CurvatureA: PUBLIC PROC [c: Curve, arcLength: REAL] RETURNS [Triple] ~ {
s: Spot ¬ SpotFromArcLength[c, arcLength];
RETURN[G3dSpline.Curvature[s.spline, s.t]];
};
CurvatureMagA: PUBLIC PROC [c: Curve, arcLength: REAL] RETURNS [REAL] ~ {
s: Spot ¬ SpotFromArcLength[c, arcLength];
RETURN[G3dSpline.CurvatureMag[s.spline, s.t]];
};
IO
Write: PUBLIC PROC [out: STREAM, c: Curve, name: ROPE ¬ NIL] ~ {
IF name # NIL THEN name ¬ Rope.Concat[name, " "];
IO.PutF[out, "Curve%g, %g sections:\n", IO.rope[name], IO.int[c.length]];
FOR n: NAT IN [0..c.length) DO
spline: Spline ¬ c[n].spline;
IO.PutRope[out, "\t"];
FOR i: NAT IN [0..4) DO
r: ARRAY [0..4) OF REAL ¬ spline[i];
IO.PutFL[out, "%g %g %g %g ",
LIST[IO.real[r[0]], IO.real[r[1]], IO.real[r[2]], IO.real[r[3]]]];
ENDLOOP;
IO.PutRope[out, "\n"];
ENDLOOP;
};
Read: PUBLIC PROC [in: STREAM] RETURNS [c: Curve] ~ {
ENABLE IO.Error, Convert.Error => GOTO Bad;
line: ROPE ¬ IO.GetLineRope[in];
p1: INT ¬ Rope.Find[line, "sections:"];
p0: INT ¬ Rope.FindBackward[line, " ", p1-2];
nSections: INT ¬ Convert.IntFromRope[Rope.Substr[line, p0, p1]];
splines: SplineSequence ¬ NEW[G3dSpline.SplineSequenceRep[nSections]];
FOR n: NAT IN [0..splines.length ¬ nSections) DO
splines[n] ¬ NEW[G3dSpline.SplineRep];
FOR i: NAT IN [0..4) DO
splines[n][i] ¬ [IO.GetReal[in], IO.GetReal[in], IO.GetReal[in], IO.GetReal[in]];
ENDLOOP;
ENDLOOP;
RETURN[CurveFromSplines[splines]];
EXITS Bad => Error[$BadInputFormat, IO.PutFR1["file position: %g", IO.int[IO.GetIndex[in]]]];
};
Sequences
CopyCurveSequence: PUBLIC PROC [curves: CurveSequence]
RETURNS [CurveSequence] ~ {
copy: CurveSequence ¬ NIL;
IF curves # NIL THEN {
copy ¬ NEW[CurveSequenceRep[curves.length]];
copy.length ¬ curves.length;
FOR n: NAT IN [0..curves.length) DO copy[n] ¬ curves[n]; ENDLOOP;
};
RETURN[copy];
};
AddToCurveSequence: PUBLIC PROC [curves: CurveSequence, curve: Curve]
RETURNS [CurveSequence]
~ {
IF curves = NIL THEN curves ¬ NEW[CurveSequenceRep[1]];
IF curves.length = curves.maxLength THEN curves ¬ LengthenCurveSequence[curves];
curves[curves.length] ¬ curve;
curves.length ¬ curves.length+1;
RETURN[curves];
};
LengthenCurveSequence: PUBLIC PROC [curves: CurveSequence, amount: REAL ¬ 1.3]
RETURNS [new: CurveSequence] ~ {
newLength: CARDINAL ¬ Real.Ceiling[amount*curves.maxLength];
newLength ¬ MAX[newLength, 3];
new ¬ NEW[CurveSequenceRep[newLength]];
FOR i: NAT IN [0..curves.length) DO new[i] ¬ curves[i]; ENDLOOP;
new.length ¬ curves.length;
};
END.