source: GTP/trunk/Lib/Vis/Preprocessing/src/KdTree.cpp @ 863

Revision 863, 24.2 KB checked in by mattausch, 18 years ago (diff)

working on preprocessor integration
added iv stuff
Line 
1#include <stack>
2#include <algorithm>
3#include <queue>
4#include "Environment.h"
5#include "Mesh.h"
6#include "KdTree.h"
7#include "ViewCell.h"
8#include "Beam.h"
9
10
11
12namespace GtpVisibilityPreprocessor {
13
14int KdNode::mailID = 1;
15
16KdNode::KdNode(KdInterior *parent):mParent(parent), mailbox(0)
17{
18  if (parent)
19    mDepth = parent->mDepth+1;
20  else
21    mDepth = 0;
22}
23
24KdTree::KdTree()
25{
26  mRoot = new KdLeaf(NULL, 0);
27  environment->GetIntValue("KdTree.Termination.maxNodes", mTermMaxNodes);
28  environment->GetIntValue("KdTree.Termination.maxDepth", mTermMaxDepth);
29  environment->GetIntValue("KdTree.Termination.minCost", mTermMinCost);
30  environment->GetFloatValue("KdTree.Termination.maxCostRatio", mMaxCostRatio);
31  environment->GetFloatValue("KdTree.Termination.ct_div_ci", mCt_div_ci);
32  environment->GetFloatValue("KdTree.splitBorder", mSplitBorder);
33
34  environment->GetBoolValue("KdTree.sahUseFaces", mSahUseFaces);
35
36  char splitType[64];
37  environment->GetStringValue("KdTree.splitMethod", splitType);
38 
39  mSplitMethod = SPLIT_SPATIAL_MEDIAN;
40  if (strcmp(splitType, "spatialMedian") == 0)
41    mSplitMethod = SPLIT_SPATIAL_MEDIAN;
42  else
43    if (strcmp(splitType, "objectMedian") == 0)
44      mSplitMethod = SPLIT_OBJECT_MEDIAN;
45    else
46      if (strcmp(splitType, "SAH") == 0)
47        mSplitMethod = SPLIT_SAH;
48      else {
49        cerr<<"Wrong kd split type "<<splitType<<endl;
50        exit(1);
51      }
52  splitCandidates = NULL;
53}
54
55bool
56KdTree::Construct()
57{
58  if (!splitCandidates)
59    splitCandidates = new vector<SortableEntry>;
60
61  // first construct a leaf that will get subdivide
62  KdLeaf *leaf = (KdLeaf *) mRoot;
63
64  mStat.nodes = 1;
65
66  mBox.Initialize();
67  ObjectContainer::const_iterator mi;
68  for ( mi = leaf->mObjects.begin();
69        mi != leaf->mObjects.end();
70        mi++) {
71        //      cout<<(*mi)->GetBox()<<endl;
72        mBox.Include((*mi)->GetBox());
73  }
74
75  cout <<"KdTree Root Box:"<<mBox<<endl;
76  mRoot = Subdivide(TraversalData(leaf, mBox, 0));
77
78  // remove the allocated array
79  delete splitCandidates;
80
81  return true;
82}
83
84KdNode *
85KdTree::Subdivide(const TraversalData &tdata)
86{
87
88  KdNode *result = NULL;
89
90  priority_queue<TraversalData> tStack;
91  //  stack<STraversalData> tStack;
92 
93  tStack.push(tdata);
94  AxisAlignedBox3 backBox, frontBox;
95
96  while (!tStack.empty()) {
97        //      cout<<mStat.Nodes()<<" "<<mTermMaxNodes<<endl;
98        if (mStat.Nodes() > mTermMaxNodes) {
99          //    if ( GetMemUsage() > maxMemory ) {
100      // count statistics on unprocessed leafs
101      while (!tStack.empty()) {
102                EvaluateLeafStats(tStack.top());
103                tStack.pop();
104      }
105      break;
106    }
107         
108       
109    TraversalData data = tStack.top();
110    tStack.pop();
111   
112    KdNode *node = SubdivideNode((KdLeaf *) data.mNode,
113                                                                 data.mBox,
114                                                                 backBox,
115                                                                 frontBox
116                                                                 );
117
118        if (result == NULL)
119      result = node;
120   
121    if (!node->IsLeaf()) {
122
123      KdInterior *interior = (KdInterior *) node;
124      // push the children on the stack
125      tStack.push(TraversalData(interior->mBack, backBox, data.mDepth+1));
126      tStack.push(TraversalData(interior->mFront, frontBox, data.mDepth+1));
127     
128    } else {
129      EvaluateLeafStats(data);
130    }
131  }
132 
133  return result;
134
135}
136
137
138
139bool
140KdTree::TerminationCriteriaMet(const KdLeaf *leaf)
141{
142  //  cerr<<"\n OBJECTS="<<leaf->mObjects.size()<<endl;
143  return
144    (leaf->mObjects.size() <= mTermMinCost) ||
145    (leaf->mDepth >= mTermMaxDepth);
146 
147}
148
149
150int
151KdTree::SelectPlane(KdLeaf *leaf,
152                                        const AxisAlignedBox3 &box,
153                                        float &position
154                                        )
155{
156  int axis = -1;
157 
158  switch (mSplitMethod)
159    {
160    case SPLIT_SPATIAL_MEDIAN: {
161      axis = box.Size().DrivingAxis();
162      position = (box.Min()[axis] + box.Max()[axis])*0.5f;
163      break;
164    }
165    case SPLIT_SAH: {
166      int objectsBack, objectsFront;
167      float costRatio;
168      bool mOnlyDrivingAxis = false;
169      if (mOnlyDrivingAxis) {
170                axis = box.Size().DrivingAxis();
171                costRatio = BestCostRatio(leaf,
172                                                                  box,
173                                                                  axis,
174                                                                  position,
175                                                                  objectsBack,
176                                                                  objectsFront);
177      } else {
178                costRatio = MAX_FLOAT;
179                for (int i=0; i < 3; i++) {
180                  float p;
181                  float r = BestCostRatio(leaf,
182                                                                  box,
183                                                                  i,
184                                                                  p,
185                                                                  objectsBack,
186                                                                  objectsFront);
187                  if (r < costRatio) {
188                        costRatio = r;
189                        axis = i;
190                        position = p;
191                  }
192                }
193      }
194     
195      if (costRatio > mMaxCostRatio) {
196                //cout<<"Too big cost ratio "<<costRatio<<endl;
197                axis = -1;
198      }
199      break;
200    }
201   
202    }
203  return axis;
204}
205
206KdNode *
207KdTree::SubdivideNode(
208                      KdLeaf *leaf,
209                      const AxisAlignedBox3 &box,
210                      AxisAlignedBox3 &backBBox,
211                      AxisAlignedBox3 &frontBBox
212                      )
213{
214 
215  if (TerminationCriteriaMet(leaf))
216          return leaf;
217
218  float position;
219 
220  // select subdivision axis
221  int axis = SelectPlane( leaf, box, position );
222
223  if (axis == -1) {
224    return leaf;
225  }
226 
227  mStat.nodes+=2;
228  mStat.splits[axis]++;
229 
230  // add the new nodes to the tree
231  KdInterior *node = new KdInterior(leaf->mParent);
232
233  node->mAxis = axis;
234  node->mPosition = position;
235  node->mBox = box;
236 
237  backBBox = box;
238  frontBBox = box;
239 
240  // first count ray sides
241  int objectsBack = 0;
242  int objectsFront = 0;
243 
244  backBBox.SetMax(axis, position);
245  frontBBox.SetMin(axis, position);
246
247  ObjectContainer::const_iterator mi;
248
249  for ( mi = leaf->mObjects.begin();
250        mi != leaf->mObjects.end();
251        mi++) {
252    // determine the side of this ray with respect to the plane
253    AxisAlignedBox3 box = (*mi)->GetBox();
254    if (box.Max(axis) > position )
255      objectsFront++;
256   
257    if (box.Min(axis) < position )
258      objectsBack++;
259  }
260
261 
262  KdLeaf *back = new KdLeaf(node, objectsBack);
263  KdLeaf *front = new KdLeaf(node, objectsFront);
264
265  // replace a link from node's parent
266  if (  leaf->mParent )
267    leaf->mParent->ReplaceChildLink(leaf, node);
268
269  // and setup child links
270  node->SetupChildLinks(back, front);
271 
272  for (mi = leaf->mObjects.begin();
273       mi != leaf->mObjects.end();
274       mi++) {
275    // determine the side of this ray with respect to the plane
276    AxisAlignedBox3 box = (*mi)->GetBox();
277
278    if (box.Max(axis) >= position )
279      front->mObjects.push_back(*mi);
280   
281    if (box.Min(axis) < position )
282      back->mObjects.push_back(*mi);
283   
284    mStat.objectRefs -= (int)leaf->mObjects.size();
285    mStat.objectRefs += objectsBack + objectsFront;
286  }
287 
288  delete leaf;
289  return node;
290}
291
292
293
294void
295KdTreeStatistics::Print(ostream &app) const
296{
297  app << "===== KdTree statistics ===============\n";
298
299  app << "#N_NODES ( Number of nodes )\n" << nodes << "\n";
300
301  app << "#N_LEAVES ( Number of leaves )\n" << Leaves() << "\n";
302
303  app << "#N_SPLITS ( Number of splits in axes x y z dx dy dz)\n";
304  for (int i=0; i<7; i++)
305    app << splits[i] <<" ";
306  app <<endl;
307
308  app << "#N_RAYREFS ( Number of rayRefs )\n" <<
309    rayRefs << "\n";
310
311  app << "#N_RAYRAYREFS  ( Number of rayRefs / ray )\n" <<
312    rayRefs/(double)rays << "\n";
313
314  app << "#N_LEAFRAYREFS  ( Number of rayRefs / leaf )\n" <<
315    rayRefs/(double)Leaves() << "\n";
316
317  app << "#N_MAXOBJECTREFS  ( Max number of rayRefs / leaf )\n" <<
318    maxObjectRefs << "\n";
319
320  app << "#N_NONEMPTYRAYREFS  ( Number of rayRefs in nonEmpty leaves / non empty leaf )\n" <<
321    rayRefsNonZeroQuery/(double)(Leaves() - zeroQueryNodes) << "\n";
322
323  app << "#N_LEAFDOMAINREFS  ( Number of query domain Refs / leaf )\n" <<
324    objectRefs/(double)Leaves() << "\n";
325
326  //  app << setprecision(4);
327
328  app << "#N_PEMPTYLEAVES  ( Percentage of leaves with zero query domains )\n"<<
329    zeroQueryNodes*100/(double)Leaves()<<endl;
330
331  app << "#N_PMAXDEPTHLEAVES ( Percentage of leaves at maxdepth )\n"<<
332    maxDepthNodes*100/(double)Leaves()<<endl;
333
334  app << "#N_PMINCOSTLEAVES  ( Percentage of leaves with minCost )\n"<<
335    minCostNodes*100/(double)Leaves()<<endl;
336
337  app << "#N_ADDED_RAYREFS  (Number of dynamically added ray references )\n"<<
338    addedRayRefs<<endl;
339
340  app << "#N_REMOVED_RAYREFS  (Number of dynamically removed ray references )\n"<<
341    removedRayRefs<<endl;
342
343  //  app << setprecision(4);
344
345  //  app << "#N_CTIME  ( Construction time [s] )\n"
346  //      << Time() << " \n";
347
348  app << "===== END OF KdTree statistics ==========\n";
349
350}
351
352
353
354void
355KdTree::EvaluateLeafStats(const TraversalData &data)
356{
357
358  // the node became a leaf -> evaluate stats for leafs
359  KdLeaf *leaf = (KdLeaf *)data.mNode;
360
361  if (data.mDepth > mTermMaxDepth)
362    mStat.maxDepthNodes++;
363 
364  if ( (int)(leaf->mObjects.size()) < mTermMinCost)
365    mStat.minCostNodes++;
366 
367 
368  if ( (int)(leaf->mObjects.size()) > mStat.maxObjectRefs)
369    mStat.maxObjectRefs = (int)leaf->mObjects.size();
370 
371}
372
373
374
375void
376KdTree::SortSplitCandidates(
377                            KdLeaf *node,
378                            const int axis
379                            )
380{
381  splitCandidates->clear();
382 
383  int requestedSize = 2*(int)node->mObjects.size();
384  // creates a sorted split candidates array
385  if (splitCandidates->capacity() > 500000 &&
386      requestedSize < (int)(splitCandidates->capacity()/10) ) {
387    delete splitCandidates;
388    splitCandidates = new vector<SortableEntry>;
389  }
390 
391  splitCandidates->reserve(requestedSize);
392 
393  // insert all queries
394  for(ObjectContainer::const_iterator mi = node->mObjects.begin();
395      mi != node->mObjects.end();
396      mi++) {
397    AxisAlignedBox3 box = (*mi)->GetBox();
398
399    splitCandidates->push_back(SortableEntry(SortableEntry::BOX_MIN,
400                                                                                         box.Min(axis),
401                                                                                         *mi)
402                                                           );
403   
404   
405    splitCandidates->push_back(SortableEntry(SortableEntry::BOX_MAX,
406                                                                                         box.Max(axis),
407                                                                                         *mi)
408                                                           );
409  }
410 
411  stable_sort(splitCandidates->begin(), splitCandidates->end());
412}
413
414
415float
416KdTree::BestCostRatio(
417                                          KdLeaf *node,
418                                          const AxisAlignedBox3 &box,
419                                          const int axis,
420                                          float &position,
421                                          int &objectsBack,
422                                          int &objectsFront
423                                          )
424{
425
426  SortSplitCandidates(node, axis);
427 
428  // go through the lists, count the number of objects left and right
429  // and evaluate the following cost funcion:
430  // C = ct_div_ci  + (ol + or)/queries
431
432  float totalIntersections = 0.0f;
433  vector<SortableEntry>::const_iterator ci;
434
435  for(ci = splitCandidates->begin();
436      ci < splitCandidates->end();
437      ci++)
438    if ((*ci).type == SortableEntry::BOX_MIN) {
439      totalIntersections += (*ci).intersectable->IntersectionComplexity();
440    }
441       
442  float intersectionsLeft = 0;
443  float intersectionsRight = totalIntersections;
444       
445  int objectsLeft = 0, objectsRight = (int)node->mObjects.size();
446 
447  float minBox = box.Min(axis);
448  float maxBox = box.Max(axis);
449  float boxArea = box.SurfaceArea();
450 
451  float minBand = minBox + mSplitBorder*(maxBox - minBox);
452  float maxBand = minBox + (1.0f - mSplitBorder)*(maxBox - minBox);
453 
454  float minSum = 1e20f;
455 
456  for(ci = splitCandidates->begin();
457      ci < splitCandidates->end();
458      ci++) {
459    switch ((*ci).type) {
460    case SortableEntry::BOX_MIN:
461      objectsLeft++;
462      intersectionsLeft += (*ci).intersectable->IntersectionComplexity();
463      break;
464    case SortableEntry::BOX_MAX:
465      objectsRight--;
466      intersectionsRight -= (*ci).intersectable->IntersectionComplexity();
467      break;
468    }
469
470    if ((*ci).value > minBand && (*ci).value < maxBand) {
471      AxisAlignedBox3 lbox = box;
472      AxisAlignedBox3 rbox = box;
473      lbox.SetMax(axis, (*ci).value);
474      rbox.SetMin(axis, (*ci).value);
475
476      float sum;
477      if (mSahUseFaces)
478                sum = intersectionsLeft*lbox.SurfaceArea() + intersectionsRight*rbox.SurfaceArea();
479      else
480                sum = objectsLeft*lbox.SurfaceArea() + objectsRight*rbox.SurfaceArea();
481     
482      //      cout<<"pos="<<(*ci).value<<"\t q=("<<ql<<","<<qr<<")\t r=("<<rl<<","<<rr<<")"<<endl;
483      //      cout<<"cost= "<<sum<<endl;
484     
485      if (sum < minSum) {
486                minSum = sum;
487                position = (*ci).value;
488               
489                objectsBack = objectsLeft;
490                objectsFront = objectsRight;
491      }
492    }
493  }
494 
495  float oldCost = mSahUseFaces ? totalIntersections : node->mObjects.size();
496  float newCost = mCt_div_ci + minSum/boxArea;
497  float ratio = newCost/oldCost;
498 
499#if 0
500  cout<<"===================="<<endl;
501  cout<<"costRatio="<<ratio<<" pos="<<position<<" t="<<(position - minBox)/(maxBox - minBox)
502      <<"\t o=("<<objectsBack<<","<<objectsFront<<")"<<endl;
503#endif
504  return ratio;
505}
506
507int
508KdTree::CastRay(
509                                Ray &ray
510                                )
511{
512  int hits = 0;
513 
514  stack<RayTraversalData> tStack;
515 
516  float maxt = 1e6;
517  float mint = 0;
518 
519  Intersectable::NewMail();
520
521  if (!mBox.GetMinMaxT(ray, &mint, &maxt))
522    return 0;
523
524  if (mint < 0)
525    mint = 0;
526 
527  maxt += Limits::Threshold;
528 
529  Vector3 entp = ray.Extrap(mint);
530  Vector3 extp = ray.Extrap(maxt);
531 
532  KdNode *node = mRoot;
533  KdNode *farChild;
534  float position;
535  int axis;
536
537 
538  while (1) {
539    if (!node->IsLeaf()) {
540      KdInterior *in = (KdInterior *) node;
541      position = in->mPosition;
542      axis = in->mAxis;
543
544      if (entp[axis] <= position) {
545                if (extp[axis] <= position) {
546                  node = in->mBack;
547                  // cases N1,N2,N3,P5,Z2,Z3
548                  continue;
549                } else {
550                  // case N4
551                  node = in->mBack;
552                  farChild = in->mFront;
553                }
554      }
555      else {
556                if (position <= extp[axis]) {
557                  node = in->mFront;
558                  // cases P1,P2,P3,N5,Z1
559                  continue;
560                } else {
561                  node = in->mFront;
562                  farChild = in->mBack;
563                  // case P4
564                }
565          }
566      // $$ modification 3.5.2004 - hints from Kamil Ghais
567      // case N4 or P4
568      float tdist = (position - ray.GetLoc(axis)) / ray.GetDir(axis);
569      tStack.push(RayTraversalData(farChild, extp, maxt));
570      extp = ray.GetLoc() + ray.GetDir()*tdist;
571      maxt = tdist;
572        } else {
573          // compute intersection with all objects in this leaf
574          KdLeaf *leaf = (KdLeaf *) node;
575          if (ray.mFlags & Ray::STORE_KDLEAVES)
576                ray.kdLeaves.push_back(leaf);
577         
578          ObjectContainer::const_iterator mi;
579          for ( mi = leaf->mObjects.begin();
580                        mi != leaf->mObjects.end();
581                        mi++) {
582                Intersectable *object = *mi;
583                if (!object->Mailed() ) {
584                  object->Mail();
585                  if (ray.mFlags & Ray::STORE_TESTED_OBJECTS)
586                        ray.testedObjects.push_back(object);
587                  hits += object->CastRay(ray);
588                }
589          }
590         
591          if (hits && ray.GetType() == Ray::LOCAL_RAY)
592                if (ray.intersections[0].mT <= maxt)
593                  break;
594         
595          // get the next node from the stack
596          if (tStack.empty())
597                break;
598         
599          entp = extp;
600          mint = maxt;
601          if (ray.GetType() == Ray::LINE_SEGMENT && mint > 1.0f)
602                break;
603         
604          RayTraversalData &s  = tStack.top();
605          node = s.mNode;
606          extp = s.mExitPoint;
607          maxt = s.mMaxT;
608          tStack.pop();
609        }
610  }
611  return hits;
612}
613
614int KdTree::CastLineSegment(const Vector3 &origin,
615                                                        const Vector3 &termination,
616                                                        vector<ViewCell *> &viewcells)
617{
618        int hits = 0;
619
620        float mint = 0.0f, maxt = 1.0f;
621        const Vector3 dir = termination - origin;
622
623        stack<RayTraversalData> tStack;
624
625        Intersectable::NewMail();
626
627        //maxt += Limits::Threshold;
628
629        Vector3 entp = origin;
630        Vector3 extp = termination;
631
632        KdNode *node = mRoot;
633        KdNode *farChild;
634
635        float position;
636        int axis;
637
638        while (1)
639        {
640                if (!node->IsLeaf())
641                {
642                        KdInterior *in = dynamic_cast<KdInterior *>(node);
643                        position = in->mPosition;
644                        axis = in->mAxis;
645
646                        if (entp[axis] <= position)
647                        {
648                                if (extp[axis] <= position)
649                                {
650                                        node = in->mBack;
651                                        // cases N1,N2,N3,P5,Z2,Z3
652                                        continue;
653                                }
654                                else
655                                {
656                                        // case N4
657                                        node = in->mBack;
658                                        farChild = in->mFront;
659                                }
660                        }
661                        else
662                        {
663                                if (position <= extp[axis])
664                                {
665                                        node = in->mFront;
666                                        // cases P1,P2,P3,N5,Z1
667                                        continue;
668                                }
669                                else
670                                {
671                                        node = in->mFront;
672                                        farChild = in->mBack;
673                                        // case P4
674                                }
675                        }
676
677                        // $$ modification 3.5.2004 - hints from Kamil Ghais
678                        // case N4 or P4
679                        float tdist = (position - origin[axis]) / dir[axis];
680                        //tStack.push(RayTraversalData(farChild, extp, maxt)); //TODO
681                        extp = origin + dir * tdist;
682                        maxt = tdist;
683                }
684                else
685                {
686                        // compute intersection with all objects in this leaf
687                        KdLeaf *leaf = dynamic_cast<KdLeaf *>(node);
688
689                        // add view cell to intersections
690                        ViewCell *vc = leaf->mViewCell;
691
692                        if (!vc->Mailed())
693                        {
694                                vc->Mail();
695                                viewcells.push_back(vc);
696                                ++ hits;
697                        }
698
699                        // get the next node from the stack
700                        if (tStack.empty())
701                                break;
702
703                        entp = extp;
704                        mint = maxt;
705                       
706                        RayTraversalData &s  = tStack.top();
707                        node = s.mNode;
708                        extp = s.mExitPoint;
709                        maxt = s.mMaxT;
710                        tStack.pop();
711                }
712        }
713
714        return hits;
715}
716
717
718void
719KdTree::CollectObjects(const AxisAlignedBox3 &box,
720                                           ObjectContainer &objects)
721{
722  stack<KdNode *> nodeStack;
723
724  nodeStack.push(mRoot);
725
726  while (!nodeStack.empty()) {
727    KdNode *node = nodeStack.top();
728    nodeStack.pop();
729    if (node->IsLeaf()) {
730      KdLeaf *leaf = (KdLeaf *)node;
731      for (int j=0; j < leaf->mObjects.size(); j++) {
732                Intersectable *object = leaf->mObjects[j];
733                if (!object->Mailed()) {
734                  object->Mail();
735                  objects.push_back(object);
736                }
737      }
738    } else {
739      KdInterior *interior = (KdInterior *)node;
740
741          if ( box.Max()[interior->mAxis] > interior->mPosition )
742                nodeStack.push(interior->mFront);
743 
744          if (box.Min()[interior->mAxis] < interior->mPosition)
745                nodeStack.push(interior->mBack);
746    }
747  }
748}
749
750void
751KdTree::CollectObjects(KdNode *n, ObjectContainer &objects)
752{
753  stack<KdNode *> nodeStack;
754
755  nodeStack.push(n);
756
757  while (!nodeStack.empty()) {
758    KdNode *node = nodeStack.top();
759    nodeStack.pop();
760    if (node->IsLeaf()) {
761      KdLeaf *leaf = (KdLeaf *)node;
762      for (int j=0; j < leaf->mObjects.size(); j++) {
763                                Intersectable *object = leaf->mObjects[j];
764                                if (!object->Mailed()) {
765                                        object->Mail();
766                                        objects.push_back(object);
767                                }
768      }
769    } else {
770      KdInterior *interior = (KdInterior *)node;
771      nodeStack.push(interior->mFront);
772      nodeStack.push(interior->mBack);
773    }
774  }
775}
776
777// Find random neighbor which was not mailed
778KdNode *
779KdTree::FindRandomNeighbor(KdNode *n,
780                                                                                                         bool onlyUnmailed
781                                                                                                         )
782{
783  stack<KdNode *> nodeStack;
784 
785  nodeStack.push(mRoot);
786
787  AxisAlignedBox3 box = GetBox(n);
788  int mask = rand();
789
790  while (!nodeStack.empty()) {
791    KdNode *node = nodeStack.top();
792    nodeStack.pop();
793    if (node->IsLeaf()) {
794      if ( node != n && (!onlyUnmailed || !node->Mailed()) )
795        return node;
796    } else {
797      KdInterior *interior = (KdInterior *)node;
798      if (interior->mPosition > box.Max(interior->mAxis))
799        nodeStack.push(interior->mBack);
800      else
801        if (interior->mPosition < box.Min(interior->mAxis))
802          nodeStack.push(interior->mFront);
803        else {
804          // random decision
805          if (mask&1)
806            nodeStack.push(interior->mBack);
807          else
808            nodeStack.push(interior->mFront);
809          mask = mask>>1;
810        }
811    }
812  }
813 
814  return NULL;
815}
816
817int
818KdTree::FindNeighbors(KdNode *n,
819                      vector<KdNode *> &neighbors,
820                      bool onlyUnmailed
821                      )
822{
823  stack<KdNode *> nodeStack;
824 
825  nodeStack.push(mRoot);
826
827  AxisAlignedBox3 box = GetBox(n);
828
829  while (!nodeStack.empty()) {
830    KdNode *node = nodeStack.top();
831    nodeStack.pop();
832    if (node->IsLeaf()) {
833      if ( node != n && (!onlyUnmailed || !node->Mailed()) )
834        neighbors.push_back(node);
835    } else {
836      KdInterior *interior = (KdInterior *)node;
837      if (interior->mPosition > box.Max(interior->mAxis))
838                                nodeStack.push(interior->mBack);
839      else
840                                if (interior->mPosition < box.Min(interior->mAxis))
841                                        nodeStack.push(interior->mFront);
842                                else {
843                                        // random decision
844                                        nodeStack.push(interior->mBack);
845                                        nodeStack.push(interior->mFront);
846                                }
847    }
848  }
849 
850  return (int)neighbors.size();
851}
852
853// Find random neighbor which was not mailed
854KdNode *
855KdTree::GetRandomLeaf(const Plane3 &plane)
856{
857  stack<KdNode *> nodeStack;
858 
859  nodeStack.push(mRoot);
860 
861  int mask = rand();
862 
863  while (!nodeStack.empty()) {
864    KdNode *node = nodeStack.top();
865    nodeStack.pop();
866    if (node->IsLeaf()) {
867      return node;
868    } else {
869      KdInterior *interior = (KdInterior *)node;
870      KdNode *next;
871        if (GetBox(interior->mBack).Side(plane) < 0)
872          next = interior->mFront;
873        else
874          if (GetBox(interior->mFront).Side(plane) < 0)
875            next = interior->mBack;
876          else {
877            // random decision
878            if (mask&1)
879              next = interior->mBack;
880            else
881              next = interior->mFront;
882            mask = mask>>1;
883          }
884        nodeStack.push(next);
885    }
886  }
887 
888 
889  return NULL;
890}
891
892void
893KdTree::CollectLeaves(vector<KdLeaf *> &leaves)
894{
895  stack<KdNode *> nodeStack;
896  nodeStack.push(mRoot);
897
898  while (!nodeStack.empty()) {
899    KdNode *node = nodeStack.top();
900    nodeStack.pop();
901    if (node->IsLeaf()) {
902      KdLeaf *leaf = (KdLeaf *)node;
903      leaves.push_back(leaf);
904    } else {
905      KdInterior *interior = (KdInterior *)node;
906      nodeStack.push(interior->mBack);
907      nodeStack.push(interior->mFront);
908    }
909  }
910}
911
912void
913KdTree::CreateAndCollectViewCells(ViewCellContainer &vc) const
914{
915  stack<KdNode *> nodeStack;
916  nodeStack.push(mRoot);
917
918  while (!nodeStack.empty()) {
919    KdNode *node = nodeStack.top();
920    nodeStack.pop();
921    if (node->IsLeaf()) {
922      KdLeaf *leaf = (KdLeaf *)node;
923          // kdtree used as view cell container => create view cell
924          KdViewCell *viewCell = new KdViewCell();
925          leaf->mViewCell = viewCell;
926          // push back pointer to this leaf
927          viewCell->mLeaf = leaf;
928      vc.push_back(viewCell);
929    } else {
930      KdInterior *interior = (KdInterior *)node;
931      nodeStack.push(interior->mBack);
932      nodeStack.push(interior->mFront);
933    }
934  }
935}
936
937int
938KdTree::CollectLeafPvs()
939{
940  int totalPvsSize = 0;
941  stack<KdNode *> nodeStack;
942 
943  nodeStack.push(mRoot);
944 
945  while (!nodeStack.empty()) {
946    KdNode *node = nodeStack.top();
947    nodeStack.pop();
948    if (node->IsLeaf()) {
949      KdLeaf *leaf = (KdLeaf *)node;
950      for (int j=0; j < leaf->mObjects.size(); j++) {
951        Intersectable *object = leaf->mObjects[j];
952        if (!object->Mailed()) {
953          object->Mail();
954          // add this node to pvs of all nodes it can see
955          KdPvsMap::iterator ni = object->mKdPvs.mEntries.begin();
956          for (; ni != object->mKdPvs.mEntries.end(); ni++) {
957            KdNode *node = (*ni).first;
958            // $$ JB TEMPORARY solution -> should add object PVS or explictly computed
959            // kd tree PVS
960                float contribution;
961                if (leaf->mKdPvs.AddSample(node, 1.0f, contribution))
962              totalPvsSize++;
963          }
964        }
965      }
966    } else {
967      KdInterior *interior = (KdInterior *)node;
968      nodeStack.push(interior->mFront);
969      nodeStack.push(interior->mBack);
970    }
971  }
972
973  return totalPvsSize;
974}
975
976
977KdNode *
978KdTree::GetRandomLeaf(const bool onlyUnmailed)
979{
980  stack<KdNode *> nodeStack;
981  nodeStack.push(mRoot);
982 
983  int mask = rand();
984       
985  while (!nodeStack.empty()) {
986    KdNode *node = nodeStack.top();
987    nodeStack.pop();
988    if (node->IsLeaf()) {
989      if ( (!onlyUnmailed || !node->Mailed()) )
990                                return node;
991    } else {
992      KdInterior *interior = (KdInterior *)node;
993                        // random decision
994                        if (mask&1)
995                                nodeStack.push(interior->mBack);
996                        else
997                                nodeStack.push(interior->mFront);
998                        mask = mask>>1;
999                }
1000        }
1001  return NULL;
1002}
1003
1004
1005int
1006KdTree::CastBeam(
1007                                 Beam &beam
1008                                 )
1009{
1010  stack<KdNode *> nodeStack;
1011  nodeStack.push(mRoot);
1012 
1013  while (!nodeStack.empty()) {
1014    KdNode *node = nodeStack.top();
1015    nodeStack.pop();
1016       
1017        int side = beam.ComputeIntersection(GetBox(node));
1018        switch (side) {
1019        case -1:
1020          beam.mKdNodes.push_back(node);
1021          break;
1022        case 0:
1023          if (node->IsLeaf())
1024                beam.mKdNodes.push_back(node);
1025          else {
1026                KdInterior *interior = (KdInterior *)node;
1027                KdNode *first = interior->mBack;
1028                KdNode *second = interior->mFront;
1029               
1030                if (interior->mAxis < 3) {
1031                  // spatial split -> decide on the order of the nodes
1032                  if (beam.mPlanes[0].mNormal[interior->mAxis] > 0)
1033                        swap(first, second);
1034                }
1035
1036                nodeStack.push(first);
1037                nodeStack.push(second);
1038          }
1039          break;
1040          // default: cull
1041        }
1042  }
1043
1044  if (beam.mFlags & Beam::STORE_OBJECTS)
1045  {
1046          vector<KdNode *>::const_iterator it, it_end = beam.mKdNodes.end();
1047       
1048          Intersectable::NewMail();
1049          for (it = beam.mKdNodes.begin(); it != it_end; ++ it)
1050          {
1051                  CollectObjects(*it, beam.mObjects);
1052          }
1053  }
1054
1055  return (int)beam.mKdNodes.size();
1056}
1057
1058}
Note: See TracBrowser for help on using the repository browser.