### Representation of Various Entities

Index is used for array subscripts.

- Index starts at zero (0).
- An array index of ‘4’ means that there are ‘5’ items in the array.
- This is standard ‘C’ practice, but is also the cause of many errors by users.

Point: is always a 3-D entity (x,y,z)

Control Points are always weighted. The following cases are distinguished

- 2-D non-rational (x, y, NL_NOZ, NL_NOW)
- 2-D rational (wx, wy, NL_NOZ, w)
- 3-D non-rational (x, y, z, NL_NOW)
- 3-D rational (wx, wy, wz, w)

where NL_NOZ and NL_NOW are predefined identifiers telling NLib that no z and/or no w component(s) are to be considered.

Knot Vector is always “clamped”, that is:

- U = {a(0),…,a(p), u(p+1),…,u(m-p-1),b(m-p),…,b(m)}

where p is the degree of the curve. There are routines to clamp or unclamp knot vectors if the need so arises.

Curve is always in a Bezier-like form, i.e. with end point interpolation. Curves are mostly created by the system. However, they can be input from a file. The form of the file is explained by the following degree-three curve example:

5 <— Highest index of control points (n)

3 <— Degree (p)

1 <— Rational (1) or non-rational (0) flag

2 <— Dimension (2 or 3)

0.100000 0.300000 1.000000 <— Pw(0)

0.700000 0.600000 2.000000 <— Pw(1)

0.350000 0.600000 1.000000 <— Pw(2)

0.650000 0.600000 1.000000 <— Pw(3)

0.975000 0.450000 1.500000 <— Pw(4)

0.900000 0.300000 1.000000 <— Pw(5)

0.000000 <— u(0)

0.000000 <— u(1)

0.000000 <— u(2)

0.000000 <— u(3)

0.300000 <— u(4)

0.700000 <— u(5)

1.000000 <— u(6)

1.000000 <— u(7)

1.000000 <— u(8)

1.000000 <— u(9)

The relationship m=n+p+1 must always hold, where m is the highest index of knots (in the above example have 9 = 5 + 3 + 1).

Surface is always in a Bezier-like form, i.e. with corner point interpolation. As in the case of curves, surfaces are mostly created by the system. However, they too can be input from a file. The form of the file is explained by the following degree 3 x 2 surface example:

4 3 <— Highest indexes in u- and v-directions (n, m)

3 2 <— Degrees in u- and v-directions (p, q)

1 <— Rational (1) or non-rational (0) flag

8.000000 -8.000000 4.000000 1.000000 <— Pw(0,0)

4.000000 -8.000000 0.000000 1.000000 <— Pw(0,1)

-4.000000 -8.000000 0.000000 1.000000 <— Pw(0,2)

-8.000000 -8.000000 6.000000 1.000000 <— Pw(0,3)

8.000000 -3.000000 0.000000 1.000000 <— Pw(1,0)

8.000000 -6.000000 0.000000 2.000000 <— Pw(1,1)

-4.000000 -3.000000 0.000000 1.000000 <— Pw(1,2)

-8.000000 -3.000000 2.000000 1.000000 <— Pw(1,3)

8.000000 0.000000 0.000000 1.000000 <— Pw(2,0)

6.000000 0.000000 3.000000 1.500000 <— Pw(2,1)

-4.000000 0.000000 2.000000 1.000000 <— Pw(2,2)

-8.000000 0.000000 4.000000 1.000000 <— Pw(2,3)

8.000000 3.000000 2.000000 1.000000 <— Pw(3,0)

4.000000 3.000000 0.000000 1.000000 <— Pw(3,1)

-8.000000 6.000000 0.000000 2.000000 <— Pw(3,2)

-8.000000 3.000000 0.000000 1.000000 <— Pw(3,3)

8.000000 8.000000 -2.000000 1.000000 <— Pw(4,0)

4.000000 8.000000 0.000000 1.000000 <— Pw(4,1)

-4.000000 8.000000 -2.000000 1.000000 <— Pw(4,2)

-8.000000 8.000000 6.000000 1.000000 <— Pw(4,3)

0.000000 <— u(0)

0.000000 <— u(1)

0.000000 <— u(2)

0.000000 <— u(3)

0.500000 <— u(4)

1.000000 <— u(5)

1.000000 <— u(6)

1.000000 <— u(7)

1.000000 <— u(8)

0.000000 <— v(0)

0.000000 <— v(1)

0.000000 <— v(2)

0.500000 <— v(3)

1.000000 <— v(4)

1.000000 <— v(5)

1.000000 <— v(6)

The relationships r=n+p+1 and s=m+q+1 must always hold, where r and s are the highest indexes of knots in the u- and v-directions, respectively (in the above example we have 8=4+3+1 and 6=3+2+1).

### Program Design

We explain the major ingredients of a typical NLib application program using a simple example below. The program reads in a NURBS curve, elevates its degree by one in three different ways, and writes the last result out to a file. The names of the input and output files are command line arguments.

#include “nurbs.h”

#include “NL_Globals.h”

main( int argc, char *argv[] )

{

NL_FLAG error;

NL_INDEX n, m, s;

NL_DEGREE p;

NL_CURVE curP, curQ, curR;

NL_KNOTVECTOR *knt

NL_STACKS S;

/* Start NURBS */

N_InitNurbs(&S);

/* Input curve */

N_CrvInitArrays(&curP);

error = N_CrvReadFromFile(&curP,argv[1],NL_YES,&S);

if ( error EQ NL_YES ) NL_PRINT_ERROR;

/***************************/

/* Elevate degree call 1 */

/***************************/

N_CrvInitArrays(&curQ);

error = N_CrvElevateDegree(&curP,1,&curQ,&S,&S);

if ( error EQ NL_YES ) NL_PRINT_ERROR;

/***************************/

/* Elevate degree call 2 */

/***************************/

/* Get knot vector, highest indexes and degree */

N_CrvGetKnotVector(&curP,&knt);

N_CrvGetArraySizes(&curP,&n,&m);

N_CrvGetDegree(&curP,&p);

/* Compute highest indexes of degree elevated curve */

N_BasisGetSpanCount(knt,p,&s);

n = n+s;

p = p+1;

m = m+s+1;

/* Allocate memory for new curve */

error = N_CrvAllocArrays(&curR,n,p,m,&S);

if ( error EQ NL_YES ) NL_PRINT_ERROR;

/* Now degree elevate */

error = N_CrvElevateDegree(&curP,1,&curR,&S,&S);

if ( error EQ NL_YES ) NL_PRINT_ERROR;

/***************************/

/* Elevate degree call 3 */

/***************************/

error = N_CrvElevateDegree(&curP,1,&curP,&S,&S);

if ( error EQ NL_YES ) NL_PRINT_ERROR;

/* Output last curve */

error = N_CrvWriteToFile(&curP,argv[2]);

if ( error EQ NL_YES ) NL_PRINT_ERROR;

/* End NURBS */

N_EndNurbs(&S);

}

There are several important points:

- Each routine that calls NLib functions must include the file nurbs.h.

This header file contains all function prototypes, and includes the headers NL_Nurbsdef.h, NL_Geom.h and NL_DataStruct.h. - The top level application routine must also include NL_Globals.h.

This file contains global parameters and constants, including important tolerances (see below). - NLib functions work on predefined data types, e.g. NL_FLAG, NL_INDEX, NL_DEGREE, NL_KNOTVECTOR, NL_CURVE and NL_STACKS. These data types act like any other standard data type such as int or double. That is, when one writes NL_CURVE curP, the system allocates memory for a curve object and tags this memory with the name curP.
- NLib routines start and end with two special functions N_InitNurbs and N_EndNurbs. The first one initializes memory stacks to NULL, and the second one traverses these stacks, deallocates all dynamically allocated memory, and kills the stacks.
- The calling mechanism of each NLib function is explained in detail in the documentation and its header.

As an example, here is the explanation of how to call N_CrvElevateDegree (as extracted from its header):

A typical calling example is:

NL_CURVE curP, curQ;

NL_INDEX t;

NL_STACKS SP, SQ;

…

(define curP, get t);

…

N_CrvInitArrays(&curQ);

N_CrvElevateDegree(&curP,t,&curQ,&SP,&SQ);

N_CrvElevateDegree(&curP,t,&curP,&SP,&SP);

…

The header explains

- how to declare variables;
- how to pass them to the routine (by value or by reference).
- how to manage memory (which stack is used for input and output

We shall now walk through the sample program

Various parameters and objects are defined. Note that three curve objects are defined, however, only a knot vector pointer is needed.

The NURBS program is initialized by N_InitNurbs

curP is initialized to NULL telling the N_CrvReadFromFile routine to allocate memory, copy the curve data into curP, and place memory pointers on to the stacks.

Call 1 curQ is initialized to NULL telling the degree elevation routine to:

- find out how much memory is needed to store the output curve, and to allocate this memory;
- degree elevate curP and to place the result in curQ; and
- place memory pointers on to the stacks S.

This is the most convenient way of calling NLib routines that return a new curve or surface.

Call 2 the new degree and the highest indexes of control points and knots are computed first. Then memory is allocated by the N_CrvAllocArrays routine. When calling the degree elevation routine with curR, the routine does the following:

- Sees that curR is not NULL and checks if sufficient memory is passed in using curR’s memory stacks (the same as curP’s in this example). If yes, it continues. If not, it sets t error flag, and returns control to the calling routine.
- Degree elevates curP and places the output curve in curR,

Call 3 when calling the degree elevation routine with the same curve pointer, the routine sees that degree elevation is to be done in-place, and does the following:

- kills curP.
- allocates new memory for curP; and
- places the degree elevated curve in the new curP.

curP is saved on a file by the N_CrvWriteToFile routine.

The program is ended by N_EndNurbs that deallocates all dynamically allocated memory.

The three calls show that there is a lot of flexibility in how to create new geometric entities. For simple operations calls 1 or 3 are recommended. If degree elevation is to be done repeatedly, call 2 is recommended, i.e. allocate enough memory, only once, for all types of degree elevations, and pass this memory to the routine.

### Array Arguments

In general, it is safe to say that we never have arrays of points ( ie point data laid out contiguously). What we have is an array of pointers to points…..so the points can be, for example, along either the row or column of an array. We gather the address of each successive point we are interested in and add the address into an array of pointers……

NL_POINT **Pnt

We do not have to rearrange the input if we switch from row to column. …just select the addresses differently.

Similarly with other objects, such as the curves in N_CrvDecomposeBez.

The argument ” ***CurQ ” means ” the address of an array of pointers to Curves”.

You have:

NL_CURVE *cur as a pointer to a curve

NL_CURVE **curQ as an array of curve pointers

( do not confuse this with a two dimensional array )

You call N_CrvDecomposeBez with (…….&curQ…..)

But N_CrvDecomposeBez receives it as (………NL_CURVE ***curQ……..) and internally builds a set list of curve pointers to newly allocated curves

NL_CURVE **curA; ….

build curA[i]; ……

*curQ= curA;

” take the address of the new curve list curA and put it into the location where curQ is pointing to (if there was no *curQ…but only curQ it would just locally change the address of curQ and not update the one that was sent).

After returning, the calling routine now has the address of the list of curve pointers ( NL_CURVE **CurQ).

## Dynamic Memory Allocation

Memory in NLib routines is allocated by special purpose routines that can be found in the space subdirectory and have the prefix S_. These routines require a memory stack to save pointers of allocated memory. The user can work with an entire hierarchy of stacks, each one of them is local to one or global to many routines. These stacks can control when and where memory is needed. Each NURBS routine has its local stack which is killed, along with all memory residing on it, upon leaving the routine. In addition to this local stack, the routine can take, as arguments, a number of global stacks that can hold memory needed when the routine in no longer active. As an example, let us look at the structure of the degree elevation routine used above.

NL_FLAG N_CrvElevateDegree( NL_CURVE *curP, NL_INDEX t, NL_CURVE *curQ, NL_STACKS *SP, NL_STACKS *SQ )

{

NL_FLAG error = NL_NO;

…

NL_STACKS SL;

N_InitNurbs(&SL);

…

alfs = N_AllocReal1dArray(p,&SL);

if ( alfs EQ NULL ) NL_QUIT;

…

N_EndNurbs(&SL);

}

N_CrvElevateDegree is given two global stacks, SP and SQ; one for curP and one for curQ. In addition to this, it has its local stack SL. If the degree is elevated in-place, the new curve curP is placed on to SP. Otherwise, the new curve curQ is placed on to SQ. Any local memory, like alfs, is stored on the local stack SL, and is released upon leaving the routine.

### Tolerances and Constants

NLib has two types of tolerances public and private. Public tolerances are declared and initialized in globals.h. Some examples follow:

- NL_MTOL models space tolerance, i.e. if two points are closer than NL_MTOL, they are considered the same. The current setting is1.0e-7.
- NL_PTOL parameter space tolerance, i.e. if two parameters are closer than NL_PTOL, they are considered the same. The current setting is 1.0e-12.
- NL_LTOL line tolerance, i.e. if the magnitude of the cross product of two lines’ unit direction vectors is less then NL_LTOL, the lines are considered parallel. The current setting is 1.0e-12.
- NL_LUDT LU-decomposition tolerance. It will not allow the computation of LU-decomposed matrices (with or without pivoting) with huge numbers even if the divisions can be performed within the available range of the floating point processor. The current setting is 1.0e-9.

Private tolerances are declared in individual routines. These tolerances are used to perform numerical computations within the scope of the individual routine. Even though they are NL_PRIVATE, they are never hard coded. They are defined on top of the routine, and adjusted in the body based on the length of the knot vector, the bounding box of the curve/surface, etc. Their initial settings are based on extensive experiments, and can be altered easily as the user sees fit.

There are a number of constants and predefined entities NLib relies on. All of them are defined in globals.h. A few examples follow:

- NL_WMIN and NL_WMAX minimum and maximum allowable weights. Current settings 1.0e-2 and 10, respectively.
- NL_DMAX highest degree. Currently 28. Warning NL_DMAX cannot be modified without considering the definition of NL_INDEX and NL_INTEGER data types. The binomial coefficients are represented for efficiency reason as integers. The current definition of NL_INDEX and NL_INTEGER is long, which may have different maximum values on different systems. The ANSI standard requires long to be at least 2,147,483,647.
- NL_ITLIM iteration limit. Currently 50 (if iterations, such as Newton, do not converge within 50 steps, they probably never will).

### Floating Point Operations

Floating point operations are performed in double precision. There are two types of error checks:

- Certain numbers are checked against predefined tolerances before any crucial operation, such as division, occurs. For example, if a number is to be divided by the distance between two points, then even if the division could be done, it is not done if the distance is less than the model space tolerance.
- When no specific tolerances are involved, operation checks are performed by the N_FloatOpIsBad routine. It’s calling mechanism is:

NL_REAL x, y;…

if ( N_FloatOpIsBad(x,y,NL_ADDITION) ) NL_ERROR(NL_NUM_ERR);

if ( N_FloatOpIsBad(x,y,NL_SUBTRACTION) ) NL_ERROR(NL_NUM_ERR);

if ( N_FloatOpIsBad(x,y,MULTIPLICATION) ) NL_ERROR(NL_NUM_ERR);

if ( N_FloatOpIsBad(x,y,NL_DIVISION) ) NL_ERROR(NL_NUM_ERR);

…

### Utilities

NLib has an extensive set of utilities to perform such tasks as making objects, input/output, converting objects, etc. One of the most important tasks performed by these functions is to hide the definition of various objects while providing access to their components. For example, to get the highest indexes of control points and knots of a curve object, one could write:

NL_INDEX n, m;

NL_CURVE *cur;

…

n = cur->pol->n;

m = cur->knt->m;

This requires knowledge of the NL_CURVE object. It also requires that the entire library be searched and modified if any change is made to its definition. In order to avoid this, NLib has an extensive set of utility routines to create objects from their individual components, and to break them into the many constituents. In the above example, one should write:

INDEX n, m;

NL_CURVE *cur;

…

N_CrvGetArraySizes(cur,&n,&m);

Each category of routines where objects are dealt with, e.g. NURBS, geometry and mathematics, has its own constructors and destructors. We strongly recommend that the application programmer use these routines. They have the following advantages:

- they are easier to understand;
- they are less prone to errors (the compiler catches the error if the are called incorrectly); and
- they allow changes to be made to the definitions of the objects.

## Error Checking

NLib has a set of routines to check errors. They are in the error directory and have the prefix E_. We strongly recommend that these routines be used to make sure that curves and surfaces are correct before entering NURBS routines. The general philosophy is: Whenever an object is passed into a function, it is assumed that its definition is correct.

### Usage Tips

Please refer to Examples/NLib for additional examples. Maybe provide FTP link.

The curve/surface evaluators use some fixed- length arrays declared with dimensions NL_MAXDEG and NL_MAXDER. The current settings (in NL_Nurbsdef.h) for these evaluators are NL_MAXDEG=28 (maximum degree) and NL_MAXDER=5 (maximum derivative order). If you need to change these, you should do so before compiling NLib (you must also change NL_DMAX in NL_Globals.h to be the same as NL_MAXDEG).

“The NURBS Book” (Springer-Verlag, ISBN 3-540-61545-8), by Les Piegl and Wayne Tiller is excellent documentation for NLib. The first 12 chapters describe the underlying mathematics and algorithms, and Chapter 13 describes the system architecture (data structures, memory management, error handling, utilities, example programs, etc). We suggest you buy it.

A few tips for using NLib (see “The NURBS Book” for more detail)

- NLib requires NURBS curves/surfaces to have “clamped”, or “nonperiodic” knot vectors; i.e. end multiplicities equal to degree plus 1. Furthermore, internal knots may not have multiplicity greater than degree.
- There are functions to “clamp” curves/surfaces (specifically: N_Iges126Crv and N_Iges128Srf), and knot removal can be used where internal knots have multiplicities greater than degree.
- When importing curves/surfaces from other systems, it is advisable to use the functions N_Iges126Crv and N_Iges128Srf. These check if the NURBs definitions are compatible with NLib.
- A main program must include the files nurbs.h and NL_Globals.h (see example programs).
- Any program should start with a call to N_InitNurbs and end with a call to N_EndNurbs (see example programs).
- Generally, if an array is passed into a function, and if an integer> (say n) is also passed, then n will usually indicate the highest valid index of the array, not the number of elements in the array. Since array indices start at 0, n+1 is the actual number of array elements. This may seem unnatural for some users, but it corresponds to the notation and algorithm descriptions of “The NURBS Book“
- All angles are to be expressed in degrees, not radians.

### References

The NLib products are based on solid research in NURBS technology. The publications below by Les Piegl and/or Wayne Tiller are examples of key papers providing a scientific basis for major NLib routines.

- Rational B-splines for curve and surface representation, IEEE Computer Graphics and Applications, Vol. 3, No. 6, 1983, pp 61-69.
- The NURBS Book, Second Revised Edition, Springer-Verlag, 1997.
- Curve and surface constructions using rational B-splines, Computer-Aided Design, Vol. 19, No. 9, 1987, 485-498.
- A menagerie of rational B-spline circles, IEEE Computer Graphics and Applications, Vol. 9, September, 1989, pp 48-56.
- Modifying the shape of rational B-splines. Part 1: curves, Computer-Aided Design, Vol. 21, No. 8, 1989, pp 509-518.
- Modifying the shape of rational B-splines. Part 2: surfaces, Computer-Aided Design, Vol. 21, No. 9, 1989, pp 538-546.
- Data reduction using cubic rational B-spline, IEEE Computer Graphics and Applications, May, pp 60-68, 1992.
- Knot-removal algorithms for NURBS curves and surfaces, Computer-Aided Design, Vol. 24, No. 8, 1992, pp 445-453.
- Delaunay triangulation using a uniform grid, IEEE Computer Graphics and Applications, May, 1993, pp 36-47.
- Software engineering approach to degree elevation of B-spline curves, Computer-Aided Design, Vol. 26, No. 1, pp 17-28 , 1994.
- Algorithm for degree reduction of B-spline curves, Computer-Aided Design, Vol. 27, No. 2, pp 101-110, 1995.
- Algorithm for approximate NURBS skinning, Computer-Aided Design, Vol. 28, No. 9, pp 699-706, 1996.
- Symbolic operators for NURBS, Computer-Aided Design, Vol. 29, No. 5, pp 361-368, 1997
- Algorithm for computing the product of two B-splines, in A. Le Mehaute, C. Rabut and L. L. Schumaker (eds.) Curves and surfaces with applications in CAGD, Vanderbilt University Press, Nashville, TN, 1997, pp 337-344.
- Geometry-based triangulation of trimmed NURBS surfaces, Computer-Aided Design, Vol. 30, No 1, pp 11-18, 1998.
- Computing the derivatives of NURBS with respect to a knot, Computer Aided Geometric Design, Vol. 15, No. 9, pp 925-934, 1998.
- Approximation of offsets of NURBS curves and surfaces, Computer-Aided Design, Vol. 31, No. 2, pp 147-156, 1999.
- Cross-sectional design with boundary constraints, Engineering with Computers, Vol. 15, 1999, pp 171-180.
- Filling n-sided regions with NURBS patches, The Visual Computer, Vol. 15, No. 2, pp 77-89, 1999.
- Reducing control points in surface interpolation, IEEE Computer Graphics and Applications, Vol. 20, No. 5, 2000, pp 70-74.
- Surface approximation to scanned data, The Visual Computer, Vol. 16, No. 7, 2000, pp 386-395.
- Curve interpolation with arbitrary end derivatives, Engineering with computers, Vol. 16, 2000, pp 73-79.
- Least-squares NURBS curve approximation with arbitrary end derivatives, Engineering with computers, Vol. 16, 2000, pp 109-116.
- Parametrization for surface fitting in reverse engineering, Computer-Aided Design, Vol. 33, No. 8, 2001, pp 593-603.