Bezier Curves and Surfaces
by Jesper Tveit

Introduction

In this tutorial I will attempt to give you an understanding of a common way to make curved lines and surfaces. I will also give you some coding examples that will show you that it actually works.

I will assume that you know how to use c++ in the coding examples. The coding examples also contain some OpenGL function calls, but they are commented and not too hard to understand.

I will also assume that you know some basic math, but I'm starting off pretty easy by introducing parametric curves so you should be able to grasp it even if you don't have the greatest math knowledge.

I have tested the coding examples, and I have looked over the text several times, so they should be fairly bug free. Still, if you find anything that looks like its wrong, let me know. You can find my email address at the bottom of the page.

Parametric curves

You are probably familiar with curves like y = x² or f(x)=x². A parametric curve is similar to those, but the x- and y-values are calculated in separate functions. To do this we need another variable (I will call it a). So instead of f(x), we say X(a) and Y(a), where X(a) gives you an x-value for each value of a, and Y(a) gives you an corresponding y-value for the same value of a. Here is an example:

X(a) = a/2
Y(a) = a+5

If we substitute a with a number we can get a point that lies on the curve/line:

a = 6

X(6) = 6/2 = 3 = x
Y(6) = 6+5 = 11 = y

We now know that (3,11) is a point on the curve.

We can easily make a curve in 3D just by defining Z(a).

Bezier curves

Simple line

The simplest form is a straight line from a control point A, to a control point B.

I will use Ax to denote the x-coordinate of control point A, and Ay to denote the y-coordinate. The same goes for B.

I am introducing another variable, b, just to make it a bit more simple to read the functions. However, the variable b is only a in disguise, b is always equal to 1-a. So if a = 0.2, then b = 1-0.2 = 0.8. So when you see the variable b, think: (1-a).

This parametric curve describes a line that goes from point A, to point B when the variable a goes from 1.0 to 0.0:

X(a) = Ax·a + Bx·b
Y(a) = Ay·a + By·b
Z(a) = Az·a + Bz·b

If you set a = 0.5 (and thus also b = 0.5, since b = 1.0-a = 1.0-0.5 = 0.5) then X(a), Y(a) and Z(a) will give you the 3D coordinate of the point on the middle of the line from point A to point B.

Note that the reason this works is because a + b is always equal to one, and when you multiply something with one, it will stay unchanged. This makes the curve behave predictably and lets you place the control points anywhere in the coordinate system, since when a = 1.0 then b = 0.0 thus making the point equal to one of the control points, and completely ignoring the other one.

Quadratic

Ok, now we have done lines, but what about curves?

This: (a+b), where b = (1.0-a) is the key to everything. We know that a polynomial function is a curve, so why not try to make (a+b) a polynomial function? Let's try:

(a+b)² = a² + 2a·b + b²

Knowing that b=1-a we see that a² + 2·a·b + b² is still equal to one, since (a+b) = 1, and 1² = 1.

We now need three control points, A, B and C, where A and C are the end points of the curve, and B decides how much and in which direction it curves. Except for that, it is the same as with the parametric line.

X(a) = Ax·a² + Bx·2·a·b + Cx·b²
Y(a) = Ay·a² + By·2·a·b + Cy·b²
Z(a) = Az·a² + Bz·2·a·b + Cz·b²

Still, if you set a = 0.5, then a² = 0.25, 2·a·b = 0.5 and b² = 0.25 giving the middle point, B, the biggest impact and making point A and C pull equally to their sides. If you set a = 1.0, then a² = 1.0, a·b = 0.0 and b² = 0.0 meaning that result is the control point A itself, just the same as setting a = 0.0 will return the coordinates of control point C.

Cubic

We can also increase the number of control points using (a+b)³ instead of (a+b)², giving you more control over how to bend the curve.

(a+b)³ = a³ + 3·a²·b + 3·a·b² + b³

Now we need four control points A, B, C and D, where A and D are the end points.

X(a) = Ax·a³ + Bx·3·a²·b + Cx·3·a·b² + Dx·b³
Y(a) = Ay·a³ + By·3·a²·b + Cy·3·a·b² + Dy·b³
Z(a) = Az·a³ + Bz·3·a²·b + Cz·3·a·b² + Dz·b³

It still works the same way as the previous ones, as a goes from 1.0 to 0.0 (and thus b from 0.0 to 1.0) the functions will return coordinates on a smooth line from control point A to control point D, curving towards B and C on the way.

You can easily use (a+b) to the power of n and get n+1 control points (where n is any integer equal to, or greater than one). But remember that the more control points you have, the more work it is, both for you and for the computer. This might not seem like much when only talking about lines, but when it comes to surfaces it definitely has something to say.

Personally I like the cubic ones best. You can easily make circles with them, and you can control the direction of the curve independently at each control point.

Example of how to draw a cubic curve using OpenGL and C++:

this is how it will look
(excluding the red lines and the letters)

// Control points (substitute these values with your own if you like)
double Ax = -2.0; double Ay = -1.0; double Az = 1.0;
double Bx = -1.0; double By = 3.0; double Bz = 1.0;
double Cx = 1.0; double Cy = -3.0; double Cz = -1.0;
double Dx = 2.0; double Dy = 1.0; double Dz = -1.0;

// Points on the curve
double X;
double Y;
double Z;

// Variable
double a = 1.0;
double b = 1.0 - a;

// Tell OGL to start drawing a line strip
glBegin(GL_LINE_STRIP);

/* We will not actually draw a curve, but we will divide the curve into small
points and draw a line between each point. If the points are close enough, it
will appear as a curved line. 20 points are plenty, and since the variable goes
from 1.0 to 0.0 we must change it by 1/20 = 0.05 each time */

for(int i = 0; i <= 20; i++)
{
  // Get a point on the curve
  X = Ax*a*a*a + Bx*3*a*a*b + Cx*3*a*b*b + Dx*b*b*b;
  Y = Ay*a*a*a + By*3*a*a*b + Cy*3*a*b*b + Dy*b*b*b;
  Z = Az*a*a*a + Bz*3*a*a*b + Cz*3*a*b*b + Dz*b*b*b;
  
  // Draw the line from point to point (assuming OGL is set up properly)
  glVertex3d(X, Y, Z);

  // Change the variable
  a -= 0.05;
  b = 1.0 - a;
}

// Tell OGL to stop drawing the line strip
glEnd();

/* Normally you will want to save the coordinates to an array for later use. And
you will probably not need to calculate the curve each frame. This code
demonstrates an easily understandable way to do it, not necessarily the most
useful way. */

Curved Surfaces

If you think of a line as 1D, and a surface as 2D, and understand that a square can be defined as length·width, then you are not far from understanding that a curved bezier surface can be defined as a bezierCurve·bezierCurve. This can be done with any (a+b) to the power of n -type curve (with n greater than or equal to one). I will demonstrate how to do this with the cubic one (n = 3).

This defines our first cubic curve (leaving out the control points):

(a+b)³ = a³ + 3·a²·b + 3·a·b² + b³

But now we need two, so here is the second one:

(c+d)³ = c³ + 3·c²·d + 3·c·d² + c³

In the last one, d=1-c, just the same way that b=1-a in the first one.

The next step is to multiply them by each other:

(a³ + 3·a²·b + 3·a·b² + b³)·(c³ + 3·c²·d + 3·c·d² + d³) =

     a³·c³     + 3·a³·c²·d   + 3·a³·c·d²   + a³·d³
   + 3·a²·b·c³ + 9·a²·b·c²·d + 9·a²·b·c·d² + 3·a²·b·d³
   + 3·a·b²·c³ + 9·a·b²·c²·d + 9·a·b²·c·d² + 3·a·b²·d³
   + b³·c³     + 3·b³·c²·d   + 3·b³·c·d²   + b³·d³

This is still equal to one, since (a+b)³ = 1, and (c+d)³ = 1, and 1·1 = 1.

The cubic curves had four control points, so it makes sense that we get sixteen control points when we multiply two of them with each other. I'm calling the control points A, B, C, D, E, F, G, H, I, J, K, L, M, N, O and P. Here is how the functions will look when we multiply in the control points:

X(a,c) = Ax·a³·c³       + Bx·3·a³·c²·d
       + Cx·3·a³·c·d²   + Dx·a³·d³
       + Ex·3·a²·b·c³   + Fx·9·a²·b·c²·d
       + Gx·9·a²·b·c·d² + Hx·3·a²·b·d³
       + Ix·3·a·b²·c³   + Jx·9·a·b²·c²·d
       + Kx·9·a·b²·c·d² + Lx·3·a·b²·d³
       + Mx·b³·c³       + Nx·3·b³·c²·d
       + Ox·3·b³·c·d²   + Px·b³·d³

Y(a,c) = Ay·a³·c³       + By·3·a³·c²·d
       + Cy·3·a³·c·d²   + Dy·a³·d³
       + Ey·3·a²·b·c³   + Fy·9·a²·b·c²·d
       + Gy·9·a²·b·c·d² + Hy·3·a²·b·d³
       + Iy·3·a·b²·c³   + Jy·9·a·b²·c²·d
       + Ky·9·a·b²·c·d² + Ly·3·a·b²·d³
       + My·b³·c³       + Ny·3·b³·c²·d
       + Oy·3·b³·c·d²   + Py·b³·d³

Z(a,c) = Az·a³·c³       + Bz·3·a³·c²·d
       + Cz·3·a³·c·d²   + Dz·a³·d³
       + Ez·3·a²·b·c³   + Fz·9·a²·b·c²·d
       + Gz·9·a²·b·c·d² + Hz·3·a²·b·d³
       + Iz·3·a·b²·c³   + Jz·9·a·b²·c²·d
       + Kz·9·a·b²·c·d² + Lz·3·a·b²·d³
       + Mz·b³·c³       + Nz·3·b³·c²·d
       + Oz·3·b³·c·d²   + Pz·b³·d³

The functions now have two variables: a and c. Passing values of a and c between 0 and 1 into the functions will give you a point on the surface. Remember that the values of a and c vary independently.

Now you can draw a curved surface, but to make it look good it needs proper lighting. To do that we need normal vectors for each point on the surface that we want to draw. To find the normal vectors we will go through four steps. First we need to find the derivative with respect to a (remember that b has to be treated as (1-a), and that c and d remain untouched):

XDA(a,c) = Ax·3·a²·c³             + Bx·9·a²·c²·d
         + Cx·9·a²·c·d²           + Dx·3·a²·d³
         + Ex·3·(2·a-3·a²)·c³     + Fx·9·(2·a-3·a²)·c²·d
         + Gx·9·(2·a-3·a²)·c·d²   + Hx·3·(2·a-3·a²)·d³
         + Ix·3·(1-4·a+3·a²)·c³   + Jx·9·(1-4·a+3·a²)·c²·d
         + Kx·9·(1-4·a+3·a²)·c·d² + Lx·3·(1-4·a+3·a²)·d³
         + Mx·3·(2·a-1-a²)·c³     + Nx·9·(2·a-1-a²)·c²·d
         + Ox·9·(2·a-1-a²)·c·d²   + Px·3·(2·a-1-a²)·d³

YDA(a,c) = Ay·3·a²·c³             + By·9·a²·c²·d
         + Cy·9·a²·c·d²           + Dy·3·a²·d³
         + Ey·3·(2·a-3·a²)·c³     + Fy·9·(2·a-3·a²)·c²·d
         + Gy·9·(2·a-3·a²)·c·d²   + Hy·3·(2·a-3·a²)·d³
         + Iy·3·(1-4·a+3·a²)·c³   + Jy·9·(1-4·a+3·a²)·c²·d
         + Ky·9·(1-4·a+3·a²)·c·d² + Ly·3·(1-4·a+3·a²)·d³
         + My·3·(2·a-1-a²)·c³     + Ny·9·(2·a-1-a²)·c²·d
         + Oy·9·(2·a-1-a²)·c·d²   + Py·3·(2·a-1-a²)·d³

ZDA(a,c) = Az·3·a²·c³             + Bz·9·a²·c²·d
         + Cz·9·a²·c·d²           + Dz·3·a²·d³
         + Ez·3·(2·a-3·a²)·c³     + Fz·9·(2·a-3·a²)·c²·d
         + Gz·9·(2·a-3·a²)·c·d²   + Hz·3·(2·a-3·a²)·d³
         + Iz·3·(1-4·a+3·a²)·c³   + Jz·9·(1-4·a+3·a²)·c²·d
         + Kz·9·(1-4·a+3·a²)·c·d² + Lz·3·(1-4·a+3·a²)·d³
         + Mz·3·(2·a-1-a²)·c³     + Nz·9·(2·a-1-a²)·c²·d
         + Oz·9·(2·a-1-a²)·c·d²   + Pz·3·(2·a-1-a²)·d³

If you are unfamiliar with derivation and think this seems very confusing then don't worry, this looks pretty confusing to anyone. The nice patterns we had have disappeared. Derivation is a complex theme so I will not explain it in depth here. Briefly though, the reason why we do this is that the derivative functions allow us to find tangent vectors on the surface, so later on we can use them to calculate the normal vectors.

The second step is to find the derivative with respect to c (this is easy, since we already have the recipe from the last one):

XDC(a,c) = Ax·3·a³·c²             + Bx·3·a³·(2·c-3·c²)
         + Cx·3·a³·(1-4·c+3·c²)   + Dx·3·a³·(-1+2·c-c²)
         + Ex·9·a²·b·c²           + Fx·9·a²·b·(2·c-3·c²)
         + Gx·9·a²·b·(1-4·c+3·c²) + Hx·9·a²·b·(-1+2·c-c²)
         + Ix·9·a·b²·c²           + Jx·9·a·b²·(2·c-3·c²)
         + Kx·9·a·b²·(1-4·c+3·c²) + Lx·9·a·b²·(-1+2·c-c²)
         + Mx·3·b³·c²             + Nx·3·b³·(2·c-3·c²)
         + Ox·3·b³·(1-4·c+3·c²)   + Px·3·b³·(-1+2·c-c²)

YDC(a,c) = Ay·3·a³·c²             + By·3·a³·(2·c-3·c²)
         + Cy·3·a³·(1-4·c+3·c²)   + Dy·3·a³·(-1+2·c-c²)
         + Ey·9·a²·b·c²           + Fy·9·a²·b·(2·c-3·c²)
         + Gy·9·a²·b·(1-4·c+3·c²) + Hy·9·a²·b·(-1+2·c-c²)
         + Iy·9·a·b²·c²           + Jy·9·a·b²·(2·c-3·c²)
         + Ky·9·a·b²·(1-4·c+3·c²) + Ly·9·a·b²·(-1+2·c-c²)
         + My·3·b³·c²             + Ny·3·b³·(2·c-3·c²)
         + Oy·3·b³·(1-4·c+3·c²)   + Py·3·b³·(-1+2·c-c²)

ZDC(a,c) = Az·3·a³·c²             + Bz·3·a³·(2·c-3·c²)
         + Cz·3·a³·(1-4·c+3·c²)   + Dz·3·a³·(-1+2·c-c²)
         + Ez·9·a²·b·c²           + Fz·9·a²·b·(2·c-3·c²)
         + Gz·9·a²·b·(1-4·c+3·c²) + Hz·9·a²·b·(-1+2·c-c²)
         + Iz·9·a·b²·c²           + Jz·9·a·b²·(2·c-3·c²)
         + Kz·9·a·b²·(1-4·c+3·c²) + Lz·9·a·b²·(-1+2·c-c²)
         + Mz·3·b³·c²             + Nz·3·b³·(2·c-3·c²)
         + Oz·3·b³·(1-4·c+3·c²)   + Pz·3·b³·(-1+2·c-c²)

Now we have a way to describe two tangent vectors on each point on the surface. This mean that we will be able to calculate the normal vector of any point on the surface.

In general, if you have two vectors, V={vx, vy, vz} and U={ux, uy, uz} then those two vectors can be 'crossed' to produce a new vector N={nx, ny, nz} that is orthogonal to both V and U (Two vectors beeing orthogonal means that the angle between them are 90°).

This is how you cross two vectors V and U and put the result in N:

nx =   (vy·uz)-(uy·vz)
ny = -((vx·uz)-(ux·vz))
nz =   (vx·uy)-(ux·vy)

This is normally refered to as V×U=N

You might have guessed what the third step is. Since we have two tangent vectors on the surface (the derived functions: {XDC, YDC, ZDC}, denoted by DC, and {XDA, YDA, ZDA} denoted by DA) we can obtain a normal vector (N={Nx, Ny, Nz}) like this:

N = DA × DC

Which is equivalent to:

Nx =   (YDA·ZDC)-(YDC·ZDA)
Ny = -((XDA·ZDC)-(XDC·ZDA))
Nz =   (XDA·YDC)-(XDC·YDA)

The fourth step is to normalize the normal vectors (to normalize a vector means to make its length equal to one; in many cases it will be assumed that a normal vector has length equal to one).

To normalize a vector we must first find the its length (L):

L = sqrt(Nx²+Ny²+Nz²)

Then we divide it by its length to obtain a new vector (N' = {Nx',Ny',Nz'}) of length one:

Nx'= Nx/L
Ny'= Ny/L
Nz'= Nz/L

N' together with X(a,c), Y(a,c) and Z(a,c) give us all we need to make a smoothly shaded curved surface.

Example of how to draw a curved surface using OpenGL and c++:

this is how it will look
(excluding the red lines and the letters)

// Sixteen control points (substitute these values with your own if you like)
double Ax = -2.0; double Ay =  2.0; double Az =  1.0;
double Bx = -1.0; double By =  3.0; double Bz = -1.0;
double Cx =  1.0; double Cy =  3.0; double Cz =  1.0;
double Dx =  2.0; double Dy =  2.0; double Dz = -1.0;

double Ex = -1.5; double Ey =  1.0; double Ez =  1.0;
double Fx = -0.5; double Fy =  1.5; double Fz = -1.0;
double Gx =  1.5; double Gy =  1.5; double Gz =  1.0;
double Hx =  2.5; double Hy =  1.0; double Hz = -1.0;

double Ix = -2.5; double Iy = -1.0; double Iz =  1.0;
double Jx = -1.5; double Jy = -0.5; double Jz = -1.0;
double Kx =  0.5; double Ky = -0.5; double Kz =  1.0;
double Lx =  1.5; double Ly = -1.0; double Lz = -1.0;

double Mx = -2.0; double My = -2.0; double Mz =  1.0;
double Nx = -1.0; double Ny = -1.0; double Nz = -1.0;
double Ox =  1.0; double Oy = -1.0; double Oz =  1.0;
double Px =  2.0; double Py = -2.0; double Pz = -1.0;


/* We will not actually draw a curved surface, but we will divide the
surface into small quads and draw them. If the quads are small enough,
it will appear as a curved surface. We will use a variable, detail, to
define how many quads to use. Since the variables goes from 1.0 to 0.0
we must change them by 1/detail from vertex to vertex. We will also
store the vertices and the normal vectors in arrays and draw them in a
separate loop */

// Detail of 10 mean that we will calculate 11·11 vertices
int    detail = 10;
double change = 1.0 / (double)detail;

// Vertices (maximum detail will now be 20·20 quads)
double Xv[21][21];
double Yv[21][21];
double Zv[21][21];

// Normal vectors
double Xn[21][21];
double Yn[21][21];
double Zn[21][21];

// Just making sure that the detail level is not set too high
if(detail > 20)
{
  detail = 20;
}

// Variables
double a = 1.0;
double b = 1.0 - a;
double c = 1.0;
double d = 1.0 - c;

// Tangent vectors
double Xta;
double Yta;
double Zta;

double Xtc;
double Ytc;
double Ztc;

/* Since we have two variables, we need two loops, we will change the
a-variable from 1.0 to 0.0 by steps of 1/detail ( = change), and for each
step we loop the c-variable from 1.0 to 0.0, thus creating a grid of
points covering the surface. Note that we could have had separate detail
 levels for the a-variable and the c-variable if we wanted to */
for(int i = 0; i <= detail; i++)
{
  for(int j = 0; j <= detail; j++)
  {
    // First get the vertices
    Xv[i][j] = Ax*a*a*a*c*c*c   + Bx*3*a*a*a*c*c*d
             + Cx*3*a*a*a*c*d*d + Dx*a*a*a*d*d*d
             + Ex*3*a*a*b*c*c*c + Fx*9*a*a*b*c*c*d
             + Gx*9*a*a*b*c*d*d + Hx*3*a*a*b*d*d*d
             + Ix*3*a*b*b*c*c*c + Jx*9*a*b*b*c*c*d
             + Kx*9*a*b*b*c*d*d + Lx*3*a*b*b*d*d*d
             + Mx*b*b*b*c*c*c   + Nx*3*b*b*b*c*c*d
             + Ox*3*b*b*b*c*d*d + Px*b*b*b*d*d*d;

    Yv[i][j] = Ay*a*a*a*c*c*c   + By*3*a*a*a*c*c*d
             + Cy*3*a*a*a*c*d*d + Dy*a*a*a*d*d*d
             + Ey*3*a*a*b*c*c*c + Fy*9*a*a*b*c*c*d
             + Gy*9*a*a*b*c*d*d + Hy*3*a*a*b*d*d*d
             + Iy*3*a*b*b*c*c*c + Jy*9*a*b*b*c*c*d
             + Ky*9*a*b*b*c*d*d + Ly*3*a*b*b*d*d*d
             + My*b*b*b*c*c*c   + Ny*3*b*b*b*c*c*d
             + Oy*3*b*b*b*c*d*d + Py*b*b*b*d*d*d;

    Zv[i][j] = Az*a*a*a*c*c*c   + Bz*3*a*a*a*c*c*d
             + Cz*3*a*a*a*c*d*d + Dz*a*a*a*d*d*d
             + Ez*3*a*a*b*c*c*c + Fz*9*a*a*b*c*c*d
             + Gz*9*a*a*b*c*d*d + Hz*3*a*a*b*d*d*d
             + Iz*3*a*b*b*c*c*c + Jz*9*a*b*b*c*c*d
             + Kz*9*a*b*b*c*d*d + Lz*3*a*b*b*d*d*d
             + Mz*b*b*b*c*c*c   + Nz*3*b*b*b*c*c*d
             + Oz*3*b*b*b*c*d*d + Pz*b*b*b*d*d*d;

    // Then use the derived functions to get the tangent vectors
    Xta = Ax*3*a*a*c*c*c           + Bx*9*a*a*c*c*d
        + Cx*9*a*a*c*d*d           + Dx*3*a*a*d*d*d
        + Ex*3*(2*a-3*a*a)*c*c*c   + Fx*9*(2*a-3*a*a)*c*c*d
        + Gx*9*(2*a-3*a*a)*c*d*d   + Hx*3*(2*a-3*a*a)*d*d*d
        + Ix*3*(1-4*a+3*a*a)*c*c*c + Jx*9*(1-4*a+3*a*a)*c*c*d
        + Kx*9*(1-4*a+3*a*a)*c*d*d + Lx*3*(1-4*a+3*a*a)*d*d*d
        + Mx*3*(2*a-1-a*a)*c*c*c   + Nx*9*(2*a-1-a*a)*c*c*d
        + Ox*9*(2*a-1-a*a)*c*d*d   + Px*3*(2*a-1-a*a)*d*d*d;

    Yta = Ay*3*a*a*c*c*c           + By*9*a*a*c*c*d
        + Cy*9*a*a*c*d*d           + Dy*3*a*a*d*d*d
        + Ey*3*(2*a-3*a*a)*c*c*c   + Fy*9*(2*a-3*a*a)*c*c*d
        + Gy*9*(2*a-3*a*a)*c*d*d   + Hy*3*(2*a-3*a*a)*d*d*d
        + Iy*3*(1-4*a+3*a*a)*c*c*c + Jy*9*(1-4*a+3*a*a)*c*c*d
        + Ky*9*(1-4*a+3*a*a)*c*d*d + Ly*3*(1-4*a+3*a*a)*d*d*d
        + My*3*(2*a-1-a*a)*c*c*c   + Ny*9*(2*a-1-a*a)*c*c*d
        + Oy*9*(2*a-1-a*a)*c*d*d   + Py*3*(2*a-1-a*a)*d*d*d;

    Zta = Az*3*a*a*c*c*c           + Bz*9*a*a*c*c*d
        + Cz*9*a*a*c*d*d           + Dz*3*a*a*d*d*d
        + Ez*3*(2*a-3*a*a)*c*c*c   + Fz*9*(2*a-3*a*a)*c*c*d
        + Gz*9*(2*a-3*a*a)*c*d*d   + Hz*3*(2*a-3*a*a)*d*d*d
        + Iz*3*(1-4*a+3*a*a)*c*c*c + Jz*9*(1-4*a+3*a*a)*c*c*d
        + Kz*9*(1-4*a+3*a*a)*c*d*d + Lz*3*(1-4*a+3*a*a)*d*d*d
        + Mz*3*(2*a-1-a*a)*c*c*c   + Nz*9*(2*a-1-a*a)*c*c*d
        + Oz*9*(2*a-1-a*a)*c*d*d   + Pz*3*(2*a-1-a*a)*d*d*d;

    Xtc = Ax*3*a*a*a*c*c           + Bx*3*a*a*a*(2*c-3*c*c)
        + Cx*3*a*a*a*(1-4*c+3*c*c) + Dx*3*a*a*a*(-1+2*c-c*c)
        + Ex*9*a*a*b*c*c           + Fx*9*a*a*b*(2*c-3*c*c)
        + Gx*9*a*a*b*(1-4*c+3*c*c) + Hx*9*a*a*b*(-1+2*c-c*c)
        + Ix*9*a*b*b*c*c           + Jx*9*a*b*b*(2*c-3*c*c)
        + Kx*9*a*b*b*(1-4*c+3*c*c) + Lx*9*a*b*b*(-1+2*c-c*c)
        + Mx*3*b*b*b*c*c           + Nx*3*b*b*b*(2*c-3*c*c)
        + Ox*3*b*b*b*(1-4*c+3*c*c) + Px*3*b*b*b*(-1+2*c-c*c);

    Ytc = Ay*3*a*a*a*c*c           + By*3*a*a*a*(2*c-3*c*c)
        + Cy*3*a*a*a*(1-4*c+3*c*c) + Dy*3*a*a*a*(-1+2*c-c*c)
        + Ey*9*a*a*b*c*c           + Fy*9*a*a*b*(2*c-3*c*c)
        + Gy*9*a*a*b*(1-4*c+3*c*c) + Hy*9*a*a*b*(-1+2*c-c*c)
        + Iy*9*a*b*b*c*c           + Jy*9*a*b*b*(2*c-3*c*c)
        + Ky*9*a*b*b*(1-4*c+3*c*c) + Ly*9*a*b*b*(-1+2*c-c*c)
        + My*3*b*b*b*c*c           + Ny*3*b*b*b*(2*c-3*c*c)
        + Oy*3*b*b*b*(1-4*c+3*c*c) + Py*3*b*b*b*(-1+2*c-c*c);

    Ztc = Az*3*a*a*a*c*c           + Bz*3*a*a*a*(2*c-3*c*c)
        + Cz*3*a*a*a*(1-4*c+3*c*c) + Dz*3*a*a*a*(-1+2*c-c*c)
        + Ez*9*a*a*b*c*c           + Fz*9*a*a*b*(2*c-3*c*c)
        + Gz*9*a*a*b*(1-4*c+3*c*c) + Hz*9*a*a*b*(-1+2*c-c*c)
        + Iz*9*a*b*b*c*c           + Jz*9*a*b*b*(2*c-3*c*c)
        + Kz*9*a*b*b*(1-4*c+3*c*c) + Lz*9*a*b*b*(-1+2*c-c*c)
        + Mz*3*b*b*b*c*c           + Nz*3*b*b*b*(2*c-3*c*c)
        + Oz*3*b*b*b*(1-4*c+3*c*c) + Pz*3*b*b*b*(-1+2*c-c*c);

    // Cross the tangent vectors, put the result to the normal vector array
    // Note: I simplified -((Xta*Ztc)-(Xtc*Zta)) to (Xtc*Zta) - (Xta*Ztc)
    Xn[i][j] = (Yta*Ztc) - (Ytc*Zta);
    Yn[i][j] = (Xtc*Zta) - (Xta*Ztc);
    Zn[i][j] = (Xta*Ytc) - (Xtc*Yta);

    // Find length of normal vector
    double length = sqrt((Xn[i][j]*Xn[i][j])+(Yn[i][j]
                         *Yn[i][j])+(Zn[i][j]*Zn[i][j]));

    // Normalize (and prevent divide by zero error)
    if(length > 0)
    {
      length = 1.0/length;
      Xn[i][j] *= length;
      Yn[i][j] *= length;
      Zn[i][j] *= length;
    }

    //change the c-variable within the inner loop
    c -= change;
    d  = 1.0 - c;
  }
  //change the a-variable outside the inner loop
  a -= change;
  b  = 1.0 - a;

  // Reset the c-variable to make it ready for the inner loop again
  c = 1.0;
  d = 1.0 - c;
}


/* Now we have two arrays, one with vertices, and one with normal vectors,
drawing them is straightforward if you know how to use a graphics API.
Following is one way to do it using openGL and triangle strips. (assuming
GL_LIGHTING etc.. has been properly set up) */
for(int m = 0; m < detail; m++)
{
  glBegin(GL_TRIANGLE_STRIP);
  for(int n = 0; n <= detail; n++)
  {
    glNormal3d(Xn[m][n],Yn[m][n],Zn[m][n]);
    glVertex3d(Xv[m][n],Yv[m][n],Zv[m][n]);

    // Note that I used real-less-than in the first loop, since I want to
    // access the m+1 entry in the array to properly draw the triangle strip
    glNormal3d(Xn[m+1][n],Yn[m+1][n],Zn[m+1][n]);
    glVertex3d(Xv[m+1][n],Yv[m+1][n],Zv[m+1][n]);
  }
  glEnd();
}

Comment

Bezier curves do not necessarily have to be curves in 3D or 2D space, they might just as well be color curves (by substituting x,y,z coordinates with r,g,b values), curves in 4D space or anything else. Personally I find it useful to use a 2D bezier patch to calculate texture coordinates that fit on a 3D bezier surface.

Sources

Various web sites and my head ;)

Questions, comments, report mistakes etc.

My name: Jesper Tveit
E-mail : 7396@online.no
Website: http://7396.home.online.no/bw/

For background info, try a google search for "Etienne Bézier"

That's it! Thanks for your attention

Discuss this article in the forums

Date this article was posted to GameDev.net: 4/27/2002
(Note that this date does not necessarily correspond to the date the article was written)