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

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