svbsp.c

   1
   2 // Shadow Volume BSP code written by Forest "LordHavoc" Hale on 2003-11-06 and placed into public domain.
   3 // Modified by LordHavoc (to make it work and other nice things like that) on 2007-01-24 and 2007-01-25
   4 // Optimized by LordHavoc on 2009-12-24 and 2009-12-25
   5
   6 #include <math.h>
   7 #include <string.h>
   8 #include "svbsp.h"
   9 #include "polygon.h"
  10
  11 #define MAX_SVBSP_POLYGONPOINTS 64
  12 #define SVBSP_CLIP_EPSILON (1.0f / 1024.0f)
  13
  14 #define SVBSP_DotProduct(a,b) ((a)[0]*(b)[0]+(a)[1]*(b)[1]+(a)[2]*(b)[2])
  15
  16 typedef struct svbsp_polygon_s
  17 {
  18         float points[MAX_SVBSP_POLYGONPOINTS][3];
  19         //unsigned char splitflags[MAX_SVBSP_POLYGONPOINTS];
  20         int facesplitflag;
  21         int numpoints;
  22 }
  23 svbsp_polygon_t;
  24
  25 static void SVBSP_PlaneFromPoints(float *plane4f, const float *p1, const float *p2, const float *p3)
  26 {
  27         float ilength;
  28         // calculate unnormalized plane
  29         plane4f[0] = (p1[1] - p2[1]) * (p3[2] - p2[2]) - (p1[2] - p2[2]) * (p3[1] - p2[1]);
  30         plane4f[1] = (p1[2] - p2[2]) * (p3[0] - p2[0]) - (p1[0] - p2[0]) * (p3[2] - p2[2]);
  31         plane4f[2] = (p1[0] - p2[0]) * (p3[1] - p2[1]) - (p1[1] - p2[1]) * (p3[0] - p2[0]);
  32         plane4f[3] = SVBSP_DotProduct(plane4f, p1);
  33         // normalize the plane normal and adjust distance accordingly
  34         ilength = (float)sqrt(SVBSP_DotProduct(plane4f, plane4f));
  35         if (ilength)
  36                 ilength = 1.0f / ilength;
  37         plane4f[0] *= ilength;
  38         plane4f[1] *= ilength;
  39         plane4f[2] *= ilength;
  40         plane4f[3] *= ilength;
  41 }
  42
  43 static void SVBSP_DividePolygon(const svbsp_polygon_t *poly, const float *plane, svbsp_polygon_t *front, svbsp_polygon_t *back, const float *dists, const int *sides)
  44 {
  45         int i, j, count = poly->numpoints, frontcount = 0, backcount = 0;
  46         float frac, ifrac, c[3], pdist, ndist;
  47         const float *nextpoint;
  48         const float *points = poly->points[0];
  49         float *outfront = front->points[0];
  50         float *outback = back->points[0];
  51         for(i = 0;i < count;i++, points += 3)
  52         {
  53                 j = i + 1;
  54                 if (j >= count)
  55                         j = 0;
  56                 if (!(sides[i] & 2))
  57                 {
  58                         outfront[frontcount*3+0] = points[0];
  59                         outfront[frontcount*3+1] = points[1];
  60                         outfront[frontcount*3+2] = points[2];
  61                         frontcount++;
  62                 }
  63                 if (!(sides[i] & 1))
  64                 {
  65                         outback[backcount*3+0] = points[0];
  66                         outback[backcount*3+1] = points[1];
  67                         outback[backcount*3+2] = points[2];
  68                         backcount++;
  69                 }
  70                 if ((sides[i] | sides[j]) == 3)
  71                 {
  72                         // don't allow splits if remaining points would overflow point buffer
  73                         if (frontcount + (count - i) > MAX_SVBSP_POLYGONPOINTS - 1)
  74                                 continue;
  75                         if (backcount + (count - i) > MAX_SVBSP_POLYGONPOINTS - 1)
  76                                 continue;
  77                         nextpoint = poly->points[j];
  78                         pdist = dists[i];
  79                         ndist = dists[j];
  80                         frac = pdist / (pdist - ndist);
  81                         ifrac = 1.0f - frac;
  82                         c[0] = points[0] * ifrac + frac * nextpoint[0];
  83                         c[1] = points[1] * ifrac + frac * nextpoint[1];
  84                         c[2] = points[2] * ifrac + frac * nextpoint[2];
  85                         outfront[frontcount*3+0] = c[0];
  86                         outfront[frontcount*3+1] = c[1];
  87                         outfront[frontcount*3+2] = c[2];
  88                         frontcount++;
  89                         outback[backcount*3+0] = c[0];
  90                         outback[backcount*3+1] = c[1];
  91                         outback[backcount*3+2] = c[2];
  92                         backcount++;
  93                 }
  94         }
  95         front->numpoints = frontcount;
  96         back->numpoints = backcount;
  97 }
  98
  99 void SVBSP_Init(svbsp_t *b, const float *origin, int maxnodes, svbsp_node_t *nodes)
 100 {
 101         memset(b, 0, sizeof(*b));
 102         b->origin[0] = origin[0];
 103         b->origin[1] = origin[1];
 104         b->origin[2] = origin[2];
 105         b->numnodes = 3;
 106         b->maxnodes = maxnodes;
 107         b->nodes = nodes;
 108         b->ranoutofnodes = 0;
 109         b->stat_occluders_rejected = 0;
 110         b->stat_occluders_accepted = 0;
 111         b->stat_occluders_fragments_accepted = 0;
 112         b->stat_occluders_fragments_rejected = 0;
 113         b->stat_queries_rejected = 0;
 114         b->stat_queries_accepted = 0;
 115         b->stat_queries_fragments_accepted = 0;
 116         b->stat_queries_fragments_rejected = 0;
 117
 118         // the bsp tree must be initialized to have two perpendicular splits axes
 119         // through origin, otherwise the polygon insertions would affect the
 120         // opposite side of the tree, which would be disasterous.
 121         //
 122         // so this code has to make 3 nodes and 4 leafs, and since the leafs are
 123         // represented by empty/solid state numbers in this system rather than
 124         // actual structs, we only need to make the 3 nodes.
 125
 126         // root node
 127         // this one splits the world into +X and -X sides
 128         b->nodes[0].plane[0] = 1;
 129         b->nodes[0].plane[1] = 0;
 130         b->nodes[0].plane[2] = 0;
 131         b->nodes[0].plane[3] = origin[0];
 132         b->nodes[0].parent = -1;
 133         b->nodes[0].children[0] = 1;
 134         b->nodes[0].children[1] = 2;
 135
 136         // +X side node
 137         // this one splits the +X half of the world into +Y and -Y
 138         b->nodes[1].plane[0] = 0;
 139         b->nodes[1].plane[1] = 1;
 140         b->nodes[1].plane[2] = 0;
 141         b->nodes[1].plane[3] = origin[1];
 142         b->nodes[1].parent = 0;
 143         b->nodes[1].children[0] = -1;
 144         b->nodes[1].children[1] = -1;
 145
 146         // -X side node
 147         // this one splits the -X half of the world into +Y and -Y
 148         b->nodes[2].plane[0] = 0;
 149         b->nodes[2].plane[1] = 1;
 150         b->nodes[2].plane[2] = 0;
 151         b->nodes[2].plane[3] = origin[1];
 152         b->nodes[2].parent = 0;
 153         b->nodes[2].children[0] = -1;
 154         b->nodes[2].children[1] = -1;
 155 }
 156
 157 static void SVBSP_InsertOccluderPolygonNodes(svbsp_t *b, int *parentnodenumpointer, int parentnodenum, const svbsp_polygon_t *poly, void (*fragmentcallback)(void *fragmentcallback_pointer1, int fragmentcallback_number1, svbsp_t *b, int numpoints, const float *points), void *fragmentcallback_pointer1, int fragmentcallback_number1)
 158 {
 159         // now we need to create up to numpoints + 1 new nodes, forming a BSP tree
 160         // describing the occluder polygon's shadow volume
 161         int i, j, p;
 162         svbsp_node_t *node;
 163
 164         // points and lines are valid testers but not occluders
 165         if (poly->numpoints < 3)
 166                 return;
 167
 168         // if there aren't enough nodes remaining, skip it
 169         if (b->numnodes + poly->numpoints + 1 >= b->maxnodes)
 170         {
 171                 b->ranoutofnodes = 1;
 172                 return;
 173         }
 174
 175         // add one node per side, then the actual occluding face node
 176
 177         // thread safety notes:
 178         // DO NOT multithread insertion, it could be made 'safe' but the results
 179         // would be inconsistent.
 180         //
 181         // it is completely safe to multithread queries in all cases.
 182         //
 183         // if an insertion is occurring the query will give intermediate results,
 184         // being blocked by some volumes but not others, which is perfectly okay
 185         // for visibility culling intended only to reduce rendering work
 186
 187         // note down the first available nodenum for the *parentnodenumpointer
 188         // line which is done last to allow multithreaded queries during an
 189         // insertion
 190         for (i = 0, p = poly->numpoints - 1;i < poly->numpoints;p = i, i++)
 191         {
 192 #if 1
 193                 // see if a parent plane describes this side
 194                 for (j = parentnodenum;j >= 0;j = b->nodes[j].parent)
 195                 {
 196                         float *parentnodeplane = b->nodes[j].plane;
 197                         if (fabs(SVBSP_DotProduct(poly->points[p], parentnodeplane) - parentnodeplane[3]) < SVBSP_CLIP_EPSILON
 198                          && fabs(SVBSP_DotProduct(poly->points[i], parentnodeplane) - parentnodeplane[3]) < SVBSP_CLIP_EPSILON
 199                          && fabs(SVBSP_DotProduct(b->origin      , parentnodeplane) - parentnodeplane[3]) < SVBSP_CLIP_EPSILON)
 200                                 break;
 201                 }
 202                 if (j >= 0)
 203                         continue; // already have a matching parent plane
 204 #endif
 205 #if 0
 206                 // skip any sides that were classified as belonging to a parent plane
 207                 if (poly->splitflags[i])
 208                         continue;
 209 #endif
 210                 // create a side plane
 211                 // anything infront of this is not inside the shadow volume
 212                 node = b->nodes + b->numnodes++;
 213                 SVBSP_PlaneFromPoints(node->plane, b->origin, poly->points[p], poly->points[i]);
 214                 // we need to flip the plane if it puts any part of the polygon on the
 215                 // wrong side
 216                 // (in this way this code treats all polygons as float sided)
 217                 //
 218                 // because speed is important this stops as soon as it finds proof
 219                 // that the orientation is right or wrong
 220                 // (we know that the plane is on one edge of the polygon, so there is
 221                 // never a case where points lie on both sides, so the first hint is
 222                 // sufficient)
 223                 for (j = 0;j < poly->numpoints;j++)
 224                 {
 225                         float d = SVBSP_DotProduct(poly->points[j], node->plane) - node->plane[3];
 226                         if (d < -SVBSP_CLIP_EPSILON)
 227                                 break;
 228                         if (d > SVBSP_CLIP_EPSILON)
 229                         {
 230                                 node->plane[0] *= -1;
 231                                 node->plane[1] *= -1;
 232                                 node->plane[2] *= -1;
 233                                 node->plane[3] *= -1;
 234                                 break;
 235                         }
 236                 }
 237                 node->parent = parentnodenum;
 238                 node->children[0] = -1; // empty
 239                 node->children[1] = -1; // empty
 240                 // link this child into the tree
 241                 *parentnodenumpointer = parentnodenum = (int)(node - b->nodes);
 242                 // now point to the child pointer for the next node to update later
 243                 parentnodenumpointer = &node->children[1];
 244         }
 245
 246 #if 1
 247         // skip the face plane if it lies on a parent plane
 248         if (!poly->facesplitflag)
 249 #endif
 250         {
 251                 // add the face-plane node
 252                 // infront is empty, behind is shadow
 253                 node = b->nodes + b->numnodes++;
 254                 SVBSP_PlaneFromPoints(node->plane, poly->points[0], poly->points[1], poly->points[2]);
 255                 // this is a flip check similar to the one above
 256                 // this one checks if the plane faces the origin, if not, flip it
 257                 if (SVBSP_DotProduct(b->origin, node->plane) - node->plane[3] < -SVBSP_CLIP_EPSILON)
 258                 {
 259                         node->plane[0] *= -1;
 260                         node->plane[1] *= -1;
 261                         node->plane[2] *= -1;
 262                         node->plane[3] *= -1;
 263                 }
 264                 node->parent = parentnodenum;
 265                 node->children[0] = -1; // empty
 266                 node->children[1] = -2; // shadow
 267                 // link this child into the tree
 268                 // (with the addition of this node, queries will now be culled by it)
 269                 *parentnodenumpointer = (int)(node - b->nodes);
 270         }
 271 }
 272
 273 static int SVBSP_AddPolygonNode(svbsp_t *b, int *parentnodenumpointer, int parentnodenum, const svbsp_polygon_t *poly, int insertoccluder, void (*fragmentcallback)(void *fragmentcallback_pointer1, int fragmentcallback_number1, svbsp_t *b, int numpoints, const float *points), void *fragmentcallback_pointer1, int fragmentcallback_number1)
 274 {
 275         int i;
 276         int s;
 277         int facesplitflag = poly->facesplitflag;
 278         int bothsides;
 279         float plane[4];
 280         float d;
 281         svbsp_polygon_t front;
 282         svbsp_polygon_t back;
 283         svbsp_node_t *node;
 284         int sides[MAX_SVBSP_POLYGONPOINTS];
 285         float dists[MAX_SVBSP_POLYGONPOINTS];
 286         if (poly->numpoints < 1)
 287                 return 0;
 288         // recurse through plane nodes
 289         while (*parentnodenumpointer >= 0)
 290         {
 291                 // get node info
 292                 parentnodenum = *parentnodenumpointer;
 293                 node = b->nodes + parentnodenum;
 294                 plane[0] = node->plane[0];
 295                 plane[1] = node->plane[1];
 296                 plane[2] = node->plane[2];
 297                 plane[3] = node->plane[3];
 298                 // calculate point dists for clipping
 299                 bothsides = 0;
 300                 for (i = 0;i < poly->numpoints;i++)
 301                 {
 302                         d = SVBSP_DotProduct(poly->points[i], plane) - plane[3];
 303                         s = 0;
 304                         if (d > SVBSP_CLIP_EPSILON)
 305                                 s = 1;
 306                         if (d < -SVBSP_CLIP_EPSILON)
 307                                 s = 2;
 308                         bothsides |= s;
 309                         dists[i] = d;
 310                         sides[i] = s;
 311                 }
 312                 // see which side the polygon is on
 313                 switch(bothsides)
 314                 {
 315                 default:
 316                 case 0:
 317                         // no need to split, this polygon is on the plane
 318                         // this case only occurs for polygons on the face plane, usually
 319                         // the same polygon (inserted twice - once as occluder, once as
 320                         // tester)
 321                         // if this is an occluder, it is redundant
 322                         if (insertoccluder)
 323                                 return 1; // occluded
 324                         // if this is a tester, test the front side, because it is
 325                         // probably the same polygon that created this node...
 326                         facesplitflag = 1;
 327                         parentnodenumpointer = &node->children[0];
 328                         continue;
 329                 case 1:
 330                         // no need to split, just go to one side
 331                         parentnodenumpointer = &node->children[0];
 332                         continue;
 333                 case 2:
 334                         // no need to split, just go to one side
 335                         parentnodenumpointer = &node->children[1];
 336                         continue;
 337                 case 3:
 338                         // lies on both sides of the plane, we need to split it
 339 #if 1
 340                         SVBSP_DividePolygon(poly, plane, &front, &back, dists, sides);
 341 #else
 342                         PolygonF_Divide(poly->numpoints, poly->points[0], plane[0], plane[1], plane[2], plane[3], SVBSP_CLIP_EPSILON, MAX_SVBSP_POLYGONPOINTS, front.points[0], &front.numpoints, MAX_SVBSP_POLYGONPOINTS, back.points[0], &back.numpoints, NULL);
 343                         if (front.numpoints > MAX_SVBSP_POLYGONPOINTS)
 344                                 front.numpoints = MAX_SVBSP_POLYGONPOINTS;
 345                         if (back.numpoints > MAX_SVBSP_POLYGONPOINTS)
 346                                 back.numpoints = MAX_SVBSP_POLYGONPOINTS;
 347 #endif
 348                         front.facesplitflag = facesplitflag;
 349                         back.facesplitflag = facesplitflag;
 350                         // recurse the sides and return the resulting occlusion flags
 351                         i  = SVBSP_AddPolygonNode(b, &node->children[0], *parentnodenumpointer, &front, insertoccluder, fragmentcallback, fragmentcallback_pointer1, fragmentcallback_number1);
 352                         i |= SVBSP_AddPolygonNode(b, &node->children[1], *parentnodenumpointer, &back , insertoccluder, fragmentcallback, fragmentcallback_pointer1, fragmentcallback_number1);
 353                         return i;
 354                 }
 355         }
 356         // leaf node
 357         if (*parentnodenumpointer == -1)
 358         {
 359                 // empty leaf node; and some geometry survived
 360                 // if inserting an occluder, replace this empty leaf with a shadow volume
 361 #if 0
 362                 for (i = 0;i < poly->numpoints-2;i++)
 363                 {
 364                         Debug_PolygonBegin(NULL, DRAWFLAG_ADDITIVE);
 365                         Debug_PolygonVertex(poly->points[  0][0], poly->points[  0][1], poly->points[  0][2], 0.0f, 0.0f, 0.25f, 0.0f, 0.0f, 1.0f);
 366                         Debug_PolygonVertex(poly->points[i+1][0], poly->points[i+1][1], poly->points[i+1][2], 0.0f, 0.0f, 0.25f, 0.0f, 0.0f, 1.0f);
 367                         Debug_PolygonVertex(poly->points[i+2][0], poly->points[i+2][1], poly->points[i+2][2], 0.0f, 0.0f, 0.25f, 0.0f, 0.0f, 1.0f);
 368                         Debug_PolygonEnd();
 369                 }
 370 #endif
 371                 if (insertoccluder)
 372                 {
 373                         b->stat_occluders_fragments_accepted++;
 374                         SVBSP_InsertOccluderPolygonNodes(b, parentnodenumpointer, parentnodenum, poly, fragmentcallback, fragmentcallback_pointer1, fragmentcallback_number1);
 375                 }
 376                 else
 377                         b->stat_queries_fragments_accepted++;
 378                 if (fragmentcallback)
 379                         fragmentcallback(fragmentcallback_pointer1, fragmentcallback_number1, b, poly->numpoints, poly->points[0]);
 380                 return 2;
 381         }
 382         else
 383         {
 384                 // otherwise it's a solid leaf which destroys all polygons inside it
 385                 if (insertoccluder)
 386                         b->stat_occluders_fragments_rejected++;
 387                 else
 388                         b->stat_queries_fragments_rejected++;
 389 #if 0
 390                 for (i = 0;i < poly->numpoints-2;i++)
 391                 {
 392                         Debug_PolygonBegin(NULL, DRAWFLAG_ADDITIVE);
 393                         Debug_PolygonVertex(poly->points[  0][0], poly->points[  0][1], poly->points[  0][2], 0.0f, 0.0f, 0.0f, 0.0f, 0.25f, 1.0f);
 394                         Debug_PolygonVertex(poly->points[i+1][0], poly->points[i+1][1], poly->points[i+1][2], 0.0f, 0.0f, 0.0f, 0.0f, 0.25f, 1.0f);
 395                         Debug_PolygonVertex(poly->points[i+2][0], poly->points[i+2][1], poly->points[i+2][2], 0.0f, 0.0f, 0.0f, 0.0f, 0.25f, 1.0f);
 396                         Debug_PolygonEnd();
 397                 }
 398 #endif
 399         }
 400         return 1;
 401 }
 402
 403 int SVBSP_AddPolygon(svbsp_t *b, int numpoints, const float *points, int insertoccluder, void (*fragmentcallback)(void *fragmentcallback_pointer1, int fragmentcallback_number1, svbsp_t *b, int numpoints, const float *points), void *fragmentcallback_pointer1, int fragmentcallback_number1)
 404 {
 405         int i;
 406         int nodenum;
 407         svbsp_polygon_t poly;
 408         // don't even consider an empty polygon
 409         // note we still allow points and lines to be tested...
 410         if (numpoints < 1)
 411                 return 0;
 412         // if the polygon has too many points, we would crash
 413         if (numpoints > MAX_SVBSP_POLYGONPOINTS)
 414                 return 0;
 415         poly.numpoints = numpoints;
 416         for (i = 0;i < numpoints;i++)
 417         {
 418                 poly.points[i][0] = points[i*3+0];
 419                 poly.points[i][1] = points[i*3+1];
 420                 poly.points[i][2] = points[i*3+2];
 421                 //poly.splitflags[i] = 0; // this edge is a valid BSP splitter - clipped edges are not (because they lie on a bsp plane)
 422                 poly.facesplitflag = 0; // this face is a valid BSP Splitter - if it lies on a bsp plane it is not
 423         }
 424 #if 0
 425 //if (insertoccluder)
 426         for (i = 0;i < poly.numpoints-2;i++)
 427         {
 428                 Debug_PolygonBegin(NULL, DRAWFLAG_ADDITIVE);
 429                 Debug_PolygonVertex(poly.points[  0][0], poly.points[  0][1], poly.points[  0][2], 0.0f, 0.0f, 0.0f, 0.25f, 0.0f, 1.0f);
 430                 Debug_PolygonVertex(poly.points[i+1][0], poly.points[i+1][1], poly.points[i+1][2], 0.0f, 0.0f, 0.0f, 0.25f, 0.0f, 1.0f);
 431                 Debug_PolygonVertex(poly.points[i+2][0], poly.points[i+2][1], poly.points[i+2][2], 0.0f, 0.0f, 0.0f, 0.25f, 0.0f, 1.0f);
 432                 Debug_PolygonEnd();
 433         }
 434 #endif
 435         nodenum = 0;
 436         i = SVBSP_AddPolygonNode(b, &nodenum, -1, &poly, insertoccluder, fragmentcallback, fragmentcallback_pointer1, fragmentcallback_number1);
 437         if (insertoccluder)
 438         {
 439                 if (i & 2)
 440                         b->stat_occluders_accepted++;
 441                 else
 442                         b->stat_occluders_rejected++;
 443         }
 444         else
 445         {
 446                 if (i & 2)
 447                         b->stat_queries_accepted++;
 448                 else
 449                         b->stat_queries_rejected++;
 450         }
 451         return i;
 452 }
 453