Two-dimensional curves show up in various applications: EPUBS, MCAD, ECAD, Postscript, and window systems, to name a few. To be useful, curves must be easy to control and easy for a computer to render. Controlling lines is simple: The user sets some control geometry (the endpoints) and knows intuitively where the line will be drawn. Controlling a curve is more difficult, since there is an infinite number of curved shapes between any two endpoints. Restricting the available curves to a fixed set of curve types, such as circles and ellipses, solves the problem, but this solution is not powerful enough for most applications. Free-form curves, whose shape the user can control without restriction, are required. Unfortunately, there is no single, obvious way to create and shape them.
Parametric curve types are distinguished by the type of control geometry describing the curve, and how it is used to generate the curve equation. A different curve type, even if it uses the same control geometry, will be interpreted into a different curve shape. The computer-graphics community makes use of a variety of curve types, trading off their different strengths and weaknesses. For example, a fast and simple Bezier curve may be ideal for representing fonts in a Post-Script printer, while a more expressive NURBS curve would better represent the complex shapes created with a solids modeling application. As a result, a sophisticated application may have to handle many different curve representations. In this article, we show how to represent a wide variety of curve representations efficiently, implementing them in C++, using a class hierarchy and an object-oriented programming style.
The formulas describing free-form curves have much in common. In fact, only two major types of curves are in common use: exterpolating and interpolating curves. They differ in how they respond to the control geometry. An application uses a set of geometric locations called "control points" to define the shape of the curve. For exterpolating curves, the control points provide a boundary called the "convex hull." The curve always remains within this boundary, sometimes not even touching the control points. Interpolative curves, however, always pass through their control points.
Both types of curves can be represented in a class/object design: The class represents free-form curves containing the control-points array and the operators that act on them. The curve object hierarchy has three features:
Figure 1 shows the class hierarchy for curves, as represented in this article. The base class is derived into two primary subclasses: the basis-matrix class and the nonuniform B-spline class. A basis-matrix curve is represented by a matrix derived using the techniques described in the accompanying text box, "Derivation of Basis Matrix." The nonuniform B-spline class is used for curves that are nonuniform in the step of the parametric variable t.
Listing One (page 60) shows the C++ header file for the curve-class hierarchy. The Basis_matrix_curve class is used for the uniform beta-spline and to derive the curve types Hermite, Bezier, uniform B-spline, and Catmull-Rom interpolating curve. The Nub_curve class defines the nonuniform, nonrational B-splines and to derive the curve type nonuniform, rational B-spline (commonly known as NURBS). Nonuniform B-splines require an extra piece of data called the "knot vector," a floating-point array containing a nondecreasing list of values that control how the curve is evaluated. Figure 2 shows the formulas used to define nonuniform B-splines.
The implementation of curves involves defining the methods that use the control points to render the curve. For the Basis_matrix_curve, the constructor defines the matrix used to compute the coefficients of the third-degree polynomial curve definition. So, to define a new type of curve, it is necessary to derive the new types from the Basis_matrix_curve class and define a constructor that computes the new basis matrix.
Note that the curves are two-dimensional; this is defined in the class Point2d. It is trivial to extend the definition to three-dimensional curves by defining a Point3d. The curves in this article are limited to fourth order or third degree. The challenge of extending the classes to other orders is left to the reader.
For a Hermite, Bezier, uniform B-spline, and Catmull-Rom spline, the constructor defines the basis matrix used to compute the coefficients of the polynomial. But for a uniform beta-spline, the bias and tension values must be initialized. The default values are bias= 1.0 and tension= 0.0. These values reduce to a basis matrix equivalent to the uniform B-spline. You can exert precise control of the curve by modifying the bias and tension parameters of the uniform beta-spline. Listing Two (page 60) shows the methods, including the constructors, for each of the various curves.
When an application creates a curve object, it must supply the control points that define the curve shape. For some curve types, the application may also have to define the know vector, or tension and bias values. The application then invokes the display_curve method to render the curve.
All curves derived from the Curve class have a display_curve method. For the Basis_matrix_curve class, the display-curve method renders a third-degree polynomial by tessellating it into vectors. But for a nonuniform B-spline class, the display_curve method is overwritten to display the curve, using the formulas in Figure 3.
For the Basis_matrix_curve class, the curve is converted to the coefficients of a third-degree polynomial. This polynomial is evaluated into a number of line segments, based on a fixed tessellation factor. The curve's tessellation into line segments is performed within the class's display_curve method.
There are two popular methods for evaluating a third-degree polynomial. The first technique is based on Horner's method. It is a method of rewriting the polynomial to reduce the number of arithmetic operations associated with its evaluation. Figure 3 shows the formulas behind Horner's method, and Listing Three (page 62) shows the C++ source code for displaying a third-degree polynomial using Horner's method. The code reduces the number of multiplies and adds so that we gain some performance, but we can do better.
A higher-performance method for evaluating polynomials surfaced in the 1970s. This method is based on evaluating the first- and second-order differences of the polynomial and performing additions within the tessellation inner loop. Again, Figure 3 shows the formulas for this technique, and Listing Four (page 62) shows the C++ source code for displaying a curve using forward differencing. The inner loop of the curve display is reduced to six additions--certainly faster than the eight multiplications and seven additions performed by Horner's method.
Listing Five (page 63) shows the main routine for displaying several different curve types. Each curve displayed defines the geometry of the control points and then invokes the display-curve method.
Using the same control points, the demonstration main routine displays several different Hermite curves. The curves' shape is modified by changing the direction vectors at the first and last control points, and modifying the vectors' magnitude and direction.
A Bezier curve is rendered for the control points (20, 20), (50, 100), (300, 50), and (100, 10). This defines a simple Bezier-style curve. Immediately after rendering the Bezier curve, a non-uniform, nonrational B-spline (NUB) is rendered. The knot vector of this curve is set to (0, 0, 0, 0, 1, 1, 1, 1). This interpolates the endpoints and extrapolates the interior control points, thus displaying the same Bezier curve. The knot vector of the NUB curve is then modified to (0, 0, 0, 1, 2, 3, 3, 3) and rendered. This curve shows the extrapolation of the control points. A simple modification to the knot vector yields a completely different curve.
The next curve displayed is a Catmull-Rom curve which interpolates the control points. It is useful in applications such as data analysis, where the control points must lie on the curve.
The demo program renders several uniform, beta-spline curves and manipulates the bias and tension parameters to demonstrate how the shape of the curve can be easily changed. Increasing the bias parameter pulls the curve to a sharp angle; while increasing the tension parameter yields a similar effect but in exactly the opposite direction.
Finally, the demonstration program renders several NURB curves through the same control points. The knot vector is manipulated to render several different curves using the same control points.
Listing Six (page 63) is the source file utilstc.c that defines all the Borland C++-specific rendering routines. Listing Seven (page 64), utilsxlib.c, defines all the rendering routines specific to machines running the X Window system. The current implementation is compiled for Borland C++ 3.0 and for UNIX/X on a Sun and IBM RS/6000.
Porting this program to another system is easy: The source code needs a C++ compiler and a display system. You write your own, system-specific version of the display-utilities source file shown in Listings Six and Seven.
To port this code, modify init_graphics_device to open and initialize your graphics device. This code should create and map windows to the display and clear the window or display to a background color. The function close_graphics_device must be modified to destroy the window or reset the hardware to the appropriate state before returning to the system.
If your hardware can scan-convert a line segment, modify the source code in the line function. Just add the appropriate code to call the system line-segment function similar to the Xlib port. If your system cannot scan-convert a line segment, you can change the set_pixel function to write the appropriate color into your display device. If your system cannot display a pixel at a given x, y position with a given color, you are out of luck.
The clear_window function sets the window or display to a background color. The VGA implementation sets the background to pixel-value 0, which is black on a standard display. The function text_output prints a string on the display. The demonstration program uses these functions to help you view the curves rendered.
Foley, J.D. et al. Computer Graphics Principles and Practice. Reading, MA: Addison-Wesley, 1990.
Hearn, Donald and Pauline M. Baker. Computer Graphics. Englewood Cliffs, NJ: Prentice Hall, 1986.
Nye, Adrian. Xlib Reference Manual for Version 11. Sebastapol, CA: O'Reilly & Associates, 1990.
Phoenix Technical Reference Series: System BIOS for IBM PC/XT/AT Computers and Compatibles. Reading, MA: Addison-Wesley, 1989.
Rodgers, David F. and J. Alan Adams. Mathematical Elements for Computer Graphics, second edition. Berkeley, CA: McGraw-Hill, 1990.
Turbo C++ Library Reference. Scotts Valley, CA: Borland International, 1990.
Wilton, Richard. Programmer's Guide to PC and PS/2: Video Systems. Redmond, WA: Microsoft Press, 1987.
In order to make them more manageable mathematically, the equations that draw curves are written parametrically. The parametric form of an equation is written as a function of one or more independent variables (in these examples, the single variable t). For example, the line equation; y= mx+b would be split into two separate equations, both depending on t: x(t)=Ax + Bxt and y(t)= Ay+ Byt. In this representation, when t=0, x and y are at one endpoint of the line; when t=1, and x and y are at the other.
Although the most common parametric curves are fourth order (containing t3, t2, t, and 1), the examples here use second order, containing only t and 1. This restricts our "curves" to straight lines.
Although a parametric equation is useful, it's hard to find coefficients needed to represent the desired curve. To make this easier, we'll break the equation up into pieces. First, we must tie the equation to some geometry that controls the type of curve drawn. The set of points that define the curve are the control points. In this example, we want the line drawn between two endpoints that we select. Assuming t ranges from 0 -- 1, we need to figure out Ax, Ay, Bx, and By in terms of the curve endpoints, which we'll call X0, Y0, and x1, y1. This leads to the equations in Example 1(a).
Now we have a parametric equation based on some understandable geometry. Example 1(b) shows this rewritten so there's only one instance of each geometry value.
Example 1(c) puts everything in matrix form. Splitting the last matrix into two parts isolates the geometry from the curve type, as in Example 1(d).
The equations are now composed of three matrices: parametric, basis, and geometry. The parametric is relatively fixed--the number of elements determines the order of the curve. The basis matrix determines how the geometry combines the geometry matrix with the equation's parameters to draw the line. The geometry part is different for every line.
To find the elements of the basis matrix, solve the matrix equation. Set the t values at 0, then 1; solving for each gives the equation in Example 1(e). The x and y values have been replaced with more general "geometry" entries. Multiplying the parameter and basis matrices together results in a set of blending functions. As the value of the parameter t changes, these equations "blend" the elements of the geometry matrix together to form curves.
Most applications contain curves made up of many connected segments. The places where the curve segments meet, as well as the ends of complete curve, are called "knots." To connect these segments smoothly, you must be able to control the connections continuity. A C0 continuity means the curve knots connect; C1 means the slopes of the curves match at the connection; and C2 means the curvatures match as well. Cn continuity is defined by the equation in Example 1(f), where the nth derivative of the curve segments match where they meet. The level of continuity possible is limited by the order of the curve used. To handle curves with multiple curve segments, we need to change our general equation to that in Example 1(g). In our example, the basis matrix didn't change. In most cases, it does. To keep the numbering the same, the is must be greater than 0. Now a curve is defined by a list of control points. Each curve segment ranges from ti to ti+1. Since this curve is uniform, ti+1-ti, the difference of t between knot values always equals 1. Note that t ranges from 0-1 in each curve segment.
Moving to order-four curves gives us two more levels of continuity to work with, allowing curve segments to be connected together smoothly. The equation for fourth-order parametric curves changes surprisingly little; see Example 1(h). Different curves are formed by defining the basis and geometry matrices. The basis matrix defines how the geometry will be selected. Since it interprets how the geometry should be mixed by the parameter t, the type of basis matrix depends on the type of geometry in the geometry matrix. The basis and geometry matrices for Hermite, Bezier, B-spline, and Catmull-Rom (an interpolating spline) are given in Example 1(i). A P in the geometry matrix indicates a control point, an R indicates a control vector.
The first three equations in Example 1(i) approximate the curves' control points; the Catmull-Rom equation approximates the curves themselves. The Beta-Spline representation is a generalization of the B-spline curve. It adds two parameters, bias (
1) and tension (2). If 1= 1 and 2= 0, then it reduces to the B-spline, as in Example 1(j). All the curves shown have been computed in 2-D (x and y) coordinates. To make 3-D curves, add a z coordinate. When rendering, convert back into nonrational by dividing each coordinate in a point by its w part, then discarding w.
--S.J. & T.M.
_IMPLEMENTING CURVES IN C++_ by Stephen P. Johnson and Tom McReynolds[LISTING ONE]
// curve.h - Base class for curves
#include <stdio.h>
class Point2d {
protected:
float x, y, w;
public:
// set x, y
void set_xy(float new_x, float new_y) {
x = new_x; y = new_y; w = 1.;
}
// set x, y and w
void set_xyw(float new_x, float new_y, float new_w) {
x = new_x; y = new_y; w = new_w;
}
void set_x(float new_x) { x = new_x; }
void set_y(float new_y) { y = new_y; }
void set_w(float new_w) { w = new_w; }
// get x, y
void get_xy(float *ret_x, float *ret_y) {
*ret_x = x; *ret_y = y;
}
// get x, y and w
void get_xy(float *ret_x, float *ret_y, float *ret_w) {
*ret_x = x; *ret_y = y; *ret_w = w;
}
float get_x(void) { return x; }
float get_y(void) { return y; }
float get_w(void) { return w; }
// print out the current x, y
void print(void) { printf("%f %f\n", x, y); }
};
class Curve {
protected:
int num_geom;
Point2d *geom;
public:
// constructor
Curve(void);
// destructor
~Curve(void);
// set the geometry vector which is also called the control points
void set_geom_vector(int count, Point2d *g);
// method to display a third degree curve
void display_curve(int n, int color, int show_geom_pts,
int show_convex_hull);
};
class Basis_matrix_curve : public Curve {
protected:
Point2d t_coeff[4];
float basis_matrix[4][4];
public:
// constructor
Basis_matrix_curve(void);
// get the coefficient matrix for the t's
void get_t_coeff(Point2d tc[4]);
// overwrite set geometry vector to add in matrix multiply
void set_geom_vector(int count, Point2d *g);
// set a user defined matrix into the basis matrix
void set_basis_matrix(float m[4][4]);
// get the current basis matrix
void get_basis_matrix(float m[4][4]);
// multiply a 4x4 by 4 points that are a 2x4 matrix
void multiply_basis_by_geometry(void);
// method to display a third degree curve
void display_curve(int n, int color, int show_geom_pts,
int show_convex_hull);
};
class Hermite_curve : public Basis_matrix_curve {
public:
// constructor for Hermite type curve
Hermite_curve(void);
};
class Bezier_curve : public Basis_matrix_curve {
public:
// constructor for Bezier type curve
Bezier_curve(void);
};
class Bspline_curve : public Basis_matrix_curve {
public:
// constructor for Bspline type curve
Bspline_curve(void);
};
class Catmull_Rom_curve : public Basis_matrix_curve {
public:
// constructor for Catmull-Rom type curve
Catmull_Rom_curve(void);
};
// uniformly shaped beta-spline
class Beta_spline_curve : public Basis_matrix_curve {
protected:
float bias;
float tension;
public:
// constructor for uniformly shaped beta-spline
Beta_spline_curve(void);
// methods for setting and getting the bias and tension
void set_bias(float new_bias) {
bias = new_bias;
update_basis_matrix();
multiply_basis_by_geometry();
}
void set_tension(float new_tension) {
tension = new_tension;
update_basis_matrix();
multiply_basis_by_geometry();
}
float get_bias(void) { return bias; }
float get_tension(void) { return tension; }
// method for updating the basis matrix
void update_basis_matrix(void);
};
class Nub_curve : public Curve {
protected:
int num_knots;
float *knots;
public:
// constructor for NUB curves
Nub_curve(void);
// destructor for NUB curves
~Nub_curve(void);
// setup the knot vector with user information
void set_knot_vector(int count, float *v);
// method to display a third degree NUB curve
void display_curve(int n, int color, int show_geom_pts,
int show_convex_hull);
};
class Nurb_curve : public Nub_curve {
public:
// method to display a third degree NURB curve
void display_curve(int n, int color, int show_geom_pts,
int show_convex_hull);
};
[LISTING TWO]
// curve.cc - Implementation of methods for base class for curves //
#include <iostream.h>
#include <stdio.h>
#include <stdlib.h>
#ifdef __TURBOC__
#include <alloc.h>
#else
#include <malloc.h>
#endif
#include <math.h>
#include "curve.h"
#define CROSS_SIZE 5
extern "C" void line(short, short, short, short, short);
// constructor for curve
Curve::Curve(void)
{
num_geom = 0;
geom = (Point2d *)NULL;
}
// destructor for curve
Curve::~Curve(void)
{
if (num_geom != 0 || geom) {
free(geom);
}
}
// set the geometry vector which is also called the control points
void Curve::set_geom_vector(int count, Point2d *g)
{
if (count != num_geom && geom) {
free(geom);
}
geom = (Point2d *)malloc(count * sizeof(Point2d));
if (geom) {
int i;
num_geom = count;
for (i = 0; i < num_geom; i++) {
geom[i] = g[i];
}
}
}
// constructor for basis matrix class
// get the coefficients
void Basis_matrix_curve::get_t_coeff(Point2d tc[4])
{
tc[0] = t_coeff[0];
tc[1] = t_coeff[1];
tc[2] = t_coeff[2];
tc[3] = t_coeff[3];
}
// overwrite set geometry vector to do basis matrix multiply
void Basis_matrix_curve::set_geom_vector(int count, Point2d *g)
{
Curve::set_geom_vector(count, g);
multiply_basis_by_geometry();
}
// constructor for basis matrix class
Basis_matrix_curve::Basis_matrix_curve(void)
{
int i, j;
Point2d g[4];
g[0].set_xy(0., 0.); g[1].set_xy(0., 0.);
g[2].set_xy(0., 0.); g[3].set_xy(0., 0.);
Curve::set_geom_vector(4, g);
t_coeff[0].set_xy(0., 0.);
t_coeff[1].set_xy(0., 0.);
t_coeff[2].set_xy(0., 0.);
t_coeff[3].set_xy(0., 0.);
for (i = 0; i < 4; i++)
for (j = 0; j < 4; j++)
if (i == j)
basis_matrix[i][j] = 1.;
else
basis_matrix[i][j] = 0.;
}
// set a user defined matrix into the basis matrix
void Basis_matrix_curve::set_basis_matrix(float m[4][4])
{
int i, j;
for (i = 0; i < 4; i++)
for (j = 0; j < 4; j++)
basis_matrix[i][j] = m[i][j];
}
// get the current basis matrix
void Basis_matrix_curve::get_basis_matrix(float m[4][4])
{
int i, j;
for (i = 0; i < 4; i++)
for (j = 0; j < 4; j++)
m[i][j] = basis_matrix[i][j];
}
// multiply a 4x4 by a 2x4 (2D geometric matrix)
void Basis_matrix_curve::multiply_basis_by_geometry()
{
t_coeff[0].set_x(
basis_matrix[0][0] * geom[0].get_x() +
basis_matrix[0][1] * geom[1].get_x() +
basis_matrix[0][2] * geom[2].get_x() +
basis_matrix[0][3] * geom[3].get_x());
t_coeff[1].set_x(
basis_matrix[1][0] * geom[0].get_x() +
basis_matrix[1][1] * geom[1].get_x() +
basis_matrix[1][2] * geom[2].get_x() +
basis_matrix[1][3] * geom[3].get_x());
t_coeff[2].set_x(
basis_matrix[2][0] * geom[0].get_x() +
basis_matrix[2][1] * geom[1].get_x() +
basis_matrix[2][2] * geom[2].get_x() +
basis_matrix[2][3] * geom[3].get_x());
t_coeff[3].set_x(
basis_matrix[3][0] * geom[0].get_x() +
basis_matrix[3][1] * geom[1].get_x() +
basis_matrix[3][2] * geom[2].get_x() +
basis_matrix[3][3] * geom[3].get_x());
t_coeff[0].set_y(
basis_matrix[0][0] * geom[0].get_y() +
basis_matrix[0][1] * geom[1].get_y() +
basis_matrix[0][2] * geom[2].get_y() +
basis_matrix[0][3] * geom[3].get_y());
t_coeff[1].set_y(
basis_matrix[1][0] * geom[0].get_y() +
basis_matrix[1][1] * geom[1].get_y() +
basis_matrix[1][2] * geom[2].get_y() +
basis_matrix[1][3] * geom[3].get_y());
t_coeff[2].set_y(
basis_matrix[2][0] * geom[0].get_y() +
basis_matrix[2][1] * geom[1].get_y() +
basis_matrix[2][2] * geom[2].get_y() +
basis_matrix[2][3] * geom[3].get_y());
t_coeff[3].set_y(
basis_matrix[3][0] * geom[0].get_y() +
basis_matrix[3][1] * geom[1].get_y() +
basis_matrix[3][2] * geom[2].get_y() +
basis_matrix[3][3] * geom[3].get_y());
}
// constructor for Hermit curve
Hermite_curve::Hermite_curve(void)
{
float m[4][4];
m[0][0] = 2.; m[0][1] = -2.; m[0][2] = 1.; m[0][3] = 1.;
m[1][0] = -3.; m[1][1] = 3.; m[1][2] = -2.; m[1][3] = -1.;
m[2][0] = 0.; m[2][1] = 0.; m[2][2] = 1.; m[2][3] = 0.;
m[3][0] = 1.; m[3][1] = 0.; m[3][2] = 0.; m[3][3] = 0.;
set_basis_matrix(m);
}
// constructor for Bezier curve
Bezier_curve::Bezier_curve(void)
{
float m[4][4];
m[0][0] = -1.; m[0][1] = 3.; m[0][2] = -3.; m[0][3] = 1.;
m[1][0] = 3.; m[1][1] = -6; m[1][2] = 3.; m[1][3] = 0.;
m[2][0] = -3.; m[2][1] = 3.; m[2][2] = 0.; m[2][3] = 0.;
m[3][0] = 1.; m[3][1] = 0.; m[3][2] = 0.; m[3][3] = 0.;
set_basis_matrix(m);
}
// constructor for Uniform Nonrational Bspline curve
Bspline_curve::Bspline_curve(void)
{
float m[4][4];
m[0][0] = -1./6.; m[0][1] = 3./6.; m[0][2] = -3./6.; m[0][3] = 1./6.;
m[1][0] = 3./6.; m[1][1] = -6./6.; m[1][2] = 3./6.; m[1][3] = 0.;
m[2][0] = -3./6.; m[2][1] = 0.; m[2][2] = 3./6.; m[2][3] = 0.;
m[3][0] = 1./6.; m[3][1] = 4./6.; m[3][2] = 1./6.; m[3][3] = 0.;
set_basis_matrix(m);
}
Catmull_Rom_curve::Catmull_Rom_curve(void)
{
float m[4][4];
m[0][0] = -1./2.; m[0][1] = 3./2.; m[0][2] = -3./2.; m[0][3] = 1./2.;
m[1][0] = 2./2.; m[1][1] = -5./2.; m[1][2] = 4./2.; m[1][3] = -1./2.;
m[2][0] = -1./2.; m[2][1] = 0.; m[2][2] = 1./2.; m[2][3] = 0.;
m[3][0] = 0.; m[3][1] = 2./2.; m[3][2] = 0.; m[3][3] = 0.;
set_basis_matrix(m);
}
Beta_spline_curve::Beta_spline_curve(void)
{
bias = 1.;
tension = 0.;
update_basis_matrix();
}
void Beta_spline_curve::update_basis_matrix(void)
{
int i, j;
float m[4][4];
float bias2 = bias*bias;
float bias3 = bias2*bias;
float delta;
m[0][0] = -2.*bias3; m[0][1] = 2.*(tension+bias3+bias2+bias);
m[0][2] = -2.*(tension+bias2+bias+1.); m[0][3] = 2.;
m[1][0] = 6.*bias3; m[1][1] = -3.*(tension+2.*bias3+2.*bias2);
m[1][2] = 3.*(tension+2.*bias2); m[1][3] = 0.;
m[2][0] = -6.*bias3; m[2][1] = 6.*(bias3-bias);
m[2][2] = 6.*bias; m[2][3] = 0.;
m[3][0] = 2.*bias3; m[3][1] = tension+4.*(bias2+bias);
m[3][2] = 2.; m[3][3] = 0.;
delta = tension + 2.*bias3 + 4*bias2 + 4*bias + 2.;
for (i = 0; i < 4; i++)
for (j = 0; j < 4; j++)
m[i][j] /= delta;
set_basis_matrix(m);
}
float b1(int i, float t, float *knots)
{
if (knots[i] <= 1. && t < knots[i+1])
return 1.;
else
return 0.;
}
float b2(int i, float t, float *knots)
{
float n, d;
float sum = 0.;
n = t - knots[i];
d = knots[i+1] - knots[i];
if (d != 0.) {
sum += (n / d) * b1(i, t, knots);
}
n = knots[i+2] - t;
d = knots[i+2] - knots[i+1];
if (d != 0.) {
sum += (n / d) * b1(i+1, t, knots);
}
return sum;
}
float b3(int i, float t, float *knots)
{
float n, d;
float sum = 0.;
n = t - knots[i];
d = knots[i+2] - knots[i];
if (d != 0.) {
sum += (n / d) * b2(i, t, knots);
}
n = knots[i+3] - t;
d = knots[i+3] - knots[i+1];
if (d != 0.) {
sum += (n / d) * b2(i+1, t, knots);
}
return sum;
}
float b4(int i, float t, float *knots)
{
float n, d;
float sum = 0.;
n = t - knots[i];
d = knots[i+3] - knots[i];
if (d != 0.) {
sum += (n / d) * b3(i, t, knots);
}
n = knots[i+4] - t;
d = knots[i+4] - knots[i+1];
if (d != 0.) {
sum += (n / d) * b3(i+1, t, knots);
}
return sum;
}
// constructor for NUB curves
Nub_curve::Nub_curve(void)
{
num_knots = 0;
knots = (float *)NULL;
}
// destructor for curve
Nub_curve::~Nub_curve(void)
{
if (num_knots != 0 || knots) {
free(knots);
}
}
// setup the knot vector with user information
void Nub_curve::set_knot_vector(int count, float *v) {
if (count != num_knots && knots) {
free(knots);
}
knots = (float *)malloc(count * sizeof(float));
if (knots) {
int i;
num_knots = count;
for (i = 0; i < num_knots; i++)
knots[i] = v[i];
}
}
// display method for a NUB curve
void Nub_curve::display_curve(int n, int color, int show_geom_pts,
int show_convex_hull)
{
int i, j, m;
float t;
float delta;
float x0, y0, x, y, x1, y1;
float bt0, bt1, bt2, bt3;
#define NUBEVAL(X, Y) \
bt0 = b4(i-3, t, knots); \
bt1 = b4(i-2, t, knots); \
bt2 = b4(i-1, t, knots); \
bt3 = b4(i, t, knots); \
X = geom[i-3].get_x() * bt0 + \
geom[i-2].get_x() * bt1 + \
geom[i-1].get_x() * bt2 + \
geom[i].get_x() * bt3; \
Y = geom[i-3].get_y() * bt0 + \
geom[i-2].get_y() * bt1 + \
geom[i-1].get_y() * bt2 + \
geom[i].get_y() * bt3;
delta = (knots[4] - knots[3]) / (float)n;
i = 3;
t = knots[3];
NUBEVAL(x0, y0);
for (j = 2; j <= n; j++) {
t += delta;
NUBEVAL(x, y);
line(
(short)x0, (short)(floor(y0)+0.5),
(short)x, (short)(floor(y)+0.5),
color);
x0 = x; y0 = y;
}
if (show_geom_pts) {
for (i = 0; i < 4; i++) {
x = geom[i].get_x();
y = floor(geom[i].get_y())+0.5;
line((short)x-CROSS_SIZE, (short)y, (short)x+CROSS_SIZE,
(short)y, color);
line((short)x, (short)y-CROSS_SIZE, (short)x,
(short)y+CROSS_SIZE, color);
}
}
if (show_convex_hull) {
for (i = 0; i < 3; i++) {
x0 = geom[i].get_x();
y0 = floor(geom[i].get_y())+0.5;
x1 = geom[i+1].get_x();
y1 = floor(geom[i+1].get_y())+0.5;
line((short)x0, (short)y0, (short)x1, (short)y1, color);
}
}
}
// display method for a NUB curve
void Nurb_curve::display_curve(int n, int color, int show_geom_pts,
int show_convex_hull)
{
int i, j;
float t;
float delta;
float x0, y0, w, x, y, x1, y1;
float bt0, bt1, bt2, bt3;
#define NURBEVAL(X, Y, W) \
bt0 = b4(i-3, t, knots); \
bt1 = b4(i-2, t, knots); \
bt2 = b4(i-1, t, knots); \
bt3 = b4(i, t, knots); \
X = geom[i-3].get_x() * bt0 + \
geom[i-2].get_x() * bt1 + \
geom[i-1].get_x() * bt2 + \
geom[i].get_x() * bt3; \
Y = geom[i-3].get_y() * bt0 + \
geom[i-2].get_y() * bt1 + \
geom[i-1].get_y() * bt2 + \
geom[i].get_y() * bt3; \
W = geom[i-3].get_w() * bt0 + \
geom[i-2].get_w() * bt1 + \
geom[i-1].get_w() * bt2 + \
geom[i].get_w() * bt3; \
if (W != 1. && W != 0.) { \
X /= W; Y /= W; \
}
delta = (knots[4] - knots[3]) / (float)n;
i = 3;
t = knots[3];
NURBEVAL(x0, y0, w);
for (j = 2; j <= n; j++) {
t += delta;
NURBEVAL(x, y, w);
line(
(short)x0, (short)(floor(y0)+0.5),
(short)x, (short)(floor(y)+0.5),
color);
x0 = x; y0 = y;
}
if (show_geom_pts) {
for (i = 0; i < 4; i++) {
x = geom[i].get_x();
y = floor(geom[i].get_y())+0.5;
line((short)x-CROSS_SIZE, (short)y, (short)x+CROSS_SIZE,
(short)y, color);
line((short)x, (short)y-CROSS_SIZE, (short)x,
(short)y+CROSS_SIZE, color);
}
}
if (show_convex_hull) {
for (i = 0; i < 3; i++) {
x0 = geom[i].get_x();
y0 = floor(geom[i].get_y())+0.5;
x1 = geom[i+1].get_x();
y1 = floor(geom[i+1].get_y())+0.5;
line((short)x0, (short)y0, (short)x1, (short)y1, color);
}
}
}
[LISTING THREE]
// Horner's method for curve display method
void Basis_matrix_curve::display_curve(int n, int color, int
show_geom_pts,
int show_convex_hull)
{
int i;
float delta;
float t, t2, t3;
float x0, y0, x, y, x1, y1;
x0 = t_coeff[3].get_x(); y0 = t_coeff[3].get_y();
delta = 1. / (float)n;
t = 0.;
for (i = 0; i <= n; i++) {
t += delta;
t2 = t * t;
t3 = t2 * t;
x = t_coeff[0].get_x() * t3 +
t_coeff[1].get_x() * t2 +
t_coeff[2].get_x() * t +
t_coeff[3].get_x();
y = t_coeff[0].get_y() * t3 +
t_coeff[1].get_y() * t2 +
t_coeff[2].get_y() * t +
t_coeff[3].get_y();
line(
(short)x0, (short)(floor(y0)+0.5),
(short)x, (short)(floor(y)+0.5),
color);
x0 = x; y0 = y;
}
if (show_geom_pts) {
for (i = 0; i < 4; i++) {
x = (short)geom[i].get_x();
y = (short)(floor(geom[i].get_y())+0.5);
line(x-CROSS_SIZE, y, x+CROSS_SIZE, y, color);
line(x, y-CROSS_SIZE, x, y+CROSS_SIZE, color);
}
}
if (show_convex_hull) {
for (i = 0; i < 3; i++) {
x0 = geom[i].get_x();
y0 = (floor(geom[i].get_y())+0.5);
x1 = geom[i+1].get_x();
y1 = (floor(geom[i+1].get_y())+0.5);
line(x0, y0, x1, y1, color);
}
}
}
[LISTING FOUR]
// Forward differencing for curve display method
void Basis_matrix_curve::display_curve(int n,
int color,
int show_geom_pts,
int show_convex_hull)
{
int i;
float d, delta, delta2, delta3;
float x0, y0, x, y, x1, y1;
float dx, d2x, d3x;
float dy, d2y, d3y;
delta = 1. / (float)n;
delta2 = delta * delta;
delta3 = delta2 * delta;
dx = t_coeff[0].get_x() * delta3 +
t_coeff[1].get_x() * delta2 +
t_coeff[2].get_x() * delta;
d2x = 6. * t_coeff[0].get_x() * delta3 +
2. * t_coeff[1].get_x() * delta2;
d3x = 6. * t_coeff[0].get_x() * delta3;
dy = t_coeff[0].get_y() * delta3 +
t_coeff[1].get_y() * delta2 +
t_coeff[2].get_y() * delta;
d2y = 6. * t_coeff[0].get_y() * delta3 +
2. * t_coeff[1].get_y() * delta2;
d3y = 6. * t_coeff[0].get_y() * delta3;
x = x0 = t_coeff[3].get_x(); y = y0 = t_coeff[3].get_y();
for (i = 1; i <= n; i++) {
x += dx; dx += d2x; d2x += d3x;
y += dy; dy += d2y; d2y += d3y;
line(
(short)x0, (short)(floor(y0)+0.5),
(short)x, (short)(floor(y)+0.5),
color);
x0 = x; y0 = y;
}
if (show_geom_pts) {
for (i = 0; i < 4; i++) {
x = (short)geom[i].get_x();
y = (short)(floor(geom[i].get_y())+0.5);
line((short)x-CROSS_SIZE, (short)y,
(short)x+CROSS_SIZE, (short)y, color);
line((short)x, (short)y-CROSS_SIZE,
(short)x, (short)y+CROSS_SIZE, color);
}
}
if (show_convex_hull) {
for (i = 0; i < 3; i++) {
x0 = geom[i].get_x();
y0 = floor(geom[i].get_y())+0.5;
x1 = geom[i+1].get_x();
y1 = floor(geom[i+1].get_y())+0.5;
line((short)x0, (short)y0,
(short)x1, (short)y1, color);
}
}
}
[LISTING FIVE]
// demo.cc - demonstrating different curves
// For Unix and X, compile with:
#include <iostream.h>
#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#include "curve.h"
extern "C" {
void init_graphics_device(void);
void close_graphics_device(void);
void clear_window(void);
void text_output(char *);
void text_output_and_wait(char *);
};
main()
{
float knot[8];
Point2d geom[4];
init_graphics_device();
Hermite_curve herm;
text_output("Hermite Curves");
geom[0].set_xy(20., 20.);
geom[1].set_xy(300., 20.);
geom[2].set_xy(0., 500.);
geom[3].set_xy(0., -500.);
herm.set_geom_vector(4, geom);
herm.display_curve(32, 14, 0, 0);
geom[0].set_xy(20., 20.);
geom[1].set_xy(300., 20.);
geom[2].set_xy(0., 200.);
geom[3].set_xy(0., -200.);
herm.set_geom_vector(4, geom);
herm.display_curve(32, 13, 0, 0);
geom[0].set_xy(20., 20.);
geom[1].set_xy(300., 20.);
geom[2].set_xy(0., 200.);
geom[3].set_xy(-100., 0.);
herm.set_geom_vector(4, geom);
herm.display_curve(32, 7, 0, 0);
geom[0].set_xy(20., 20.);
geom[1].set_xy(300., 20.);
geom[2].set_xy(200., 200.);
geom[3].set_xy(-200., -200.);
herm.set_geom_vector(4, geom);
herm.display_curve(32, 6, 0, 0);
text_output_and_wait("Press return:");
Bezier_curve bez;
Nub_curve nub;
clear_window();
text_output("Bezier and NUB Curves");
geom[0].set_xy(20., 20.);
geom[1].set_xy(50., 180.);
geom[2].set_xy(300., 50.);
geom[3].set_xy(100., 10.);
bez.set_geom_vector(4, geom);
bez.display_curve(32, 15, 1, 1);
text_output_and_wait("Press return:");
knot[0] = knot[1] = knot[2] = knot[3] = 0.;
knot[4] = knot[5] = knot[6] = knot[7] = 1.;
nub.set_knot_vector(8, knot);
nub.set_geom_vector(4, geom);
nub.display_curve(32, 6, 0, 0);
knot[0] = knot[1] = knot[2] = 0.;
knot[3] = 1.; knot[4] = 2.;
knot[5] = knot[6] = knot[7] = 3.;
nub.set_knot_vector(8, knot);
nub.set_geom_vector(4, geom);
nub.display_curve(32, 5, 0, 0);
text_output_and_wait("Press return:");
Catmull_Rom_curve cm;
clear_window();
text_output("Catmull-Rom Curve");
geom[0].set_xy(20., 20.);
geom[1].set_xy(50., 180.);
geom[2].set_xy(300., 50.);
geom[3].set_xy(100., 10.);
cm.set_geom_vector(4, geom);
cm.display_curve(32, 15, 1, 1);
text_output_and_wait("Press return:");
Beta_spline_curve beta;
clear_window();
text_output("Beta-splines");
geom[0].set_xy(20., 20.);
geom[1].set_xy(170., 180.);
geom[2].set_xy(150., 180.);
geom[3].set_xy(300., 20.);
beta.set_geom_vector(4, geom);
beta.display_curve(32, 7, 1, 1);
text_output_and_wait("Press return:");
Nurb_curve nurb;
clear_window();
text_output("NURB Curves with varying Knots");
geom[0].set_xy(20., 20.);
geom[1].set_xy(50., 180.);
geom[2].set_xy(300., 50.);
geom[3].set_xy(100., 10.);
knot[0] = knot[1] = knot[2] = knot[3] = 0.;
knot[4] = knot[5] = knot[6] = knot[7] = 1.;
nurb.set_knot_vector(8, knot);
nurb.set_geom_vector(4, geom);
nurb.display_curve(32, 15, 1, 1);
knot[0] = knot[1] = knot[2] = 0.; knot[3] = 1.;
knot[4] = 2.; knot[5] = knot[6] = knot[7] = 3.;
nurb.set_knot_vector(8, knot);
nurb.set_geom_vector(4, geom);
nurb.display_curve(32, 14, 0, 0);
knot[0] = 0.; knot[1] = 1.; knot[2] = 2.; knot[3] = 3.;
knot[4] = 4.; knot[5] = 5.; knot[6] = 6.; knot[7] = 7.;
nurb.set_knot_vector(8, knot);
nurb.set_geom_vector(4, geom);
nurb.display_curve(32, 13, 0, 0);
knot[0] = 0.; knot[1] = 0.; knot[2] = 1.; knot[3] = 1.;
knot[4] = 2.; knot[5] = 2.; knot[6] = 3.; knot[7] = 3.;
nurb.set_knot_vector(8, knot);
nurb.set_geom_vector(4, geom);
nurb.display_curve(32, 12, 0, 0);
text_output_and_wait("Press return:");
close_graphics_device();
return 0;
}
[LISTING SIX]
// utilTC.cc - display utilities for Turbo C
// For Borland C++, compile with:
// bcc -c -mh -P -v -Ie:\bc\include util.cc
#include <iostream.h>
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <conio.h>
#include <alloc.h>
#include <dos.h>
#include <math.h>
#include "curve.h"
#define WIDTH 320
#define HEIGHT 200
#define ABS(x) ((x) < 0) ? -(x) : (x)
void init_graphics_device(void);
void close_graphics_device(void);
void set_pixel(short, short, short, short, short, short, short);
void line(short, short, short, short, short);
void clear_window(void);
void text_output(char *);
void text_output_and_wait(char *);
void display_string(int, char *);
// Global variables associated with the VGA display
unsigned char far *vga_fb = (unsigned char far *)MK_FP(0xA000,0);
void init_graphics_device(void)
{
union REGS inregs, outregs;
inregs.h.ah = 0; // set video mode command
inregs.h.al = 0x13; // 320x200x8 mode
int86(0x10, &inregs, &outregs);
inregs.h.ah = 0xf; // check video mode
int86(0x10, &inregs, &outregs);
if (inregs.h.al != 0x13) {
printf("cannot get 320x200x8 VGA mode : got %x\n",
inregs.h.al);
exit(1);
}
}
void close_graphics_device(void)
{
union REGS inregs, outregs;
inregs.h.ah = 0; // set video mode command
inregs.h.al = 0x2; // 80x25 mode
int86(0x10, &inregs, &outregs);
}
void line(short x1, short y1, short x2, short y2, short color)
{
short x, y;
short deltax, deltay;
short temp;
short err;
short i;
short swap;
short s1, s2;
x = x1; y = y1;
deltax = (short)ABS(x2 - x1);
deltay = (short)ABS(y2 - y1);
if ((x2 - x1) < 0.) s1 = -1; else s1 = 1;
if ((y2 - y1) < 0.) s2 = -1; else s2 = 1;
if (deltay > deltax) {
temp = deltax;
deltax = deltay;
deltay = temp;
swap = 1;
}
else
swap = 0;
err = 2 * deltay - deltax;
for (i = 1; i <= deltax; i++) {
set_pixel(0, 0, WIDTH-1, HEIGHT-1, x, y, color);
while (err >= 0) {
if (swap) x += s1; else y += s2;
err -= 2 * deltax;
}
if (swap) y += s2; else x += s1;
err += 2 * deltay;
}
}
void set_pixel(short xmin, short ymin, short xmax, short ymax,
short x, short y, short color)
{
unsigned char far *pixel;
unsigned long offset;
if (x >= xmin && x <= xmax && y >= ymin && y <= ymax) {
offset = (unsigned long)(((unsigned long)WIDTH * (unsigned long)y) +
(unsigned long)x);
pixel = (unsigned char far *)(unsigned long)vga_fb +
(unsigned long)offset;
*pixel = (unsigned char)color;
}
}
void clear_window()
{
unsigned long far *pixel32;
unsigned long count;
// Clear 4 8-bit pixels at once by writing a 32 bit value
// with a 32-bit plane mask.
pixel32 = (unsigned long far *)vga_fb;
// compute number of 8-bit pixels and then divide by 4 for words
count = (unsigned long)WIDTH * (unsigned long)HEIGHT;
count = count >> 2;
while (count > 0) {
*pixel32 = (unsigned long)0;
pixel32++;
count--;
}
}
void text_output(char *str)
{
display_string(23, str);
}
void text_output_and_wait(char *str)
{
display_string(24, str);
getch();
}
void display_string(int row, char *str)
{
int i;
union REGS inregs, outregs;
for (i = 0; i < strlen(str); i++) {
inregs.h.ah = 0x2; // set cursor position
inregs.h.bh = 0x0; // display page
inregs.h.dh = row; // near bottom of display
inregs.h.dl = i; // column
int86(0x10, &inregs, &outregs);
inregs.h.ah = 0x9; // write character
inregs.h.bh = 0x0; // display page
inregs.h.bl = 0x7; // black BG/white FG
inregs.h.al = str[i]; // character
inregs.x.cx = 1; // repeat count
int86(0x10, &inregs, &outregs);
}
}
[LISTING SEVEN]
/* utilXlib.cc - utilities for Xlib */
#include <X11/Xlib.h>
#include <X11/X.h>
#include <math.h>
#define WIDTH 320
#define HEIGHT 220
void init_graphics_device(void);
void close_graphics_device(void);
void line(short, short, short, short, short);
void clear_window(void);
void text_output(char *);
void text_output_and_wait(char *);
Display *dpy;
Window win;
GC gc;
void init_graphics_device(void)
{
XEvent event;
int screen;
int done;
/* Open the X display and get the screen */
dpy = XOpenDisplay(NULL);
if (!dpy) {
printf("could not open display\n");
exit(1);
}
screen = DefaultScreen(dpy);
/* Create a simple window which will be a child to the root window of
* the display. */
win = XCreateSimpleWindow(dpy, RootWindow(dpy, screen),
0, 0, WIDTH, HEIGHT, 2, BlackPixel(dpy, screen),
WhitePixel(dpy, screen));
/* Display the window */
XMapWindow(dpy, win);
XSelectInput(dpy, win, ExposureMask | KeyPressMask |
ButtonPressMask);
done = 0;
while(!done) {
XNextEvent(dpy, &event);
switch(event.type) {
case Expose:
done = 1;
break;
}
}
/* Create an X graphics context for rendering vectors that compromise
* curves. */
gc = XCreateGC(dpy, win, 0, NULL);
XSetForeground(dpy, gc, BlackPixel(dpy, screen));
}
void close_graphics_device(void)
{
exit(0);
}
void line(short x1, short y1, short x2, short y2, short color)
{
/* Draw a 2D line using X calls and then flush the pixels to the server. */
XDrawLine(dpy, win, gc, x1, y1, x2, y2);
XFlush(dpy);
}
void clear_window()
{
/* Clear the X window to the background color. */
XClearWindow(dpy, win);
XFlush(dpy);
}
void text_output(char *str)
{
/* Put a string into the window. */
XDrawString(dpy, win, gc, 0, 192, str, strlen(str));
XFlush(dpy);
}
void text_output_and_wait(char *str)
{
/* Put a string into the window and wait for user input. */
XDrawString(dpy, win, gc, 0, 210, str, strlen(str));
XFlush(dpy);
getchar();
}
Copyright © 1992, Dr. Dobb's Journal