Main Page   Class Hierarchy   Compound List   File List   Compound Members  

PLib::NurbsCurve Class Template Reference

A NURBS curve class. More...

#include <nurbs.hh>

Inheritance diagram for PLib::NurbsCurve::

PLib::ParaCurve PLib::NurbsCurveSP List of all members.

Public Methods

 NurbsCurve ()
 NurbsCurve (const NurbsCurve< T, N > &nurb)
 NurbsCurve (const Vector< HPoint_nD< T, N > > &P1, const Vector< T > &U1, int deg=3)
 NurbsCurve (const Vector< Point_nD< T, N > > &P1, const Vector< T > &W, const Vector< T > &U1, int deg=3)
virtual ~NurbsCurve ()
int degree () const
const Vector< HPoint_nD< T,
N > > & 
ctrlPnts () const
const HPoint_nD< T, N > ctrlPnts (int i) const
const Vector< T > & knot () const
knot (int i) const
void resize (int n, int Deg)
virtual void reset (const Vector< HPoint_nD< T, N > > &P1, const Vector< T > &U1, int deg)
virtual NurbsCurve & operator= (const NurbsCurve< T, N > &)
virtual HPoint_nD< T, N > operator() (T u) const
HPoint_nD< T, N > hpointAt (T u) const
HPoint_nD< T, N > hpointAt (T u, int span) const
void deriveAtH (T u, int, Vector< HPoint_nD< T, N > > &) const
void deriveAt (T u, int, Vector< Point_nD< T, N > > &) const
void deriveAtH (T u, int, int, Vector< HPoint_nD< T, N > > &) const
void deriveAt (T u, int, int, Vector< Point_nD< T, N > > &) const
Point_nD< T, N > derive3D (T u, int d) const
HPoint_nD< T, N > derive (T u, int d) const
Point_nD< T, N > normal (T u, const Point_nD< T, N > &v) const
HPoint_nD< T, N > firstD (T u) const
HPoint_nD< T, N > firstD (T u, int span) const
Point_nD< T, N > firstDn (T u) const
Point_nD< T, N > firstDn (T u, int span) const
basisFun (T u, int i, int p=-1) const
void basisFuns (T u, int span, Vector< T > &N) const
void dersBasisFuns (int n, T u, int span, Matrix< T > &N) const
minKnot () const
maxKnot () const
int findSpan (T u) const
void findMultSpan (T u, int &r, int &s) const
int findMult (int r) const
int findKnot (T u) const
getRemovalBnd (int r, int s) const
void removeKnot (int r, int s, int num)
void removeKnotsBound (const Vector< T > &ub, Vector< T > &ek, T E)
int knotInsertion (T u, int r, NurbsCurve< T, N > &nc)
void refineKnotVector (const Vector< T > &X)
void refineKnotVectorClosed (const Vector< T > &X)
void mergeKnotVector (const Vector< T > &Um)
void clamp ()
void unclamp ()
int leastSquares (const Vector< Point_nD< T, N > > &Q, int degC, int n)
int leastSquares (const Vector< Point_nD< T, N > > &Q, int degC, int n, const Vector< T > &ub)
int leastSquaresH (const Vector< HPoint_nD< T, N > > &Q, int degC, int n, const Vector< T > &ub)
int leastSquares (const Vector< Point_nD< T, N > > &Q, int degC, int n, const Vector< T > &ub, const Vector< T > &knot)
int leastSquaresH (const Vector< HPoint_nD< T, N > > &Q, int degC, int n, const Vector< T > &ub, const Vector< T > &knot)
int leastSquaresClosed (const Vector< Point_nD< T, N > > &Q, int degC, int n)
int leastSquaresClosed (const Vector< Point_nD< T, N > > &Q, int degC, int n, const Vector< T > &ub)
int leastSquaresClosedH (const Vector< HPoint_nD< T, N > > &Q, int degC, int n, const Vector< T > &ub)
int leastSquaresClosed (const Vector< Point_nD< T, N > > &Q, int degC, int n, const Vector< T > &ub, const Vector< T > &knot)
int leastSquaresClosedH (const Vector< HPoint_nD< T, N > > &Q, int degC, int n, const Vector< T > &ub, const Vector< T > &knot)
void globalApproxErrBnd (Vector< Point_nD< T, N > > &Q, int deg, T E)
void globalApproxErrBnd (Vector< Point_nD< T, N > > &Q, Vector< T > &ub, int deg, T E)
void globalApproxErrBnd2 (Vector< Point_nD< T, N > > &Q, int degC, T E)
void globalApproxErrBnd3 (Vector< Point_nD< T, N > > &Q, int degC, T E)
void globalApproxErrBnd3 (Vector< Point_nD< T, N > > &Q, const Vector< T > &ub, int degC, T E)
void globalInterp (const Vector< Point_nD< T, N > > &Q, int d)
void globalInterp (const Vector< Point_nD< T, N > > &Q, const Vector< T > &ub, int d)
void globalInterpH (const Vector< HPoint_nD< T, N > > &Q, int d)
void globalInterpH (const Vector< HPoint_nD< T, N > > &Q, const Vector< T > &U, int d)
void globalInterpH (const Vector< HPoint_nD< T, N > > &Q, const Vector< T > &ub, const Vector< T > &U, int d)
void globalInterpClosed (const Vector< Point_nD< T, N > > &Qw, int d)
void globalInterpClosed (const Vector< Point_nD< T, N > > &Qw, const Vector< T > &ub, int d)
void globalInterpClosedH (const Vector< HPoint_nD< T, N > > &Qw, int d)
void globalInterpClosedH (const Vector< HPoint_nD< T, N > > &Qw, const Vector< T > &U, int d)
void globalInterpClosedH (const Vector< HPoint_nD< T, N > > &Qw, const Vector< T > &ub, const Vector< T > &U, int d)
void globalInterpClosed (const Vector< Point_nD< T, N > > &Qw, const Vector< T > &ub, const Vector< T > &Uc, int d)
void globalInterpD (const Vector< Point_nD< T, N > > &Q, const Vector< Point_nD< T, N > > &D, int d, int unitD, T a=1.0)
void projectTo (const Point_nD< T, N > &p, T guess, T &u, Point_nD< T, N > &r, T e1=0.001, T e2=0.001, int maxTry=100) const
length (T eps=0.001, int n=100) const
lengthIn (T us, T ue, T eps=0.001, int n=100) const
lengthF (T) const
lengthF (T, int) const
void makeCircle (const Point_nD< T, N > &O, const Point_nD< T, N > &X, const Point_nD< T, N > &Y, T r, double as, double ae)
void makeCircle (const Point_nD< T, N > &O, T r, double as, double ae)
void makeCircle (const Point_nD< T, N > &O, T r)
void makeLine (const Point_nD< T, N > &P0, const Point_nD< T, N > &P1, int d)
virtual void degreeElevate (int t)
void decompose (NurbsCurveArray< T, N > &c) const
void decomposeClosed (NurbsCurveArray< T, N > &c) const
int splitAt (T u, NurbsCurve< T, N > &cl, NurbsCurve< T, N > &cu) const
int mergeOf (const NurbsCurve< T, N > &cl, const NurbsCurve< T, N > &cu)
void transform (const MatrixRT< T > &A)
void modCP (int i, const HPoint_nD< T, N > &a)
void modCPby (int i, const HPoint_nD< T, N > &a)
virtual void modKnot (const Vector< T > &knotU)
int movePoint (T u, const Point_nD< T, N > &delta)
int movePoint (T u, const BasicArray< Point_nD< T, N > > &delta)
int movePoint (const BasicArray< T > &ur, const BasicArray< Point_nD< T, N > > &D)
int movePoint (const BasicArray< T > &ur, const BasicArray< Point_nD< T, N > > &D, const BasicArray_INT &Dr, const BasicArray_INT &Dk)
int movePoint (const BasicArray< T > &ur, const BasicArray< Point_nD< T, N > > &D, const BasicArray_INT &Dr, const BasicArray_INT &Dk, const BasicArray_INT &fixCP)
int setTangent (T u, const Point_nD< T, N > &T0)
int setTangentAtEnd (const Point_nD< T, N > &T0, const Point_nD< T, N > &T1)
int read (const char *)
int write (const char *) const
virtual int read (ifstream &fin)
int write (ofstream &fout) const
int writePS (const char *, int cp=0, T magFact=T(-1), T dash=T(5), bool bOpen=true) const
int writePSp (const char *, const Vector< Point_nD< T, N > > &, const Vector< Point_nD< T, N > > &, int cp=0, T magFact=0.0, T dash=5.0, bool bOpen=true) const
int writeVRML (ostream &fout, T radius, int K, const Color &color, int Nu, int Nv, T u_s, T u_e) const
int writeVRML (const char *filename, T radius, int K, const Color &color, int Nu, int Nv, T u_s, T u_e) const
int writeVRML (const char *filename, T radius=1, int K=5, const Color &color=whiteColor, int Nu=20, int Nv=20) const
int writeVRML (ostream &fout, T radius=1, int K=5, const Color &color=whiteColor, int Nu=20, int Nv=20) const
int writeVRML97 (const char *filename, T radius, int K, const Color &color, int Nu, int Nv, T u_s, T u_e) const
int writeVRML97 (ostream &fout, T radius, int K, const Color &color, int Nu, int Nv, T u_s, T u_e) const
int writeVRML97 (const char *filename, T radius=1, int K=5, const Color &color=whiteColor, int Nu=20, int Nv=20) const
int writeVRML97 (ostream &fout, T radius=1, int K=5, const Color &color=whiteColor, int Nu=20, int Nv=20) const
int writeDisplayLINE (const char *filename, int iNu, const Color &color=blueColor, T fA=1) const
int writeDisplayLINE (const char *filename, const Color &color, int iNu, T u_s, T u_e) const
void drawImg (Image_UBYTE &Img, unsigned char color=255, T step=0.01)
void drawImg (Image_Color &Img, const Color &color, T step=0.01)
void drawAaImg (Image_Color &Img, const Color &color, int precision=3, int alpha=1)
void drawAaImg (Image_Color &Img, const Color &color, const NurbsCurve< T, 3 > &profile, int precision=3, int alpha=1)
NurbsSurface< T, 3 > drawAaImg (Image_Color &Img, const Color &color, const NurbsCurve< T, 3 > &profile, const NurbsCurve< T, 3 > &scaling, int precision=3, int alpha=1)
BasicList< Point_nD< T, N > > tesselate (T tolerance, BasicList< T > *uk) const

Protected Attributes

Vector< HPoint_nD< T, N > > P
Vector< T > U
int deg_

Friends

HPoint_nD< T, N > C (T u, const NurbsCurve< T, N > &nurb)
Point_nD< T, N > Cp (T u, const NurbsCurve< T, N > &nurb)
void generateCompatibleCurves (NurbsCurveArray< T, N > &ca)

Related Functions

(Note that these are not member functions.)

chordLengthParam (const Vector< Point_nD< T, N > > &Q, Vector< T > &ub)

Detailed Description

template<class T, int N>
class PLib::NurbsCurve< T, N >

A NURBS curve class.

This class is used to represent and manipulate NURBS curve. The curves are composed of points in 4D. They can have any degree and have any number of control points.

Author:
Philippe Lavoie
Date:
4 Oct. 1996


Constructor & Destructor Documentation

template<class T, int N>
PLib::NurbsCurve< T, N >::NurbsCurve  
 

default constructor.

Author:
Philippe Lavoie
Date:
24 January 1997

template<class T, int N>
PLib::NurbsCurve< T, N >::NurbsCurve const NurbsCurve< T, N > &    nurb
 

A copy constructor.

Parameters:
nurb  the NURBS curve to copy
Author:
Philippe Lavoie
Date:
24 January 1997

template<class T, int N>
PLib::NurbsCurve< T, N >::NurbsCurve const Vector< HPoint_nD< T, N > > &    P1,
const Vector< T > &    U1,
int    Degree = 3
 

Constructor with control points in 4D.

Parameters:
P1  the control points
U1  the knot vector
Degree  the degree of the curve
Warning:
The size of P1,U1 and Degree must agree: P.n()+degree+1=U.n()
Author:
Philippe Lavoie
Date:
24 January 1997

template<class T, int N>
PLib::NurbsCurve< T, N >::NurbsCurve const Vector< Point_nD< T, N > > &    P1,
const Vector< T > &    W,
const Vector< T > &    U1,
int    Degree = 3
 

Constructor with control points in 3D.

Parameters:
P1  --> the control point vector
W  --> the weight for each control points
U1  --> the knot vector
Degree  --> the degree of the curve
Warning:
The size of P1,U1 and Degree must agree: P.n()+degree+1=U.n()
Author:
Philippe Lavoie
Date:
24 January 1997


Member Function Documentation

template<class T, int D>
T PLib::NurbsCurve< T, D >::basisFun   u,
int    i,
int    p = -1
const
 

Computes the basis function of the curve.

Computes the i basis function of degree p of the curve at parameter u.

You can have more information about this function in the LaTeX version.

Parameters:
u  the parametric variable
i  specifies which basis function to compute
p  the degree to which the basis function is computed
Returns:
the value of N_{ip}(u)
Author:
Philippe Lavoie
Date:
24 January 1997

template<class T, int D>
void PLib::NurbsCurve< T, D >::basisFuns   u,
int    i,
Vector< T > &    N
const
 

computes the non-zero basis functions of the curve.

Computes the non-zero basis functions and puts the result into N. N has a size of deg+1. To relate N to the basis functions, Basis[span -deg +i] = N[i] for i=0...deg.

You can find more information in the LaTeX version.

Parameters:
u  the parametric value
i  the non-zero span of the basis functions
N  the non-zero basis functions
Warning:
u and i must be valid values
Author:
Philippe Lavoie
Date:
24 January 1997

template<class T, int N>
void PLib::NurbsCurve< T, N >::clamp  
 

clamp a NURBS curve.

A clamped NURBS curve has degree+1 equal knots at both ends of the knot vector.

Author:
Philippe Lavoie
Date:
27 April 1999

template<class T, int N>
void PLib::NurbsCurve< T, N >::decompose NurbsCurveArray< T, N > &    c const
 

Decompose the curve into Bézier segments.

This function decomposes the curve into an array of 4D Bézier segments.

Parameters:
c  an array of Bézier segments
Warning:
The end Bézier segments will not be valid if the NURBS curve is not clamped.
Author:
Philippe Lavoie
Date:
16 February 1997

template<class T, int D>
void PLib::NurbsCurve< T, D >::decomposeClosed NurbsCurveArray< T, D > &    c const
 

decompose the closed curve into Bézier segments.

This function decomposes a closed curve into an array of Bézier segments.

Parameters:
c  an array of Bézier segments
Author:
Alejandro Frangi
Date:
30 July 1998

template<class T, int N>
void PLib::NurbsCurve< T, N >::degreeElevate int    t [virtual]
 

degree elevate a curve a number of times.

For more information, see A5.9 on p 206 of the NURBS book

Parameters:
t  the number of times to increase the degree of the curve
Author:
Philippe Lavoie
Date:
24 January, 1997

Reimplemented in PLib::NurbsCurveSP, and PLib::NurbsCurveSP< float, 3 >.

template<class T, int N>
HPoint_nD< T, N > PLib::NurbsCurve< T, N >::derive   u,
int    d
const
 

Computes the derivative of degree of the curve at parameter u.

For more information on the algorithm used, see A3.2 on p 93 of the NurbsBook.

Parameters:
u  the parametric value to evaluate at
d  the degree of the derivative
Returns:
The derivative d in 4D at the parameter u
Warning:
u and d must be in a valid range.
Author:
Philippe Lavoie
Date:
24 January, 1997

template<class T, int N>
Point_nD< T, N > PLib::NurbsCurve< T, N >::derive3D   u,
int    d
const
 

Computes the derivative of degree d of the curve at parameter u.

For more information on the algorithm used, see A3.2 on p 93 of the NurbsBook.

Parameters:
u  the parametric value to evaluate at
d  the degree of the derivative
Returns:
The derivative d in norma space at the parameter u
Warning:
u and d must be in a valid range.
Author:
Philippe Lavoie
Date:
24 January, 1997

template<class T, int N>
void PLib::NurbsCurve< T, N >::deriveAt   u,
int    d,
int    span,
Vector< Point_nD< T, N > > &    ders
const
 

Computes the derivative of the curve at the parameter u.

Parameters:
u  the parameter at which the derivative is computed
d  the degree of derivation
span  the span of u.
ders  the vector containing the derivatives of the point at u.
Warning:
u and $d$ must be in a valid range.
Author:
Philippe Lavoie
Date:
9 October 1998

template<class T, int N>
void PLib::NurbsCurve< T, N >::deriveAt   u,
int    d,
Vector< Point_nD< T, N > > &    ders
const [virtual]
 

Computes the derivative at the parameter u.

Parameters:
u  the parameter at which the derivative is computed
d  the degree of derivation
ders  the vector containing the derivatives of the point at u.
\wanring u and d must be in a valid range.
Author:
Philippe Lavoie
Date:
24 January, 1997

Reimplemented from PLib::ParaCurve.

template<class T, int N>
void PLib::NurbsCurve< T, N >::deriveAtH   u,
int    d,
int    span,
Vector< HPoint_nD< T, N > > &    ders
const
 

Computes the derivative of degree d of the curve at parameter u.

For more information on the algorithm used, see A3.2 on p 93 of the NurbsBook.

Parameters:
u  the parametric value to evaluate at
d  the degree of the derivative
span  the span of u
ders  a vector containing the derivatives of the curve at u.
Warning:
u and d must be in a valid range.
Author:
Philippe Lavoie
Date:
9 October, 1998

template<class T, int N>
void PLib::NurbsCurve< T, N >::deriveAtH   u,
int    d,
Vector< HPoint_nD< T, N > > &    ders
const [virtual]
 

Computes the derivative of degree d of the curve at parameter u in the homonegeous domain.

For more information on the algorithm used, see A3.2 on p 93 of the NurbsBook.

Parameters:
u  the parametric value to evaluate at
d  the degree of the derivative
ders  a vector containing the derivatives of the curve at u.
Warning:
u and d must be in a valid range.
Author:
Philippe Lavoie
Date:
24 January, 1997

Reimplemented from PLib::ParaCurve.

template<class T, int N>
void PLib::NurbsCurve< T, N >::dersBasisFuns int    n,
  u,
int    span,
Matrix< T > &    ders
const
 

Compute the derivatives functions at u of the NURBS curve.

For information on the algorithm, see A2.3 on p72 of the NURBS book.

The result is stored in the ders matrix, where ders is of size (n+1,deg+1) and the derivative N'_i(u) = ders(1,i=span-deg+j) where j=0...deg+1.

Parameters:
n  the degree of the derivation
u  the parametric value
span  the span for the basis functions
ders  A matrix containing the derivatives of the curve.
Warning:
n, u and span must be valid values.
Author:
Philippe Lavoie
Date:
24 January 1997

template<class T, int N>
NurbsSurface< T, 3 > PLib::NurbsCurve< T, N >::drawAaImg Image_Color &    Img,
const Color   color,
const NurbsCurve< T, 3 > &    profile,
const NurbsCurve< T, 3 > &    scaling,
int    precision = 3,
int    alpha = 1
 

Draws an anti-aliased NURBS curve on an image.

This will draw the NURBS by using a brush profile. The drawing is performed by averaging the intensity of the profile at the pixels.

This function generates a sweep surface by using the profile given in its argument. The sweep is always performed by following the y-axis of the profile. A scaling function is also used when sweeping. This is used to vary the shape of the profile while it's being swept (see the sweep member function of NurbsSurface<T,N> for more details).

Parameters:
Img  draws the nurbs curve to this Image
color  the line is drawn in this color
profile  the profile of the NURBS curve to draw
scaling  the scaling to give the profile while drawing the curve
precision  this number influences the number of points used for averaging purposes.
alpha  a flag indicating if the profile is used as an alpha chanel. If so, the line doesn't overwrite, it blends the line with the image already present in Img.
Warning:
This routine is very slow; use normal drawing for speed or lower the precision factor.
Author:
Philippe Lavoie
Date:
25 July 1997

template<class T, int N>
void PLib::NurbsCurve< T, N >::drawAaImg Image_Color &    Img,
const Color   color,
const NurbsCurve< T, 3 > &    profile,
int    precision = 3,
int    alpha = 1
 

draws an anti-aliased NURBS curve on an image.

This will draw the NURBS by using a user-defined brush profile. The drawing is performed by averaging the intensity of the profile at the pixels.

Parameters:
Img  draws the nurbs curve to this Image
color  the line is drawn in this color
profile  the profile of the NURBS curve to draw
precision  this number influences the number of points used for averaging purposes.
alpha  a flag indicating if the profile is used as an alpha chanel. If so, the line doesn't overwrite, it blends the line with the image already present in Img.
Warning:
This routine is very slow; use normal drawing for speed.
Author:
Philippe Lavoie
Date:
22 August 1997

template<class T, int N>
void PLib::NurbsCurve< T, N >::drawAaImg Image_Color &    Img,
const Color   color,
int    precision = 3,
int    alpha = 1
 

Draws an anti-aliased NURBS curve on an image.

This will draw the NURBS by using a circular brush profile. The drawing is performed by averaging the intensity of the profile at the pixels.

Parameters:
Img  draws the nurbs curve to this Image
color  the line is drawn in this color
precision  this number influences the number of points used for averaging purposes.
alpha  a flag indicating if the profile is used as an alpha chanel. If so, the line doesn't overwrite, it blends the line with the image already present in Img.
Warning:
This routine is very slow; use normal drawing for speed.
Author:
Philippe Lavoie
Date:
25 July 1997

template<class T, int N>
void PLib::NurbsCurve< T, N >::drawImg Image_Color &    Img,
const Color   color,
  step = 0.01
 

Draws a NURBS curve on an image.

This will draw very primitively the NURBS curve on an image. The drawing assumes the line is only in the xy plane (the z is not used for now).

The algorithm finds the points on the curve at a step parametric intervall between them and join them by a line. No fancy stuff.

Parameters:
Img  draws the nurbs curve to this Image
color  the line is drawn in this color
step  the parametric distance between two computed points.
Author:
Philippe Lavoie
Date:
24 January 1997

template<class T, int N>
void PLib::NurbsCurve< T, N >::drawImg Image_UBYTE &    Img,
unsigned char    color = 255,
  step = 0.01
 

draws a NURBS curve on an image.

This will draw very primitively the NURBS curve on an image. The drawing assumes the line is only in the xy plane (the z is not used for now).

The algorithm finds the points on the curve at a step parametric intervall between them and join them by a line. No fancy stuff.

Parameters:
Img  <-- draws the nurbs curve to this Image
color  --> the line is drawn in this color
step  --> the parametric distance between two computed points.
Author:
Philippe Lavoie
Date:
24 January 1997

template<class T, int N>
int PLib::NurbsCurve< T, N >::findKnot   u const
 

Finds the knot k for which u is in the range [u_k,u_{k+1}).

Parameters:
u  parametric value
Returns:
the index k
Warning:
u must be in a valid range.
Author:
Philippe Lavoie
Date:
24 January, 1997

template<class T, int N>
int PLib::NurbsCurve< T, N >::findMult int    r const
 

Finds the multiplicity of a knot.

Parameters:
r  the knot to observe
Returns:
the multiplicity of the knot
Warning:
r must be a valid knot index
Author:
Philippe Lavoie
Date:
24 January, 1997

template<class T, int N>
void PLib::NurbsCurve< T, N >::findMultSpan   u,
int &    r,
int &    s
const
 

Finds the multiplicity of a knot at a parametric value.

Finds the index of the knot at parametric value u and returns its multiplicity.

Parameters:
u  the parametric value
r  the knot of interest
s  the multiplicity of this knot
Warning:
u must be in a valid range.
Author:
Philippe Lavoie
Date:
24 January, 1997

template<class T, int N>
int PLib::NurbsCurve< T, N >::findSpan   u const
 

Determines the knot span index.

Determines the knot span for which their exists non-zero basis functions. The span is the index k for which the parameter u is valid in the [u_k,u_{k+1}] range.

Parameters:
u  the parametric value
Returns:
the span index at u.
Warning:
u must be in a valid range
Author:
Philippe Lavoie
Date:
24 January 1997 \modified 20 January, 1999 (Alejandro Frangi)

template<class T, int D>
HPoint_nD< T, D > PLib::NurbsCurve< T, D >::firstD   u,
int    span
const
 

Computes the first derivative.

Computes the first derivative in the honogenous space with the span given.

Parameters:
u  compute the derivative at this parameter
span  the span of u
Returns:
The first derivative of the point in the homogoneous space
Warning:
u and span must be in a valid range
Author:
Philippe Lavoie
Date:
13 October 1998

template<class T, int D>
HPoint_nD< T, D > PLib::NurbsCurve< T, D >::firstD   u const
 

Computes the first derivative Computes the first derivative in the 4D homogenous space.

Parameters:
u  --> compute the derivative at this parameter
Returns:
The first derivative in homogenous space
Warning:
u must be in the valid range
Author:
Philippe Lavoie
Date:
13 October 1998

template<class T, int N>
Point_nD< T, N > PLib::NurbsCurve< T, N >::firstDn   u,
int    span
const
 

Computes the first derivative.

Computes the first derivative in the normal space (3D or 2D).

Parameters:
u  --> compute the derivative at this parameter
span  --> the span of u
Warning:
u and span must be in a valid range
Author:
Philippe Lavoie
Date:
13 October 1998

template<class T, int N>
Point_nD< T, N > PLib::NurbsCurve< T, N >::firstDn   u const
 

Computes the first derivative.

Computes the first derivative in the normal space.

Parameters:
u  compute the derivative at this parameter
span  the span of u
Returns:
The first derivative in normal space
Warning:
u and span must be in a valid range
Author:
Philippe Lavoie
Date:
13 October 1998

template<class T, int N>
T PLib::NurbsCurve< T, N >::getRemovalBnd int    r,
int    s
const
 

Get the knot removal error bound for an internal knot.

Get the knot removal error bound for an internal knot r (non-rational). For more information on the algorithm, see A9.8 from the Nurbs book on page 428.

Parameters:
curve  a NURBS curve
r  the index of the internal knot to check
s  the multiplicity of that knot
Returns:
The maximum distance between the new curve and the old one
Author:
Philippe Lavoie
Date:
24 January, 1997

template<class T, int N>
void PLib::NurbsCurve< T, N >::globalApproxErrBnd Vector< Point_nD< T, N > > &    Q,
Vector< T > &    ub,
int    degC,
  E
 

Approximation of a curve bounded to a certain error.

It is a type II approximation: it starts with a lot of control points then tries to eliminate as much as it can as long as the curve stays within a certain error bound.

The method uses least squares fitting along with knot removal techniques. It is the algorithm A9.10 on p 431 of the NURBS book.

Parameters:
Q  the points to approximate
ub  the vector of parameters where the points are located
degree  the degree of the approximation curve
E  the maximum error allowed
Warning:
ub and Q must be of the same size
Author:
Philippe Lavoie
Date:
24 January 1997

template<class T, int N>
void PLib::NurbsCurve< T, N >::globalApproxErrBnd Vector< Point_nD< T, N > > &    Q,
int    degC,
  E
 

Approximation of a curve bounded to a certain error.

It is a type II approximation: it starts with a lot of control points then tries to eliminate as much as it can as long as the curve stays within a certain error bound.

The method uses least squares fitting along with knot removal techniques. It is the algorithm A9.10 on p 431 of the NURBS book.

Parameters:
Q  the points to approximate
degree  the degree of the approximation curve
E  the maximum error allowed
Author:
Philippe Lavoie
Date:
24 January 1997

template<class T, int N>
void PLib::NurbsCurve< T, N >::globalApproxErrBnd2 Vector< Point_nD< T, N > > &    Q,
int    degC,
  E
 

Approximation of a curve bounded to a certain error.

The algorithm is quite simplistic but in some cases gives better result than globalApproxErrBnd (when the error allowed is low and the data is very close to each other).

This algorithm generates a first degree interpolation of the data points, degree elevates the curve by 1 degree recomputes the error around each points, removes the knots which are within a certain error range then repeats the process until the desired degree is reached.

Parameters:
Q  the points to approximate
degC  the degree of the approximation curve
E  the maximum error allowed
Author:
Philippe Lavoie
Date:
24 January 1997

template<class T, int N>
void PLib::NurbsCurve< T, N >::globalApproxErrBnd3 Vector< Point_nD< T, N > > &    Q,
const Vector< T > &    ub,
int    degC,
  E
 

Approximation of a curve bounded to a certain error.

The algorithm is quite simplistic but in some cases gives better result than globalApproxErrBnd (when the error allowed is low and the data is very close to each other).

This algorithm generates a first degree interpolation of the data points, degree elevates it to the degree requested then removes all the control points which are within the error bound.

Parameters:
Q  the points to approximate
ub  the parameter vector
degC  the degree of the approximation curve
E  the maximum error allowed
Author:
Philippe Lavoie
Date:
24 January 1997

template<class T, int N>
void PLib::NurbsCurve< T, N >::globalApproxErrBnd3 Vector< Point_nD< T, N > > &    Q,
int    degC,
  E
 

Approximation of a curve bounded to a certain error.

The algorithm is quite simplistic but in some cases gives better result than globalApproxErrBnd (when the error allowed is low and the data is very close to each other).

This algorithm generates a first degree interpolation of the data points, degree elevates it to the degree requested then removes all the control points which are within the error bound.

Parameters:
Q  the points to approximate
degC  the degree of the approximation curve
E  the maximum error allowed
Author:
Philippe Lavoie
Date:
24 January 1997

template<class T, int D>
void PLib::NurbsCurve< T, D >::globalInterp const Vector< Point_nD< T, D > > &    Q,
const Vector< T > &    ub,
int    d
 

global curve interpolation with points in 3D.

Global curve interpolation with points in 3D and with the parametric values specified.

Parameters:
Q  the 3D points to interpolate
ub  the parametric values
d  the degree of the interpolation
Warning:
The number of points to interpolate must be greater than the degree specified for the curve.
Author:
Philippe Lavoie
Date:
24 January, 1997

template<class T, int N>
void PLib::NurbsCurve< T, N >::globalInterp const Vector< Point_nD< T, N > > &    Q,
int    d
 

global curve interpolation with points in 3D.

Parameters:
Q  the 3D points to interpolate
d  the degree of the interpolation
Warning:
The number of points to interpolate must be greater than the degree specified for the curve.
Author:
Philippe Lavoie
Date:
24 January, 1997

template<class T, int D>
void PLib::NurbsCurve< T, D >::globalInterpClosed const Vector< Point_nD< T, D > > &    Qw,
const Vector< T > &    ub,
const Vector< T > &    Uc,
int    d
 

global curve interpolation with homogenous points.

Global curve interpolation with 4D points, a knot vector defined and the parametric value vector defined.The curve will have C(d-1) continuity at the point u=0 and u=1.

Parameters:
Qw  the 3D points to interpolate (wrapped around)
ub  the parametric values vector
Uc  the knot vector computed using knotAveragingC
d  the degree of the closed curve
Warning:
The number of points to interpolate must be greater than the degree specified for the curve. Uc must be compatible with the values given for Q.n(), ub.n().
Author:
Alejandro Frangi
Date:
13 July, 1998

template<class T, int D>
void PLib::NurbsCurve< T, D >::globalInterpClosed const Vector< Point_nD< T, D > > &    Qw,
const Vector< T > &    ub,
int    d
 

global curve interpolation with homogenous points.

Global curve interpolation with 4D points, a knot vector defined and the parametric value vector defined.

Parameters:
Q  the 3D points to interpolate
ub  the parametric values vector
d  the degree of the closed curve
Warning:
The number of points to interpolate must be greater than the degree specified for the curve. Uc must be compatible with the values given for Q.n(), ub.n() and d.
Author:
Alejandro Frangi
Date:
13 July, 1998

template<class T, int N>
void PLib::NurbsCurve< T, N >::globalInterpClosed const Vector< Point_nD< T, N > > &    Qw,
int    d
 

global closed curve interpolation with a list of points.

Global curve interpolation with points in 3D. This function will generate a closed curve with C(d-1) continuity between the parameters u=0 and u=1

Parameters:
Q  the 3D points to interpolate
d  the degree of the interpolation
Warning:
The number of points to interpolate must be greater than the degree specified for the curve.
Author:
Alejandro Frangi
Date:
13 July, 1998

template<class T, int D>
void PLib::NurbsCurve< T, D >::globalInterpClosedH const Vector< HPoint_nD< T, D > > &    Qw,
const Vector< T > &    ub,
const Vector< T > &    Uc,
int    d
 

global curve interpolation with homogenous points.

Global curve interpolation with 4D points, a knot vector defined and the parametric value vector defined.The curve will have C(d-1) continuity at the point u=0 and u=1.

Parameters:
Qw  the 3D points to interpolate (wrapped around)
ub  the parametric values vector
Uc  the knot vector computed using knotAveragingClosed
d  the degree of the closed curve
Warning:
The number of points to interpolate must be greater than the degree specified for the curve. Uc must be compatible with the values given for Q.n(), ub.n().
Author:
Alejandro Frangi
Date:
13 July, 1998

template<class T, int D>
void PLib::NurbsCurve< T, D >::globalInterpClosedH const Vector< HPoint_nD< T, D > > &    Qw,
const Vector< T > &    ub,
int    d
 

global curve interpolation with homogenous points.

Global curve interpolation with 4D points, a knot vector defined and the parametric value vector defined.

Parameters:
Q  the 3D points to interpolate
ub  the parametric values vector
d  the degree of the closed curve
Warning:
The number of points to interpolate must be greater than the degree specified for the curve. Uc must be compatible with the values given for Q.n(), ub.n() and d.
Author:
Alejandro Frangi
Date:
13 July, 1998

template<class T, int D>
void PLib::NurbsCurve< T, D >::globalInterpClosedH const Vector< HPoint_nD< T, D > > &    Qw,
int    d
 

global close curve interpolation with points in homogenous space.

Global curve interpolation with points in homogenouse space with C(d-1) continuity in the wrap-around point.

Parameters:
Qw  the points in 4D to interpolate
d  the degree of the closed curve
Warning:
The number of points to interpolate must be greater than the degree specified for the curve. The interpolation degree is only 3. The first and last interpolation points should be equal.
Author:
Alejandro Frangi
Date:
13 July, 1998

template<class T, int nD>
void PLib::NurbsCurve< T, nD >::globalInterpD const Vector< Point_nD< T, nD > > &    Q,
const Vector< Point_nD< T, nD > > &    D,
int    d,
int    unitD,
  a = 1.0
 

global curve interpolation with the 1st derivatives of the points specified.

Global curve interpolation with the 1st degree derivative specified for each of the points to interpolate. This will generate a number of control points 2 times greater than the number of interpolation points.

If the derivative specified is of unit length, i.e. the tangent vectors are given, then the derivative at each point will be multiplied by the chord length. A second multiplicative factor can be specified to make the derivative even greater. This second number should be close to 1.0.

For more information about the algorithm, refer to section 9.2.4 on page 373 of the NURBS book.

Parameters:
Q  the points to interpolate
D  the first derivative at these points
d  the degree of the curve
unitD  a flag specifying if the derivatives are unit vectors
a  a multiplicative factor for the tangent (must be greater than 0 or the function aborts).
Warning:
The number of points to interpolate must be greater than the degree specified for the curve.
Author:
Philippe Lavoie
Date:
3 September, 1997

template<class T, int D>
void PLib::NurbsCurve< T, D >::globalInterpH const Vector< HPoint_nD< T, D > > &    Q,
const Vector< T > &    ub,
const Vector< T > &    Uc,
int    d
 

global curve interpolation with 4D points, a knot vector defined and the parametric value vector defined.

Global curve interpolation with 4D points, a knot vector defined and the parametric value vector defined.

Parameters:
Q  the 3D points to interpolate
ub  the parametric values vector
Uc  the knot vector to set the curve to
d  the degree of the interpolation
Warning:
The number of points to interpolate must be greater than the degree specified for the curve. Uc must be compatible with the values given for Q.n(), ub.n() and d.
Author:
Philippe Lavoie
Date:
24 January, 1997

template<class T, int D>
void PLib::NurbsCurve< T, D >::globalInterpH const Vector< HPoint_nD< T, D > > &    Q,
const Vector< T > &    Uc,
int    d
 

global curve interpolation with 4D points and a knot vector defined.

Global curve interpolation with 4D points and a knot vector defined.

Parameters:
Q  the 3D points to interpolate
Uc  the knot vector to set the curve to
d  the degree of the interpolation
Warning:
The number of points to interpolate must be greater than the degree specified for the curve.
Author:
Philippe Lavoie
Date:
24 January, 1997

template<class T, int D>
void PLib::NurbsCurve< T, D >::globalInterpH const Vector< HPoint_nD< T, D > > &    Q,
int    d
 

global curve interpolation with points in 4D.

Global curve interpolation with points in 4D

Parameters:
Q  the points in 4D to interpolate
d  the degree of the curve interpolation
Warning:
The number of points to interpolate must be greater than the degree specified for the curve.
Author:
Philippe Lavoie
Date:
24 January, 1997

template<class T, int N>
HPoint_nD< T, N > PLib::NurbsCurve< T, N >::hpointAt   u,
int    span
const [virtual]
 

Evaluates the curve in homogenous space at parameter u.

For more details on the algorithm, see A4.1 on page 124 of the Nurbs Book.

Parameters:
u  the parametric value at which the curve is evaluated
span  the span of u
Returns:
the 4D point at C(u)
Warning:
the parametric value must be in a valid range
Author:
Philippe Lavoie
Date:
24 January, 1997

Reimplemented from PLib::ParaCurve.

template<class T, int N>
int PLib::NurbsCurve< T, N >::knotInsertion   u,
int    r,
NurbsCurve< T, N > &    nc
 

It inserts a knot a number of times.

It inserts the knot u, r times and generates the curve nc. For more information, see A5.1 on page 151 of the NURBS book

Parameters:
u  the knot to insert
r  the number of times to insert u.
nc  the resulting NURBS curve
Returns:
the number of knots inserted (0 if none)
Warning:
the routine throws a NurbsError if the u value is not inside the clamped range.
Author:
Philippe Lavoie
Date:
24 January 1997

template<class T, int D>
int PLib::NurbsCurve< T, D >::leastSquares const Vector< Point_nD< T, D > > &    Q,
int    degC,
int    n,
const Vector< T > &    ub,
const Vector< T > &    knot
 

A least squares curve approximation.

This routines generates a curve that approrimates the points in the least square sense, you can find more details in the LaTeX version.

For more details, see section 9.4.1 on page 491 of the NURBS book.

Parameters:
Q  the vector of 3D points
degC  the degree of the curve
n  the number of control points in the new curve
ub  the knot coefficients
knot  the knot vector to use for the curve
Returns:
1 if succesfull, 0 it the number of points to approximate the curve with is too big compared to the number of points.
Warning:
the variable curve must contain a valid knot vector.
Author:
Philippe Lavoie
Date:
24 January 1997

template<class T, int N>
int PLib::NurbsCurve< T, N >::leastSquares const Vector< Point_nD< T, N > > &    Q,
int    degC,
int    n,
const Vector< T > &    ub
 

A least squares curve approximation.

This routines generates a curve that approrimates the points in the least square sense, you can find more details in the LaTeX version.

For more details, see section 9.4.1 on page 491 of the NURBS book.

Parameters:
Q  the vector of 3D points
degC  the degree of the curve
n  the number of control points in the new curve
ub  the knot coefficients
Warning:
the variable curve must contain a valid knot vector.
Author:
Philippe Lavoie
Date:
24 January 1997

template<class T, int N>
int PLib::NurbsCurve< T, N >::leastSquares const Vector< Point_nD< T, N > > &    Q,
int    degC,
int    n
 

A least squares curve approximation.

This routines generates a curve that approrimates the points in the least square sense, you can find more details in the LaTeX version.

For more details, see section 9.4.1 on page 491 of the NURBS book.

Parameters:
Q  the vector of 3D points
degC  the degree of the curve
n  the number of control points in the new curve.
Warning:
deg must be smaller than Q.n().
Author:
Philippe Lavoie
Date:
24 January, 1997

template<class T, int D>
int PLib::NurbsCurve< T, D >::leastSquaresClosed const Vector< Point_nD< T, D > > &    Qw,
int    degC,
int    nCP,
const Vector< T > &    ub,
const Vector< T > &    knots
 

A least squares curve approximation for closed curves.

This routine finds a closed NURBS curve that satisfy a least square criteria. The resulting curve will generally not pass through the input points (except the first point).

For more details, see section 9.4.1 on page 491 of the NURBS book.

Parameters:
Qw  the vector of 4D points (wrapped around)
degC  the degree of the curve
nCP  the number of (distinct) control points in the new curve
ub  the parameter values of Qw obtained via chordLengthParamClosed()
knots  the knots for the control points obtained via knotApproximationClosed()
Author:
Alejandro Frangi
Date:
30 July 1998

template<class T, int N>
int PLib::NurbsCurve< T, N >::leastSquaresClosed const Vector< Point_nD< T, N > > &    Qw,
int    degC,
int    nCP,
const Vector< T > &    ub
 

A least squares curve approximation for closed curves.

This routine finds a closed NURBS curve that satisfy a least square criteria. The resulting curve will generally not pass through the input points (except the first point).

For more details, see section 9.4.1 on page 491 of the NURBS book.

Parameters:
Qw  the vector of 3D points (wrapped around)
degC  the degree of the curve
nCP  the number of (distinct) control points in the new curve
ub  the parameter values of Qw obtained via chordLengthParamClosed()
Author:
Alejandro Frangi
Date:
30 July 1998

template<class T, int N>
int PLib::NurbsCurve< T, N >::leastSquaresClosed const Vector< Point_nD< T, N > > &    Qw,
int    degC,
int    nCP
 

A least squares curve approximation for closed curves.

This routine finds a closed NURBS curve that satisfy a least square criteria. The resulting curve will generally not pass through the input points (except the first point).

For more details, see section 9.4.1 on page 491 of the NURBS book.

Parameters:
Qw  the vector of 3D points (wrapped around)
degC  the degree of the curve
nCP  the number of (distinct) control points in the new curve
Author:
Alejandro Frangi
Date:
30 July 1998

template<class T, int D>
int PLib::NurbsCurve< T, D >::leastSquaresClosedH const Vector< HPoint_nD< T, D > > &    Qw,
int    degC,
int    nCP,
const Vector< T > &    ub,
const Vector< T > &    knots
 

A least squares curve approximation for closed curves.

This routine finds a closed NURBS curve that satisfy a least square criteria. The resulting curve will generally not pass through the input points (except the first point).

For more details, see section 9.4.1 on page 491 of the NURBS book.

Parameters:
Qw  the vector of 4D points (wrapped around)
degC  the degree of the curve
nCP  the number of (distinct) control points in the new curve
ub  the parameter values of Qw obtained via chordLengthParamClosed()
knots  the knots for the control points obtained via knotApproximationC()
Author:
Alejandro Frangi
Date:
30 July 1998

template<class T, int N>
int PLib::NurbsCurve< T, N >::leastSquaresClosedH const Vector< HPoint_nD< T, N > > &    Qw,
int    degC,
int    nCP,
const Vector< T > &    ub
 

A least squares curve approximation for closed curves.

This routine finds a closed NURBS curve that satisfy a least square criteria. The resulting curve will generally not pass through the input points (except the first point).

For more details, see section 9.4.1 on page 491 of the NURBS book.

Parameters:
Qw  the vector of 4D points (wrapped around)
degC  the degree of the curve
nCP  the number of (distinct) control points in thenew curve
ub  the parameter values of Qw obtained via chordLengthParamClosed()
Author:
Alejandro Frangi
Date:
30 July 1998

template<class T, int D>
int PLib::NurbsCurve< T, D >::leastSquaresH const Vector< HPoint_nD< T, D > > &    Q,
int    degC,
int    n,
const Vector< T > &    ub,
const Vector< T > &    knot
 

A least squares curve approximation.

This routines generates a curve that approrimates the points in the least square sense, you can find more details in the LaTeX version.

For more details, see section 9.4.1 on page 491 of the NURBS book.

Parameters:
Q  the vector of 4D points
degC  the degree of the curve
n  the number of control points in the new curve
ub  the knot coefficients
knot  the knot vector to use for the curve
Returns:
1 if succesfull, 0 it the number of points to approximate the curve with is too big compared to the number of points.
Warning:
the variable curve must contain a valid knot vector.
Author:
Philippe Lavoie
Date:
24 January 1997

template<class T, int N>
int PLib::NurbsCurve< T, N >::leastSquaresH const Vector< HPoint_nD< T, N > > &    Q,
int    degC,
int    n,
const Vector< T > &    ub
 

A least squares curve approximation.

This routines generates a curve that approrimates the points in the least square sense, you can find more details in the LaTeX version.

For more details, see section 9.4.1 on page 491 of the NURBS book.

Parameters:
Q  the vector of 4D points
degC  the degree of the curve
n  the number of control points in the new curve
ub  the knot coefficients
Warning:
the variable curve must contain a valid knot vector.
Author:
Philippe Lavoie
Date:
24 January 1997

template<class T, int N>
T PLib::NurbsCurve< T, N >::length   eps = 0.001,
int    n = 100
const
 

Computes the length of the curve.

Computes an approximation of the length of the curve using a numerical automatic integrator.

That integrator uses a Chebyshev Series Expansion to perform its approximation. This is why you can change the value $n$ which sets the number of elements in the series.

The method is simple, integrate between each span. This is necessary in case the tangant of a point at u_i is undefined. Add the result and return this as the approximation.

Parameters:
n  the number of element in the Chebyshev series
eps  the accepted relative error
Returns:
the length of the NURBS curve.
Author:
Philippe Lavoie
Date:
22 September 1998

template<class T, int N>
T PLib::NurbsCurve< T, N >::lengthF   u,
int    span
const
 

length needs to integrate a function over an interval to determine the length of the NURBS curve. Well, this is the function.

Parameters:
u  the parameter
span  the span of the parameter
Returns:
square root of the square of the x,y and z value
Author:
Philippe Lavoie
Date:
22 September 1998

template<class T, int N>
T PLib::NurbsCurve< T, N >::lengthF   u const
 

The function used by length length needs to integrate a function over an interval to determine the length of the NURBS curve. Well, this is the function.

Parameters:
u  --> the parameter
Returns:
square root of the square of the x,y and z value
Author:
Philippe Lavoie
Date:
22 September 1998

template<class T, int N>
T PLib::NurbsCurve< T, N >::lengthIn   us,
  ue,
  eps = 0.001,
int    n = 100
const
 

Computes the length of the curve inside [u_s,u_e].

Computes an approximation of the length of the curve using a numerical automatic integrator. The length is computed for the range [u_s,u_e]

That integrator uses a Chebyshev Series Expansion to perform its approximation. This is why you can change the value n which sets the number of elements in the series.

The method is similar to the one used by length excepted that it needs to check for the range.

Parameters:
us  the starting range
ue  the end of the range
n  the number of element in the Chebyshev series
eps  the accepted relative error
Returns:
the length of the NURBS curve.
Warning:
ue must be greater than us and both must be in a valid range.
Author:
Philippe Lavoie
Date:
22 September 1998

template<class T, int N>
void PLib::NurbsCurve< T, N >::makeCircle const Point_nD< T, N > &    O,
  r
 

generates a circular curve.

Generates a circle of radius $r$ at origin $O$. The curve is drawn in the $xy$-axis.

Parameters:
O  the center of the circle
r  the radius of the circle
Author:
Philippe Lavoie
Date:
3 May, 1999

template<class T, int N>
void PLib::NurbsCurve< T, N >::makeCircle const Point_nD< T, N > &    O,
  r,
double    as,
double    ae
 

generates a circular curve.

Generates a circular curve of radius $r$ at origin $O$. The curve is drawn in the $xy$-axis.

Parameters:
O  the center of the circle
r  the radius of the circle
as  start angle measured with respect to the x-axis
ae  end angle measured with respect to the y-axis
Author:
Philippe Lavoie
Date:
25 July, 1997

template<class T, int N>
void PLib::NurbsCurve< T, N >::makeCircle const Point_nD< T, N > &    O,
const Point_nD< T, N > &    X,
const Point_nD< T, N > &    Y,
  r,
double    as,
double    ae
 

generates a circular curve.

Generates parts of a circle, starting at angle as and finishing at ae with a radius r and having the origin located at O. The X and Y vector describe the local x-axis and the local y-axis of the circle.

The degrees are specified in radians.

Parameters:
O  the center of the circle
X  unit length vector lying in the plane of the circle
Y  unit length vector lying in the plane of the circle
r  the radius of the circle
as  start angle in radians measured with respect to $X$
ae  end angle in radians measured with respect to $X$
Author:
Philippe Lavoie
Date:
25 July, 1997

template<class T, int N>
void PLib::NurbsCurve< T, N >::makeLine const Point_nD< T, N > &    P0,
const Point_nD< T, N > &    P1,
int    d
 

Generate a straight line.

Generate a straight line going from point P0 to point P1 of degree d.

Parameters:
P0  the beginning of the line
P1  the end of the line
d  the degree of the curve
Warning:
d must be greater or equal to 1
Author:
Philippe Lavoie
Date:
22 September 1998

template<class T, int N>
void PLib::NurbsCurve< T, N >::mergeKnotVector const Vector< T > &    Um
 

Merges the knot vector of a curve with another knot vector.

Will merge the Knot vector U with the one from the curve and it will refine the curve appropriately.

Parameters:
Um  the knot vector to merge with
Warning:
the knot U must be common with the one from the curve c
Author:
Philippe Lavoie
Date:
24 January, 1997

Reimplemented in PLib::NurbsCurveSP.

template<class T, int N>
int PLib::NurbsCurve< T, N >::mergeOf const NurbsCurve< T, N > &    cl,
const NurbsCurve< T, N > &    cu
 

The curve is the result of mergin two curves.

Parameters:
cl  the lower curve
cu  the upper curve
Returns:
1 if the operation is succesfull, 0 otherwise.
Warning:
You must make sure the the knot vectors are compatible, i.e. that the end knots of the lower curve are the same as the first knots of the upper curve. The curves must also be of the same degree.
Author:
Philippe Lavoie
Date:
16 October 1997

template<class T, int N>
int PLib::NurbsCurve< T, N >::movePoint const BasicArray< T > &    ur,
const BasicArray< Point_nD< T, N > > &    D,
const BasicArray_INT &    Dr,
const BasicArray_INT &    Dk,
const BasicArray_INT &    fixCP
 

Moves a point with some constraint.

This will modify the NURBS curve by respecting a certain number of constraints. u_r specifies the parameters on which the constraints should be applied. The constraint are defined by D_r^{(k)} which requires 3 vectors to fully qualify. D specifies the value of the constraint and D_r and D_k are used to specify on which parameter the constraint is applied and of what degree.

A second constraint fixCP consists of specifying which control points can not be moved by the routine.

For example, if you want to move the point C(0.5) by (10,0,10) and fix the point C(0.6) on the current curve (a move of (0,0,0)) but change its 1st derivative by (0,20,0). Doing this without modifying control point 4 . Then the following values must be inputed to the routine. u_r = [0.5, 0.6], D = [(10,0,10), (0,0,0), (0,20,0)], D_r = [0, 1, 1], D_k = [0, 0, 1] and fixCP= 4.

The values in D should be ordered in respect with r and k. i.e. for D_i=D_{r_i}^{(k_i)}, then i < j implies that r_i < r_j and that either r_i < r_j or k_i < k_j.

See section 11.5.1 of the NURBS book for an explanation of the algorithm.

Parameters:
ur  the vector of parameters on which a constraint is applied
D  a vector of the value of D_r^{(k)}
Dr  a vector specifying the value of r for D
Dk  a vector specifying the value of k for D
fixCP  a vector specifying which control points {can not} be modified.
Returns:
1 if the operation is possible, 0 if the problem is ill defined i.e. there isn't enough information to find a unique solution (the system is overdetermined) or that the system has non-independant components.
Warning:
The values of ur should not repeat. D,Dr and Dk must be of the same size.
Author:
Philippe Lavoie
Date:
24 January 1997

template<class T, int N>
int PLib::NurbsCurve< T, N >::movePoint const BasicArray< T > &    ur,
const BasicArray< Point_nD< T, N > > &    D,
const BasicArray_INT &    Dr,
const BasicArray_INT &    Dk
 

Moves a point with some constraint.

This will modify the NURBS curve by respecting a certain number of constraints. u_r specifies the parameters on which the constraints should be applied. The constraint are defined by $D_r^{(k)}$ which requires 3 vectors to fully qualify. D specifies the value of the constraint and D_r and D_k are used to specify on which parameter the constraint is applied and of what degree.

For example, if you want to move the point C(0.5) by (10,0,10) and fix the point C(0.6) on the current curve (a move of (0,0,0)) but change its 1st derivative by (0,20,0). Then the following values must be inputed to the routine. u_r = [0.5,0.6], D = [(10,0,10), (0,0,0), (0,20,0)], D_r = [0, 1, 1] and D_k = [0, 0, 1].

The values in D should be ordered in respect with r and k. {i.e.} for D_i=D_{r_i}^{(k_i)}, then i < j implies that r_i < r_j and that either r_i < r_j or k_i < k_j.

See section 11.5.1 of the NURBS book for an explanation of the algorithm.

Parameters:
ur  the vector of parameters on which a constraint is applied
D  a vector of the value of D_r^{(k)}
Dr  a vector specifying the value of r for D
Dk  a vector specifying the value of k for D
Returns:
1 if the operation is possible, 0 if the problem is ill defined i.e. there isn't enough information to find a unique solution (the system is overdetermined) or that the system has non-independant components.
Warning:
The values inside ur should not repeat. D,Dr and Dk must be of the same size.
Author:
Philippe Lavoie
Date:
24 January 1997

template<class T, int N>
int PLib::NurbsCurve< T, N >::movePoint const BasicArray< T > &    ur,
const BasicArray< Point_nD< T, N > > &    D
 

Moves a point with some constraint.

This will modify the NURBS curve by respecting a certain number of constraints. $u_r$ specifies the parameters on which the constraints should be applied. The constraint are defined by $D$ which specifies the vector by which the points should move.

For example, if you want to move the point C(0.5) by (10,0,10) and fix the point C(0.6) on the current curve (a move of (0,0,0)). u_r = 0.5, 0.6 and D = (10,0,10), (0,0,0)

The u_r vector should be in an increasing order.

See section 11.5.1 of the NURBS book for an explanation of the algorithm.

Parameters:
ur  the vector of parameters on which a constraint is applied
D  a vector of the value of D_r^{(0)}
Returns:
1 if the operation is possible, 0 if the problem is ill defined i.e. there isn't enough information to find a unique solution (the system is overdetermined) or that the system has non-independant components.
Warning:
ur and D must have the same size and the values inside ur should not repeat and they should be ordered
Author:
Philippe Lavoie
Date:
24 January 1997

template<class T, int N>
int PLib::NurbsCurve< T, N >::movePoint   u,
const BasicArray< Point_nD< T, N > > &    delta
 

Moves a point in the NURBS curve.

This modifies the curve such that the point C(u) is moved by delta. Delta is a vector containing the movement as D^{(k)} where (k) specifies the derivative. Thus at D[0], this specifies the 0th derivative movement, at D[1] it specifies the 1st derivative movement of the point. i.e. Suppose that C(u) = (10,20,3) then a D[0] = (10,10,10) will move the point to C(u) = (20,30,13)

See section 11.5.1 of the NURBS book for an explanation of the algorithm.

Parameters:
u  the point to move
delta  the vector of movement
Author:
Philippe Lavoie
Date:
24 January 1997

template<class T, int N>
int PLib::NurbsCurve< T, N >::movePoint   u,
const Point_nD< T, N > &    delta
 

Moves a point on the NURBS curve.

This modifies the curve such that the point $C(u)$ is moved by delta.

See section 11.5.1 of the NURBS book for an explanation of the algorithm.

Parameters:
u  the point to move
delta  the movement in 3D of the point at $C(u)$
Author:
Philippe Lavoie
Date:
24 January 1997

template<class T, int N>
Point_nD< T, N > PLib::NurbsCurve< T, N >::normal   u,
const Point_nD< T, N > &    v
const
 

Computes the normal of the curve at u from a vector.

Computes the normal of the curve at u from a vector. If the curve lies only in the xy-plane, then calling the function with the vector v = (0,0,1) (the z$axis) will yield a proper normal for this curve.

Parameters:
u  the parameter at which the normal is computed
v  the vector to compute the normal with
Returns:
the normal vector in 3D.
Warning:
u must be in a valid range.
Author:
Philippe Lavoie
Date:
2 September, 1997

template<class T, int N>
HPoint_nD< T, N > PLib::NurbsCurve< T, N >::operator()   u const [virtual]
 

Evaluates the curve in 4D at parameter u.

For more details on the algorithm, see A4.1 on page 124 of the Nurbs Book.

Parameters:
u  the parametric value at which the curve is evaluated
Returns:
the 4D point at C(u)
Warning:
the parametric value must be in a valid range
Author:
Philippe Lavoie
Date:
24 January, 1997

Reimplemented from PLib::ParaCurve.

template<class T, int N>
NurbsCurve< T, N > & PLib::NurbsCurve< T, N >::operator= const NurbsCurve< T, N > &    a [virtual]
 

The assignment operator for a NURBS curve.

Parameters:
curve  the NURBS curve to copy
Returns:
A reference to itself
Warning:
The curve being copied must be valid, otherwise strange results might occur.
Author:
Philippe Lavoie
Date:
24 January 1997

Reimplemented in PLib::NurbsCurveSP.

template<class T, int N>
void PLib::NurbsCurve< T, N >::projectTo const Point_nD< T, N > &    p,
  guess,
T &    u,
Point_nD< T, N > &    r,
  e1 = 0.001,
  e2 = 0.001,
int    maxTry = 100
const
 

project a point onto the curve.

It finds the closest point in the curve to a point $p$. For more information, see A6.4 and A6.5 on p 231 of the NURBS book

Parameters:
p  the point p being projected
guess  initial guess
u  the parametric value for which C(u) is closest to p.
r  the point at C(u)
e1  the minimal error
e2  the maximal error
maxTry  the maximum number of times to try before returning from the function
Author:
Philippe Lavoie
Date:
24 January 1997

template<class T, int N>
int PLib::NurbsCurve< T, N >::read ifstream &    fin [virtual]
 

reads a NurbsCurve<T,N> from a file.

Parameters:
fin  an input file stream
Returns:
0 if an error occurs, 1 otherwise
Author:
Philippe Lavoie
Date:
24 January 1997

Reimplemented in PLib::NurbsCurveSP, PLib::NurbsCurveGL, and PLib::NurbsCurveSP< float, 3 >.

template<class T, int N>
int PLib::NurbsCurve< T, N >::read const char *    filename
 

Reads a NurbsCurve<T,N> from a file.

Parameters:
filename  the filename to read the curve from
Returns:
0 if an error occurs, 1 otherwise
Author:
Philippe Lavoie
Date:
24 January 1997

template<class T, int N>
void PLib::NurbsCurve< T, N >::refineKnotVector const Vector< T > &    X
 

Refine the curve knot vector.

For more information, see A5.4 on page 164 of the NURBS book

Parameters:
X  the new knots to insert in the knot vector
Author:
Philippe Lavoie
Date:
24 January, 1997

Reimplemented in PLib::NurbsCurveSP.

template<class T, int N>
void PLib::NurbsCurve< T, N >::refineKnotVectorClosed const Vector< T > &    X
 

refine the closed curve knot vector.

Adapted from algorithm A5.4 on page 164 of the NURBS book.

Parameters:
X  the knot vector to refine with
Author:
Alejandro F Frangi
Date:
17 July, 1998

template<class T, int N>
void PLib::NurbsCurve< T, N >::removeKnot int    r,
int    s,
int    num
 

Removes an internal knot from a curve. This is A5.8 on p185 from the NURB book modified to not check for tolerance before removing the knot.

Parameters:
r  the knot to remove
s  the multiplicity of the knot
num  the number of times to try to remove the knot
Warning:
r must be an internal knot.
Author:
Philippe Lavoie
Date:
24 January 1997

Reimplemented in PLib::NurbsCurveSP, and PLib::NurbsCurveSP< float, 3 >.

template<class T, int N>
void PLib::NurbsCurve< T, N >::removeKnotsBound const Vector< T > &    ub,
Vector< T > &    ek,
  E
 

Remove knots from a curve without exceeding an error bound.

For more information about the algorithm, see A9.9 on p429 of the NURB book.

Parameters:
ub  the knot coefficients
ek  the error after removing
Author:
Philippe Lavoie
Date:
24 January 1997

Reimplemented in PLib::NurbsCurveSP.

template<class T, int N>
void PLib::NurbsCurve< T, N >::reset const Vector< HPoint_nD< T, N > > &    P1,
const Vector< T > &    U1,
int    degree
[virtual]
 

Resets a NURBS curve to new values.

Parameters:
P1  the new values for the control points
U1  the new values for the knot vector
Degree  the new degree of the curve
Warning:
The size of P1,U1 and Degree must agree: P.n()+degree+1=U.n()
Author:
Philippe Lavoie
Date:
24 January 1997

Reimplemented in PLib::NurbsCurveSP.

template<class T, int N>
void PLib::NurbsCurve< T, N >::resize int    n,
int    Deg
 

Resizes a NURBS curve.

Resizes a NURBS curve. The old values are lost and new ones have to be created.

Parameters:
n  the new number of control points for the curve
Deg  the new degree for the curve
Author:
Philippe Lavoie
Date:
24 January 1997

template<class T, int N>
int PLib::NurbsCurve< T, N >::setTangent   u,
const Point_nD< T, N > &    T0
 

set the tangent at a point.

Parameters:
u  the parametric value of the point
T  the tangent value to set to
Warning:
the tangent must be a unit length vector or an odd behavior might occur
Author:
Philippe Lavoie
Date:
2 March 1999

template<class T, int N>
int PLib::NurbsCurve< T, N >::setTangentAtEnd const Point_nD< T, N > &    T0,
const Point_nD< T, N > &    T1
 

set the tangent at the end points.

Parameters:
T0  the tangent at the begining of the curve
T1  the tangent at the end of the curve
Warning:
the tangent must be a unit length vector or and odd behavior might occur
Author:
Philippe Lavoie
Date:
2 March 1999

template<class T, int N>
int PLib::NurbsCurve< T, N >::splitAt   u,
NurbsCurve< T, N > &    cl,
NurbsCurve< T, N > &    cu
const
 

Splits the curve into two curves.

Parameters:
u  splits at this parametric value
cl  the lower curve
cu  the upper curve
Returns:
1 if the operation, 0 otherwise.
Warning:
You must make sure that you split at a valid parametric value. You can't split a curve at its end points.
Author:
Philippe Lavoie
Date:
16 October 1997

template<class T, int N>
BasicList< Point_nD< T, N > > PLib::NurbsCurve< T, N >::tesselate   tolerance,
BasicList< T > *    uk
const
 

Generates a list of points from the curve.

Generates a list of points from the curve. The list is generated within a user specified tolerance.

Parameters:
tolerance  --> the tolerance for the tesselation.
Returns:
The list of points.
Author:
Philippe Lavoie
Date:
17 October 1997

template<class T, int N>
void PLib::NurbsCurve< T, N >::transform const MatrixRT< T > &    A
 

Performs geometrical modifications.

Each control points will be modified by a rotation-translation matrix.

Parameters:
A  the rotation-translation matrix
Author:
Philippe Lavoie
Date:
22 August 1997

template<class T, int N>
void PLib::NurbsCurve< T, N >::unclamp  
 

unclamp a NURBS curve.

An unclamped NURBS curve doesn't have any equal knots at both ends of the knot vector.

Author:
Philippe Lavoie
Date:
27 April 1999

template<class T, int N>
int PLib::NurbsCurve< T, N >::write ofstream &    fout const
 

Writes a NurbsCurve<T,N> to an output stream.

Parameters:
fout  the output stream
Returns:
0 if an error occurs, 1 otherwise
Author:
Philippe Lavoie
Date:
24 January 1997

Reimplemented in PLib::NurbsCurveGL.

template<class T, int N>
int PLib::NurbsCurve< T, N >::write const char *    filename const
 

Writes a NurbsCurve<T,N> to a file.

Parameters:
filename  the filename to write to.
Returns:
0 if an error occurs, 1 otherwise
Author:
Philippe Lavoie
Date:
24 January 1997

template<class T, int N>
int PLib::NurbsCurve< T, N >::writeDisplayLINE const char *    filename,
const Color   color,
int    iNu,
  u_s,
  u_e
const
 

Writes the curve to a file in Display format as a line object.

This function writes a surface in LINE ascii format to interface with Display (Copyright 1993,1994,1995 David MacDonald, McConnell Brain Imaging Centre), Montreal Neurological Institute, McGill University.

Parameters:
filename  the name of the file to save to
color  the color of the curve
Nu  the number of points along the path
u_s  the starting parametric value for u
u_e  the end parametric value for u
Returns:
returns 1 on success, 0 otherwise.
Author:
Alejandro Frangi
Date:
15 Juanuary 1999

template<class T, int N>
int PLib::NurbsCurve< T, N >::writeDisplayLINE const char *    filename,
int    iNu,
const Color   color = blueColor,
  fA = 1
const
 

Writes the curve to in Display format as a LINE object.

This function writes a surface in LINE ascii format to interface with Display (Copyright 1993,1994,1995 David MacDonald, McConnell Brain Imaging Centre), Montreal Neurological Institute, McGill University.

Parameters:
filename  the name of the file to save to
iNu  the number of straight segments in the line
color  the color of the line
fA  alpha blending parameter of the line
Returns:
returns 1 on success, 0 otherwise.
Author:
Alejandro Frangi
Date:
21 January 1999

template<class T, int N>
int PLib::NurbsCurve< T, N >::writePS const char *    filename,
int    cp = 0,
  magFact = T(-1),
  dash = T(5),
bool    bOpen = true
const
 

Writes the curve in the postscript format to a file.

Parameters:
filename  the file to write the postscript file to
cp  a flag indicating if the control points should be drawn, 0 = no and 1 = yes
magFact  a magnification factor, the 2D point of the control points will be magnified by this value. The size is measured in postscript points. If the magFact is set to a value smaller or equal to 0, than the program will try to guess a magnification factor such that the curve is large enough to fill the page.
dash  the size of the dash in postscript points . A size smaller or equal to 0 indicates that the line joining the control points is plain.
Returns:
0 if an error occurs, 1 otherwise
Warning:
If the weights of the curve are not all at 1, the result might not be representative of the true NURBS curve.
Author:
Philippe Lavoie
Date:
16 February 1997

template<class T, int N>
int PLib::NurbsCurve< T, N >::writePSp const char *    filename,
const Vector< Point_nD< T, N > > &    points,
const Vector< Point_nD< T, N > > &    vectors,
int    cp = 0,
  magFact = 0.0,
  dash = 5.0,
bool    bOpen = true
const
 

Writes a post-script file representing the curve.

Writes the curve in the postscript format to a file, it also draws the points defined in $points$ with their associated vectors if $vector$ is used.

Parameters:
filename  the file to write the postscript file to
points  draws these additional points as empty circles
vectors  specify a vector associated with the points (this can be an empty vector)
cp  a flag indicating if the control points should be drawn, 0 = no and 1 = yes
magFact  a magnification factor, the 2D point of the control points will be magnified by this value. The size is measured in postscript points. If the magFact is set to a value smaller or equal to 0, than the program will try to guess a magnification factor such that the curve is large enough to fill the page.
dash  the size of the dash in postscript points . A size smaller or equal to 0 indicates that the line joining the control points is plain.
Returns:
0 if an error occurs, 1 otherwise
Warning:
If the weights of the curve are not all at 1, the result might not be representative of the true NURBS curve. If vectors is used, then it must be of the same size as points. If a vector element is (0,0,0) it will not be drawn.
Author:
Philippe Lavoie
Date:
16 February 1997

template<class T, int N>
int PLib::NurbsCurve< T, N >::writeVRML const char *    filename,
  radius,
int    K,
const Color   color,
int    Nu,
int    Nv,
  u_s,
  u_e
const
 

Writes the curve to a VRML file.

A circle is swept around the trajectory made by the curve. The resulting surface is saved as a VRML file.

Parameters:
filename  the name of the VRML file to save to
radius  the radius of the line
K  the minimum number of interpolation
color  the color of the line
Nu  the number of points for the circle
Nv  the number of points along the path
u_s  the starting parametric value for u
u_e  the end parametric value for u
Returns:
returns 1 on success, 0 otherwise.
Author:
Philippe Lavoie
Date:
25 July 1997

template<class T, int N>
int PLib::NurbsCurve< T, N >::writeVRML ostream &    fout,
  radius,
int    K,
const Color   color,
int    Nu,
int    Nv,
  u_s,
  u_e
const
 

Writes the curve to a VRML file.

A circle is swept around the trajectory made by the curve. The resulting surface is saved as a VRML file.

Parameters:
filename  the name of the ostream to write to
radius  the radius of the line
K  the minimum number of interpolation
color  the color of the line
Nu  the number of points for the circle
Nv  the number of points along the path
u_s  the starting parametric value for u
u_e  the end parametric value for u
Returns:
returns 1 on success, 0 otherwise.
Author:
Philippe Lavoie
Date:
25 July 1997

template<class T, int N>
int PLib::NurbsCurve< T, N >::writeVRML97 ostream &    fout,
  radius,
int    K,
const Color   color,
int    Nu,
int    Nv,
  u_s,
  u_e
const
 

Writes the curve to a VRML file.

A circle is swept around the trajectory made by the curve. The resulting surface is saved as a VRML file.

Parameters:
filename  the name of the ostream to write to
radius  the radius of the line
K  the minimum number of interpolation
color  the color of the line
Nu  the number of points for the circle
Nv  the number of points along the path
u_s  the starting parametric value for u
u_e  the end parametric value for u
Returns:
returns 1 on success, 0 otherwise.
Author:
Philippe Lavoie
Date:
4 May 1999

template<class T, int N>
int PLib::NurbsCurve< T, N >::writeVRML97 const char *    filename,
  radius,
int    K,
const Color   color,
int    Nu,
int    Nv,
  u_s,
  u_e
const
 

Writes the curve to a VRML file.

A circle is swept around the trajectory made by the curve. The resulting surface is saved as a VRML file.

Parameters:
filename  the name of the VRML file to save to
radius  the radius of the line
K  the minimum number of interpolation
color  the color of the line
Nu  the number of points for the circle
Nv  the number of points along the path
u_s  the starting parametric value for u
u_e  the end parametric value for u
Returns:
returns 1 on success, 0 otherwise.
Author:
Philippe Lavoie
Date:
4 May 1999


Friends And Related Function Documentation

template<class T, int N>
T chordLengthParam const Vector< Point_nD< T, N > > &    Q,
Vector< T > &    ub
[related]
 

chord length parameterization.

Performs chord length parameterization from a vector of points. There are more details about this function in the LaTeX version.

Parameters:
Q  a vector of 3D points
ub  the result of chord length parameterization
Returns:
the total chord length of the points.
Author:
Philippe Lavoie
Date:
24 January, 1997

template<class T, int N>
void generateCompatibleCurves NurbsCurveArray< T, N > &    ca [friend]
 

Generate compatible curves from an array of curves \relates NurbsCurveArray.

This routine will put to the same degree all the curves in the array and it will ensure that they have the same knot vector.

Parameters:
ca  the array of Nurbs curves
Warning:
the knot vector of all the curves must be in the range [0,1]
Author:
Philippe Lavoie
Date:
24 January, 1997


The documentation for this class was generated from the following files:
Generated on Sat Jan 26 15:25:28 2002 for NURBS++ by doxygen1.2.12 written by Dimitri van Heesch, © 1997-2001