source: trunk/VUT/GtpVisibilityPreprocessor/src/VspKdTree.cpp @ 418

Revision 418, 36.4 KB checked in by mattausch, 19 years ago (diff)
Line 
1
2// ================================================================
3// $Id: lsds_kdtree.cpp,v 1.18 2005/04/16 09:34:21 bittner Exp $
4// ****************************************************************
5/**
6   The KD tree based LSDS
7*/
8// Initial coding by
9/**
10   @author Jiri Bittner
11*/
12
13// Standard headers
14#include <stack>
15#include <queue>
16#include <algorithm>
17#include <fstream>
18#include <string>
19
20#include "VspKdTree.h"
21
22#include "Environment.h"
23#include "VssRay.h"
24#include "Intersectable.h"
25#include "Ray.h"
26
27// Static variables
28int VspKdTreeLeaf::mailID = 0;
29int VspKdTree::sTermMaxDepth = 10;
30float VspKdTree::sTermMinSize = 0.1;
31int VspKdTree::sTermMinPvs = 10;
32int VspKdTree::sTermMinRays = 20;
33float VspKdTree::sTermMaxCostRatio = 0.1;
34float VspKdTree::sTermMaxRayContribution = 0.1;
35
36
37
38inline void AddObject2Pvs(Intersectable *object,
39                                                  const int side,
40                                                  int &pvsBack,
41                                                  int &pvsFront)
42{
43        if (!object)
44                return;
45               
46        if (side <= 0)
47        {
48                if (!object->Mailed() && !object->Mailed(2))
49                {
50                        ++ pvsBack;
51
52                        if (object->Mailed(1))
53                                object->Mail(2);
54                        else
55                                object->Mail();
56                }
57        }
58       
59        if (side >= 0)
60        {
61                if (!object->Mailed(1) && !object->Mailed(2))
62                {
63                        ++ pvsFront;
64                        if (object->Mailed())
65                                object->Mail(2);
66                        else
67                                object->Mail(1);
68                }
69        }
70}
71
72/**************************************************************/
73 *                VspKdTreeInterior implementation            *
74/**************************************************************/
75
76VspKdTreeInterior::VspKdTreeInterior(VspKdTreeInterior *p):
77VspKdTreeNode(p), mBack(NULL), mFront(NULL), mAccesses(0), mLastAccessTime(-1)
78{
79}
80
81int VspKdTreeInterior::GetAccessTime()
82{
83        return lastAccessTime;
84}
85
86void VspKdTreeInterior::SetupChildLinks(VspKdTreeNode *b, VspKdTreeNode *f)
87{
88        mBack = b;
89        mFront = f;
90        b->mParent = f->mParent = this;
91}
92
93void VspKdTreeInterior::ReplaceChildLink(VspKdTreeNode *oldChild, VspKdTreeNode *newChild)
94{
95        if (mBack == oldChild)
96                mBack = newChild;
97        else
98                mFront = newChild;
99}
100
101int VspKdTreeInterior::Type() const 
102{
103        return VspKdInterior;
104}
105 
106VspKdTreeInterior::~VspKdTreeInterior()
107{
108        DEL_PTR(mBack);
109        DEL_PTR(mFront);
110}
111 
112void VspKdTreeInterior::Print(ostream &s) const
113{
114        if (axis == 0)
115                s << "x ";
116        else
117                if (axis == 1)
118                        s << "y ";
119                else
120                        s << "z ";
121        s << position << " ";
122
123        mBack->Print(s);
124        mFront->Print(s);
125}
126       
127int VspKdTreeInterior::ComputeRayIntersection(const RayInfo &rayData, float &t)
128{
129        return rayData.ComputeRayIntersection(axis, position, t);
130}
131
132
133// Constructor
134VspKdTree::VspKdTree()
135{
136        environment->GetIntValue("VspKdTree.maxDepth", sTermMaxDepth);
137        environment->GetIntValue("VspKdTree.minPvs", sTermMinPvs);
138        environment->GetIntValue("VspKdTree.minRays", sTermMinRays);
139        environment->GetFloatValue("VspKdTree.maxRayContribution", sTermMaxRayContribution);
140        environment->GetFloatValue("VspKdTree.maxCostRatio", sTermMaxCostRatio);
141
142        environment->GetFloatValue("VspKdTree.minSize", sTermMinSize);
143        termMinSize = sqr(termMinSize);
144         
145        environment->GetFloatValue("VspKdTree.epsilon", epsilon);
146        environment->GetFloatValue("VspKdTree.ct_div_ci", ct_div_ci);
147       
148        environment->GetFloatValue("VspKdTree.maxTotalMemory", maxTotalMemory);
149        environment->GetFloatValue("VspKdTree.maxStaticMemory", maxStaticMemory);
150 
151 
152        float refDirAngle;
153        environment->GetIntValue("VspKdTree.accessTimeThreshold", accessTimeThreshold);
154        //= 1000;
155        environment->GetIntValue("VspKdTree.minCollapseDepth", minCollapseDepth);
156 
157
158        // split type
159        char sname[128];
160        environment->GetStringValue("VspKdTree.splitType", sname);
161        string name(sname);
162       
163        if (name.compare("regular") == 0)
164                splitType = ESplitRegular;
165        else
166        {
167                if (name.compare("heuristic") == 0)
168                        splitType = ESplitHeuristic;
169                else
170                {
171                        cerr << "Invalid VspKdTree split type " << name << endl;
172                        exit(1);
173                }
174        }
175
176        environment->GetBoolValue("VspKdTree.randomize", randomize);
177
178        mRoot = NULL;
179
180        splitCandidates = new vector<SortableEntry>;
181}
182
183
184VspKdTree::~VspKdTree()
185{
186        if (mRoot)
187                delete mRoot;
188}
189
190void VspKdStatistics::Print(ostream &app) const
191{
192  app << "===== VspKdTree statistics ===============\n";
193
194  app << "#N_RAYS ( Number of rays )\n"
195      << rays <<endl;
196 
197  app << "#N_INITPVS ( Initial PVS size )\n"
198      << initialPvsSize <<endl;
199 
200  app << "#N_NODES ( Number of nodes )\n" << nodes << "\n";
201
202  app << "#N_LEAVES ( Number of leaves )\n" << Leaves() << "\n";
203
204  app << "#N_SPLITS ( Number of splits in axes x y z dx dy dz \n";
205  for (int i=0; i<7; i++)
206    app << splits[i] <<" ";
207  app <<endl;
208
209  app << "#N_RAYREFS ( Number of rayRefs )\n" <<
210    rayRefs << "\n";
211
212  app << "#N_RAYRAYREFS  ( Number of rayRefs / ray )\n" <<
213    rayRefs/(double)rays << "\n";
214
215  app << "#N_LEAFRAYREFS  ( Number of rayRefs / leaf )\n" <<
216    rayRefs/(double)Leaves() << "\n";
217
218  app << "#N_MAXRAYREFS  ( Max number of rayRefs / leaf )\n" <<
219    maxRayRefs << "\n";
220
221
222  //  app << setprecision(4);
223
224  app << "#N_PMAXDEPTHLEAVES ( Percentage of leaves at maxdepth )\n"<<
225    maxDepthNodes*100/(double)Leaves()<<endl;
226
227  app << "#N_PMINPVSLEAVES  ( Percentage of leaves with minCost )\n"<<
228    minPvsNodes*100/(double)Leaves()<<endl;
229
230  app << "#N_PMINRAYSLEAVES  ( Percentage of leaves with minCost )\n"<<
231    minRaysNodes*100/(double)Leaves()<<endl;
232       
233        app << "#N_PMINSIZELEAVES  ( Percentage of leaves with minSize )\n"<<
234    minSizeNodes*100/(double)Leaves()<<endl;
235
236        app << "#N_PMAXRAYCONTRIBLEAVES  ( Percentage of leaves with maximal ray contribution )\n"<<
237    maxRayContribNodes*100/(double)Leaves()<<endl;
238
239  app << "#N_ADDED_RAYREFS  (Number of dynamically added ray references )\n"<<
240    addedRayRefs<<endl;
241
242  app << "#N_REMOVED_RAYREFS  (Number of dynamically removed ray references )\n"<<
243    removedRayRefs<<endl;
244
245  //  app << setprecision(4);
246
247  app << "#N_CTIME  ( Construction time [s] )\n"
248      << Time() << " \n";
249
250  app << "===== END OF VspKdTree statistics ==========\n";
251
252}
253
254
255void
256VspKdTreeLeaf::UpdatePvsSize()
257{
258        if (!mValidPvs)
259        {
260                Intersectable::NewMail();
261                int pvsSize = 0;
262                for(VspKdTreeNode::RayInfoContainer::iterator ri = rays.begin();
263                        ri != rays.end(); ++ ri)
264                {
265                        if ((*ri).mRay->IsActive())
266                        {
267                                Intersectable *object;
268#if BIDIRECTIONAL_RAY
269                                object = (*ri).mRay->mOriginObject;
270                       
271                                if (object && !object->Mailed())
272                                {
273                                        ++ pvsSize;
274                                        object->Mail();
275                                }
276#endif
277                                object = (*ri).mRay->mTerminationObject;
278                                if (object && !object->Mailed())
279                                {
280                                        ++ pvsSize;
281                                        object->Mail();
282                                }
283                        }
284                }
285
286                mPvsSize = pvsSize;
287                mValidPvs = true;
288        }
289}
290
291
292void
293VspKdTree::Construct(VssRayContainer &rays,
294                                         AxisAlignedBox3 *forcedBoundingBox)
295{
296        mStat.Start();
297 
298        maxMemory = maxStaticMemory;
299
300        if (root)
301                delete root;
302
303        root = new VspKdTreeLeaf(NULL, rays.size());
304
305        // first construct a leaf that will get subdivide
306        VspKdTreeLeaf *leaf = (VspKdTreeLeaf *) root;
307
308       
309        mStat.nodes = 1;
310       
311        bbox.Initialize();
312        dirBBox.Initialize();
313 
314       
315        for (VssRayContainer::const_iterator ri = rays.begin(); ri != rays.end(); ++ ri)
316        {
317                leaf->AddRay(VspKdTreeNode::RayInfo(*ri));
318
319                bbox.Include((*ri)->GetOrigin());
320                bbox.Include((*ri)->GetTermination());   
321               
322                dirBBox.Include(Vector3((*ri)->GetDirParametrization(0), (*ri)->GetDirParametrization(1), 0));
323        }
324       
325        if (forcedBoundingBox)
326                bbox = *forcedBoundingBox;
327
328        cout<<"Bbox = "<<bbox<<endl;
329        cout<<"Dirr Bbox = "<<dirBBox<<endl;
330
331        mStat.rays = leaf->rays.size();
332        leaf->UpdatePvsSize();
333       
334        mStat.initialPvsSize = leaf->GetPvsSize();
335       
336        // Subdivide();
337        root = Subdivide(TraversalData(leaf, bbox, 0));
338
339        if (splitCandidates)
340        {
341                // force realease of this vector
342                delete splitCandidates;
343                splitCandidates = new vector<SortableEntry>;
344        }
345 
346        mStat.Stop();
347
348        mStat.Print(cout);
349        cout<<"#Total memory=" << GetMemUsage() << endl;
350}
351
352
353
354VspKdTreeNode *VspKdTree::Subdivide(const TraversalData &tdata)
355{
356        VspKdTreeNode *result = NULL;
357
358        priority_queue<TraversalData> tStack;
359        //  stack<TraversalData> tStack;
360 
361        tStack.push(tdata);
362
363        AxisAlignedBox3 backBox;
364        AxisAlignedBox3 frontBox;
365 
366        int lastMem = 0;
367        while (!tStack.empty())
368        {
369                float mem = GetMemUsage();
370               
371                if ( lastMem/10 != ((int)mem)/10)
372                {
373                        cout << mem << " MB" << endl;
374                }
375                lastMem = (int)mem;
376               
377                if (  mem > maxMemory )
378                {
379                        // count statistics on unprocessed leafs
380                        while (!tStack.empty())
381                        {
382                                EvaluateLeafStats(tStack.top());
383                                tStack.pop();
384                        }
385                        break;
386                }
387   
388                TraversalData data = tStack.top();
389                tStack.pop();
390   
391               
392                VspKdTreeNode *node = SubdivideNode((VspKdTreeLeaf *) data.node,
393                                                                                        data.bbox, backBox,     frontBox);
394                if (result == NULL)
395                        result = node;
396   
397                if (!node->IsLeaf())
398                {
399                        VspKdTreeInterior *interior = (VspKdTreeInterior *) node;
400
401                        // push the children on the stack
402                        tStack.push(TraversalData(interior->back, backBox, data.depth+1));
403                        tStack.push(TraversalData(interior->front, frontBox, data.depth+1));
404                }
405                else
406                {
407                        EvaluateLeafStats(data);
408                }
409        }
410
411        return result;
412}
413
414
415// returns selected plane for subdivision
416int
417VspKdTree::SelectPlane(VspKdTreeLeaf *leaf,
418                                           const AxisAlignedBox3 &box,
419                                           float &position,
420                                           int &raysBack,
421                                           int &raysFront,
422                                           int &pvsBack,
423                                           int &pvsFront)
424{
425        int minDirDepth = 6;
426        int axis;
427        float costRatio;
428       
429        if (splitType == ESplitRegular) {
430                costRatio = BestCostRatioRegular(leaf,
431                                                                                 axis,
432                                                                                 position,
433                                                                                 raysBack,
434                                                                                 raysFront,
435                                                                                 pvsBack,
436                                                                                 pvsFront);
437
438        }
439        else
440        {
441                if (splitType == ESplitHeuristic)
442                        costRatio = BestCostRatioHeuristic(leaf,
443                                                                                           axis,       
444                                                                                           position,
445                                                                                           raysBack,
446                                                                                           raysFront,
447                                                                                           pvsBack,
448                                                                                           pvsFront);
449                else {
450                        cerr << "VspKdTree: Unknown split heuristics\n";
451                        exit(1);
452                }
453  }
454
455        if (costRatio > termMaxCostRatio)
456        {
457                //              cout<<"Too big cost ratio "<<costRatio<<endl;
458                return -1;
459        }
460
461#if 0   
462        cout<<
463                "pvs="<<leaf->mPvsSize<<
464                " rays="<<leaf->rays.size()<<
465                " rc="<<leaf->GetAvgRayContribution()<<
466                " axis="<<axis<<endl;
467#endif
468       
469        return axis;
470}
471
472
473                                                       
474
475float
476VspKdTree::EvalCostRatio(VspKdTreeLeaf *leaf,
477                                                 const int axis,
478                                                 const float position,
479                                                 int &raysBack,
480                                                 int &raysFront,
481                                                 int &pvsBack,
482                                                 int &pvsFront)
483{
484        raysBack = 0;
485        raysFront = 0;
486        pvsFront = 0;
487        pvsBack = 0;
488
489        float newCost;
490
491        Intersectable::NewMail(3);
492       
493        // eval pvs size
494        int pvsSize = leaf->GetPvsSize();
495
496        Intersectable::NewMail(3);
497
498        if (axis <= VspKdTreeNode::SPLIT_Z) {
499    // this is the main ray classification loop!
500    for(VspKdTreeNode::RayInfoContainer::iterator ri = leaf->rays.begin();
501                ri != leaf->rays.end(); ri++)
502                if ((*ri).mRay->IsActive())
503                {
504                        // determine the side of this ray with respect to the plane
505                        int side = (*ri).ComputeRayIntersection(axis, position, (*ri).mRay->mT);
506                        //                              (*ri).mRay->mSide = side;
507                       
508                        if (side <= 0)
509                                raysBack ++;
510                               
511                        if (side >= 0)
512                                raysFront++;
513                               
514                        AddObject2Pvs((*ri).mRay->mTerminationObject, side, pvsBack, pvsFront);
515                }
516               
517                AxisAlignedBox3 box = GetBBox(leaf);
518               
519                float minBox = box.Min(axis);
520                float maxBox = box.Max(axis);
521                float sizeBox = maxBox - minBox;
522               
523                //              float sum = raysBack*(position - minBox) + raysFront*(maxBox - position);
524                float sum = pvsBack*(position - minBox) + pvsFront*(maxBox - position);
525
526                newCost = ct_div_ci + sum/sizeBox;
527        }
528        else
529        {
530                // directional split
531                for(VspKdTreeNode::RayInfoContainer::iterator ri = leaf->rays.begin();
532                        ri != leaf->rays.end(); ++ ri)
533                        if ((*ri).mRay->IsActive())
534                        {
535                                // determine the side of this ray with respect to the plane
536                                int side;
537                                if ((*ri).mRay->GetDirParametrization(axis - 3) > position)
538                                        side = 1;
539                                else
540                                        side = -1;
541
542                                if (side <= 0)
543                                        raysBack++;
544                               
545                                if (side >= 0)
546                                        raysFront++;
547
548                                //                              (*ri).mRay->mSide = side;
549                                AddObject2Pvs((*ri).mRay->mTerminationObject, side, pvsBack, pvsFront);
550                        }
551       
552                        AxisAlignedBox3 box = GetDirBBox(leaf);
553               
554                        float minBox = box.Min(axis - 3);
555                        float maxBox = box.Max(axis - 3);
556                        float sizeBox = maxBox - minBox;
557               
558                        // float sum = raysBack*(position - minBox) + raysFront*(maxBox - position);
559                        float sum =     pvsBack*(position - minBox) + pvsFront*(maxBox - position);
560                //              float sum = pvsBack + pvsFront;
561                newCost = ct_div_ci + sum/sizeBox;
562        }
563
564        //      cout<<axis<<" "<<pvsSize<<" "<<pvsBack<<" "<<pvsFront<<endl;
565        //  float oldCost = leaf->rays.size();
566        float oldCost = pvsSize;
567       
568        float ratio = newCost / oldCost;
569        //      cout<<"ratio="<<ratio<<endl;
570        return ratio;
571}
572
573float
574VspKdTree::BestCostRatioRegular(VspKdTreeLeaf *leaf,
575                                                                int &axis,
576                                                                float &position,
577                                                                int &raysBack,
578                                                                int &raysFront,
579                                                                int &pvsBack,
580                                                                int &pvsFront)
581{
582        int nRaysBack[6], nRaysFront[6];
583        int nPvsBack[6], nPvsFront[6];
584        float nPosition[6];
585        float nCostRatio[6];
586        int bestAxis = -1;
587       
588        AxisAlignedBox3 sBox = GetBBox(leaf);
589        AxisAlignedBox3 dBox = GetDirBBox(leaf);
590
591        // int sAxis = box.Size().DrivingAxis();
592        int sAxis = sBox.Size().DrivingAxis();
593        int dAxis = dBox.Size().DrivingAxis() + 3;
594
595        bool onlyDrivingAxis = true;
596
597        for (axis = 0; axis < 6; ++ axis)
598        {
599                if (!onlyDrivingAxis || axis == sAxis || axis == dAxis)
600                {
601                        if (axis < 3)
602                                nPosition[axis] = (sBox.Min()[axis] + sBox.Max()[axis])*0.5f;
603                        else
604                                nPosition[axis] = (dBox.Min()[axis-3] + dBox.Max()[axis-3])*0.5f;
605                       
606                        nCostRatio[axis] = EvalCostRatio(leaf,
607                                                                                         axis,
608                                                                                         nPosition[axis],
609                                                                                         nRaysBack[axis],
610                                                                                         nRaysFront[axis],
611                                                                                         nPvsBack[axis],
612                                                                                         nPvsFront[axis]);
613                       
614                        if (bestAxis == -1)
615                                bestAxis = axis;
616                        else
617                                if (nCostRatio[axis] < nCostRatio[bestAxis])
618                                        bestAxis = axis;
619                }
620        }
621
622        axis = bestAxis;
623        position = nPosition[bestAxis];
624
625        raysBack = nRaysBack[bestAxis];
626        raysFront = nRaysFront[bestAxis];
627
628        pvsBack = nPvsBack[bestAxis];
629        pvsFront = nPvsFront[bestAxis];
630       
631        return nCostRatio[bestAxis];
632}
633
634float VspKdTree::BestCostRatioHeuristic(VspKdTreeLeaf *leaf,
635                                                                                int &axis,
636                                                                                float &position,
637                                                                                int &raysBack,
638                                                                                int &raysFront,
639                                                                                int &pvsBack,
640                                                                                int &pvsFront)
641{
642        AxisAlignedBox3 box = GetBBox(leaf);
643        //      AxisAlignedBox3 dirBox = GetDirBBox(node);
644       
645        axis = box.Size().DrivingAxis();
646               
647        SortSplitCandidates(leaf, axis);
648 
649        // go through the lists, count the number of objects left and right
650        // and evaluate the following cost funcion:
651        // C = ct_div_ci  + (ql*rl + qr*rr)/queries
652       
653        int rl=0, rr = leaf->rays.size();
654        int pl=0, pr = leaf->GetPvsSize();
655        float minBox = box.Min(axis);
656        float maxBox = box.Max(axis);
657        float sizeBox = maxBox - minBox;
658       
659        float minBand = minBox + 0.1*(maxBox - minBox);
660        float maxBand = minBox + 0.9*(maxBox - minBox);
661       
662        float sum = rr*sizeBox;
663        float minSum = 1e20f;
664       
665        Intersectable::NewMail();
666        // set all object as belonging to the fron pvs
667        for(VspKdTreeNode::RayInfoContainer::iterator ri = leaf->rays.begin();
668                ri != leaf->rays.end(); ++ ri)
669        {
670                if ((*ri).mRay->IsActive())
671                {
672                        Intersectable *object = (*ri).mRay->mTerminationObject;
673                       
674                        if (object)
675                                if (!object->Mailed())
676                                {
677                                        object->Mail();
678                                        object->mCounter = 1;
679                                }
680                                else
681                                        ++ object->mCounter;
682                }
683        }
684       
685        Intersectable::NewMail();
686       
687        for (vector<SortableEntry>::const_iterator ci = splitCandidates->begin();
688                ci < splitCandidates->end(); ++ ci)
689        {
690                VssRay *ray;
691                       
692                switch ((*ci).type)
693                {
694                        case SortableEntry::ERayMin:
695                                {
696                                        ++ rl;
697                                        ray = (VssRay *) (*ci).data;
698                                        Intersectable *object = ray->mTerminationObject;
699                                        if (object && !object->Mailed())
700                                        {
701                                                object->Mail();
702                                                ++ pl;
703                                        }
704                                        break;
705                                }
706                        case SortableEntry::ERayMax:
707                                {
708                                        -- rr;
709                                       
710                                        ray = (VssRay *) (*ci).data;
711                                        Intersectable *object = ray->mTerminationObject;
712                                       
713                                        if (object)
714                                        {
715                                                if (-- object->mCounter == 0)
716                                                -- pr;
717                                        }
718                                               
719                                        break;
720                                }
721                }
722                       
723                if ((*ci).value > minBand && (*ci).value < maxBand)
724                {
725                        sum = pl*((*ci).value - minBox) + pr*(maxBox - (*ci).value);
726               
727                        //  cout<<"pos="<<(*ci).value<<"\t q=("<<ql<<","<<qr<<")\t r=("<<rl<<","<<rr<<")"<<endl;
728                        // cout<<"cost= "<<sum<<endl;
729                       
730                        if (sum < minSum)
731                        {
732                                minSum = sum;
733                                position = (*ci).value;
734                       
735                                raysBack = rl;
736                                raysFront = rr;
737                       
738                                pvsBack = pl;
739                                pvsFront = pr;
740                               
741                        }
742                }
743        }
744
745        float oldCost = leaf->GetPvsSize();
746        float newCost = ct_div_ci + minSum / sizeBox;
747        float ratio = newCost / oldCost;
748 
749        //  cout<<"===================="<<endl;
750        //  cout<<"costRatio="<<ratio<<" pos="<<position<<" t="<<(position - minBox)/(maxBox - minBox)
751        //      <<"\t q=("<<queriesBack<<","<<queriesFront<<")\t r=("<<raysBack<<","<<raysFront<<")"<<endl;
752       
753        return ratio;
754}
755
756void VspKdTree::SortSplitCandidates(VspKdTreeLeaf *node,
757                                                                        const int axis)
758{
759        splitCandidates->clear();
760 
761        int requestedSize = 2 * (int)(node->rays.size());
762        // creates a sorted split candidates array
763        if (splitCandidates->capacity() > 500000 &&
764                requestedSize < (int)(splitCandidates->capacity()/10) )
765        {
766        delete splitCandidates;
767                splitCandidates = new vector<SortableEntry>;
768        }
769 
770        splitCandidates->reserve(requestedSize);
771
772        // insert all queries
773        for(VspKdTreeNode::RayInfoContainer::const_iterator ri = node->rays.begin();
774                ri < node->rays.end(); ++ ri)
775        {
776                bool positive = (*ri).mRay->HasPosDir(axis);
777                splitCandidates->push_back(SortableEntry(positive ? SortableEntry::ERayMin : SortableEntry::ERayMax,
778                                                                                                 (*ri).ExtrapOrigin(axis), (void *)&*ri));
779               
780                splitCandidates->push_back(SortableEntry(positive ? SortableEntry::ERayMax : SortableEntry::ERayMin,
781                                                                                                 (*ri).ExtrapTermination(axis), (void *)&*ri));
782        }
783       
784        stable_sort(splitCandidates->begin(), splitCandidates->end());
785}
786
787
788void VspKdTree::EvaluateLeafStats(const TraversalData &data)
789{
790        // the node became a leaf -> evaluate stats for leafs
791        VspKdTreeLeaf *leaf = (VspKdTreeLeaf *)data.node;
792
793        if (data.depth >= termMaxDepth)
794                ++ mStat.maxDepthNodes;
795 
796        //  if ( (int)(leaf->rays.size()) < termMinCost)
797        //    stat.minCostNodes++;
798        if (leaf->GetPvsSize() < termMinPvs)
799                ++ mStat.minPvsNodes;
800
801        if (leaf->GetPvsSize() < termMinRays)
802                ++ mStat.minRaysNodes;
803
804        if (0 && leaf->GetAvgRayContribution() > termMaxRayContribution)
805                ++ mStat.maxRayContribNodes;
806       
807        if (SqrMagnitude(data.bbox.Size()) <= termMinSize)
808                ++ mStat.minSizeNodes;
809
810//      if ((int)(leaf->rays.size()) > stat.maxRayRefs)
811//              mStat.maxRayRefs = (int)leaf->rays.size();
812}
813
814
815
816VspKdTreeNode *VspKdTree::SubdivideNode(VspKdTreeLeaf *leaf,
817                                                                                const AxisAlignedBox3 &box,
818                                                                                AxisAlignedBox3 &backBBox,
819                                                                                AxisAlignedBox3 &frontBBox)
820{
821        if ( (leaf->GetPvsSize() < termMinPvs) || (leaf->rays.size() < termMinRays) ||
822        // (leaf->GetAvgRayContribution() > termMaxRayContribution ) ||
823       (leaf->depth >= termMaxDepth) || SqrMagnitude(box.Size()) <= termMinSize )
824        {
825
826#if 0
827                if (leaf->depth >= termMaxDepth) {
828                        cout<<"Warning: max depth reached depth="<<(int)leaf->depth<<" rays="<<leaf->rays.size()<<endl;
829                        cout<<"Bbox: "<<GetBBox(leaf)<<" dirbbox:"<<GetDirBBox(leaf)<<endl;
830                }
831#endif
832               
833                return leaf;
834        }
835       
836        float position;
837        // first count ray sides
838        int raysBack;
839        int raysFront;
840        int pvsBack;
841        int pvsFront;
842       
843        // select subdivision axis
844        int axis = SelectPlane( leaf, box, position, raysBack, raysFront, pvsBack, pvsFront);
845
846        //      cout<<"rays back="<<raysBack<<" rays front="<<raysFront<<" pvs back="<<pvsBack<<" pvs front="<<
847        //              pvsFront<<endl;
848
849        if (axis == -1)
850        {
851                return leaf;
852        }
853 
854        mStat.nodes += 2;
855        //++ mStat.splits[axis];
856
857        // add the new nodes to the tree
858        VspKdTreeInterior *node = new VspKdTreeInterior(leaf->parent);
859
860        node->axis = axis;
861        node->position = position;
862        node->bbox = box;
863        node->dirBBox = GetDirBBox(leaf);
864 
865        backBBox = box;
866        frontBBox = box;
867
868        VspKdTreeLeaf *back = new VspKdTreeLeaf(node, raysBack);
869        back->SetPvsSize(pvsBack);
870        VspKdTreeLeaf *front = new VspKdTreeLeaf(node, raysFront);
871        front->SetPvsSize(pvsFront);
872
873        // replace a link from node's parent
874        if (leaf->parent)
875                leaf->parent->ReplaceChildLink(leaf, node);
876        // and setup child links
877        node->SetupChildLinks(back, front);
878       
879        if (axis <= VspKdTreeNode::SPLIT_Z)
880        {
881                backBBox.SetMax(axis, position);
882                frontBBox.SetMin(axis, position);
883               
884                for(VspKdTreeNode::RayInfoContainer::iterator ri = leaf->rays.begin();
885                                ri != leaf->rays.end(); ri++)
886                {
887                        if ((*ri).mRay->IsActive())
888                        {
889                                // first unref ray from the former leaf
890                                (*ri).mRay->Unref();
891
892                                // determine the side of this ray with respect to the plane
893                                int side = node->ComputeRayIntersection(*ri, (*ri).mRay->mT);
894
895                                if (side == 0)
896                                {
897                                        if ((*ri).mRay->HasPosDir(axis))
898                                        {
899                                                back->AddRay(VspKdTreeNode::RayInfo((*ri).mRay,
900                                                                                                                        (*ri).mMinT,
901                                                                                                                        (*ri).mRay->mT));
902                                                front->AddRay(VspKdTreeNode::RayInfo((*ri).mRay,
903                                                                                                                         (*ri).mRay->mT,
904                                                                                                                         (*ri).mMaxT));
905                                        }
906                                        else
907                                        {
908                                                back->AddRay(VspKdTreeNode::RayInfo((*ri).mRay,
909                                                                                                                        (*ri).mRay->mT,
910                                                                                                                        (*ri).mMaxT));
911                                                front->AddRay(VspKdTreeNode::RayInfo((*ri).mRay,
912                                                                                                                         (*ri).mMinT,
913                                                                                                                         (*ri).mRay->mT));
914                                        }
915                                }
916                                else
917                                        if (side == 1)
918                                                front->AddRay(*ri);
919                                        else
920                                                back->AddRay(*ri);
921                        } else
922                                (*ri).mRay->Unref();
923                }
924        }
925        else
926        {
927                // rays front/back
928   
929        for(VspKdTreeNode::RayInfoContainer::iterator ri = leaf->rays.begin();
930                        ri != leaf->rays.end(); ++ ri)
931                {
932                        if ((*ri).mRay->IsActive())
933                        {
934                                // first unref ray from the former leaf
935                                (*ri).mRay->Unref();
936
937                                int side;
938                                if ((*ri).mRay->GetDirParametrization(axis - 3) > position)
939                                        side = 1;
940                                else
941                                        side = -1;
942                               
943                                if (side == 1)
944                                        front->AddRay(*ri);
945                                else
946                                        back->AddRay(*ri);             
947                        }
948                        else
949                                (*ri).mRay->Unref();
950                }
951        }
952       
953        // update stats
954        mStat.rayRefs -= (int)leaf->rays.size();
955        mStat.rayRefs += raysBack + raysFront;
956
957        delete leaf;
958        return node;
959}
960
961int VspKdTree::ReleaseMemory(const int time)
962{
963        stack<VspKdTreeNode *> tstack;
964
965        // find a node in the tree which subtree will be collapsed
966        int maxAccessTime = time - accessTimeThreshold;
967        int released;
968
969        tstack.push(root);
970
971        while (!tstack.empty())
972        {
973                VspKdTreeNode *node = tstack.top();
974                tstack.pop();
975   
976                if (!node->IsLeaf())
977                {
978                        VspKdTreeInterior *in = (VspKdTreeInterior *)node;
979                        // cout<<"depth="<<(int)in->depth<<" time="<<in->lastAccessTime<<endl;
980                        if (in->depth >= minCollapseDepth && in->lastAccessTime <= maxAccessTime)
981                        {
982                                released = CollapseSubtree(node, time);
983                                break;
984                        }
985                        if (in->back->GetAccessTime() <  in->front->GetAccessTime())
986                        {
987                                tstack.push(in->front);
988                                tstack.push(in->back);
989                        }
990                        else
991                        {
992                                tstack.push(in->back);
993                                tstack.push(in->front);
994                        }
995                }
996        }
997
998        while (tstack.empty())
999        {
1000                // could find node to collaps...
1001                // cout<<"Could not find a node to release "<<endl;
1002                break;
1003        }
1004 
1005        return released;
1006}
1007
1008
1009
1010
1011VspKdTreeNode *VspKdTree::SubdivideLeaf(VspKdTreeLeaf *leaf,
1012                                                                                const float sizeThreshold)
1013{
1014        VspKdTreeNode *node = leaf;
1015
1016        AxisAlignedBox3 leafBBox = GetBBox(leaf);
1017
1018        static int pass = 0;
1019        ++ pass;
1020       
1021        // check if we should perform a dynamic subdivision of the leaf
1022        if (// leaf->rays.size() > (unsigned)termMinCost &&
1023                (leaf->GetPvsSize() >= termMinPvs) &&
1024                (SqrMagnitude(leafBBox.Size()) > sizeThreshold) )
1025        {
1026        // memory check and realese...
1027                if (GetMemUsage() > maxTotalMemory)
1028                        ReleaseMemory(pass);
1029               
1030                AxisAlignedBox3 backBBox, frontBBox;
1031
1032                // subdivide the node
1033                node = SubdivideNode(leaf, leafBBox, backBBox, frontBBox);
1034        }
1035        return node;
1036}
1037
1038
1039
1040void
1041VspKdTree::UpdateRays(VssRayContainer &remove,
1042                                          VssRayContainer &add)
1043{
1044        VspKdTreeLeaf::NewMail();
1045
1046        // schedule rays for removal
1047        for(VssRayContainer::const_iterator ri = remove.begin();
1048                ri != remove.end(); ++ ri)
1049        {
1050                (*ri)->ScheduleForRemoval(); 
1051        }
1052
1053        int inactive = 0;
1054
1055        for (VssRayContainer::const_iterator ri = remove.begin(); ri != remove.end(); ++ ri)
1056        {
1057                if ((*ri)->ScheduledForRemoval())
1058                        //      RemoveRay(*ri, NULL, false);
1059                        // !!! BUG - with true it does not work correctly - aggreated delete
1060                        RemoveRay(*ri, NULL, true);
1061                else
1062                        ++ inactive;
1063        }
1064
1065        //  cout<<"all/inactive"<<remove.size()<<"/"<<inactive<<endl;
1066 
1067        for(VssRayContainer::const_iterator ri = add.begin(); ri != add.end(); ++ ri)
1068        {
1069                AddRay(*ri);
1070        }
1071}
1072
1073
1074void VspKdTree::RemoveRay(VssRay *ray,
1075                                                  vector<VspKdTreeLeaf *> *affectedLeaves,
1076                                                  const bool removeAllScheduledRays)
1077{
1078        stack<RayTraversalData> tstack;
1079
1080        tstack.push(RayTraversalData(root, VspKdTreeNode::RayInfo(ray)));
1081 
1082        RayTraversalData data;
1083
1084        // cout<<"Number of ray refs = "<<ray->RefCount()<<endl;
1085        while (!tstack.empty())
1086        {
1087                data = tstack.top();
1088                tstack.pop();
1089
1090                if (!data.node->IsLeaf())
1091                {
1092                        // split the set of rays in two groups intersecting the
1093                        // two subtrees
1094
1095                        TraverseInternalNode(data, tstack);
1096        }
1097                else
1098                {
1099                        // remove the ray from the leaf
1100                        // find the ray in the leaf and swap it with the last ray...
1101                        VspKdTreeLeaf *leaf = (VspKdTreeLeaf *)data.node;
1102     
1103                        if (!leaf->Mailed())
1104                        {
1105                                leaf->Mail();
1106                               
1107                                if (affectedLeaves)
1108                                        affectedLeaves->push_back(leaf);
1109       
1110                                if (removeAllScheduledRays)
1111                                {
1112                                        int tail = (int)leaf->rays.size() - 1;
1113
1114                                        for (int i=0; i < (int)(leaf->rays.size()); ++ i)
1115                                        {
1116                                                if (leaf->rays[i].mRay->ScheduledForRemoval())
1117                                                {
1118                                                        // find a ray to replace it with
1119                                                        while (tail >= i && leaf->rays[tail].mRay->ScheduledForRemoval())
1120                                                        {
1121                                                                ++ mStat.removedRayRefs;
1122                                                                leaf->rays[tail].mRay->Unref();
1123                                                                leaf->rays.pop_back();
1124                                                               
1125                                                                -- tail;
1126                                                        }
1127                                                       
1128                                                        if (tail < i)
1129                                                                break;
1130             
1131                                                        ++ mStat.removedRayRefs;
1132             
1133                                                        leaf->rays[i].mRay->Unref();
1134                                                        leaf->rays[i] = leaf->rays[tail];
1135                                                        leaf->rays.pop_back();
1136                                                        -- tail;
1137                                                }
1138                                        }
1139                                }
1140                        }
1141     
1142                        if (!removeAllScheduledRays)
1143                                for (int i=0; i < (int)leaf->rays.size(); i++)
1144                                {
1145                                        if (leaf->rays[i].mRay == ray)
1146                                        {
1147                                                ++ mStat.removedRayRefs;
1148                                                ray->Unref();
1149                                                leaf->rays[i] = leaf->rays[leaf->rays.size() - 1];
1150                                                leaf->rays.pop_back();
1151                                            // check this ray again
1152                                                break;
1153                                        }
1154                                }
1155                }
1156        }
1157
1158        if (ray->RefCount() != 0)
1159        {
1160                cerr << "Error: Number of remaining refs = " << ray->RefCount() << endl;
1161                exit(1);
1162        }
1163
1164}
1165
1166void VspKdTree::AddRay(VssRay *ray)
1167{
1168        stack<RayTraversalData> tstack;
1169 
1170        tstack.push(RayTraversalData(root, VspKdTreeNode::RayInfo(ray)));
1171 
1172        RayTraversalData data;
1173 
1174        while (!tstack.empty())
1175        {
1176                data = tstack.top();
1177                tstack.pop();
1178
1179                if (!data.node->IsLeaf())
1180                {
1181                        TraverseInternalNode(data, tstack);
1182                }
1183                else
1184                {
1185                        // remove the ray from the leaf
1186                        // find the ray in the leaf and swap it with the last ray...
1187                        VspKdTreeLeaf *leaf = (VspKdTreeLeaf *)data.node;
1188                        leaf->AddRay(data.rayData);
1189                        ++ mStat.addedRayRefs;
1190                }
1191        }
1192}
1193
1194void VspKdTree::TraverseInternalNode(RayTraversalData &data,
1195                                                                         stack<RayTraversalData> &tstack)
1196{
1197        VspKdTreeInterior *in =  (VspKdTreeInterior *) data.node;
1198 
1199        if (in->axis <= VspKdTreeNode::SPLIT_Z)
1200        {
1201            // determine the side of this ray with respect to the plane
1202                int side = in->ComputeRayIntersection(data.rayData, data.rayData.mRay->mT);
1203   
1204      if (side == 0)
1205          {
1206                if (data.rayData.mRay->HasPosDir(in->axis))
1207                {
1208                        tstack.push(RayTraversalData(in->back,
1209                                                VspKdTreeNode::RayInfo(data.rayData.mRay, data.rayData.mMinT,
1210                                                                                           data.rayData.mRay->mT)));
1211                               
1212                        tstack.push(RayTraversalData(in->front,
1213                                                VspKdTreeNode::RayInfo(data.rayData.mRay, data.rayData.mRay->mT,
1214                                                                                           data.rayData.mMaxT)));
1215       
1216                }
1217                else
1218                {
1219                        tstack.push(RayTraversalData(in->back,
1220                                                                                 VspKdTreeNode::RayInfo(data.rayData.mRay,
1221                                                                                                                                data.rayData.mRay->mT,
1222                                                                                                                                data.rayData.mMaxT)));
1223                        tstack.push(RayTraversalData(in->front,
1224                                                                                 VspKdTreeNode::RayInfo(data.rayData.mRay,
1225                                                                                                                                data.rayData.mMinT,
1226                                                                                                                                data.rayData.mRay->mT)));
1227                }
1228          }
1229          else
1230          if (side == 1)
1231                          tstack.push(RayTraversalData(in->front, data.rayData));
1232                  else
1233                          tstack.push(RayTraversalData(in->back, data.rayData));
1234        }
1235        else
1236        {
1237                // directional split
1238                if (data.rayData.mRay->GetDirParametrization(in->axis - 3) > in->position)
1239                        tstack.push(RayTraversalData(in->front, data.rayData));
1240                else
1241                        tstack.push(RayTraversalData(in->back, data.rayData));
1242        }
1243}
1244
1245
1246int VspKdTree::CollapseSubtree(VspKdTreeNode *sroot, const int time)
1247{
1248        // first count all rays in the subtree
1249        // use mail 1 for this purpose
1250        stack<VspKdTreeNode *> tstack;
1251       
1252        int rayCount = 0;
1253        int totalRayCount = 0;
1254        int collapsedNodes = 0;
1255
1256#if DEBUG_COLLAPSE
1257        cout<<"Collapsing subtree"<<endl;
1258        cout<<"acessTime="<<sroot->GetAccessTime()<<endl;
1259        cout<<"depth="<<(int)sroot->depth<<endl;
1260#endif
1261
1262        // tstat.collapsedSubtrees++;
1263        // tstat.collapseDepths += (int)sroot->depth;
1264        // tstat.collapseAccessTimes += time - sroot->GetAccessTime();
1265        tstack.push(sroot);
1266        VssRay::NewMail();
1267       
1268        while (!tstack.empty())
1269        {
1270                ++ collapsedNodes;
1271               
1272                VspKdTreeNode *node = tstack.top();
1273                tstack.pop();
1274                if (node->IsLeaf())
1275                {
1276                        VspKdTreeLeaf *leaf = (VspKdTreeLeaf *) node;
1277                       
1278                        for(VspKdTreeNode::RayInfoContainer::iterator ri = leaf->rays.begin();
1279                                        ri != leaf->rays.end(); ++ ri)
1280                        {
1281                                ++ totalRayCount;
1282                               
1283                                if ((*ri).mRay->IsActive() && !(*ri).mRay->Mailed())
1284                                {
1285                                        (*ri).mRay->Mail();
1286                                        ++ rayCount;
1287                                }
1288                        }
1289                }
1290                else
1291                {
1292                        tstack.push(((VspKdTreeInterior *)node)->back);
1293                        tstack.push(((VspKdTreeInterior *)node)->front);
1294                }
1295        }
1296 
1297        VssRay::NewMail();
1298
1299        // create a new node that will hold the rays
1300        VspKdTreeLeaf *newLeaf = new VspKdTreeLeaf(sroot->parent, rayCount);
1301       
1302        if (newLeaf->parent)
1303                newLeaf->parent->ReplaceChildLink(sroot, newLeaf);
1304
1305        tstack.push(sroot);
1306
1307        while (!tstack.empty())
1308        {
1309                VspKdTreeNode *node = tstack.top();
1310                tstack.pop();
1311
1312                if (node->IsLeaf())
1313                {
1314                        VspKdTreeLeaf *leaf = (VspKdTreeLeaf *) node;
1315
1316                        for(VspKdTreeNode::RayInfoContainer::iterator ri = leaf->rays.begin();
1317                                        ri != leaf->rays.end(); ++ ri)
1318                        {
1319                                // unref this ray from the old node             
1320                                if ((*ri).mRay->IsActive())
1321                                {
1322                                        (*ri).mRay->Unref();
1323                                        if (!(*ri).mRay->Mailed())
1324                                        {
1325                                                (*ri).mRay->Mail();
1326                                                newLeaf->AddRay(*ri);
1327                                        }
1328                                }
1329                                else
1330                                        (*ri).mRay->Unref();
1331                        }
1332                }
1333                else
1334                {
1335                        tstack.push(((VspKdTreeInterior *)node)->back);
1336                        tstack.push(((VspKdTreeInterior *)node)->front);
1337                }
1338        }
1339 
1340        // delete the node and all its children
1341        delete sroot;
1342 
1343        // for(VspKdTreeNode::SRayContainer::iterator ri = newLeaf->rays.begin();
1344    //      ri != newLeaf->rays.end(); ++ ri)
1345        //     (*ri).ray->UnMail(2);
1346
1347#if DEBUG_COLLAPSE
1348        cout<<"Total memory before="<<GetMemUsage()<<endl;
1349#endif
1350
1351        mStat.nodes -= collapsedNodes - 1;
1352        mStat.rayRefs -= totalRayCount - rayCount;
1353 
1354#if DEBUG_COLLAPSE
1355        cout << "collapsed nodes" << collapsedNodes << endl;
1356        cout << "collapsed rays" << totalRayCount - rayCount << endl;
1357        cout << "Total memory after=" << GetMemUsage() << endl;
1358        cout << "================================" << endl;
1359#endif
1360
1361        //  tstat.collapsedNodes += collapsedNodes;
1362        //  tstat.collapsedRays += totalRayCount - rayCount;
1363    return totalRayCount - rayCount;
1364}
1365
1366
1367int VspKdTree::GetPvsSize(VspKdTreeNode *node, const AxisAlignedBox3 &box) const
1368{
1369        stack<VspKdTreeNode *> tstack;
1370        tstack.push(root);
1371
1372        Intersectable::NewMail();
1373        int pvsSize = 0;
1374       
1375        while (!tstack.empty())
1376        {
1377                VspKdTreeNode *node = tstack.top();
1378                tstack.pop();
1379   
1380          if (node->IsLeaf())
1381                {
1382                        VspKdTreeLeaf *leaf = (VspKdTreeLeaf *)node;
1383                        for (VspKdTreeNode::RayInfoContainer::iterator ri = leaf->rays.begin();
1384                                ri != leaf->rays.end(); ++ ri)
1385                        {
1386                                if ((*ri).mRay->IsActive())
1387                                {
1388                                        Intersectable *object;
1389#if BIDIRECTIONAL_RAY
1390                                        object = (*ri).mRay->mOriginObject;
1391                                       
1392                                        if (object && !object->Mailed())
1393                                        {
1394                                                ++ pvsSize;
1395                                                object->Mail();
1396                                        }
1397#endif
1398                                        object = (*ri).mRay->mTerminationObject;
1399                                        if (object && !object->Mailed())
1400                                        {
1401                                                ++ pvsSize;
1402                                                object->Mail();
1403                                        }
1404                                }
1405                        }
1406                }
1407                else
1408                {
1409                        VspKdTreeInterior *in = (VspKdTreeInterior *)node;
1410       
1411                        if (in->axis < 3)
1412                        {
1413                                if (box.Max(in->axis) >= in->position)
1414                                        tstack.push(in->front);
1415                               
1416                                if (box.Min(in->axis) <= in->position)
1417                                        tstack.push(in->back);
1418                        }
1419                        else
1420                        {
1421                                // both nodes for directional splits
1422                                tstack.push(in->front);
1423                                tstack.push(in->back);
1424                        }
1425                }
1426        }
1427
1428        return pvsSize;
1429}
1430
1431void VspKdTree::GetRayContributionStatistics(float &minRayContribution,
1432                                                                                         float &maxRayContribution,
1433                                                                                         float &avgRayContribution)
1434{
1435        stack<VspKdTreeNode *> tstack;
1436        tstack.push(root);
1437       
1438        minRayContribution = 1.0f;
1439        maxRayContribution = 0.0f;
1440
1441        float sumRayContribution = 0.0f;
1442        int leaves = 0;
1443
1444        while (!tstack.empty())
1445        {
1446                VspKdTreeNode *node = tstack.top();
1447                tstack.pop();
1448                if (node->IsLeaf())
1449                {
1450                        leaves++;
1451                        VspKdTreeLeaf *leaf = (VspKdTreeLeaf *)node;
1452                        float c = leaf->GetAvgRayContribution();
1453
1454                        if (c > maxRayContribution)
1455                                maxRayContribution = c;
1456                        if (c < minRayContribution)
1457                                minRayContribution = c;
1458                        sumRayContribution += c;
1459                       
1460                }
1461                else {
1462                        VspKdTreeInterior *in = (VspKdTreeInterior *)node;
1463                        // both nodes for directional splits
1464                        tstack.push(in->front);
1465                        tstack.push(in->back);
1466                }
1467        }
1468       
1469        cout<<"sum="<<sumRayContribution<<endl;
1470        cout<<"leaves="<<leaves<<endl;
1471        avgRayContribution = sumRayContribution/(float)leaves;
1472}
1473
1474
1475int VspKdTree::GenerateRays(const float ratioPerLeaf,
1476                                                        SimpleRayContainer &rays)
1477{
1478        stack<VspKdTreeNode *> tstack;
1479        tstack.push(root);
1480       
1481        while (!tstack.empty())
1482        {
1483                VspKdTreeNode *node = tstack.top();
1484                tstack.pop();
1485               
1486                if (node->IsLeaf())
1487                {
1488                        VspKdTreeLeaf *leaf = (VspKdTreeLeaf *)node;
1489                        float c = leaf->GetAvgRayContribution();
1490                        int num = (int)(c*ratioPerLeaf + 0.5);
1491                        //                      cout<<num<<" ";
1492
1493                        for (int i=0; i < num; i++)
1494                        {
1495                                Vector3 origin = GetBBox(leaf).GetRandomPoint();
1496                                Vector3 dirVector = GetDirBBox(leaf).GetRandomPoint();
1497                                Vector3 direction = Vector3(sin(dirVector.x), sin(dirVector.y), cos(dirVector.x));
1498                                //cout<<"dir vector.x="<<dirVector.x<<"direction'.x="<<atan2(direction.x, direction.y)<<endl;
1499                                rays.push_back(SimpleRay(origin, direction));
1500                        }
1501                       
1502                }
1503                else
1504                {
1505                        VspKdTreeInterior *in = (VspKdTreeInterior *)node;
1506                        // both nodes for directional splits
1507                        tstack.push(in->front);
1508                        tstack.push(in->back);
1509                }
1510        }
1511
1512        return (int)rays.size();
1513}
1514
1515
1516float VspKdTree::GetAvgPvsSize()
1517{
1518        stack<VspKdTreeNode *> tstack;
1519        tstack.push(root);
1520
1521        int sumPvs = 0;
1522        int leaves = 0;
1523       
1524        while (!tstack.empty())
1525        {
1526                VspKdTreeNode *node = tstack.top();
1527                tstack.pop();
1528               
1529                if (node->IsLeaf())
1530                {
1531                        VspKdTreeLeaf *leaf = (VspKdTreeLeaf *)node;
1532                        // update pvs size
1533                        leaf->UpdatePvsSize();
1534                        sumPvs += leaf->GetPvsSize();
1535                        leaves++;
1536                }
1537                else
1538                {
1539                        VspKdTreeInterior *in = (VspKdTreeInterior *)node;
1540                        // both nodes for directional splits
1541                        tstack.push(in->front);
1542                        tstack.push(in->back);
1543                }
1544        }
1545
1546        return sumPvs / (float)leaves;
1547}
1548
1549VspKdTreeNode *VspKdTree::GetRoot() const
1550{
1551        return mRoot;
1552}
1553
1554AxisAlignedBox3 VspKdTree::GetBBox(const VspKdTreeNode *node)
1555{
1556        if (node->parent == NULL)
1557                return mBbox;
1558
1559        if (!node->IsLeaf())
1560                return ((VspKdTreeInterior *)node)->mBbox;
1561
1562        if (node->parent->axis >= 3)
1563                return node->parent->mBbox;
1564   
1565        AxisAlignedBox3 box(node->parent->mBbox);
1566        if (node->parent->front == node)
1567                box.SetMin(node->parent->axis, node->parent->position);
1568        else
1569                box.SetMax(node->parent->axis, node->parent->position);
1570        return mBbox;
1571}
1572
1573int     VspKdTree::GetRootPvsSize() const
1574{
1575        return GetPvsSize(mRoot, mBox);
1576}
Note: See TracBrowser for help on using the repository browser.