tools/quake2/qdata_heretic2/common/polylib.c

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