tools/quake2/common/polylib.c

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