File indexing completed on 2024-03-24 15:17:33
0001 /* 0002 ** Author: Eric Veach, July 1994. 0003 ** 0004 */ 0005 0006 #include "gluos.h" 0007 #include "mesh.h" 0008 #include "tess.h" 0009 #include "normal.h" 0010 #include <math.h> 0011 #include <assert.h> 0012 0013 #ifndef TRUE 0014 #define TRUE 1 0015 #endif 0016 #ifndef FALSE 0017 #define FALSE 0 0018 #endif 0019 0020 #define Dot(u, v) (u[0] * v[0] + u[1] * v[1] + u[2] * v[2]) 0021 0022 #if 0 0023 static void Normalize( GLdouble v[3] ) 0024 { 0025 GLdouble len = v[0]*v[0] + v[1]*v[1] + v[2]*v[2]; 0026 0027 assert( len > 0 ); 0028 len = sqrt( len ); 0029 v[0] /= len; 0030 v[1] /= len; 0031 v[2] /= len; 0032 } 0033 #endif 0034 0035 #undef ABS 0036 #define ABS(x) ((x) < 0 ? -(x) : (x)) 0037 0038 static int LongAxis(GLdouble v[3]) 0039 { 0040 int i = 0; 0041 0042 if (ABS(v[1]) > ABS(v[0])) 0043 { 0044 i = 1; 0045 } 0046 if (ABS(v[2]) > ABS(v[i])) 0047 { 0048 i = 2; 0049 } 0050 return i; 0051 } 0052 0053 static void ComputeNormal(GLUtesselator *tess, GLdouble norm[3]) 0054 { 0055 GLUvertex *v, *v1, *v2; 0056 GLdouble c, tLen2, maxLen2; 0057 GLdouble maxVal[3], minVal[3], d1[3], d2[3], tNorm[3]; 0058 GLUvertex *maxVert[3], *minVert[3]; 0059 GLUvertex *vHead = &tess->mesh->vHead; 0060 int i; 0061 0062 maxVal[0] = maxVal[1] = maxVal[2] = -2 * GLU_TESS_MAX_COORD; 0063 minVal[0] = minVal[1] = minVal[2] = 2 * GLU_TESS_MAX_COORD; 0064 0065 for (v = vHead->next; v != vHead; v = v->next) 0066 { 0067 for (i = 0; i < 3; ++i) 0068 { 0069 c = v->coords[i]; 0070 if (c < minVal[i]) 0071 { 0072 minVal[i] = c; 0073 minVert[i] = v; 0074 } 0075 if (c > maxVal[i]) 0076 { 0077 maxVal[i] = c; 0078 maxVert[i] = v; 0079 } 0080 } 0081 } 0082 0083 /* Find two vertices separated by at least 1/sqrt(3) of the maximum 0084 * distance between any two vertices 0085 */ 0086 i = 0; 0087 if (maxVal[1] - minVal[1] > maxVal[0] - minVal[0]) 0088 { 0089 i = 1; 0090 } 0091 if (maxVal[2] - minVal[2] > maxVal[i] - minVal[i]) 0092 { 0093 i = 2; 0094 } 0095 if (minVal[i] >= maxVal[i]) 0096 { 0097 /* All vertices are the same -- normal doesn't matter */ 0098 norm[0] = 0; 0099 norm[1] = 0; 0100 norm[2] = 1; 0101 return; 0102 } 0103 0104 /* Look for a third vertex which forms the triangle with maximum area 0105 * (Length of normal == twice the triangle area) 0106 */ 0107 maxLen2 = 0; 0108 v1 = minVert[i]; 0109 v2 = maxVert[i]; 0110 d1[0] = v1->coords[0] - v2->coords[0]; 0111 d1[1] = v1->coords[1] - v2->coords[1]; 0112 d1[2] = v1->coords[2] - v2->coords[2]; 0113 for (v = vHead->next; v != vHead; v = v->next) 0114 { 0115 d2[0] = v->coords[0] - v2->coords[0]; 0116 d2[1] = v->coords[1] - v2->coords[1]; 0117 d2[2] = v->coords[2] - v2->coords[2]; 0118 tNorm[0] = d1[1] * d2[2] - d1[2] * d2[1]; 0119 tNorm[1] = d1[2] * d2[0] - d1[0] * d2[2]; 0120 tNorm[2] = d1[0] * d2[1] - d1[1] * d2[0]; 0121 tLen2 = tNorm[0] * tNorm[0] + tNorm[1] * tNorm[1] + tNorm[2] * tNorm[2]; 0122 if (tLen2 > maxLen2) 0123 { 0124 maxLen2 = tLen2; 0125 norm[0] = tNorm[0]; 0126 norm[1] = tNorm[1]; 0127 norm[2] = tNorm[2]; 0128 } 0129 } 0130 0131 if (maxLen2 <= 0) 0132 { 0133 /* All points lie on a single line -- any decent normal will do */ 0134 norm[0] = norm[1] = norm[2] = 0; 0135 norm[LongAxis(d1)] = 1; 0136 } 0137 } 0138 0139 static void CheckOrientation(GLUtesselator *tess) 0140 { 0141 GLdouble area; 0142 GLUface *f, *fHead = &tess->mesh->fHead; 0143 GLUvertex *v, *vHead = &tess->mesh->vHead; 0144 GLUhalfEdge *e; 0145 0146 /* When we compute the normal automatically, we choose the orientation 0147 * so that the sum of the signed areas of all contours is non-negative. 0148 */ 0149 area = 0; 0150 for (f = fHead->next; f != fHead; f = f->next) 0151 { 0152 e = f->anEdge; 0153 if (e->winding <= 0) 0154 continue; 0155 do 0156 { 0157 area += (e->Org->s - e->Dst->s) * (e->Org->t + e->Dst->t); 0158 e = e->Lnext; 0159 } while (e != f->anEdge); 0160 } 0161 if (area < 0) 0162 { 0163 /* Reverse the orientation by flipping all the t-coordinates */ 0164 for (v = vHead->next; v != vHead; v = v->next) 0165 { 0166 v->t = -v->t; 0167 } 0168 tess->tUnit[0] = -tess->tUnit[0]; 0169 tess->tUnit[1] = -tess->tUnit[1]; 0170 tess->tUnit[2] = -tess->tUnit[2]; 0171 } 0172 } 0173 0174 #ifdef FOR_TRITE_TEST_PROGRAM 0175 #include <stdlib.h> 0176 extern int RandomSweep; 0177 #define S_UNIT_X (RandomSweep ? (2 * drand48() - 1) : 1.0) 0178 #define S_UNIT_Y (RandomSweep ? (2 * drand48() - 1) : 0.0) 0179 #else 0180 #if defined(SLANTED_SWEEP) 0181 /* The "feature merging" is not intended to be complete. There are 0182 * special cases where edges are nearly parallel to the sweep line 0183 * which are not implemented. The algorithm should still behave 0184 * robustly (ie. produce a reasonable tesselation) in the presence 0185 * of such edges, however it may miss features which could have been 0186 * merged. We could minimize this effect by choosing the sweep line 0187 * direction to be something unusual (ie. not parallel to one of the 0188 * coordinate axes). 0189 */ 0190 #define S_UNIT_X 0.50941539564955385 /* Pre-normalized */ 0191 #define S_UNIT_Y 0.86052074622010633 0192 #else 0193 #define S_UNIT_X 1.0 0194 #define S_UNIT_Y 0.0 0195 #endif 0196 #endif 0197 0198 /* Determine the polygon normal and project vertices onto the plane 0199 * of the polygon. 0200 */ 0201 void __gl_projectPolygon(GLUtesselator *tess) 0202 { 0203 GLUvertex *v, *vHead = &tess->mesh->vHead; 0204 GLdouble norm[3]; 0205 GLdouble *sUnit, *tUnit; 0206 int i, computedNormal = FALSE; 0207 0208 norm[0] = tess->normal[0]; 0209 norm[1] = tess->normal[1]; 0210 norm[2] = tess->normal[2]; 0211 if (norm[0] == 0 && norm[1] == 0 && norm[2] == 0) 0212 { 0213 ComputeNormal(tess, norm); 0214 computedNormal = TRUE; 0215 } 0216 sUnit = tess->sUnit; 0217 tUnit = tess->tUnit; 0218 i = LongAxis(norm); 0219 0220 #if defined(FOR_TRITE_TEST_PROGRAM) || defined(TRUE_PROJECT) 0221 /* Choose the initial sUnit vector to be approximately perpendicular 0222 * to the normal. 0223 */ 0224 Normalize(norm); 0225 0226 sUnit[i] = 0; 0227 sUnit[(i + 1) % 3] = S_UNIT_X; 0228 sUnit[(i + 2) % 3] = S_UNIT_Y; 0229 0230 /* Now make it exactly perpendicular */ 0231 w = Dot(sUnit, norm); 0232 sUnit[0] -= w * norm[0]; 0233 sUnit[1] -= w * norm[1]; 0234 sUnit[2] -= w * norm[2]; 0235 Normalize(sUnit); 0236 0237 /* Choose tUnit so that (sUnit,tUnit,norm) form a right-handed frame */ 0238 tUnit[0] = norm[1] * sUnit[2] - norm[2] * sUnit[1]; 0239 tUnit[1] = norm[2] * sUnit[0] - norm[0] * sUnit[2]; 0240 tUnit[2] = norm[0] * sUnit[1] - norm[1] * sUnit[0]; 0241 Normalize(tUnit); 0242 #else 0243 /* Project perpendicular to a coordinate axis -- better numerically */ 0244 sUnit[i] = 0; 0245 sUnit[(i + 1) % 3] = S_UNIT_X; 0246 sUnit[(i + 2) % 3] = S_UNIT_Y; 0247 0248 tUnit[i] = 0; 0249 tUnit[(i + 1) % 3] = (norm[i] > 0) ? -S_UNIT_Y : S_UNIT_Y; 0250 tUnit[(i + 2) % 3] = (norm[i] > 0) ? S_UNIT_X : -S_UNIT_X; 0251 #endif 0252 0253 /* Project the vertices onto the sweep plane */ 0254 for (v = vHead->next; v != vHead; v = v->next) 0255 { 0256 v->s = Dot(v->coords, sUnit); 0257 v->t = Dot(v->coords, tUnit); 0258 } 0259 if (computedNormal) 0260 { 0261 CheckOrientation(tess); 0262 } 0263 }