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

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