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

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