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

Revision 1758, 64.8 KB checked in by mattausch, 18 years ago (diff)

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