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

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