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

Revision 1122, 24.8 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    if (box.Max(axis) >= position )
299      front->mObjects.push_back(*mi);
300   
301    if (box.Min(axis) < position )
302      back->mObjects.push_back(*mi);
303   
304    mStat.objectRefs -= (int)leaf->mObjects.size();
305    mStat.objectRefs += objectsBack + objectsFront;
306  }
307
308  delete leaf;
309  return node;
310}
311
312
313void
314KdTreeStatistics::Print(ostream &app) const
315{
316  app << "===== KdTree statistics ===============\n";
317
318  app << "#N_NODES ( Number of nodes )\n" << nodes << "\n";
319
320  app << "#N_LEAVES ( Number of leaves )\n" << Leaves() << "\n";
321
322  app << "#N_SPLITS ( Number of splits in axes x y z dx dy dz)\n";
323  for (int i=0; i<7; i++)
324    app << splits[i] <<" ";
325  app <<endl;
326
327  app << "#N_RAYREFS ( Number of rayRefs )\n" <<
328    rayRefs << "\n";
329
330  app << "#N_RAYRAYREFS  ( Number of rayRefs / ray )\n" <<
331    rayRefs/(double)rays << "\n";
332
333  app << "#N_LEAFRAYREFS  ( Number of rayRefs / leaf )\n" <<
334    rayRefs/(double)Leaves() << "\n";
335
336  app << "#N_MAXOBJECTREFS  ( Max number of rayRefs / leaf )\n" <<
337    maxObjectRefs << "\n";
338
339  app << "#N_NONEMPTYRAYREFS  ( Number of rayRefs in nonEmpty leaves / non empty leaf )\n" <<
340    rayRefsNonZeroQuery/(double)(Leaves() - zeroQueryNodes) << "\n";
341
342  app << "#N_LEAFDOMAINREFS  ( Number of query domain Refs / leaf )\n" <<
343    objectRefs/(double)Leaves() << "\n";
344
345  //  app << setprecision(4);
346
347  app << "#N_PEMPTYLEAVES  ( Percentage of leaves with zero query domains )\n"<<
348    zeroQueryNodes*100/(double)Leaves()<<endl;
349
350  app << "#N_PMAXDEPTHLEAVES ( Percentage of leaves at maxdepth )\n"<<
351    maxDepthNodes*100/(double)Leaves()<<endl;
352
353  app << "#N_PMINCOSTLEAVES  ( Percentage of leaves with minCost )\n"<<
354    minCostNodes*100/(double)Leaves()<<endl;
355
356  app << "#N_ADDED_RAYREFS  (Number of dynamically added ray references )\n"<<
357    addedRayRefs<<endl;
358
359  app << "#N_REMOVED_RAYREFS  (Number of dynamically removed ray references )\n"<<
360    removedRayRefs<<endl;
361
362  //  app << setprecision(4);
363
364  //  app << "#N_CTIME  ( Construction time [s] )\n"
365  //      << Time() << " \n";
366
367  app << "===== END OF KdTree statistics ==========\n";
368
369}
370
371
372
373void
374KdTree::EvaluateLeafStats(const TraversalData &data)
375{
376
377  // the node became a leaf -> evaluate stats for leafs
378  KdLeaf *leaf = (KdLeaf *)data.mNode;
379
380  if (data.mDepth > mTermMaxDepth)
381    mStat.maxDepthNodes++;
382 
383  if ( (int)(leaf->mObjects.size()) < mTermMinCost)
384    mStat.minCostNodes++;
385 
386 
387  if ( (int)(leaf->mObjects.size()) > mStat.maxObjectRefs)
388    mStat.maxObjectRefs = (int)leaf->mObjects.size();
389 
390}
391
392
393
394void
395KdTree::SortSplitCandidates(
396                            KdLeaf *node,
397                            const int axis
398                            )
399{
400  splitCandidates->clear();
401 
402  int requestedSize = 2*(int)node->mObjects.size();
403  // creates a sorted split candidates array
404  if (splitCandidates->capacity() > 500000 &&
405      requestedSize < (int)(splitCandidates->capacity()/10) ) {
406    delete splitCandidates;
407    splitCandidates = new vector<SortableEntry>;
408  }
409 
410  splitCandidates->reserve(requestedSize);
411 
412  // insert all queries
413  for(ObjectContainer::const_iterator mi = node->mObjects.begin();
414      mi != node->mObjects.end();
415      mi++) {
416    AxisAlignedBox3 box = (*mi)->GetBox();
417
418    splitCandidates->push_back(SortableEntry(SortableEntry::BOX_MIN,
419                                                                                         box.Min(axis),
420                                                                                         *mi)
421                                                           );
422   
423   
424    splitCandidates->push_back(SortableEntry(SortableEntry::BOX_MAX,
425                                                                                         box.Max(axis),
426                                                                                         *mi)
427                                                           );
428  }
429 
430  stable_sort(splitCandidates->begin(), splitCandidates->end());
431}
432
433
434float
435KdTree::BestCostRatio(
436                                          KdLeaf *node,
437                                          const AxisAlignedBox3 &box,
438                                          const int axis,
439                                          float &position,
440                                          int &objectsBack,
441                                          int &objectsFront
442                                          )
443{
444
445  SortSplitCandidates(node, axis);
446 
447  // go through the lists, count the number of objects left and right
448  // and evaluate the following cost funcion:
449  // C = ct_div_ci  + (ol + or)/queries
450
451  float totalIntersections = 0.0f;
452  vector<SortableEntry>::const_iterator ci;
453
454  for(ci = splitCandidates->begin();
455      ci < splitCandidates->end();
456      ci++)
457    if ((*ci).type == SortableEntry::BOX_MIN) {
458      totalIntersections += (*ci).intersectable->IntersectionComplexity();
459    }
460       
461  float intersectionsLeft = 0;
462  float intersectionsRight = totalIntersections;
463       
464  int objectsLeft = 0, objectsRight = (int)node->mObjects.size();
465 
466  float minBox = box.Min(axis);
467  float maxBox = box.Max(axis);
468  float boxArea = box.SurfaceArea();
469 
470  float minBand = minBox + mSplitBorder*(maxBox - minBox);
471  float maxBand = minBox + (1.0f - mSplitBorder)*(maxBox - minBox);
472 
473  float minSum = 1e20f;
474 
475  for(ci = splitCandidates->begin();
476      ci < splitCandidates->end();
477      ci++) {
478    switch ((*ci).type) {
479    case SortableEntry::BOX_MIN:
480      objectsLeft++;
481      intersectionsLeft += (*ci).intersectable->IntersectionComplexity();
482      break;
483    case SortableEntry::BOX_MAX:
484      objectsRight--;
485      intersectionsRight -= (*ci).intersectable->IntersectionComplexity();
486      break;
487    }
488
489    if ((*ci).value > minBand && (*ci).value < maxBand) {
490      AxisAlignedBox3 lbox = box;
491      AxisAlignedBox3 rbox = box;
492      lbox.SetMax(axis, (*ci).value);
493      rbox.SetMin(axis, (*ci).value);
494
495      float sum;
496      if (mSahUseFaces)
497                sum = intersectionsLeft*lbox.SurfaceArea() + intersectionsRight*rbox.SurfaceArea();
498      else
499                sum = objectsLeft*lbox.SurfaceArea() + objectsRight*rbox.SurfaceArea();
500     
501      //      cout<<"pos="<<(*ci).value<<"\t q=("<<ql<<","<<qr<<")\t r=("<<rl<<","<<rr<<")"<<endl;
502      //      cout<<"cost= "<<sum<<endl;
503     
504      if (sum < minSum) {
505                minSum = sum;
506                position = (*ci).value;
507               
508                objectsBack = objectsLeft;
509                objectsFront = objectsRight;
510      }
511    }
512  }
513 
514  float oldCost = mSahUseFaces ? totalIntersections : node->mObjects.size();
515  float newCost = mCt_div_ci + minSum/boxArea;
516  float ratio = newCost/oldCost;
517 
518#if 0
519  cout<<"===================="<<endl;
520  cout<<"costRatio="<<ratio<<" pos="<<position<<" t="<<(position - minBox)/(maxBox - minBox)
521      <<"\t o=("<<objectsBack<<","<<objectsFront<<")"<<endl;
522#endif
523  return ratio;
524}
525
526int
527KdTree::CastRay(
528                                Ray &ray
529                                )
530{
531  int hits = 0;
532 
533  stack<RayTraversalData> tStack;
534 
535  float maxt = 1e6;
536  float mint = 0;
537 
538  Intersectable::NewMail();
539
540  if (!mBox.GetMinMaxT(ray, &mint, &maxt))
541    return 0;
542
543  if (mint < 0)
544    mint = 0;
545 
546  maxt += Limits::Threshold;
547 
548  Vector3 entp = ray.Extrap(mint);
549  Vector3 extp = ray.Extrap(maxt);
550 
551  KdNode *node = mRoot;
552  KdNode *farChild;
553  float position;
554  int axis;
555
556 
557  while (1) {
558    if (!node->IsLeaf()) {
559      KdInterior *in = (KdInterior *) node;
560      position = in->mPosition;
561      axis = in->mAxis;
562
563      if (entp[axis] <= position) {
564                if (extp[axis] <= position) {
565                  node = in->mBack;
566                  // cases N1,N2,N3,P5,Z2,Z3
567                  continue;
568                } else {
569                  // case N4
570                  node = in->mBack;
571                  farChild = in->mFront;
572                }
573      }
574      else {
575                if (position <= extp[axis]) {
576                  node = in->mFront;
577                  // cases P1,P2,P3,N5,Z1
578                  continue;
579                } else {
580                  node = in->mFront;
581                  farChild = in->mBack;
582                  // case P4
583                }
584          }
585      // $$ modification 3.5.2004 - hints from Kamil Ghais
586      // case N4 or P4
587      float tdist = (position - ray.GetLoc(axis)) / ray.GetDir(axis);
588      tStack.push(RayTraversalData(farChild, extp, maxt));
589      extp = ray.GetLoc() + ray.GetDir()*tdist;
590      maxt = tdist;
591        } else {
592          // compute intersection with all objects in this leaf
593          KdLeaf *leaf = (KdLeaf *) node;
594          if (ray.mFlags & Ray::STORE_KDLEAVES)
595                ray.kdLeaves.push_back(leaf);
596         
597          ObjectContainer::const_iterator mi;
598          for ( mi = leaf->mObjects.begin();
599                        mi != leaf->mObjects.end();
600                        mi++) {
601                Intersectable *object = *mi;
602                if (!object->Mailed() ) {
603                  object->Mail();
604                  if (ray.mFlags & Ray::STORE_TESTED_OBJECTS)
605                        ray.testedObjects.push_back(object);
606
607                  static int oi=1;
608                  if (MeshDebug)
609                        cout<<"Object "<<oi++;
610                 
611                  hits += object->CastRay(ray);
612
613                  if (MeshDebug) {
614                        if (ray.intersections.size())
615                        cout<<"nearest t="<<ray.intersections[0].mT<<endl;
616                  else
617                        cout<<"nearest t=-INF"<<endl;
618                  }
619                 
620                }
621          }
622         
623          if (hits && ray.GetType() == Ray::LOCAL_RAY)
624                if (ray.intersections[0].mT <= maxt)
625                  break;
626         
627          // get the next node from the stack
628          if (tStack.empty())
629                break;
630         
631          entp = extp;
632          mint = maxt;
633          if (ray.GetType() == Ray::LINE_SEGMENT && mint > 1.0f)
634                break;
635         
636          RayTraversalData &s  = tStack.top();
637          node = s.mNode;
638          extp = s.mExitPoint;
639          maxt = s.mMaxT;
640          tStack.pop();
641        }
642  }
643  return hits;
644}
645
646int KdTree::CastLineSegment(const Vector3 &origin,
647                                                        const Vector3 &termination,
648                                                        ViewCellContainer &viewcells)
649{
650        int hits = 0;
651
652        float mint = 0.0f, maxt = 1.0f;
653        const Vector3 dir = termination - origin;
654
655        stack<RayTraversalData> tStack;
656
657        Intersectable::NewMail();
658
659        //maxt += Limits::Threshold;
660
661        Vector3 entp = origin;
662        Vector3 extp = termination;
663
664        KdNode *node = mRoot;
665        KdNode *farChild;
666
667        float position;
668        int axis;
669
670        while (1)
671        {
672                if (!node->IsLeaf())
673                {
674                        KdInterior *in = dynamic_cast<KdInterior *>(node);
675                        position = in->mPosition;
676                        axis = in->mAxis;
677
678                        if (entp[axis] <= position)
679                        {
680                                if (extp[axis] <= position)
681                                {
682                                        node = in->mBack;
683                                        // cases N1,N2,N3,P5,Z2,Z3
684                                        continue;
685                                }
686                                else
687                                {
688                                        // case N4
689                                        node = in->mBack;
690                                        farChild = in->mFront;
691                                }
692                        }
693                        else
694                        {
695                                if (position <= extp[axis])
696                                {
697                                        node = in->mFront;
698                                        // cases P1,P2,P3,N5,Z1
699                                        continue;
700                                }
701                                else
702                                {
703                                        node = in->mFront;
704                                        farChild = in->mBack;
705                                        // case P4
706                                }
707                        }
708
709                        // $$ modification 3.5.2004 - hints from Kamil Ghais
710                        // case N4 or P4
711                        float tdist = (position - origin[axis]) / dir[axis];
712                        //tStack.push(RayTraversalData(farChild, extp, maxt)); //TODO
713                        extp = origin + dir * tdist;
714                        maxt = tdist;
715                }
716                else
717                {
718                        // compute intersection with all objects in this leaf
719                        KdLeaf *leaf = dynamic_cast<KdLeaf *>(node);
720
721                        // add view cell to intersections
722                        ViewCell *vc = leaf->mViewCell;
723
724                        if (!vc->Mailed())
725                        {
726                                vc->Mail();
727                                viewcells.push_back(vc);
728                                ++ hits;
729                        }
730
731                        // get the next node from the stack
732                        if (tStack.empty())
733                                break;
734
735                        entp = extp;
736                        mint = maxt;
737                       
738                        RayTraversalData &s  = tStack.top();
739                        node = s.mNode;
740                        extp = s.mExitPoint;
741                        maxt = s.mMaxT;
742                        tStack.pop();
743                }
744        }
745
746        return hits;
747}
748
749
750void
751KdTree::CollectObjects(const AxisAlignedBox3 &box,
752                                           ObjectContainer &objects)
753{
754  stack<KdNode *> nodeStack;
755
756  nodeStack.push(mRoot);
757
758  while (!nodeStack.empty()) {
759    KdNode *node = nodeStack.top();
760    nodeStack.pop();
761    if (node->IsLeaf()) {
762      KdLeaf *leaf = (KdLeaf *)node;
763      for (int j=0; j < leaf->mObjects.size(); j++) {
764                Intersectable *object = leaf->mObjects[j];
765                if (!object->Mailed() && Overlap(box, object->GetBox())) {
766                  object->Mail();
767                  objects.push_back(object);
768                }
769      }
770    } else {
771      KdInterior *interior = (KdInterior *)node;
772
773          if ( box.Max()[interior->mAxis] > interior->mPosition )
774                nodeStack.push(interior->mFront);
775 
776          if (box.Min()[interior->mAxis] < interior->mPosition)
777                nodeStack.push(interior->mBack);
778    }
779  }
780}
781
782void
783KdTree::CollectObjects(KdNode *n, ObjectContainer &objects)
784{
785  stack<KdNode *> nodeStack;
786
787  nodeStack.push(n);
788
789  while (!nodeStack.empty()) {
790    KdNode *node = nodeStack.top();
791    nodeStack.pop();
792    if (node->IsLeaf()) {
793      KdLeaf *leaf = (KdLeaf *)node;
794      for (int j=0; j < leaf->mObjects.size(); j++) {
795                                Intersectable *object = leaf->mObjects[j];
796                                if (!object->Mailed()) {
797                                        object->Mail();
798                                        objects.push_back(object);
799                                }
800      }
801    } else {
802      KdInterior *interior = (KdInterior *)node;
803      nodeStack.push(interior->mFront);
804      nodeStack.push(interior->mBack);
805    }
806  }
807}
808
809// Find random neighbor which was not mailed
810KdNode *
811KdTree::FindRandomNeighbor(KdNode *n,
812                                                                                                         bool onlyUnmailed
813                                                                                                         )
814{
815  stack<KdNode *> nodeStack;
816 
817  nodeStack.push(mRoot);
818
819  AxisAlignedBox3 box = GetBox(n);
820  int mask = rand();
821
822  while (!nodeStack.empty()) {
823    KdNode *node = nodeStack.top();
824    nodeStack.pop();
825    if (node->IsLeaf()) {
826      if ( node != n && (!onlyUnmailed || !node->Mailed()) )
827        return node;
828    } else {
829      KdInterior *interior = (KdInterior *)node;
830      if (interior->mPosition > box.Max(interior->mAxis))
831        nodeStack.push(interior->mBack);
832      else
833        if (interior->mPosition < box.Min(interior->mAxis))
834          nodeStack.push(interior->mFront);
835        else {
836          // random decision
837          if (mask&1)
838            nodeStack.push(interior->mBack);
839          else
840            nodeStack.push(interior->mFront);
841          mask = mask>>1;
842        }
843    }
844  }
845 
846  return NULL;
847}
848
849int
850KdTree::FindNeighbors(KdNode *n,
851                      vector<KdNode *> &neighbors,
852                      bool onlyUnmailed
853                      )
854{
855  stack<KdNode *> nodeStack;
856 
857  nodeStack.push(mRoot);
858
859  AxisAlignedBox3 box = GetBox(n);
860
861  while (!nodeStack.empty()) {
862    KdNode *node = nodeStack.top();
863    nodeStack.pop();
864    if (node->IsLeaf()) {
865      if ( node != n && (!onlyUnmailed || !node->Mailed()) )
866        neighbors.push_back(node);
867    } else {
868      KdInterior *interior = (KdInterior *)node;
869      if (interior->mPosition > box.Max(interior->mAxis))
870                                nodeStack.push(interior->mBack);
871      else
872                                if (interior->mPosition < box.Min(interior->mAxis))
873                                        nodeStack.push(interior->mFront);
874                                else {
875                                        // random decision
876                                        nodeStack.push(interior->mBack);
877                                        nodeStack.push(interior->mFront);
878                                }
879    }
880  }
881 
882  return (int)neighbors.size();
883}
884
885// Find random neighbor which was not mailed
886KdNode *
887KdTree::GetRandomLeaf(const Plane3 &plane)
888{
889  stack<KdNode *> nodeStack;
890 
891  nodeStack.push(mRoot);
892 
893  int mask = rand();
894 
895  while (!nodeStack.empty()) {
896    KdNode *node = nodeStack.top();
897    nodeStack.pop();
898    if (node->IsLeaf()) {
899      return node;
900    } else {
901      KdInterior *interior = (KdInterior *)node;
902      KdNode *next;
903        if (GetBox(interior->mBack).Side(plane) < 0)
904          next = interior->mFront;
905        else
906          if (GetBox(interior->mFront).Side(plane) < 0)
907            next = interior->mBack;
908          else {
909            // random decision
910            if (mask&1)
911              next = interior->mBack;
912            else
913              next = interior->mFront;
914            mask = mask>>1;
915          }
916        nodeStack.push(next);
917    }
918  }
919 
920 
921  return NULL;
922}
923
924void
925KdTree::CollectLeaves(vector<KdLeaf *> &leaves)
926{
927  stack<KdNode *> nodeStack;
928  nodeStack.push(mRoot);
929
930  while (!nodeStack.empty()) {
931    KdNode *node = nodeStack.top();
932    nodeStack.pop();
933    if (node->IsLeaf()) {
934      KdLeaf *leaf = (KdLeaf *)node;
935      leaves.push_back(leaf);
936    } else {
937      KdInterior *interior = (KdInterior *)node;
938      nodeStack.push(interior->mBack);
939      nodeStack.push(interior->mFront);
940    }
941  }
942}
943
944void
945KdTree::CreateAndCollectViewCells(ViewCellContainer &vc) const
946{
947  stack<KdNode *> nodeStack;
948  nodeStack.push(mRoot);
949
950  while (!nodeStack.empty()) {
951    KdNode *node = nodeStack.top();
952    nodeStack.pop();
953    if (node->IsLeaf()) {
954      KdLeaf *leaf = (KdLeaf *)node;
955          // kdtree used as view cell container => create view cell
956          KdViewCell *viewCell = new KdViewCell();
957          leaf->mViewCell = viewCell;
958          // push back pointer to this leaf
959          viewCell->mLeaf = leaf;
960      vc.push_back(viewCell);
961    } else {
962      KdInterior *interior = (KdInterior *)node;
963      nodeStack.push(interior->mBack);
964      nodeStack.push(interior->mFront);
965    }
966  }
967}
968
969int
970KdTree::CollectLeafPvs()
971{
972  int totalPvsSize = 0;
973  stack<KdNode *> nodeStack;
974 
975  nodeStack.push(mRoot);
976 
977  while (!nodeStack.empty()) {
978    KdNode *node = nodeStack.top();
979    nodeStack.pop();
980    if (node->IsLeaf()) {
981      KdLeaf *leaf = (KdLeaf *)node;
982      for (int j=0; j < leaf->mObjects.size(); j++) {
983        Intersectable *object = leaf->mObjects[j];
984        if (!object->Mailed()) {
985          object->Mail();
986          // add this node to pvs of all nodes it can see
987          KdPvsMap::iterator ni = object->mKdPvs.mEntries.begin();
988          for (; ni != object->mKdPvs.mEntries.end(); ni++) {
989            KdNode *node = (*ni).first;
990            // $$ JB TEMPORARY solution -> should add object PVS or explictly computed
991            // kd tree PVS
992                float contribution;
993                if (leaf->mKdPvs.AddSample(node, 1.0f, contribution))
994              totalPvsSize++;
995          }
996        }
997      }
998    } else {
999      KdInterior *interior = (KdInterior *)node;
1000      nodeStack.push(interior->mFront);
1001      nodeStack.push(interior->mBack);
1002    }
1003  }
1004
1005  return totalPvsSize;
1006}
1007
1008
1009KdNode *
1010KdTree::GetRandomLeaf(const bool onlyUnmailed)
1011{
1012  stack<KdNode *> nodeStack;
1013  nodeStack.push(mRoot);
1014 
1015  int mask = rand();
1016       
1017  while (!nodeStack.empty()) {
1018    KdNode *node = nodeStack.top();
1019    nodeStack.pop();
1020    if (node->IsLeaf()) {
1021      if ( (!onlyUnmailed || !node->Mailed()) )
1022                                return node;
1023    } else {
1024      KdInterior *interior = (KdInterior *)node;
1025                        // random decision
1026                        if (mask&1)
1027                                nodeStack.push(interior->mBack);
1028                        else
1029                                nodeStack.push(interior->mFront);
1030                        mask = mask>>1;
1031                }
1032        }
1033  return NULL;
1034}
1035
1036
1037int
1038KdTree::CastBeam(
1039                                 Beam &beam
1040                                 )
1041{
1042  stack<KdNode *> nodeStack;
1043  nodeStack.push(mRoot);
1044 
1045  while (!nodeStack.empty()) {
1046    KdNode *node = nodeStack.top();
1047    nodeStack.pop();
1048       
1049        int side = beam.ComputeIntersection(GetBox(node));
1050        switch (side) {
1051        case -1:
1052          beam.mKdNodes.push_back(node);
1053          break;
1054        case 0:
1055          if (node->IsLeaf())
1056                beam.mKdNodes.push_back(node);
1057          else {
1058                KdInterior *interior = (KdInterior *)node;
1059                KdNode *first = interior->mBack;
1060                KdNode *second = interior->mFront;
1061               
1062                if (interior->mAxis < 3) {
1063                  // spatial split -> decide on the order of the nodes
1064                  if (beam.mPlanes[0].mNormal[interior->mAxis] > 0)
1065                        swap(first, second);
1066                }
1067
1068                nodeStack.push(first);
1069                nodeStack.push(second);
1070          }
1071          break;
1072          // default: cull
1073        }
1074  }
1075
1076  if (beam.mFlags & Beam::STORE_OBJECTS)
1077  {
1078          vector<KdNode *>::const_iterator it, it_end = beam.mKdNodes.end();
1079       
1080          Intersectable::NewMail();
1081          for (it = beam.mKdNodes.begin(); it != it_end; ++ it)
1082          {
1083                  CollectObjects(*it, beam.mObjects);
1084          }
1085  }
1086
1087  return (int)beam.mKdNodes.size();
1088}
1089
1090
1091}
Note: See TracBrowser for help on using the repository browser.