Chapter 7
Surface Modeling

7.1  Introduction

Although most mechanical parts can be modeled simply using planes, spheres, cones, and toruses, many other geometric entities such as telephone handsets, car and airplane bodies, and sculptures may not be represented by these primitive surfaces. In CAD systems, Bézier and B-spline surfaces are used as the standard in industry for representing free-form or sculptured surfaces. For this reason, Bézier and B-spline surfaces are sometimes called free-form surfaces. In this chapter, we shall extend our knowledge on Bézier and B-spline curves to surfaces. Although it is very involved in processing a surface in terms of data arrangement and computation, the fundamental mathematics on surfaces is directly derived from curves. If you have read through the previous chapter on curve modeling and completed all examples, you will find that the underlying mathematics on parametric surfaces is fairly straightforward. Virtually, you do not need any new algorithms to evaluate Bézier and B-spline surfaces.

Bézier and B-spline surfaces are used in computer graphics and CAD systems for either iterative design of geometric shapes or representation of prescribed surfaces. The latter usually involves interpolation and approximation techniques. For example, to represent a clay model of a car by B-spline surfaces we may collect a set of digitized data from the clay model, then fit these digitized data by B-spline surfaces using interpolation or least square approximation method. In this chapter, we cannot cover all surface construction methods used for geometric shape design. Instead, we shall restrict ourselves to a number of essential methods to construct ruled surfaces, surfaces of revolution, and swept surfaces and show you how to use these techniques for interactive design of surfaces.

7.2  Primitive surfaces

Plane: An infinite plane is uniquely defined by a point p on the plane and a unit normal e₃ of the plane as shown in Figure 7.1.

To represent a bounded (or a finite) plane, we need two additional unit vectors that furnish the coordinate system, known as the frame, so that we can represent the plane parametrically over some finite intervals. Let e₁ and e₂ be two unit vectors such that e₁, e₂, and e₃ are mutually orthogonal. Then, a finite plane may be defined as

• Find the coordinate index i such that l = min{|e₃¹|, |e₃²|, |e₃³|}, where |e₃¹|, |e₃²|, and |e₃³| are the three coordinates of e₃.
• Let e₁ⁱ = 1 and other coordinates of e₁ be zero.
• Set e₁ = e₁-le₃.
• Normalize e₁ to a unit vector, i.e., e₁/Ö{1-l²}. Since l is the minimum absolute coordinate of e₃, it is then numerically stable to normalize e₁.

Implementing the above algorithm in Java gives

lambda = Math.abs(e_3.x);
e_1.x = 1.0;
e_1.y = e_1.z = 0.0;
if (lambda > Math.abs(e_3.y))
{
   lambda = e_3.y;
   e_1.y = 1.0;
   e_1.x = 0.0;
}
if (lambda < Math.abs(e_3.z))
{
   lambda = e_3.z;
   e_1.z = 1.0;
   e_1.x = e_1.y = 0.0;
}
temp = Math.sqrt(1 - lambda * lambda);
e_1.x = (e_1.x - lambda * e_3.x) / temp;
e_1.y = (e_1.y - lambda * e_3.y) / temp;
e_1.z = (e_1.x - lambda * e_3.z) / temp;

The above algorithm will be used throughout this book to compute an orthogonal unit vector with respect to another unit vector e₃ whenever there is a need.

We now compute the first order partial derivatives. Differentiating S(u,v) with respect to u and v gives:

Sphere: A sphere is uniquely defined by the center C and radius r of the sphere. In order to represent a sphere in a vector space, we usually require the specification of a unit vector e₃ pointing from the center to the ``north pole'' (actually, an arbitrary point on the sphere) as shown in Figure 7.2.

Accordingly, we may construct two unit vectors e₁ and e₂ such that e₁, e₂, and e₃ are mutually orthogonal and hence represent the sphere in the local frame (i.e., e₁, e₂, and e₃) with origin coincides with the center of the sphere. Let the angle parameter u Î [0,p] be oriented from the south to the north pole and v Î [0,2p] be, by default, oriented in accordance with the right-hand rule with respect to the vector e₃ (see Figure 7.2). Then, projecting an arbitrary point p onto the e₃ axis and the e₁-e₂ plane gives -rcosu and rsinu respectively. If rsinu is further projected onto the e₁ and e₂ axes, we obtain respectively the coordinates on these two axes as rsinucosv and rsinusinv. Consequently, we can write p in the local coordinate system as

The above code is a direct implementation of the math equations discussed previous. Save the above code under the name GSphere.java in the directory GeomLib and compile it to get GSphere.class. Let us now create the applet ShowSphere to display a predefine sphere. The code is listed below.

Save the above code in Ch6 and compile it to get ShowSphere.class. We can then create a HTML file to invoke the applet.

Cone/Cylinder: A cone is defined by the base ellipse r(v) and the major half angle aÎ [0,p/2) as shown in Figure 7.3.

Let u be the arc length parameter along the the direction of the normal of the base ellipse. Then, the scale factor for an ellipse at u to the base ellipse is

Differentiating S(u,v) with respect to u and v gives the first order partial derivatives:

The above code is a direct implementation of the math equations discussed previously. Let's save it in the directory GeomLib and compile it to get GCone.class. We design the applet ShowCone that allows you to display the predefined cone. You can also click the scroll bar on the applet to change the angle a so as to obtain different cones. The source code of the ShowCone is given below.

Save the code in the directory Ch6 and compile it to get ShowCone.class. Then, write a HTML file to invoke the applet.

Torus: A torus is generated by moving a minor (or outside) circle along a major (or inner circle) with the center of the minor circle on the major circle as shown in Figure 7.4.

We denote the center, radius, and unit normal of the major circle respectively by C, R, and e₃. If this major circle is parametrized by v in radian, we assume that the flow direction of v agrees with the right-hand rule with respect to the normal e₃. We also denote the radius and parameter of the minor circle by r and u. In order to represent the torus in a vector-valued parametric form, we choose two unit vectors e₁ and e₂ such that e₁, e₂, and e₃ are mutually orthogonal. Accordingly, the three unit vector and the center of the major circle furnish the local frame.

For v = 0, an arbitrary point on the torus in the local frame is given by

The above code is a straightforward implementation of the math equations discussed previously. Save this code in the directory GeomLib and compile it to get GTorus.class. We then design the ShowTorus class that allows you to display a torus and change the radius of the major circle by clicking the scroll bar. The source code is listed as follows.

7.3  Mathematical background on Bézier Surface

We start our discussion with a straight line. Representing the straight line by a Bézier curve of degree one (i.e., order 2) defined over v Î [0,1] gives

• For a fixed u parameter, å_{i = 0}ⁿ p_i,jB_i,n(u) are constant and denoted by p_j. Hence, S(u,v) = å_{j = 0}^m p_j B_j,m(v) is a Bézier curve in v and is called the iso-u curve.
• We rewrite the surface as

  S(u,v) = n å i = 0  æ è m å j = 0  pi,j Bi,m(v) ö ø Bj,n(u).

  Then, for a fixed v parameter, å_{j = 0}^m p_i,jB_i,m(v) are constant and denoted by p_i. Therefore, S(u,v) = å_{i = 0}ⁿ p_i B_i,n(u) is a Bézier curve in u and is called the iso-v curve.

We now consider the evaluation of Bézier surfaces. Mathematically speaking, we do not need a new algorithm to evaluate a point on a Bézier surface with respect to the given u and v parameters. We can use the bases evaluation method discussed in the section of Mathematical background on Bézier curve to compute the basis functions B_i,n(u) and B_j,m(v). Let BFunct_u[i] and BFunct_v[j] store B_i,n(u) and B_j,m(v) respectively. Then, the point on the Bézier surface is computed as follows:

   rx = ry = rz = 0.0;
   for (j=0; j<=m; j++)
   {
      for (i=0; i<=n; i++)
      {
         temp = BFunct_u[i] * BFunct_v[j];
         rx += x[j][i] * temp;
         ry += y[j][i] * temp;
         rz += z[j][i] * temp;
      }
   }

It is also true that we do not need a new algorithm to subdivide a Bézier surface. For example, if we want to subdivide the surface given in (7-2) at u = 0.5, we simply call the SubdivideCurve() method discussed in the previous chapter to halve

    S(u,v) = å_{j = 0}^m(å_{i = 0}ⁿp_i,j¹ B_i,n(u)+ å_{i = 0}ⁿp_i,j² B_i,n(u))B_j,m(v)

             = å_{j = 0}^m(å_{i = 0}ⁿp_i,j¹ B_i,n(u))B_j,m(v)+ å_{j = 0}^m(å_{i = 0}ⁿp_i,j² B_i,n(u))B_j,m(v).

As it is seen, the surface is now subdivided into two.

Analogous to a rational Bézier curve, a rational Bézier surface is given by

We now move to the next section to design an applet that allows you to display and manipulate a Bézier surface.

7.4  Creating Bézier surface applet

Let us first implement the GBezierSurface class that defines a Bézier surface and provides member methods for data setting and geometric processing. Similar to the GBezierCurve class, an one-dimensional array is used in GBezierSurface class to store the control vertices of a Bézier surface as shown below:

poles[]={x₀₀,x₁₀,¼,x_n0,¼,x_0m,x_1m,¼, x_nm, y₀₀,y₁₀,¼,y_n0,¼}

Therefore, there are (m+1)×(n+1) elements for each dimension. The source code of GBezierSurface is listed as follows.

As you may note, we did not implement a method to evaluate the Bézier basis functions. This is because the computation of basis functions for a surface is the same as that for a curve. Therefore, when there is a need to evaluate the Bézier basis functions, we can create a Bézier curve object to access its member method EvalBases() defined in the GBezierCurve class. Since this curve object does not hold any control vertices, we may consider it as a ``dummy'' object. Save the above code under the name GBezierSurface.java in the GeomLib directory and compile it to get GBezierSurface.class.

Let us now design an applet that allows you to display and manipulate a predefined Bézier surface with option to change the density of grids. We call this applet the ShowBezSurf and type the following code.

What we need to emphasize in the above code is the way to display the control polyhedron and the resulting surface. The polyhedron is shown by plotting the polygons first in u-direction, then in v direction. The surface is displayed by plotting a number of iso-v curves and then iso-u curves. Therefore, the surface is actually seen through iso-curves that form the grids (or meshes) on the surface. This is an inexpensive method to visualize the surface. Save the above code under the name ShowBezSurf.java in the directory Ch6 and compile it to get ShowBezSurf.class. Then, we write a HTML file to invoke the applet. If everything works fine, you should see an applet simialr to that in Figure 7.5.

In CAGD applications, we often encounter a case where the shape of a single Bézier surface is not rich enough to describe a complicated geometric object (e.g., a arc body). One solution to such a problem is to blend a number of Bézier surfaces together with some constraints. The surface that is made of several Bézier surfaces is called a composite Bézier surface. Although a composite Bézier surface can meet the requirement of describing a complicated geometric shape, it is often not preferred in applications since it usually requires more control vertices than a B-spline surface does. Therefore, let us move to the next section to study the mathematics on B-spline surfaces.

7.5  Mathematical background on B-spline Surface

Let a B-spline curve of order k_v with q control points be given by

• For a fixed u parameter, å_{i = 0}^p-1 d_i,j N_i,ku(u) are constant and denoted by d_j. Hence, S(u,v) = å_{j = 0}^q-1 d_j N_j,kv(v) is a B-spline curve in v and is called the iso-u curve.
• We rewrite the surface as

  S(u,v) = p-1 å i = 0  æ è q-1 å j = 0  di,j Ni,kv(v) ö ø Nj,ku(u).

  Then, for a fixed v parameter, å_{j = 0}^q-1 d_i,jN_i,kv(v) are constant and denoted by d_i. Therefore, S(u,v) = å_{i = 0}^p-1 d_i N_i,ku(u) is a B-spline curve in u and is called the iso-v curve.

We now consider the evaluation of B-spline surfaces. Mathematically speaking, we do not need a new algorithm to evaluate a point on a B-spline surface. Assume that we want to evaluate a B-spline surface with respect to the parameters u Î [u_ju, u_ju+1) and v Î [v_jv, v_jv+1). We then use the basis evaluation method discussed in the section of Mathematical background on B-spline curve to compute the basis functions N_i,ku(u) and N_j,kv(v). Due to the local support property of B-splines, the non-vanishing basis functions are N_ju-ku+1(u),¼, N_ju(u) and N_jv-kv+1(v),¼, N_jv(v). Therefore, the point on the surface is given by

It is also true that we do not need a new algorithm to insert a knot into the existing u and v knot vectors of a B-spline surface. Suppose that we want to insert a new knot, say [`u], into the u knot vector of the surface given in (7-3) such that there are m [`u] in the existing u knot vector. We then simply apply the same knot insertion algorithm discussed in the section of Mathematical background on B-spline curve to each B-spline curve

Analogous to a rational B-spline curve, a rational B-spline surface is given by

Enough for the math. Let us move to the next section to design two applets that allow you to display and manipulate B-spline surfaces.

7.6  Creating B-spline surface applet

poles[]={x_0,0,x_1,0,¼,x_p-1,0,¼,x_0,q-1, x_1,q-1, ¼, x_p-1,q-1,  y_0,0,y_1,0,¼, y_p-1,0, ¼}

Therefore, there are p×q elements for each dimension. The actual source code of GBsplineSurface is listed below.

Similar to the way we implemented GBezierSurface, we also create ``dummy'' B-spline curve objects to access member methods defined in the GBsplineCurve class, which avoids duplication of the same code. The EvalPoint() and EvalIsoCurve() methods are simply the implementation of the method discussed in the previous section. Save this code under the name GBsplineSurface.java in the GeomLib directory and compile it to get GBsplineSurface.class.

We now create the first applet ShowBspSurf that displays a B-spline surface with option to choose the density of grids. It should be pointed out that the GBsplineSurface class can handle generic polynomial B-spline surfaces. For simplicity, however, the surface in ShowBspSurf applet is predefined. When you understand the implementation detail, you can either modify or create a new applet that allows you to generate a sophisticated B-spline surface. Let us type the following code.

Save the code under the name ShowBspSurf.java in the directory Ch6 and compile it to get ShowBspSurf.class. Then, we write a HTML web page to invoke the applet. If everything is okay, you should be able to see an applet similar to that in Figure 7.6.

The B-spline surfaces in the applets are created by specifying the control vertices that form the control polyhedrons of the surfaces. In subsequent sections, we shall show you some important surface construction methods that start by designing curves then extend them to generate ruled surfaces, surfaces of revolution, and swept surfaces.

7.7  Ruled surface

Differentiating S(u,v) with respect to u and v gives the first order partial derivatives:

In what follows, we shall design an applet that allows you to create different ruled surfaces by manipulating two B-spline curves r₁(u) and r₂(u). Let us first define the GRuledSurface class of which the source code is given below.

The GRuledSurface class has two member variables bspcv1 and bspcv2 that define the ruled surface. The member method UpdateCurve is used to update one of the base curves. The EvalPoint method evaluate a point on the surface with respect to the given u and v parameters. Finally, the EvalIsovCurve method computes a set of points with respect to given v. It is noted that we do not need a method to evaluate iso-u curves since they are straight lines. Save the above code under the name GRuledSurface.java in the directory GeomLib and compile it to get GRuledSurface.class.

We now design the applet ShowRuleSurf that has two buttons Design Curve and Create Surface. When you click on the Design Curve button, you can manipulate the predefined B-spline curves. When you are satisfied with the shapes of two curves, you then click the Create Surface to construct a ruled surface. The source code is given as follows.

We save the above code under the name ShowRuleSurf.java in the directory Ch6 and compile it to get ShowRuleSurf.class. Then, write a HTML web page ShowRuleSurf.html to invoke the applet which looks similar to Figure 7.7.

7.8  Surface of revolution

Let C(u) = r(u)-s. Then, the three coordinates of C(u) in the local frame are given by

    S(u,v) = s + E₁(u)(cosve₁+ sinve₂) + E₂(u)(-sinve₁+cosve₂ ) + E₃(u)e₃

             = s+(E₁(u)cosv-E₂(u)sinv)e₁ + (E₁(u)sinv+E₂(u)cosv)e₂ + E₃(u)e₃

If we differentiate S(u,v) with respect to u and v, we then obtain the first order partial derivatives:

       [(¶S)/(¶u)] = (E₁^¢(u)cosv-E₂^¢(u) sinv )e₁ +(E₁^¢(u)sinv + E₂^¢(u)cosv)e₂ + E₃^¢(u)e₃,

       [(¶S)/(¶v)] = -(E₁(u)sinv+E₂(u)cosv)e₁+(E₁(u)cosv -E₂(u)sinv)e₂,

where E₁^¢(u) = C^¢(u)·e₁ = r^¢(u)·e₁, E₂^¢(u) = r^¢(u)·e₂, and E₃^¢ (u) = r^¢(u) ·e₃.

Having understood the mathematics on a surface of revolution, we can now design the GRevolSurface class of which the code is given below:

A surface of revolution is defined by a profile curve and an axis which in turn is defined by a base point and the direction vector e₃. Therefore, the GRevolSurface class has respectively three member variables bspcv, basePoint, and e_3. We also declare three additional member variables e_1, e_2 and sigma that represent e₁, e₂, and s. They are declared as private since there is no need for users to be involved with them. Similar to other geometric surfaces classes, we have some member methods that deal with updating and evaluating. Save the above code under the name GRevolSurface.java in GeomLib and compile it to get GRevolSurface.class.

Let's create the applet ShowRevSurf that allows us to generate and manipulate surfaces of revolution. The source code of ShowRevSurf is given below.

We save the above code under the name ShowRevSurf.java in the directory GeomLib and compile it to get ShowRevSurf.class. We then write a HTML web page named ShowRevSurf.html to invoke the applet. If verything works well, we should see an applet similar to that in Figure 7.8.

7.9  Swept surface