source: GTP/trunk/Lib/Vis/Preprocessing/src/BvHierarchy.cpp @ 1718

Revision 1718, 64.1 KB checked in by mattausch, 18 years ago (diff)
Line 
1#include <stack>
2#include <time.h>
3#include <iomanip>
4
5#include "BvHierarchy.h"
6#include "ViewCell.h"
7#include "Plane3.h"
8#include "Mesh.h"
9#include "common.h"
10#include "Environment.h"
11#include "Polygon3.h"
12#include "Ray.h"
13#include "AxisAlignedBox3.h"
14#include "Exporter.h"
15#include "Plane3.h"
16#include "ViewCellsManager.h"
17#include "Beam.h"
18#include "VspTree.h"
19#include "HierarchyManager.h"
20
21
22namespace GtpVisibilityPreprocessor {
23
24
25#define PROBABILIY_IS_BV_VOLUME 1
26#define USE_FIXEDPOINT_T 0
27#define USE_VOLUMES_FOR_HEURISTICS 1
28
29int BvhNode::sMailId = 10000; //2147483647;
30int BvhNode::sReservedMailboxes = 1;
31
32BvHierarchy *BvHierarchy::BvhSubdivisionCandidate::sBvHierarchy = NULL;
33
34
35/// sorting operator
36inline static bool ilt(Intersectable *obj1, Intersectable *obj2)
37{
38        return obj1->mId < obj2->mId;
39}
40
41
42/***************************************************************/
43/*              class BvhNode implementation                   */
44/***************************************************************/
45
46BvhNode::BvhNode():
47mParent(NULL),
48mMailbox(0),
49mRenderCostDecr(0),
50mMemoryIncr(0),
51mPvsEntriesIncr(0),
52mTimeStamp(0)
53{
54       
55}
56
57BvhNode::BvhNode(const AxisAlignedBox3 &bbox):
58mParent(NULL),
59mBoundingBox(bbox),
60mMailbox(0),
61mMemoryIncr(0),
62mRenderCostDecr(0),
63mPvsEntriesIncr(0),
64mTimeStamp(0)
65{
66}
67
68
69BvhNode::BvhNode(const AxisAlignedBox3 &bbox, BvhInterior *parent):
70mBoundingBox(bbox),
71mParent(parent),
72mMailbox(0),
73mMemoryIncr(0),
74mRenderCostDecr(0),
75mPvsEntriesIncr(0),
76mTimeStamp(0)
77{
78}
79
80
81bool BvhNode::IsRoot() const
82{
83        return mParent == NULL;
84}
85
86
87BvhInterior *BvhNode::GetParent()
88{
89        return mParent;
90}
91
92
93void BvhNode::SetParent(BvhInterior *parent)
94{
95        mParent = parent;
96}
97
98
99
100/******************************************************************/
101/*              class BvhInterior implementation                  */
102/******************************************************************/
103
104
105BvhLeaf::BvhLeaf(const AxisAlignedBox3 &bbox):
106BvhNode(bbox),
107mSubdivisionCandidate(NULL)
108{
109        mActiveNode = this;
110}
111
112
113BvhLeaf::BvhLeaf(const AxisAlignedBox3 &bbox, BvhInterior *parent):
114BvhNode(bbox, parent)
115{
116        mActiveNode = this;
117}
118
119
120BvhLeaf::BvhLeaf(const AxisAlignedBox3 &bbox,
121                                 BvhInterior *parent,
122                                 const int numObjects):
123BvhNode(bbox, parent)
124{
125        mObjects.reserve(numObjects);
126        mActiveNode = this;
127}
128
129
130bool BvhLeaf::IsLeaf() const
131{
132        return true;
133}
134
135
136BvhLeaf::~BvhLeaf()
137{
138}
139
140
141void BvhLeaf::CollectObjects(ObjectContainer &objects)
142{
143        ObjectContainer::const_iterator oit, oit_end = mObjects.end();
144        for (oit = mObjects.begin(); oit != oit_end; ++ oit)
145        {
146                objects.push_back(*oit);
147        }
148}
149
150/******************************************************************/
151/*              class BvhInterior implementation                  */
152/******************************************************************/
153
154
155BvhInterior::BvhInterior(const AxisAlignedBox3 &bbox):
156BvhNode(bbox), mFront(NULL), mBack(NULL)
157{
158}
159
160
161BvhInterior::BvhInterior(const AxisAlignedBox3 &bbox, BvhInterior *parent):
162BvhNode(bbox, parent), mFront(NULL), mBack(NULL)
163{
164}
165
166
167void BvhInterior::ReplaceChildLink(BvhNode *oldChild, BvhNode *newChild)
168{
169        if (mBack == oldChild)
170                mBack = newChild;
171        else
172                mFront = newChild;
173}
174
175
176bool BvhInterior::IsLeaf() const
177{
178        return false;
179}
180
181
182BvhInterior::~BvhInterior()
183{
184        DEL_PTR(mFront);
185        DEL_PTR(mBack);
186}
187
188
189void BvhInterior::SetupChildLinks(BvhNode *front, BvhNode *back)
190{
191    mBack = back;
192    mFront = front;
193}
194
195
196void BvhInterior::CollectObjects(ObjectContainer &objects)
197{
198        mFront->CollectObjects(objects);
199        mBack->CollectObjects(objects);
200}
201
202
203/*******************************************************************/
204/*                  class BvHierarchy implementation               */
205/*******************************************************************/
206
207
208BvHierarchy::BvHierarchy():
209mRoot(NULL),
210mTimeStamp(1)
211{
212        ReadEnvironment();
213        mSubdivisionCandidates = new SortableEntryContainer;
214        for (int i = 0; i < 3; ++ i)
215                mSortedObjects[i] = NULL;
216}
217
218
219BvHierarchy::~BvHierarchy()
220{
221        // delete bvh intersectables
222        BvhIntersectableMap::iterator it, it_end = mBvhIntersectables.end();
223
224        for (it = mBvhIntersectables.begin(); it != mBvhIntersectables.end(); ++ it)
225        {
226                DEL_PTR((*it).second);
227        }
228
229        // delete the local subdivision candidates
230        DEL_PTR(mSubdivisionCandidates);
231
232        // delete the presorted objects
233        for (int i = 0; i < 3; ++ i)
234        {
235                DEL_PTR(mSortedObjects[i]);
236        }
237       
238        // delete the tree
239        DEL_PTR(mRoot);
240}
241
242
243void BvHierarchy::ReadEnvironment()
244{
245        bool randomize = false;
246        Environment::GetSingleton()->GetBoolValue("VspTree.Construction.randomize", randomize);
247         
248        // initialise random generator for heuristics
249        if (randomize)
250                Randomize();
251
252        //////////////////////////////
253        //-- termination criteria for autopartition
254
255        Environment::GetSingleton()->GetIntValue("BvHierarchy.Termination.maxDepth", mTermMaxDepth);
256        Environment::GetSingleton()->GetIntValue("BvHierarchy.Termination.maxLeaves", mTermMaxLeaves);
257        Environment::GetSingleton()->GetIntValue("BvHierarchy.Termination.minObjects", mTermMinObjects);
258        Environment::GetSingleton()->GetIntValue("BvHierarchy.Termination.minRays", mTermMinRays);
259        Environment::GetSingleton()->GetFloatValue(
260                "BvHierarchy.Termination.minProbability", mTermMinProbability);
261    Environment::GetSingleton()->GetIntValue("BvHierarchy.Termination.missTolerance", mTermMissTolerance);
262
263
264        //////////////////////////////
265        //-- max cost ratio for early tree termination
266
267        Environment::GetSingleton()->GetFloatValue("BvHierarchy.Termination.maxCostRatio", mTermMaxCostRatio);
268        Environment::GetSingleton()->GetFloatValue("BvHierarchy.Termination.minGlobalCostRatio",
269                mTermMinGlobalCostRatio);
270        Environment::GetSingleton()->GetIntValue("BvHierarchy.Termination.globalCostMissTolerance",
271                mTermGlobalCostMissTolerance);
272
273
274        //////////////////////////////
275        //-- factors for subdivision heuristics
276
277        // if only the driving axis is used for splits
278        Environment::GetSingleton()->GetBoolValue("BvHierarchy.splitUseOnlyDrivingAxis", mOnlyDrivingAxis);
279        Environment::GetSingleton()->GetFloatValue("BvHierarchy.maxStaticMemory", mMaxMemory);
280        Environment::GetSingleton()->GetBoolValue("BvHierarchy.useCostHeuristics", mUseCostHeuristics);
281        Environment::GetSingleton()->GetBoolValue("BvHierarchy.useSah", mUseSah);
282
283    char subdivisionStatsLog[100];
284        Environment::GetSingleton()->GetStringValue("BvHierarchy.subdivisionStats", subdivisionStatsLog);
285        mSubdivisionStats.open(subdivisionStatsLog);
286
287        Environment::GetSingleton()->GetFloatValue(
288                "BvHierarchy.Construction.renderCostDecreaseWeight", mRenderCostDecreaseWeight);
289       
290        Environment::GetSingleton()->GetBoolValue(
291                "BvHierarchy.Construction.useGlobalSorting", mUseGlobalSorting);
292
293        Environment::GetSingleton()->GetIntValue("BvHierarchy.minRaysForVisibility", mMinRaysForVisibility);
294       
295        //mUseBboxAreaForSah = false;
296        mUseBboxAreaForSah = true;
297
298        /////////////
299        //-- debug output
300
301        Debug << "******* Bvh hierarchy options ******** " << endl;
302    Debug << "max depth: " << mTermMaxDepth << endl;
303        Debug << "min probabiliy: " << mTermMinProbability<< endl;
304        Debug << "min objects: " << mTermMinObjects << endl;
305        Debug << "max cost ratio: " << mTermMaxCostRatio << endl;
306        Debug << "miss tolerance: " << mTermMissTolerance << endl;
307        Debug << "max leaves: " << mTermMaxLeaves << endl;
308        Debug << "randomize: " << randomize << endl;
309        Debug << "min global cost ratio: " << mTermMinGlobalCostRatio << endl;
310        Debug << "global cost miss tolerance: " << mTermGlobalCostMissTolerance << endl;
311        Debug << "only driving axis: " << mOnlyDrivingAxis << endl;
312        Debug << "max memory: " << mMaxMemory << endl;
313        Debug << "use cost heuristics: " << mUseCostHeuristics << endl;
314        Debug << "use sah (surface area heuristics: " << mUseSah << endl;
315        Debug << "subdivision stats log: " << subdivisionStatsLog << endl;
316        Debug << "split borders: " << mSplitBorder << endl;
317        Debug << "render cost decrease weight: " << mRenderCostDecreaseWeight << endl;
318        Debug << "use global sort: " << mUseGlobalSorting << endl;
319        Debug << "minimal rays for visibility: " << mMinRaysForVisibility << endl;
320
321        Debug << endl;
322}
323
324
325void BvHierarchy::AssociateObjectsWithLeaf(BvhLeaf *leaf)
326{
327        ObjectContainer::const_iterator oit, oit_end = leaf->mObjects.end();
328
329        for (oit = leaf->mObjects.begin(); oit != oit_end; ++ oit)
330        {
331                (*oit)->mBvhLeaf = leaf;
332        }
333}
334
335
336static int CountRays(const ObjectContainer &objects)
337{
338        int nRays = 0;
339
340        ObjectContainer::const_iterator oit, oit_end = objects.end();
341
342        for (oit = objects.begin(); oit != oit_end; ++ oit)
343        {
344                nRays += (int)(*oit)->GetOrCreateRays()->size();
345        }
346
347        return nRays;
348}
349                                                                       
350
351BvhInterior *BvHierarchy::SubdivideNode(const BvhSubdivisionCandidate &sc,
352                                                                                BvhTraversalData &frontData,
353                                                                                BvhTraversalData &backData)
354{
355        const BvhTraversalData &tData = sc.mParentData;
356        BvhLeaf *leaf = tData.mNode;
357        AxisAlignedBox3 parentBox = leaf->GetBoundingBox();
358
359        // update stats: we have two new leaves
360        mBvhStats.nodes += 2;
361
362        if (tData.mDepth > mBvhStats.maxDepth)
363        {
364                mBvhStats.maxDepth = tData.mDepth;
365        }
366
367        // add the new nodes to the tree
368        BvhInterior *node = new BvhInterior(parentBox, leaf->GetParent());
369       
370
371        //////////////////
372        //-- create front and back leaf
373
374        AxisAlignedBox3 fbox = EvalBoundingBox(sc.mFrontObjects, &parentBox);
375        AxisAlignedBox3 bbox = EvalBoundingBox(sc.mBackObjects, &parentBox);
376
377        BvhLeaf *back = new BvhLeaf(bbox, node, (int)sc.mBackObjects.size());
378        BvhLeaf *front = new BvhLeaf(fbox, node, (int)sc.mFrontObjects.size());
379
380        BvhInterior *parent = leaf->GetParent();
381
382        // replace a link from node's parent
383        if (parent)
384        {
385                parent->ReplaceChildLink(leaf, node);
386                node->SetParent(parent);
387        }
388        else // no parent => this node is the root
389        {
390                mRoot = node;
391        }
392
393        // and setup child links
394        node->SetupChildLinks(front, back);
395
396        ++ mBvhStats.splits;
397
398
399        ////////////////////////////////////////
400        //-- fill  front and back traversal data with the new values
401
402        frontData.mDepth = backData.mDepth = tData.mDepth + 1;
403
404        frontData.mNode = front;
405        backData.mNode = back;
406
407        back->mObjects = sc.mBackObjects;
408        front->mObjects = sc.mFrontObjects;
409
410        // if the number of rays is too low, no assumptions can be made
411        // (=> switch to surface area heuristics?)
412        frontData.mNumRays = CountRays(sc.mFrontObjects);
413        backData.mNumRays = CountRays(sc.mBackObjects);
414
415        AssociateObjectsWithLeaf(back);
416        AssociateObjectsWithLeaf(front);
417   
418#if PROBABILIY_IS_BV_VOLUME
419        // volume of bvh (= probability that this bvh can be seen)
420        frontData.mProbability = fbox.GetVolume();
421        backData.mProbability = bbox.GetVolume();
422#else
423        // compute probability of this node being visible,
424        // i.e., volume of the view cells that can see this node
425        frontData.mProbability = EvalViewCellsVolume(sc.mFrontObjects);
426        backData.mProbability = EvalViewCellsVolume(sc.mBackObjects);
427#endif
428
429    // how often was max cost ratio missed in this branch?
430        frontData.mMaxCostMisses = sc.GetMaxCostMisses();
431        backData.mMaxCostMisses = sc.GetMaxCostMisses();
432       
433        // set the time stamp so the order of traversal can be reconstructed
434        node->mTimeStamp = mHierarchyManager->mTimeStamp ++;
435               
436        // assign the objects in sorted order
437        if (mUseGlobalSorting)
438        {
439                AssignSortedObjects(sc, frontData, backData);
440        }
441       
442        // return the new interior node
443        return node;
444}
445
446
447BvhNode *BvHierarchy::Subdivide(SplitQueue &tQueue,
448                                                                SubdivisionCandidate *splitCandidate,
449                                                                const bool globalCriteriaMet)
450{
451        BvhSubdivisionCandidate *sc = dynamic_cast<BvhSubdivisionCandidate *>(splitCandidate);
452        BvhTraversalData &tData = sc->mParentData;
453
454        BvhNode *currentNode = tData.mNode;
455
456        if (!LocalTerminationCriteriaMet(tData) && !globalCriteriaMet)
457        {       
458                //////////////
459                //-- continue subdivision
460
461                BvhTraversalData tFrontData;
462                BvhTraversalData tBackData;
463                       
464                // create new interior node and two leaf node
465                currentNode = SubdivideNode(*sc, tFrontData, tBackData);
466       
467                // decrease the weighted average cost of the subdivisoin
468                mTotalCost -= sc->GetRenderCostDecrease();
469                mPvsEntries += sc->GetPvsEntriesIncr();
470
471                // subdivision statistics
472                if (1) PrintSubdivisionStats(*sc);
473
474
475                ///////////////////////////
476                //-- push the new split candidates on the queue
477               
478                BvhSubdivisionCandidate *frontCandidate = new BvhSubdivisionCandidate(tFrontData);
479                BvhSubdivisionCandidate *backCandidate = new BvhSubdivisionCandidate(tBackData);
480
481                EvalSubdivisionCandidate(*frontCandidate);
482                EvalSubdivisionCandidate(*backCandidate);
483       
484                // cross reference
485                tFrontData.mNode->SetSubdivisionCandidate(frontCandidate);
486                tBackData.mNode->SetSubdivisionCandidate(backCandidate);
487
488                //cout << "f: " << frontCandidate->GetPriority() << " b: " << backCandidate->GetPriority() << endl;
489                tQueue.Push(frontCandidate);
490                tQueue.Push(backCandidate);
491        }
492
493        /////////////////////////////////
494        //-- node is a leaf => terminate traversal
495
496        if (currentNode->IsLeaf())
497        {
498                /////////////////////
499                //-- store additional info
500                EvaluateLeafStats(tData);
501       
502                const bool mStoreRays = true;
503                if (mStoreRays)
504                {
505                        BvhLeaf *leaf = dynamic_cast<BvhLeaf *>(currentNode);
506                        CollectRays(leaf->mObjects, leaf->mVssRays);
507                }
508               
509                //////////////////////////////////////
510               
511                // this leaf is no candidate for splitting anymore
512                // => detach subdivision candidate
513                tData.mNode->SetSubdivisionCandidate(NULL);
514                // detach node so we don't delete it with the traversal data
515                tData.mNode = NULL;
516        }
517       
518        return currentNode;
519}
520
521
522void BvHierarchy::EvalSubdivisionCandidate(BvhSubdivisionCandidate &splitCandidate,
523                                                                                   bool computeSplitPlane)
524{
525        if (computeSplitPlane)
526        {
527                const bool sufficientSamples =
528                        splitCandidate.mParentData.mNumRays > mMinRaysForVisibility;
529
530                const bool useVisibiliyBasedHeuristics =
531                        !mUseSah &&
532                        (mHierarchyManager->GetViewSpaceSubdivisionType() ==
533                        HierarchyManager::KD_BASED_VIEWSPACE_SUBDIV) &&
534                        sufficientSamples;
535
536                // compute best object partition
537                const float ratio =     SelectObjectPartition(splitCandidate.mParentData,
538                                                                                                  splitCandidate.mFrontObjects,
539                                                                                                  splitCandidate.mBackObjects,
540                                                                                                  useVisibiliyBasedHeuristics);
541       
542                // cost ratio violated?
543                const bool maxCostRatioViolated = mTermMaxCostRatio < ratio;
544                const int previousMisses = splitCandidate.mParentData.mMaxCostMisses;
545
546                splitCandidate.SetMaxCostMisses(maxCostRatioViolated ? previousMisses + 1 : previousMisses);
547
548        }
549
550        BvhLeaf *leaf = splitCandidate.mParentData.mNode;
551
552        const float oldProp = EvalViewCellsVolume(leaf->mObjects);
553        const float oldRenderCost = EvalRenderCost(leaf->mObjects);
554               
555        // compute global decrease in render cost
556        const float newRenderCost = EvalRenderCost(splitCandidate.mFrontObjects) +
557                                                                EvalRenderCost(splitCandidate.mBackObjects);
558
559        const float renderCostDecr = oldRenderCost - newRenderCost;
560       
561        splitCandidate.SetRenderCostDecrease(renderCostDecr);
562
563        // increase in pvs entries
564        const int pvsEntriesIncr = EvalPvsEntriesIncr(splitCandidate);
565        splitCandidate.SetPvsEntriesIncr(pvsEntriesIncr);
566
567#ifdef GTP_DEBUG
568        Debug << "old render cost: " << oldRenderCost << endl;
569        Debug << "new render cost: " << newRenderCost << endl;
570        Debug << "render cost decrease: " << renderCostDecr << endl;
571#endif
572
573        float priority;
574
575        // surface area heuristics is used when there is no view space subdivision available.
576        // In order to have some prioritized traversal, use this formula instead
577        if (mHierarchyManager->GetViewSpaceSubdivisionType() ==
578                HierarchyManager::NO_VIEWSPACE_SUBDIV)
579        {
580                priority = EvalSahCost(leaf);
581        }
582        else
583        {
584                // take render cost of node into account
585                // otherwise danger of being stuck in a local minimum!
586                const float factor = mRenderCostDecreaseWeight;
587                priority = factor * renderCostDecr + (1.0f - factor) * oldRenderCost;
588
589                if (mHierarchyManager->mConsiderMemory2)
590                {
591                        priority /= ((float)splitCandidate.GetPvsEntriesIncr() + mHierarchyManager->mMemoryConst);
592                }
593        }
594
595        // compute global decrease in render cost
596        splitCandidate.SetPriority(priority);
597}
598
599
600int BvHierarchy::EvalPvsEntriesIncr(BvhSubdivisionCandidate &splitCandidate) const
601{
602        const int oldPvsSize = CountViewCells(splitCandidate.mParentData.mNode->mObjects);
603       
604        const int fPvsSize = CountViewCells(splitCandidate.mFrontObjects);
605        const int bPvsSize = CountViewCells(splitCandidate.mBackObjects);
606
607        return fPvsSize + bPvsSize - oldPvsSize;
608}
609
610
611inline bool BvHierarchy::LocalTerminationCriteriaMet(const BvhTraversalData &data) const
612{
613        const bool terminationCriteriaMet =
614                        (0
615                        || ((int)data.mNode->mObjects.size() <= 1)//mTermMinObjects)
616                        //|| (data.mProbability <= mTermMinProbability)
617                        //|| (data.mNumRays <= mTermMinRays)
618                 );
619
620#ifdef _DEBUG
621        if (terminationCriteriaMet)
622        {
623                cout << "bvh local termination criteria met:" << endl;
624                cout << "objects: " << data.mNode->mObjects.size() << " " << mTermMinObjects << endl;
625        }
626#endif
627        return terminationCriteriaMet;
628}
629
630
631inline bool BvHierarchy::GlobalTerminationCriteriaMet(const BvhTraversalData &data) const
632{
633        // note: tracking for global cost termination
634        // does not make much sense for interleaved vsp / osp partition
635        // as it is the responsibility of the hierarchy manager
636
637        const bool terminationCriteriaMet =
638                (0
639                || (mBvhStats.Leaves() >= mTermMaxLeaves)
640                //|| (mBvhStats.mGlobalCostMisses >= mTermGlobalCostMissTolerance)
641                //|| mOutOfMemory
642                );
643
644#ifdef GTP_DEBUG
645        if (terminationCriteriaMet)
646        {
647                Debug << "bvh global termination criteria met:" << endl;
648                Debug << "cost misses: " << mBvhStats.mGlobalCostMisses << " " << mTermGlobalCostMissTolerance << endl;
649                Debug << "leaves: " << mBvhStats.Leaves() << " " << mTermMaxLeaves << endl;
650        }
651#endif
652        return terminationCriteriaMet;
653}
654
655
656void BvHierarchy::EvaluateLeafStats(const BvhTraversalData &data)
657{
658        // the node became a leaf -> evaluate stats for leafs
659        BvhLeaf *leaf = data.mNode;
660       
661        ++ mCreatedLeaves;
662
663       
664        if (data.mProbability <= mTermMinProbability)
665        {
666                ++ mBvhStats.minProbabilityNodes;
667        }
668
669        ////////////////////////////////////////////
670        // depth related stuff
671
672        if (data.mDepth < mBvhStats.minDepth)
673        {
674                mBvhStats.minDepth = data.mDepth;
675        }
676
677        if (data.mDepth >= mTermMaxDepth)
678        {
679        ++ mBvhStats.maxDepthNodes;
680        }
681
682        // accumulate depth to compute average depth
683        mBvhStats.accumDepth += data.mDepth;
684
685
686        ////////////////////////////////////////////
687        // objects related stuff
688
689        // note: the sum should alwaysbe total number of objects for bvh
690        mBvhStats.objectRefs += (int)leaf->mObjects.size();
691
692        if ((int)leaf->mObjects.size() <= mTermMinObjects)
693        {
694             ++ mBvhStats.minObjectsNodes;
695        }
696
697        if (leaf->mObjects.empty())
698        {
699                ++ mBvhStats.emptyNodes;
700        }
701
702        if ((int)leaf->mObjects.size() > mBvhStats.maxObjectRefs)
703        {
704                mBvhStats.maxObjectRefs = (int)leaf->mObjects.size();
705        }
706
707        if ((int)leaf->mObjects.size() < mBvhStats.minObjectRefs)
708        {
709                mBvhStats.minObjectRefs = (int)leaf->mObjects.size();
710        }
711
712        ////////////////////////////////////////////
713        // ray related stuff
714
715        // note: this number should always accumulate to the total number of rays
716        mBvhStats.rayRefs += data.mNumRays;
717       
718        if (data.mNumRays <= mTermMinRays)
719        {
720             ++ mBvhStats.minRaysNodes;
721        }
722
723        if (data.mNumRays > mBvhStats.maxRayRefs)
724        {
725                mBvhStats.maxRayRefs = data.mNumRays;
726        }
727
728        if (data.mNumRays < mBvhStats.minRayRefs)
729        {
730                mBvhStats.minRayRefs = data.mNumRays;
731        }
732
733#ifdef _DEBUG
734        cout << "depth: " << data.mDepth << " objects: " << (int)leaf->mObjects.size()
735                 << " rays: " << data.mNumRays << " rays / objects "
736                 << (float)data.mNumRays / (float)leaf->mObjects.size() << endl;
737#endif
738}
739
740
741#if 0
742
743/// compute object boundaries using spatial mid split
744float BvHierarchy::EvalLocalObjectPartition(const BvhTraversalData &tData,
745                                                                                        const int axis,
746                                                                                        ObjectContainer &objectsFront,
747                                                                                        ObjectContainer &objectsBack)
748{
749        const float maxBox = tData.mBoundingBox.Max(axis);
750        const float minBox = tData.mBoundingBox.Min(axis);
751
752        float midPoint = (maxBox + minBox) * 0.5f;
753
754        ObjectContainer::const_iterator oit, oit_end = tData.mNode->mObjects.end();
755       
756        for (oit = tData.mNode->mObjects.begin(); oit != oit_end; ++ oit)
757        {
758                Intersectable *obj = *oit;
759                const AxisAlignedBox3 box = obj->GetBox();
760
761                const float objMid = (box.Max(axis) + box.Min(axis)) * 0.5f;
762
763                // object mailed => belongs to back objects
764                if (objMid < midPoint)
765                {
766                        objectsBack.push_back(obj);
767                }
768                else
769                {
770                        objectsFront.push_back(obj);
771                }
772        }
773
774        const float oldRenderCost = EvalRenderCost(tData.mNode->mObjects);
775        const float newRenderCost = EvalRenderCost(objectsFront) * EvalRenderCost(objectsBack);
776
777        const float ratio = newRenderCost / oldRenderCost;
778        return ratio;
779}
780
781#else
782
783/// compute object partition by getting balanced objects on the left and right side
784float BvHierarchy::EvalLocalObjectPartition(const BvhTraversalData &tData,
785                                                                                        const int axis,
786                                                                                        ObjectContainer &objectsFront,
787                                                                                        ObjectContainer &objectsBack)
788{
789        PrepareLocalSubdivisionCandidates(tData, axis);
790       
791        SortableEntryContainer::const_iterator cit, cit_end = mSubdivisionCandidates->end();
792
793        int i = 0;
794        const int border = (int)tData.mNode->mObjects.size() / 2;
795
796    for (cit = mSubdivisionCandidates->begin(); cit != cit_end; ++ cit, ++ i)
797        {
798                Intersectable *obj = (*cit).mObject;
799
800                // object mailed => belongs to back objects
801                if (i < border)
802                {
803                        objectsBack.push_back(obj);
804                }
805                else
806                {
807                        objectsFront.push_back(obj);
808                }
809        }
810
811#if 1
812        const float cost = (tData.mNode->GetBoundingBox().Size().DrivingAxis() == axis) ? -1.0f : 0.0f;
813#else
814        const float oldRenderCost = EvalRenderCost(tData.mNode->mObjects);
815        const float newRenderCost = EvalRenderCost(objectsFront) * EvalRenderCost(objectsBack);
816
817        const float cost = newRenderCost / oldRenderCost;
818#endif
819
820        return cost;
821}
822#endif
823
824#if 1
825
826float BvHierarchy::EvalSah(const BvhTraversalData &tData,
827                                                   const int axis,
828                                                   ObjectContainer &objectsFront,
829                                                   ObjectContainer &objectsBack)
830{
831        // go through the lists, count the number of objects left and right
832        // and evaluate the following cost funcion:
833        // C = ct_div_ci  + (ol + or) / queries
834        PrepareLocalSubdivisionCandidates(tData, axis);
835
836        const float totalRenderCost = EvalAbsCost(tData.mNode->mObjects);
837        float objectsLeft = 0, objectsRight = totalRenderCost;
838 
839        const AxisAlignedBox3 nodeBbox = tData.mNode->GetBoundingBox();
840        const float boxArea = nodeBbox.SurfaceArea();
841
842        float minSum = 1e20f;
843 
844        float minBorder = nodeBbox.Max(axis);
845        float maxBorder = nodeBbox.Min(axis);
846        float areaLeft = 0, areaRight = 0;
847
848        SortableEntryContainer::const_iterator currentPos =
849                mSubdivisionCandidates->begin();
850       
851        vector<float> bordersRight;
852
853        // we keep track of both borders of the bounding boxes =>
854        // store the events in descending order
855
856        bordersRight.resize(mSubdivisionCandidates->size());
857
858        SortableEntryContainer::reverse_iterator rcit =
859                mSubdivisionCandidates->rbegin(), rcit_end =
860                mSubdivisionCandidates->rend();
861
862        vector<float>::reverse_iterator rbit = bordersRight.rbegin();
863
864        for (; rcit != rcit_end; ++ rcit, ++ rbit)
865        {
866                Intersectable *obj = (*rcit).mObject;
867                const AxisAlignedBox3 obox = obj->GetBox();
868
869                if (obox.Min(axis) < minBorder)
870                {
871                        minBorder = obox.Min(axis);
872                }
873
874                (*rbit) = minBorder;
875        }
876
877        // temporary surface areas
878        float al = 0;
879        float ar = boxArea;
880
881        vector<float>::const_iterator bit = bordersRight.begin();
882        SortableEntryContainer::const_iterator cit, cit_end = mSubdivisionCandidates->end();
883
884        for (cit = mSubdivisionCandidates->begin(); cit != cit_end; ++ cit, ++ bit)
885        {
886                Intersectable *obj = (*cit).mObject;
887
888                const float renderCost = mViewCellsManager->EvalRenderCost(obj);
889               
890                objectsLeft += renderCost;
891                objectsRight -= renderCost;
892
893                const AxisAlignedBox3 obox = obj->GetBox();
894
895                // the borders of the bounding boxes have changed
896                if (obox.Max(axis) > maxBorder)
897                {
898                        maxBorder = obox.Max(axis);
899                }
900
901                minBorder = (*bit);
902
903                AxisAlignedBox3 lbox = nodeBbox;
904                AxisAlignedBox3 rbox = nodeBbox;
905
906                lbox.SetMax(axis, maxBorder);
907                rbox.SetMin(axis, minBorder);
908
909                al = lbox.SurfaceArea();
910                ar = rbox.SurfaceArea();
911
912                const bool noValidSplit = ((objectsLeft <= Limits::Small) || (objectsRight <= Limits::Small));
913                const float sum =  noValidSplit ? 1e25 : objectsLeft * al + objectsRight * ar;
914     
915                /*cout << "pos=" << (*cit).mPos << "\t q=(" << objectsLeft << "," << objectsRight <<")\t r=("
916                         << lbox.SurfaceArea() << "," << rbox.SurfaceArea() << ")" << endl;
917                cout << "minborder: " << minBorder << " maxborder: " << maxBorder << endl;
918            cout << "cost= " << sum << endl;*/
919       
920                if (sum < minSum)
921                {       
922                        minSum = sum;
923                        areaLeft = al;
924                        areaRight = ar;
925
926                        // objects belong to left side now
927                        for (; currentPos != (cit + 1); ++ currentPos);
928                }
929        }
930
931        ////////////
932        //-- assign object to front and back volume
933
934        // belongs to back bv
935        for (cit = mSubdivisionCandidates->begin(); cit != currentPos; ++ cit)
936                objectsBack.push_back((*cit).mObject);
937       
938        // belongs to front bv
939        for (cit = currentPos; cit != cit_end; ++ cit)
940                objectsFront.push_back((*cit).mObject);
941
942        float newCost = minSum / boxArea;
943        float ratio = newCost / totalRenderCost;
944 
945#ifdef GTP_DEBUG
946        cout << "\n\nobjects=(" << (int)objectsBack.size() << "," << (int)objectsFront.size() << " of "
947                 << (int)tData.mNode->mObjects.size() << ")\t area=("
948                 << areaLeft << ", " << areaRight << ", " << boxArea << ")" << endl
949                 << "cost= " << newCost << " oldCost=" << totalRenderCost / boxArea << endl;
950#endif
951
952        return ratio;
953}
954
955#else
956
957float BvHierarchy::EvalSah(const BvhTraversalData &tData,
958                                                   const int axis,
959                                                   ObjectContainer &objectsFront,
960                                                   ObjectContainer &objectsBack)
961{
962        // go through the lists, count the number of objects left and right
963        // and evaluate the following cost funcion:
964        // C = ct_div_ci  + (ol + or) / queries
965        PrepareLocalSubdivisionCandidates(tData, axis);
966
967        const float totalRenderCost = EvalAbsCost(tData.mNode->mObjects);
968        float objectsLeft = 0, objectsRight = totalRenderCost;
969 
970        const AxisAlignedBox3 nodeBbox = tData.mNode->GetBoundingBox();
971
972        const float minBox = nodeBbox.Min(axis);
973        const float maxBox = nodeBbox.Max(axis);
974        const float boxArea = nodeBbox.SurfaceArea();
975
976        float minSum = 1e20f;
977 
978        Vector3 minBorder = nodeBbox.Max();
979        Vector3 maxBorder = nodeBbox.Min();
980
981        float areaLeft = 0, areaRight = 0;
982
983        SortableEntryContainer::const_iterator currentPos =
984                mSubdivisionCandidates->begin();
985       
986        vector<Vector3> bordersRight;
987
988        // we keep track of both borders of the bounding boxes =>
989        // store the events in descending order
990        bordersRight.resize(mSubdivisionCandidates->size());
991
992        SortableEntryContainer::reverse_iterator rcit =
993                mSubdivisionCandidates->rbegin(), rcit_end =
994                mSubdivisionCandidates->rend();
995
996        vector<Vector3>::reverse_iterator rbit = bordersRight.rbegin();
997
998        for (; rcit != rcit_end; ++ rcit, ++ rbit)
999        {
1000                Intersectable *obj = (*rcit).mObject;
1001                const AxisAlignedBox3 obox = obj->GetBox();
1002
1003                for (int i = 0; i < 3; ++ i)
1004                {
1005                        if (obox.Min(i) < minBorder[i])
1006                        {
1007                                minBorder[i] = obox.Min(i);
1008                        }
1009                }
1010
1011                (*rbit) = minBorder;
1012        }
1013
1014        // temporary surface areas
1015        float al = 0;
1016        float ar = boxArea;
1017
1018        vector<Vector3>::const_iterator bit = bordersRight.begin();
1019        SortableEntryContainer::const_iterator cit, cit_end =
1020                mSubdivisionCandidates->end();
1021
1022        for (cit = mSubdivisionCandidates->begin(); cit != cit_end; ++ cit, ++ bit)
1023        {
1024                Intersectable *obj = (*cit).mObject;
1025
1026                const float renderCost = mViewCellsManager->EvalRenderCost(obj);
1027               
1028                objectsLeft += renderCost;
1029                objectsRight -= renderCost;
1030
1031                const AxisAlignedBox3 obox = obj->GetBox();
1032
1033                AxisAlignedBox3 lbox = nodeBbox;
1034                AxisAlignedBox3 rbox = nodeBbox;
1035       
1036                // the borders of the left bounding box have changed
1037                for (int i = 0; i < 3; ++ i)
1038                {
1039                        if (obox.Max(i) > maxBorder[i])
1040                        {
1041                                maxBorder[i] = obox.Max(i);
1042                        }
1043                }
1044
1045                minBorder = (*bit);
1046
1047                lbox.SetMax(maxBorder);
1048                rbox.SetMin(minBorder);
1049
1050                al = lbox.SurfaceArea();
1051                ar = rbox.SurfaceArea();
1052       
1053                const bool noValidSplit = ((objectsLeft <= Limits::Small) || (objectsRight <= Limits::Small));
1054                const float sum =  noValidSplit ? 1e25 : objectsLeft * al + objectsRight * ar;
1055     
1056                /*cout << "pos=" << (*cit).mPos << "\t q=(" << objectsLeft << "," << objectsRight <<")\t r=("
1057                         << lbox.SurfaceArea() << "," << rbox.SurfaceArea() << ")" << endl;
1058                cout << "minborder: " << minBorder << " maxborder: " << maxBorder << endl;
1059            cout << "cost= " << sum << endl;*/
1060       
1061                if (sum < minSum)
1062                {       
1063                        minSum = sum;
1064                        areaLeft = al;
1065                        areaRight = ar;
1066
1067                        // objects belong to left side now
1068                        for (; currentPos != (cit + 1); ++ currentPos);
1069                }
1070        }
1071
1072        /////////////
1073        //-- assign object to front and back volume
1074
1075        // belongs to back bv
1076        for (cit = mSubdivisionCandidates->begin(); cit != currentPos; ++ cit)
1077                objectsBack.push_back((*cit).mObject);
1078       
1079        // belongs to front bv
1080        for (cit = currentPos; cit != cit_end; ++ cit)
1081                objectsFront.push_back((*cit).mObject);
1082
1083        float newCost = minSum / boxArea;
1084        float ratio = newCost / totalRenderCost;
1085 
1086#ifdef GTP_DEBUG
1087        cout << "\n\nobjects=(" << (int)objectsBack.size() << "," << (int)objectsFront.size() << " of "
1088                 << (int)tData.mNode->mObjects.size() << ")\t area=("
1089                 << areaLeft << ", " << areaRight << ", " << boxArea << ")" << endl
1090                 << "cost= " << newCost << " oldCost=" << totalRenderCost / boxArea << endl;
1091#endif
1092
1093        return ratio;
1094}
1095
1096#endif
1097
1098static bool PrepareOutput(const int axis,
1099                                                  const int leaves,
1100                                                  ofstream &sumStats,
1101                                                  ofstream &vollStats,
1102                                                  ofstream &volrStats)
1103{
1104        if ((axis == 0) && (leaves > 0) && (leaves < 90))
1105        {
1106                char str[64];   
1107                sprintf(str, "tmp/bvh_heur_sum-%04d.log", leaves);
1108                sumStats.open(str);
1109                sprintf(str, "tmp/bvh_heur_voll-%04d.log", leaves);
1110                vollStats.open(str);
1111                sprintf(str, "tmp/bvh_heur_volr-%04d.log", leaves);
1112                volrStats.open(str);
1113        }
1114
1115        return sumStats.is_open() && vollStats.is_open() && volrStats.is_open();
1116}
1117
1118
1119static void PrintHeuristics(const float objectsRight,
1120                                                        const float sum,
1121                                                        const float volLeft,
1122                                                        const float volRight,
1123                                                        const float viewSpaceVol,
1124                                                        ofstream &sumStats,
1125                                                        ofstream &vollStats,
1126                                                        ofstream &volrStats)
1127{
1128        sumStats
1129                << "#Position\n" << objectsRight << endl
1130                << "#Sum\n" << sum / viewSpaceVol << endl
1131                << "#Vol\n" << (volLeft +  volRight) / viewSpaceVol << endl;
1132
1133        vollStats
1134                << "#Position\n" << objectsRight << endl
1135                << "#Vol\n" << volLeft / viewSpaceVol << endl;
1136
1137        volrStats
1138                << "#Position\n" << objectsRight << endl
1139                << "#Vol\n" << volRight / viewSpaceVol << endl;
1140}
1141
1142
1143float BvHierarchy::EvalLocalCostHeuristics(const BvhTraversalData &tData,
1144                                                                                   const int axis,
1145                                                                                   ObjectContainer &objectsFront,
1146                                                                                   ObjectContainer &objectsBack)
1147{
1148        ///////////////////////////////////////////////////////
1149        //-- go through the lists, count the number of objects left
1150        //-- and right and evaluate the cost funcion
1151
1152        // prepare the heuristics by setting mailboxes and counters.
1153        const float totalVol = PrepareHeuristics(tData, axis);
1154       
1155        // local helper variables
1156        float volLeft = 0;
1157        float volRight = totalVol;
1158       
1159        const float nTotalObjects = EvalAbsCost(tData.mNode->mObjects);
1160        float nObjectsLeft = 0;
1161        float nObjectsRight = nTotalObjects;
1162
1163        const float viewSpaceVol = mViewCellsManager->GetViewSpaceBox().GetVolume();
1164
1165        SortableEntryContainer::const_iterator backObjectsStart =
1166                mSubdivisionCandidates->begin();
1167
1168        /////////////////////////////////
1169        //-- the parameters for the current optimum
1170
1171        float volBack = volLeft;
1172        float volFront = volRight;
1173        float newRenderCost = nTotalObjects * totalVol;
1174
1175#ifdef GTP_DEBUG
1176        ofstream sumStats;
1177        ofstream vollStats;
1178        ofstream volrStats;
1179
1180        const bool printStats =
1181                PrepareOutput(axis, mBvhStats.Leaves(), sumStats, vollStats, volrStats);
1182#endif
1183
1184        ///////////////////////////////////////////////////
1185        //-- the sweep heuristics
1186        //-- traverse through events and find best split plane
1187
1188        SortableEntryContainer::const_iterator cit,
1189                cit_end = cit_end = mSubdivisionCandidates->end();
1190
1191        for (cit = mSubdivisionCandidates->begin(); cit != cit_end; ++ cit)
1192        {
1193                Intersectable *object = (*cit).mObject;
1194       
1195                // evaluate change in l and r volume
1196                // voll = view cells that see only left node (i.e., left pvs)
1197                // volr = view cells that see only right node (i.e., right pvs)
1198                EvalHeuristicsContribution(object, volLeft, volRight);
1199
1200                const float rc = mViewCellsManager->EvalRenderCost(object);
1201
1202                nObjectsLeft += rc;
1203                nObjectsRight -= rc;
1204               
1205                const bool noValidSplit = ((nObjectsLeft <= Limits::Small) || (nObjectsRight <= Limits::Small));
1206
1207                // the heuristics
1208            const float sum = noValidSplit ?
1209                        1e25 : volLeft * (float)nObjectsLeft + volRight * (float)nObjectsRight;
1210
1211#ifdef GTP_DEBUG
1212                if (printStats)
1213                {
1214                        PrintHeuristics(nObjectsRight, sum, volLeft, volRight, viewSpaceVol,
1215                                                        sumStats, vollStats, volrStats);
1216                }
1217#endif
1218
1219                if (sum < newRenderCost)
1220                {
1221                        newRenderCost = sum;
1222
1223                        volBack = volLeft;
1224                        volFront = volRight;
1225
1226                        // objects belongs to left side now
1227                        for (; backObjectsStart != (cit + 1); ++ backObjectsStart);
1228                }
1229        }
1230
1231        ////////////////////////////////////////////
1232        //-- assign object to front and back volume
1233
1234        // belongs to back bv
1235        for (cit = mSubdivisionCandidates->begin(); cit != backObjectsStart; ++ cit)
1236        {
1237                objectsBack.push_back((*cit).mObject);
1238        }
1239        // belongs to front bv
1240        for (cit = backObjectsStart; cit != cit_end; ++ cit)
1241        {
1242                objectsFront.push_back((*cit).mObject);
1243        }
1244
1245        // render cost of the old parent
1246        const float oldRenderCost = (float)nTotalObjects * totalVol + Limits::Small;
1247        // the relative cost ratio
1248        const float ratio = newRenderCost / oldRenderCost;
1249
1250#ifdef GTP_DEBUG
1251        Debug << "\n§§§§ bvh eval const decrease §§§§" << endl
1252                  << "back pvs: " << (int)objectsBack.size() << " front pvs: "
1253                  << (int)objectsFront.size() << " total pvs: " << nTotalObjects << endl
1254                  << "back p: " << volBack / viewSpaceVol << " front p "
1255                  << volFront / viewSpaceVol << " p: " << totalVol / viewSpaceVol << endl
1256                  << "old rc: " << oldRenderCost / viewSpaceVol << " new rc: "
1257                  << newRenderCost / viewSpaceVol << endl
1258                  << "render cost decrease: "
1259                  << oldRenderCost / viewSpaceVol - newRenderCost / viewSpaceVol << endl;
1260#endif
1261
1262        return ratio;
1263}
1264
1265
1266void BvHierarchy::PrepareLocalSubdivisionCandidates(const BvhTraversalData &tData,
1267                                                                                                        const int axis)                                                                                 
1268{
1269        //-- insert object queries
1270        ObjectContainer *objects = mUseGlobalSorting ?
1271                tData.mSortedObjects[axis] : &tData.mNode->mObjects;
1272
1273        CreateLocalSubdivisionCandidates(*objects, &mSubdivisionCandidates, !mUseGlobalSorting, axis);
1274}
1275
1276
1277void BvHierarchy::CreateLocalSubdivisionCandidates(const ObjectContainer &objects,
1278                                                                                                  SortableEntryContainer **subdivisionCandidates,
1279                                                                                                  const bool sort,
1280                                                                                                  const int axis)
1281{
1282        (*subdivisionCandidates)->clear();
1283
1284        // compute requested size and look if subdivision candidate has to be recomputed
1285        const int requestedSize = (int)objects.size() * 2;
1286       
1287        // creates a sorted split candidates array
1288        if ((*subdivisionCandidates)->capacity() > 500000 &&
1289                requestedSize < (int)((*subdivisionCandidates)->capacity() / 10) )
1290        {
1291        delete (*subdivisionCandidates);
1292                (*subdivisionCandidates) = new SortableEntryContainer;
1293        }
1294
1295        (*subdivisionCandidates)->reserve(requestedSize);
1296
1297        ObjectContainer::const_iterator oit, oit_end = objects.end();
1298
1299        for (oit = objects.begin(); oit < oit_end; ++ oit)
1300        {
1301                Intersectable *object = *oit;
1302                const AxisAlignedBox3 &box = object->GetBox();
1303                const float midPt = (box.Min(axis) + box.Max(axis)) * 0.5f;
1304
1305                (*subdivisionCandidates)->push_back(SortableEntry(object, midPt));
1306        }
1307
1308        if (sort)
1309        {       // no presorted candidate list
1310                stable_sort((*subdivisionCandidates)->begin(), (*subdivisionCandidates)->end());
1311        }
1312}
1313
1314
1315const BvhStatistics &BvHierarchy::GetStatistics() const
1316{
1317        return mBvhStats;
1318}
1319
1320
1321float BvHierarchy::PrepareHeuristics(const BvhTraversalData &tData, const int axis)
1322{       
1323        BvhLeaf *leaf = tData.mNode;
1324        float vol = 0;
1325
1326    // sort so we can use a sweep from right to left
1327        PrepareLocalSubdivisionCandidates(tData, axis);
1328       
1329        // collect and mark the view cells as belonging to front pvs
1330        ViewCellContainer viewCells;
1331        CollectViewCells(tData.mNode->mObjects, viewCells, true);
1332                       
1333        ViewCellContainer::const_iterator vit, vit_end = viewCells.end();
1334        for (vit = viewCells.begin(); vit != vit_end; ++ vit)
1335        {
1336#if USE_VOLUMES_FOR_HEURISTICS
1337                const float volIncr = (*vit)->GetVolume();
1338#else
1339                const float volIncr = 1.0f;
1340#endif
1341                vol += volIncr;
1342        }
1343
1344        // we will mail view cells switching to the back side
1345        ViewCell::NewMail();
1346       
1347        return vol;
1348}
1349
1350///////////////////////////////////////////////////////////
1351
1352
1353void BvHierarchy::EvalHeuristicsContribution(Intersectable *obj,
1354                                                                                         float &volLeft,
1355                                                                                         float &volRight)
1356{
1357        // collect all view cells associated with this objects
1358        // (also multiple times, if they are pierced by several rays)
1359        ViewCellContainer viewCells;
1360        const bool useMailboxing = false;
1361
1362        CollectViewCells(obj, viewCells, useMailboxing);
1363
1364        // classify view cells and compute volume contri accordingly
1365        // possible view cell classifications:
1366        // view cell mailed => view cell can be seen from left child node
1367        // view cell counter > 0 view cell can be seen from right child node
1368        // combined: view cell volume belongs to both nodes
1369        ViewCellContainer::const_iterator vit, vit_end = viewCells.end();
1370       
1371        for (vit = viewCells.begin(); vit != vit_end; ++ vit)
1372        {
1373                // view cells can also be seen from left child node
1374                ViewCell *viewCell = *vit;
1375#if USE_VOLUMES_FOR_HEURISTICS
1376                const float vol = viewCell->GetVolume();
1377#else
1378                const float vol = 1.0f;
1379#endif
1380                if (!viewCell->Mailed())
1381                {
1382                        viewCell->Mail();
1383                        // we now see view cell from both nodes
1384                        // => add volume to left node
1385                        volLeft += vol;
1386                }
1387
1388                // last reference into the right node
1389                if (-- viewCell->mCounter == 0)
1390                {       
1391                        // view cell was previously seen from both nodes  =>
1392                        // remove volume from right node
1393                        volRight -= vol;
1394                }
1395        }
1396}
1397
1398
1399void BvHierarchy::SetViewCellsManager(ViewCellsManager *vcm)
1400{
1401        mViewCellsManager = vcm;
1402}
1403
1404
1405AxisAlignedBox3 BvHierarchy::GetBoundingBox() const
1406{
1407        return mBoundingBox;
1408}
1409
1410
1411float BvHierarchy::SelectObjectPartition(const BvhTraversalData &tData,
1412                                                                                 ObjectContainer &frontObjects,
1413                                                                                 ObjectContainer &backObjects,
1414                                                                                 bool useVisibilityBasedHeuristics)
1415{
1416        ObjectContainer nFrontObjects[3];
1417        ObjectContainer nBackObjects[3];
1418        float nCostRatio[3];
1419
1420        int sAxis = 0;
1421        int bestAxis = -1;
1422
1423        if (mOnlyDrivingAxis)
1424        {
1425                const AxisAlignedBox3 box = tData.mNode->GetBoundingBox();
1426                sAxis = box.Size().DrivingAxis();
1427        }
1428       
1429        ////////////////////////////////////
1430        //-- evaluate split cost for all three axis
1431       
1432        for (int axis = 0; axis < 3; ++ axis)
1433        {
1434                if (!mOnlyDrivingAxis || (axis == sAxis))
1435                {
1436                        if (mUseCostHeuristics)
1437                        {
1438                                //////////////////////////////////
1439                //-- split objects using heuristics
1440                               
1441                                if (useVisibilityBasedHeuristics)
1442                                {
1443                                        //cout << "v rays: " << tData.mNumRays << " " << endl;
1444                                        ///////////
1445                                        //-- heuristics using objects weighted by view cells volume
1446                                        nCostRatio[axis] =
1447                                                EvalLocalCostHeuristics(tData,
1448                                                                                                axis,
1449                                                                                                nFrontObjects[axis],
1450                                                                                                nBackObjects[axis]);
1451                                }
1452                                else
1453                                {       //cout << "e rays: " << tData.mNumRays << " ";
1454                                        //////////////////
1455                                        //-- view cells not constructed yet     => use surface area heuristic                   
1456                                        nCostRatio[axis] = EvalSah(tData,
1457                                                                                           axis,
1458                                                                                           nFrontObjects[axis],
1459                                                                                           nBackObjects[axis]);
1460                                }
1461                        }
1462                        else
1463                        {
1464                                //-- split objects using some simple criteria
1465                                nCostRatio[axis] =
1466                                        EvalLocalObjectPartition(tData, axis, nFrontObjects[axis], nBackObjects[axis]);
1467                        }
1468
1469                        if ((bestAxis == -1) || (nCostRatio[axis] < nCostRatio[bestAxis]))
1470                        {
1471                                bestAxis = axis;
1472                        }
1473                }
1474        }
1475
1476    ////////////////
1477        //-- assign values
1478
1479        frontObjects = nFrontObjects[bestAxis];
1480        backObjects = nBackObjects[bestAxis];
1481
1482        //cout << "val: " << nCostRatio[bestAxis] << " axis: " << bestAxis << endl;
1483        return nCostRatio[bestAxis];
1484}
1485
1486
1487int BvHierarchy::AssociateObjectsWithRays(const VssRayContainer &rays) const
1488{
1489        int nRays = 0;
1490        VssRayContainer::const_iterator rit, rit_end = rays.end();
1491
1492        VssRay::NewMail();
1493
1494    for (rit = rays.begin(); rit != rays.end(); ++ rit)
1495        {
1496                VssRay *ray = (*rit);
1497
1498                if (ray->mTerminationObject)
1499                {
1500                        ray->mTerminationObject->GetOrCreateRays()->push_back(ray);
1501                        if (!ray->Mailed())
1502                        {
1503                                ray->Mail();
1504                                ++ nRays;
1505                        }
1506                }
1507#if COUNT_ORIGIN_OBJECTS
1508                if (ray->mOriginObject)
1509                {
1510                        ray->mOriginObject->GetOrCreateRays()->push_back(ray);
1511
1512                        if (!ray->Mailed())
1513                        {
1514                                ray->Mail();
1515                                ++ nRays;
1516                        }
1517                }
1518#endif
1519        }
1520
1521        return nRays;
1522}
1523
1524
1525void BvHierarchy::PrintSubdivisionStats(const SubdivisionCandidate &sc)
1526{
1527        const float costDecr = sc.GetRenderCostDecrease();     
1528
1529        mSubdivisionStats
1530                        << "#Leaves\n" << mBvhStats.Leaves() << endl
1531                        << "#RenderCostDecrease\n" << costDecr << endl
1532                        << "#TotalRenderCost\n" << mTotalCost << endl
1533                        << "#EntriesInPvs\n" << mPvsEntries << endl;
1534}
1535
1536
1537void BvHierarchy::CollectRays(const ObjectContainer &objects,
1538                                                          VssRayContainer &rays) const
1539{
1540        VssRay::NewMail();
1541        ObjectContainer::const_iterator oit, oit_end = objects.end();
1542
1543        // evaluate reverse pvs and view cell volume on left and right cell
1544        // note: should I take all leaf objects or rather the objects hit by rays?
1545        for (oit = objects.begin(); oit != oit_end; ++ oit)
1546        {
1547                Intersectable *obj = *oit;
1548                VssRayContainer::const_iterator rit, rit_end = obj->GetOrCreateRays()->end();
1549
1550                for (rit = obj->GetOrCreateRays()->begin(); rit < rit_end; ++ rit)
1551                {
1552                        VssRay *ray = (*rit);
1553
1554                        if (!ray->Mailed())
1555                        {
1556                                ray->Mail();
1557                                rays.push_back(ray);
1558                        }
1559                }
1560        }
1561}
1562
1563
1564float BvHierarchy::EvalAbsCost(const ObjectContainer &objects)// const
1565{
1566#if USE_BETTER_RENDERCOST_EST
1567        ObjectContainer::const_iterator oit, oit_end = objects.end();
1568
1569        for (oit = objects.begin(); oit != oit_end; ++ oit)
1570        {
1571                objRenderCost += ViewCellsManager::GetRendercost(*oit);
1572        }
1573#else
1574        return (float)objects.size();
1575#endif
1576}
1577
1578
1579float BvHierarchy::EvalSahCost(BvhLeaf *leaf)
1580{
1581        ////////////////
1582        //-- surface area heuristics
1583        if (leaf->mObjects.empty())
1584                return 0.0f;
1585
1586        const AxisAlignedBox3 box = GetBoundingBox(leaf);
1587        const float area = box.SurfaceArea();
1588        const float viewSpaceArea = mViewCellsManager->GetViewSpaceBox().SurfaceArea();
1589
1590        return EvalAbsCost(leaf->mObjects) * area / viewSpaceArea;
1591}
1592
1593
1594float BvHierarchy::EvalRenderCost(const ObjectContainer &objects) const
1595{       
1596        ///////////////
1597        //-- render cost heuristics
1598
1599        const float viewSpaceVol = mViewCellsManager->GetViewSpaceBox().GetVolume();
1600
1601        // probability that view point lies in a view cell which sees this node
1602        const float p = EvalViewCellsVolume(objects) / viewSpaceVol;
1603    const float objRenderCost = EvalAbsCost(objects);
1604       
1605        return objRenderCost * p;
1606}
1607
1608
1609AxisAlignedBox3 BvHierarchy::EvalBoundingBox(const ObjectContainer &objects,
1610                                                                                         const AxisAlignedBox3 *parentBox) const
1611{
1612        // if there are no objects in this box, box size is set to parent box size.
1613        // Question: Invalidate box instead?
1614        if (parentBox && objects.empty())
1615                return *parentBox;
1616
1617        AxisAlignedBox3 box;
1618        box.Initialize();
1619
1620        ObjectContainer::const_iterator oit, oit_end = objects.end();
1621
1622        for (oit = objects.begin(); oit != oit_end; ++ oit)
1623        {
1624                Intersectable *obj = *oit;
1625                // grow bounding box to include all objects
1626                box.Include(obj->GetBox());
1627        }
1628
1629        return box;
1630}
1631
1632
1633void BvHierarchy::CollectLeaves(BvhNode *root, vector<BvhLeaf *> &leaves) const
1634{
1635        stack<BvhNode *> nodeStack;
1636        nodeStack.push(root);
1637
1638        while (!nodeStack.empty())
1639        {
1640                BvhNode *node = nodeStack.top();
1641                nodeStack.pop();
1642
1643                if (node->IsLeaf())
1644                {
1645                        BvhLeaf *leaf = (BvhLeaf *)node;
1646                        leaves.push_back(leaf);
1647                }
1648                else
1649                {
1650                        BvhInterior *interior = (BvhInterior *)node;
1651
1652                        nodeStack.push(interior->GetBack());
1653                        nodeStack.push(interior->GetFront());
1654                }
1655        }
1656}
1657
1658
1659AxisAlignedBox3 BvHierarchy::GetBoundingBox(BvhNode *node) const
1660{
1661        return node->GetBoundingBox();
1662}
1663
1664
1665void BvHierarchy::CollectViewCells(const ObjectContainer &objects,
1666                                                                   ViewCellContainer &viewCells,
1667                                                                   const bool setCounter) const
1668{
1669        // no view cells yet
1670        if (0 && mHierarchyManager->GetViewSpaceSubdivisionType() ==
1671                HierarchyManager::NO_VIEWSPACE_SUBDIV)
1672                return;
1673
1674        ViewCell::NewMail();
1675        ObjectContainer::const_iterator oit, oit_end = objects.end();
1676
1677        // loop through all object and collect view cell pvs of this node
1678        for (oit = objects.begin(); oit != oit_end; ++ oit)
1679        {
1680                CollectViewCells(*oit, viewCells, true, setCounter);
1681        }
1682}
1683
1684
1685void BvHierarchy::CollectViewCells(Intersectable *obj,
1686                                                                   ViewCellContainer &viewCells,
1687                                                                   const bool useMailBoxing,
1688                                                                   const bool setCounter) const
1689{
1690        VssRayContainer::const_iterator rit, rit_end = obj->GetOrCreateRays()->end();
1691
1692        for (rit = obj->GetOrCreateRays()->begin(); rit < rit_end; ++ rit)
1693        {
1694                VssRay *ray = (*rit);
1695                ViewCellContainer tmpViewCells;
1696       
1697                mHierarchyManager->mVspTree->GetViewCells(*ray, tmpViewCells);
1698
1699                // matt: probably slow to allocate memory for view cells every time
1700                ViewCellContainer::const_iterator vit, vit_end = tmpViewCells.end();
1701
1702                for (vit = tmpViewCells.begin(); vit != vit_end; ++ vit)
1703                {
1704                        ViewCell *vc = *vit;
1705
1706                        // store view cells
1707                        if (!useMailBoxing || !vc->Mailed())
1708                        {
1709                                if (useMailBoxing)
1710                                {
1711                                        vc->Mail();
1712                                        if (setCounter)
1713                                        {
1714                                                vc->mCounter = 0;
1715                                        }
1716                                }
1717                                viewCells.push_back(vc);
1718                        }
1719                       
1720                        if (setCounter)
1721                        {
1722                                ++ vc->mCounter;
1723                        }
1724                }
1725        }
1726}
1727
1728
1729int BvHierarchy::CountViewCells(Intersectable *obj) const
1730{
1731        int result = 0;
1732       
1733        VssRayContainer::const_iterator rit, rit_end = obj->GetOrCreateRays()->end();
1734
1735        for (rit = obj->GetOrCreateRays()->begin(); rit < rit_end; ++ rit)
1736        {
1737                VssRay *ray = (*rit);
1738                ViewCellContainer tmpViewCells;
1739       
1740                mHierarchyManager->mVspTree->GetViewCells(*ray, tmpViewCells);
1741               
1742                ViewCellContainer::const_iterator vit, vit_end = tmpViewCells.end();
1743                for (vit = tmpViewCells.begin(); vit != vit_end; ++ vit)
1744                {
1745                        ViewCell *vc = *vit;
1746
1747                        // store view cells
1748                        if (!vc->Mailed())
1749                        {
1750                                vc->Mail();
1751                                ++ result;
1752                        }
1753                }
1754        }
1755
1756        return result;
1757}
1758
1759
1760int BvHierarchy::CountViewCells(const ObjectContainer &objects) const
1761{
1762        int nViewCells = 0;
1763        ViewCell::NewMail();
1764        ObjectContainer::const_iterator oit, oit_end = objects.end();
1765
1766        // loop through all object and collect view cell pvs of this node
1767        for (oit = objects.begin(); oit != oit_end; ++ oit)
1768        {
1769                nViewCells += CountViewCells(*oit);
1770        }
1771
1772        return nViewCells;
1773}
1774
1775
1776void BvHierarchy::CollectDirtyCandidates(BvhSubdivisionCandidate *sc,
1777                                                                                 vector<SubdivisionCandidate *> &dirtyList,
1778                                                                                 const bool onlyUnmailed)
1779{
1780        BvhTraversalData &tData = sc->mParentData;
1781        BvhLeaf *node = tData.mNode;
1782       
1783        ViewCellContainer viewCells;
1784        ViewCell::NewMail();
1785        CollectViewCells(node->mObjects, viewCells, true);
1786
1787        if (0) cout << "collected " << (int)viewCells.size() << " dirty candidates" << endl;
1788       
1789        // split candidates handling
1790        // these view cells  are thrown into dirty list
1791        ViewCellContainer::const_iterator vit, vit_end = viewCells.end();
1792
1793        for (vit = viewCells.begin(); vit != vit_end; ++ vit)
1794        {
1795        VspViewCell *vc = dynamic_cast<VspViewCell *>(*vit);
1796                VspLeaf *leaf = vc->mLeaves[0];
1797       
1798                SubdivisionCandidate *candidate = leaf->GetSubdivisionCandidate();
1799               
1800                // is this leaf still a split candidate?
1801                if (candidate && (!onlyUnmailed || !candidate->Mailed()))
1802                {
1803                        candidate->Mail();
1804                        dirtyList.push_back(candidate);
1805                }
1806        }
1807}
1808
1809
1810BvhNode *BvHierarchy::GetRoot() const
1811{
1812        return mRoot;
1813}
1814
1815
1816bool BvHierarchy::IsObjectInLeaf(BvhLeaf *leaf, Intersectable *object) const
1817{
1818        ObjectContainer::const_iterator oit =
1819                lower_bound(leaf->mObjects.begin(), leaf->mObjects.end(), object, ilt);
1820                               
1821        // objects sorted by id
1822        if ((oit != leaf->mObjects.end()) && ((*oit)->GetId() == object->GetId()))
1823        {
1824                return true;
1825        }
1826        else
1827        {
1828                return false;
1829        }
1830}
1831
1832
1833BvhLeaf *BvHierarchy::GetLeaf(Intersectable *object, BvhNode *node) const
1834{
1835        // rather use the simple version
1836        if (!object)
1837                return NULL;
1838       
1839        return object->mBvhLeaf;
1840       
1841        ///////////////////////////////////////
1842        // start from root of tree
1843       
1844        if (node == NULL)
1845                node = mRoot;
1846       
1847        vector<BvhLeaf *> leaves;
1848
1849        stack<BvhNode *> nodeStack;
1850        nodeStack.push(node);
1851 
1852        BvhLeaf *leaf = NULL;
1853 
1854        while (!nodeStack.empty()) 
1855        {
1856                BvhNode *node = nodeStack.top();
1857                nodeStack.pop();
1858       
1859                if (node->IsLeaf())
1860                {
1861                        leaf = dynamic_cast<BvhLeaf *>(node);
1862
1863                        if (IsObjectInLeaf(leaf, object))
1864                        {
1865                                return leaf;
1866                        }
1867                }
1868                else   
1869                {       
1870                        // find point
1871                        BvhInterior *interior = dynamic_cast<BvhInterior *>(node);
1872       
1873                        if (interior->GetBack()->GetBoundingBox().Includes(object->GetBox()))
1874                        {
1875                                nodeStack.push(interior->GetBack());
1876                        }
1877                       
1878                        // search both sides as we are using bounding volumes
1879                        if (interior->GetFront()->GetBoundingBox().Includes(object->GetBox()))
1880                        {
1881                                nodeStack.push(interior->GetFront());
1882                        }
1883                }
1884        }
1885 
1886        return leaf;
1887}
1888
1889
1890BvhIntersectable *BvHierarchy::GetOrCreateBvhIntersectable(BvhLeaf *node)
1891{
1892        // search nodes
1893        std::map<BvhLeaf *, BvhIntersectable *>::const_iterator
1894                it = mBvhIntersectables.find(node);
1895
1896        if (it != mBvhIntersectables.end())
1897        {
1898                return (*it).second;
1899        }
1900
1901        // not in map => create new entry
1902        BvhIntersectable *bvhObj = new BvhIntersectable(node);
1903        mBvhIntersectables[node] = bvhObj;
1904
1905        return bvhObj;
1906}
1907
1908
1909bool BvHierarchy::Export(OUT_STREAM &stream)
1910{
1911        ExportNode(mRoot, stream);
1912
1913        return true;
1914}
1915
1916
1917void BvHierarchy::ExportObjects(BvhLeaf *leaf, OUT_STREAM &stream)
1918{
1919        ObjectContainer::const_iterator oit, oit_end = leaf->mObjects.end();
1920        for (oit = leaf->mObjects.begin(); oit != oit_end; ++ oit)
1921        {
1922                stream << (*oit)->GetId() << " ";
1923        }
1924}
1925
1926
1927void BvHierarchy::ExportNode(BvhNode *node, OUT_STREAM &stream)
1928{
1929        if (node->IsLeaf())
1930        {
1931                BvhLeaf *leaf = dynamic_cast<BvhLeaf *>(node);
1932                const AxisAlignedBox3 box = leaf->GetBoundingBox();
1933                stream << "<Leaf"
1934                           << " min=\"" << box.Min().x << " " << box.Min().y << " " << box.Min().z << "\""
1935                           << " max=\"" << box.Max().x << " " << box.Max().y << " " << box.Max().z << "\""
1936                           << " objects=\"";
1937               
1938                //-- export objects
1939                ExportObjects(leaf, stream);
1940               
1941                stream << "\" />" << endl;
1942        }
1943        else
1944        {       
1945                BvhInterior *interior = dynamic_cast<BvhInterior *>(node);
1946                const AxisAlignedBox3 box = interior->GetBoundingBox();
1947
1948                stream << "<Interior"
1949                           << " min=\"" << box.Min().x << " " << box.Min().y << " " << box.Min().z << "\""
1950                           << " max=\"" << box.Max().x << " " << box.Max().y << " " << box.Max().z
1951                           << "\">" << endl;
1952
1953                ExportNode(interior->GetBack(), stream);
1954                ExportNode(interior->GetFront(), stream);
1955
1956                stream << "</Interior>" << endl;
1957        }
1958}
1959
1960
1961float BvHierarchy::EvalViewCellsVolume(const ObjectContainer &objects) const
1962{
1963        float vol = 0;
1964
1965        ViewCellContainer viewCells;
1966        CollectViewCells(objects, viewCells);
1967
1968        ViewCellContainer::const_iterator vit, vit_end = viewCells.end();
1969
1970        for (vit = viewCells.begin(); vit != vit_end; ++ vit)
1971        {
1972                vol += (*vit)->GetVolume();
1973        }
1974
1975        return vol;
1976}
1977
1978
1979void BvHierarchy::Initialise(const ObjectContainer &objects)
1980{
1981        AxisAlignedBox3 box = EvalBoundingBox(objects);
1982
1983        ///////
1984        //-- create new root
1985
1986        BvhLeaf *bvhleaf = new BvhLeaf(box, NULL, (int)objects.size());
1987        bvhleaf->mObjects = objects;
1988        mRoot = bvhleaf;
1989
1990        // compute bounding box from objects
1991        mBoundingBox = mRoot->GetBoundingBox();
1992
1993        // associate root with current objects
1994        AssociateObjectsWithLeaf(bvhleaf);
1995}
1996
1997
1998/*
1999Mesh *BvHierarchy::MergeLeafToMesh()
2000{
2001        vector<BvhLeaf *> leaves;
2002        CollectLeaves(leaves);
2003
2004        vector<BvhLeaf *>::const_iterator lit, lit_end = leaves.end();
2005
2006        for (lit = leaves.begin(); lit != lit_end; ++ lit)
2007        {
2008                Mesh *mesh = MergeLeafToMesh(*lit);
2009        }
2010}*/
2011
2012
2013SubdivisionCandidate *BvHierarchy::PrepareConstruction(const VssRayContainer &sampleRays,
2014                                                                                                           const ObjectContainer &objects)
2015{
2016        ///////////////////////////////////////
2017        //-- we assume that we have objects sorted by their id =>
2018        //-- we don't have to sort them here and an binary search
2019        //-- for identifying if a object is in a leaf.
2020       
2021        mBvhStats.Reset();
2022        mBvhStats.Start();
2023        mBvhStats.nodes = 1;
2024               
2025        // store pointer to this tree
2026        BvhSubdivisionCandidate::sBvHierarchy = this;
2027       
2028        // root and bounding box was already constructed
2029        BvhLeaf *bvhLeaf = dynamic_cast<BvhLeaf *>(mRoot);
2030
2031        // multiply termination criterium for comparison,
2032        // so it can be set between zero and one and
2033        // no division is necessary during traversal
2034
2035#if PROBABILIY_IS_BV_VOLUME
2036        mTermMinProbability *= mBoundingBox.GetVolume();
2037        // probability that bounding volume is seen
2038        const float prop = GetBoundingBox().GetVolume();
2039#else
2040        mTermMinProbability *= mVspTree->GetBoundingBox().GetVolume();
2041        // probability that volume is "seen" from the view cells
2042        const float prop = EvalViewCellsVolume(objects);
2043#endif
2044
2045        // only rays intersecting objects in node are interesting
2046        const int nRays = AssociateObjectsWithRays(sampleRays);
2047        //Debug << "using " << nRays << " of " << (int)sampleRays.size() << " rays" << endl;
2048
2049        // create bvh traversal data
2050        BvhTraversalData oData(bvhLeaf, 0, prop, nRays);
2051
2052        // create sorted object lists for the first data
2053        if (mUseGlobalSorting)
2054        {
2055                AssignInitialSortedObjectList(oData);
2056        }
2057       
2058
2059        ///////////////////
2060        //-- add first candidate for object space partition     
2061
2062        BvhSubdivisionCandidate *oSubdivisionCandidate = new BvhSubdivisionCandidate(oData);
2063
2064        // evaluate priority
2065        EvalSubdivisionCandidate(*oSubdivisionCandidate);
2066        bvhLeaf->SetSubdivisionCandidate(oSubdivisionCandidate);
2067
2068        mTotalCost = EvalRenderCost(objects);
2069        mPvsEntries = CountViewCells(objects);
2070
2071        PrintSubdivisionStats(*oSubdivisionCandidate);
2072       
2073        return oSubdivisionCandidate;
2074}
2075
2076
2077void BvHierarchy::AssignInitialSortedObjectList(BvhTraversalData &tData)
2078{
2079        // we sort the objects as a preprocess so they don't have
2080        // to be sorted for each split
2081        for (int i = 0; i < 3; ++ i)
2082        {
2083                // create new objects
2084                if (!mSortedObjects[i])
2085                {
2086                        mSortedObjects[i] = new SortableEntryContainer();
2087                        CreateLocalSubdivisionCandidates(tData.mNode->mObjects, &mSortedObjects[i], true, i);
2088                }
2089
2090                // copy list into traversal data list
2091                tData.mSortedObjects[i] = new ObjectContainer();
2092                tData.mSortedObjects[i]->reserve((int)mSortedObjects[i]->size());
2093
2094                SortableEntryContainer::const_iterator oit, oit_end = mSortedObjects[i]->end();
2095
2096                for (oit = mSortedObjects[i]->begin(); oit != oit_end; ++ oit)
2097                {
2098                        tData.mSortedObjects[i]->push_back((*oit).mObject);
2099                }
2100        }
2101}
2102
2103
2104void BvHierarchy::AssignSortedObjects(const BvhSubdivisionCandidate &sc,
2105                                                                          BvhTraversalData &frontData,
2106                                                                          BvhTraversalData &backData)
2107{
2108        Intersectable::NewMail();
2109
2110        // we sorted the objects as a preprocess so they don't have
2111        // to be sorted for each split
2112        ObjectContainer::const_iterator fit, fit_end = sc.mFrontObjects.end();
2113
2114        for (fit = sc.mFrontObjects.begin(); fit != fit_end; ++ fit)
2115        {
2116                (*fit)->Mail();
2117        }
2118
2119        for (int i = 0; i < 3; ++ i)
2120        {
2121                frontData.mSortedObjects[i] = new ObjectContainer();
2122                backData.mSortedObjects[i] = new ObjectContainer();
2123
2124                frontData.mSortedObjects[i]->reserve((int)sc.mFrontObjects.size());
2125                backData.mSortedObjects[i]->reserve((int)sc.mFrontObjects.size());
2126
2127                ObjectContainer::const_iterator oit, oit_end = sc.mParentData.mSortedObjects[i]->end();
2128
2129                for (oit = sc.mParentData.mSortedObjects[i]->begin(); oit != oit_end; ++ oit)
2130                {
2131                        if ((*oit)->Mailed())
2132                        {
2133                                frontData.mSortedObjects[i]->push_back(*oit);
2134                        }
2135                        else
2136                        {
2137                                backData.mSortedObjects[i]->push_back(*oit);
2138                        }
2139                }
2140        }
2141}
2142
2143
2144SubdivisionCandidate *BvHierarchy::Reset(const VssRayContainer &sampleRays,
2145                                                                                 const ObjectContainer &objects)
2146{
2147        // reset stats
2148        mBvhStats.Reset();
2149        mBvhStats.Start();
2150        mBvhStats.nodes = 1;
2151
2152        // reset root
2153        DEL_PTR(mRoot);
2154       
2155        BvhLeaf *bvhleaf = new BvhLeaf(mBoundingBox, NULL, (int)objects.size());
2156        bvhleaf->mObjects = objects;
2157        mRoot = bvhleaf;
2158       
2159#if PROBABILIY_IS_BV_VOLUME
2160        mTermMinProbability *= mBoundingBox.GetVolume();
2161        // probability that bounding volume is seen
2162        const float prop = GetBoundingBox().GetVolume();
2163#else
2164        mTermMinProbability *= mVspTree->GetBoundingBox().GetVolume();
2165        // probability that volume is "seen" from the view cells
2166        const float prop = EvalViewCellsVolume(objects);
2167#endif
2168
2169        const int nRays = CountRays(objects);
2170        BvhLeaf *bvhLeaf = dynamic_cast<BvhLeaf *>(mRoot);
2171
2172        // create bvh traversal data
2173        BvhTraversalData oData(bvhLeaf, 0, prop, nRays);
2174
2175        AssignInitialSortedObjectList(oData);
2176       
2177
2178        ///////////////////
2179        //-- add first candidate for object space partition     
2180
2181        BvhSubdivisionCandidate *oSubdivisionCandidate =
2182                new BvhSubdivisionCandidate(oData);
2183
2184        EvalSubdivisionCandidate(*oSubdivisionCandidate);
2185        bvhLeaf->SetSubdivisionCandidate(oSubdivisionCandidate);
2186
2187        const float viewSpaceVol = mViewCellsManager->GetViewSpaceBox().GetVolume();
2188        mTotalCost = (float)objects.size() * prop / viewSpaceVol;
2189
2190        PrintSubdivisionStats(*oSubdivisionCandidate);
2191
2192        return oSubdivisionCandidate;
2193}
2194
2195
2196void BvhStatistics::Print(ostream &app) const
2197{
2198        app << "=========== BvHierarchy statistics ===============\n";
2199
2200        app << setprecision(4);
2201
2202        app << "#N_CTIME  ( Construction time [s] )\n" << Time() << " \n";
2203
2204        app << "#N_NODES ( Number of nodes )\n" << nodes << "\n";
2205
2206        app << "#N_INTERIORS ( Number of interior nodes )\n" << Interior() << "\n";
2207
2208        app << "#N_LEAVES ( Number of leaves )\n" << Leaves() << "\n";
2209
2210        app << "#AXIS_ALIGNED_SPLITS (number of axis aligned splits)\n" << splits << endl;
2211
2212        app << "#N_MAXCOSTNODES  ( Percentage of leaves with terminated because of max cost ratio )\n"
2213                << maxCostNodes * 100 / (double)Leaves() << endl;
2214
2215        app << "#N_PMINPROBABILITYLEAVES  ( Percentage of leaves with mininum probability )\n"
2216                << minProbabilityNodes * 100 / (double)Leaves() << endl;
2217
2218
2219        //////////////////////////////////////////////////
2220       
2221        app << "#N_PMAXDEPTHLEAVES ( Percentage of leaves at maximum depth )\n"
2222                <<      maxDepthNodes * 100 / (double)Leaves() << endl;
2223       
2224        app << "#N_PMAXDEPTH ( Maximal reached depth )\n" << maxDepth << endl;
2225
2226        app << "#N_PMINDEPTH ( Minimal reached depth )\n" << minDepth << endl;
2227
2228        app << "#AVGDEPTH ( average depth )\n" << AvgDepth() << endl;
2229
2230       
2231        ////////////////////////////////////////////////////////
2232       
2233        app << "#N_PMINOBJECTSLEAVES  ( Percentage of leaves with mininum objects )\n"
2234                << minObjectsNodes * 100 / (double)Leaves() << endl;
2235
2236        app << "#N_MAXOBJECTREFS  ( Max number of object refs / leaf )\n" << maxObjectRefs << "\n";
2237
2238        app << "#N_MINOBJECTREFS  ( Min number of object refs / leaf )\n" << minObjectRefs << "\n";
2239
2240        app << "#N_EMPTYLEAFS ( Empty leafs )\n" << emptyNodes << "\n";
2241       
2242        app << "#N_PAVGOBJECTSLEAVES  ( average object refs / leaf)\n" << AvgObjectRefs() << endl;
2243
2244
2245        ////////////////////////////////////////////////////////
2246       
2247        app << "#N_PMINRAYSLEAVES  ( Percentage of leaves with mininum rays )\n"
2248                << minRaysNodes * 100 / (double)Leaves() << endl;
2249
2250        app << "#N_MAXRAYREFS  ( Max number of ray refs / leaf )\n" << maxRayRefs << "\n";
2251
2252        app << "#N_MINRAYREFS  ( Min number of ray refs / leaf )\n" << minRayRefs << "\n";
2253       
2254        app << "#N_PAVGRAYLEAVES  ( average ray refs / leaf )\n" << AvgRayRefs() << endl;
2255       
2256        app << "#N_PAVGRAYCONTRIBLEAVES  ( Average ray contribution)\n" <<
2257                rayRefs / (double)objectRefs << endl;
2258
2259        app << "#N_PMAXRAYCONTRIBLEAVES  ( Percentage of leaves with maximal ray contribution )\n"<<
2260                maxRayContriNodes * 100 / (double)Leaves() << endl;
2261
2262        app << "#N_PGLOBALCOSTMISSES ( Global cost misses )\n" << mGlobalCostMisses << endl;
2263
2264        app << "========== END OF BvHierarchy statistics ==========\n";
2265}
2266
2267
2268// TODO: return memory usage in MB
2269float BvHierarchy::GetMemUsage() const
2270{
2271        return (float)(sizeof(BvHierarchy)
2272                                   + mBvhStats.Leaves() * sizeof(BvhLeaf)
2273                                   + mBvhStats.Interior() * sizeof(BvhInterior)
2274                                   ) / float(1024 * 1024);
2275}
2276
2277
2278void BvHierarchy::SetActive(BvhNode *node) const
2279{
2280        vector<BvhLeaf *> leaves;
2281
2282        // sets the pointers to the currently active view cells
2283        CollectLeaves(node, leaves);
2284        vector<BvhLeaf *>::const_iterator lit, lit_end = leaves.end();
2285
2286        for (lit = leaves.begin(); lit != lit_end; ++ lit)
2287        {
2288                (*lit)->SetActiveNode(node);
2289        }
2290}
2291
2292
2293BvhNode *BvHierarchy::SubdivideAndCopy(SplitQueue &tQueue,
2294                                                                           SubdivisionCandidate *splitCandidate)
2295{
2296        BvhSubdivisionCandidate *sc =
2297                dynamic_cast<BvhSubdivisionCandidate *>(splitCandidate);
2298        BvhTraversalData &tData = sc->mParentData;
2299
2300        BvhNode *currentNode = tData.mNode;
2301        BvhNode *oldNode = (BvhNode *)splitCandidate->mEvaluationHack;
2302
2303        if (!oldNode->IsLeaf())
2304        {       
2305                //////////////
2306                //-- continue subdivision
2307
2308                BvhTraversalData tFrontData;
2309                BvhTraversalData tBackData;
2310                       
2311                BvhInterior *oldInterior = dynamic_cast<BvhInterior *>(oldNode);
2312               
2313                sc->mFrontObjects.clear();
2314                sc->mBackObjects.clear();
2315
2316                oldInterior->GetFront()->CollectObjects(sc->mFrontObjects);
2317                oldInterior->GetBack()->CollectObjects(sc->mBackObjects);
2318               
2319                // evaluate the changes in render cost and pvs entries
2320                EvalSubdivisionCandidate(*sc, false);
2321
2322                // create new interior node and two leaf node
2323                currentNode = SubdivideNode(*sc, tFrontData, tBackData);
2324       
2325                //oldNode->mRenderCostDecr += sc->GetRenderCostDecrease();
2326                //oldNode->mPvsEntriesIncr += sc->GetPvsEntriesIncr();
2327               
2328                oldNode->mRenderCostDecr = sc->GetRenderCostDecrease();
2329                oldNode->mPvsEntriesIncr = sc->GetPvsEntriesIncr();
2330               
2331                ///////////////////////////
2332                //-- push the new split candidates on the queue
2333               
2334                BvhSubdivisionCandidate *frontCandidate = new BvhSubdivisionCandidate(tFrontData);
2335                BvhSubdivisionCandidate *backCandidate = new BvhSubdivisionCandidate(tBackData);
2336
2337                frontCandidate->SetPriority((float)-oldInterior->GetFront()->mTimeStamp);
2338                backCandidate->SetPriority((float)-oldInterior->GetBack()->mTimeStamp);
2339
2340                frontCandidate->mEvaluationHack = oldInterior->GetFront();
2341                backCandidate->mEvaluationHack = oldInterior->GetBack();
2342
2343                // cross reference
2344                tFrontData.mNode->SetSubdivisionCandidate(frontCandidate);
2345                tBackData.mNode->SetSubdivisionCandidate(backCandidate);
2346
2347                //cout << "f: " << frontCandidate->GetPriority() << " b: " << backCandidate->GetPriority() << endl;
2348                tQueue.Push(frontCandidate);
2349                tQueue.Push(backCandidate);
2350        }
2351
2352        /////////////////////////////////
2353        //-- node is a leaf => terminate traversal
2354
2355        if (currentNode->IsLeaf())
2356        {
2357                // this leaf is no candidate for splitting anymore
2358                // => detach subdivision candidate
2359                tData.mNode->SetSubdivisionCandidate(NULL);
2360                // detach node so we don't delete it with the traversal data
2361                tData.mNode = NULL;
2362        }
2363       
2364        return currentNode;
2365}
2366
2367
2368void BvHierarchy::CollectObjects(const AxisAlignedBox3 &box, ObjectContainer &objects)
2369{
2370        stack<BvhNode *> nodeStack;
2371
2372        nodeStack.push(mRoot);
2373
2374        while (!nodeStack.empty())
2375        {
2376                BvhNode *node = nodeStack.top();
2377
2378                nodeStack.pop();
2379
2380                if (node->IsLeaf())
2381                {
2382                        BvhLeaf *leaf = (BvhLeaf *)node;
2383
2384                        ObjectContainer::const_iterator oit, oit_end = leaf->mObjects.end();
2385
2386                        for (oit = leaf->mObjects.begin(); oit != oit_end; ++oit)
2387                        {
2388                                Intersectable *object = *oit;
2389                                if (Overlap(box, object->GetBox()))
2390                                {
2391                                        object->Mail();
2392                                        objects.push_back(object);
2393                                }
2394                        }
2395                }
2396                else
2397                {
2398                        BvhInterior *interior = (BvhInterior *)node;
2399
2400                        if (Overlap(box, interior->GetBoundingBox()))
2401                                nodeStack.push(interior->GetFront());
2402
2403                        if (Overlap(box, interior->GetBoundingBox()))
2404                                nodeStack.push(interior->GetBack());
2405                }
2406        }
2407}
2408
2409}
Note: See TracBrowser for help on using the repository browser.