source: GTP/trunk/Lib/Vis/Preprocessing/src/VspOspTree.cpp @ 1123

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