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

Revision 1076, 25.5 KB checked in by mattausch, 18 years ago (diff)

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