source: trunk/VUT/GtpVisibilityPreprocessor/src/KdTree.cpp @ 580

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

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