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

Revision 2643, 40.5 KB checked in by mattausch, 17 years ago (diff)

compiling under release internal

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#include "ViewCell.h"
10#include "IntersectableWrapper.h"
11
12
13
14using namespace std;
15
16
17#define TYPE_INTERIOR -2
18#define TYPE_LEAF -3
19
20
21namespace GtpVisibilityPreprocessor {
22
23int KdNode::sMailId = 1;
24int KdNode::sReservedMailboxes = 1;
25int KdNode::sMailId2 = 1;
26int KdNode::sReservedMailboxes2 = 1;
27
28
29inline static bool ilt(Intersectable *obj1, Intersectable *obj2)
30{
31  return obj1->mId < obj2->mId;
32}
33
34
35KdNode::KdNode(KdInterior *parent):
36mParent(parent), mMailbox(0), mIntersectable(NULL), mMailbox2(0)
37{
38  if (parent)
39    mDepth = parent->mDepth+1;
40  else
41    mDepth = 0;
42}
43
44
45KdInterior::~KdInterior()
46{
47  // recursivly destroy children
48  DEL_PTR(mFront);
49  DEL_PTR(mBack);
50}
51
52
53KdLeaf::~KdLeaf()
54{
55  DEL_PTR(mViewCell); 
56}
57
58
59KdTree::KdTree()
60{
61  mRoot = new KdLeaf(NULL, 0);
62  Environment::GetSingleton()->GetIntValue("KdTree.Termination.maxNodes",
63                                           mTermMaxNodes);
64  Environment::GetSingleton()->GetIntValue("KdTree.Termination.maxDepth",
65                                           mTermMaxDepth);
66  Environment::GetSingleton()->GetIntValue("KdTree.Termination.minCost",
67                                           mTermMinCost);
68  Environment::GetSingleton()->GetFloatValue("KdTree.Termination.maxCostRatio",
69                                             mMaxCostRatio);
70  Environment::GetSingleton()->GetFloatValue("KdTree.Termination.ct_div_ci",
71                                             mCt_div_ci);
72  Environment::GetSingleton()->GetFloatValue("KdTree.splitBorder",
73                                             mSplitBorder);
74  Environment::GetSingleton()->GetFloatValue("KdTree.pvsArea",
75                                             mKdPvsArea);
76
77  Environment::GetSingleton()->GetBoolValue("KdTree.sahUseFaces",
78                                            mSahUseFaces);
79
80  char splitType[64];
81  Environment::GetSingleton()->GetStringValue("KdTree.splitMethod", splitType);
82 
83  mSplitMethod = SPLIT_SPATIAL_MEDIAN;
84  if (strcmp(splitType, "spatialMedian") == 0)
85    mSplitMethod = SPLIT_SPATIAL_MEDIAN;
86  else
87    if (strcmp(splitType, "objectMedian") == 0)
88      mSplitMethod = SPLIT_OBJECT_MEDIAN;
89    else
90      if (strcmp(splitType, "SAH") == 0)
91        mSplitMethod = SPLIT_SAH;
92      else {
93        cerr<<"Wrong kd split type "<<splitType<<endl;
94        exit(1);
95      }
96  splitCandidates = NULL;
97}
98
99
100KdTree::~KdTree()
101{
102    DEL_PTR(mRoot);
103        CLEAR_CONTAINER(mKdIntersectables);
104}
105
106
107bool
108KdTree::Construct()
109{
110  if (!splitCandidates)
111    splitCandidates = new vector<SortableEntry *>;
112
113  // first construct a leaf that will get subdivide
114  KdLeaf *leaf = (KdLeaf *) mRoot;
115
116  mStat.nodes = 1;
117
118  mBox.Initialize();
119  ObjectContainer::const_iterator mi;
120  for ( mi = leaf->mObjects.begin();
121        mi != leaf->mObjects.end();
122        mi++) {
123        //      cout<<(*mi)->GetBox()<<endl;
124        mBox.Include((*mi)->GetBox());
125  }
126
127  cout <<"KdTree Root Box:"<<mBox<<endl;
128  mRoot = Subdivide(TraversalData(leaf, mBox, 0));
129
130  // remove the allocated array
131  CLEAR_CONTAINER(*splitCandidates);
132  delete splitCandidates;
133
134  float area = GetPvsArea();
135 
136  SetPvsTerminationNodes(area);
137 
138  return true;
139}
140
141float
142KdTree::GetPvsArea() const {
143  return GetBox().SurfaceArea() * mKdPvsArea;
144}
145
146KdNode *
147KdTree::Subdivide(const TraversalData &tdata)
148{
149
150  KdNode *result = NULL;
151
152  priority_queue<TraversalData> tStack;
153  //  stack<STraversalData> tStack;
154 
155  tStack.push(tdata);
156  AxisAlignedBox3 backBox, frontBox;
157 
158  while (!tStack.empty()) {
159        //      cout<<mStat.Nodes()<<" "<<mTermMaxNodes<<endl;
160        if (mStat.Nodes() > mTermMaxNodes) {
161          //    if ( GetMemUsage() > maxMemory ) {
162      // count statistics on unprocessed leafs
163      while (!tStack.empty()) {
164                EvaluateLeafStats(tStack.top());
165                tStack.pop();
166      }
167      break;
168    }
169         
170       
171    TraversalData data = tStack.top();
172    tStack.pop();
173   
174    KdNode *node = SubdivideNode((KdLeaf *) data.mNode,
175                                 data.mBox,
176                                 backBox,
177                                 frontBox
178                                 );
179
180    if (result == NULL)
181      result = node;
182   
183    if (!node->IsLeaf()) {
184      KdInterior *interior = (KdInterior *) node;
185      // push the children on the stack
186      tStack.push(TraversalData(interior->mBack, backBox, data.mDepth+1));
187      tStack.push(TraversalData(interior->mFront, frontBox, data.mDepth+1));
188     
189    }
190    else {
191      EvaluateLeafStats(data);
192    }
193  }
194 
195  return result;
196
197}
198
199
200bool
201KdTree::TerminationCriteriaMet(const KdLeaf *leaf)
202{
203  const bool criteriaMet =
204    ((int)leaf->mObjects.size() <= mTermMinCost) ||
205    (leaf->mDepth >= mTermMaxDepth);
206 
207  if (0 && criteriaMet)
208    cerr<<"\n OBJECTS="<<(int)leaf->mObjects.size()<<endl;
209 
210  return criteriaMet;
211}
212
213
214int
215KdTree::SelectPlane(KdLeaf *leaf,
216                    const AxisAlignedBox3 &box,
217                    float &position
218                    )
219{
220  int axis = -1;
221 
222  switch (mSplitMethod)
223    {
224    case SPLIT_SPATIAL_MEDIAN: {
225      axis = box.Size().DrivingAxis();
226      position = (box.Min()[axis] + box.Max()[axis])*0.5f;
227      break;
228    }
229    case SPLIT_SAH: {
230      int objectsBack, objectsFront;
231      float costRatio;
232      bool mOnlyDrivingAxis = true;
233
234      if (mOnlyDrivingAxis) {
235        axis = box.Size().DrivingAxis();
236        costRatio = BestCostRatio(leaf,
237                                  box,
238                                  axis,
239                                  position,
240                                  objectsBack,
241                                  objectsFront);
242      }
243      else {
244        // for all 3 axes
245        costRatio = MAX_FLOAT;
246        for (int i=0; i < 3; i++) {
247          float p;
248          float r = BestCostRatio(leaf,
249                                  box,
250                                  i,
251                                  p,
252                                  objectsBack,
253                                  objectsFront);
254          if (r < costRatio) {
255            costRatio = r;
256            axis = i;
257            position = p;
258          }
259        }
260      }
261     
262      if (costRatio > mMaxCostRatio) {
263        //cout<<"Too big cost ratio "<<costRatio<<endl;
264        axis = -1;
265      }
266      break;
267    }
268   
269    }
270  return axis;
271}
272
273KdNode*
274KdTree::SubdivideNode(
275                      KdLeaf *leaf,
276                      const AxisAlignedBox3 &box,
277                      AxisAlignedBox3 &backBBox,
278                      AxisAlignedBox3 &frontBBox
279                      )
280
281  if (TerminationCriteriaMet(leaf))
282    return leaf;
283
284  float position;
285 
286  // select subdivision axis
287  int axis = SelectPlane( leaf, box, position );
288
289  if (axis == -1) {
290    return leaf;
291  }
292 
293  mStat.nodes += 2;
294  mStat.splits[axis]++;
295 
296  // add the new nodes to the tree
297  KdInterior *node = new KdInterior(leaf->mParent);
298
299  node->mAxis = axis;
300  node->mPosition = position;
301  node->mBox = box;
302 
303  backBBox = box;
304  frontBBox = box;
305 
306  // first count ray sides
307  int objectsBack = 0;
308  int objectsFront = 0;
309 
310  backBBox.SetMax(axis, position);
311  frontBBox.SetMin(axis, position);
312
313  ObjectContainer::const_iterator mi;
314
315  for ( mi = leaf->mObjects.begin();
316        mi != leaf->mObjects.end();
317        mi++) {
318    // determine the side of this ray with respect to the plane
319    AxisAlignedBox3 box = (*mi)->GetBox();
320    if (box.Max(axis) > position )
321      objectsFront++;
322   
323    if (box.Min(axis) < position )
324      objectsBack++;
325  }
326
327  KdLeaf *back = new KdLeaf(node, objectsBack);
328  KdLeaf *front = new KdLeaf(node, objectsFront);
329 
330  // replace a link from node's parent
331  if (  leaf->mParent )
332    leaf->mParent->ReplaceChildLink(leaf, node);
333
334  // and setup child links
335  node->SetupChildLinks(back, front);
336 
337  for (mi = leaf->mObjects.begin();
338       mi != leaf->mObjects.end();
339       mi++) {
340    // determine the side of this ray with respect to the plane
341    AxisAlignedBox3 box = (*mi)->GetBox();
342       
343        // matt: no more ref
344        // for handling multiple objects: keep track of references
345        //if (leaf->IsRoot())
346        //      (*mi)->mReferences = 1; // initialise references at root
347       
348        // matt: no more ref
349        //-- (*mi)->mReferences; // remove parent ref
350
351
352    if (box.Max(axis) >= position )
353    {
354      front->mObjects.push_back(*mi);
355      //++ (*mi)->mReferences;
356    }
357       
358    if (box.Min(axis) < position )
359    {
360      back->mObjects.push_back(*mi);
361      // matt: no more ref
362      //  ++ (*mi)->mReferences;
363    }
364
365    mStat.objectRefs -= (int)leaf->mObjects.size();
366    mStat.objectRefs += objectsBack + objectsFront;
367  }
368
369  // store objects referenced in more than one leaf
370  // for easy access
371  ProcessMultipleRefs(back);
372  ProcessMultipleRefs(front);
373 
374  delete leaf;
375  return node;
376}
377
378
379void KdTree::ProcessMultipleRefs(KdLeaf *leaf) const
380{
381  // find objects from multiple kd-leaves
382  ObjectContainer::const_iterator oit, oit_end = leaf->mObjects.end();
383
384  for (oit = leaf->mObjects.begin(); oit != oit_end; ++ oit)
385  {
386    Intersectable *object = *oit;
387               
388    // matt: no more ref
389    /*
390      if (object->mReferences > 1) {
391        leaf->mMultipleObjects.push_back(object);
392      }
393    */
394  }
395}
396
397
398void
399KdTreeStatistics::Print(ostream &app) const
400{
401  app << "===== KdTree statistics ===============\n";
402
403  app << "#N_NODES ( Number of nodes )\n" << nodes << "\n";
404
405  app << "#N_LEAVES ( Number of leaves )\n" << Leaves() << "\n";
406
407  app << "#N_SPLITS ( Number of splits in axes x y z dx dy dz)\n";
408  for (int i=0; i<7; i++)
409    app << splits[i] <<" ";
410  app <<endl;
411
412  app << "#N_RAYREFS ( Number of rayRefs )\n" <<
413    rayRefs << "\n";
414
415  app << "#N_RAYRAYREFS  ( Number of rayRefs / ray )\n" <<
416    rayRefs/(double)rays << "\n";
417
418  app << "#N_LEAFRAYREFS  ( Number of rayRefs / leaf )\n" <<
419    rayRefs/(double)Leaves() << "\n";
420
421  app << "#N_MAXOBJECTREFS  ( Max number of object refs / leaf )\n" <<
422    maxObjectRefs << "\n";
423
424  app << "#N_NONEMPTYRAYREFS  ( Number of rayRefs in nonEmpty leaves / non empty leaf )\n" <<
425    rayRefsNonZeroQuery/(double)(Leaves() - zeroQueryNodes) << "\n";
426
427  app << "#N_LEAFDOMAINREFS  ( Number of query domain Refs / leaf )\n" <<
428    objectRefs/(double)Leaves() << "\n";
429
430  //  app << setprecision(4);
431
432  app << "#N_PEMPTYLEAVES  ( Percentage of leaves with zero query domains )\n"<<
433    zeroQueryNodes*100/(double)Leaves()<<endl;
434
435  app << "#N_PMAXDEPTHLEAVES ( Percentage of leaves at maxdepth )\n"<<
436    maxDepthNodes*100/(double)Leaves()<<endl;
437
438  app << "#N_PMINCOSTLEAVES  ( Percentage of leaves with minCost )\n"<<
439    minCostNodes*100/(double)Leaves()<<endl;
440
441  app << "#N_ADDED_RAYREFS  (Number of dynamically added ray references )\n"<<
442    addedRayRefs<<endl;
443
444  app << "#N_REMOVED_RAYREFS  (Number of dynamically removed ray references )\n"<<
445    removedRayRefs<<endl;
446
447  //  app << setprecision(4);
448
449  //  app << "#N_CTIME  ( Construction time [s] )\n"
450  //      << Time() << " \n";
451
452  app << "===== END OF KdTree statistics ==========\n";
453
454}
455
456
457
458void
459KdTree::EvaluateLeafStats(const TraversalData &data)
460{
461
462  // the node became a leaf -> evaluate stats for leafs
463  KdLeaf *leaf = (KdLeaf *)data.mNode;
464
465  if (data.mDepth > mTermMaxDepth)
466    mStat.maxDepthNodes++;
467 
468  if ( (int)(leaf->mObjects.size()) < mTermMinCost)
469    mStat.minCostNodes++;
470 
471 
472  if ( (int)(leaf->mObjects.size()) > mStat.maxObjectRefs)
473    mStat.maxObjectRefs = (int)leaf->mObjects.size();
474 
475}
476
477
478
479void
480KdTree::SortSubdivisionCandidates(
481                                  KdLeaf *node,
482                                  const int axis
483                                  )
484{
485  CLEAR_CONTAINER(*splitCandidates);
486  //splitCandidates->clear();
487 
488  int requestedSize = 2*(int)node->mObjects.size();
489 
490  // creates a sorted split candidates array
491  if (splitCandidates->capacity() > 500000 &&
492      requestedSize < (int)(splitCandidates->capacity()/10) ) {         
493    delete splitCandidates;
494    splitCandidates = new vector<SortableEntry *>;
495  }
496 
497  splitCandidates->reserve(requestedSize);
498 
499  // insert all queries
500  for(ObjectContainer::const_iterator mi = node->mObjects.begin();
501      mi != node->mObjects.end();
502      mi++)
503  {
504    AxisAlignedBox3 box = (*mi)->GetBox();
505
506    splitCandidates->push_back(new SortableEntry(SortableEntry::BOX_MIN,
507                                                 box.Min(axis),
508                                                 *mi)
509                               );
510   
511    splitCandidates->push_back(new SortableEntry(SortableEntry::BOX_MAX,
512                                                 box.Max(axis),
513                                                 *mi)
514                               );
515  }
516 
517  stable_sort(splitCandidates->begin(), splitCandidates->end(), iltS);
518}
519
520
521float
522KdTree::BestCostRatio(
523                      KdLeaf *node,
524                      const AxisAlignedBox3 &box,
525                      const int axis,
526                      float &position,
527                      int &objectsBack,
528                      int &objectsFront
529                      )
530{
531
532#define DEBUG_COST 0
533
534#if DEBUG_COST
535  static int nodeId = -1;
536  char filename[256];
537 
538  static int lastAxis = 100;
539  if (axis <= lastAxis)
540        nodeId++;
541
542  lastAxis = axis;
543 
544  sprintf(filename, "sah-cost%d-%d.log", nodeId, axis);
545  ofstream costStream;
546 
547  if (nodeId < 100)
548    costStream.open(filename);
549
550#endif
551 
552  SortSubdivisionCandidates(node, axis);
553 
554  // go through the lists, count the number of objects left and right
555  // and evaluate the following cost funcion:
556  // C = ct_div_ci  + (ol + or)/queries
557
558  float totalIntersections = 0.0f;
559  vector<SortableEntry *>::const_iterator ci;
560
561  for(ci = splitCandidates->begin();
562      ci < splitCandidates->end();
563      ci++)
564    if ((*ci)->type == SortableEntry::BOX_MIN) {
565      totalIntersections += (*ci)->intersectable->IntersectionComplexity();
566    }
567       
568  float intersectionsLeft = 0;
569  float intersectionsRight = totalIntersections;
570       
571  int objectsLeft = 0, objectsRight = (int)node->mObjects.size();
572 
573  float minBox = box.Min(axis);
574  float maxBox = box.Max(axis);
575  float boxArea = box.SurfaceArea();
576 
577  float minBand = minBox + mSplitBorder*(maxBox - minBox);
578  float maxBand = minBox + (1.0f - mSplitBorder)*(maxBox - minBox);
579 
580  float minSum = 1e20f;
581 
582  for(ci = splitCandidates->begin();
583      ci < splitCandidates->end();
584      ci++) {
585    switch ((*ci)->type) {
586    case SortableEntry::BOX_MIN:
587      objectsLeft++;
588      intersectionsLeft += (*ci)->intersectable->IntersectionComplexity();
589      break;
590    case SortableEntry::BOX_MAX:
591      objectsRight--;
592      intersectionsRight -= (*ci)->intersectable->IntersectionComplexity();
593      break;
594    } // switch
595
596    if ((*ci)->value > minBand && (*ci)->value < maxBand) {
597      AxisAlignedBox3 lbox = box;
598      AxisAlignedBox3 rbox = box;
599      lbox.SetMax(axis, (*ci)->value);
600      rbox.SetMin(axis, (*ci)->value);
601
602      float sum;
603      if (mSahUseFaces)
604        sum = intersectionsLeft*lbox.SurfaceArea() + intersectionsRight*rbox.SurfaceArea();
605      else
606        sum = objectsLeft*lbox.SurfaceArea() + objectsRight*rbox.SurfaceArea();
607     
608      //      cout<<"pos="<<(*ci).value<<"\t q=("<<ql<<","<<qr<<")\t r=("<<rl<<","<<rr<<")"<<endl;
609      //      cout<<"cost= "<<sum<<endl;
610
611#if DEBUG_COST
612      if (nodeId < 100) {
613        float oldCost = mSahUseFaces ? totalIntersections : node->mObjects.size();
614        float newCost = mCt_div_ci + sum/boxArea;
615        float ratio = newCost/oldCost;
616        costStream<<(*ci)->value<<" "<<ratio<<endl;
617      }
618#endif
619     
620      if (sum < minSum) {
621        minSum = sum;
622        position = (*ci)->value;
623       
624        objectsBack = objectsLeft;
625        objectsFront = objectsRight;
626      }
627    }
628  } // for ci
629 
630  float oldCost = mSahUseFaces ? totalIntersections : node->mObjects.size();
631  float newCost = mCt_div_ci + minSum/boxArea;
632  float ratio = newCost/oldCost;
633 
634#if 0
635  cout<<"===================="<<endl;
636  cout<<"costRatio="<<ratio<<" pos="<<position<<" t="<<(position - minBox)/(maxBox - minBox)
637      <<"\t o=("<<objectsBack<<","<<objectsFront<<")"<<endl;
638#endif
639  return ratio;
640}
641
642int
643KdTree::CastRay(
644                Ray &ray
645                )
646{
647 
648  int hits = 0;
649 
650  stack<RayTraversalData> tStack;
651 
652  float maxt = 1e6;
653  float mint = 0;
654
655  //  ray.ComputeInvertedDir();
656  Intersectable::NewMail();
657
658  if (!mBox.GetMinMaxT(ray, &mint, &maxt))
659    return 0;
660
661  if (mint < 0)
662    mint = 0;
663
664 
665  maxt += Limits::Threshold;
666 
667  Vector3 entp = ray.Extrap(mint);
668  Vector3 extp = ray.Extrap(maxt);
669 
670  KdNode *node = mRoot;
671  KdNode *farChild;
672  float position;
673  int axis;
674
675 
676  while (1) {
677    if (!node->IsLeaf()) {
678      KdInterior *in = (KdInterior *) node;
679      position = in->mPosition;
680      axis = in->mAxis;
681
682      if (entp[axis] <= position) {
683        if (extp[axis] <= position) {
684          node = in->mBack;
685          // cases N1,N2,N3,P5,Z2,Z3
686          continue;
687        }
688        else {
689          // case N4
690          node = in->mBack;
691          farChild = in->mFront;
692        }
693      }
694      else {
695        if (position <= extp[axis]) {
696          node = in->mFront;
697          // cases P1,P2,P3,N5,Z1
698          continue;
699        }
700        else {
701          node = in->mFront;
702          farChild = in->mBack;
703          // case P4
704        }
705      }
706      // $$ modification 3.5.2004 - hints from Kamil Ghais
707      // case N4 or P4
708      float tdist = (position - ray.GetLoc(axis)) / ray.GetDir(axis);
709      tStack.push(RayTraversalData(farChild, extp, maxt));
710      extp = ray.GetLoc() + ray.GetDir()*tdist;
711      maxt = tdist;
712    }
713    else {
714      // compute intersection with all objects in this leaf
715      KdLeaf *leaf = (KdLeaf *) node;
716      if (ray.mFlags & Ray::STORE_KDLEAVES)
717        ray.kdLeaves.push_back(leaf);
718         
719      ObjectContainer::const_iterator mi;
720      for ( mi = leaf->mObjects.begin();
721            mi != leaf->mObjects.end();
722            mi++) {
723        Intersectable *object = *mi;
724        if (!object->Mailed() ) {
725          object->Mail();
726          if (ray.mFlags & Ray::STORE_TESTED_OBJECTS)
727                        ray.testedObjects.push_back(object);
728                 
729          static int oi=1;
730          if (MeshDebug)
731            cout<<"Object "<<oi++;
732         
733          hits += object->CastRay(ray);
734         
735          if (MeshDebug) {
736            if (!ray.intersections.empty())
737              cout<<"nearest t="<<ray.intersections[0].mT<<endl;
738            else
739              cout<<"nearest t=-INF"<<endl;
740          }       
741        }
742      }
743         
744      if (hits && ray.GetType() == Ray::LOCAL_RAY)
745        if (ray.intersections[0].mT <= maxt)
746          break;
747     
748      // get the next node from the stack
749      if (tStack.empty())
750        break;
751     
752      entp = extp;
753      mint = maxt;
754      if (ray.GetType() == Ray::LINE_SEGMENT && mint > 1.0f)
755        break;
756     
757      RayTraversalData &s  = tStack.top();
758      node = s.mNode;
759      extp = s.mExitPoint;
760      maxt = s.mMaxT;
761      tStack.pop();
762    }
763  }
764  return hits;
765}
766
767int KdTree::CastLineSegment2(const Vector3 &origin,
768                            const Vector3 &termination,
769                            ViewCellContainer &viewcells)
770{
771  int hits = 0;
772
773  float mint = 0.0f, maxt = 1.0f;
774  const Vector3 dir = termination - origin;
775 
776  stack<RayTraversalData> tStack;
777 
778  Intersectable::NewMail();
779 
780  //maxt += Limits::Threshold;
781 
782  Vector3 entp = origin;
783  Vector3 extp = termination;
784 
785  KdNode *node = mRoot;
786  KdNode *farChild;
787 
788  float position;
789  int axis;
790
791  while (1)
792  {
793    if (!node->IsLeaf())
794    {
795      KdInterior *in = static_cast<KdInterior *>(node);
796      position = in->mPosition;
797      axis = in->mAxis;
798     
799      if (entp[axis] <= position)
800      {
801        if (extp[axis] <= position)
802        {
803          node = in->mBack;
804          // cases N1,N2,N3,P5,Z2,Z3
805          continue;
806        }
807        else
808        {
809          // case N4
810          node = in->mBack;
811          farChild = in->mFront;
812        }
813      }
814      else
815      {
816        if (position <= extp[axis])
817        {
818          node = in->mFront;
819          // cases P1,P2,P3,N5,Z1
820          continue;
821        }
822        else
823        {
824          node = in->mFront;
825          farChild = in->mBack;
826          // case P4
827        }
828      }
829     
830      // $$ modification 3.5.2004 - hints from Kamil Ghais
831      // case N4 or P4
832      float tdist = (position - origin[axis]) / dir[axis];
833      //tStack.push(RayTraversalData(farChild, extp, maxt)); //TODO
834      extp = origin + dir * tdist;
835      maxt = tdist;
836    }
837    else
838    {
839      // compute intersection with all objects in this leaf
840      KdLeaf *leaf = static_cast<KdLeaf *>(node);
841     
842      // add view cell to intersections
843      ViewCell *vc = leaf->mViewCell;
844     
845      if (!vc->Mailed())
846      {
847        vc->Mail();
848        viewcells.push_back(vc);
849        ++ hits;
850      }
851     
852      // get the next node from the stack
853      if (tStack.empty())
854        break;
855     
856      entp = extp;
857      mint = maxt;
858     
859      RayTraversalData &s  = tStack.top();
860      node = s.mNode;
861      extp = s.mExitPoint;
862      maxt = s.mMaxT;
863      tStack.pop();
864    }
865  }
866  return hits;
867}
868
869int KdTree::CastLineSegment(const Vector3 &origin,
870                            const Vector3 &termination,
871                            ViewCellContainer &viewcells)
872{
873  int hits = 0;
874
875  float mint = 0.0f, maxt = 1.0f;
876  const Vector3 dir = termination - origin;
877 
878  stack<RayTraversalData> tStack;
879 
880  Intersectable::NewMail();
881 
882  //maxt += Limits::Threshold;
883 
884  Vector3 entp = origin;
885  Vector3 extp = termination;
886 
887  KdNode *node = mRoot;
888  KdNode *farChild;
889 
890  float position;
891  float tdist;
892
893  while (1)
894  {
895    if (!node->IsLeaf())
896    {
897      KdInterior *in = static_cast<KdInterior *>(node);
898      position = in->mPosition;
899      switch (in->mAxis) {
900          case 0:
901                if (entp.x <= position)
902                  {
903                        if (extp.x <= position)
904                          {
905                                node = in->mBack;
906                                // cases N1,N2,N3,P5,Z2,Z3
907                                continue;
908                          }
909                        else
910                          {
911                                // case N4
912                                node = in->mBack;
913                                farChild = in->mFront;
914                          }
915                  }
916                else
917                  {
918                        if (position <= extp.x)
919                          {
920                                node = in->mFront;
921                                // cases P1,P2,P3,N5,Z1
922                                continue;
923                          }
924                        else
925                          {
926                                node = in->mFront;
927                                farChild = in->mBack;
928                                // case P4
929                          }
930                  }
931               
932                // $$ modification 3.5.2004 - hints from Kamil Ghais
933                // case N4 or P4
934                tdist = (position - origin.x) / dir.x;
935                //tStack.push(RayTraversalData(farChild, extp, maxt)); //TODO
936                extp = origin + dir * tdist;
937                maxt = tdist;
938         
939          break;
940        case 1:
941       
942                if (entp.y <= position)
943                  {
944                        if (extp.y <= position)
945                          {
946                                node = in->mBack;
947                                // cases N1,N2,N3,P5,Z2,Z3
948                                continue;
949                          }
950                        else
951                          {
952                                // case N4
953                                node = in->mBack;
954                                farChild = in->mFront;
955                          }
956                  }
957                else
958                  {
959                        if (position <= extp.y)
960                          {
961                                node = in->mFront;
962                                // cases P1,P2,P3,N5,Z1
963                                continue;
964                          }
965                        else
966                          {
967                                node = in->mFront;
968                                farChild = in->mBack;
969                                // case P4
970                          }
971                  }
972               
973                // $$ modification 3.5.2004 - hints from Kamil Ghais
974                // case N4 or P4
975                tdist = (position - origin.y) / dir.y;
976                //tStack.push(RayTraversalData(farChild, extp, maxt)); //TODO
977                extp = origin + dir * tdist;
978                maxt = tdist;
979               
980                break;
981          case 2:
982                if (entp.z <= position)
983                  {
984                        if (extp.z <= position)
985                          {
986                                node = in->mBack;
987                                // cases N1,N2,N3,P5,Z2,Z3
988                                continue;
989                          }
990                        else
991                          {
992                                // case N4
993                                node = in->mBack;
994                                farChild = in->mFront;
995                          }
996                  }
997                else
998                  {
999                        if (position <= extp.z)
1000                          {
1001                                node = in->mFront;
1002                                // cases P1,P2,P3,N5,Z1
1003                                continue;
1004                          }
1005                        else
1006                          {
1007                                node = in->mFront;
1008                                farChild = in->mBack;
1009                                // case P4
1010                          }
1011                  }
1012               
1013                // $$ modification 3.5.2004 - hints from Kamil Ghais
1014                // case N4 or P4
1015                tdist = (position - origin.z) / dir.z;
1016                //tStack.push(RayTraversalData(farChild, extp, maxt)); //TODO
1017                extp = origin + dir * tdist;
1018                maxt = tdist;
1019         
1020          break;
1021         
1022          }
1023         
1024        }
1025    else
1026          {
1027      // compute intersection with all objects in this leaf
1028      KdLeaf *leaf = static_cast<KdLeaf *>(node);
1029     
1030      // add view cell to intersections
1031      ViewCell *vc = leaf->mViewCell;
1032     
1033      if (!vc->Mailed())
1034      {
1035        vc->Mail();
1036        viewcells.push_back(vc);
1037        ++ hits;
1038      }
1039     
1040      // get the next node from the stack
1041      if (tStack.empty())
1042                break;
1043     
1044      entp = extp;
1045      mint = maxt;
1046     
1047      RayTraversalData &s  = tStack.top();
1048      node = s.mNode;
1049      extp = s.mExitPoint;
1050      maxt = s.mMaxT;
1051      tStack.pop();
1052    }
1053  }
1054  return hits;
1055}
1056
1057void
1058KdTree::CollectKdObjects(const AxisAlignedBox3 &box,
1059                         ObjectContainer &objects
1060                        )
1061{
1062  stack<KdNode *> nodeStack;
1063 
1064  nodeStack.push(mRoot);
1065
1066  while (!nodeStack.empty()) {
1067    KdNode *node = nodeStack.top();
1068    nodeStack.pop();
1069    if (node->IsLeaf() || node->mPvsTermination == 1)  {
1070      if (!node->Mailed()) {
1071                Intersectable *object = GetOrCreateKdIntersectable(node);
1072                node->Mail();
1073                objects.push_back(object);
1074      }
1075    }
1076    else {
1077      KdInterior *interior = (KdInterior *)node;
1078         
1079      if ( box.Max()[interior->mAxis] > interior->mPosition )
1080        nodeStack.push(interior->mFront);
1081     
1082      if (box.Min()[interior->mAxis] < interior->mPosition)
1083        nodeStack.push(interior->mBack);
1084    }
1085  } // while
1086}
1087
1088void
1089KdTree::CollectSmallKdObjects(const AxisAlignedBox3 &box,
1090                                                          ObjectContainer &objects,
1091                                                          const float maxAreaFrac
1092                                                          )
1093{
1094  stack<KdNode *> nodeStack;
1095 
1096  nodeStack.push(mRoot);
1097  float maxArea = GetPvsArea()*maxAreaFrac;
1098 
1099  while (!nodeStack.empty()) {
1100    KdNode *node = nodeStack.top();
1101    nodeStack.pop();
1102        if (!node->Mailed()) {
1103          if (node->IsLeaf() || GetSurfaceArea(node) <= maxArea)  {
1104                Intersectable *object = GetOrCreateKdIntersectable(node);
1105                if (!node->Mailed()) {
1106                  node->Mail();
1107                  objects.push_back(object);
1108                }
1109          }
1110          else {
1111                KdInterior *interior = (KdInterior *)node;
1112               
1113                if ( box.Max()[interior->mAxis] > interior->mPosition )
1114                  nodeStack.push(interior->mFront);
1115               
1116                if (box.Min()[interior->mAxis] < interior->mPosition)
1117                  nodeStack.push(interior->mBack);
1118          }
1119        }
1120  } // while
1121}
1122
1123
1124void
1125KdTree::CollectObjects(const AxisAlignedBox3 &box,
1126                       ObjectContainer &objects)
1127{
1128  stack<KdNode *> nodeStack;
1129
1130  nodeStack.push(mRoot);
1131
1132  while (!nodeStack.empty()) {
1133    KdNode *node = nodeStack.top();
1134    nodeStack.pop();
1135    if (node->IsLeaf()) {
1136      KdLeaf *leaf = (KdLeaf *)node;
1137      for (int j=0; j < leaf->mObjects.size(); j++) {
1138        Intersectable *object = leaf->mObjects[j];
1139        if (!object->Mailed() && Overlap(box, object->GetBox())) {
1140          object->Mail();
1141          objects.push_back(object);
1142        }
1143      }
1144    } else {
1145      KdInterior *interior = (KdInterior *)node;
1146
1147      if ( box.Max()[interior->mAxis] > interior->mPosition )
1148        nodeStack.push(interior->mFront);
1149     
1150      if (box.Min()[interior->mAxis] < interior->mPosition)
1151        nodeStack.push(interior->mBack);
1152    }
1153  }
1154}
1155
1156
1157void KdTree::CollectObjectsWithoutMail(KdNode *n, ObjectContainer &objects)
1158{
1159        stack<KdNode *> nodeStack;
1160
1161        nodeStack.push(n);
1162
1163        while (!nodeStack.empty())
1164        {
1165                KdNode *node = nodeStack.top();
1166                nodeStack.pop();
1167
1168                if (node->IsLeaf())
1169                {
1170                        KdLeaf *leaf = (KdLeaf *)node;
1171                        for (int j=0; j < leaf->mObjects.size(); j++)
1172                        {
1173                                Intersectable *object = leaf->mObjects[j];
1174                                objects.push_back(object);
1175                        }
1176                }
1177                else
1178                {
1179                        KdInterior *interior = (KdInterior *)node;
1180
1181                        nodeStack.push(interior->mFront);
1182                        nodeStack.push(interior->mBack);
1183                }
1184        }
1185}
1186
1187
1188void
1189KdTree::CollectObjects(KdNode *n, ObjectContainer &objects)
1190{
1191  stack<KdNode *> nodeStack;
1192
1193  nodeStack.push(n);
1194
1195  while (!nodeStack.empty()) {
1196    KdNode *node = nodeStack.top();
1197    nodeStack.pop();
1198    if (node->IsLeaf()) {
1199      KdLeaf *leaf = (KdLeaf *)node;
1200      for (int j=0; j < leaf->mObjects.size(); j++) {
1201        Intersectable *object = leaf->mObjects[j];
1202        if (!object->Mailed()) {
1203          object->Mail();
1204          objects.push_back(object);
1205        }
1206      }
1207    } else {
1208      KdInterior *interior = (KdInterior *)node;
1209      nodeStack.push(interior->mFront);
1210      nodeStack.push(interior->mBack);
1211    }
1212  }
1213}
1214
1215// Find random neighbor which was not mailed
1216KdNode *
1217KdTree::FindRandomNeighbor(KdNode *n,
1218                           bool onlyUnmailed
1219                           )
1220{
1221  stack<KdNode *> nodeStack;
1222 
1223  nodeStack.push(mRoot);
1224
1225  AxisAlignedBox3 box = GetBox(n);
1226  int mask = rand();
1227
1228  while (!nodeStack.empty()) {
1229    KdNode *node = nodeStack.top();
1230    nodeStack.pop();
1231    if (node->IsLeaf()) {
1232      if ( node != n && (!onlyUnmailed || !node->Mailed()) )
1233        return node;
1234    } else {
1235      KdInterior *interior = (KdInterior *)node;
1236      if (interior->mPosition > box.Max(interior->mAxis))
1237        nodeStack.push(interior->mBack);
1238      else
1239        if (interior->mPosition < box.Min(interior->mAxis))
1240          nodeStack.push(interior->mFront);
1241        else {
1242          // random decision
1243          if (mask&1)
1244            nodeStack.push(interior->mBack);
1245          else
1246            nodeStack.push(interior->mFront);
1247          mask = mask>>1;
1248        }
1249    }
1250  }
1251 
1252  return NULL;
1253}
1254
1255int
1256KdTree::FindNeighbors(KdNode *n,
1257                      vector<KdNode *> &neighbors,
1258                      bool onlyUnmailed
1259                      )
1260{
1261  stack<KdNode *> nodeStack;
1262 
1263  nodeStack.push(mRoot);
1264
1265  AxisAlignedBox3 box = GetBox(n);
1266
1267  while (!nodeStack.empty()) {
1268    KdNode *node = nodeStack.top();
1269    nodeStack.pop();
1270    if (node->IsLeaf()) {
1271      if ( node != n && (!onlyUnmailed || !node->Mailed()) )
1272        neighbors.push_back(node);
1273    }
1274    else {
1275      KdInterior *interior = (KdInterior *)node;
1276      if (interior->mPosition > box.Max(interior->mAxis))
1277        nodeStack.push(interior->mBack);
1278      else {
1279        if (interior->mPosition < box.Min(interior->mAxis))
1280          nodeStack.push(interior->mFront);
1281        else {
1282          // random decision
1283          nodeStack.push(interior->mBack);
1284          nodeStack.push(interior->mFront);
1285        }
1286      }
1287    }
1288  }
1289 
1290  return (int)neighbors.size();
1291}
1292
1293// Find random neighbor which was not mailed
1294KdNode *
1295KdTree::GetRandomLeaf(const Plane3 &plane)
1296{
1297  stack<KdNode *> nodeStack;
1298 
1299  nodeStack.push(mRoot);
1300 
1301  int mask = rand();
1302 
1303  while (!nodeStack.empty()) {
1304    KdNode *node = nodeStack.top();
1305    nodeStack.pop();
1306    if (node->IsLeaf()) {
1307      return node;
1308    } else {
1309      KdInterior *interior = (KdInterior *)node;
1310      KdNode *next;
1311        if (GetBox(interior->mBack).Side(plane) < 0)
1312          next = interior->mFront;
1313        else
1314          if (GetBox(interior->mFront).Side(plane) < 0)
1315            next = interior->mBack;
1316          else {
1317            // random decision
1318            if (mask&1)
1319              next = interior->mBack;
1320            else
1321              next = interior->mFront;
1322            mask = mask>>1;
1323          }
1324        nodeStack.push(next);
1325    }
1326  }
1327 
1328 
1329  return NULL;
1330}
1331
1332void
1333KdTree::CollectLeaves(vector<KdLeaf *> &leaves)
1334{
1335  stack<KdNode *> nodeStack;
1336  nodeStack.push(mRoot);
1337
1338  while (!nodeStack.empty()) {
1339    KdNode *node = nodeStack.top();
1340    nodeStack.pop();
1341    if (node->IsLeaf()) {
1342      KdLeaf *leaf = (KdLeaf *)node;
1343      leaves.push_back(leaf);
1344    } else {
1345      KdInterior *interior = (KdInterior *)node;
1346      nodeStack.push(interior->mBack);
1347      nodeStack.push(interior->mFront);
1348    }
1349  }
1350}
1351
1352void
1353KdTree::CreateAndCollectViewCells(ViewCellContainer &vc) const
1354{
1355  stack<KdNode *> nodeStack;
1356  nodeStack.push(mRoot);
1357
1358  while (!nodeStack.empty()) {
1359    KdNode *node = nodeStack.top();
1360    nodeStack.pop();
1361    if (node->IsLeaf()) {
1362      KdLeaf *leaf = (KdLeaf *)node;
1363          // kdtree used as view cell container => create view cell
1364          KdViewCell *viewCell = new KdViewCell();
1365          leaf->mViewCell = viewCell;
1366          // push back pointer to this leaf
1367          viewCell->mLeaves[0] = leaf;
1368      vc.push_back(viewCell);
1369    } else {
1370      KdInterior *interior = (KdInterior *)node;
1371      nodeStack.push(interior->mBack);
1372      nodeStack.push(interior->mFront);
1373    }
1374  }
1375}
1376
1377
1378KdNode *
1379KdTree::GetRandomLeaf(const bool onlyUnmailed)
1380{
1381  stack<KdNode *> nodeStack;
1382  nodeStack.push(mRoot);
1383 
1384  int mask = rand();
1385       
1386  while (!nodeStack.empty()) {
1387    KdNode *node = nodeStack.top();
1388    nodeStack.pop();
1389    if (node->IsLeaf()) {
1390      if ( (!onlyUnmailed || !node->Mailed()) )
1391                                return node;
1392    } else {
1393      KdInterior *interior = (KdInterior *)node;
1394                        // random decision
1395                        if (mask&1)
1396                                nodeStack.push(interior->mBack);
1397                        else
1398                                nodeStack.push(interior->mFront);
1399                        mask = mask>>1;
1400                }
1401        }
1402  return NULL;
1403}
1404
1405
1406int KdTree::CastBeam(Beam &beam)
1407{
1408  stack<KdNode *> nodeStack;
1409  nodeStack.push(mRoot);
1410 
1411  while (!nodeStack.empty()) {
1412    KdNode *node = nodeStack.top();
1413    nodeStack.pop();
1414       
1415        int side = beam.ComputeIntersection(GetBox(node));
1416        switch (side) {
1417        case -1:
1418          beam.mKdNodes.push_back(node);
1419          break;
1420        case 0:
1421          if (node->IsLeaf())
1422                beam.mKdNodes.push_back(node);
1423          else {
1424                KdInterior *interior = (KdInterior *)node;
1425                KdNode *first = interior->mBack;
1426                KdNode *second = interior->mFront;
1427               
1428                if (interior->mAxis < 3) {
1429                  // spatial split -> decide on the order of the nodes
1430                  if (beam.mPlanes[0].mNormal[interior->mAxis] > 0)
1431                        swap(first, second);
1432                }
1433
1434                nodeStack.push(first);
1435                nodeStack.push(second);
1436          }
1437          break;
1438          // default: cull
1439        }
1440  }
1441
1442  if (beam.mFlags & Beam::STORE_OBJECTS)
1443  {
1444          vector<KdNode *>::const_iterator it, it_end = beam.mKdNodes.end();
1445       
1446          Intersectable::NewMail();
1447          for (it = beam.mKdNodes.begin(); it != it_end; ++ it)
1448          {
1449                  CollectObjects(*it, beam.mObjects);
1450          }
1451  }
1452
1453  return (int)beam.mKdNodes.size();
1454}
1455
1456
1457void KdTree::ExportBinLeaf(OUT_STREAM &stream, KdLeaf *leaf)
1458{
1459  ObjectContainer::const_iterator it, it_end = leaf->mObjects.end();
1460 
1461  int type = TYPE_LEAF;
1462  int size = (int)leaf->mObjects.size();
1463 
1464  stream.write(reinterpret_cast<char *>(&type), sizeof(int));
1465  stream.write(reinterpret_cast<char *>(&size), sizeof(int));
1466 
1467  for (it = leaf->mObjects.begin(); it != it_end; ++ it)
1468    {   
1469      Intersectable *obj = *it;
1470      int id = obj->mId;               
1471     
1472      //stream.write(reinterpret_cast<char *>(&origin), sizeof(Vector3));
1473      stream.write(reinterpret_cast<char *>(&id), sizeof(int));
1474    }
1475}
1476
1477
1478KdLeaf *KdTree::ImportBinLeaf(IN_STREAM &stream,
1479                              KdInterior *parent,
1480                              const ObjectContainer &objects)
1481{
1482  int leafId = TYPE_LEAF;
1483  int objId = leafId;
1484  int size;
1485 
1486  stream.read(reinterpret_cast<char *>(&size), sizeof(int));
1487  KdLeaf *leaf = new KdLeaf(parent, size);
1488
1489  MeshInstance dummyInst(NULL);
1490       
1491  // read object ids
1492  // note: this can also be done geometrically
1493  for (int i = 0; i < size; ++ i)
1494    {   
1495      stream.read(reinterpret_cast<char *>(&objId), sizeof(int));
1496      dummyInst.SetId(objId);
1497     
1498      ObjectContainer::const_iterator oit =
1499        lower_bound(objects.begin(), objects.end(), (Intersectable *)&dummyInst, ilt);
1500     
1501      if ((oit != objects.end()) && ((*oit)->GetId() == objId))
1502        leaf->mObjects.push_back(*oit);
1503      else
1504        Debug << "error: object with id " << objId << " does not exist" << endl;
1505    }
1506 
1507  return leaf;
1508}
1509
1510
1511void KdTree::ExportBinInterior(OUT_STREAM &stream, KdInterior *interior)
1512{
1513        int interiorid = TYPE_INTERIOR;
1514        stream.write(reinterpret_cast<char *>(&interiorid), sizeof(int));
1515
1516        int axis = interior->mAxis;
1517        float pos = interior->mPosition;
1518
1519        stream.write(reinterpret_cast<char *>(&axis), sizeof(int));
1520        stream.write(reinterpret_cast<char *>(&pos), sizeof(float));
1521}
1522
1523
1524KdInterior *KdTree::ImportBinInterior(IN_STREAM  &stream, KdInterior *parent)
1525{
1526        KdInterior *interior = new KdInterior(parent);
1527
1528        int axis;
1529        float pos;
1530
1531        stream.read(reinterpret_cast<char *>(&axis), sizeof(int));
1532        stream.read(reinterpret_cast<char *>(&pos), sizeof(float));
1533
1534        interior->mAxis = axis;
1535        interior->mPosition = pos;
1536
1537        return interior;
1538}
1539
1540
1541bool KdTree::ExportBinTree(const string &filename)
1542{
1543        OUT_STREAM stream(filename.c_str(), OUT_BIN_MODE);
1544       
1545        if (!stream.is_open()) return false;
1546
1547        // export binary version of mesh
1548        queue<KdNode *> tStack;
1549        tStack.push(mRoot);
1550
1551        while(!tStack.empty())
1552        {
1553                KdNode *node = tStack.front();
1554                tStack.pop();
1555               
1556                if (node->IsLeaf())
1557                {
1558                        //Debug << "l";
1559                        ExportBinLeaf(stream, static_cast<KdLeaf *>(node));
1560                }
1561                else
1562                {
1563                        //Debug << "i";
1564                        KdInterior *interior = static_cast<KdInterior *>(node);
1565
1566                        ExportBinInterior(stream, interior);
1567                       
1568                        tStack.push(interior->mFront);
1569                        tStack.push(interior->mBack);
1570                }
1571        }
1572
1573        return true;
1574}
1575
1576
1577KdNode *KdTree::ImportNextNode(IN_STREAM &stream,
1578                                                           KdInterior *parent,
1579                                                           const ObjectContainer &objects)
1580{
1581        int nodeType;
1582        stream.read(reinterpret_cast<char *>(&nodeType), sizeof(int));
1583
1584        if (nodeType == TYPE_LEAF)
1585                return ImportBinLeaf(stream, static_cast<KdInterior *>(parent), objects);
1586
1587        if (nodeType == TYPE_INTERIOR)
1588                return ImportBinInterior(stream, static_cast<KdInterior *>(parent));
1589
1590        Debug << "error! loading failed!" << endl;
1591        return NULL;
1592}
1593
1594
1595bool KdTree::ImportBinTree(const string &filename, ObjectContainer &objects)
1596{
1597        // export binary version of mesh
1598        queue<TraversalData> tStack;
1599        IN_STREAM stream(filename.c_str(), IN_BIN_MODE);
1600
1601        if (!stream.is_open()) return false;
1602
1603        // sort objects by their id
1604        //      if (!is_sorted(objects.begin(), objects.end(), ilt))
1605        sort(objects.begin(), objects.end(), ilt);
1606
1607        mBox.Initialize();
1608        ObjectContainer::const_iterator oit, oit_end = objects.end();
1609
1610        ///////////////////////////
1611        //-- compute bounding box of object space
1612
1613        for (oit = objects.begin(); oit != oit_end; ++ oit)
1614        {
1615                const AxisAlignedBox3 box = (*oit)->GetBox();
1616                mBox.Include(box);
1617        }
1618
1619        // hack: we make a new root
1620        DEL_PTR(mRoot);
1621
1622        mRoot = ImportNextNode(stream, NULL, objects);
1623
1624        tStack.push(TraversalData(mRoot, mBox, 0));
1625        mStat.Reset();
1626        mStat.nodes = 1;
1627
1628        while(!tStack.empty())
1629        {
1630                TraversalData tData = tStack.front();
1631                tStack.pop();
1632
1633                KdNode *node = tData.mNode;
1634
1635                if (!node->IsLeaf())
1636                {
1637                        mStat.nodes += 2;
1638
1639                        //Debug << "i" ;
1640                        KdInterior *interior = static_cast<KdInterior *>(node);
1641                        interior->mBox = tData.mBox;
1642
1643                        KdNode *front = ImportNextNode(stream, interior, objects);
1644                        KdNode *back = ImportNextNode(stream, interior, objects);
1645
1646                        interior->SetupChildLinks(back, front);
1647
1648                        ++ mStat.splits[interior->mAxis];
1649
1650                        // compute new bounding box
1651                        AxisAlignedBox3 frontBox, backBox;
1652
1653                        tData.mBox.Split(interior->mAxis,
1654                                interior->mPosition,
1655                                frontBox,
1656                                backBox);
1657
1658                        tStack.push(TraversalData(front, frontBox, tData.mDepth + 1));                 
1659                        tStack.push(TraversalData(back, backBox, tData.mDepth + 1));
1660                }
1661                else
1662                {
1663                        EvaluateLeafStats(tData);
1664                        //cout << "l";
1665                }
1666        }
1667
1668        float area = GetBox().SurfaceArea()*mKdPvsArea;
1669
1670        SetPvsTerminationNodes(area);
1671
1672        Debug << mStat << endl;
1673
1674        return true;
1675}
1676
1677
1678KdIntersectable *
1679KdTree::GetOrCreateKdIntersectable(KdNode *node)
1680{
1681        if (node == NULL)
1682                return NULL;
1683
1684        if (node->mIntersectable == NULL)
1685        {
1686                // not in map => create new entry
1687                const int id = (int)mKdIntersectables.size();
1688                KdIntersectable *kdObj = new KdIntersectable(node, GetBox(node));
1689                node->mIntersectable = kdObj;
1690
1691                mKdIntersectables.push_back(kdObj);
1692                kdObj->SetId(id);
1693#ifdef USE_BIT_PVS
1694                // hack: for kd pvs the kd intersecables are the pvs objects
1695                ObjectPvsIterator::sObjects.push_back(kdObj);
1696#endif
1697        }
1698
1699        return node->mIntersectable;
1700}
1701
1702
1703void
1704KdTree::SetPvsTerminationNodes(const float maxArea)
1705{
1706  stack<KdNode *> nodeStack;
1707 
1708  nodeStack.push(mRoot);
1709
1710  float area = 0.0f;
1711  float radius = 0.0f;
1712  int nodes = 0;
1713 
1714  while (!nodeStack.empty()) {
1715    KdNode *node = nodeStack.top();
1716    nodeStack.pop();
1717       
1718    node->mPvsTermination = 0;
1719    if (node->IsLeaf() || (GetSurfaceArea(node) <= maxArea) ) {
1720      area += GetSurfaceArea(node);
1721      radius += GetBox(node).Radius();
1722      nodes++;
1723      node->mPvsTermination = 1;
1724      // create dummy kd intersectable
1725      Intersectable *object = GetOrCreateKdIntersectable(node);
1726    } else {
1727      KdInterior *interior = (KdInterior *)node;
1728      nodeStack.push(interior->mFront);
1729      nodeStack.push(interior->mBack);
1730    }
1731  }
1732 
1733  if (nodes) {
1734    area /= nodes;
1735    radius /= nodes;
1736    cout<<"Number of nodes for storing in the PVS = "<<nodes<<endl;
1737    cout<<"Average rel. node area = "<<area/GetSurfaceArea(mRoot)<<endl;
1738    cout<<"Average rel. node radius = "<<radius/GetBox(mRoot).Radius()<<endl;
1739    cout<<"Avg node radius = "<<radius<<endl;
1740  }
1741 
1742}
1743
1744
1745/*KdNode *
1746KdTree::GetPvsNode(const Triangle3 &tri) const
1747{
1748  KdNode *node = mRoot;
1749 
1750  while (node->mPvsTermination == 0 )
1751  {
1752          KdInterior *inter = (KdInterior *)node;
1753          if (point[inter->mAxis] < inter->mPosition)
1754                  node = inter->mBack;
1755          else
1756                  node = inter->mFront;
1757  }
1758
1759  return node;
1760}*/
1761
1762
1763
1764KdNode *
1765KdTree::GetPvsNode(const Vector3 &point) const
1766{
1767  KdNode *node = mRoot;
1768 
1769  while (node->mPvsTermination == 0 ) {
1770    KdInterior *inter = (KdInterior *)node;
1771    if (point[inter->mAxis] < inter->mPosition)
1772      node = inter->mBack;
1773    else
1774      node = inter->mFront;
1775  }
1776 
1777  return node;
1778}
1779
1780KdNode *
1781KdTree::GetNode(const Vector3 &point,
1782                                const float maxArea) const
1783{
1784  KdNode *node = mRoot;
1785 
1786  while (!node->IsLeaf() && (GetSurfaceArea(node) > maxArea) ) {
1787    KdInterior *inter = (KdInterior *)node;
1788    if (point[inter->mAxis] < inter->mPosition)
1789      node = inter->mBack;
1790    else
1791      node = inter->mFront;
1792  }
1793 
1794  return node;
1795}
1796
1797
1798void KdTree::GetBoxIntersections(const AxisAlignedBox3 &box,
1799                                 vector<KdLeaf *> &leaves)
1800{
1801  stack<KdNode *> tStack;
1802 
1803  tStack.push(mRoot);
1804 
1805  while (!tStack.empty())
1806  {
1807          KdNode *node = tStack.top();
1808          tStack.pop();
1809
1810          if (node->IsLeaf())
1811          {
1812                  leaves.push_back(static_cast<KdLeaf *>(node));
1813          }
1814          else // interior
1815          {
1816                  KdInterior *interior = static_cast<KdInterior *>(node);
1817
1818                  if (box.Max(interior->mAxis) >= interior->mPosition)
1819                  {
1820                          tStack.push(interior->mFront);
1821                  }
1822
1823                  if (box.Min(interior->mAxis) < interior->mPosition)
1824                  {
1825                          tStack.push(interior->mBack);
1826                  }
1827          }
1828  }
1829}
1830
1831
1832
1833}
Note: See TracBrowser for help on using the repository browser.