Creating a spline surface involves taking the product of the same spline blending functions used for spline curves as follows
where the control points form a 2D array Pij. Most of the properties of the spline curve also apply to spline surfaces. For example
As an example the following show the surface given a 4x5 matrix of control points, degree 3 in both directions, the surface is sampled on a 30x40 grid.
Zero order continuity between spline surfaces is ensured if the control points along the joining edge are the same. First order continuity can be obtained by matching the control point slope across the boundary. This also applies to attempts to form a closed spline surface.
Forming a closed surface is easier if third order splines are being used because the surface is tangential to the end planes. In the following image the surface on the right starts and ends at the origin. The surface on the right has an extra control point so that the start and end planes have the same tangent.
/*
Create an example of a spline surfaces
*/
#define NI 3
#define NJ 4
XYZ inp[NI+1][NJ+1];
#define TI 3
#define TJ 3
int knotsI[NI+TI+1];
int knotsJ[NJ+TJ+1];
#define RESOLUTIONI 30
#define RESOLUTIONJ 40
XYZ outp[RESOLUTIONI][RESOLUTIONJ];
int main(argc,argv)
int argc;
char **argv;
{
int i,j,ki,kj;
double intervalI,incrementI;
double intervalJ,incrementJ;
double bi,bj;
/* Create a random surface */
srandom(1111);
for (i=0;i<=NI;i++) {
for (j=0;j<=NJ;j++) {
inp[i][j].x = i;
inp[i][j].y = j;
inp[i][j].z = (random() % 10000) / 5000.0 - 1;
}
}
/* Step size along the curve */
incrementI = (NI - TI + 2) / ((double)RESOLUTIONI - 1);
incrementJ = (NJ - TJ + 2) / ((double)RESOLUTIONJ - 1);
/* Calculate the knots */
SplineKnots(knotsI,NI,TI);
SplineKnots(knotsJ,NJ,TJ);
intervalI = 0;
for (i=0;i<RESOLUTIONI-1;i++) {
intervalJ = 0;
for (j=0;j<RESOLUTIONJ-1;j++) {
outp[i][j].x = 0;
outp[i][j].y = 0;
outp[i][j].z = 0;
for (ki=0;ki<=NI;ki++) {
for (kj=0;kj<=NJ;kj++) {
bi = SplineBlend(ki,TI,knotsI,intervalI);
bj = SplineBlend(kj,TJ,knotsJ,intervalJ);
outp[i][j].x += (inp[ki][kj].x * bi * bj);
outp[i][j].y += (inp[ki][kj].y * bi * bj);
outp[i][j].z += (inp[ki][kj].z * bi * bj);
}
}
intervalJ += incrementJ;
}
intervalI += incrementI;
}
intervalI = 0;
for (i=0;i<RESOLUTIONI-1;i++) {
outp[i][RESOLUTIONJ-1].x = 0;
outp[i][RESOLUTIONJ-1].y = 0;
outp[i][RESOLUTIONJ-1].z = 0;
for (ki=0;ki<=NI;ki++) {
bi = SplineBlend(ki,TI,knotsI,intervalI);
outp[i][RESOLUTIONJ-1].x += (inp[ki][NJ].x * bi);
outp[i][RESOLUTIONJ-1].y += (inp[ki][NJ].y * bi);
outp[i][RESOLUTIONJ-1].z += (inp[ki][NJ].z * bi);
}
intervalI += incrementI;
}
outp[i][RESOLUTIONJ-1] = inp[NI][NJ];
intervalJ = 0;
for (j=0;j<RESOLUTIONJ-1;j++) {
outp[RESOLUTIONI-1][j].x = 0;
outp[RESOLUTIONI-1][j].y = 0;
outp[RESOLUTIONI-1][j].z = 0;
for (kj=0;kj<=NJ;kj++) {
bj = SplineBlend(kj,TJ,knotsJ,intervalJ);
outp[RESOLUTIONI-1][j].x += (inp[NI][kj].x * bj);
outp[RESOLUTIONI-1][j].y += (inp[NI][kj].y * bj);
outp[RESOLUTIONI-1][j].z += (inp[NI][kj].z * bj);
}
intervalJ += incrementJ;
}
outp[RESOLUTIONI-1][j] = inp[NI][NJ];
printf("LIST\n");
/* Display the surface, in this case in OOGL format for GeomView */
printf("{ = CQUAD\n");
for (i=0;i<RESOLUTIONI-1;i++) {
for (j=0;j<RESOLUTIONJ-1;j++) {
printf("%g %g %g 1 1 1 1\n",
outp[i][j].x, outp[i][j].y, outp[i][j].z);
printf("%g %g %g 1 1 1 1\n",
outp[i][j+1].x, outp[i][j+1].y, outp[i][j+1].z);
printf("%g %g %g 1 1 1 1\n",
outp[i+1][j+1].x,outp[i+1][j+1].y,outp[i+1][j+1].z);
printf("%g %g %g 1 1 1 1\n",
outp[i+1][j].x, outp[i+1][j].y, outp[i+1][j].z);
}
}
printf("}\n");
/* Control point polygon */
for (i=0;i<NI;i++) {
for (j=0;j<NJ;j++) {
printf("{ = SKEL 4 1 \n");
printf("%g %g %g \n",inp[i][j].x,inp[i][j].y,inp[i][j].z);
printf("%g %g %g \n",inp[i][j+1].x,inp[i][j+1].y,inp[i][j+1].z);
printf("%g %g %g \n",inp[i+1][j+1].x,inp[i+1][j+1].y,inp[i+1][j+1].z);
printf("%g %g %g \n",inp[i+1][j].x,inp[i+1][j].y,inp[i+1][j].z);
printf("5 0 1 2 3 0\n");
printf("}\n");
}
}
}