source: GTP/trunk/Lib/Vis/Preprocessing/src/VspTree.cpp @ 1308

Revision 1308, 69.9 KB checked in by mattausch, 18 years ago (diff)

changed hierarchy construction: up to a certain level view space split is taken

Line 
1#include <stack>
2#include <time.h>
3#include <iomanip>
4
5#include "ViewCell.h"
6#include "Plane3.h"
7#include "VspTree.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 "KdTree.h"
19#include "KdIntersectable.h"
20#include "HierarchyManager.h"
21#include "BvHierarchy.h"
22#include "OspTree.h"
23
24
25
26namespace GtpVisibilityPreprocessor {
27
28
29#define USE_FIXEDPOINT_T 0
30
31
32//-- static members
33
34VspTree *VspTree::VspSubdivisionCandidate::sVspTree = NULL;
35int VspNode::sMailId = 1;
36
37// variable for debugging volume contribution for heuristics
38static float debugVol;
39
40
41// pvs penalty can be different from pvs size
42inline static float EvalPvsPenalty(const int pvs,
43                                                                   const int lower,
44                                                                   const int upper)
45{
46        // clamp to minmax values
47        if (pvs < lower)
48        {
49                return (float)lower;
50        }
51        else if (pvs > upper)
52        {
53                return (float)upper;
54        }
55        return (float)pvs;
56}
57
58
59static bool ViewCellHasMultipleReferences(Intersectable *obj,
60                                                                                  ViewCell *vc,
61                                                                                  bool checkOnlyMailed)
62{
63        MailablePvsData *vdata = obj->mViewCellPvs.Find(vc);
64//return false;
65        if (vdata)
66        {
67                // more than one view cell sees this object inside different kd cells
68                if (!checkOnlyMailed || !vdata->Mailed())
69                {
70                        if (checkOnlyMailed)
71                                vdata->Mail();
72                        //Debug << "sumpdf: " << vdata->mSumPdf << endl;
73                        if (vdata->mSumPdf > 1.5f)
74                                return true;
75                }
76        }
77
78        return false;
79}
80
81
82void VspTreeStatistics::Print(ostream &app) const
83{
84        app << "=========== VspTree statistics ===============\n";
85
86        app << setprecision(4);
87
88        app << "#N_CTIME  ( Construction time [s] )\n" << Time() << " \n";
89
90        app << "#N_NODES ( Number of nodes )\n" << nodes << "\n";
91
92        app << "#N_INTERIORS ( Number of interior nodes )\n" << Interior() << "\n";
93
94        app << "#N_LEAVES ( Number of leaves )\n" << Leaves() << "\n";
95
96        app << "#AXIS_ALIGNED_SPLITS (number of axis aligned splits)\n" << splits[0] + splits[1] + splits[2] << endl;
97
98        app << "#N_SPLITS ( Number of splits in axes x y z)\n";
99
100        for (int i = 0; i < 3; ++ i)
101                app << splits[i] << " ";
102        app << endl;
103
104        app << "#N_PMAXDEPTHLEAVES ( Percentage of leaves at maximum depth )\n"
105                <<      maxDepthNodes * 100 / (double)Leaves() << endl;
106
107        app << "#N_PMINPVSLEAVES  ( Percentage of leaves with mininimal PVS )\n"
108                << minPvsNodes * 100 / (double)Leaves() << endl;
109
110        app << "#N_PMINRAYSLEAVES  ( Percentage of leaves with minimal number of rays)\n"
111                << minRaysNodes * 100 / (double)Leaves() << endl;
112
113        app << "#N_MAXCOSTNODES  ( Percentage of leaves with terminated because of max cost ratio )\n"
114                << maxCostNodes * 100 / (double)Leaves() << endl;
115
116        app << "#N_PMINPROBABILITYLEAVES  ( Percentage of leaves with mininum probability )\n"
117                << minProbabilityNodes * 100 / (double)Leaves() << endl;
118
119        app << "#N_PMAXRAYCONTRIBLEAVES  ( Percentage of leaves with maximal ray contribution )\n"
120                <<      maxRayContribNodes * 100 / (double)Leaves() << endl;
121
122        app << "#N_PMAXDEPTH ( Maximal reached depth )\n" << maxDepth << endl;
123
124        app << "#N_PMINDEPTH ( Minimal reached depth )\n" << minDepth << endl;
125
126        app << "#AVGDEPTH ( average depth )\n" << AvgDepth() << endl;
127
128        app << "#N_INVALIDLEAVES (number of invalid leaves )\n" << invalidLeaves << endl;
129
130        app << "#N_RAYS (number of rays / leaf)\n" << AvgRays() << endl;
131        //app << "#N_PVS: " << pvs << endl;
132
133        //app << "#N_MAXOBJECTREFS  ( Max number of object refs / leaf )\n" << maxObjectRefs << "\n";
134
135        app << "========== END OF VspTree statistics ==========\n";
136}
137
138
139
140/******************************************************************/
141/*                  class VspNode implementation                  */
142/******************************************************************/
143
144
145VspNode::VspNode():
146mParent(NULL), mTreeValid(true), mTimeStamp(0)
147{}
148
149
150VspNode::VspNode(VspInterior *parent):
151mParent(parent), mTreeValid(true)
152{}
153
154
155bool VspNode::IsRoot() const
156{
157        return mParent == NULL;
158}
159
160
161VspInterior *VspNode::GetParent()
162{
163        return mParent;
164}
165
166
167void VspNode::SetParent(VspInterior *parent)
168{
169        mParent = parent;
170}
171
172
173bool VspNode::IsSibling(VspNode *n) const
174{
175        return  ((this != n) && mParent &&
176                         (mParent->GetFront() == n) || (mParent->GetBack() == n));
177}
178
179
180int VspNode::GetDepth() const
181{
182        int depth = 0;
183        VspNode *p = mParent;
184       
185        while (p)
186        {
187                p = p->mParent;
188                ++ depth;
189        }
190
191        return depth;
192}
193
194
195bool VspNode::TreeValid() const
196{
197        return mTreeValid;
198}
199
200
201void VspNode::SetTreeValid(const bool v)
202{
203        mTreeValid = v;
204}
205
206
207/****************************************************************/
208/*              class VspInterior implementation                */
209/****************************************************************/
210
211
212VspInterior::VspInterior(const AxisAlignedPlane &plane):
213mPlane(plane), mFront(NULL), mBack(NULL)
214{}
215
216
217VspInterior::~VspInterior()
218{
219        DEL_PTR(mFront);
220        DEL_PTR(mBack);
221}
222
223
224bool VspInterior::IsLeaf() const
225{
226        return false;
227}
228
229
230VspNode *VspInterior::GetBack()
231{
232        return mBack;
233}
234
235
236VspNode *VspInterior::GetFront()
237{
238        return mFront;
239}
240
241
242AxisAlignedPlane VspInterior::GetPlane() const
243{
244        return mPlane;
245}
246
247
248float VspInterior::GetPosition() const
249{
250        return mPlane.mPosition;
251}
252
253
254int VspInterior::GetAxis() const
255{
256        return mPlane.mAxis;
257}
258
259
260void VspInterior::ReplaceChildLink(VspNode *oldChild, VspNode *newChild)
261{
262        if (mBack == oldChild)
263                mBack = newChild;
264        else
265                mFront = newChild;
266}
267
268
269void VspInterior::SetupChildLinks(VspNode *front, VspNode *back)
270{
271    mBack = back;
272    mFront = front;
273}
274
275
276AxisAlignedBox3 VspInterior::GetBoundingBox() const
277{
278        return mBoundingBox;
279}
280
281
282void VspInterior::SetBoundingBox(const AxisAlignedBox3 &box)
283{
284        mBoundingBox = box;
285}
286
287
288int VspInterior::Type() const
289{
290        return Interior;
291}
292
293
294
295/****************************************************************/
296/*                  class VspLeaf implementation                */
297/****************************************************************/
298
299
300VspLeaf::VspLeaf(): mViewCell(NULL), mPvs(NULL), mSubdivisionCandidate(NULL)
301{
302}
303
304
305VspLeaf::~VspLeaf()
306{
307        DEL_PTR(mPvs);
308        CLEAR_CONTAINER(mVssRays);
309}
310
311
312int VspLeaf::Type() const
313{
314        return Leaf;
315}
316
317
318VspLeaf::VspLeaf(ViewCellLeaf *viewCell):
319mViewCell(viewCell)
320{
321}
322
323
324VspLeaf::VspLeaf(VspInterior *parent):
325VspNode(parent), mViewCell(NULL), mPvs(NULL)
326{}
327
328
329
330VspLeaf::VspLeaf(VspInterior *parent, ViewCellLeaf *viewCell):
331VspNode(parent), mViewCell(viewCell), mPvs(NULL)
332{
333}
334
335ViewCellLeaf *VspLeaf::GetViewCell() const
336{
337        return mViewCell;
338}
339
340void VspLeaf::SetViewCell(ViewCellLeaf *viewCell)
341{
342        mViewCell = viewCell;
343}
344
345
346bool VspLeaf::IsLeaf() const
347{
348        return true;
349}
350
351
352/*************************************************************************/
353/*                       class VspTree implementation                    */
354/*************************************************************************/
355
356
357VspTree::VspTree():
358mRoot(NULL),
359mOutOfBoundsCell(NULL),
360mStoreRays(false),
361mTimeStamp(1),
362mHierarchyManager(NULL)
363{
364        bool randomize = false;
365        Environment::GetSingleton()->GetBoolValue("VspTree.Construction.randomize", randomize);
366        if (randomize)
367                Randomize(); // initialise random generator for heuristics
368
369        //-- termination criteria for autopartition
370        Environment::GetSingleton()->GetIntValue("VspTree.Termination.maxDepth", mTermMaxDepth);
371        Environment::GetSingleton()->GetIntValue("VspTree.Termination.minPvs", mTermMinPvs);
372        Environment::GetSingleton()->GetIntValue("VspTree.Termination.minRays", mTermMinRays);
373        Environment::GetSingleton()->GetFloatValue("VspTree.Termination.minProbability", mTermMinProbability);
374        Environment::GetSingleton()->GetFloatValue("VspTree.Termination.maxRayContribution", mTermMaxRayContribution);
375       
376        Environment::GetSingleton()->GetIntValue("VspTree.Termination.missTolerance", mTermMissTolerance);
377        Environment::GetSingleton()->GetIntValue("VspTree.Termination.maxViewCells", mMaxViewCells);
378
379        //-- max cost ratio for early tree termination
380        Environment::GetSingleton()->GetFloatValue("VspTree.Termination.maxCostRatio", mTermMaxCostRatio);
381
382        Environment::GetSingleton()->GetFloatValue("VspTree.Termination.minGlobalCostRatio", mTermMinGlobalCostRatio);
383        Environment::GetSingleton()->GetIntValue("VspTree.Termination.globalCostMissTolerance", mTermGlobalCostMissTolerance);
384
385        //-- factors for bsp tree split plane heuristics
386        Environment::GetSingleton()->GetFloatValue("VspTree.Termination.ct_div_ci", mCtDivCi);
387
388        //-- partition criteria
389        Environment::GetSingleton()->GetFloatValue("VspTree.Construction.epsilon", mEpsilon);
390        Environment::GetSingleton()->GetFloatValue("VspTree.Construction.renderCostDecreaseWeight", mRenderCostDecreaseWeight);
391
392        // if only the driving axis is used for axis aligned split
393        Environment::GetSingleton()->GetBoolValue("VspTree.splitUseOnlyDrivingAxis", mOnlyDrivingAxis);
394       
395        Environment::GetSingleton()->GetIntValue("VspTree.maxTests", mMaxTests);
396        Environment::GetSingleton()->GetFloatValue("VspTree.maxStaticMemory", mMaxMemory);
397
398        Environment::GetSingleton()->GetBoolValue("VspTree.useCostHeuristics", mUseCostHeuristics);
399        Environment::GetSingleton()->GetBoolValue("VspTree.simulateOctree", mCirculatingAxis);
400       
401        //Environment::GetSingleton()->GetBoolValue("VspTree.useKdPvsForHeuristics", mUseKdPvsForHeuristics);
402
403        char subdivisionStatsLog[100];
404        Environment::GetSingleton()->GetStringValue("VspTree.subdivisionStats", subdivisionStatsLog);
405        mSubdivisionStats.open(subdivisionStatsLog);
406
407        Environment::GetSingleton()->GetFloatValue("VspTree.Construction.minBand", mMinBand);
408        Environment::GetSingleton()->GetFloatValue("VspTree.Construction.maxBand", mMaxBand);
409       
410
411        //-- debug output
412
413        Debug << "******* VSP options ******** " << endl;
414
415    Debug << "max depth: " << mTermMaxDepth << endl;
416        Debug << "min PVS: " << mTermMinPvs << endl;
417        Debug << "min probabiliy: " << mTermMinProbability << endl;
418        Debug << "min rays: " << mTermMinRays << endl;
419        Debug << "max ray contri: " << mTermMaxRayContribution << endl;
420        Debug << "max cost ratio: " << mTermMaxCostRatio << endl;
421        Debug << "miss tolerance: " << mTermMissTolerance << endl;
422        Debug << "max view cells: " << mMaxViewCells << endl;
423        Debug << "randomize: " << randomize << endl;
424
425        Debug << "min global cost ratio: " << mTermMinGlobalCostRatio << endl;
426        Debug << "global cost miss tolerance: " << mTermGlobalCostMissTolerance << endl;
427        Debug << "only driving axis: " << mOnlyDrivingAxis << endl;
428        Debug << "max memory: " << mMaxMemory << endl;
429        Debug << "use cost heuristics: " << mUseCostHeuristics << endl;
430        Debug << "subdivision stats log: " << subdivisionStatsLog << endl;
431        Debug << "render cost decrease weight: " << mRenderCostDecreaseWeight << endl;
432
433        Debug << "circulating axis: " << mCirculatingAxis << endl;
434        Debug << "minband: " << mMinBand << endl;
435        Debug << "maxband: " << mMaxBand << endl;
436
437        mLocalSubdivisionCandidates = new vector<SortableEntry>;
438
439        Debug << endl;
440}
441
442
443VspViewCell *VspTree::GetOutOfBoundsCell()
444{
445        return mOutOfBoundsCell;
446}
447
448
449VspViewCell *VspTree::GetOrCreateOutOfBoundsCell()
450{
451        if (!mOutOfBoundsCell)
452        {
453                mOutOfBoundsCell = new VspViewCell();
454                mOutOfBoundsCell->SetId(-1);
455                mOutOfBoundsCell->SetValid(false);
456        }
457
458        return mOutOfBoundsCell;
459}
460
461
462const VspTreeStatistics &VspTree::GetStatistics() const
463{
464        return mVspStats;
465}
466
467
468VspTree::~VspTree()
469{
470        DEL_PTR(mRoot);
471        DEL_PTR(mLocalSubdivisionCandidates);
472}
473
474
475void VspTree::ComputeBoundingBox(const VssRayContainer &rays,
476                                                                 AxisAlignedBox3 *forcedBoundingBox)
477{       
478        if (forcedBoundingBox)
479        {
480                mBoundingBox = *forcedBoundingBox;
481        }
482        else // compute vsp tree bounding box
483        {
484                mBoundingBox.Initialize();
485                VssRayContainer::const_iterator rit, rit_end = rays.end();
486
487                //-- compute bounding box
488        for (rit = rays.begin(); rit != rit_end; ++ rit)
489                {
490                        VssRay *ray = *rit;
491
492                        // compute bounding box of view space
493                        mBoundingBox.Include(ray->GetTermination());
494                        mBoundingBox.Include(ray->GetOrigin());
495                }
496        }
497}
498
499
500void VspTree::AddSubdivisionStats(const int viewCells,
501                                                                  const float renderCostDecr,
502                                                                  const float totalRenderCost,
503                                                                  const float avgRenderCost)
504{
505        mSubdivisionStats
506                        << "#ViewCells\n" << viewCells << endl
507                        << "#RenderCostDecrease\n" << renderCostDecr << endl
508                        << "#TotalRenderCost\n" << totalRenderCost << endl
509                        << "#AvgRenderCost\n" << avgRenderCost << endl;
510}
511
512
513// TODO: return memory usage in MB
514float VspTree::GetMemUsage() const
515{
516        return (float)
517                 (sizeof(VspTree) +
518                  mVspStats.Leaves() * sizeof(VspLeaf) +
519                  mCreatedViewCells * sizeof(VspViewCell) +
520                  mVspStats.pvs * sizeof(PvsData) +
521                  mVspStats.Interior() * sizeof(VspInterior) +
522                  mVspStats.accumRays * sizeof(RayInfo)) / (1024.0f * 1024.0f);
523}
524
525
526inline bool VspTree::LocalTerminationCriteriaMet(const VspTraversalData &data) const
527{
528        const bool localTerminationCriteriaMet = (0
529                //|| ((int)data.mRays->size() <= mTermMinRays)
530                //|| (data.mPvs <= mTermMinPvs)
531                //|| (data.mProbability <= mTermMinProbability)
532                //|| (data.GetAvgRayContribution() > mTermMaxRayContribution)
533                //|| (data.mDepth >= mTermMaxDepth)
534                );
535
536        if (1 && localTerminationCriteriaMet)
537        {
538                Debug << "********local termination *********" << endl;
539                Debug << "rays: " << (int)data.mRays->size() << "  " << mTermMinRays << endl;
540                Debug << "pvs: " << data.mPvs << " " << mTermMinPvs << endl;
541                Debug << "p: " <<  data.mProbability << " " << mTermMinProbability << endl;
542                Debug << "avg contri: " << data.GetAvgRayContribution() << " " << mTermMaxRayContribution << endl;
543                Debug << "depth " << data.mDepth << " " << mTermMaxDepth << endl;
544        }
545
546        return localTerminationCriteriaMet;             
547}
548
549
550inline bool VspTree::GlobalTerminationCriteriaMet(const VspTraversalData &data) const
551{
552        const bool terminationCriteriaMet = (0
553                // || mOutOfMemory
554                || (mVspStats.Leaves() >= mMaxViewCells)
555        //|| (mGlobalCostMisses >= mTermGlobalCostMissTolerance)
556                );
557
558        if (0 && terminationCriteriaMet)
559        {
560                Debug << "********* terminationCriteriaMet *********" << endl;
561                Debug << "cost misses: " << mGlobalCostMisses << " " << mTermGlobalCostMissTolerance << endl;
562                Debug << "leaves: " << mVspStats.Leaves() << " " <<  mMaxViewCells << endl;
563        }
564        return terminationCriteriaMet;
565}
566
567
568void VspTree::CreateViewCell(VspTraversalData &tData, const bool updatePvs)
569{
570        //-- create new view cell
571        VspLeaf *leaf = dynamic_cast<VspLeaf *>(tData.mNode);
572       
573        VspViewCell *viewCell = new VspViewCell();
574    leaf->SetViewCell(viewCell);
575       
576        int conSamp = 0;
577        float sampCon = 0.0f;
578
579        if (updatePvs)
580        {
581                //-- update pvs of view cell
582                AddSamplesToPvs(leaf, *tData.mRays, sampCon, conSamp);
583
584                // update scalar pvs size value
585                ObjectPvs &pvs = viewCell->GetPvs();
586                mViewCellsManager->UpdateScalarPvsSize(viewCell, pvs.CountObjectsInPvs(), pvs.GetSize());
587
588                mVspStats.contributingSamples += conSamp;
589                mVspStats.sampleContributions += (int)sampCon;
590        }
591
592
593        //-- store additional info
594        if (mStoreRays)
595        {
596                RayInfoContainer::const_iterator it, it_end = tData.mRays->end();
597
598                for (it = tData.mRays->begin(); it != it_end; ++ it)
599                {
600                        (*it).mRay->Ref();                     
601                        leaf->mVssRays.push_back((*it).mRay);
602                }
603        }
604               
605        // set view cell values
606        viewCell->mLeaf = leaf;
607
608        viewCell->SetVolume(tData.mProbability);
609    leaf->mProbability = tData.mProbability;
610}
611
612
613void VspTree::EvalSubdivisionStats(const SubdivisionCandidate &sc)
614{
615        const float costDecr = sc.GetRenderCostDecrease();
616       
617        AddSubdivisionStats(mVspStats.Leaves(),
618                                                costDecr,
619                                                mTotalCost,
620                                                (float)mTotalPvsSize / (float)mVspStats.Leaves());
621}
622
623
624VspNode *VspTree::Subdivide(SplitQueue &tQueue,
625                                                        SubdivisionCandidate *splitCandidate,
626                                                        const bool globalCriteriaMet)
627{
628        // todo remove dynamic cast
629        VspSubdivisionCandidate *sc =
630                dynamic_cast<VspSubdivisionCandidate *>(splitCandidate);
631
632        VspTraversalData &tData = sc->mParentData;
633        VspNode *newNode = tData.mNode;
634
635        if (!LocalTerminationCriteriaMet(tData) && !globalCriteriaMet)
636        {       
637                //-- continue subdivision
638                VspTraversalData tFrontData;
639                VspTraversalData tBackData;
640               
641                // create new interior node and two leaf node
642                const AxisAlignedPlane splitPlane = sc->mSplitPlane;
643                const int maxCostMisses = sc->mMaxCostMisses;
644
645                newNode = SubdivideNode(splitPlane, tData, tFrontData, tBackData);
646       
647                // how often was max cost ratio missed in this branch?
648                tFrontData.mMaxCostMisses = maxCostMisses;
649                tBackData.mMaxCostMisses = maxCostMisses;
650                       
651                mTotalCost -= sc->GetRenderCostDecrease();
652                mTotalPvsSize += tFrontData.mPvs + tBackData.mPvs - tData.mPvs;
653
654                // subdivision statistics
655                if (1) EvalSubdivisionStats(*sc);
656               
657                //-- evaluate new split candidates for global greedy cost heuristics
658
659                VspSubdivisionCandidate *frontCandidate = new VspSubdivisionCandidate(tFrontData);
660                VspSubdivisionCandidate *backCandidate = new VspSubdivisionCandidate(tBackData);
661
662                EvalSubdivisionCandidate(*frontCandidate);
663                EvalSubdivisionCandidate(*backCandidate);
664
665                // cross reference
666                tFrontData.mNode->SetSubdivisionCandidate(frontCandidate);
667                tBackData.mNode->SetSubdivisionCandidate(backCandidate);
668
669                tQueue.Push(frontCandidate);
670                tQueue.Push(backCandidate);
671               
672                // delete old leaf node
673                //DEL_PTR(tData.mNode);
674        }
675
676        if (newNode->IsLeaf()) // subdivision terminated
677        {
678                // view cell is created during subdivision
679                //CreateViewCell(tData);
680
681                VspLeaf *leaf = dynamic_cast<VspLeaf *>(newNode);
682                ViewCell *viewCell = leaf->GetViewCell();
683
684                int conSamp = 0;
685                float sampCon = 0.0f;
686
687#if 0
688                //-- store pvs optained from rays
689                AddSamplesToPvs(leaf, *tData.mRays, sampCon, conSamp);
690
691                // update scalar pvs size value
692                ObjectPvs &pvs = viewCell->GetPvs();
693                mViewCellsManager->UpdateScalarPvsSize(viewCell, pvs.CountObjectsInPvs(), pvs.GetSize());
694
695                mVspStats.contributingSamples += conSamp;
696                mVspStats.sampleContributions += (int)sampCon;
697#endif
698                //-- store additional info
699                if (mStoreRays)
700                {
701                        RayInfoContainer::const_iterator it, it_end = tData.mRays->end();
702                        for (it = tData.mRays->begin(); it != it_end; ++ it)
703                        {
704                                (*it).mRay->Ref();                     
705                                leaf->mVssRays.push_back((*it).mRay);
706                        }
707                }
708
709                // finally evaluate statistics for this leaf
710                EvaluateLeafStats(tData);
711                // detach subdivision candidate: this leaf is no candidate for
712                // splitting anymore
713                tData.mNode->SetSubdivisionCandidate(NULL);
714                // detach node so it won't get deleted
715                tData.mNode = NULL;
716        }
717
718        return newNode;
719}
720
721
722void VspTree::EvalSubdivisionCandidate(VspSubdivisionCandidate &splitCandidate)
723{
724        float frontProb;
725        float backProb;
726       
727        VspLeaf *leaf = dynamic_cast<VspLeaf *>(splitCandidate.mParentData.mNode);
728       
729        // compute locally best split plane
730        const float ratio = SelectSplitPlane(splitCandidate.mParentData,
731                                                                                 splitCandidate.mSplitPlane,
732                                                                                 frontProb,
733                                                                                 backProb);
734
735        const bool maxCostRatioViolated = mTermMaxCostRatio < ratio;
736
737        // max cost threshold violated?
738        splitCandidate.mMaxCostMisses = maxCostRatioViolated  ?
739                        splitCandidate.mParentData.mMaxCostMisses + 1:
740                        splitCandidate.mParentData.mMaxCostMisses;
741       
742        float oldRenderCost;
743
744        // compute global decrease in render cost
745        const float renderCostDecr = EvalRenderCostDecrease(splitCandidate.mSplitPlane,
746                                                                                                                splitCandidate.mParentData,
747                                                                                                                oldRenderCost);
748
749        splitCandidate.SetRenderCostDecrease(renderCostDecr);
750
751#if 1
752        const float priority = (float)-splitCandidate.mParentData.mDepth;
753#else
754        // take render cost of node into account
755        // otherwise danger of being stuck in a local minimum!!
756        const float factor = mRenderCostDecreaseWeight;
757        const float priority = factor * renderCostDecr + (1.0f - factor) * oldRenderCost;
758#endif
759       
760        splitCandidate.SetPriority(priority);
761}
762
763
764VspInterior *VspTree::SubdivideNode(const AxisAlignedPlane &splitPlane,
765                                                                        VspTraversalData &tData,
766                                                                        VspTraversalData &frontData,
767                                                                        VspTraversalData &backData)
768{
769        VspLeaf *leaf = dynamic_cast<VspLeaf *>(tData.mNode);
770       
771        //-- the front and back traversal data is filled with the new values
772
773        frontData.mDepth = tData.mDepth + 1;
774        backData.mDepth = tData.mDepth + 1;
775
776        frontData.mRays = new RayInfoContainer();
777        backData.mRays = new RayInfoContainer();
778
779        //-- subdivide rays
780        SplitRays(splitPlane,
781                          *tData.mRays,
782                          *frontData.mRays,
783                          *backData.mRays);
784
785        //-- compute pvs
786
787        frontData.mPvs = EvalPvsSize(*frontData.mRays);
788        backData.mPvs = EvalPvsSize(*backData.mRays);
789
790        Debug << "f pvs: " << frontData.mPvs << " b pvs: " << backData.mPvs << " pvs " << tData.mPvs << endl;
791
792        // split front and back node geometry and compute area
793        tData.mBoundingBox.Split(splitPlane.mAxis,
794                                                         splitPlane.mPosition,
795                                                         frontData.mBoundingBox,
796                                                         backData.mBoundingBox);
797
798        frontData.mProbability = frontData.mBoundingBox.GetVolume();
799        backData.mProbability = tData.mProbability - frontData.mProbability;
800
801
802    ///////////////////////////////////////////
803        // subdivide further
804
805        // store maximal and minimal depth
806        if (tData.mDepth > mVspStats.maxDepth)
807        {
808                Debug << "max depth increases to " << tData.mDepth
809                          << " at " << mVspStats.Leaves() << " leaves" << endl;
810                mVspStats.maxDepth = tData.mDepth;
811        }
812
813        // two more leaves
814        mVspStats.nodes += 2;
815   
816        VspInterior *interior = new VspInterior(splitPlane);
817
818
819        //-- create front and back leaf
820
821        VspInterior *parent = leaf->GetParent();
822
823        // replace a link from node's parent
824        if (parent)
825        {
826                parent->ReplaceChildLink(leaf, interior);
827                interior->SetParent(parent);
828
829                // remove "parent" view cell from pvs of all objects (traverse trough rays)
830                RemoveParentViewCellReferences(tData.mNode->GetViewCell());
831        }
832        else // new root
833        {
834                mRoot = interior;
835        }
836
837        VspLeaf *frontLeaf = new VspLeaf(interior);
838        VspLeaf *backLeaf = new VspLeaf(interior);
839
840        // and setup child links
841        interior->SetupChildLinks(frontLeaf, backLeaf);
842       
843        // add bounding box
844        interior->SetBoundingBox(tData.mBoundingBox);
845
846        // set front and back leaf
847        frontData.mNode = frontLeaf;
848        backData.mNode = backLeaf;
849
850        // explicitely create front and back view cell
851        CreateViewCell(frontData, false);
852        CreateViewCell(backData, false);
853
854
855#if WORK_WITH_VIEWCELL_PVS
856        // create front and back view cell
857        // add front and back view cell to "Potentially Visbilie View Cells"
858        // of the objects in front and back pvs
859
860        AddViewCellReferences(frontLeaf->GetViewCell());
861        AddViewCellReferences(backLeaf->GetViewCell());
862#endif
863
864        interior->mTimeStamp = mTimeStamp ++;
865       
866        return interior;
867}
868
869
870void VspTree::RemoveParentViewCellReferences(ViewCell *parent) const
871{
872        KdLeaf::NewMail();
873
874        // remove the parents from the object pvss
875        ObjectPvsMap::const_iterator oit, oit_end = parent->GetPvs().mEntries.end();
876
877        for (oit = parent->GetPvs().mEntries.begin(); oit != oit_end; ++ oit)
878        {
879                Intersectable *object = (*oit).first;
880                // HACK: make sure that the view cell is removed from the pvs
881                const float high_contri = 9999999;
882
883                // remove reference count of view cells
884                object->mViewCellPvs.RemoveSample(parent, high_contri);
885        }
886}
887
888
889void VspTree::AddViewCellReferences(ViewCell *vc) const
890{
891        KdLeaf::NewMail();
892
893        // Add front view cell to the object pvsss
894        ObjectPvsMap::const_iterator oit, oit_end = vc->GetPvs().mEntries.end();
895
896        for (oit = vc->GetPvs().mEntries.begin(); oit != oit_end; ++ oit)
897        {
898                Intersectable *object = (*oit).first;
899
900                // increase reference count of view cells
901                object->mViewCellPvs.AddSample(vc, 1);
902        }
903}
904
905
906void VspTree::AddSamplesToPvs(VspLeaf *leaf,
907                                                          const RayInfoContainer &rays,
908                                                          float &sampleContributions,
909                                                          int &contributingSamples)
910{
911        sampleContributions = 0;
912        contributingSamples = 0;
913 
914        RayInfoContainer::const_iterator it, it_end = rays.end();
915 
916        ViewCellLeaf *vc = leaf->GetViewCell();
917
918        // add contributions from samples to the PVS
919        for (it = rays.begin(); it != it_end; ++ it)
920        {
921                float sc = 0.0f;
922                VssRay *ray = (*it).mRay;
923
924                bool madeContrib = false;
925                float contribution;
926
927                Intersectable *obj = ray->mTerminationObject;
928
929                if (obj)
930                {
931                        madeContrib =
932                                mViewCellsManager->AddSampleToPvs(
933                                        obj,
934                                        ray->mTermination,
935                                        vc,
936                                        ray->mPdf,
937                                        contribution);
938
939                        sc += contribution;
940                }
941
942                obj = ray->mOriginObject;
943
944                if (obj)
945                {
946                        madeContrib =
947                                mViewCellsManager->AddSampleToPvs(
948                                        obj,
949                                        ray->mOrigin,
950                                        vc,
951                                        ray->mPdf,
952                                        contribution);
953
954                        sc += contribution;
955                }
956
957                if (madeContrib)
958                {
959                        ++ contributingSamples;
960                }
961
962                // store rays for visualization
963                if (0) leaf->mVssRays.push_back(new VssRay(*ray));
964        }
965}
966
967
968void VspTree::SortSubdivisionCandidates(const RayInfoContainer &rays,
969                                                                  const int axis,
970                                                                  float minBand,
971                                                                  float maxBand)
972{
973        mLocalSubdivisionCandidates->clear();
974
975        int requestedSize = 2 * (int)(rays.size());
976
977        // creates a sorted split candidates array
978        if (mLocalSubdivisionCandidates->capacity() > 500000 &&
979                requestedSize < (int)(mLocalSubdivisionCandidates->capacity() / 10) )
980        {
981        delete mLocalSubdivisionCandidates;
982                mLocalSubdivisionCandidates = new vector<SortableEntry>;
983        }
984
985        mLocalSubdivisionCandidates->reserve(requestedSize);
986
987        float pos;
988        RayInfoContainer::const_iterator rit, rit_end = rays.end();
989
990        //-- insert all queries
991        for (rit = rays.begin(); rit != rit_end; ++ rit)
992        {
993                const bool positive = (*rit).mRay->HasPosDir(axis);
994                               
995                pos = (*rit).ExtrapOrigin(axis);
996               
997                mLocalSubdivisionCandidates->push_back(SortableEntry(positive ? SortableEntry::ERayMin : SortableEntry::ERayMax,
998                                                                        pos, (*rit).mRay));
999
1000                pos = (*rit).ExtrapTermination(axis);
1001
1002                mLocalSubdivisionCandidates->push_back(SortableEntry(positive ? SortableEntry::ERayMax : SortableEntry::ERayMin,
1003                                                                        pos, (*rit).mRay));
1004        }
1005
1006        stable_sort(mLocalSubdivisionCandidates->begin(), mLocalSubdivisionCandidates->end());
1007}
1008
1009
1010int VspTree::PrepareHeuristics(KdLeaf *leaf)
1011{       
1012        int pvsSize = 0;
1013       
1014        if (!leaf->Mailed())
1015        {
1016                leaf->Mail();
1017                leaf->mCounter = 1;
1018                // add objects without the objects which are in several kd leaves
1019                pvsSize += (int)(leaf->mObjects.size() - leaf->mMultipleObjects.size());
1020        }
1021        else
1022        {
1023                ++ leaf->mCounter;
1024        }
1025
1026        //-- the objects belonging to several leaves must be handled seperately
1027        ObjectContainer::const_iterator oit, oit_end = leaf->mMultipleObjects.end();
1028
1029        for (oit = leaf->mMultipleObjects.begin(); oit != oit_end; ++ oit)
1030        {
1031                Intersectable *object = *oit;
1032                                               
1033                if (!object->Mailed())
1034                {
1035                        object->Mail();
1036                        object->mCounter = 1;
1037
1038                        ++ pvsSize;
1039                }
1040                else
1041                {
1042                        ++ object->mCounter;
1043                }
1044        }
1045       
1046        return pvsSize;
1047}
1048
1049
1050int VspTree::PrepareHeuristics(const RayInfoContainer &rays)
1051{       
1052        Intersectable::NewMail();
1053        KdNode::NewMail();
1054        BvhLeaf::NewMail();
1055
1056        int pvsSize = 0;
1057
1058        RayInfoContainer::const_iterator ri, ri_end = rays.end();
1059
1060    // set all kd nodes / objects as belonging to the front pvs
1061        for (ri = rays.begin(); ri != ri_end; ++ ri)
1062        {
1063                VssRay *ray = (*ri).mRay;
1064               
1065                pvsSize += PrepareHeuristics(*ray, true);
1066                pvsSize += PrepareHeuristics(*ray, false);
1067        }
1068
1069        return pvsSize;
1070}
1071
1072
1073int VspTree::EvalMaxEventContribution(KdLeaf *leaf) const
1074{
1075        int pvs = 0;
1076
1077        // leaf falls out of right pvs
1078        if (-- leaf->mCounter == 0)
1079        {
1080                pvs -= ((int)leaf->mObjects.size() - (int)leaf->mMultipleObjects.size());
1081        }
1082
1083        //-- separately handle objects which are in several kd leaves
1084
1085        ObjectContainer::const_iterator oit, oit_end = leaf->mMultipleObjects.end();
1086
1087        for (oit = leaf->mMultipleObjects.begin(); oit != oit_end; ++ oit)
1088        {
1089                Intersectable *object = *oit;
1090
1091                if (-- object->mCounter == 0)
1092                {
1093                        ++ pvs;
1094                }
1095        }
1096
1097        return pvs;
1098}
1099
1100
1101int VspTree::EvalMinEventContribution(KdLeaf *leaf) const
1102{
1103        if (leaf->Mailed())
1104                return 0;
1105       
1106        leaf->Mail();
1107
1108        // add objects without those which are part of several kd leaves
1109        int pvs = ((int)leaf->mObjects.size() - (int)leaf->mMultipleObjects.size());
1110
1111        // separately handle objects which are part of several kd leaves
1112        ObjectContainer::const_iterator oit, oit_end = leaf->mMultipleObjects.end();
1113
1114        for (oit = leaf->mMultipleObjects.begin(); oit != oit_end; ++ oit)
1115        {
1116                Intersectable *object = *oit;
1117
1118                // object not previously in pvs
1119                if (!object->Mailed())
1120                {
1121                        object->Mail();
1122                        ++ pvs;
1123                }
1124        }       
1125
1126        return pvs;
1127}
1128
1129
1130void VspTree::EvalHeuristics(const SortableEntry &ci,
1131                                                         int &pvsLeft,
1132                                                         int &pvsRight) const
1133{
1134        VssRay *ray = ci.ray;
1135
1136        // eval changes in pvs causes by min event
1137        if (ci.type == SortableEntry::ERayMin)
1138        {
1139                pvsLeft += EvalMinEventContribution(*ray, true);
1140                pvsLeft += EvalMinEventContribution(*ray, false);
1141        }
1142        else // eval changes in pvs causes by max event
1143        {
1144                pvsRight -= EvalMaxEventContribution(*ray, true);
1145                pvsRight -= EvalMaxEventContribution(*ray, false);
1146        }
1147}
1148
1149
1150float VspTree::EvalLocalCostHeuristics(const VspTraversalData &tData,
1151                                                                           const AxisAlignedBox3 &box,
1152                                                                           const int axis,
1153                                                                           float &position)
1154{
1155        // get subset of rays
1156        RayInfoContainer usedRays;
1157
1158        if (mMaxTests < (int)tData.mRays->size())
1159        {
1160                GetRayInfoSets(*tData.mRays, mMaxTests, usedRays);
1161        }
1162        else
1163        {
1164                usedRays = *tData.mRays;
1165        }
1166
1167        const float minBox = box.Min(axis);
1168        const float maxBox = box.Max(axis);
1169
1170        const float sizeBox = maxBox - minBox;
1171
1172        const float minBand = minBox + mMinBand * sizeBox;
1173        const float maxBand = minBox + mMaxBand * sizeBox;
1174
1175        SortSubdivisionCandidates(usedRays, axis, minBand, maxBand);
1176
1177        // prepare the sweep
1178        // note: returns pvs size => no need t give pvs size as function parameter
1179        const int pvsSize = PrepareHeuristics(usedRays);
1180
1181        // go through the lists, count the number of objects left and right
1182        // and evaluate the following cost funcion:
1183        // C = ct_div_ci  + (ql*rl + qr*rr)/queries
1184
1185        int pvsl = 0;
1186        int pvsr = pvsSize;
1187
1188        int pvsBack = pvsl;
1189        int pvsFront = pvsr;
1190
1191        float sum = (float)pvsSize * sizeBox;
1192        float minSum = 1e20f;
1193
1194        // if no good split can be found, take mid split
1195        position = minBox + 0.5f * sizeBox;
1196       
1197        // the relative cost ratio
1198        float ratio = 99999999.0f;
1199        bool splitPlaneFound = false;
1200
1201        Intersectable::NewMail();
1202        KdLeaf::NewMail();
1203        BvhLeaf::NewMail();
1204
1205        //-- traverse through visibility events
1206        vector<SortableEntry>::const_iterator ci, ci_end = mLocalSubdivisionCandidates->end();
1207
1208        for (ci = mLocalSubdivisionCandidates->begin(); ci != ci_end; ++ ci)
1209        {
1210                // compute changes to front and back pvs
1211                EvalHeuristics(*ci, pvsl, pvsr);
1212
1213                // Note: sufficient to compare size of bounding boxes of front and back side?
1214                if (((*ci).value >= minBand) && ((*ci).value <= maxBand))
1215                {
1216                        float currentPos;
1217                       
1218                        // HACK: current positition is BETWEEN visibility events
1219                        if (0 && ((ci + 1) != ci_end))
1220                                currentPos = ((*ci).value + (*(ci + 1)).value) * 0.5f;
1221                        else
1222                                currentPos = (*ci).value;                       
1223
1224                        sum = pvsl * ((*ci).value - minBox) + pvsr * (maxBox - (*ci).value);
1225                       
1226
1227                        if (sum < minSum)
1228                        {
1229                                splitPlaneFound = true;
1230
1231                                minSum = sum;
1232                                position = currentPos;
1233                               
1234                                pvsBack = pvsl;
1235                                pvsFront = pvsr;
1236                        }
1237                }
1238        }
1239       
1240       
1241        //-- compute cost
1242
1243        const int lowerPvsLimit = mViewCellsManager->GetMinPvsSize();
1244        const int upperPvsLimit = mViewCellsManager->GetMaxPvsSize();
1245
1246        const float pOverall = sizeBox;
1247        const float pBack = position - minBox;
1248        const float pFront = maxBox - position;
1249
1250        const float penaltyOld = EvalPvsPenalty(pvsSize, lowerPvsLimit, upperPvsLimit);
1251    const float penaltyFront = EvalPvsPenalty(pvsFront, lowerPvsLimit, upperPvsLimit);
1252        const float penaltyBack = EvalPvsPenalty(pvsBack, lowerPvsLimit, upperPvsLimit);
1253       
1254        const float oldRenderCost = penaltyOld * pOverall + Limits::Small;
1255        const float newRenderCost = penaltyFront * pFront + penaltyBack * pBack;
1256
1257        if (splitPlaneFound)
1258        {
1259                ratio = newRenderCost / oldRenderCost;
1260        }
1261       
1262        const float volRatio = tData.mBoundingBox.GetVolume() / (sizeBox * mBoundingBox.GetVolume());
1263
1264        Debug << "\n§§§§ eval local cost §§§§" << endl
1265                  << "back pvs: " << penaltyBack << " front pvs: " << penaltyFront << " total pvs: " << penaltyOld << endl
1266                  << "back p: " << pBack * volRatio << " front p " << pFront * volRatio << " p: " << pOverall * volRatio << endl
1267                  << "old rc: " << oldRenderCost * volRatio << " new rc: " << newRenderCost * volRatio << endl
1268                  << "render cost decrease: " << oldRenderCost * volRatio - newRenderCost * volRatio << endl;
1269
1270        return ratio;
1271}
1272
1273
1274float VspTree::SelectSplitPlane(const VspTraversalData &tData,
1275                                                                AxisAlignedPlane &plane,
1276                                                                float &pFront,
1277                                                                float &pBack)
1278{
1279        float nPosition[3];
1280        float nCostRatio[3];
1281        float nProbFront[3];
1282        float nProbBack[3];
1283
1284        // create bounding box of node geometry
1285        AxisAlignedBox3 box = tData.mBoundingBox;
1286               
1287        int sAxis = 0;
1288        int bestAxis = -1;
1289
1290        // if we use some kind of specialised fixed axis
1291    const bool useSpecialAxis =
1292                mOnlyDrivingAxis || mCirculatingAxis;
1293       
1294        if (mCirculatingAxis)
1295        {
1296                int parentAxis = 0;
1297                VspNode *parent = tData.mNode->GetParent();
1298
1299                if (parent)
1300                        parentAxis = dynamic_cast<VspInterior *>(parent)->GetAxis();
1301
1302                sAxis = (parentAxis + 1) % 3;
1303        }
1304        else if (mOnlyDrivingAxis)
1305        {
1306                sAxis = box.Size().DrivingAxis();
1307        }
1308       
1309        for (int axis = 0; axis < 3; ++ axis)
1310        {
1311                if (!useSpecialAxis || (axis == sAxis))
1312                {
1313                        if (mUseCostHeuristics)
1314                        {
1315                                //-- place split plane using heuristics
1316                                nCostRatio[axis] =
1317                                        EvalLocalCostHeuristics(tData,
1318                                                                                        box,
1319                                                                                        axis,
1320                                                                                        nPosition[axis]);                       
1321                        }
1322                        else
1323                        {
1324                                //-- split plane position is spatial median
1325                                nPosition[axis] = (box.Min()[axis] + box.Max()[axis]) * 0.5f;
1326                                nCostRatio[axis] = EvalLocalSplitCost(tData,
1327                                                                                                          box,
1328                                                                                                          axis,
1329                                                                                                          nPosition[axis],
1330                                                                                                          nProbFront[axis],
1331                                                                                                          nProbBack[axis]);
1332                        }
1333                                               
1334                        if (bestAxis == -1)
1335                        {
1336                                bestAxis = axis;
1337                        }
1338                        else if (nCostRatio[axis] < nCostRatio[bestAxis])
1339                        {
1340                                bestAxis = axis;
1341                        }
1342                }
1343        }
1344
1345        //-- assign values of best split
1346        plane.mAxis = bestAxis;
1347        plane.mPosition = nPosition[bestAxis]; // split plane position
1348
1349        pFront = nProbFront[bestAxis];
1350        pBack = nProbBack[bestAxis];
1351
1352        return nCostRatio[bestAxis];
1353}
1354
1355
1356float VspTree::EvalRenderCostDecrease(const AxisAlignedPlane &candidatePlane,
1357                                                                          const VspTraversalData &data,
1358                                                                          float &normalizedOldRenderCost) const
1359{
1360        float pvsFront = 0;
1361        float pvsBack = 0;
1362        float totalPvs = 0;
1363
1364        const float viewSpaceVol = mBoundingBox.GetVolume();
1365
1366        // create unique ids for pvs heuristics
1367        Intersectable::NewMail(3);
1368        KdLeaf::NewMail(3);
1369        BvhLeaf::NewMail(3);
1370
1371        RayInfoContainer::const_iterator rit, rit_end = data.mRays->end();
1372
1373        for (rit = data.mRays->begin(); rit != rit_end; ++ rit)
1374        {
1375                RayInfo rayInf = *rit;
1376
1377                float t;
1378                VssRay *ray = rayInf.mRay;
1379
1380                // classify ray
1381                const int cf =
1382                        rayInf.ComputeRayIntersection(candidatePlane.mAxis,
1383                                                                                  candidatePlane.mPosition, t);
1384
1385                // evaluate contribution of ray endpoint to front and back pvs
1386                // with respect to the classification
1387                UpdateContributionsToPvs(*ray, true, cf, pvsFront, pvsBack, totalPvs); 
1388                UpdateContributionsToPvs(*ray, false, cf, pvsFront, pvsBack, totalPvs);
1389        }
1390
1391        AxisAlignedBox3 frontBox;
1392        AxisAlignedBox3 backBox;
1393
1394        data.mBoundingBox.Split(candidatePlane.mAxis, candidatePlane.mPosition, frontBox, backBox);
1395
1396        // probability that view point lies in back / front node
1397        float pOverall = data.mProbability;
1398        float pFront = pFront = frontBox.GetVolume();
1399        float pBack = pOverall - pFront;
1400
1401
1402        //-- pvs rendering heuristics
1403
1404        const int lowerPvsLimit = mViewCellsManager->GetMinPvsSize();
1405        const int upperPvsLimit = mViewCellsManager->GetMaxPvsSize();
1406
1407        // only render cost heuristics or combined with standard deviation
1408        const float penaltyOld = EvalPvsPenalty((int)totalPvs, lowerPvsLimit, upperPvsLimit);
1409    const float penaltyFront = EvalPvsPenalty((int)pvsFront, lowerPvsLimit, upperPvsLimit);
1410        const float penaltyBack = EvalPvsPenalty((int)pvsBack, lowerPvsLimit, upperPvsLimit);
1411                       
1412        const float oldRenderCost = pOverall * penaltyOld;
1413        const float newRenderCost = penaltyFront * pFront + penaltyBack * pBack;
1414
1415        normalizedOldRenderCost = oldRenderCost / viewSpaceVol;
1416
1417        const float renderCostDecrease = (oldRenderCost - newRenderCost) / viewSpaceVol;
1418       
1419        Debug << "\neval vsp render cost decrease" << endl
1420                  << "back pvs: " << pvsBack << " front pvs " << pvsFront << " total pvs: " << totalPvs << endl
1421                  << "back p: " << pBack / viewSpaceVol << " front p " << pFront / viewSpaceVol << " p: " << pOverall / viewSpaceVol << endl
1422                  << "old rc: " << normalizedOldRenderCost << " new rc: " << newRenderCost / viewSpaceVol << endl
1423                  << "render cost decrease: " << renderCostDecrease << endl;
1424
1425        return renderCostDecrease;
1426}
1427
1428
1429
1430float VspTree::EvalLocalSplitCost(const VspTraversalData &data,
1431                                                                  const AxisAlignedBox3 &box,
1432                                                                  const int axis,
1433                                                                  const float &position,
1434                                                                  float &pFront,
1435                                                                  float &pBack) const
1436{
1437        float pvsTotal = 0;
1438        float pvsFront = 0;
1439        float pvsBack = 0;
1440       
1441        // create unique ids for pvs heuristics
1442        Intersectable::NewMail(3);
1443        BvhLeaf::NewMail(3);
1444        KdLeaf::NewMail(3);
1445
1446        const int pvsSize = data.mPvs;
1447
1448        RayInfoContainer::const_iterator rit, rit_end = data.mRays->end();
1449
1450        // this is the main ray classification loop!
1451        for(rit = data.mRays->begin(); rit != rit_end; ++ rit)
1452        {
1453                VssRay *ray = (*rit).mRay;
1454
1455                // determine the side of this ray with respect to the plane
1456                float t;
1457                const int side = (*rit).ComputeRayIntersection(axis, position, t);
1458       
1459                UpdateContributionsToPvs(*ray, true, side, pvsFront, pvsBack, pvsTotal);
1460                UpdateContributionsToPvs(*ray, false, side, pvsFront, pvsBack, pvsTotal);
1461        }
1462
1463
1464        //-- evaluate cost heuristics
1465
1466        float pOverall = data.mProbability;
1467
1468        // we use spatial mid split => simplified computation
1469        pBack = pFront = pOverall * 0.5f;
1470       
1471        const float newCost = pvsBack * pBack + pvsFront * pFront;
1472        const float oldCost = (float)pvsSize * pOverall + Limits::Small;
1473       
1474#ifdef _DEBUG
1475        Debug << "axis: " << axis << " " << pvsSize << " " << pvsBack << " " << pvsFront << endl;
1476        Debug << "p: " << pFront << " " << pBack << " " << pOverall << endl;
1477#endif
1478
1479        return  (mCtDivCi + newCost) / oldCost;
1480}
1481
1482
1483void VspTree::UpdateContributionsToPvs(Intersectable *obj,
1484                                                                           const int cf,
1485                                                                           float &frontPvs,
1486                                                                           float &backPvs,
1487                                                                           float &totalPvs) const
1488{
1489        if (!obj) return;
1490
1491        //const float renderCost = mViewCellsManager->SimpleRay &raynderCost(obj);
1492        const int renderCost = 1;
1493
1494        // object in no pvs => new
1495        if (!obj->Mailed() && !obj->Mailed(1) && !obj->Mailed(2))
1496        {
1497                totalPvs += renderCost;
1498        }
1499
1500        // QUESTION matt: is it safe to assume that
1501        // the object belongs to no pvs in this case?
1502        //if (cf == Ray::COINCIDENT) return;
1503
1504        if (cf >= 0) // front pvs
1505        {
1506                if (!obj->Mailed() && !obj->Mailed(2))
1507                {
1508                        frontPvs += renderCost;
1509               
1510                        // already in back pvs => in both pvss
1511                        if (obj->Mailed(1))
1512                                obj->Mail(2);
1513                        else
1514                                obj->Mail();
1515                }
1516        }
1517
1518        if (cf <= 0) // back pvs
1519        {
1520                if (!obj->Mailed(1) && !obj->Mailed(2))
1521                {
1522                        backPvs += renderCost;
1523               
1524                        // already in front pvs => in both pvss
1525                        if (obj->Mailed())
1526                                obj->Mail(2);
1527                        else
1528                                obj->Mail(1);
1529                }
1530        }
1531}
1532
1533
1534void VspTree::UpdateContributionsToPvs(BvhLeaf *leaf,
1535                                                                           const int cf,
1536                                                                           float &frontPvs,
1537                                                                           float &backPvs,
1538                                                                           float &totalPvs) const
1539{
1540        if (!leaf) return;
1541
1542        const int renderCost = (int)leaf->mObjects.size();
1543
1544        // leaf in no pvs => new
1545        if (!leaf->Mailed() && !leaf->Mailed(1) && !leaf->Mailed(2))
1546        {
1547                totalPvs += renderCost;
1548        }
1549
1550        // QUESTION matt: is it safe to assume that
1551        // the leaf belongs to no pvs in this case?
1552        //if (cf == Ray::COINCIDENT) return;
1553
1554        if (cf >= 0) // front pvs
1555        {
1556                if (!leaf->Mailed() && !leaf->Mailed(2))
1557                {
1558                        frontPvs += renderCost;
1559       
1560                        // already in back pvs => in both pvss
1561                        if (leaf->Mailed(1))
1562                                leaf->Mail(2);
1563                        else
1564                                leaf->Mail();
1565                }
1566        }
1567
1568        if (cf <= 0) // back pvs
1569        {
1570                if (!leaf->Mailed(1) && !leaf->Mailed(2))
1571                {
1572                        backPvs += renderCost;
1573               
1574                        // already in front pvs => in both pvss
1575                        if (leaf->Mailed())
1576                        {
1577                                leaf->Mail(2);
1578                        }
1579                        else
1580                                leaf->Mail(1);
1581                }
1582        }
1583}
1584
1585
1586
1587void VspTree::UpdateContributionsToPvs(KdLeaf *leaf,
1588                                                                           const int cf,
1589                                                                           float &frontPvs,
1590                                                                           float &backPvs,
1591                                                                           float &totalPvs) const
1592{
1593        if (!leaf) return;
1594
1595        // the objects which are referenced in this and only this leaf
1596        const int contri = (int)(leaf->mObjects.size() - leaf->mMultipleObjects.size());
1597       
1598        // newly found leaf
1599        if (!leaf->Mailed() && !leaf->Mailed(1) && !leaf->Mailed(2))
1600        {
1601                totalPvs += contri;
1602        }
1603
1604        // recursivly update contributions of yet unclassified objects
1605        ObjectContainer::const_iterator oit, oit_end = leaf->mMultipleObjects.end();
1606
1607        for (oit = leaf->mMultipleObjects.begin(); oit != oit_end; ++ oit)
1608        {       
1609                UpdateContributionsToPvs(*oit, cf, frontPvs, backPvs, totalPvs);
1610    }   
1611       
1612        // QUESTION matt: is it safe to assume that
1613        // the object belongs to no pvs in this case?
1614        //if (cf == Ray::COINCIDENT) return;
1615
1616        if (cf >= 0) // front pvs
1617        {
1618                if (!leaf->Mailed() && !leaf->Mailed(2))
1619                {
1620                        frontPvs += contri;
1621               
1622                        // already in back pvs => in both pvss
1623                        if (leaf->Mailed(1))
1624                                leaf->Mail(2);
1625                        else
1626                                leaf->Mail();
1627                }
1628        }
1629
1630        if (cf <= 0) // back pvs
1631        {
1632                if (!leaf->Mailed(1) && !leaf->Mailed(2))
1633                {
1634                        backPvs += contri;
1635               
1636                        // already in front pvs => in both pvss
1637                        if (leaf->Mailed())
1638                                leaf->Mail(2);
1639                        else
1640                                leaf->Mail(1);
1641                }
1642        }
1643}
1644
1645
1646void VspTree::CollectLeaves(vector<VspLeaf *> &leaves,
1647                                                        const bool onlyUnmailed,
1648                                                        const int maxPvsSize) const
1649{
1650        stack<VspNode *> nodeStack;
1651        nodeStack.push(mRoot);
1652
1653        while (!nodeStack.empty())
1654        {
1655                VspNode *node = nodeStack.top();
1656                nodeStack.pop();
1657               
1658                if (node->IsLeaf())
1659                {
1660                        // test if this leaf is in valid view space
1661                        VspLeaf *leaf = dynamic_cast<VspLeaf *>(node);
1662                        if (leaf->TreeValid() &&
1663                                (!onlyUnmailed || !leaf->Mailed()) &&
1664                                ((maxPvsSize < 0) || (leaf->GetViewCell()->GetPvs().CountObjectsInPvs() <= maxPvsSize)))
1665                        {
1666                                leaves.push_back(leaf);
1667                        }
1668                }
1669                else
1670                {
1671                        VspInterior *interior = dynamic_cast<VspInterior *>(node);
1672
1673                        nodeStack.push(interior->GetBack());
1674                        nodeStack.push(interior->GetFront());
1675                }
1676        }
1677}
1678
1679
1680AxisAlignedBox3 VspTree::GetBoundingBox() const
1681{
1682        return mBoundingBox;
1683}
1684
1685
1686VspNode *VspTree::GetRoot() const
1687{
1688        return mRoot;
1689}
1690
1691
1692void VspTree::EvaluateLeafStats(const VspTraversalData &data)
1693{
1694        // the node became a leaf -> evaluate stats for leafs
1695        VspLeaf *leaf = dynamic_cast<VspLeaf *>(data.mNode);
1696
1697
1698        if (data.mPvs > mVspStats.maxPvs)
1699        {
1700                mVspStats.maxPvs = data.mPvs;
1701        }
1702
1703        mVspStats.pvs += data.mPvs;
1704
1705        if (data.mDepth < mVspStats.minDepth)
1706        {
1707                mVspStats.minDepth = data.mDepth;
1708        }
1709       
1710        if (data.mDepth >= mTermMaxDepth)
1711        {
1712        ++ mVspStats.maxDepthNodes;
1713                //Debug << "new max depth: " << mVspStats.maxDepthNodes << endl;
1714        }
1715
1716        // accumulate rays to compute rays /  leaf
1717        mVspStats.accumRays += (int)data.mRays->size();
1718
1719        if (data.mPvs < mTermMinPvs)
1720                ++ mVspStats.minPvsNodes;
1721
1722        if ((int)data.mRays->size() < mTermMinRays)
1723                ++ mVspStats.minRaysNodes;
1724
1725        if (data.GetAvgRayContribution() > mTermMaxRayContribution)
1726                ++ mVspStats.maxRayContribNodes;
1727
1728        if (data.mProbability <= mTermMinProbability)
1729                ++ mVspStats.minProbabilityNodes;
1730       
1731        // accumulate depth to compute average depth
1732        mVspStats.accumDepth += data.mDepth;
1733
1734        ++ mCreatedViewCells;
1735
1736#ifdef _DEBUG
1737        Debug << "BSP stats: "
1738                  << "Depth: " << data.mDepth << " (max: " << mTermMaxDepth << "), "
1739                  << "PVS: " << data.mPvs << " (min: " << mTermMinPvs << "), "
1740                  << "#rays: " << (int)data.mRays->size() << " (max: " << mTermMinRays << "), "
1741                  << "#pvs: " << leaf->GetViewCell()->GetPvs().CountObjectsInPvs() << "), "
1742                  << "#avg ray contrib (pvs): " << (float)data.mPvs / (float)data.mRays->size() << endl;
1743#endif
1744}
1745
1746
1747void VspTree::CollectViewCells(ViewCellContainer &viewCells, bool onlyValid) const
1748{
1749        ViewCell::NewMail();
1750        CollectViewCells(mRoot, onlyValid, viewCells, true);
1751}
1752
1753
1754void VspTree::CollapseViewCells()
1755{
1756// TODO matt
1757#if HAS_TO_BE_REDONE
1758        stack<VspNode *> nodeStack;
1759
1760        if (!mRoot)
1761                return;
1762
1763        nodeStack.push(mRoot);
1764       
1765        while (!nodeStack.empty())
1766        {
1767                VspNode *node = nodeStack.top();
1768                nodeStack.pop();
1769               
1770                if (node->IsLeaf())
1771        {
1772                        BspViewCell *viewCell = dynamic_cast<VspLeaf *>(node)->GetViewCell();
1773
1774                        if (!viewCell->GetValid())
1775                        {
1776                                BspViewCell *viewCell = dynamic_cast<VspLeaf *>(node)->GetViewCell();
1777       
1778                                ViewCellContainer leaves;
1779                                mViewCellsTree->CollectLeaves(viewCell, leaves);
1780
1781                                ViewCellContainer::const_iterator it, it_end = leaves.end();
1782
1783                                for (it = leaves.begin(); it != it_end; ++ it)
1784                                {
1785                                        VspLeaf *l = dynamic_cast<BspViewCell *>(*it)->mLeaf;
1786                                        l->SetViewCell(GetOrCreateOutOfBoundsCell());
1787                                        ++ mVspStats.invalidLeaves;
1788                                }
1789
1790                                // add to unbounded view cell
1791                                GetOrCreateOutOfBoundsCell()->GetPvs().AddPvs(viewCell->GetPvs());
1792                                DEL_PTR(viewCell);
1793                        }
1794                }
1795                else
1796                {
1797                        VspInterior *interior = dynamic_cast<VspInterior *>(node);
1798               
1799                        nodeStack.push(interior->GetFront());
1800                        nodeStack.push(interior->GetBack());
1801                }
1802        }
1803
1804        Debug << "invalid leaves: " << mVspStats.invalidLeaves << endl;
1805#endif
1806}
1807
1808
1809void VspTree::CollectRays(VssRayContainer &rays)
1810{
1811        vector<VspLeaf *> leaves;
1812        CollectLeaves(leaves);
1813
1814        vector<VspLeaf *>::const_iterator lit, lit_end = leaves.end();
1815
1816        for (lit = leaves.begin(); lit != lit_end; ++ lit)
1817        {
1818                VspLeaf *leaf = *lit;
1819                VssRayContainer::const_iterator rit, rit_end = leaf->mVssRays.end();
1820
1821                for (rit = leaf->mVssRays.begin(); rit != rit_end; ++ rit)
1822                        rays.push_back(*rit);
1823        }
1824}
1825
1826
1827void VspTree::SetViewCellsManager(ViewCellsManager *vcm)
1828{
1829        mViewCellsManager = vcm;
1830}
1831
1832
1833void VspTree::ValidateTree()
1834{
1835        mVspStats.invalidLeaves = 0;
1836
1837        stack<VspNode *> nodeStack;
1838
1839        if (!mRoot)
1840                return;
1841
1842        nodeStack.push(mRoot);
1843
1844        while (!nodeStack.empty())
1845        {
1846                VspNode *node = nodeStack.top();
1847                nodeStack.pop();
1848               
1849                if (node->IsLeaf())
1850                {
1851                        VspLeaf *leaf = dynamic_cast<VspLeaf *>(node);
1852
1853                        if (!leaf->GetViewCell()->GetValid())
1854                                ++ mVspStats.invalidLeaves;
1855
1856                        // validity flags don't match => repair
1857                        if (leaf->GetViewCell()->GetValid() != leaf->TreeValid())
1858                        {
1859                                leaf->SetTreeValid(leaf->GetViewCell()->GetValid());
1860                                PropagateUpValidity(leaf);
1861                        }
1862                }
1863                else
1864                {
1865                        VspInterior *interior = dynamic_cast<VspInterior *>(node);
1866               
1867                        nodeStack.push(interior->GetFront());
1868                        nodeStack.push(interior->GetBack());
1869                }
1870        }
1871
1872        Debug << "invalid leaves: " << mVspStats.invalidLeaves << endl;
1873}
1874
1875
1876
1877void VspTree::CollectViewCells(VspNode *root,
1878                                                                  bool onlyValid,
1879                                                                  ViewCellContainer &viewCells,
1880                                                                  bool onlyUnmailed) const
1881{
1882        stack<VspNode *> nodeStack;
1883
1884        if (!root)
1885                return;
1886
1887        nodeStack.push(root);
1888       
1889        while (!nodeStack.empty())
1890        {
1891                VspNode *node = nodeStack.top();
1892                nodeStack.pop();
1893               
1894                if (node->IsLeaf())
1895                {
1896                        if (!onlyValid || node->TreeValid())
1897                        {
1898                                ViewCellLeaf *leafVc = dynamic_cast<VspLeaf *>(node)->GetViewCell();
1899
1900                                ViewCell *viewCell = mViewCellsTree->GetActiveViewCell(leafVc);
1901                                               
1902                                if (!onlyUnmailed || !viewCell->Mailed())
1903                                {
1904                                        viewCell->Mail();
1905                                        viewCells.push_back(viewCell);
1906                                }
1907                        }
1908                }
1909                else
1910                {
1911                        VspInterior *interior = dynamic_cast<VspInterior *>(node);
1912               
1913                        nodeStack.push(interior->GetFront());
1914                        nodeStack.push(interior->GetBack());
1915                }
1916        }
1917}
1918
1919
1920int VspTree::FindNeighbors(VspLeaf *n,
1921                                                   vector<VspLeaf *> &neighbors,
1922                                                   const bool onlyUnmailed) const
1923{
1924        stack<VspNode *> nodeStack;
1925        nodeStack.push(mRoot);
1926
1927        const AxisAlignedBox3 box = GetBoundingBox(n);
1928
1929        while (!nodeStack.empty())
1930        {
1931                VspNode *node = nodeStack.top();
1932                nodeStack.pop();
1933
1934                if (node->IsLeaf())
1935                {
1936                        VspLeaf *leaf = dynamic_cast<VspLeaf *>(node);
1937
1938                        if (leaf != n && (!onlyUnmailed || !leaf->Mailed()))
1939                                neighbors.push_back(leaf);
1940                }
1941                else
1942                {
1943                        VspInterior *interior = dynamic_cast<VspInterior *>(node);
1944                       
1945                        if (interior->GetPosition() > box.Max(interior->GetAxis()))
1946                                nodeStack.push(interior->GetBack());
1947                        else
1948                        {
1949                                if (interior->GetPosition() < box.Min(interior->GetAxis()))
1950                                        nodeStack.push(interior->GetFront());
1951                                else
1952                                {
1953                                        // random decision
1954                                        nodeStack.push(interior->GetBack());
1955                                        nodeStack.push(interior->GetFront());
1956                                }
1957                        }
1958                }
1959        }
1960
1961        return (int)neighbors.size();
1962}
1963
1964
1965// Find random neighbor which was not mailed
1966VspLeaf *VspTree::GetRandomLeaf(const Plane3 &plane)
1967{
1968        stack<VspNode *> nodeStack;
1969        nodeStack.push(mRoot);
1970 
1971        int mask = rand();
1972 
1973        while (!nodeStack.empty())
1974        {
1975                VspNode *node = nodeStack.top();
1976               
1977                nodeStack.pop();
1978               
1979                if (node->IsLeaf())
1980                {
1981                        return dynamic_cast<VspLeaf *>(node);
1982                }
1983                else
1984                {
1985                        VspInterior *interior = dynamic_cast<VspInterior *>(node);
1986                        VspNode *next;
1987                       
1988                        if (GetBoundingBox(interior->GetBack()).Side(plane) < 0)
1989                        {
1990                                next = interior->GetFront();
1991                        }
1992            else
1993                        {
1994                                if (GetBoundingBox(interior->GetFront()).Side(plane) < 0)
1995                                {
1996                                        next = interior->GetBack();
1997                                }
1998                                else
1999                                {
2000                                        // random decision
2001                                        if (mask & 1)
2002                                                next = interior->GetBack();
2003                                        else
2004                                                next = interior->GetFront();
2005                                        mask = mask >> 1;
2006                                }
2007                        }
2008                       
2009                        nodeStack.push(next);
2010                }
2011        }
2012 
2013        return NULL;
2014}
2015
2016
2017VspLeaf *VspTree::GetRandomLeaf(const bool onlyUnmailed)
2018{
2019        stack<VspNode *> nodeStack;
2020
2021        nodeStack.push(mRoot);
2022
2023        int mask = rand();
2024
2025        while (!nodeStack.empty())
2026        {
2027                VspNode *node = nodeStack.top();
2028                nodeStack.pop();
2029
2030                if (node->IsLeaf())
2031                {
2032                        if ( (!onlyUnmailed || !node->Mailed()) )
2033                                return dynamic_cast<VspLeaf *>(node);
2034                }
2035                else
2036                {
2037                        VspInterior *interior = dynamic_cast<VspInterior *>(node);
2038
2039                        // random decision
2040                        if (mask & 1)
2041                                nodeStack.push(interior->GetBack());
2042                        else
2043                                nodeStack.push(interior->GetFront());
2044
2045                        mask = mask >> 1;
2046                }
2047        }
2048
2049        return NULL;
2050}
2051
2052
2053int VspTree::EvalPvsSize(const RayInfoContainer &rays) const
2054{
2055        int pvsSize = 0;
2056       
2057        Intersectable::NewMail();
2058        KdNode::NewMail();
2059        BvhLeaf::NewMail();
2060
2061        RayInfoContainer::const_iterator rit, rit_end = rays.end();
2062
2063        for (rit = rays.begin(); rit != rays.end(); ++ rit)
2064        {
2065                VssRay *ray = (*rit).mRay;
2066
2067                pvsSize += EvalContributionToPvs(*ray, true);
2068                pvsSize += EvalContributionToPvs(*ray, false);
2069        }
2070       
2071        return pvsSize;
2072}
2073
2074
2075float VspTree::GetEpsilon() const
2076{
2077        return mEpsilon;
2078}
2079
2080
2081int VspTree::CastLineSegment(const Vector3 &origin,
2082                                                         const Vector3 &termination,
2083                             ViewCellContainer &viewcells,
2084                                                         const bool useMailboxing)
2085{
2086        int hits = 0;
2087
2088        float mint = 0.0f, maxt = 1.0f;
2089        const Vector3 dir = termination - origin;
2090
2091        stack<LineTraversalData> tStack;
2092
2093        Vector3 entp = origin;
2094        Vector3 extp = termination;
2095
2096        VspNode *node = mRoot;
2097        VspNode *farChild;
2098
2099        float position;
2100        int axis;
2101
2102        while (1)
2103        {
2104                if (!node->IsLeaf())
2105                {
2106                        VspInterior *in = dynamic_cast<VspInterior *>(node);
2107                        position = in->GetPosition();
2108                        axis = in->GetAxis();
2109
2110                        if (entp[axis] <= position)
2111                        {
2112                                if (extp[axis] <= position)
2113                                {
2114                                        node = in->GetBack();
2115                                        // cases N1,N2,N3,P5,Z2,Z3
2116                                        continue;
2117                                } else
2118                                {
2119                                        // case N4
2120                                        node = in->GetBack();
2121                                        farChild = in->GetFront();
2122                                }
2123                        }
2124                        else
2125                        {
2126                                if (position <= extp[axis])
2127                                {
2128                                        node = in->GetFront();
2129                                        // cases P1,P2,P3,N5,Z1
2130                                        continue;
2131                                }
2132                                else
2133                                {
2134                                        node = in->GetFront();
2135                                        farChild = in->GetBack();
2136                                        // case P4
2137                                }
2138                        }
2139
2140                        // $$ modification 3.5.2004 - hints from Kamil Ghais
2141                        // case N4 or P4
2142                        const float tdist = (position - origin[axis]) / dir[axis];
2143                        tStack.push(LineTraversalData(farChild, extp, maxt)); //TODO
2144
2145                        extp = origin + dir * tdist;
2146                        maxt = tdist;
2147                }
2148                else
2149                {
2150                        // compute intersection with all objects in this leaf
2151                        VspLeaf *leaf = dynamic_cast<VspLeaf *>(node);
2152                        ViewCell *vc = leaf->GetViewCell();
2153
2154                        // don't have to mail if each view cell belongs to exactly one leaf
2155                        if (!useMailboxing || !vc->Mailed())
2156                        {
2157                                if (useMailboxing)
2158                                        vc->Mail();
2159
2160                                viewcells.push_back(vc);
2161                                ++ hits;
2162                        }
2163#if 0
2164                        leaf->mRays.push_back(RayInfo(new VssRay(origin, termination, NULL, NULL, 0)));
2165#endif
2166                        // get the next node from the stack
2167                        if (tStack.empty())
2168                                break;
2169
2170                        entp = extp;
2171                        mint = maxt;
2172                       
2173                        LineTraversalData &s  = tStack.top();
2174                        node = s.mNode;
2175                        extp = s.mExitPoint;
2176                        maxt = s.mMaxT;
2177                        tStack.pop();
2178                }
2179        }
2180
2181        return hits;
2182}
2183
2184
2185int VspTree::CastRay(Ray &ray)
2186{
2187        int hits = 0;
2188
2189        stack<LineTraversalData> tStack;
2190        const Vector3 dir = ray.GetDir();
2191
2192        float maxt, mint;
2193
2194        if (!mBoundingBox.GetRaySegment(ray, mint, maxt))
2195                return 0;
2196
2197        Intersectable::NewMail();
2198        ViewCell::NewMail();
2199
2200        Vector3 entp = ray.Extrap(mint);
2201        Vector3 extp = ray.Extrap(maxt);
2202
2203        const Vector3 origin = entp;
2204
2205        VspNode *node = mRoot;
2206        VspNode *farChild = NULL;
2207
2208        float position;
2209        int axis;
2210
2211        while (1)
2212        {
2213                if (!node->IsLeaf())
2214                {
2215                        VspInterior *in = dynamic_cast<VspInterior *>(node);
2216                        position = in->GetPosition();
2217                        axis = in->GetAxis();
2218
2219                        if (entp[axis] <= position)
2220                        {
2221                                if (extp[axis] <= position)
2222                                {
2223                                        node = in->GetBack();
2224                                        // cases N1,N2,N3,P5,Z2,Z3
2225                                        continue;
2226                                }
2227                                else
2228                                {
2229                                        // case N4
2230                                        node = in->GetBack();
2231                                        farChild = in->GetFront();
2232                                }
2233                        }
2234                        else
2235                        {
2236                                if (position <= extp[axis])
2237                                {
2238                                        node = in->GetFront();
2239                                        // cases P1,P2,P3,N5,Z1
2240                                        continue;
2241                                }
2242                                else
2243                                {
2244                                        node = in->GetFront();
2245                                        farChild = in->GetBack();
2246                                        // case P4
2247                                }
2248                        }
2249
2250                        // $$ modification 3.5.2004 - hints from Kamil Ghais
2251                        // case N4 or P4
2252                        const float tdist = (position - origin[axis]) / dir[axis];
2253                        tStack.push(LineTraversalData(farChild, extp, maxt)); //TODO
2254                        extp = origin + dir * tdist;
2255                        maxt = tdist;
2256                }
2257                else
2258                {
2259                        // compute intersection with all objects in this leaf
2260                        VspLeaf *leaf = dynamic_cast<VspLeaf *>(node);
2261                        ViewCell *vc = leaf->GetViewCell();
2262
2263                        if (!vc->Mailed())
2264                        {
2265                                vc->Mail();
2266                                // todo: add view cells to ray
2267                                ++ hits;
2268                        }
2269#if 0
2270                        leaf->mRays.push_back(RayInfo(new VssRay(origin, termination, NULL, NULL, 0)));
2271#endif
2272                        // get the next node from the stack
2273                        if (tStack.empty())
2274                                break;
2275
2276                        entp = extp;
2277                        mint = maxt;
2278                       
2279                        LineTraversalData &s  = tStack.top();
2280                        node = s.mNode;
2281                        extp = s.mExitPoint;
2282                        maxt = s.mMaxT;
2283                        tStack.pop();
2284                }
2285        }
2286
2287        return hits;
2288}
2289
2290
2291ViewCell *VspTree::GetViewCell(const Vector3 &point, const bool active)
2292{
2293        if (mRoot == NULL)
2294                return NULL;
2295
2296        stack<VspNode *> nodeStack;
2297        nodeStack.push(mRoot);
2298 
2299        ViewCellLeaf *viewcell = NULL;
2300 
2301        while (!nodeStack.empty()) 
2302        {
2303                VspNode *node = nodeStack.top();
2304                nodeStack.pop();
2305       
2306                if (node->IsLeaf())
2307                {
2308                        viewcell = dynamic_cast<VspLeaf *>(node)->GetViewCell();
2309                        break;
2310                }
2311                else   
2312                {       
2313                        VspInterior *interior = dynamic_cast<VspInterior *>(node);
2314     
2315                        // random decision
2316                        if (interior->GetPosition() - point[interior->GetAxis()] < 0)
2317                                nodeStack.push(interior->GetBack());
2318                        else
2319                                nodeStack.push(interior->GetFront());
2320                }
2321        }
2322 
2323        if (active)
2324                return mViewCellsTree->GetActiveViewCell(viewcell);
2325        else
2326                return viewcell;
2327}
2328
2329
2330bool VspTree::ViewPointValid(const Vector3 &viewPoint) const
2331{
2332        VspNode *node = mRoot;
2333
2334        while (1)
2335        {
2336                // early exit
2337                if (node->TreeValid())
2338                        return true;
2339
2340                if (node->IsLeaf())
2341                        return false;
2342                       
2343                VspInterior *in = dynamic_cast<VspInterior *>(node);
2344                                       
2345                if (in->GetPosition() - viewPoint[in->GetAxis()] <= 0)
2346                {
2347                        node = in->GetBack();
2348                }
2349                else
2350                {
2351                        node = in->GetFront();
2352                }
2353        }
2354
2355        // should never come here
2356        return false;
2357}
2358
2359
2360void VspTree::PropagateUpValidity(VspNode *node)
2361{
2362        const bool isValid = node->TreeValid();
2363
2364        // propagative up invalid flag until only invalid nodes exist over this node
2365        if (!isValid)
2366        {
2367                while (!node->IsRoot() && node->GetParent()->TreeValid())
2368                {
2369                        node = node->GetParent();
2370                        node->SetTreeValid(false);
2371                }
2372        }
2373        else
2374        {
2375                // propagative up valid flag until one of the subtrees is invalid
2376                while (!node->IsRoot() && !node->TreeValid())
2377                {
2378            node = node->GetParent();
2379                        VspInterior *interior = dynamic_cast<VspInterior *>(node);
2380                       
2381                        // the parent is valid iff both leaves are valid
2382                        node->SetTreeValid(interior->GetBack()->TreeValid() &&
2383                                                           interior->GetFront()->TreeValid());
2384                }
2385        }
2386}
2387
2388
2389bool VspTree::Export(OUT_STREAM &stream)
2390{
2391        ExportNode(mRoot, stream);
2392
2393        return true;
2394}
2395
2396
2397void VspTree::ExportNode(VspNode *node, OUT_STREAM &stream)
2398{
2399        if (node->IsLeaf())
2400        {
2401                VspLeaf *leaf = dynamic_cast<VspLeaf *>(node);
2402                ViewCell *viewCell = mViewCellsTree->GetActiveViewCell(leaf->GetViewCell());
2403
2404                int id = -1;
2405                if (viewCell != mOutOfBoundsCell)
2406                        id = viewCell->GetId();
2407
2408                stream << "<Leaf viewCellId=\"" << id << "\" />" << endl;
2409        }
2410        else
2411        {       
2412                VspInterior *interior = dynamic_cast<VspInterior *>(node);
2413       
2414                AxisAlignedPlane plane = interior->GetPlane();
2415                stream << "<Interior plane=\"" << plane.mPosition << " "
2416                           << plane.mAxis << "\">" << endl;
2417
2418                ExportNode(interior->GetBack(), stream);
2419                ExportNode(interior->GetFront(), stream);
2420
2421                stream << "</Interior>" << endl;
2422        }
2423}
2424
2425
2426int VspTree::SplitRays(const AxisAlignedPlane &plane,
2427                                           RayInfoContainer &rays,
2428                                           RayInfoContainer &frontRays,
2429                                           RayInfoContainer &backRays) const
2430{
2431        int splits = 0;
2432
2433        RayInfoContainer::const_iterator rit, rit_end = rays.end();
2434
2435        for (rit = rays.begin(); rit != rit_end; ++ rit)
2436        {
2437                RayInfo bRay = *rit;
2438               
2439                VssRay *ray = bRay.mRay;
2440                float t;
2441
2442                // get classification and receive new t
2443                // test if start point behind or in front of plane
2444                const int side = bRay.ComputeRayIntersection(plane.mAxis, plane.mPosition, t);
2445                       
2446                if (side == 0)
2447                {
2448                        ++ splits;
2449
2450                        if (ray->HasPosDir(plane.mAxis))
2451                        {
2452                                backRays.push_back(RayInfo(ray, bRay.GetMinT(), t));
2453                                frontRays.push_back(RayInfo(ray, t, bRay.GetMaxT()));
2454                        }
2455                        else
2456                        {
2457                                frontRays.push_back(RayInfo(ray, bRay.GetMinT(), t));
2458                                backRays.push_back(RayInfo(ray, t, bRay.GetMaxT()));
2459                        }
2460                }
2461                else if (side == 1)
2462                {
2463                        frontRays.push_back(bRay);
2464                }
2465                else
2466                {
2467                        backRays.push_back(bRay);
2468                }
2469        }
2470
2471        return splits;
2472}
2473
2474
2475AxisAlignedBox3 VspTree::GetBoundingBox(VspNode *node) const
2476{
2477        if (!node->GetParent())
2478                return mBoundingBox;
2479
2480        if (!node->IsLeaf())
2481        {
2482                return (dynamic_cast<VspInterior *>(node))->GetBoundingBox();           
2483        }
2484
2485        VspInterior *parent = dynamic_cast<VspInterior *>(node->GetParent());
2486
2487        AxisAlignedBox3 box(parent->GetBoundingBox());
2488
2489        if (parent->GetFront() == node)
2490                box.SetMin(parent->GetAxis(), parent->GetPosition());
2491    else
2492                box.SetMax(parent->GetAxis(), parent->GetPosition());
2493
2494        return box;
2495}
2496
2497
2498int VspTree::ComputeBoxIntersections(const AxisAlignedBox3 &box,
2499                                                                         ViewCellContainer &viewCells) const
2500{
2501        stack<VspNode *> nodeStack;
2502 
2503        ViewCell::NewMail();
2504
2505        while (!nodeStack.empty())
2506        {
2507                VspNode *node = nodeStack.top();
2508                nodeStack.pop();
2509
2510                const AxisAlignedBox3 bbox = GetBoundingBox(node);
2511
2512                if (bbox.Includes(box))
2513                {
2514                        // node geometry is contained in box
2515                        CollectViewCells(node, true, viewCells, true);
2516                }
2517                else if (Overlap(bbox, box))
2518                {
2519                        if (node->IsLeaf())
2520                        {
2521                                BspLeaf *leaf = dynamic_cast<BspLeaf *>(node);
2522                       
2523                                if (!leaf->GetViewCell()->Mailed() && leaf->TreeValid())
2524                                {
2525                                        leaf->GetViewCell()->Mail();
2526                                        viewCells.push_back(leaf->GetViewCell());
2527                                }
2528                        }
2529                        else
2530                        {
2531                                VspInterior *interior = dynamic_cast<VspInterior *>(node);
2532                       
2533                                VspNode *first = interior->GetFront();
2534                                VspNode *second = interior->GetBack();
2535           
2536                                nodeStack.push(first);
2537                                nodeStack.push(second);
2538                        }
2539                }       
2540                // default: cull
2541        }
2542
2543        return (int)viewCells.size();
2544}
2545
2546
2547void VspTree::PreprocessRays(const VssRayContainer &sampleRays,
2548                                                         RayInfoContainer &rays)
2549{
2550        VssRayContainer::const_iterator rit, rit_end = sampleRays.end();
2551
2552        long startTime = GetTime();
2553
2554        cout << "storing rays ... ";
2555
2556        Intersectable::NewMail();
2557
2558        //-- store rays and objects
2559        for (rit = sampleRays.begin(); rit != rit_end; ++ rit)
2560        {
2561                VssRay *ray = *rit;
2562                float minT, maxT;
2563                static Ray hray;
2564
2565                hray.Init(*ray);
2566               
2567                // TODO: not very efficient to implictly cast between rays types
2568                if (GetBoundingBox().GetRaySegment(hray, minT, maxT))
2569                {
2570                        float len = ray->Length();
2571
2572                        if (!len)
2573                                len = Limits::Small;
2574
2575                        rays.push_back(RayInfo(ray, minT / len, maxT / len));
2576                }
2577        }
2578
2579        cout << "finished in " << TimeDiff(startTime, GetTime()) * 1e-3 << " secs" << endl;
2580}
2581
2582
2583void VspTree::GetViewCells(const VssRay &ray, ViewCellContainer &viewCells)
2584{
2585        static Ray hray;
2586        hray.Init(ray);
2587       
2588        float tmin = 0, tmax = 1.0;
2589
2590        if (!mBoundingBox.GetRaySegment(hray, tmin, tmax) || (tmin > tmax))
2591                return;
2592
2593        const Vector3 origin = hray.Extrap(tmin);
2594        const Vector3 termination = hray.Extrap(tmax);
2595
2596        // view cells were not precomputed
2597        // don't mail because we need mailboxing for something else
2598        CastLineSegment(origin, termination, viewCells, false);
2599}
2600
2601
2602SubdivisionCandidate *VspTree::PrepareConstruction(const VssRayContainer &sampleRays,
2603                                                                                                   AxisAlignedBox3 *forcedViewSpace,
2604                                                                                                   RayInfoContainer &rays)
2605{       
2606        mVspStats.Reset();
2607        mVspStats.Start();
2608        mVspStats.nodes = 1;
2609
2610        // store pointer to this tree
2611        VspSubdivisionCandidate::sVspTree = this;
2612       
2613        // compute view space bounding box
2614        ComputeBoundingBox(sampleRays, forcedViewSpace);
2615
2616        // initialise termination criteria
2617        mTermMinProbability *= mBoundingBox.GetVolume();
2618        mGlobalCostMisses = 0;
2619
2620        // get clipped rays
2621        PreprocessRays(sampleRays, rays);
2622
2623        const int pvsSize = EvalPvsSize(rays);
2624       
2625        Debug <<  "pvs size: " << (int)pvsSize << endl;
2626        Debug <<  "rays size: " << (int)rays.size() << endl;
2627
2628        //-- prepare view space partition
2629        const float prop = mBoundingBox.GetVolume();
2630       
2631        // we assume that leaf was already created
2632        VspLeaf *leaf = new VspLeaf();
2633        mRoot = leaf;
2634
2635        // first vsp traversal data
2636        VspTraversalData vData(leaf, 0, &rays, pvsSize, prop, mBoundingBox);
2637
2638        // create first view cell
2639        CreateViewCell(vData, false);
2640               
2641#if WORK_WITH_VIEWCELL_PVS
2642        // add first view cell to all the objects view cell pvs
2643        ObjectPvsMap::const_iterator oit,
2644                oit_end = leaf->GetViewCell()->GetPvs().mEntries.end();
2645
2646        for (oit = leaf->GetViewCell()->GetPvs().mEntries.begin(); oit != oit_end; ++ oit)
2647        {
2648                Intersectable *obj = (*oit).first;
2649                obj->mViewCellPvs.AddSample(leaf->GetViewCell(), 1);
2650        }
2651#endif
2652
2653        //-- compute first split candidate
2654
2655        VspSubdivisionCandidate *splitCandidate = new VspSubdivisionCandidate(vData);
2656    EvalSubdivisionCandidate(*splitCandidate);
2657        leaf->SetSubdivisionCandidate(splitCandidate);
2658
2659        mTotalCost = (float)pvsSize;
2660        EvalSubdivisionStats(*splitCandidate);
2661
2662        return splitCandidate;
2663}
2664
2665
2666void VspTree::CollectDirtyCandidate(const VssRay &ray,
2667                                                                        const bool isTermination,
2668                                                                        vector<SubdivisionCandidate *> &dirtyList) const
2669{
2670
2671        Intersectable *obj;
2672        Vector3 pt;
2673        KdNode *node;
2674
2675        ray.GetSampleData(isTermination, pt, &obj, &node);
2676       
2677        if (!obj) return;
2678
2679        switch (mHierarchyManager->GetObjectSpaceSubdivisonType())
2680        {
2681        case HierarchyManager::KD_BASED_OBJ_SUBDIV:
2682                {
2683                        KdLeaf *leaf = mHierarchyManager->mOspTree->GetLeaf(pt, node);
2684
2685                        if (!leaf->Mailed())
2686                        {
2687                                leaf->Mail();
2688                                dirtyList.push_back(leaf->mSubdivisionCandidate);
2689                        }
2690                        break;
2691                }
2692        case HierarchyManager::BV_BASED_OBJ_SUBDIV:
2693                {
2694                        BvhLeaf *leaf = mHierarchyManager->mBvHierarchy->GetLeaf(obj);
2695
2696                        if (!leaf->Mailed())
2697                        {
2698                                leaf->Mail();
2699                                if (leaf->GetSubdivisionCandidate()) // a candidate still attached to this node
2700                                {
2701                                        dirtyList.push_back(leaf->GetSubdivisionCandidate());
2702                                }
2703                        }
2704                        break;
2705                }
2706        default:
2707                break;
2708        }
2709}
2710
2711
2712void VspTree::CollectDirtyCandidates(VspSubdivisionCandidate *sc,
2713                                                                         vector<SubdivisionCandidate *> &dirtyList)
2714{
2715        VspTraversalData &tData = sc->mParentData;
2716        VspLeaf *node = tData.mNode;
2717       
2718        KdLeaf::NewMail();
2719        BvhLeaf::NewMail();
2720       
2721        RayInfoContainer::const_iterator rit, rit_end = tData.mRays->end();
2722
2723        // add all kd nodes seen by the rays
2724        for (rit = tData.mRays->begin(); rit != rit_end; ++ rit)
2725        {
2726                VssRay *ray = (*rit).mRay;
2727               
2728                CollectDirtyCandidate(*ray, true, dirtyList);
2729                CollectDirtyCandidate(*ray, false, dirtyList);
2730        }
2731}
2732
2733
2734int VspTree::EvalMaxEventContribution(const VssRay &ray,
2735                                                                          const bool isTermination) const
2736{
2737        Intersectable *obj;
2738        Vector3 pt;
2739        KdNode *node;
2740
2741        ray.GetSampleData(isTermination, pt, &obj, &node);
2742
2743        if (!obj)
2744                return 0;
2745
2746        int pvs = 0;
2747
2748        switch (mHierarchyManager->GetObjectSpaceSubdivisonType())
2749        {
2750        case HierarchyManager::NO_OBJ_SUBDIV:
2751                {
2752                        if (-- obj->mCounter == 0)
2753                                ++ pvs;
2754                        break;
2755                }
2756        case HierarchyManager::KD_BASED_OBJ_SUBDIV:
2757                {
2758                        KdLeaf *leaf = mHierarchyManager->mOspTree->GetLeaf(pt, node);
2759
2760                        // add contributions of the kd nodes
2761                        pvs += EvalMaxEventContribution(leaf);
2762                        break;
2763                }
2764        case HierarchyManager::BV_BASED_OBJ_SUBDIV:
2765                {
2766                        BvhLeaf *leaf = mHierarchyManager->mBvHierarchy->GetLeaf(obj);
2767
2768                        if (-- leaf->mCounter == 0)
2769                                pvs += (int)leaf->mObjects.size();
2770                        break;
2771                }
2772        default:
2773                break;
2774        }
2775
2776        return pvs;
2777}
2778
2779
2780int VspTree::PrepareHeuristics(const VssRay &ray, const bool isTermination)
2781{
2782        int pvsSize = 0;
2783       
2784        Intersectable *obj;
2785        Vector3 pt;
2786        KdNode *node;
2787
2788        ray.GetSampleData(isTermination, pt, &obj, &node);
2789
2790        if (!obj)
2791                return 0;
2792
2793        switch (mHierarchyManager->GetObjectSpaceSubdivisonType())
2794        {
2795        case HierarchyManager::NO_OBJ_SUBDIV:
2796                {
2797                        if (!obj->Mailed())
2798                        {
2799                                obj->Mail();
2800                                obj->mCounter = 0;
2801                                ++ pvsSize;
2802                        }
2803                        ++ obj->mCounter;       
2804                        break;
2805                }
2806        case HierarchyManager::KD_BASED_OBJ_SUBDIV:
2807                {
2808                        KdLeaf *leaf = mHierarchyManager->mOspTree->GetLeaf(pt, node);
2809                        pvsSize += PrepareHeuristics(leaf);     
2810                        break;
2811                }
2812        case HierarchyManager::BV_BASED_OBJ_SUBDIV:
2813                {
2814                        BvhLeaf *leaf = mHierarchyManager->mBvHierarchy->GetLeaf(obj);
2815
2816                        if (!leaf->Mailed())
2817                        {
2818                                leaf->Mail();
2819                                leaf->mCounter = 0;
2820                                pvsSize += (int)leaf->mObjects.size();
2821                        }
2822                        ++ leaf->mCounter;     
2823                        break;
2824                }
2825        default:
2826                break;
2827        }
2828
2829        return pvsSize;
2830}
2831
2832
2833int VspTree::EvalMinEventContribution(const VssRay &ray,
2834                                                                          const bool isTermination) const
2835{
2836        Intersectable *obj;
2837        Vector3 pt;
2838        KdNode *node;
2839
2840        ray.GetSampleData(isTermination, pt, &obj, &node);
2841
2842        if (!obj) return 0;
2843
2844        int pvs = 0;
2845
2846        switch (mHierarchyManager->GetObjectSpaceSubdivisonType())
2847        {
2848        case HierarchyManager::NO_OBJ_SUBDIV:
2849                {
2850                        if (!obj->Mailed())
2851                        {
2852                                obj->Mail();
2853                                ++ pvs;
2854                        }
2855                        break;
2856                }
2857        case HierarchyManager::KD_BASED_OBJ_SUBDIV:
2858                {
2859                        KdLeaf *leaf = mHierarchyManager->mOspTree->GetLeaf(pt, node);
2860                        // add contributions of the kd nodes
2861                        pvs += EvalMinEventContribution(leaf);                         
2862                        break;
2863                }
2864        case HierarchyManager::BV_BASED_OBJ_SUBDIV:
2865                {
2866                        BvhLeaf *leaf = mHierarchyManager->mBvHierarchy->GetLeaf(obj);
2867                        if (!leaf->Mailed())
2868                        {
2869                                leaf->Mail();
2870                                pvs += (int)leaf->mObjects.size();
2871                        }
2872                        break;
2873                }
2874        default:
2875                break;
2876        }
2877
2878        return pvs;
2879}
2880
2881
2882void VspTree::UpdateContributionsToPvs(const VssRay &ray,
2883                                                                           const bool isTermination,
2884                                                                           const int cf,
2885                                                                           float &frontPvs,
2886                                                                           float &backPvs,
2887                                                                           float &totalPvs) const
2888{
2889        Intersectable *obj;
2890        Vector3 pt;
2891        KdNode *node;
2892
2893        ray.GetSampleData(isTermination, pt, &obj, &node);
2894
2895        if (!obj) return;
2896
2897        switch (mHierarchyManager->GetObjectSpaceSubdivisonType())
2898        {
2899                case HierarchyManager::NO_OBJ_SUBDIV:
2900                {
2901                        // find front and back pvs for origing and termination object
2902                        UpdateContributionsToPvs(obj, cf, frontPvs, backPvs, totalPvs);
2903                        break;
2904                }
2905                case HierarchyManager::KD_BASED_OBJ_SUBDIV:
2906                {
2907                        KdLeaf *leaf = mHierarchyManager->mOspTree->GetLeaf(pt, node);
2908                        UpdateContributionsToPvs(leaf, cf, frontPvs, backPvs, totalPvs);
2909                        break;
2910                }
2911                case HierarchyManager::BV_BASED_OBJ_SUBDIV:
2912                {
2913                        BvhLeaf *leaf = mHierarchyManager->mBvHierarchy->GetLeaf(obj);
2914                        UpdateContributionsToPvs(leaf, cf, frontPvs, backPvs, totalPvs);
2915                        break;
2916                }
2917        }
2918}
2919
2920
2921int VspTree::EvalContributionToPvs(const VssRay &ray,
2922                                                                   const bool isTermination) const
2923{       
2924        Intersectable *obj;
2925        Vector3 pt;
2926        KdNode *node;
2927
2928        ray.GetSampleData(isTermination, pt, &obj, &node);
2929
2930        if (!obj) return 0;
2931
2932        int pvs = 0;
2933
2934        switch(mHierarchyManager->GetObjectSpaceSubdivisonType())
2935        {
2936        case HierarchyManager::NO_OBJ_SUBDIV:
2937                {
2938                        if (!obj->Mailed())
2939                        {
2940                                obj->Mail();
2941                                ++ pvs;
2942                        }
2943                        break;
2944                }
2945        case HierarchyManager::KD_BASED_OBJ_SUBDIV:
2946                {
2947                        KdLeaf *leaf = mHierarchyManager->mOspTree->GetLeaf(pt, node);
2948
2949                        pvs += EvalContributionToPvs(leaf);
2950                        break;
2951                }
2952        case HierarchyManager::BV_BASED_OBJ_SUBDIV:
2953                {
2954                        BvhLeaf *bvhleaf = mHierarchyManager->mBvHierarchy->GetLeaf(obj);
2955
2956                        if (!bvhleaf->Mailed())
2957                        {
2958                                bvhleaf->Mail();
2959                                pvs += (int)bvhleaf->mObjects.size();
2960                        }
2961                        break;
2962                }
2963        default:
2964                break;
2965        }
2966
2967        return pvs;
2968}
2969
2970
2971int VspTree::EvalContributionToPvs(KdLeaf *leaf) const
2972{
2973        if (leaf->Mailed()) // leaf already mailed
2974                return 0;
2975       
2976        leaf->Mail();
2977
2978        int pvs = (int)(leaf->mObjects.size() - leaf->mMultipleObjects.size());
2979
2980        ObjectContainer::const_iterator oit, oit_end = leaf->mMultipleObjects.end();
2981
2982        for (oit = leaf->mMultipleObjects.begin(); oit != oit_end; ++ oit)
2983        {
2984                Intersectable *obj = *oit;
2985                if (!obj->Mailed())
2986                {
2987                        obj->Mail();
2988                        ++ pvs;
2989                }
2990        }
2991
2992        return pvs;
2993}
2994
2995
2996}
Note: See TracBrowser for help on using the repository browser.