home *** CD-ROM | disk | FTP | other *** search
/ gondwana.ecr.mu.oz.au/pub/ / Graphics.tar / Graphics / atomart.tar.gz / atomart.tar / stree.c < prev    next >
C/C++ Source or Header  |  1990-06-21  |  12KB  |  551 lines
  1. #include <stdio.h>
  2. #include <math.h>
  3. #include "atomart.h"
  4. #include "gram.h"
  5. #include "macro.h"
  6.  
  7. extern object    *oblist;
  8. extern light    *lights;
  9. extern hlist    *fhlist;
  10. extern int    maxhitlevel;
  11.  
  12. extern pixel    backcol;
  13. extern colour    ambient;
  14.  
  15. extern double    power();
  16.  
  17. extern float    tolerance;
  18.  
  19. #define    U    0
  20. #define    V    1
  21.  
  22. extern object    *oblist;
  23. extern light    *lights;
  24. extern hlist    *fhlist;
  25. extern int    maxhitlevel;
  26.  
  27. extern pixel    backcol;
  28. extern colour    ambient;
  29.  
  30. extern double    power();
  31.  
  32. extern float    tolerance;
  33.  
  34. /*
  35.  * shrinkbbox
  36.  *
  37.  *    reduce the size of the bounding box for the leaf (if possible).
  38.  */
  39. shrinkbbox(leaf)
  40.     stree    *leaf;
  41. {
  42.     olist    *ol;
  43.     sphere    *s;
  44.     float    val;
  45.     float    minx, miny, minz;
  46.     float    maxx, maxy, maxz;
  47.  
  48.     if (leaf->u.leaf.oblist == (olist *)NULL)
  49.         return;
  50.  
  51.     s = leaf->u.leaf.oblist->obj->u.sph;
  52.  
  53.     minx = s->cent.x - s->rad;
  54.     miny = s->cent.y - s->rad;
  55.     minz = s->cent.z - s->rad;
  56.  
  57.     maxx = s->cent.x + s->rad;
  58.     maxy = s->cent.y + s->rad;
  59.     maxz = s->cent.z + s->rad;
  60.  
  61.     for (ol = leaf->u.leaf.oblist->nxt; ol != (olist *)NULL; ol = ol->nxt) {
  62.         s = ol->obj->u.sph;
  63.  
  64.         val = s->cent.x - s->rad;
  65.         if (val < minx)
  66.             minx = val;
  67.         val = s->cent.x + s->rad;
  68.         if (val > maxx)
  69.             maxx = val;
  70.  
  71.         val = s->cent.y - s->rad;
  72.         if (val < miny)
  73.             miny = val;
  74.         val = s->cent.y + s->rad;
  75.         if (val > maxy)
  76.             maxy = val;
  77.  
  78.         val = s->cent.z - s->rad;
  79.         if (val < minz)
  80.             minz = val;
  81.         val = s->cent.z + s->rad;
  82.         if (val > maxz)
  83.             maxz = val;
  84.     }
  85.  
  86.     if (minx > leaf->bb.min[X])
  87.         leaf->bb.min[X] = minx;
  88.     if (miny > leaf->bb.min[Y])
  89.         leaf->bb.min[Y] = miny;
  90.     if (minz > leaf->bb.min[Z])
  91.         leaf->bb.min[Z] = minz;
  92.  
  93.     if (maxx < leaf->bb.max[X])
  94.         leaf->bb.max[X] = maxx;
  95.     if (maxy < leaf->bb.max[Y])
  96.         leaf->bb.max[Y] = maxy;
  97.     if (maxz < leaf->bb.max[Z])
  98.         leaf->bb.max[Z] = maxz;
  99. }
  100.  
  101. /*
  102.  * splitleaf
  103.  *
  104.  *    split a leaf along the middle of an axis.
  105.  */
  106. splitleaf(axis, leaf, val)
  107.     int    axis;
  108.     stree    *leaf;
  109.     float    val;
  110. {
  111.     olist    *objs, *nxtobj, *tmpobj;
  112.     stree    *left, *right;
  113.     int    leftcount, rightcount;
  114.     float    min, max;
  115.     sphere    *s;
  116.  
  117.     if (leaf->type != SPLITABLELEAF)
  118.         return;
  119.  
  120.     left = (stree *)smalloc(sizeof(stree));
  121.     right = (stree *)smalloc(sizeof(stree));
  122.  
  123.     leftcount = rightcount = 0;
  124.  
  125.     left->u.leaf.oblist = right->u.leaf.oblist = (olist *)NULL;
  126.      
  127.     left->bb = right->bb = leaf->bb;
  128.     left->bb.max[axis] = val;
  129.     right->bb.min[axis] = val;
  130.                   
  131.     for (objs = leaf->u.leaf.oblist; objs != (olist *)NULL; objs = nxtobj) {
  132.         nxtobj = objs->nxt;
  133.         s = objs->obj->u.sph;
  134.         switch (axis) {
  135.         case X:
  136.             min = s->cent.x - s->rad;
  137.             max = s->cent.x + s->rad;
  138.             break;
  139.         case Y:
  140.             min = s->cent.y - s->rad;
  141.             max = s->cent.y + s->rad;
  142.             break;
  143.         case Z:
  144.             min = s->cent.z - s->rad;
  145.             max = s->cent.z + s->rad;
  146.             break;
  147.         default:
  148.             fatal("atomart: bad axis in splitleaf.\n");
  149.         }
  150.  
  151.         if (max < val) {            /* to the left */
  152.             leftcount++;
  153.             objs->nxt = left->u.leaf.oblist;
  154.             left->u.leaf.oblist = objs;
  155.         } else if (min > val) {            /* to the right */
  156.             rightcount++;
  157.             objs->nxt = right->u.leaf.oblist;
  158.             right->u.leaf.oblist = objs;
  159.         } else {                /* in both (sigh) */
  160.             tmpobj = (olist *)smalloc(sizeof(olist));
  161.             tmpobj->obj = objs->obj;
  162.  
  163.             leftcount++;
  164.             tmpobj->nxt = left->u.leaf.oblist;
  165.             left->u.leaf.oblist = tmpobj;
  166.  
  167.             rightcount++;
  168.             objs->nxt = right->u.leaf.oblist;
  169.             right->u.leaf.oblist = objs;
  170.         }
  171.     }
  172.  
  173.     shrinkbbox(left);
  174.     shrinkbbox(right);
  175.  
  176.     left->u.leaf.depth = right->u.leaf.depth = leaf->u.leaf.depth + 1;
  177.     left->type = (leftcount < MAXOBJS || leaf->u.leaf.depth == 13) ? LEAF : SPLITABLELEAF;
  178.     right->type = (rightcount < MAXOBJS || leaf->u.leaf.depth == 13) ? LEAF : SPLITABLELEAF;
  179.  
  180.     leaf->type = axis;
  181.     leaf->splitval = val;
  182.     leaf->u.branch.left = left;
  183.     leaf->u.branch.right = right;
  184.  
  185.     /*
  186.     printf("%f %f %f %f %f\n", leaf->bb.min[X], leaf->bb.max[X], leaf->bb.min[Y], leaf->bb.max[Y], leaf->bb.min[Z], leaf->bb.max[Z]);
  187.     printf("%d %d\n", leftcount, rightcount);
  188.     printf("%f %f %f %f %f\n", left->bb.min[X], left->bb.max[X], left->bb.min[Y], left->bb.max[Y], left->bb.min[Z], left->bb.max[Z]);
  189.     printf("%f %f %f %f %f\n", right->bb.min[X], right->bb.max[X], right->bb.min[Y], right->bb.max[Y], right->bb.min[Z], right->bb.max[Z]);
  190.     */
  191. }
  192.  
  193. /*
  194.  * split
  195.  *
  196.  *    split a splittable leaf into two nodes.
  197.  */
  198. void
  199. split(root)
  200.     stree    *root;
  201. {
  202.     int    total, left[3], right[3];
  203.     float    vals[3], min, max;
  204.     olist    *ol;
  205.     sphere    *s;
  206.  
  207.     total = 0;
  208.  
  209.     vals[X] = (root->bb.max[X] + root->bb.min[X]) / 2;
  210.     vals[Y] = (root->bb.max[Y] + root->bb.min[Y]) / 2;
  211.     vals[Z] = (root->bb.max[Z] + root->bb.min[Z]) / 2;
  212.  
  213.     left[X] = right[X] = 0;
  214.     left[Y] = right[Y] = 0;
  215.     left[Z] = right[Z] = 0;
  216.  
  217.     for (ol = root->u.leaf.oblist; ol != (olist *)NULL; ol = ol->nxt) {
  218.         s = ol->obj->u.sph;
  219.  
  220.         total++;
  221.  
  222.         min = s->cent.x - s->rad;    /* check for x */
  223.         max = s->cent.x + s->rad;
  224.         if (max < vals[X])        /* to the left */
  225.             left[X]++;
  226.         else if (min > vals[X])    /* to the right */
  227.             right[X]++;
  228.         else {            /* in both (sigh) */
  229.             left[X]++;
  230.             right[X]++;
  231.         }
  232.  
  233.         min = s->cent.y - s->rad;    /* check for y */
  234.         max = s->cent.y + s->rad;
  235.         if (max < vals[Y])        /* to the left */
  236.             left[Y]++;
  237.         else if (min > vals[Y])    /* to the right */
  238.             right[Y]++;
  239.         else {            /* in both (sigh) */
  240.             left[Y]++;
  241.             right[Y]++;
  242.         }
  243.  
  244.         min = s->cent.z - s->rad;    /* check for z */
  245.         max = s->cent.z + s->rad;
  246.         if (max < vals[Z])        /* to the left */
  247.             left[Z]++;
  248.         else if (min > vals[Z])    /* to the right */
  249.             right[Z]++;
  250.         else {            /* in both (sigh) */
  251.             left[Z]++;
  252.             right[Z]++;
  253.         }
  254.     }
  255.  
  256.     total += total / 3;
  257.  
  258.     if ((left[X] + right[X]) < (left[Y] + right[Y]))
  259.         if ((left[X] + right[X]) < (left[Z] + right[Z])) {
  260.             if (left[X] + right[X] <= total)
  261.                 splitleaf(X, root, vals[X]);
  262.             else
  263.                 root->type = LEAF;
  264.         } else {
  265.             if (left[Z] + right[Z] <= total)
  266.                 splitleaf(Z, root, vals[Z]);
  267.             else
  268.                 root->type = LEAF;
  269.         }
  270.     else
  271.         if ((left[Y] + right[Y]) < (left[Z] + right[Z])) {
  272.             if (left[Y] + right[Y] <= total)
  273.                 splitleaf(Y, root, vals[Y]);
  274.             else
  275.                 root->type = LEAF;
  276.         } else {
  277.             if (left[Z] + right[Z] <= total)
  278.                 splitleaf(Z, root, vals[Z]);
  279.             else
  280.                 root->type = LEAF;
  281.         }
  282. }
  283.  
  284. /*
  285.  * treei
  286.  *
  287.  *    returns the first hit on the bsp tree tr.
  288.  */
  289. hlist *
  290. treei(r, root)
  291.     register ray    *r;
  292.     stree        *root;
  293. {
  294.     register olist    *ol;
  295.     register hlist    *hit, *prevhits, *lp, *hp;
  296.     object        *o;
  297.     stree        *nxtbranch;
  298.     register int    p, u, v;
  299.     register float    val, t1, t2, othert1, othert2;
  300.     float        orgs[DIMS], dirs[DIMS];
  301.  
  302.     prevhits = (hlist *)NULL;
  303.  
  304.     hit = (hlist *)NULL;
  305.  
  306.     if (root->type == SPLITABLELEAF) {
  307.         split(root);
  308.         hit = treei(r, root);
  309.     } else if (root->type == LEAF) {
  310.  
  311.         for (ol = root->u.leaf.oblist; ol != (olist *)NULL; ol = ol->nxt) {
  312.             o = ol->obj;
  313.             if ((hit = o->intersects(r, o)) != (hlist *)NULL) {
  314.                 if (prevhits == (hlist *)NULL) {
  315.                     prevhits = hit;
  316.                 } else {
  317.                     if (prevhits->t > hit->t) {
  318.                         for (hp = prevhits; hp != (hlist *)NULL; hp = lp) {
  319.                             lp = hp->nxt;
  320.                             release(hp);
  321.                         }
  322.  
  323.                         prevhits = hit;
  324.                     } else {
  325.                         for (hp = hit; hp != (hlist *)NULL; hp = lp) {
  326.                             lp = hp->nxt;
  327.                             release(hp);
  328.                         }
  329.  
  330.                     }
  331.                 }
  332.             }
  333.         }
  334.  
  335.         if (prevhits != (hlist *)NULL) {
  336.             val = r->org.x + r->dir.x * prevhits->t;
  337.  
  338.             if (val > root->bb.max[X] || val < root->bb.min[X]) {
  339.                 for (hp = prevhits; hp != (hlist *)NULL; hp = lp) {
  340.                     lp = hp->nxt;
  341.                     release(hp);
  342.                 }
  343.                 return((hlist *)NULL);
  344.             }
  345.  
  346.             val = r->org.y + r->dir.y * prevhits->t;
  347.             if (val > root->bb.max[Y] || val < root->bb.min[Y]) {
  348.                 for (hp = prevhits; hp != (hlist *)NULL; hp = lp) {
  349.                     lp = hp->nxt;
  350.                     release(hp);
  351.                 }
  352.                 return((hlist *)NULL);
  353.             }
  354.  
  355.             val = r->org.z + r->dir.z * prevhits->t;
  356.             if (val > root->bb.max[Z] || val < root->bb.min[Z]) {
  357.                 for (hp = prevhits; hp != (hlist *)NULL; hp = lp) {
  358.                     lp = hp->nxt;
  359.                     release(hp);
  360.                 }
  361.                 return((hlist *)NULL);
  362.             }
  363.  
  364.         }
  365.  
  366.         return(prevhits);
  367.     } else {
  368.         orgs[X] = r->org.x;
  369.         orgs[Y] = r->org.y;
  370.         orgs[Z] = r->org.z;
  371.  
  372.         dirs[X] = r->dir.x;
  373.         dirs[Y] = r->dir.y;
  374.         dirs[Z] = r->dir.z;
  375.  
  376.         p = root->type;
  377.  
  378.         if (orgs[p] < root->splitval) {
  379.  
  380.             if (orgs[p] < root->bb.min[p] && dirs[p] < 0.0)
  381.                 return((hlist *)NULL);
  382.  
  383.             hit = treei(r, root->u.branch.left);
  384.  
  385.             if (hit != (hlist *)NULL) {
  386.                 return(hit);
  387.             } else if (dirs[p] <= 0.0)
  388.                 return((hlist *)NULL);
  389.  
  390.             nxtbranch = root->u.branch.right;
  391.             t1 = (nxtbranch->bb.min[p] - orgs[p]) / dirs[p];
  392.             t2 = (nxtbranch->bb.max[p] - orgs[p]) / dirs[p];
  393.  
  394.         } else {
  395.  
  396.             if (orgs[p] > root->bb.max[p] && dirs[p] > 0.0)
  397.                 return((hlist *)NULL);
  398.  
  399.             hit = treei(r, root->u.branch.right);
  400.  
  401.             if (hit != (hlist *)NULL) {
  402.                 return(hit);
  403.             } else if (dirs[p] >= 0.0)
  404.                 return((hlist *)NULL);
  405.  
  406.             nxtbranch = root->u.branch.left;
  407.             t1 = (nxtbranch->bb.max[p] - orgs[p]) / dirs[p];
  408.             t2 = (nxtbranch->bb.min[p] - orgs[p]) / dirs[p];
  409.         } 
  410.  
  411.         if (p == X) {
  412.             u = Y;
  413.             v = Z;
  414.         } else if (p == Y) {
  415.             u = Z;
  416.             v = X;
  417.         } else {
  418.             u = X;
  419.             v = Y;
  420.         }
  421.  
  422.         if (dirs[u] != 0.0) {
  423.             othert1 = (nxtbranch->bb.min[u] - orgs[u]) / dirs[u];
  424.             othert2 = (nxtbranch->bb.max[u] - orgs[u]) / dirs[u];
  425.             if (othert1 > othert2) {
  426.                 if (othert1 < t1 || othert2 > t2)
  427.                     return((hlist *)NULL);
  428.             } else {
  429.                 if (othert2 < t1 || othert1 > t2)
  430.                     return((hlist *)NULL);
  431.             }
  432.         } else {
  433.             if (nxtbranch->bb.min[u] > orgs[u])
  434.                 return((hlist *)NULL);
  435.             if (nxtbranch->bb.max[u] < orgs[u])
  436.                 return((hlist *)NULL);
  437.         }
  438.  
  439.         if (dirs[v] != 0.0) {
  440.             othert1 = (nxtbranch->bb.min[v] - orgs[v]) / dirs[v];
  441.             othert2 = (nxtbranch->bb.max[v] - orgs[v]) / dirs[v];
  442.             if (othert1 > othert2) {
  443.                 if (othert1 < t1 || othert2 > t2)
  444.                     return((hlist *)NULL);
  445.             } else {
  446.                 if (othert2 < t1 || othert1 > t2)
  447.                     return((hlist *)NULL);
  448.             }
  449.         } else {
  450.             if (nxtbranch->bb.min[v] > orgs[v])
  451.                 return((hlist *)NULL);
  452.             if (nxtbranch->bb.max[v] < orgs[v])
  453.                 return((hlist *)NULL);
  454.         }
  455.  
  456.         hit = treei(r, nxtbranch);
  457.     }
  458.  
  459.     return(hit);
  460. }
  461.  
  462. /*
  463.  * streei
  464.  *
  465.  *    returns the first hit on the bsp tree in o
  466.  */
  467. hlist *
  468. streei(r, o)
  469.     ray    *r;
  470.     object    *o;
  471. {
  472.     return(treei(r, o->u.tree));
  473. }
  474.  
  475. /*
  476.  * adjustbbox
  477.  *
  478.  *    adjust the bounding box for a tree for the object o
  479.  */
  480. adjustbbox(tree, o)
  481.     stree    *tree;
  482.     object    *o;
  483. {
  484.     register float    val;
  485.  
  486.     val = o->u.sph->cent.x - o->u.sph->rad;
  487.     if (val < tree->bb.min[X])
  488.         tree->bb.min[X] = val;
  489.     val = o->u.sph->cent.x + o->u.sph->rad;
  490.     if (val > tree->bb.max[X])
  491.         tree->bb.max[X] = val;
  492.  
  493.     val = o->u.sph->cent.y - o->u.sph->rad;
  494.     if (val < tree->bb.min[Y])
  495.         tree->bb.min[Y] = val;
  496.     val = o->u.sph->cent.y + o->u.sph->rad;
  497.     if (val > tree->bb.max[Y])
  498.         tree->bb.max[Y] = val;
  499.  
  500.     val = o->u.sph->cent.z - o->u.sph->rad;
  501.     if (val < tree->bb.min[Z])
  502.         tree->bb.min[Z] = val;
  503.     val = o->u.sph->cent.z + o->u.sph->rad;
  504.     if (val > tree->bb.max[Z])
  505.         tree->bb.max[Z] = val;
  506. }
  507.  
  508. /*
  509.  * treeinit
  510.  *
  511.  *    set up a bsp tree
  512.  */
  513. object *
  514. treeinit(obs)
  515.     object    *obs;
  516. {
  517.     object        *obj, *o;
  518.     olist        *ol;
  519.     stree        *tree;
  520.  
  521.     obj = (object *)smalloc(sizeof(object));
  522.  
  523.     obj->type = STREE;
  524.     obj->intersects = streei;
  525.     obj->u.tree = tree = (stree *)smalloc(sizeof(stree));
  526.  
  527.     tree->type = SPLITABLELEAF;
  528.  
  529.     tree->u.leaf.oblist = (olist *)NULL;
  530.     tree->u.leaf.depth = 0;
  531.  
  532.     tree->bb.min[X] = obs->u.sph->cent.x - obs->u.sph->rad;
  533.     tree->bb.max[X] = obs->u.sph->cent.x + obs->u.sph->rad;
  534.  
  535.     tree->bb.min[Y] = obs->u.sph->cent.y - obs->u.sph->rad;
  536.     tree->bb.max[Y] = obs->u.sph->cent.y + obs->u.sph->rad;
  537.  
  538.     tree->bb.min[Z] = obs->u.sph->cent.z - obs->u.sph->rad;
  539.     tree->bb.max[Z] = obs->u.sph->cent.z + obs->u.sph->rad;
  540.  
  541.     for (o = obs; o != (object *)NULL; o = o->nxt) {
  542.         ol = (olist *)smalloc(sizeof(olist));
  543.         ol->nxt = tree->u.leaf.oblist;
  544.         ol->obj = o;
  545.         adjustbbox(tree, o);
  546.         tree->u.leaf.oblist = ol;
  547.     }
  548.  
  549.     return(obj);
  550. }
  551.