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

RevLine 
[372]1#include <stack>
2#include <algorithm>
3#include <queue>
4#include "Environment.h"
5#include "Mesh.h"
6#include "KdTree.h"
[469]7#include "ViewCell.h"
[505]8#include "Beam.h"
[2116]9#include "ViewCell.h"
10#include "IntersectableWrapper.h"
[372]11
12
[2638]13
[2176]14using namespace std;
15
[372]16
[2332]17#define TYPE_INTERIOR -2
18#define TYPE_LEAF -3
19
20
[863]21namespace GtpVisibilityPreprocessor {
[860]22
[1176]23int KdNode::sMailId = 1;
24int KdNode::sReservedMailboxes = 1;
[2643]25int KdNode::sMailId2 = 1;
26int KdNode::sReservedMailboxes2 = 1;
[863]27
[2066]28
[1197]29inline static bool ilt(Intersectable *obj1, Intersectable *obj2)
30{
[2575]31  return obj1->mId < obj2->mId;
[1197]32}
33
34
[2539]35KdNode::KdNode(KdInterior *parent):
[2643]36mParent(parent), mMailbox(0), mIntersectable(NULL), mMailbox2(0)
[372]37{
[2575]38  if (parent)
39    mDepth = parent->mDepth+1;
40  else
41    mDepth = 0;
[372]42}
43
[1002]44
45KdInterior::~KdInterior()
46{
[2575]47  // recursivly destroy children
48  DEL_PTR(mFront);
49  DEL_PTR(mBack);
[1002]50}
51
52
[1999]53KdLeaf::~KdLeaf()
54{
[2575]55  DEL_PTR(mViewCell); 
[1999]56}
57
58
[372]59KdTree::KdTree()
[2539]60{
[372]61  mRoot = new KdLeaf(NULL, 0);
[2575]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);
[372]76
[2575]77  Environment::GetSingleton()->GetBoolValue("KdTree.sahUseFaces",
78                                            mSahUseFaces);
[372]79
80  char splitType[64];
[1004]81  Environment::GetSingleton()->GetStringValue("KdTree.splitMethod", splitType);
[372]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
[1002]99
100KdTree::~KdTree()
101{
[1999]102    DEL_PTR(mRoot);
103        CLEAR_CONTAINER(mKdIntersectables);
[1002]104}
105
106
[372]107bool
108KdTree::Construct()
109{
110  if (!splitCandidates)
[2003]111    splitCandidates = new vector<SortableEntry *>;
[372]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++) {
[752]123        //      cout<<(*mi)->GetBox()<<endl;
124        mBox.Include((*mi)->GetBox());
[372]125  }
126
[752]127  cout <<"KdTree Root Box:"<<mBox<<endl;
[372]128  mRoot = Subdivide(TraversalData(leaf, mBox, 0));
129
130  // remove the allocated array
[2003]131  CLEAR_CONTAINER(*splitCandidates);
[372]132  delete splitCandidates;
133
[2606]134  float area = GetPvsArea();
135 
[1989]136  SetPvsTerminationNodes(area);
137 
[372]138  return true;
139}
140
[2606]141float
142KdTree::GetPvsArea() const {
143  return GetBox().SurfaceArea() * mKdPvsArea;
144}
145
[372]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;
[2636]157 
[372]158  while (!tStack.empty()) {
[752]159        //      cout<<mStat.Nodes()<<" "<<mTermMaxNodes<<endl;
160        if (mStat.Nodes() > mTermMaxNodes) {
161          //    if ( GetMemUsage() > maxMemory ) {
[372]162      // count statistics on unprocessed leafs
163      while (!tStack.empty()) {
[752]164                EvaluateLeafStats(tStack.top());
165                tStack.pop();
[372]166      }
167      break;
168    }
[752]169         
170       
[372]171    TraversalData data = tStack.top();
172    tStack.pop();
173   
174    KdNode *node = SubdivideNode((KdLeaf *) data.mNode,
[2575]175                                 data.mBox,
176                                 backBox,
177                                 frontBox
178                                 );
[752]179
[2575]180    if (result == NULL)
[372]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     
[2575]189    }
190    else {
[372]191      EvaluateLeafStats(data);
192    }
193  }
194 
195  return result;
196
197}
198
199
200bool
201KdTree::TerminationCriteriaMet(const KdLeaf *leaf)
202{
[2575]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;
[372]211}
212
213
214int
215KdTree::SelectPlane(KdLeaf *leaf,
[2575]216                    const AxisAlignedBox3 &box,
217                    float &position
218                    )
[372]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;
[1584]232      bool mOnlyDrivingAxis = true;
[1387]233
[2575]234      if (mOnlyDrivingAxis) {
235        axis = box.Size().DrivingAxis();
236        costRatio = BestCostRatio(leaf,
237                                  box,
238                                  axis,
239                                  position,
240                                  objectsBack,
241                                  objectsFront);
[372]242      }
[2575]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      }
[372]261     
262      if (costRatio > mMaxCostRatio) {
[2575]263        //cout<<"Too big cost ratio "<<costRatio<<endl;
264        axis = -1;
[372]265      }
266      break;
267    }
268   
269    }
270  return axis;
271}
272
[2575]273KdNode*
[372]274KdTree::SubdivideNode(
275                      KdLeaf *leaf,
276                      const AxisAlignedBox3 &box,
277                      AxisAlignedBox3 &backBBox,
278                      AxisAlignedBox3 &frontBBox
279                      )
[2575]280
[372]281  if (TerminationCriteriaMet(leaf))
[2575]282    return leaf;
[469]283
[372]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 
[2575]293  mStat.nodes += 2;
[372]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);
[1074]329 
[372]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();
[1736]342       
343        // matt: no more ref
[1143]344        // for handling multiple objects: keep track of references
[1736]345        //if (leaf->IsRoot())
346        //      (*mi)->mReferences = 1; // initialise references at root
[1143]347       
[1736]348        // matt: no more ref
349        //-- (*mi)->mReferences; // remove parent ref
[372]350
[1143]351
[2575]352    if (box.Max(axis) >= position )
353    {
[372]354      front->mObjects.push_back(*mi);
[2575]355      //++ (*mi)->mReferences;
356    }
[1143]357       
[372]358    if (box.Min(axis) < position )
[2575]359    {
[372]360      back->mObjects.push_back(*mi);
[2575]361      // matt: no more ref
362      //  ++ (*mi)->mReferences;
363    }
[1143]364
[1144]365    mStat.objectRefs -= (int)leaf->mObjects.size();
366    mStat.objectRefs += objectsBack + objectsFront;
367  }
368
[2575]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;
[372]376}
377
378
[1143]379void KdTree::ProcessMultipleRefs(KdLeaf *leaf) const
380{
[2575]381  // find objects from multiple kd-leaves
382  ObjectContainer::const_iterator oit, oit_end = leaf->mObjects.end();
[1143]383
[2575]384  for (oit = leaf->mObjects.begin(); oit != oit_end; ++ oit)
385  {
386    Intersectable *object = *oit;
[1143]387               
[2575]388    // matt: no more ref
389    /*
390      if (object->mReferences > 1) {
391        leaf->mMultipleObjects.push_back(object);
392      }
393    */
394  }
[1143]395}
396
397
[372]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
[667]407  app << "#N_SPLITS ( Number of splits in axes x y z dx dy dz)\n";
[372]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
[1184]421  app << "#N_MAXOBJECTREFS  ( Max number of object refs / leaf )\n" <<
[372]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)
[469]473    mStat.maxObjectRefs = (int)leaf->mObjects.size();
[372]474 
475}
476
477
478
479void
[1233]480KdTree::SortSubdivisionCandidates(
[2575]481                                  KdLeaf *node,
482                                  const int axis
483                                  )
[372]484{
[2575]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 *>;
[372]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();
[2066]502      mi++)
503  {
[2575]504    AxisAlignedBox3 box = (*mi)->GetBox();
[372]505
[2575]506    splitCandidates->push_back(new SortableEntry(SortableEntry::BOX_MIN,
507                                                 box.Min(axis),
508                                                 *mi)
509                               );
[372]510   
[2575]511    splitCandidates->push_back(new SortableEntry(SortableEntry::BOX_MAX,
512                                                 box.Max(axis),
513                                                 *mi)
514                               );
[372]515  }
516 
[2066]517  stable_sort(splitCandidates->begin(), splitCandidates->end(), iltS);
[372]518}
519
520
521float
522KdTree::BestCostRatio(
[2575]523                      KdLeaf *node,
524                      const AxisAlignedBox3 &box,
525                      const int axis,
526                      float &position,
527                      int &objectsBack,
528                      int &objectsFront
529                      )
[372]530{
531
[1387]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)
[2575]548    costStream.open(filename);
[1387]549
550#endif
551 
[1233]552  SortSubdivisionCandidates(node, axis);
[372]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;
[2003]559  vector<SortableEntry *>::const_iterator ci;
[372]560
561  for(ci = splitCandidates->begin();
562      ci < splitCandidates->end();
563      ci++)
[2003]564    if ((*ci)->type == SortableEntry::BOX_MIN) {
565      totalIntersections += (*ci)->intersectable->IntersectionComplexity();
[372]566    }
567       
568  float intersectionsLeft = 0;
569  float intersectionsRight = totalIntersections;
570       
[469]571  int objectsLeft = 0, objectsRight = (int)node->mObjects.size();
[372]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 
[469]580  float minSum = 1e20f;
[372]581 
582  for(ci = splitCandidates->begin();
583      ci < splitCandidates->end();
584      ci++) {
[2003]585    switch ((*ci)->type) {
[372]586    case SortableEntry::BOX_MIN:
587      objectsLeft++;
[2003]588      intersectionsLeft += (*ci)->intersectable->IntersectionComplexity();
[372]589      break;
590    case SortableEntry::BOX_MAX:
591      objectsRight--;
[2003]592      intersectionsRight -= (*ci)->intersectable->IntersectionComplexity();
[372]593      break;
[2575]594    } // switch
[372]595
[2003]596    if ((*ci)->value > minBand && (*ci)->value < maxBand) {
[372]597      AxisAlignedBox3 lbox = box;
598      AxisAlignedBox3 rbox = box;
[2003]599      lbox.SetMax(axis, (*ci)->value);
600      rbox.SetMin(axis, (*ci)->value);
[372]601
602      float sum;
603      if (mSahUseFaces)
[2575]604        sum = intersectionsLeft*lbox.SurfaceArea() + intersectionsRight*rbox.SurfaceArea();
[372]605      else
[2575]606        sum = objectsLeft*lbox.SurfaceArea() + objectsRight*rbox.SurfaceArea();
[372]607     
608      //      cout<<"pos="<<(*ci).value<<"\t q=("<<ql<<","<<qr<<")\t r=("<<rl<<","<<rr<<")"<<endl;
609      //      cout<<"cost= "<<sum<<endl;
[1387]610
611#if DEBUG_COST
[2575]612      if (nodeId < 100) {
[1387]613        float oldCost = mSahUseFaces ? totalIntersections : node->mObjects.size();
614        float newCost = mCt_div_ci + sum/boxArea;
615        float ratio = newCost/oldCost;
[2003]616        costStream<<(*ci)->value<<" "<<ratio<<endl;
[2575]617      }
[1387]618#endif
[2575]619     
[372]620      if (sum < minSum) {
[2575]621        minSum = sum;
622        position = (*ci)->value;
623       
624        objectsBack = objectsLeft;
625        objectsFront = objectsRight;
[372]626      }
627    }
[2575]628  } // for ci
[372]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(
[2575]644                Ray &ray
645                )
[372]646{
[1984]647 
[372]648  int hits = 0;
649 
650  stack<RayTraversalData> tStack;
651 
652  float maxt = 1e6;
653  float mint = 0;
[1583]654
[1584]655  //  ray.ComputeInvertedDir();
[372]656  Intersectable::NewMail();
657
658  if (!mBox.GetMinMaxT(ray, &mint, &maxt))
659    return 0;
660
661  if (mint < 0)
662    mint = 0;
[1584]663
[372]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;
[752]674
[372]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) {
[2575]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        }
[372]693      }
694      else {
[2575]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      }
[372]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;
[2575]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);
[752]718         
[2575]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)
[1867]727                        ray.testedObjects.push_back(object);
728                 
[2575]729          static int oi=1;
730          if (MeshDebug)
731            cout<<"Object "<<oi++;
[752]732         
[2575]733          hits += object->CastRay(ray);
[752]734         
[2575]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      }
[752]743         
[2575]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    }
[372]763  }
764  return hits;
765}
766
[2618]767int KdTree::CastLineSegment2(const Vector3 &origin,
[2575]768                            const Vector3 &termination,
769                            ViewCellContainer &viewcells)
[469]770{
[2575]771  int hits = 0;
[2017]772
[2575]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;
[469]790
[2575]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)
[469]802        {
[2575]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;
[469]812        }
[2575]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;
[469]867}
868
[2618]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
[1974]1057void
1058KdTree::CollectKdObjects(const AxisAlignedBox3 &box,
[2575]1059                         ObjectContainer &objects
1060                        )
[1974]1061{
1062  stack<KdNode *> nodeStack;
1063 
1064  nodeStack.push(mRoot);
[469]1065
[1974]1066  while (!nodeStack.empty()) {
1067    KdNode *node = nodeStack.top();
1068    nodeStack.pop();
[2575]1069    if (node->IsLeaf() || node->mPvsTermination == 1)  {
1070      if (!node->Mailed()) {
[2638]1071                Intersectable *object = GetOrCreateKdIntersectable(node);
1072                node->Mail();
1073                objects.push_back(object);
[2575]1074      }
1075    }
1076    else {
[1974]1077      KdInterior *interior = (KdInterior *)node;
1078         
[2575]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);
[1974]1084    }
[2575]1085  } // while
[1974]1086}
1087
[372]1088void
[2606]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()) {
[2634]1106                  node->Mail();
[2606]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
[859]1125KdTree::CollectObjects(const AxisAlignedBox3 &box,
[2575]1126                       ObjectContainer &objects)
[859]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++) {
[2575]1138        Intersectable *object = leaf->mObjects[j];
1139        if (!object->Mailed() && Overlap(box, object->GetBox())) {
1140          object->Mail();
1141          objects.push_back(object);
1142        }
[859]1143      }
1144    } else {
1145      KdInterior *interior = (KdInterior *)node;
1146
[2575]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);
[859]1152    }
1153  }
1154}
1155
[2643]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
[859]1188void
[372]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++) {
[2575]1201        Intersectable *object = leaf->mObjects[j];
1202        if (!object->Mailed()) {
1203          object->Mail();
1204          objects.push_back(object);
1205        }
[372]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,
[2575]1218                           bool onlyUnmailed
1219                           )
[372]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);
[2575]1273    }
1274    else {
[372]1275      KdInterior *interior = (KdInterior *)node;
1276      if (interior->mPosition > box.Max(interior->mAxis))
[2575]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      }
[372]1287    }
1288  }
1289 
[469]1290  return (int)neighbors.size();
[372]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
[469]1352void
1353KdTree::CreateAndCollectViewCells(ViewCellContainer &vc) const
1354{
1355  stack<KdNode *> nodeStack;
1356  nodeStack.push(mRoot);
[372]1357
[469]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
[471]1364          KdViewCell *viewCell = new KdViewCell();
1365          leaf->mViewCell = viewCell;
1366          // push back pointer to this leaf
[1551]1367          viewCell->mLeaves[0] = leaf;
[471]1368      vc.push_back(viewCell);
[469]1369    } else {
1370      KdInterior *interior = (KdInterior *)node;
1371      nodeStack.push(interior->mBack);
1372      nodeStack.push(interior->mFront);
1373    }
1374  }
1375}
1376
[1736]1377
[372]1378KdNode *
1379KdTree::GetRandomLeaf(const bool onlyUnmailed)
1380{
[507]1381  stack<KdNode *> nodeStack;
[372]1382  nodeStack.push(mRoot);
[507]1383 
1384  int mask = rand();
[372]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}
[504]1404
1405
[2539]1406int KdTree::CastBeam(Beam &beam)
[504]1407{
[507]1408  stack<KdNode *> nodeStack;
1409  nodeStack.push(mRoot);
[504]1410 
[507]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                }
[504]1433
[507]1434                nodeStack.push(first);
1435                nodeStack.push(second);
1436          }
1437          break;
1438          // default: cull
1439        }
1440  }
1441
[532]1442  if (beam.mFlags & Beam::STORE_OBJECTS)
1443  {
1444          vector<KdNode *>::const_iterator it, it_end = beam.mKdNodes.end();
[535]1445       
[532]1446          Intersectable::NewMail();
1447          for (it = beam.mKdNodes.begin(); it != it_end; ++ it)
1448          {
1449                  CollectObjects(*it, beam.mObjects);
1450          }
1451  }
1452
[837]1453  return (int)beam.mKdNodes.size();
[504]1454}
[860]1455
[1074]1456
[1201]1457void KdTree::ExportBinLeaf(OUT_STREAM &stream, KdLeaf *leaf)
[1194]1458{
[2575]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));
[1194]1474    }
[878]1475}
[1194]1476
1477
[1201]1478KdLeaf *KdTree::ImportBinLeaf(IN_STREAM &stream,
[2575]1479                              KdInterior *parent,
1480                              const ObjectContainer &objects)
[1194]1481{
[2575]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);
[1196]1488
[2575]1489  MeshInstance dummyInst(NULL);
[1633]1490       
[2575]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;
[1194]1508}
1509
1510
[1201]1511void KdTree::ExportBinInterior(OUT_STREAM &stream, KdInterior *interior)
[1194]1512{
[1196]1513        int interiorid = TYPE_INTERIOR;
[1194]1514        stream.write(reinterpret_cast<char *>(&interiorid), sizeof(int));
1515
1516        int axis = interior->mAxis;
[1196]1517        float pos = interior->mPosition;
[1194]1518
1519        stream.write(reinterpret_cast<char *>(&axis), sizeof(int));
[1196]1520        stream.write(reinterpret_cast<char *>(&pos), sizeof(float));
[1194]1521}
1522
1523
[1201]1524KdInterior *KdTree::ImportBinInterior(IN_STREAM  &stream, KdInterior *parent)
[1194]1525{
[1196]1526        KdInterior *interior = new KdInterior(parent);
[1194]1527
[2539]1528        int axis;
1529        float pos;
[1194]1530
1531        stream.read(reinterpret_cast<char *>(&axis), sizeof(int));
[1196]1532        stream.read(reinterpret_cast<char *>(&pos), sizeof(float));
1533
1534        interior->mAxis = axis;
1535        interior->mPosition = pos;
1536
1537        return interior;
[1194]1538}
1539
1540
1541bool KdTree::ExportBinTree(const string &filename)
1542{
[1201]1543        OUT_STREAM stream(filename.c_str(), OUT_BIN_MODE);
[1194]1544       
[2332]1545        if (!stream.is_open()) return false;
[1194]1546
1547        // export binary version of mesh
[1197]1548        queue<KdNode *> tStack;
[1194]1549        tStack.push(mRoot);
1550
1551        while(!tStack.empty())
1552        {
[1197]1553                KdNode *node = tStack.front();
[1196]1554                tStack.pop();
1555               
1556                if (node->IsLeaf())
1557                {
[1197]1558                        //Debug << "l";
[2017]1559                        ExportBinLeaf(stream, static_cast<KdLeaf *>(node));
[1196]1560                }
1561                else
1562                {
[1197]1563                        //Debug << "i";
[2017]1564                        KdInterior *interior = static_cast<KdInterior *>(node);
[1194]1565
[1196]1566                        ExportBinInterior(stream, interior);
1567                       
1568                        tStack.push(interior->mFront);
1569                        tStack.push(interior->mBack);
[1194]1570                }
1571        }
1572
1573        return true;
1574}
1575
1576
[2539]1577KdNode *KdTree::ImportNextNode(IN_STREAM &stream,
1578                                                           KdInterior *parent,
1579                                                           const ObjectContainer &objects)
[1194]1580{
[1196]1581        int nodeType;
1582        stream.read(reinterpret_cast<char *>(&nodeType), sizeof(int));
1583
1584        if (nodeType == TYPE_LEAF)
[2017]1585                return ImportBinLeaf(stream, static_cast<KdInterior *>(parent), objects);
[1196]1586
1587        if (nodeType == TYPE_INTERIOR)
[2017]1588                return ImportBinInterior(stream, static_cast<KdInterior *>(parent));
[1196]1589
1590        Debug << "error! loading failed!" << endl;
1591        return NULL;
1592}
1593
1594
[2539]1595bool KdTree::ImportBinTree(const string &filename, ObjectContainer &objects)
[1196]1596{
[2615]1597        // export binary version of mesh
1598        queue<TraversalData> tStack;
1599        IN_STREAM stream(filename.c_str(), IN_BIN_MODE);
[1197]1600
[2615]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;
[1194]1675}
1676
[1999]1677
[1594]1678KdIntersectable *
1679KdTree::GetOrCreateKdIntersectable(KdNode *node)
1680{
[2615]1681        if (node == NULL)
1682                return NULL;
[1194]1683
[2116]1684        if (node->mIntersectable == NULL)
1685        {
1686                // not in map => create new entry
[2117]1687                const int id = (int)mKdIntersectables.size();
[2571]1688                KdIntersectable *kdObj = new KdIntersectable(node, GetBox(node));
[2117]1689                node->mIntersectable = kdObj;
[2615]1690
[2117]1691                mKdIntersectables.push_back(kdObj);
1692                kdObj->SetId(id);
1693#ifdef USE_BIT_PVS
[2615]1694                // hack: for kd pvs the kd intersecables are the pvs objects
1695                ObjectPvsIterator::sObjects.push_back(kdObj);
[2117]1696#endif
[2615]1697        }
[1761]1698
[2615]1699        return node->mIntersectable;
[1387]1700}
[1594]1701
[1989]1702
1703void
[2575]1704KdTree::SetPvsTerminationNodes(const float maxArea)
[1989]1705{
1706  stack<KdNode *> nodeStack;
1707 
1708  nodeStack.push(mRoot);
1709
[1995]1710  float area = 0.0f;
1711  float radius = 0.0f;
1712  int nodes = 0;
1713 
[1989]1714  while (!nodeStack.empty()) {
1715    KdNode *node = nodeStack.top();
[2575]1716    nodeStack.pop();
[2606]1717       
[2575]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);
[1989]1730    }
1731  }
[2575]1732 
[1995]1733  if (nodes) {
[2575]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;
[1995]1740  }
1741 
[1989]1742}
1743
[2615]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
[1594]1764KdNode *
[1989]1765KdTree::GetPvsNode(const Vector3 &point) const
1766{
1767  KdNode *node = mRoot;
1768 
1769  while (node->mPvsTermination == 0 ) {
[2575]1770    KdInterior *inter = (KdInterior *)node;
1771    if (point[inter->mAxis] < inter->mPosition)
1772      node = inter->mBack;
1773    else
1774      node = inter->mFront;
[1989]1775  }
1776 
1777  return node;
1778}
1779
1780KdNode *
[1594]1781KdTree::GetNode(const Vector3 &point,
1782                                const float maxArea) const
1783{
1784  KdNode *node = mRoot;
1785 
1786  while (!node->IsLeaf() && (GetSurfaceArea(node) > maxArea) ) {
[2575]1787    KdInterior *inter = (KdInterior *)node;
1788    if (point[inter->mAxis] < inter->mPosition)
1789      node = inter->mBack;
1790    else
1791      node = inter->mFront;
[1594]1792  }
1793 
1794  return node;
1795}
1796
1797
[1999]1798void KdTree::GetBoxIntersections(const AxisAlignedBox3 &box,
[2575]1799                                 vector<KdLeaf *> &leaves)
[1999]1800{
[2575]1801  stack<KdNode *> tStack;
1802 
1803  tStack.push(mRoot);
1804 
1805  while (!tStack.empty())
[2615]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  }
[1594]1829}
[1633]1830
[1999]1831
1832
[1633]1833}
Note: See TracBrowser for help on using the repository browser.