hawley@adobe.COM (Steve Hawley) (02/21/91)
Here is the spline code (originally written by Theo Pavlidis at Bell Labs)
optimized to a degree, but with a reasonable amount of readability (read:
it can still be better).
/*
* spline(p, n)
* draws an approximate bicubic spline through the control points in the
* array p. Don't use the first and last slots -- they're used to weight
* the endpoints.
* The previous version called for an offset. None has been included here.
* Also, the previous version had a kind of screwy notion as to what n was.
* Here, n refers to the total number of control points. If n is 4, the
* array p must have 6 entries in it. This 0th and 5th slots are used to
* weight control points, slots 1-4 inclusive are the 4 control points.
*/
/*
* Scaling factor and the number of subdivisions per curve. For best
* performance, SCALE should be power of 2 and STEPS should evenly
* divide SCALE. Also SCALE should always be greater than STEPS.
*/
#define SCALE 1024L
#define STEPS 8
spline(p, n)
register Point *p;
short n;
{
register long w, t1, t2, t3;
register short i, j;
p[0] = p[1];
MoveTo(p->h, p->v);
p[n+1] = p[n];
for (i=0; i < n; i++, p++) {
for (j=0; j < STEPS; j++) {
/*
* calculate the three terms for the polynomial
*/
w = (SCALE/STEPS) * j;
t1 = w*w / (2 * SCALE);
w -= SCALE/2;
t2 = 3*SCALE/4 - w*w / SCALE;
w -= SCALE/2;
t3 = w*w / (2 * SCALE);
LineTo((t1 * (p+2)->h + t2 * (p+1)->h +
t3 * p->h + SCALE/2) / SCALE,
(t1 * (p+2)->v + t2 * (p+1)->v +
t3 * p->v + SCALE/2) / SCALE);
/* The following call to Line() will redraw the
* last point in the line segment. This allows
* spline() to be used with patXor modes. If you're
* not going to use xor drawing modes, take this
* call out.
*/
Line(0, 0);
}
}
}
typical usage:
static Point k[7] = {
{0, 0}, {50, 50}, {75, 100}, {75, 120}, {200, 30}, {85, 150},
{0, 0} };
...
PenMode(patXor);
spline(k, 5);
...
Steve Hawley
hawley@adobe.com
--
"Did you know that a cow was *MURDERED* to make that jacket?"
"Yes. I didn't think there were any witnesses, so I guess I'll have to kill
you too." -Jake Johansenspencer@eecs.umich.edu (Spencer W. Thomas) (02/23/91)
Yes, you can do this, but it's more expensive than just computing a
priori how many points you need for a given spline segment and then
just blasting them out. A sketch of the code appears below.
The second method involves looking at the curvature (second
derivative) and using that to estimate the number of segments needed.
Sederberg and Parry ("Comparison of three curve intersection
algorithms," CAD, 1986 (sorry, I don't have vol/no reference), pp
58-63) attribute the following method to G Wang (published in Zhejiang
U Journal). Given a Bezier curve (you can always convert a polynomial
spline to a collection of Bezier curves by subdividing at the knots)
with control points (x[i],y[i]) i=0..n, compute
L0 = max( | x[i] - 2x[i+1]+ x[i+2] |, | y[i] - 2y[i+1] + y[i+2] | )
i = 0 .. n-2
then
r0 = log4(sqrt(2)*n*(n-1)*L0 / 8*epsilon)
is the number of segments to draw to stay within epsilon distance of
the curve.
/* Draw a spline (or Bezier curve) with recursive adaptive subdivision. */
draw_spline(s)
spline *s;
{
spline *s1, *s2;
if ( is_straight(s) )
draw_seg( first_pt(s), last_pt(s) );
else
{
subdivide( s, &s1, &s2 )
draw_spline( s1 );
draw_spline( s2 );
}
}
is_straight(s)
spline *s;
{
l = line( first_pt(s), last_pt(s) );
for (p = each control point in s)
{
/* Actually, you want to make sure the curve doesn't
"double back", either. The usual method is to
project each of the control points onto the line,
and make sure that they appear monotonically.
(I.e., the 't' values for the projections should
be in increasing or decreasing sequence.) I leave
this as an exercise.
*/
if ( distance( p, l ) > EPSILON )
return FALSE;
}
return TRUE;
}
subdivide( s, s1, s2 )
spline *s;
spline **s1, **s2;
{
/* Use a spline subdivision algorithm to split s into two
halves, usually at the midpoint of its parametric range.
If using Bezier curves instead of splines, just split it
in the middle. The first half is returned via s1 and the
second half via s2 (the storage is dynamically allocated
in this "implementation".)
*/
}
--
=Spencer W. Thomas EECS Dept, U of Michigan, Ann Arbor, MI 48109
spencer@eecs.umich.edu 313-936-2616 (8-6 E[SD]T M-F)