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

Revision 752, 23.4 KB checked in by mattausch, 18 years ago (diff)

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