polylib.c (14599B)
1 /* 2 =========================================================================== 3 Copyright (C) 1999-2005 Id Software, Inc. 4 5 This file is part of Quake III Arena source code. 6 7 Quake III Arena source code is free software; you can redistribute it 8 and/or modify it under the terms of the GNU General Public License as 9 published by the Free Software Foundation; either version 2 of the License, 10 or (at your option) any later version. 11 12 Quake III Arena source code is distributed in the hope that it will be 13 useful, but WITHOUT ANY WARRANTY; without even the implied warranty of 14 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 15 GNU General Public License for more details. 16 17 You should have received a copy of the GNU General Public License 18 along with Foobar; if not, write to the Free Software 19 Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA 20 =========================================================================== 21 */ 22 23 #include "cmdlib.h" 24 #include "mathlib.h" 25 #include "polylib.h" 26 #include "qfiles.h" 27 28 29 extern int numthreads; 30 31 // counters are only bumped when running single threaded, 32 // because they are an awefull coherence problem 33 int c_active_windings; 34 int c_peak_windings; 35 int c_winding_allocs; 36 int c_winding_points; 37 38 #define BOGUS_RANGE WORLD_SIZE 39 40 void pw(winding_t *w) 41 { 42 int i; 43 for (i=0 ; i<w->numpoints ; i++) 44 printf ("(%5.1f, %5.1f, %5.1f)\n",w->p[i][0], w->p[i][1],w->p[i][2]); 45 } 46 47 48 /* 49 ============= 50 AllocWinding 51 ============= 52 */ 53 winding_t *AllocWinding (int points) 54 { 55 winding_t *w; 56 int s; 57 58 if (numthreads == 1) 59 { 60 c_winding_allocs++; 61 c_winding_points += points; 62 c_active_windings++; 63 if (c_active_windings > c_peak_windings) 64 c_peak_windings = c_active_windings; 65 } 66 s = sizeof(vec_t)*3*points + sizeof(int); 67 w = malloc (s); 68 memset (w, 0, s); 69 return w; 70 } 71 72 void FreeWinding (winding_t *w) 73 { 74 if (*(unsigned *)w == 0xdeaddead) 75 Error ("FreeWinding: freed a freed winding"); 76 *(unsigned *)w = 0xdeaddead; 77 78 if (numthreads == 1) 79 c_active_windings--; 80 free (w); 81 } 82 83 /* 84 ============ 85 RemoveColinearPoints 86 ============ 87 */ 88 int c_removed; 89 90 void RemoveColinearPoints (winding_t *w) 91 { 92 int i, j, k; 93 vec3_t v1, v2; 94 int nump; 95 vec3_t p[MAX_POINTS_ON_WINDING]; 96 97 nump = 0; 98 for (i=0 ; i<w->numpoints ; i++) 99 { 100 j = (i+1)%w->numpoints; 101 k = (i+w->numpoints-1)%w->numpoints; 102 VectorSubtract (w->p[j], w->p[i], v1); 103 VectorSubtract (w->p[i], w->p[k], v2); 104 VectorNormalize(v1,v1); 105 VectorNormalize(v2,v2); 106 if (DotProduct(v1, v2) < 0.999) 107 { 108 VectorCopy (w->p[i], p[nump]); 109 nump++; 110 } 111 } 112 113 if (nump == w->numpoints) 114 return; 115 116 if (numthreads == 1) 117 c_removed += w->numpoints - nump; 118 w->numpoints = nump; 119 memcpy (w->p, p, nump*sizeof(p[0])); 120 } 121 122 /* 123 ============ 124 WindingPlane 125 ============ 126 */ 127 void WindingPlane (winding_t *w, vec3_t normal, vec_t *dist) 128 { 129 vec3_t v1, v2; 130 131 VectorSubtract (w->p[1], w->p[0], v1); 132 VectorSubtract (w->p[2], w->p[0], v2); 133 CrossProduct (v2, v1, normal); 134 VectorNormalize (normal, normal); 135 *dist = DotProduct (w->p[0], normal); 136 137 } 138 139 /* 140 ============= 141 WindingArea 142 ============= 143 */ 144 vec_t WindingArea (winding_t *w) 145 { 146 int i; 147 vec3_t d1, d2, cross; 148 vec_t total; 149 150 total = 0; 151 for (i=2 ; i<w->numpoints ; i++) 152 { 153 VectorSubtract (w->p[i-1], w->p[0], d1); 154 VectorSubtract (w->p[i], w->p[0], d2); 155 CrossProduct (d1, d2, cross); 156 total += 0.5 * VectorLength ( cross ); 157 } 158 return total; 159 } 160 161 void WindingBounds (winding_t *w, vec3_t mins, vec3_t maxs) 162 { 163 vec_t v; 164 int i,j; 165 166 mins[0] = mins[1] = mins[2] = 99999; 167 maxs[0] = maxs[1] = maxs[2] = -99999; 168 169 for (i=0 ; i<w->numpoints ; i++) 170 { 171 for (j=0 ; j<3 ; j++) 172 { 173 v = w->p[i][j]; 174 if (v < mins[j]) 175 mins[j] = v; 176 if (v > maxs[j]) 177 maxs[j] = v; 178 } 179 } 180 } 181 182 /* 183 ============= 184 WindingCenter 185 ============= 186 */ 187 void WindingCenter (winding_t *w, vec3_t center) 188 { 189 int i; 190 float scale; 191 192 VectorCopy (vec3_origin, center); 193 for (i=0 ; i<w->numpoints ; i++) 194 VectorAdd (w->p[i], center, center); 195 196 scale = 1.0/w->numpoints; 197 VectorScale (center, scale, center); 198 } 199 200 /* 201 ================= 202 BaseWindingForPlane 203 ================= 204 */ 205 winding_t *BaseWindingForPlane (vec3_t normal, vec_t dist) 206 { 207 int i, x; 208 vec_t max, v; 209 vec3_t org, vright, vup; 210 winding_t *w; 211 212 // find the major axis 213 214 max = -BOGUS_RANGE; 215 x = -1; 216 for (i=0 ; i<3; i++) 217 { 218 v = fabs(normal[i]); 219 if (v > max) 220 { 221 x = i; 222 max = v; 223 } 224 } 225 if (x==-1) 226 Error ("BaseWindingForPlane: no axis found"); 227 228 VectorCopy (vec3_origin, vup); 229 switch (x) 230 { 231 case 0: 232 case 1: 233 vup[2] = 1; 234 break; 235 case 2: 236 vup[0] = 1; 237 break; 238 } 239 240 v = DotProduct (vup, normal); 241 VectorMA (vup, -v, normal, vup); 242 VectorNormalize (vup, vup); 243 244 VectorScale (normal, dist, org); 245 246 CrossProduct (vup, normal, vright); 247 248 VectorScale (vup, MAX_WORLD_COORD, vup); 249 VectorScale (vright, MAX_WORLD_COORD, vright); 250 251 // project a really big axis aligned box onto the plane 252 w = AllocWinding (4); 253 254 VectorSubtract (org, vright, w->p[0]); 255 VectorAdd (w->p[0], vup, w->p[0]); 256 257 VectorAdd (org, vright, w->p[1]); 258 VectorAdd (w->p[1], vup, w->p[1]); 259 260 VectorAdd (org, vright, w->p[2]); 261 VectorSubtract (w->p[2], vup, w->p[2]); 262 263 VectorSubtract (org, vright, w->p[3]); 264 VectorSubtract (w->p[3], vup, w->p[3]); 265 266 w->numpoints = 4; 267 268 return w; 269 } 270 271 /* 272 ================== 273 CopyWinding 274 ================== 275 */ 276 winding_t *CopyWinding (winding_t *w) 277 { 278 int size; 279 winding_t *c; 280 281 c = AllocWinding (w->numpoints); 282 size = (int)((winding_t *)0)->p[w->numpoints]; 283 memcpy (c, w, size); 284 return c; 285 } 286 287 /* 288 ================== 289 ReverseWinding 290 ================== 291 */ 292 winding_t *ReverseWinding (winding_t *w) 293 { 294 int i; 295 winding_t *c; 296 297 c = AllocWinding (w->numpoints); 298 for (i=0 ; i<w->numpoints ; i++) 299 { 300 VectorCopy (w->p[w->numpoints-1-i], c->p[i]); 301 } 302 c->numpoints = w->numpoints; 303 return c; 304 } 305 306 307 /* 308 ============= 309 ClipWindingEpsilon 310 ============= 311 */ 312 void ClipWindingEpsilon (winding_t *in, vec3_t normal, vec_t dist, 313 vec_t epsilon, winding_t **front, winding_t **back) 314 { 315 vec_t dists[MAX_POINTS_ON_WINDING+4]; 316 int sides[MAX_POINTS_ON_WINDING+4]; 317 int counts[3]; 318 static vec_t dot; // VC 4.2 optimizer bug if not static 319 int i, j; 320 vec_t *p1, *p2; 321 vec3_t mid; 322 winding_t *f, *b; 323 int maxpts; 324 325 counts[0] = counts[1] = counts[2] = 0; 326 327 // determine sides for each point 328 for (i=0 ; i<in->numpoints ; i++) 329 { 330 dot = DotProduct (in->p[i], normal); 331 dot -= dist; 332 dists[i] = dot; 333 if (dot > epsilon) 334 sides[i] = SIDE_FRONT; 335 else if (dot < -epsilon) 336 sides[i] = SIDE_BACK; 337 else 338 { 339 sides[i] = SIDE_ON; 340 } 341 counts[sides[i]]++; 342 } 343 sides[i] = sides[0]; 344 dists[i] = dists[0]; 345 346 *front = *back = NULL; 347 348 if (!counts[0]) 349 { 350 *back = CopyWinding (in); 351 return; 352 } 353 if (!counts[1]) 354 { 355 *front = CopyWinding (in); 356 return; 357 } 358 359 maxpts = in->numpoints+4; // cant use counts[0]+2 because 360 // of fp grouping errors 361 362 *front = f = AllocWinding (maxpts); 363 *back = b = AllocWinding (maxpts); 364 365 for (i=0 ; i<in->numpoints ; i++) 366 { 367 p1 = in->p[i]; 368 369 if (sides[i] == SIDE_ON) 370 { 371 VectorCopy (p1, f->p[f->numpoints]); 372 f->numpoints++; 373 VectorCopy (p1, b->p[b->numpoints]); 374 b->numpoints++; 375 continue; 376 } 377 378 if (sides[i] == SIDE_FRONT) 379 { 380 VectorCopy (p1, f->p[f->numpoints]); 381 f->numpoints++; 382 } 383 if (sides[i] == SIDE_BACK) 384 { 385 VectorCopy (p1, b->p[b->numpoints]); 386 b->numpoints++; 387 } 388 389 if (sides[i+1] == SIDE_ON || sides[i+1] == sides[i]) 390 continue; 391 392 // generate a split point 393 p2 = in->p[(i+1)%in->numpoints]; 394 395 dot = dists[i] / (dists[i]-dists[i+1]); 396 for (j=0 ; j<3 ; j++) 397 { // avoid round off error when possible 398 if (normal[j] == 1) 399 mid[j] = dist; 400 else if (normal[j] == -1) 401 mid[j] = -dist; 402 else 403 mid[j] = p1[j] + dot*(p2[j]-p1[j]); 404 } 405 406 VectorCopy (mid, f->p[f->numpoints]); 407 f->numpoints++; 408 VectorCopy (mid, b->p[b->numpoints]); 409 b->numpoints++; 410 } 411 412 if (f->numpoints > maxpts || b->numpoints > maxpts) 413 Error ("ClipWinding: points exceeded estimate"); 414 if (f->numpoints > MAX_POINTS_ON_WINDING || b->numpoints > MAX_POINTS_ON_WINDING) 415 Error ("ClipWinding: MAX_POINTS_ON_WINDING"); 416 } 417 418 419 /* 420 ============= 421 ChopWindingInPlace 422 ============= 423 */ 424 void ChopWindingInPlace (winding_t **inout, vec3_t normal, vec_t dist, vec_t epsilon) 425 { 426 winding_t *in; 427 vec_t dists[MAX_POINTS_ON_WINDING+4]; 428 int sides[MAX_POINTS_ON_WINDING+4]; 429 int counts[3]; 430 static vec_t dot; // VC 4.2 optimizer bug if not static 431 int i, j; 432 vec_t *p1, *p2; 433 vec3_t mid; 434 winding_t *f; 435 int maxpts; 436 437 in = *inout; 438 counts[0] = counts[1] = counts[2] = 0; 439 440 // determine sides for each point 441 for (i=0 ; i<in->numpoints ; i++) 442 { 443 dot = DotProduct (in->p[i], normal); 444 dot -= dist; 445 dists[i] = dot; 446 if (dot > epsilon) 447 sides[i] = SIDE_FRONT; 448 else if (dot < -epsilon) 449 sides[i] = SIDE_BACK; 450 else 451 { 452 sides[i] = SIDE_ON; 453 } 454 counts[sides[i]]++; 455 } 456 sides[i] = sides[0]; 457 dists[i] = dists[0]; 458 459 if (!counts[0]) 460 { 461 FreeWinding (in); 462 *inout = NULL; 463 return; 464 } 465 if (!counts[1]) 466 return; // inout stays the same 467 468 maxpts = in->numpoints+4; // cant use counts[0]+2 because 469 // of fp grouping errors 470 471 f = AllocWinding (maxpts); 472 473 for (i=0 ; i<in->numpoints ; i++) 474 { 475 p1 = in->p[i]; 476 477 if (sides[i] == SIDE_ON) 478 { 479 VectorCopy (p1, f->p[f->numpoints]); 480 f->numpoints++; 481 continue; 482 } 483 484 if (sides[i] == SIDE_FRONT) 485 { 486 VectorCopy (p1, f->p[f->numpoints]); 487 f->numpoints++; 488 } 489 490 if (sides[i+1] == SIDE_ON || sides[i+1] == sides[i]) 491 continue; 492 493 // generate a split point 494 p2 = in->p[(i+1)%in->numpoints]; 495 496 dot = dists[i] / (dists[i]-dists[i+1]); 497 for (j=0 ; j<3 ; j++) 498 { // avoid round off error when possible 499 if (normal[j] == 1) 500 mid[j] = dist; 501 else if (normal[j] == -1) 502 mid[j] = -dist; 503 else 504 mid[j] = p1[j] + dot*(p2[j]-p1[j]); 505 } 506 507 VectorCopy (mid, f->p[f->numpoints]); 508 f->numpoints++; 509 } 510 511 if (f->numpoints > maxpts) 512 Error ("ClipWinding: points exceeded estimate"); 513 if (f->numpoints > MAX_POINTS_ON_WINDING) 514 Error ("ClipWinding: MAX_POINTS_ON_WINDING"); 515 516 FreeWinding (in); 517 *inout = f; 518 } 519 520 521 /* 522 ================= 523 ChopWinding 524 525 Returns the fragment of in that is on the front side 526 of the cliping plane. The original is freed. 527 ================= 528 */ 529 winding_t *ChopWinding (winding_t *in, vec3_t normal, vec_t dist) 530 { 531 winding_t *f, *b; 532 533 ClipWindingEpsilon (in, normal, dist, ON_EPSILON, &f, &b); 534 FreeWinding (in); 535 if (b) 536 FreeWinding (b); 537 return f; 538 } 539 540 541 /* 542 ================= 543 CheckWinding 544 545 ================= 546 */ 547 void CheckWinding (winding_t *w) 548 { 549 int i, j; 550 vec_t *p1, *p2; 551 vec_t d, edgedist; 552 vec3_t dir, edgenormal, facenormal; 553 vec_t area; 554 vec_t facedist; 555 556 if (w->numpoints < 3) 557 Error ("CheckWinding: %i points",w->numpoints); 558 559 area = WindingArea(w); 560 if (area < 1) 561 Error ("CheckWinding: %f area", area); 562 563 WindingPlane (w, facenormal, &facedist); 564 565 for (i=0 ; i<w->numpoints ; i++) 566 { 567 p1 = w->p[i]; 568 569 for (j=0 ; j<3 ; j++) 570 if (p1[j] > MAX_WORLD_COORD || p1[j] < MIN_WORLD_COORD) 571 Error ("CheckFace: BUGUS_RANGE: %f",p1[j]); 572 573 j = i+1 == w->numpoints ? 0 : i+1; 574 575 // check the point is on the face plane 576 d = DotProduct (p1, facenormal) - facedist; 577 if (d < -ON_EPSILON || d > ON_EPSILON) 578 Error ("CheckWinding: point off plane"); 579 580 // check the edge isnt degenerate 581 p2 = w->p[j]; 582 VectorSubtract (p2, p1, dir); 583 584 if (VectorLength (dir) < ON_EPSILON) 585 Error ("CheckWinding: degenerate edge"); 586 587 CrossProduct (facenormal, dir, edgenormal); 588 VectorNormalize (edgenormal, edgenormal); 589 edgedist = DotProduct (p1, edgenormal); 590 edgedist += ON_EPSILON; 591 592 // all other points must be on front side 593 for (j=0 ; j<w->numpoints ; j++) 594 { 595 if (j == i) 596 continue; 597 d = DotProduct (w->p[j], edgenormal); 598 if (d > edgedist) 599 Error ("CheckWinding: non-convex"); 600 } 601 } 602 } 603 604 605 /* 606 ============ 607 WindingOnPlaneSide 608 ============ 609 */ 610 int WindingOnPlaneSide (winding_t *w, vec3_t normal, vec_t dist) 611 { 612 qboolean front, back; 613 int i; 614 vec_t d; 615 616 front = qfalse; 617 back = qfalse; 618 for (i=0 ; i<w->numpoints ; i++) 619 { 620 d = DotProduct (w->p[i], normal) - dist; 621 if (d < -ON_EPSILON) 622 { 623 if (front) 624 return SIDE_CROSS; 625 back = qtrue; 626 continue; 627 } 628 if (d > ON_EPSILON) 629 { 630 if (back) 631 return SIDE_CROSS; 632 front = qtrue; 633 continue; 634 } 635 } 636 637 if (back) 638 return SIDE_BACK; 639 if (front) 640 return SIDE_FRONT; 641 return SIDE_ON; 642 } 643 644 645 /* 646 ================= 647 AddWindingToConvexHull 648 649 Both w and *hull are on the same plane 650 ================= 651 */ 652 #define MAX_HULL_POINTS 128 653 void AddWindingToConvexHull( winding_t *w, winding_t **hull, vec3_t normal ) { 654 int i, j, k; 655 float *p, *copy; 656 vec3_t dir; 657 float d; 658 int numHullPoints, numNew; 659 vec3_t hullPoints[MAX_HULL_POINTS]; 660 vec3_t newHullPoints[MAX_HULL_POINTS]; 661 vec3_t hullDirs[MAX_HULL_POINTS]; 662 qboolean hullSide[MAX_HULL_POINTS]; 663 qboolean outside; 664 665 if ( !*hull ) { 666 *hull = CopyWinding( w ); 667 return; 668 } 669 670 numHullPoints = (*hull)->numpoints; 671 memcpy( hullPoints, (*hull)->p, numHullPoints * sizeof(vec3_t) ); 672 673 for ( i = 0 ; i < w->numpoints ; i++ ) { 674 p = w->p[i]; 675 676 // calculate hull side vectors 677 for ( j = 0 ; j < numHullPoints ; j++ ) { 678 k = ( j + 1 ) % numHullPoints; 679 680 VectorSubtract( hullPoints[k], hullPoints[j], dir ); 681 VectorNormalize( dir, dir ); 682 CrossProduct( normal, dir, hullDirs[j] ); 683 } 684 685 outside = qfalse; 686 for ( j = 0 ; j < numHullPoints ; j++ ) { 687 VectorSubtract( p, hullPoints[j], dir ); 688 d = DotProduct( dir, hullDirs[j] ); 689 if ( d >= ON_EPSILON ) { 690 outside = qtrue; 691 } 692 if ( d >= -ON_EPSILON ) { 693 hullSide[j] = qtrue; 694 } else { 695 hullSide[j] = qfalse; 696 } 697 } 698 699 // if the point is effectively inside, do nothing 700 if ( !outside ) { 701 continue; 702 } 703 704 // find the back side to front side transition 705 for ( j = 0 ; j < numHullPoints ; j++ ) { 706 if ( !hullSide[ j % numHullPoints ] && hullSide[ (j + 1) % numHullPoints ] ) { 707 break; 708 } 709 } 710 if ( j == numHullPoints ) { 711 continue; 712 } 713 714 // insert the point here 715 VectorCopy( p, newHullPoints[0] ); 716 numNew = 1; 717 718 // copy over all points that aren't double fronts 719 j = (j+1)%numHullPoints; 720 for ( k = 0 ; k < numHullPoints ; k++ ) { 721 if ( hullSide[ (j+k) % numHullPoints ] && hullSide[ (j+k+1) % numHullPoints ] ) { 722 continue; 723 } 724 copy = hullPoints[ (j+k+1) % numHullPoints ]; 725 VectorCopy( copy, newHullPoints[numNew] ); 726 numNew++; 727 } 728 729 numHullPoints = numNew; 730 memcpy( hullPoints, newHullPoints, numHullPoints * sizeof(vec3_t) ); 731 } 732 733 FreeWinding( *hull ); 734 w = AllocWinding( numHullPoints ); 735 w->numpoints = numHullPoints; 736 *hull = w; 737 memcpy( w->p, hullPoints, numHullPoints * sizeof(vec3_t) ); 738 } 739 740