[_CDMiwok_01_]<cedarcommon2.0>Fit>NodesImpl.mesa!1

NodeImpl.mesa

last edited by Michael Plass August 30, 1982 2:06 pm

last edited by Maureen Stone September 21, 1984 11:54:41 am PDT

Doug Wyatt, September 5, 1985 2:34:29 pm PDT

Maureen Stone, September 28, 1987 4:49:51 pm PDT

DIRECTORY

Complex USING [Abs, Add, SqrAbs, Sub, VEC],

Cubic2 USING [Bezier, BezierToCoeffs, Coeffs],

DynFit USING [FitSegments, Segments],

LSPiece USING [FitPiece, Metrics],

Nodes USING [Progress],

RealFns USING [SinDeg, SqRt],

Seq USING [ComplexSequence, ComplexSequenceRec, NatSequence, NatSequenceRec, RealSequence, RealSequenceRec],

Vector USING [Cross, Dot, Mul];

NodesImpl: CEDAR PROGRAM

IMPORTS Complex, Cubic2, DynFit, LSPiece, Vector, RealFns

EXPORTS Nodes =

BEGIN

DynNodes: PUBLIC PROC [samples: Seq.ComplexSequence, closed: BOOLEAN, penalty: REAL] RETURNS [nodes: Seq.NatSequence] = {

finds nodes by fitting a polygon using DynFit

endpoints: DynFit.Segments;

[segments: endpoints] ← DynFit.FitSegments[samples, closed, penalty];

nodes ← NEW[Seq.NatSequenceRec[endpoints.length]];

nodes[0] ← 0;

FOR i:NAT IN [0..endpoints.length) DO

nodes[i] ← endpoints[i];

ENDLOOP;

};

CubicTangents: PUBLIC PROC [samples: Seq.ComplexSequence, closed: BOOLEAN, metrics: LSPiece.Metrics, nodes: Seq.NatSequence ← NIL] RETURNS [tangents: Seq.ComplexSequence] = {

Sets tangents at nodes by fitting a cubic between neighboring nodes

Calls ICubicTangents and packs the results in a sequence

i: NAT ← 0;

progress: Nodes.Progress = {

tangents[i] ← tangent;

i ← i+1;

RETURN[FALSE];

};

tangents ← IF nodes#NIL

THEN NEW[Seq.ComplexSequenceRec[nodes.length]]

ELSE NEW[Seq.ComplexSequenceRec[samples.length]];

ICubicTangents[progress, samples, closed, metrics, nodes];

};

Progress: PROC[tangent: Complex.VEC, cubic: Cubic2.Bezier] RETURNS [stop: BOOLEAN];

ICubicTangents: PUBLIC PROC [progress: Nodes.Progress, samples: Seq.ComplexSequence, closed: BOOLEAN, metrics: LSPiece.Metrics, nodes: Seq.NatSequence ← NIL] = {

Sets tangents at nodes by fitting a cubic between neighboring nodes

n: NAT ← samples.length;

t: Seq.RealSequence ← NEW[Seq.RealSequenceRec[n]];

lastTangent: Complex.VEC;

bezier: Cubic2.Bezier;

prev: PROC [index: NAT] RETURNS [NAT] = {

IF closed THEN RETURN[IF index=0 THEN nodes.length-1 ELSE index-1]

ELSE RETURN[index];

};

next: PROC [index: NAT] RETURNS [NAT] = {

IF closed THEN RETURN[IF index=nodes.length-1 THEN 1 ELSE index+1]

ELSE RETURN[index];

};

pointsBetween: PROC [from, to: NAT] RETURNS [NAT] = {

IF closed AND to<from THEN RETURN[samples.length-from+to-1]

ELSE RETURN[to-from-1];

};

FOR i: NAT IN [0..nodes.length) DO

p: NAT ← nodes[prev[i]];

q: NAT ← nodes[i];

r: NAT ← nodes[next[i]];

nPoints: NAT ← pointsBetween[p,r];

e: REAL;

IF nPoints < minPoints THEN

lastTangent ← CheapTangent[samples[p], samples[q], samples[r], maxSin]

ELSE {

[b:bezier, maxDev: e] ← LSPiece.FitPiece[z: samples, t: t, from: p, thru: r,

metrics: metrics, initFree: freeEnds, finalFree: freeEnds];

lastTangent ←IF e > metrics.sumErr THEN [0,0] ELSE TangentVector[bezier,t[q]];

};

IF progress[lastTangent,bezier] THEN EXIT;

ENDLOOP;

};

freeEnds: BOOLEAN ← TRUE;

minPoints: NAT ← 4;

maxSin: REAL ← 0.5;

TangentVector: PROC [b: Cubic2.Bezier, t: REAL] RETURNS [tang: Complex.VEC] = {

c: Cubic2.Coeffs ← Cubic2.BezierToCoeffs[b];

tang.x ← (3*c.c3.x*t+2*c.c2.x)*t + c.c1.x;

tang.y ← (3*c.c3.y*t+2*c.c2.y)*t + c.c1.y;

};

QuickTangents: PUBLIC PROC [samples: Seq.ComplexSequence, closed: BOOLEAN, maxAngle: REAL, nodes: Seq.NatSequence ← NIL] RETURNS [tangents: Seq.ComplexSequence] = {

computes tangents by differencing neighbors

tangents ← DiffTangents[samples, closed, maxAngle, nodes, CheapTangent];

};

SquareTangents: PUBLIC PROC [samples: Seq.ComplexSequence, closed: BOOLEAN, maxAngle: REAL, nodes: Seq.NatSequence ← NIL] RETURNS [tangents: Seq.ComplexSequence] = {

computes tangents by differencing neighbors

tangents ← DiffTangents[samples, closed, maxAngle, nodes, SquareTangent];

};

TangentProc: TYPE = PROC[a, b, c: Complex.VEC, maxSin: REAL] RETURNS [t: Complex.VEC];

DiffTangents: PROC [samples: Seq.ComplexSequence, closed: BOOLEAN, maxAngle: REAL, nodes: Seq.NatSequence, tangentProc: TangentProc] RETURNS [tangents: Seq.ComplexSequence]= {

maxSin: REAL ← RealFns.SinDeg[maxAngle];

tangents ← NEW[Seq.ComplexSequenceRec[nodes.length]];

tangents[0] ← IF closed THEN

CheapTangent[samples[nodes[nodes.length-1]], samples[nodes[0]], samples[nodes[1]], maxSin]

ELSE [0,0];

FOR i: NAT IN [1..nodes.length-1) DO

tangents[i] ← CheapTangent[samples[nodes[i-1]], samples[nodes[i]], samples[nodes[i+1]], maxSin];

ENDLOOP;

tangents[nodes.length-1] ← IF closed THEN

CheapTangent[samples[nodes[nodes.length-2]], samples[nodes[nodes.length-1]], samples[nodes[0]], maxSin]

ELSE [0,0];

};

CheapTangent: PROC [a, b, c: Complex.VEC, maxSin: REAL] RETURNS [t: Complex.VEC] = {

d1: Complex.VEC ← Complex.Sub[b,a];

d2: Complex.VEC ← Complex.Sub[c,b];

absd1d2: REAL ← RealFns.SqRt[Complex.SqrAbs[d1]*Complex.SqrAbs[d2]];

IF Vector.Dot[d1, d2] < 0 OR ABS[Vector.Cross[d1, d2]] >= maxSin * absd1d2 THEN t ← [0,0]

ELSE t ← Complex.Sub[c,a];

};

SquareTangent: PROC [a, b, c: Complex.VEC, maxSin: REAL] RETURNS [t: Complex.VEC] = {

d1: Complex.VEC ← Complex.Sub[b,a];

d2: Complex.VEC ← Complex.Sub[c,b];

absd1d2: REAL ← RealFns.SqRt[Complex.SqrAbs[d1]*Complex.SqrAbs[d2]];

IF Vector.Dot[d1, d2] < 0 OR ABS[Vector.Cross[d1, d2]] >= maxSin * absd1d2 THEN t ← [0,0]

ELSE t ← Complex.Add[Vector.Mul[d1, Complex.Abs[d1]], Vector.Mul[d2, Complex.Abs[d2]]];

};

END.

Michael Plass, August 30, 1982 2:06 pm. Added CubicTangents.

SetTangents: PROC = {-- range err maxit => . Sets tangents at nodes by locally fitting 2*range + 1 samples

maxit: NAT ← PopInteger[];

err: REAL ← GetReal[];

range: NAT ← PopInteger[];

node,cusp: Seq.NatSequence;

z: Seq.ComplexSequence ← Curve.CurrentSamples[Curve.defaultHandle];

n: NAT ← z.length;

t: Seq.RealSequence ← NEW[Seq.RealSequenceRec[n]];

bezier: Cubic2.Bezier;

nextCusp,cindex: INT ← 0;

forward: PROC[pt: NAT, range: INTEGER] RETURNS[INTEGER] = {

np,nc: INTEGER;

np ← (pt+range) MOD z.length;

IF cusp.length = 0 THEN RETURN[np]

ELSE {

i: INTEGER ← IF pt=nextCusp THEN (cindex+1) MOD cusp.length ELSE cindex;

nc ← cusp[i];

RETURN[MIN[np,nc]];

};

back: PROC[pt: NAT, range: INTEGER] RETURNS[INTEGER] = {

np,nc: INTEGER;

np ← pt-range;

IF np<0 THEN np ← z.length+np;

IF cusp.length=0 THEN RETURN[np]

ELSE {

nc ← IF cindex-1 <0 THEN cusp[cusp.length-1] ELSE cusp[cindex-1];

RETURN[MAX[np,nc]];

};

[nodes:node] ← Curve.CurrentNodes[Curve.defaultHandle];

[cusps: cusp] ← Curve.CurrentCusps[Curve.defaultHandle];

IF cusp.length=0 THEN nextCusp ← -1 ELSE nextCusp ← cusp[0];

FOR i: NAT IN [0..node.length) DO

q: NAT ← node[i];

in: Complex.VEC;

IF q=nextCusp THEN {

p: INTEGER ← back[q,2*range];

r: INTEGER ← forward[q,2*range];

[b:bezier] ← LSPiece.FitPiece[z: z, t: t, from: p, thru: q,

eps: .0001, maxd: err, maxit: maxit, initFree: TRUE, finalFree: TRUE];

in ← TangentVector[bezier,t[q]];

[b:bezier] ← LSPiece.FitPiece[z: z, t: t, from: q, thru: r,

eps: .0001, maxd: err, maxit: maxit, initFree: TRUE, finalFree: TRUE];

Curve.AddCusp[Curve.defaultHandle, q,in,TangentVector[bezier,t[q]]];

cindex ← cindex+1;

nextCusp ← IF cindex < cusp.length THEN cusp[cindex] ELSE -1;

}

ELSE {

p: INTEGER ← back[q,range];

r: INTEGER ← forward[q,range];

[b:bezier] ← LSPiece.FitPiece[z: z, t: t, from: p, thru: r,

eps: .0001, maxd: err, maxit: maxit, initFree: TRUE, finalFree: TRUE];

Curve.AddNode[Curve.defaultHandle, q,TangentVector[bezier,t[q]]];

};

IF GetJaMBreak[] THEN EXIT;

ENDLOOP;

};