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

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