/* Date: Wed, 10 Dec 86 22:49:14 pst From: michael zyda Message-Id: <8612110649.AA26225@nps-cs> To: info-iris@sumex-aim One of the frustrating things about the IRIS is its lack of support for filled surface patches. One of my students has put together a package that looks like the IRIS patch routines but produces a filled surface instead. Also included with this is a simple lighting model... Michael J. Zyda Naval Postgraduate School Code 52, Dept. of Computer Science Monterey, California 93943 (408) 646-2305 zyda@nps-cs.arpa ucbvax!dual!lll-crg!nps-cs!zyda Files: Makefile patch3.c (main program) npatch.c (patch replacement library) lightpoly.c (light a polygon routine) npoly_orient.c (compute a polygon's vertex orientation) */ "Makefile" ************************************************** CFLAGS = -Zf ALL = npatch \ npoly_orient \ lightpoly \ patch3 all: $(ALL) clean: rm -f *.o delete: rm -f *.o $(ALL) npatch: cc -c npatch.c $(CFLAGS) npoly_orient: cc -c npoly_orient.c $(CFLAGS) lightpoly: cc -c lightpoly.c $(CFLAGS) patch3: cc -o patch3 patch3.c npatch.o npoly_orient.o lightpoly.o $(CFLAGS) -Zg patch3: npatch.o npoly_orient.o lightpoly.o patch3.c ***************************************************** /* this is file patch3.c It draws the Bezier patch of program patch.c and changes one of the values in the control point array each time through the loop if MOUSE3 or MOUSE1 is pressed. MOUSE3 increases the Y coordinate of the control point. MOUSE1 decreases the Y coordinate of the control point. MOUSE2 allows the change of which control point is selected. All three mouse buttons pressed exits. Pressing the outer two mouse buttons (MOUSE1 and MOUSE3) toggles on/off the filling of the surface patch. */ #define DASHED 1 #include "gl.h" #include "device.h" Matrix bezier = { -1.0, 3.0, -3.0, 1.0, 3.0, -6.0, 3.0, 0.0, -3.0, 3.0, 0.0, 0.0, 1.0, 0.0, 0.0, 0.0 }; #define BEZIER 1 /* Control point arrays... These arrays are separate for each coordinate, i.e. geomx holds all the x coordinates, geomy holds all the y coordinates, and geomz holds all the z coordinates. */ Coord geomx[4][4] = { 0.0, 100.0, 200.0, 300.0, 0.0, 100.0, 200.0, 300.0, 700.0, 600.0, 500.0, 400.0, 700.0, 600.0, 500.0, 400.0 }; Coord geomy[4][4] = { 400.0, 500.0, 600.0, 700.0, 0.0, 100.0, 200.0, 300.0, 0.0, 100.0, 200.0, 300.0, 400.0, 500.0, 600.0, 700.0 }; Coord geomz[4][4] = { 100.0, 200.0, 300.0, 400.0, 100.0, 200.0, 300.0, 400.0, 100.0, 200.0, 300.0, 400.0, 100.0, 200.0, 300.0, 400.0 }; int toggle; /* a toggle switch to turn the user's routine on/off */ main(argc, argv) int argc; char *argv[]; { int light_poly(); /* what user routine we will use */ static long row = 1; /* currently selected row */ static long col = 2; /* currently selected col */ Colorindex wmask; /* writemask value */ int i; /* temp loop variable */ /* init the graphics package */ ginit(); /* put the iris into double buffer mode */ doublebuffer(); /* configure the iris */ gconfig(); /* enable all bit planes for writing */ wmask = (1 << getplanes()) -1; writemask(wmask); /* make sure we can see the cursor */ setcursor(0, RED, wmask); /* full screen viewport */ viewport(0,1023,0,767); /* orthographic viewport */ ortho(0.0,1023.0,0.0,767.0,1023.0,-1023.0); /* translate the picture */ translate(250.0,0.0,0.0); /* make the picture smaller */ scale(0.75,0.75,0.75); /* turn the cursor on */ curson(); /* set the picksize for the middle mouse button hits */ picksize(30,30); /* define a dashed linestyle */ deflinestyle(DASHED,0xf0f0); /* setup the Bezier patch stuff */ setupbezier(argc); /* set up a color ramp... */ for(i = 0 ; i < 256 ; i ++) { mapcolor(i + 8, i , i , i ); } while(TRUE) { if(getbutton(MOUSE1) && getbutton(MOUSE2) && getbutton(MOUSE3)) { break; } /* turn on patch filling */ if(getbutton(MOUSE1) && !getbutton(MOUSE2) && getbutton(MOUSE3)) { User_Routine(toggle = !toggle); ringbell(); while(getbutton(MOUSE1) && getbutton(MOUSE3)); } /* should we increase the Y coordinate? */ if(getbutton(MOUSE3) && !getbutton(MOUSE2) && !getbutton(MOUSE1)) { /* increase the y coordinate */ geomy[row][col] = geomy[row][col] * 1.10; } /* should we decrease the Y coordinate? */ if(getbutton(MOUSE1) && !getbutton(MOUSE2) && !getbutton(MOUSE3)) { /* decrease the Y coordinate */ geomy[row][col] = geomy[row][col] * 0.90; } /* should we change which control point we are using ? */ if(getbutton(MOUSE2) && !getbutton(MOUSE1) && !getbutton(MOUSE3)) { resetcontrolpt(&row,&col,geomx,geomy,geomz); } /* redraw the picture */ color(CYAN); clear(); /* draw the Bezier patch */ linewidth(1); color(BLACK); npatch(geomx,geomy,geomz); /* highlight the changing control point's loc */ controlptoutlines(geomx,geomy,geomz); /* put up some text */ linewidth(1); color(BLACK); cmov2i(200,10); charstr("The Bezier Patch"); swapbuffers(); } /* clear the screen reset */ color(BLACK); clear(); textinit(); swapbuffers(); gflush(); gexit(); } /* this is function setupbezier */ setupbezier(argc) { /* define a basis matrix called BEZIER */ ndefbasis(BEZIER,bezier); /* specify a basis for each direction in the patch */ npatchbasis(BEZIER,BEZIER); /* curve segments will be drawn in the s direction. curve segments will be drawn in the t direction */ npatchcurves(10,10); /* specify the minimum number of line segments for each direction. The actual number of line segments is a multiple of the number of curve segments being drawn in the opposite direction. */ npatchprecision(20,20); Set_User_Routine_for_npatch(light_poly); User_Routine(argc > 1 ? (toggle = 1) : (toggle = 0)); } /* this is controlptoutlines This routine draws the outlines of the input control point arrays. */ controlptoutlines(geomx,geomy,geomz) float geomx[4][4], geomy[4][4], geomz[4][4]; /* control pt arrays */ { long i,j; /* temp locations for the drawing of the control pts */ /* make the control pt lines red */ color(RED); setlinestyle(1); linewidth(3); /* loop in one direction */ for(i=0; i < 4; i=i+1) { for(j=0; j < 4; j=j+1) { if(j == 0) { move(geomx[i][j],geomy[i][j],geomz[i][j]); } else { draw(geomx[i][j],geomy[i][j],geomz[i][j]); } } } /* loop in the other direction */ for(j=0; j < 4; j=j+1) { for(i=0; i < 4; i=i+1) { if(i == 0) { move(geomx[i][j],geomy[i][j],geomz[i][j]); } else { draw(geomx[i][j],geomy[i][j],geomz[i][j]); } } } setlinestyle(0); } /* this is resetcontrolpt.c This routine uses the IRIS pick mechanism to select which control pt should be modified. */ resetcontrolpt(row,col,geomx,geomy,geomz) long *row, *col; /* currently selected row/col of the control pt array */ float geomx[4][4], geomy[4][4], geomz[4][4]; /* control pt arrays */ { short i,j; /* loop temps */ long nhits; /* number of menu pick hit sets */ #define MAXNAMES 50 /* max numberof names spaces */ short names[MAXNAMES]; /* array to hold returned names */ short k,l; /* temp loop variables */ short namecount; /* name counter for pick mode */ /* we have a hit on mouse2, this means we should do iris pick mode to find out what we hit. */ /* save the current view matrix or else... */ pushmatrix(); /* turn on pick mode */ pick(names,MAXNAMES); /* re-do the ortho viewing matrix */ ortho(0.0,1023.0,0.0,767.0,1023.0,-1023.0); /* translate the picture */ translate(250.0,0.0,0.0); /* make the picture larger */ scale(0.75,0.75,0.75); /* draw the pick picture with the associated names. These objects will not actually be displayed. */ namecount = -1; /* initialname for the stack */ for(i=0; i < 4; i=i+1) { for(j=0; j < 4; j=j+1) { namecount = namecount +1; pushname(namecount); /* put a name on the stack */ /* draw a point at each contrl pt as a hit target. */ pnt(geomx[i][j],geomy[i][j],geomz[i][j]); /* get rid of the last name put on */ popname(); } } /* endfor i */ /* turn off pick mode */ nhits = endpick(names); /* restore the matrix */ popmatrix(); /* did we get any hits? */ if(nhits > 0) { /* we had some hits */ j = -1; /* set up ptr for namesbuffer decoding */ /* for each hit set... */ for(i=0; i < nhits; i=i+1) { j = j + 1; /* j points at count of names in names for this hit set */ k = names[j]; /* k is the count of names this hit set */ /* loop for all names this set */ for(l=0; l < k; l=l+1) { j = j + 1; /* j points at the next names entry */ /* compute this row and column */ *row = names[j]/4; *col = names[j]%4; } /* endfor l ... */ } /* endfor i... (endfor all hit sets) */ } /* endif nhits > 0 */ } /* this is routine light_poly The purpose of this routine is to light the created triangles if called. */ int light_poly(Triangle) Coord Triangle[3][3]; { lightthepoly(Triangle, 3, 350.0, -1150.0, 250.0, 350.0, 1150.0, 0.0, 50, 264); } npatch.c **************************************************** /* this is file npatch.c It is a collection of routines that shadow the IRIS patch routines. The purpose of this set of routines is to provide filled surface patches instead of just wireframes. */ #include "gl.h" #include "stdio.h" #define OFF 0 typedef struct { Coord x; Coord y; Coord z; } Point; /* how points are stored */ static struct list_elem { long nid; /* matrix id number */ float *nmatrix; /* pointer to the matrix */ struct list_elem *next_elem_ptr; }; /* struct used to track user supplied matrices */ static Matrix Precision_Matrix_U ; /* used for forward difference computation in the U direction */ static Matrix Precision_Matrix_V ; /* used for forward difference computation in the V direction */ static int _Default_system_routine(); /* what the system does if you do not specify */ static float *current_U_basis; /* pointer to the currently defined U basis matrix */ static float *current_V_basis; /* pointer to the currently defined V basis matrix */ static long UCURVES; /* how many curves in the U direction */ static long VCURVES; /* how many curves in the V direction */ static long USEGMENTS; /* how may curve segments in the U direction */ static long VSEGMENTS; /* how may curve segments in the V direction */ static struct list_elem *head_of_list = NULL; /* pointer to linked list of user supplied basis matrices */ static int (*_User_routine)() = _Default_system_routine; /* set initial user routine to default */ static int _User_routine_is_on = OFF; /* initially do not call user's function */ /******************************************************************************* This routine is what the system automatically does if the user does not supply a particular intercept function. *******************************************************************************/ static int _Default_system_routine(Triangle) Coord Triangle[3][3]; { poly(3, Triangle); } /******************************************************************************* This function allows the user to supply an intercept function that will be used to manipulate the triangles that are generated in the the decomposition of the surface patch. *******************************************************************************/ void Set_User_Routine_for_npatch(routine) int (*routine)(); { _User_routine = routine; } /******************************************************************************* This function allows the user to turn his intercept function on or off at will. *******************************************************************************/ void User_Routine(boolean) int boolean; { _User_routine_is_on = boolean; } /******************************************************************************* This function allows the user to reset his intercept routine to the same routine used as the system default. *******************************************************************************/ void Set_Default_Routine_for_npatch() { _User_routine = _Default_system_routine; } /******************************************************************************* This function is equivalent to the IRIS defbasis function in that it allows the user to define a basis matrix for use in the generation of patches. matrix is saved and is associated with id. id may then be used in subsequent calls to npatchbasis. *******************************************************************************/ void ndefbasis(id, matrix) long id; Matrix matrix; { static int has_been_called = FALSE; struct list_elem *new_elem_to_add; struct list_elem *walking_ptr; float *pmatrix; int row; int column; /* get memory for the data structure */ new_elem_to_add = (struct list_elem *) malloc(sizeof(struct list_elem)); pmatrix = (float *) calloc(sizeof(Matrix), sizeof(float)); /* make a copy of the basis matrix passed in by the user */ for(row = 0; row < 4; row++) for(column = 0; column < 4; column++) { *(pmatrix + (4*row) + column) = matrix[row][column]; } /* put in the information into the data structure */ new_elem_to_add -> nid = id; new_elem_to_add -> nmatrix = pmatrix; /* determine how to add this information into the linked list of basis matrices */ switch(has_been_called) /* does a list already exist ? */ { case TRUE : /* Yes */ walking_ptr = head_of_list; /* point to the beginning of the list */ /* traverse the list */ while((walking_ptr -> nid != id) && (walking_ptr -> next_elem_ptr != head_of_list)) { walking_ptr = walking_ptr -> next_elem_ptr; } if(walking_ptr -> nid == id) /* this id already exists so we can re-use the data structure already in the list */ { /* get rid of the old basis matrix */ cfree(walking_ptr -> nmatrix, sizeof(Matrix), sizeof(float)); /* point to the replacement matrix */ walking_ptr -> nmatrix = pmatrix; /* get rid of the un-needed data structure */ cfree(new_elem_to_add, 1, sizeof(struct list_elem)); } else /* the id does not exist in the current linked list of basis matrices so we must add it into it */ { /* manipulate the pointers to add it in */ new_elem_to_add -> next_elem_ptr = head_of_list -> next_elem_ptr; head_of_list -> next_elem_ptr = new_elem_to_add; } break; case FALSE : /* No */ /* no linked list of basis matrices exists so we start up one */ head_of_list = (struct list_elem *) malloc(sizeof(struct list_elem)); head_of_list -> next_elem_ptr = head_of_list; new_elem_to_add -> next_elem_ptr = head_of_list -> next_elem_ptr; head_of_list -> next_elem_ptr = new_elem_to_add; has_been_called = TRUE; /* make sure we do not do this again */ break; default : break; } } /******************************************************************************* This function is equivalent to the IRIS patchbasis function in that it sets the current basis matrices for both the U and V parametric directions of a surface patch. The current U and V bases are used when the npatch command is issued. *******************************************************************************/ int npatchbasis(uid, vid) long uid, vid; { struct list_elem *walking_ptr; if(head_of_list == NULL) /* no linked list of basis matrices */ { fprintf(stderr, "\nnpatchbasis: no basis matrices defined\n\n"); exit(-1); } /* traverse the list looking for the desired U basis matrix */ walking_ptr = head_of_list; while((walking_ptr -> nid != uid) && (walking_ptr -> next_elem_ptr != head_of_list)) { walking_ptr = walking_ptr -> next_elem_ptr; } if((walking_ptr -> nid != uid) && (walking_ptr -> next_elem_ptr == head_of_list)) /* the desired U basis matrix does not exist in the linked list */ { fprintf(stderr, "\nnpatchbasis: undefined U basis matrix %d\n\n", uid); exit(-1); } current_U_basis = walking_ptr -> nmatrix; /* traverse the list looking for the desired V basis matrix */ walking_ptr = head_of_list; while((walking_ptr -> nid != vid) && (walking_ptr -> next_elem_ptr != head_of_list)) { walking_ptr = walking_ptr -> next_elem_ptr; } if((walking_ptr -> nid != vid) && (walking_ptr -> next_elem_ptr == head_of_list)) /* the desired V basis matrix does not exist in the linked list */ { fprintf(stderr, "\nnpatchbasis: undefined V basis matrix %d\n\n", vid); exit(-1); } current_V_basis = walking_ptr -> nmatrix; } /******************************************************************************* This function is similar to the IRIS patchcurves command. ucurves and vcurves set the subdivision parameters used in decomposing the surface patch. An individual patch will be decomposed into a ucurves * vcurves grid with each grid generating two triangles. *******************************************************************************/ void npatchcurves(ucurves, vcurves) long ucurves, vcurves; { /* prevent a negative number of curves */ UCURVES = (ucurves < 1 ? 10 : ucurves); VCURVES = (vcurves < 1 ? 10 : vcurves); /* set up the Precision_Matrix_U used in the forward difference along the U direction */ Precision_Matrix_U[0][0] = 6.0/(float)(UCURVES * UCURVES * UCURVES); Precision_Matrix_U[1][0] = 6.0/(float)(UCURVES * UCURVES * UCURVES); Precision_Matrix_U[1][1] = 2.0/(float)(UCURVES * UCURVES); Precision_Matrix_U[2][0] = 1.0/(float)(UCURVES * UCURVES * UCURVES); Precision_Matrix_U[2][1] = 1.0/(float)(UCURVES * UCURVES); Precision_Matrix_U[2][2] = 1.0/(float)(UCURVES); Precision_Matrix_U[3][3] = 1.0; /* set up the Precision_Matrix_V used in the forward difference along the V direction */ Precision_Matrix_V[0][0] = 6.0/(float)(VCURVES * VCURVES * VCURVES); Precision_Matrix_V[1][0] = 6.0/(float)(VCURVES * VCURVES * VCURVES); Precision_Matrix_V[1][1] = 2.0/(float)(VCURVES * VCURVES); Precision_Matrix_V[2][0] = 1.0/(float)(VCURVES * VCURVES * VCURVES); Precision_Matrix_V[2][1] = 1.0/(float)(VCURVES * VCURVES); Precision_Matrix_V[2][2] = 1.0/(float)(VCURVES); Precision_Matrix_V[3][3] = 1.0; } /******************************************************************************* This function is similar to the IRIS patchprecision routine but is only used to maintain consistency with the IRIS routines. It's results are not used by any other function in the suite. *******************************************************************************/ void npatchprecision(usegments, vsegments) long usegments, vsegments; { /* prevent a negative number of segments */ USEGMENTS = (usegments < 1 ? 10 : usegments); VSEGMENTS = (vsegments < 1 ? 10 : vsegments); } /******************************************************************************* This function does the actual display for the surface patch. *******************************************************************************/ static void nspeckle(COORD_ARRAY_POINTER) Point *COORD_ARRAY_POINTER[4]; { register Point *Patch_array; /* get enough memory to hold all the points that will be generated */ Patch_array = (Point *) calloc((UCURVES+1)*(VCURVES+1), sizeof(Point)); { register Point *WHERE_IT_IS_AT; register unsigned int total_points; register unsigned int point_count; register unsigned int t_count; register unsigned int count; register unsigned int Row; register unsigned int Column; Matrix control_matrix; Matrix inter1; /* for every curve in the the U parametric direction */ for(point_count = 0, total_points = 0 ; point_count <= UCURVES ; point_count++) { /* build a control matrix for the current curve */ for(count = 0; count < 4; count++) { WHERE_IT_IS_AT = (Point *) (COORD_ARRAY_POINTER[count] + point_count); control_matrix[count][0] = WHERE_IT_IS_AT -> x; control_matrix[count][1] = WHERE_IT_IS_AT -> y; control_matrix[count][2] = WHERE_IT_IS_AT -> z; } /* generate the matrix to compute the forward difference on */ pushmatrix(); loadmatrix(control_matrix); multmatrix(current_V_basis); multmatrix(Precision_Matrix_V); getmatrix(inter1); popmatrix(); /* save all the points generated for that curve */ for(t_count = 0 ; t_count <= VCURVES ; t_count++, total_points++) { (Patch_array + total_points) -> x = inter1[3][0] ; (Patch_array + total_points) -> y = inter1[3][1] ; (Patch_array + total_points) -> z = inter1[3][2] ; /* do the forward difference */ for(Row = 3; Row > 0; Row--) for(Column = 0; Column < 4; Column++) inter1[Row][Column] = inter1[Row][Column] + inter1[Row - 1][Column]; } } } { Coord Triangle_1[4][3]; register Point *WHERE_IT_IS_AT; register unsigned int Row; register unsigned int Column; /* decompose the patch into its individual triangles */ for(Row = 0; Row < UCURVES; Row++) for(Column = 0 ; Column < VCURVES; Column++) { WHERE_IT_IS_AT = (Patch_array + ((VCURVES + 1) * Row) + Column); Triangle_1[0][0] = WHERE_IT_IS_AT -> x; Triangle_1[0][1] = WHERE_IT_IS_AT -> y; Triangle_1[0][2] = WHERE_IT_IS_AT -> z; Triangle_1[1][0] = (WHERE_IT_IS_AT + 1) -> x; Triangle_1[1][1] = (WHERE_IT_IS_AT + 1) -> y; Triangle_1[1][2] = (WHERE_IT_IS_AT + 1) -> z; Triangle_1[2][0] = (WHERE_IT_IS_AT + VCURVES + 1) -> x; Triangle_1[2][1] = (WHERE_IT_IS_AT + VCURVES + 1) -> y; Triangle_1[2][2] = (WHERE_IT_IS_AT + VCURVES + 1) -> z; Triangle_1[3][0] = (WHERE_IT_IS_AT + VCURVES + 2) -> x; Triangle_1[3][1] = (WHERE_IT_IS_AT + VCURVES + 2) -> y; Triangle_1[3][2] = (WHERE_IT_IS_AT + VCURVES + 2) -> z; if(_User_routine_is_on) /* does the user have an intercept routine ? */ { /* Yes */ (*_User_routine) (&Triangle_1[0][0]); (*_User_routine) (&Triangle_1[1][0]); } else /* no user routine, so we use the default */ { _Default_system_routine(&Triangle_1[0][0]); _Default_system_routine(&Triangle_1[1][0]); } } } /* return the memory used to the system */ cfree(Patch_array, (UCURVES+1)*(VCURVES+1), sizeof(Point)); } /******************************************************************************* This function rearranges the input matrices into a form readily used in computing points along the four curves defined by those matrices. It then computes those points for each curve using the technique of forwards differencing of a matrix. *******************************************************************************/ void npatch(geomx, geomy, geomz) Coord geomx[4][4], geomy[4][4], geomz[4][4]; { register Point *COORD_ARRAY_POINTER[4]; /* one control matrix for each curve */ Matrix ctrl_pts1; Matrix ctrl_pts2; Matrix ctrl_pts3; Matrix ctrl_pts4; /* load the appropriate control matrix for each curve */ ctrl_pts1[0][0] = geomx[0][0]; ctrl_pts1[0][1] = geomy[0][0]; ctrl_pts1[0][2] = geomz[0][0]; ctrl_pts1[1][0] = geomx[1][0]; ctrl_pts1[1][1] = geomy[1][0]; ctrl_pts1[1][2] = geomz[1][0]; ctrl_pts1[2][0] = geomx[2][0]; ctrl_pts1[2][1] = geomy[2][0]; ctrl_pts1[2][2] = geomz[2][0]; ctrl_pts1[3][0] = geomx[3][0]; ctrl_pts1[3][1] = geomy[3][0]; ctrl_pts1[3][2] = geomz[3][0]; ctrl_pts2[0][0] = geomx[0][1]; ctrl_pts2[0][1] = geomy[0][1]; ctrl_pts2[0][2] = geomz[0][1]; ctrl_pts2[1][0] = geomx[1][1]; ctrl_pts2[1][1] = geomy[1][1]; ctrl_pts2[1][2] = geomz[1][1]; ctrl_pts2[2][0] = geomx[2][1]; ctrl_pts2[2][1] = geomy[2][1]; ctrl_pts2[2][2] = geomz[2][1]; ctrl_pts2[3][0] = geomx[3][1]; ctrl_pts2[3][1] = geomy[3][1]; ctrl_pts2[3][2] = geomz[3][1]; ctrl_pts3[0][0] = geomx[0][2]; ctrl_pts3[0][1] = geomy[0][2]; ctrl_pts3[0][2] = geomz[0][2]; ctrl_pts3[1][0] = geomx[1][2]; ctrl_pts3[1][1] = geomy[1][2]; ctrl_pts3[1][2] = geomz[1][2]; ctrl_pts3[2][0] = geomx[2][2]; ctrl_pts3[2][1] = geomy[2][2]; ctrl_pts3[2][2] = geomz[2][2]; ctrl_pts3[3][0] = geomx[3][2]; ctrl_pts3[3][1] = geomy[3][2]; ctrl_pts3[3][2] = geomz[3][2]; ctrl_pts4[0][0] = geomx[0][3]; ctrl_pts4[0][1] = geomy[0][3]; ctrl_pts4[0][2] = geomz[0][3]; ctrl_pts4[1][0] = geomx[1][3]; ctrl_pts4[1][1] = geomy[1][3]; ctrl_pts4[1][2] = geomz[1][3]; ctrl_pts4[2][0] = geomx[2][3]; ctrl_pts4[2][1] = geomy[2][3]; ctrl_pts4[2][2] = geomz[2][3]; ctrl_pts4[3][0] = geomx[3][3]; ctrl_pts4[3][1] = geomy[3][3]; ctrl_pts4[3][2] = geomz[3][3]; { register Point *WHERE_TO_SAVE; register unsigned int Row; register unsigned int Column; register unsigned int point_count; register unsigned int count; float *matrix_pointer[4]; Matrix inter1; matrix_pointer[0] = (float *) ctrl_pts1; matrix_pointer[1] = (float *) ctrl_pts2; matrix_pointer[2] = (float *) ctrl_pts3; matrix_pointer[3] = (float *) ctrl_pts4; /* get enough memory to hold the points */ COORD_ARRAY_POINTER[0] = (Point *) calloc(UCURVES + 1, sizeof(Point)); COORD_ARRAY_POINTER[1] = (Point *) calloc(UCURVES + 1, sizeof(Point)); COORD_ARRAY_POINTER[2] = (Point *) calloc(UCURVES + 1, sizeof(Point)); COORD_ARRAY_POINTER[3] = (Point *) calloc(UCURVES + 1, sizeof(Point)); /* for each curve */ for(count = 0; count < 4; count++) { /* generate the matrix used in the forward difference for this curve */ pushmatrix(); loadmatrix(matrix_pointer[count]); multmatrix(current_U_basis); multmatrix(Precision_Matrix_U); getmatrix(inter1); popmatrix(); /* for each curve generate points in the U parametric direction */ for(point_count = 0 ; point_count <= UCURVES ; point_count++) { WHERE_TO_SAVE = (Point *) (COORD_ARRAY_POINTER[count] + point_count); WHERE_TO_SAVE -> x = inter1[3][0]; WHERE_TO_SAVE -> y = inter1[3][1]; WHERE_TO_SAVE -> z = inter1[3][2]; /* compute the forward difference */ for(Row = 3; Row > 0; Row--) for(Column = 0; Column < 4; Column++) inter1[Row][Column] = inter1[Row][Column] + inter1[Row - 1][Column]; } } } nspeckle(COORD_ARRAY_POINTER); /* call function to finish computations and display */ /* return the memory used back to the system */ cfree(COORD_ARRAY_POINTER[0],(UCURVES + 1), sizeof(Point)); cfree(COORD_ARRAY_POINTER[1],(UCURVES + 1), sizeof(Point)); cfree(COORD_ARRAY_POINTER[2],(UCURVES + 1), sizeof(Point)); cfree(COORD_ARRAY_POINTER[3],(UCURVES + 1), sizeof(Point)); } lightpoly.c ************************************************ /* this is file lightpoly.c It is a routine that computes lighting for a polygon based upon the angle between the Normal vector of the polygon and the direction to the light source. lightthepoly(xyz,ncoords,ax,ay,az,lx,ly,lz,colormin,colormax) xyz[][3] = floating coords of the polygon. ncoords = number of coordinates. ax,ay,az = interior point of the whole object. Used to determine outward facing normal of the polygon. This is the same point of reference that would be used for backface polygon removal. lx,ly,lz = vector pointing in direction of the light source. colormin, colormax = indices used for the colors assigned to this polygon. The user is responsible for setting up the color ramp. Note: the routine also puts the polygons out ordered counterclockwise with respect to the interior point for ease of backface polygon removal. */ #include #include #define PIDIV2 1.570796327 #define CLOCKWISE 1 #define COLUMN 3 lightthepoly(xyz,ncoords,ax,ay,az,lx,ly,lz,colormin,colormax) Coord xyz[][3]; unsigned int ncoords; Coord ax,ay,az; /* interior point of the whole object. */ Coord lx,ly,lz; /* direction to the light source */ int colormin,colormax; /* color min/max indices */ { Coord *txyz; /* temp coord hold */ register unsigned short int i,j; /* loop temps */ int npoly_orient(); /* direction test function */ Coord v1[3],v2[3]; /* vectors used to compute the polygon's normal */ Coord normal[3]; /* the polygon's normal */ Coord normalmag; /* normal's magnitude */ Coord lightmag; /* light's magnitude */ float dotprod; /* dot product of N and L */ float radians; /* angle between N and L */ unsigned short int colortouse; /* color to use in drawing the polygon */ /* allocate memory for the temporary array */ /* this is a two dimensional array of ncoords by 3 */ txyz = (Coord *) calloc((ncoords * 3), sizeof(Coord)); /* orient the polygon so that its counterclockwise with respect to the interior point */ if(npoly_orient(ncoords,xyz,ax,ay,az) == CLOCKWISE) { /* the polygon is clockwise, reverse it. */ for(i=0; i < ncoords; i=i+1) { for(j=0; j < COLUMN; j=j+1) { *(txyz + (i * COLUMN) + j) = xyz[ncoords-i-1][j]; } } } else { /* no need to reverse */ for(i=0; i < ncoords; i=i+1) { for(j=0; j < COLUMN; j=j+1) { *(txyz + (COLUMN * i) + j) = xyz[i][j]; } } } /* the coordinates are ordered counterclockwise in array txyz */ /* compute the normal vector for the polygon using the first three vertices...*/ /* compute the first vector to use in the computation */ v1[0] = *(txyz + 6) - *(txyz + 3); /* txyz[2][0] - txyz[1][0] */ v1[1] = *(txyz + 7) - *(txyz + 4); /* txyz[2][1] - txyz[1][1] */ v1[2] = *(txyz + 8) - *(txyz + 5); /* txyz[2][2] - txyz[1][2] */ /* compute the second vector to use in computing the normal */ v2[0] = *(txyz ) - *(txyz + 3); /* txyz[0][0] - txyz[1][0] */ v2[1] = *(txyz + 1) - *(txyz + 4); /* txyz[0][1] - txyz[1][1] */ v2[2] = *(txyz + 2) - *(txyz + 5); /* txyz[0][2] - txyz[1][2] */ /* the normal is v1 x v2 */ normal[0] = v1[1]*v2[2] - v1[2]*v2[1]; normal[1] = v1[2]*v2[0] - v1[0]*v2[2]; normal[2] = v1[0]*v2[1] - v1[1]*v2[0]; /* compute the magnitude of the normal */ normalmag = sqrt((normal[0]*normal[0]) + (normal[1]*normal[1])+ (normal[2]*normal[2])); /* compute the magnitude of the light */ lightmag = sqrt((lx * lx) + (ly * ly) + (lz * lz)); /* compute N . L (normal dot product with the light source direction) */ dotprod = ((normal[0] * lx) + (normal[1] * ly) + (normal[2] * lz)); /* compute the unit normal */ dotprod = ((dotprod/(normalmag * lightmag))); /* dotprod = cos(theta) of the angle between N and L. Convert to angle in radians */ radians = acos(dotprod); /* compute the color we should use */ if(-PIDIV2 <= radians && radians <= PIDIV2) { /* if the angle is negative, set to positive */ if(radians < 0.0) { radians = -radians; } colortouse = ((colormax-colormin)/PIDIV2)*(PIDIV2-radians)+colormin; } else { colortouse = colormin; } /* set the color */ color(colortouse); /* draw the poly */ polf(ncoords,txyz); /* free up memory allocated for the temporary array */ cfree(txyz, (ncoords * 3), sizeof(Coord)); } npoly_orient.c ********************************************* /* this is file npoly_orient.c Its purpose is to determine the orientation of a polygon with respect to its interior point. */ #include #include int npoly_orient(ncoords,xyz,xinside,yinside,zinside) unsigned int ncoords; Coord xyz[][3]; Coord xinside, yinside, zinside; { register unsigned short int i,j; /* loop temps */ Coord center[3]; /* center coordinate of the polygon */ Coord a[3], b[3]; /* vector hold locations for the vectors that run from the center coordinate to the points of the polygon */ Coord xn[3], xmn[3]; /* points on line containing normal that are on opposite sides of the plane containing the polygon. */ float distton; /* distance to point n from pt inside. */ float disttomn; /* distance to point -n from pt inside. */ Coord normal[3]; /* the normal vector computed from a x b */ /* compute the center coordinate of the polygon */ center[0] = 0.0; center[1] = 0.0; center[2] = 0.0; for(i=0; i < ncoords; i++) { for(j=0; j < 3; j++) { center[j] += xyz[i][j]; } } /* divide out by the number of coordinates */ for(j=0; j < 3; j++) { center[j] = center[j]/(float)ncoords; } /* check the first 2 coordinates of the polygon for their direction */ /* compute vector a. It runs from the center coordinate to coordinate 0 */ for(j=0; j < 3; j++) { a[j] = xyz[0][j] - center[j]; } /* compute vector b. It runs from the center coordinate to coordinate 1 */ for(j=0; j <3; j++) { b[j] = xyz[1][j] - center[j]; } /* compute a x b to get the normal vector */ normal[0] = a[1]*b[2] - a[2]*b[1]; normal[1] = a[2]*b[0] - a[0]*b[2]; normal[2] = a[0]*b[1] - a[1]*b[0]; /* compute point n, offset pt from center in direction of normal */ for(j=0; j < 3; j++) { xn[j] = center[j] + normal[j]; } /* compute point -n, offset pt from center in opposite direction from normal. */ for(j=0; j < 3; j++) { xmn[j] = center[j] - normal[j]; } /* compute the distance the inside pt is from point n */ distton = sqrt((xn[0] - xinside) * (xn[0] - xinside) + (xn[1] - yinside) * (xn[1] - yinside) + (xn[2] - zinside) * (xn[2] - zinside)); /* compute the distance the inside pt is from point -n */ disttomn = sqrt((xmn[0] - xinside) * (xmn[0] - xinside) + (xmn[1] - yinside) * (xmn[1] - yinside) + (xmn[2] - zinside) * (xmn[2] - zinside)); /* if the dist(n) < dist(-n), then n points back towards the inside point and is on the same side of the plane as inside. a x b is then clockwise. */ if(distton < disttomn) { return(1); /* clockwise */ } else { return(0); /* counterclockwise */ } }