/* * FUNCTION: * This file contains a number of utilities useful to 3D graphics in * general, and to the generation of tubing and extrusions in particular * * HISTORY: * Written by Linas Vepstas, August 1991 * Updated to correctly handle degenerate cases, Linas, February 1993 */ #include #include "port.h" #include "vvector.h" #define BACKWARDS_INTERSECT (2) /* ========================================================== */ /* * the Degenerate_Tolerance token represents the greatest amount by * which different scales in a graphics environment can differ before * they should be considered "degenerate". That is, when one vector is * a million times longer than another, changces are that the second will * be less than a pixel long, and therefore was probably meant to be * degenerate (by the CAD package, etc.) But what should this tolerance * be? At least 1 in onethousand (since screen sizes are 1K pixels), but * les than 1 in 4 million (since this is the limit of single-precision * floating point accuracy). Of course, if double precision were used, * then the tolerance could be increased. * * Potentially, this naive assumption could cause problems if the CAD * package attempts to zoom in on small details, and turns out, certain * points should not hvae been degenerate. The problem presented here * is that the tolerance could run out before single-precision ran * out, and so the CAD packages would perceive this as a "bug". * One alternative is to fiddle around & try to tighten the tolerance. * However, the right alternative is to code the graphics pipeline in * double-precision (and tighten the tolerance). * * By the way, note that Degernate Tolerance is a "dimensionless" * quantitiy -- it has no units -- it does not measure feet, inches, * millimeters or pixels. It is used only in the computations of ratios * and relative lengths. */ /* * Right now, the tolerance is set to 2 parts in a million, which * corresponds to a 19-bit distinction of mantissas. Note that * single-precsion numbers have 24 bit mantissas. */ #define DEGENERATE_TOLERANCE (0.000002) /* ========================================================== */ /* * The macro and subroutine INTERSECT are designed to compute the * intersection of a line (defined by the points v1 and v2) and a plane * (defined as plane which is normal to the vector n, and contains the * point p). Both return the point sect, which is the point of * interesection. * * This MACRO attemps to be fairly robust by checking for a divide by * zero. */ /* ========================================================== */ /* * HACK ALERT * The intersection parameter t has the nice property that if t>1, * then the intersection is "in front of" p1, and if t<0, then the * intersection is "behind" p2. Unfortunately, as the intersecting plane * and the line become parallel, t wraps through infinity -- i.e. t can * become so large that t becomes "greater than infinity" and comes back * as a negative number (i.e. winding number hopped by one unit). We * have no way of detecting this situation without adding gazzillions * of lines of code of topological algebra to detect the winding number; * and this would be incredibly difficult, and ruin performance. * * Thus, we've installed a cheap hack for use by the "cut style" drawing * routines. If t proves to be a large negative number (more negative * than -5), then we assume that t was positive and wound through * infinity. This makes most cuts look good, without introducing bogus * cuts at infinity. */ /* ========================================================== */ #define INTERSECT(sect,p,n,v1,v2) \ { \ gleDouble deno, numer, t, omt; \ \ deno = (v1[0] - v2[0]) * n[0]; \ deno += (v1[1] - v2[1]) * n[1]; \ deno += (v1[2] - v2[2]) * n[2]; \ \ if (deno == 0.0) { \ VEC_COPY (n, v1); \ /* printf ("Intersect: Warning: line is coplanar with plane \n"); */ \ } else { \ \ numer = (p[0] - v2[0]) * n[0]; \ numer += (p[1] - v2[1]) * n[1]; \ numer += (p[2] - v2[2]) * n[2]; \ \ t = numer / deno; \ omt = 1.0 - t; \ \ sect[0] = t * v1[0] + omt * v2[0]; \ sect[1] = t * v1[1] + omt * v2[1]; \ sect[2] = t * v1[2] + omt * v2[2]; \ } \ } /* ========================================================== */ /* * The macro and subroutine BISECTING_PLANE compute a normal vector that * describes the bisecting plane between three points (v1, v2 and v3). * This bisecting plane has the following properties: * 1) it contains the point v2 * 2) the angle it makes with v21 == v2 - v1 is equal to the angle it * makes with v32 == v3 - v2 * 3) it is perpendicular to the plane defined by v1, v2, v3. * * Having input v1, v2, and v3, it returns a unit vector n. * * In some cases, the user may specify degenerate points, and still * expect "reasonable" or "obvious" behaviour. The "expected" * behaviour for these degenerate cases is: * * 1) if v1 == v2 == v3, then return n=0 * 2) if v1 == v2, then return v32 (normalized). * 3) if v2 == v3, then return v21 (normalized). * 4) if v1, v2 and v3 co-linear, then return v21 (normalized). * * Mathematically, these special cases "make sense" -- we just have to * code around potential divide-by-zero's in the code below. */ /* ========================================================== */ #define BISECTING_PLANE(valid,n,v1,v2,v3) \ { \ double v21[3], v32[3]; \ double len21, len32; \ double dot; \ \ VEC_DIFF (v21, v2, v1); \ VEC_DIFF (v32, v3, v2); \ \ VEC_LENGTH (len21, v21); \ VEC_LENGTH (len32, v32); \ \ if (len21 <= DEGENERATE_TOLERANCE * len32) { \ \ if (len32 == 0.0) { \ /* all three points lie ontop of one-another */ \ VEC_ZERO (n); \ valid = FALSE; \ } else { \ /* return a normalized copy of v32 as bisector */ \ len32 = 1.0 / len32; \ VEC_SCALE (n, len32, v32); \ valid = TRUE; \ } \ \ } else { \ \ valid = TRUE; \ \ if (len32 <= DEGENERATE_TOLERANCE * len21) { \ /* return a normalized copy of v21 as bisector */ \ len21 = 1.0 / len21; \ VEC_SCALE (n, len21, v21); \ \ } else { \ \ /* normalize v21 to be of unit length */ \ len21 = 1.0 / len21; \ VEC_SCALE (v21, len21, v21); \ \ /* normalize v32 to be of unit length */ \ len32 = 1.0 / len32; \ VEC_SCALE (v32, len32, v32); \ \ VEC_DOT_PRODUCT (dot, v32, v21); \ \ /* if dot == 1 or -1, then points are colinear */ \ if ((dot >= (1.0-DEGENERATE_TOLERANCE)) || \ (dot <= (-1.0+DEGENERATE_TOLERANCE))) { \ VEC_COPY (n, v21); \ } else { \ \ /* go do the full computation */ \ n[0] = dot * (v32[0] + v21[0]) - v32[0] - v21[0]; \ n[1] = dot * (v32[1] + v21[1]) - v32[1] - v21[1]; \ n[2] = dot * (v32[2] + v21[2]) - v32[2] - v21[2]; \ \ /* if above if-test's passed, \ * n should NEVER be of zero length */ \ VEC_NORMALIZE (n); \ } \ } \ } \ } /* ========================================================== */ /* * The block of code below is ifdef'd out, and is here for reference * purposes only. It performs the "mathematically right thing" for * computing a bisecting plane, but is, unfortunately, subject ot noise * in the presence of near degenerate points. Since computer graphics, * due to sloppy coding, laziness, or correctness, is filled with * degenerate points, we can't really use this version. The code above * is far more appropriate for graphics. */ #ifdef MATHEMATICALLY_EXACT_GRAPHICALLY_A_KILLER #define BISECTING_PLANE(n,v1,v2,v3) \ { \ double v21[3], v32[3]; \ double len21, len32; \ double dot; \ \ VEC_DIFF (v21, v2, v1); \ VEC_DIFF (v32, v3, v2); \ \ VEC_LENGTH (len21, v21); \ VEC_LENGTH (len32, v32); \ \ if (len21 == 0.0) { \ \ if (len32 == 0.0) { \ /* all three points lie ontop of one-another */ \ VEC_ZERO (n); \ valid = FALSE; \ } else { \ /* return a normalized copy of v32 as bisector */ \ len32 = 1.0 / len32; \ VEC_SCALE (n, len32, v32); \ } \ \ } else { \ \ /* normalize v21 to be of unit length */ \ len21 = 1.0 / len21; \ VEC_SCALE (v21, len21, v21); \ \ if (len32 == 0.0) { \ /* return a normalized copy of v21 as bisector */ \ VEC_COPY (n, v21); \ } else { \ \ /* normalize v32 to be of unit length */ \ len32 = 1.0 / len32; \ VEC_SCALE (v32, len32, v32); \ \ VEC_DOT_PRODUCT (dot, v32, v21); \ \ /* if dot == 1 or -1, then points are colinear */ \ if ((dot == 1.0) || (dot == -1.0)) { \ VEC_COPY (n, v21); \ } else { \ \ /* go do the full computation */ \ n[0] = dot * (v32[0] + v21[0]) - v32[0] - v21[0]; \ n[1] = dot * (v32[1] + v21[1]) - v32[1] - v21[1]; \ n[2] = dot * (v32[2] + v21[2]) - v32[2] - v21[2]; \ \ /* if above if-test's passed, \ * n should NEVER be of zero length */ \ VEC_NORMALIZE (n); \ } \ } \ } \ } #endif /* ========================================================== */ /* * This macro computes the plane perpendicular to the the plane * defined by three points, and whose normal vector is givven as the * difference between the two vectors ... * * (See way below for the "math" model if you want to understand this. * The comments about relative errors above apply here.) */ #define CUTTING_PLANE(valid,n,v1,v2,v3) \ { \ double v21[3], v32[3]; \ double len21, len32; \ double lendiff; \ \ VEC_DIFF (v21, v2, v1); \ VEC_DIFF (v32, v3, v2); \ \ VEC_LENGTH (len21, v21); \ VEC_LENGTH (len32, v32); \ \ if (len21 <= DEGENERATE_TOLERANCE * len32) { \ \ if (len32 == 0.0) { \ /* all three points lie ontop of one-another */ \ VEC_ZERO (n); \ valid = FALSE; \ } else { \ /* return a normalized copy of v32 as cut-vector */ \ len32 = 1.0 / len32; \ VEC_SCALE (n, len32, v32); \ valid = TRUE; \ } \ \ } else { \ \ valid = TRUE; \ \ if (len32 <= DEGENERATE_TOLERANCE * len21) { \ /* return a normalized copy of v21 as cut vector */ \ len21 = 1.0 / len21; \ VEC_SCALE (n, len21, v21); \ } else { \ \ /* normalize v21 to be of unit length */ \ len21 = 1.0 / len21; \ VEC_SCALE (v21, len21, v21); \ \ /* normalize v32 to be of unit length */ \ len32 = 1.0 / len32; \ VEC_SCALE (v32, len32, v32); \ \ VEC_DIFF (n, v21, v32); \ VEC_LENGTH (lendiff, n); \ \ /* if the perp vector is very small, then the two \ * vectors are darn near collinear, and the cut \ * vector is probably poorly defined. */ \ if (lendiff < DEGENERATE_TOLERANCE) { \ VEC_ZERO (n); \ valid = FALSE; \ } else { \ lendiff = 1.0 / lendiff; \ VEC_SCALE (n, lendiff, n); \ } \ } \ } \ } /* ========================================================== */ #ifdef MATHEMATICALLY_EXACT_GRAPHICALLY_A_KILLER #define CUTTING_PLANE(n,v1,v2,v3) \ { \ double v21[3], v32[3]; \ \ VEC_DIFF (v21, v2, v1); \ VEC_DIFF (v32, v3, v2); \ \ VEC_NORMALIZE (v21); \ VEC_NORMALIZE (v32); \ \ VEC_DIFF (n, v21, v32); \ VEC_NORMALIZE (n); \ } #endif /* ============================================================ */ /* This macro is used in several places to cycle through a series of * points to find the next non-degenerate point in a series */ #define FIND_NON_DEGENERATE_POINT(inext,npoints,len,diff,point_array) \ { \ gleDouble slen; \ gleDouble summa[3]; \ \ do { \ /* get distance to next point */ \ VEC_DIFF (diff, point_array[inext+1], point_array[inext]); \ VEC_LENGTH (len, diff); \ VEC_SUM (summa, point_array[inext+1], point_array[inext]); \ VEC_LENGTH (slen, summa); \ slen *= DEGENERATE_TOLERANCE; \ inext ++; \ } while ((len <= slen) && (inext < npoints-1)); \ } /* ========================================================== */ extern int bisecting_plane (gleDouble n[3], /* returned */ gleDouble v1[3], /* input */ gleDouble v2[3], /* input */ gleDouble v3[3]); /* input */