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

Revision 1015, 51.3 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
19namespace GtpVisibilityPreprocessor {
20
21#define USE_FIXEDPOINT_T 0
22
23
24//-- static members
25
26int VspOspTree::sFrontId = 0;
27int VspOspTree::sBackId = 0;
28int VspOspTree::sFrontAndBackId = 0;
29
30
31
32// pvs penalty can be different from pvs size
33inline static float EvalPvsPenalty(const int pvs,
34                                                                   const int lower,
35                                                                   const int upper)
36{
37        // clamp to minmax values
38        if (pvs < lower)
39                return (float)lower;
40        if (pvs > upper)
41                return (float)upper;
42
43        return (float)pvs;
44}
45
46
47int VspNode::sMailId = 1;
48
49
50/******************************************************************/
51/*                  class VspNode implementation                  */
52/******************************************************************/
53
54
55VspNode::VspNode():
56mParent(NULL), mTreeValid(true), mTimeStamp(0)
57{}
58
59
60VspNode::VspNode(VspInterior *parent):
61mParent(parent), mTreeValid(true)
62{}
63
64
65bool VspNode::IsRoot() const
66{
67        return mParent == NULL;
68}
69
70
71VspInterior *VspNode::GetParent()
72{
73        return mParent;
74}
75
76
77void VspNode::SetParent(VspInterior *parent)
78{
79        mParent = parent;
80}
81
82
83bool VspNode::IsSibling(VspNode *n) const
84{
85        return  ((this != n) && mParent &&
86                         (mParent->GetFront() == n) || (mParent->GetBack() == n));
87}
88
89
90int VspNode::GetDepth() const
91{
92        int depth = 0;
93        VspNode *p = mParent;
94       
95        while (p)
96        {
97                p = p->mParent;
98                ++ depth;
99        }
100
101        return depth;
102}
103
104
105bool VspNode::TreeValid() const
106{
107        return mTreeValid;
108}
109
110
111void VspNode::SetTreeValid(const bool v)
112{
113        mTreeValid = v;
114}
115
116
117/****************************************************************/
118/*              class VspInterior implementation                */
119/****************************************************************/
120
121
122VspInterior::VspInterior(const AxisAlignedPlane &plane):
123mPlane(plane), mFront(NULL), mBack(NULL)
124{}
125
126
127VspInterior::~VspInterior()
128{
129        DEL_PTR(mFront);
130        DEL_PTR(mBack);
131}
132
133
134bool VspInterior::IsLeaf() const
135{
136        return false;
137}
138
139
140VspNode *VspInterior::GetBack()
141{
142        return mBack;
143}
144
145
146VspNode *VspInterior::GetFront()
147{
148        return mFront;
149}
150
151
152AxisAlignedPlane VspInterior::GetPlane() const
153{
154        return mPlane;
155}
156
157
158float VspInterior::GetPosition() const
159{
160        return mPlane.mPosition;
161}
162
163
164int VspInterior::GetAxis() const
165{
166        return mPlane.mAxis;
167}
168
169
170void VspInterior::ReplaceChildLink(VspNode *oldChild, VspNode *newChild)
171{
172        if (mBack == oldChild)
173                mBack = newChild;
174        else
175                mFront = newChild;
176}
177
178
179void VspInterior::SetupChildLinks(VspNode *b, VspNode *f)
180{
181    mBack = b;
182    mFront = f;
183}
184
185
186AxisAlignedBox3 VspInterior::GetBox() const
187{
188        return mBox;
189}
190
191
192void VspInterior::SetBox(const AxisAlignedBox3 &box)
193{
194        mBox = box;
195}
196
197
198int VspInterior::Type() const
199{
200        return Interior;
201}
202
203/****************************************************************/
204/*                  class VspLeaf implementation                */
205/****************************************************************/
206
207
208VspLeaf::VspLeaf(): mViewCell(NULL), mPvs(NULL)
209{
210}
211
212
213VspLeaf::~VspLeaf()
214{
215        DEL_PTR(mPvs);
216}
217
218int VspLeaf::Type() const
219{
220        return Leaf;
221}
222
223
224
225VspLeaf::VspLeaf(ViewCellLeaf *viewCell):
226mViewCell(viewCell)
227{
228}
229
230
231VspLeaf::VspLeaf(VspInterior *parent):
232VspNode(parent), mViewCell(NULL), mPvs(NULL)
233{}
234
235
236
237VspLeaf::VspLeaf(VspInterior *parent, ViewCellLeaf *viewCell):
238VspNode(parent), mViewCell(viewCell), mPvs(NULL)
239{
240}
241
242ViewCellLeaf *VspLeaf::GetViewCell() const
243{
244        return mViewCell;
245}
246
247void VspLeaf::SetViewCell(ViewCellLeaf *viewCell)
248{
249        mViewCell = viewCell;
250}
251
252
253bool VspLeaf::IsLeaf() const
254{
255        return true;
256}
257
258
259/******************************************************************************/
260/*                       class VspOspTree implementation                      */
261/******************************************************************************/
262
263
264VspOspTree::VspOspTree():
265mRoot(NULL),
266mOutOfBoundsCell(NULL),
267mStoreRays(false),
268mUseRandomAxis(false),
269mTimeStamp(1)
270{
271        bool randomize = false;
272        Environment::GetSingleton()->GetBoolValue("VspOspTree.Construction.randomize", randomize);
273        if (randomize)
274                Randomize(); // initialise random generator for heuristics
275
276        //-- termination criteria for autopartition
277        Environment::GetSingleton()->GetIntValue("VspOspTree.Termination.maxDepth", mTermMaxDepth);
278        Environment::GetSingleton()->GetIntValue("VspOspTree.Termination.minPvs", mTermMinPvs);
279        Environment::GetSingleton()->GetIntValue("VspOspTree.Termination.minRays", mTermMinRays);
280        Environment::GetSingleton()->GetFloatValue("VspOspTree.Termination.minProbability", mTermMinProbability);
281        Environment::GetSingleton()->GetFloatValue("VspOspTree.Termination.maxRayContribution", mTermMaxRayContribution);
282        Environment::GetSingleton()->GetFloatValue("VspOspTree.Termination.maxCostRatio", mTermMaxCostRatio);
283        Environment::GetSingleton()->GetIntValue("VspOspTree.Termination.missTolerance", mTermMissTolerance);
284        Environment::GetSingleton()->GetIntValue("VspOspTree.Termination.maxViewCells", mMaxViewCells);
285
286        //-- max cost ratio for early tree termination
287        Environment::GetSingleton()->GetFloatValue("VspOspTree.Termination.maxCostRatio", mTermMaxCostRatio);
288
289        Environment::GetSingleton()->GetFloatValue("VspOspTree.Termination.minGlobalCostRatio", mTermMinGlobalCostRatio);
290        Environment::GetSingleton()->GetIntValue("VspOspTree.Termination.globalCostMissTolerance", mTermGlobalCostMissTolerance);
291
292        // HACK//mTermMinPolygons = 25;
293
294        //-- factors for bsp tree split plane heuristics
295        Environment::GetSingleton()->GetFloatValue("VspOspTree.Termination.ct_div_ci", mCtDivCi);
296
297        //-- partition criteria
298        Environment::GetSingleton()->GetFloatValue("VspOspTree.Construction.epsilon", mEpsilon);
299
300        // if only the driving axis is used for axis aligned split
301        Environment::GetSingleton()->GetBoolValue("VspOspTree.splitUseOnlyDrivingAxis", mOnlyDrivingAxis);
302       
303        //Environment::GetSingleton()->GetFloatValue("VspOspTree.maxTotalMemory", mMaxTotalMemory);
304        Environment::GetSingleton()->GetFloatValue("VspOspTree.maxStaticMemory", mMaxMemory);
305
306        Environment::GetSingleton()->GetBoolValue("VspOspTree.useCostHeuristics", mUseCostHeuristics);
307        Environment::GetSingleton()->GetBoolValue("VspOspTree.simulateOctree", mCirculatingAxis);
308        Environment::GetSingleton()->GetBoolValue("VspOspTree.useRandomAxis", mUseRandomAxis);
309        Environment::GetSingleton()->GetIntValue("VspOspTree.nodePriorityQueueType", mNodePriorityQueueType);
310       
311        char subdivisionStatsLog[100];
312        Environment::GetSingleton()->GetStringValue("VspOspTree.subdivisionStats", subdivisionStatsLog);
313        mSubdivisionStats.open(subdivisionStatsLog);
314
315        Environment::GetSingleton()->GetFloatValue("VspOspTree.Construction.minBand", mMinBand);
316        Environment::GetSingleton()->GetFloatValue("VspOspTree.Construction.maxBand", mMaxBand);
317       
318
319        //-- debug output
320
321        Debug << "******* VSP BSP options ******** " << endl;
322    Debug << "max depth: " << mTermMaxDepth << endl;
323        Debug << "min PVS: " << mTermMinPvs << endl;
324        Debug << "min probabiliy: " << mTermMinProbability << endl;
325        Debug << "min rays: " << mTermMinRays << endl;
326        Debug << "max ray contri: " << mTermMaxRayContribution << endl;
327        Debug << "max cost ratio: " << mTermMaxCostRatio << endl;
328        Debug << "miss tolerance: " << mTermMissTolerance << endl;
329        Debug << "max view cells: " << mMaxViewCells << endl;
330        Debug << "randomize: " << randomize << endl;
331
332        Debug << "min global cost ratio: " << mTermMinGlobalCostRatio << endl;
333        Debug << "global cost miss tolerance: " << mTermGlobalCostMissTolerance << endl;
334        Debug << "only driving axis: " << mOnlyDrivingAxis << endl;
335        Debug << "max memory: " << mMaxMemory << endl;
336        Debug << "use cost heuristics: " << mUseCostHeuristics << endl;
337        Debug << "subdivision stats log: " << subdivisionStatsLog << endl;
338        Debug << "use random axis: " << mUseRandomAxis << endl;
339        Debug << "priority queue type: " << mNodePriorityQueueType << endl;
340        Debug << "circulating axis: " << mCirculatingAxis << endl;
341        Debug << "minband: " << mMinBand << endl;
342        Debug << "maxband: " << mMaxBand << endl;
343       
344
345        mSplitCandidates = new vector<SortableEntry>;
346
347        Debug << endl;
348}
349
350
351VspViewCell *VspOspTree::GetOutOfBoundsCell()
352{
353        return mOutOfBoundsCell;
354}
355
356
357VspViewCell *VspOspTree::GetOrCreateOutOfBoundsCell()
358{
359        if (!mOutOfBoundsCell)
360        {
361                mOutOfBoundsCell = new VspViewCell();
362                mOutOfBoundsCell->SetId(-1);
363                mOutOfBoundsCell->SetValid(false);
364        }
365
366        return mOutOfBoundsCell;
367}
368
369
370const VspTreeStatistics &VspOspTree::GetStatistics() const
371{
372        return mVspStats;
373}
374
375
376VspOspTree::~VspOspTree()
377{
378        DEL_PTR(mRoot);
379        DEL_PTR(mSplitCandidates);
380}
381
382
383void VspOspTree::Construct(const VssRayContainer &sampleRays,
384                                                   AxisAlignedBox3 *forcedBoundingBox)
385{
386        mVspStats.nodes = 1;
387        mBox.Initialize();      // initialise BSP tree bounding box
388
389        if (forcedBoundingBox)
390                mBox = *forcedBoundingBox;
391       
392        PolygonContainer polys;
393        RayInfoContainer *rays = new RayInfoContainer();
394
395        VssRayContainer::const_iterator rit, rit_end = sampleRays.end();
396
397        long startTime = GetTime();
398
399        cout << "Extracting polygons from rays ... ";
400
401        Intersectable::NewMail();
402
403        int numObj = 0;
404
405        //-- store rays
406        for (rit = sampleRays.begin(); rit != rit_end; ++ rit)
407        {
408                VssRay *ray = *rit;
409
410                float minT, maxT;
411
412                static Ray hray;
413                hray.Init(*ray);
414
415                // TODO: not very efficient to implictly cast between rays types
416                if (mBox.GetRaySegment(hray, minT, maxT))
417                {
418                        float len = ray->Length();
419
420                        if (!len)
421                                len = Limits::Small;
422
423                        rays->push_back(RayInfo(ray, minT / len, maxT / len));
424                }
425        }
426
427        mTermMinProbability *= mBox.GetVolume();
428
429        mGlobalCostMisses = 0;
430
431        cout << "finished" << endl;
432
433        // use split cost priority queue
434        Construct(rays);
435}
436
437
438// TODO: return memory usage in MB
439float VspOspTree::GetMemUsage() const
440{
441        return (float)
442                 (sizeof(VspOspTree) +
443                  mVspStats.Leaves() * sizeof(VspLeaf) +
444                  mCreatedViewCells * sizeof(VspViewCell) +
445                  mVspStats.pvs * sizeof(ObjectPvsData) +
446                  mVspStats.Interior() * sizeof(VspInterior) +
447                  mVspStats.accumRays * sizeof(RayInfo)) / (1024.0f * 1024.0f);
448}
449
450
451void VspOspTree::Construct(RayInfoContainer *rays)
452{
453        VspOspSplitQueue tQueue;
454
455        mRoot = new VspLeaf();
456
457        const float prop = mBox.GetVolume();
458
459        VspOspTraversalData tData(mRoot,
460                                                          0,
461                                                          rays,
462                              ComputePvsSize(*rays),
463                                                          prop,
464                                                          mBox);
465
466
467        // compute first split candidate
468        VspOspSplitCandidate splitCandidate;
469        EvalSplitCandidate(tData, splitCandidate);
470
471        tQueue.push(splitCandidate);
472
473        mTotalCost = tData.mPvs * tData.mProbability / mBox.GetVolume();
474        mTotalPvsSize = tData.mPvs;
475       
476        mSubdivisionStats
477                        << "#ViewCells\n1\n" <<  endl
478                        << "#RenderCostDecrease\n0\n" << endl
479                        << "#SplitCandidateCost\n0\n" << endl
480                        << "#TotalRenderCost\n" << mTotalCost << endl
481                        << "#AvgRenderCost\n" << mTotalPvsSize << endl;
482
483        Debug << "total cost: " << mTotalCost << endl;
484       
485       
486        mVspStats.Start();
487        cout << "Constructing vsp bsp tree ... \n";
488
489        long startTime = GetTime();     
490        int nLeaves = 500;
491        int nViewCells = 500;
492
493        // used for intermediate time measurements and progress
494        long interTime = GetTime();     
495
496        mOutOfMemory = false;
497
498        mCreatedViewCells = 0;
499       
500        while (!tQueue.empty())
501        {
502                splitCandidate = tQueue.top();
503            tQueue.pop();               
504
505                // cost ratio of cost decrease / totalCost
506                float costRatio = splitCandidate.GetCost() / mTotalCost;
507
508                //Debug << "cost ratio: " << costRatio << endl;
509
510                if (costRatio < mTermMinGlobalCostRatio)
511                        ++ mGlobalCostMisses;
512               
513                if (0 && !mOutOfMemory)
514                {
515                        float mem = GetMemUsage();
516
517                        if (mem > mMaxMemory)
518                        {
519                                mOutOfMemory = true;
520                                Debug << "memory limit reached: " << mem << endl;
521                        }
522                }
523
524                // subdivide leaf node
525                VspNode *r = Subdivide(tQueue, splitCandidate);
526
527                if (r == mRoot)
528                        Debug << "VSP BSP tree construction time spent at root: "
529                                  << TimeDiff(startTime, GetTime())*1e-3 << "s" << endl;
530
531                if (mVspStats.Leaves() >= nLeaves)
532                {
533                        nLeaves += 500;
534
535                        cout << "leaves=" << mVspStats.Leaves() << endl;
536                        Debug << "needed "
537                                  << TimeDiff(interTime, GetTime())*1e-3
538                                  << " secs to create 500 view cells" << endl;
539                        interTime = GetTime();
540                }
541
542                if (mCreatedViewCells == nViewCells)
543                {
544                        nViewCells += 500;
545
546                        cout << "generated " << mCreatedViewCells << " viewcells" << endl;
547                }
548        }
549
550        Debug << "Used Memory: " << GetMemUsage() << " MB" << endl << endl;
551        cout << "finished\n";
552
553        mVspStats.Stop();
554}
555
556
557bool VspOspTree::LocalTerminationCriteriaMet(const VspOspTraversalData &data) const
558{
559        return
560                (((int)data.mRays->size() <= mTermMinRays) ||
561                 (data.mPvs <= mTermMinPvs)   ||
562                 (data.mProbability <= mTermMinProbability) ||
563                 (data.GetAvgRayContribution() > mTermMaxRayContribution) ||
564                 (data.mDepth >= mTermMaxDepth));
565}
566
567
568bool VspOspTree::GlobalTerminationCriteriaMet(const VspOspTraversalData &data) const
569{
570        return
571                (mOutOfMemory ||
572                (mVspStats.Leaves() >= mMaxViewCells) ||
573                (mGlobalCostMisses >= mTermGlobalCostMissTolerance));
574}
575
576
577// subdivide using a split plane queue
578VspNode *VspOspTree::Subdivide(VspOspSplitQueue &tQueue,
579                                                           VspOspSplitCandidate &splitCandidate)
580{
581        VspOspTraversalData &tData = splitCandidate.mParentData;
582
583        VspNode *newNode = tData.mNode;
584
585        if (!LocalTerminationCriteriaMet(tData) && !GlobalTerminationCriteriaMet(tData))
586        {       
587                VspOspTraversalData tFrontData;
588                VspOspTraversalData tBackData;
589
590                //-- continue subdivision
591               
592                // create new interior node and two leaf node
593                const AxisAlignedPlane splitPlane = splitCandidate.mSplitPlane;
594
595                newNode = SubdivideNode(splitPlane, tData, tFrontData, tBackData);
596       
597                const int maxCostMisses = splitCandidate.mMaxCostMisses;
598
599
600                // how often was max cost ratio missed in this branch?
601                tFrontData.mMaxCostMisses = maxCostMisses;
602                tBackData.mMaxCostMisses = maxCostMisses;
603                       
604               
605                if (1)
606                {
607                        float cFront = (float)tFrontData.mPvs * tFrontData.mProbability;
608                        float cBack = (float)tBackData.mPvs * tBackData.mProbability;
609                        float cData = (float)tData.mPvs * tData.mProbability;
610
611                       
612                        float costDecr =
613                                (cFront + cBack - cData) / mBox.GetVolume();
614
615                        mTotalCost += costDecr;
616                        mTotalPvsSize += tFrontData.mPvs + tBackData.mPvs - tData.mPvs;
617
618                        mSubdivisionStats
619                                        << "#ViewCells\n" << mVspStats.Leaves() << endl
620                                        << "#RenderCostDecrease\n" << -costDecr << endl
621                                        << "#SplitCandidateCost\n" << splitCandidate.GetCost() << endl
622                                        << "#TotalRenderCost\n" << mTotalCost << endl
623                                        << "#AvgRenderCost\n" << (float)mTotalPvsSize / (float)mVspStats.Leaves() << endl;
624                }
625
626       
627                //-- push the new split candidates on the stack
628                VspOspSplitCandidate frontCandidate;
629                VspOspSplitCandidate backCandidate;
630
631                EvalSplitCandidate(tFrontData, frontCandidate);
632                EvalSplitCandidate(tBackData, backCandidate);
633       
634                tQueue.push(frontCandidate);
635                tQueue.push(backCandidate);
636       
637                // delete old leaf node
638                DEL_PTR(tData.mNode);
639        }
640
641
642        //-- terminate traversal and create new view cell
643        if (newNode->IsLeaf())
644        {
645                VspLeaf *leaf = dynamic_cast<VspLeaf *>(newNode);
646                VspViewCell *viewCell = new VspViewCell();
647
648        leaf->SetViewCell(viewCell);
649               
650                //-- update pvs
651                int conSamp = 0;
652                float sampCon = 0.0f;
653                AddToPvs(leaf, *tData.mRays, sampCon, conSamp);
654
655                // update scalar pvs size value
656                mViewCellsManager->SetScalarPvsSize(viewCell, viewCell->GetPvs().GetSize());
657
658                mVspStats.contributingSamples += conSamp;
659                mVspStats.sampleContributions +=(int) sampCon;
660
661                //-- store additional info
662                if (mStoreRays)
663                {
664                        RayInfoContainer::const_iterator it, it_end = tData.mRays->end();
665                        for (it = tData.mRays->begin(); it != it_end; ++ it)
666                        {
667                                (*it).mRay->Ref();                     
668                                leaf->mVssRays.push_back((*it).mRay);
669                        }
670                }
671
672                // should I check here?
673                viewCell->mLeaf = leaf;
674
675                viewCell->SetVolume(tData.mProbability);
676        leaf->mProbability = tData.mProbability;
677
678                EvaluateLeafStats(tData);               
679        }
680
681        //-- cleanup
682        tData.Clear();
683
684        return newNode;
685}
686
687
688void VspOspTree::EvalSplitCandidate(VspOspTraversalData &tData,
689                                                                        VspOspSplitCandidate &splitData)
690{
691        VspOspTraversalData frontData;
692        VspOspTraversalData backData;
693
694        VspLeaf *leaf = dynamic_cast<VspLeaf *>(tData.mNode);
695
696        // compute locally best split plane
697        const bool success =
698                SelectPlane(tData, splitData.mSplitPlane,
699                                        frontData.mProbability, backData.mProbability);
700
701        //TODO
702        // compute global decrease in render cost
703        splitData.mRenderCost = EvalRenderCostDecrease(splitData.mSplitPlane, tData);
704        splitData.mParentData = tData;
705        splitData.mMaxCostMisses = success ? tData.mMaxCostMisses : tData.mMaxCostMisses + 1;
706}
707
708
709VspInterior *VspOspTree::SubdivideNode(const AxisAlignedPlane &splitPlane,
710                                                                           VspOspTraversalData &tData,
711                                                                           VspOspTraversalData &frontData,
712                                                                           VspOspTraversalData &backData)
713{
714        VspLeaf *leaf = dynamic_cast<VspLeaf *>(tData.mNode);
715       
716        //-- the front and back traversal data is filled with the new values
717        frontData.mDepth = tData.mDepth + 1;
718        frontData.mRays = new RayInfoContainer();
719       
720        backData.mDepth = tData.mDepth + 1;
721        backData.mRays = new RayInfoContainer();
722       
723        //-- subdivide rays
724        SplitRays(splitPlane,
725                          *tData.mRays,
726                          *frontData.mRays,
727                          *backData.mRays);
728
729        //-- compute pvs
730        frontData.mPvs = ComputePvsSize(*frontData.mRays);
731        backData.mPvs = ComputePvsSize(*backData.mRays);
732
733        // split front and back node geometry and compute area
734        tData.mBoundingBox.Split(splitPlane.mAxis, splitPlane.mPosition,
735                                                         frontData.mBoundingBox, backData.mBoundingBox);
736
737
738        frontData.mProbability = frontData.mBoundingBox.GetVolume();
739        backData.mProbability = tData.mProbability - frontData.mProbability;
740
741       
742    ///////////////////////////////////////////
743        // subdivide further
744
745        // store maximal and minimal depth
746        if (tData.mDepth > mVspStats.maxDepth)
747        {
748                Debug << "max depth increases to " << tData.mDepth << " at " << mVspStats.Leaves() << " leaves" << endl;
749                mVspStats.maxDepth = tData.mDepth;
750        }
751
752        mVspStats.nodes += 2;
753
754   
755        VspInterior *interior = new VspInterior(splitPlane);
756
757#ifdef _DEBUG
758        Debug << interior << endl;
759#endif
760
761
762        //-- create front and back leaf
763
764        VspInterior *parent = leaf->GetParent();
765
766        // replace a link from node's parent
767        if (parent)
768        {
769                parent->ReplaceChildLink(leaf, interior);
770                interior->SetParent(parent);
771        }
772        else // new root
773        {
774                mRoot = interior;
775        }
776
777        // and setup child links
778        interior->SetupChildLinks(new VspLeaf(interior), new VspLeaf(interior));
779
780        frontData.mNode = interior->GetFront();
781        backData.mNode = interior->GetBack();
782
783        interior->mTimeStamp = mTimeStamp ++;
784       
785        return interior;
786}
787
788
789void VspOspTree::AddToPvs(VspLeaf *leaf,
790                                                  const RayInfoContainer &rays,
791                                                  float &sampleContributions,
792                                                  int &contributingSamples)
793{
794        sampleContributions = 0;
795        contributingSamples = 0;
796 
797        RayInfoContainer::const_iterator it, it_end = rays.end();
798 
799        ViewCellLeaf *vc = leaf->GetViewCell();
800 
801        // add contributions from samples to the PVS
802        for (it = rays.begin(); it != it_end; ++ it)
803        {
804                float sc = 0.0f;
805                VssRay *ray = (*it).mRay;
806
807                bool madeContrib = false;
808                float contribution;
809
810                if (ray->mTerminationObject)
811                {
812                        if (vc->AddPvsSample(ray->mTerminationObject, ray->mPdf, contribution))
813                                madeContrib = true;
814                        sc += contribution;
815                }
816         
817                if (ray->mOriginObject)
818                {
819                        if (vc->AddPvsSample(ray->mOriginObject, ray->mPdf, contribution))
820                                madeContrib = true;
821                        sc += contribution;
822                }
823         
824          sampleContributions += sc;
825
826          if (madeContrib)
827                  ++ contributingSamples;
828               
829          //leaf->mVssRays.push_back(ray);
830        }
831}
832
833
834void VspOspTree::SortSplitCandidates(const RayInfoContainer &rays,
835                                                                         const int axis,
836                                                                         float minBand,
837                                                                         float maxBand)
838{
839        mSplitCandidates->clear();
840
841        int requestedSize = 2 * (int)(rays.size());
842        // creates a sorted split candidates array
843        if (mSplitCandidates->capacity() > 500000 &&
844                requestedSize < (int)(mSplitCandidates->capacity() / 10) )
845        {
846        delete mSplitCandidates;
847                mSplitCandidates = new vector<SortableEntry>;
848        }
849
850        mSplitCandidates->reserve(requestedSize);
851
852        float pos;
853
854        // float values => don't compare with exact values
855        if (0)
856        {
857                minBand += Limits::Small;
858                maxBand -= Limits::Small;
859        }
860
861        // insert all queries
862        for (RayInfoContainer::const_iterator ri = rays.begin(); ri < rays.end(); ++ ri)
863        {
864                const bool positive = (*ri).mRay->HasPosDir(axis);
865                               
866                pos = (*ri).ExtrapOrigin(axis);
867                // clamp to min / max band
868                if (0) ClipValue(pos, minBand, maxBand);
869               
870                mSplitCandidates->push_back(SortableEntry(positive ? SortableEntry::ERayMin : SortableEntry::ERayMax,
871                                                                        pos, (*ri).mRay));
872
873                pos = (*ri).ExtrapTermination(axis);
874                // clamp to min / max band
875                if (0) ClipValue(pos, minBand, maxBand);
876
877                mSplitCandidates->push_back(SortableEntry(positive ? SortableEntry::ERayMax : SortableEntry::ERayMin,
878                                                                        pos, (*ri).mRay));
879        }
880
881        stable_sort(mSplitCandidates->begin(), mSplitCandidates->end());
882}
883
884
885float VspOspTree::BestCostRatioHeuristics(const RayInfoContainer &rays,
886                                                                                  const AxisAlignedBox3 &box,
887                                                                                  const int pvsSize,
888                                                                                  const int axis,
889                                          float &position)
890{
891        const float minBox = box.Min(axis);
892        const float maxBox = box.Max(axis);
893
894        const float sizeBox = maxBox - minBox;
895
896        const float minBand = minBox + mMinBand * sizeBox;
897        const float maxBand = minBox + mMaxBand * sizeBox;
898
899        SortSplitCandidates(rays, axis, minBand, maxBand);
900
901        // go through the lists, count the number of objects left and right
902        // and evaluate the following cost funcion:
903        // C = ct_div_ci  + (ql*rl + qr*rr)/queries
904
905        int pvsl = 0;
906        int pvsr = pvsSize;
907
908        int pvsBack = pvsl;
909        int pvsFront = pvsr;
910
911        float sum = (float)pvsSize * sizeBox;
912        float minSum = 1e20f;
913
914       
915        // if no border can be found, take mid split
916        position = minBox + 0.5f * sizeBox;
917       
918        // the relative cost ratio
919        float ratio = /*Limits::Infinity;*/99999999.0f;
920        bool splitPlaneFound = false;
921
922        Intersectable::NewMail();
923
924        RayInfoContainer::const_iterator ri, ri_end = rays.end();
925
926        // set all object as belonging to the front pvs
927        for(ri = rays.begin(); ri != ri_end; ++ ri)
928        {
929                Intersectable *oObject = (*ri).mRay->mOriginObject;
930                Intersectable *tObject = (*ri).mRay->mTerminationObject;
931
932                if (oObject)
933                {
934                        if (!oObject->Mailed())
935                        {
936                                oObject->Mail();
937                                oObject->mCounter = 1;
938                        }
939                        else
940                        {
941                                ++ oObject->mCounter;
942                        }
943                }
944                if (tObject)
945                {
946                        if (!tObject->Mailed())
947                        {
948                                tObject->Mail();
949                                tObject->mCounter = 1;
950                        }
951                        else
952                        {
953                                ++ tObject->mCounter;
954                        }
955                }
956        }
957
958        Intersectable::NewMail();
959
960        vector<SortableEntry>::const_iterator ci, ci_end = mSplitCandidates->end();
961
962        for (ci = mSplitCandidates->begin(); ci != ci_end; ++ ci)
963        {
964                VssRay *ray;
965                ray = (*ci).ray;
966               
967                Intersectable *oObject = ray->mOriginObject;
968                Intersectable *tObject = ray->mTerminationObject;
969               
970
971                switch ((*ci).type)
972                {
973                        case SortableEntry::ERayMin:
974                                {
975                                        if (oObject && !oObject->Mailed())
976                                        {
977                                                oObject->Mail();
978                                                ++ pvsl;
979                                        }
980                                        if (tObject && !tObject->Mailed())
981                                        {
982                                                tObject->Mail();
983                                                ++ pvsl;
984                                        }
985                                        break;
986                                }
987                        case SortableEntry::ERayMax:
988                                {
989                                        if (oObject)
990                                        {
991                                                if (-- oObject->mCounter == 0)
992                                                        -- pvsr;
993                                        }
994
995                                        if (tObject)
996                                        {
997                                                if (-- tObject->mCounter == 0)
998                                                        -- pvsr;
999                                        }
1000
1001                                        break;
1002                                }
1003                }
1004
1005               
1006               
1007                // Note: sufficient to compare size of bounding boxes of front and back side?
1008                if (((*ci).value >= minBand) && ((*ci).value <= maxBand))
1009                {
1010                        sum = pvsl * ((*ci).value - minBox) + pvsr * (maxBox - (*ci).value);
1011
1012                        //Debug  << "pos=" << (*ci).value << "\t pvs=(" <<  pvsl << "," << pvsr << ")" << endl;
1013                        //Debug << "cost= " << sum << endl;
1014
1015                        if (sum < minSum)
1016                        {
1017                                splitPlaneFound = true;
1018
1019                                minSum = sum;
1020                                position = (*ci).value;
1021                               
1022                                pvsBack = pvsl;
1023                                pvsFront = pvsr;
1024                        }
1025                }
1026        }
1027       
1028       
1029        // -- compute cost
1030        const int lowerPvsLimit = mViewCellsManager->GetMinPvsSize();
1031        const int upperPvsLimit = mViewCellsManager->GetMaxPvsSize();
1032
1033        const float pOverall = sizeBox;
1034
1035        const float pBack = position - minBox;
1036        const float pFront = maxBox - position;
1037       
1038        const float penaltyOld = EvalPvsPenalty(pvsSize, lowerPvsLimit, upperPvsLimit);
1039    const float penaltyFront = EvalPvsPenalty(pvsFront, lowerPvsLimit, upperPvsLimit);
1040        const float penaltyBack = EvalPvsPenalty(pvsBack, lowerPvsLimit, upperPvsLimit);
1041       
1042        const float oldRenderCost = penaltyOld * pOverall;
1043        const float newRenderCost = penaltyFront * pFront + penaltyBack * pBack;
1044
1045        if (splitPlaneFound)
1046        {
1047                ratio = newRenderCost / (oldRenderCost + Limits::Small);
1048        }
1049        //if (axis != 1)
1050        //Debug << "axis=" << axis << " costRatio=" << ratio << " pos=" << position << " t=" << (position - minBox) / (maxBox - minBox)
1051         //    <<"\t pb=(" << pvsBack << ")\t pf=(" << pvsFront << ")" << endl;
1052
1053        return ratio;
1054}
1055
1056
1057float VspOspTree::SelectPlane(const VspOspTraversalData &tData,
1058                                                          AxisAlignedPlane &plane,
1059                                                          float &pFront,
1060                                                          float &pBack)
1061{
1062        float nPosition[3];
1063        float nCostRatio[3];
1064        float nProbFront[3];
1065        float nProbBack[3];
1066
1067        // create bounding box of node geometry
1068        AxisAlignedBox3 box = tData.mBoundingBox;
1069               
1070        int sAxis = 0;
1071        int bestAxis = -1;
1072
1073
1074        // if we use some kind of specialised fixed axis
1075    const bool useSpecialAxis =
1076                mOnlyDrivingAxis || mUseRandomAxis || mCirculatingAxis;
1077
1078        if (mUseRandomAxis)
1079                sAxis = Random(3);
1080        else if (mCirculatingAxis)
1081                sAxis = (tData.mAxis + 1) % 3;
1082        else if (mOnlyDrivingAxis)
1083                sAxis = box.Size().DrivingAxis();
1084
1085
1086        for (int axis = 0; axis < 3 ; ++ axis)
1087        {
1088                if (!useSpecialAxis || (axis == sAxis))
1089                {
1090                        //-- place split plane using heuristics
1091
1092                        if (mUseCostHeuristics)
1093                        {
1094                                nCostRatio[axis] =
1095                                        BestCostRatioHeuristics(*tData.mRays,
1096                                                                                    box,
1097                                                                                        tData.mPvs,
1098                                                                                        axis,
1099                                                                                        nPosition[axis]);                       
1100                        }
1101                        else //-- split plane position is spatial median
1102                        {
1103                                nPosition[axis] = (box.Min()[axis] + box.Max()[axis]) * 0.5f;
1104
1105                                nCostRatio[axis] = EvalSplitCost(tData,
1106                                                                                                 box,
1107                                                                                                 axis,
1108                                                                                                 nPosition[axis],
1109                                                                                                 nProbFront[axis],
1110                                                                                                 nProbBack[axis]);
1111                        }
1112                                               
1113                        if (bestAxis == -1)
1114                        {
1115                                bestAxis = axis;
1116                        }
1117                        else if (nCostRatio[axis] < nCostRatio[bestAxis])
1118                        {
1119                                bestAxis = axis;
1120                        }
1121                }
1122        }
1123
1124        //-- assign values
1125        plane.mAxis = bestAxis;
1126        pFront = nProbFront[bestAxis];
1127        pBack = nProbBack[bestAxis];
1128
1129        //-- split plane position
1130    plane.mPosition = nPosition[bestAxis];
1131
1132        return nCostRatio[bestAxis];
1133}
1134
1135
1136inline void VspOspTree::GenerateUniqueIdsForPvs()
1137{
1138        Intersectable::NewMail(); sBackId = Intersectable::sMailId;
1139        Intersectable::NewMail(); sFrontId = Intersectable::sMailId;
1140        Intersectable::NewMail(); sFrontAndBackId = Intersectable::sMailId;
1141}
1142
1143
1144float VspOspTree::EvalRenderCostDecrease(const AxisAlignedPlane &candidatePlane,
1145                                                                                 const VspOspTraversalData &data) const
1146{
1147        float pvsFront = 0;
1148        float pvsBack = 0;
1149        float totalPvs = 0;
1150
1151        // probability that view point lies in back / front node
1152        float pOverall = data.mProbability;
1153        float pFront = 0;
1154        float pBack = 0;
1155
1156
1157        // create unique ids for pvs heuristics
1158        GenerateUniqueIdsForPvs();
1159       
1160        RayInfoContainer::const_iterator rit, rit_end = data.mRays.end();
1161
1162        for (rit = data.mRays.begin(); rit != rit_end; ++ rit)
1163        {
1164                RayInfo rayInf = *rit;
1165
1166                float t;
1167                VssRay *ray = rayInf.mRay;
1168                const int cf =
1169                        rayInf.ComputeRayIntersection(candidatePlane.mAxis,
1170                                                                                  candidatePlane.mPosition, t);
1171
1172                // find front and back pvs for origing and termination object
1173                AddObjToPvs(ray->mTerminationObject, cf, pvsFront, pvsBack, totalPvs);
1174                AddObjToPvs(ray->mOriginObject, cf, pvsFront, pvsBack, totalPvs);
1175        }
1176
1177        AxisAlignedBox3 frontBox;
1178        AxisAlignedBox3 backBox;
1179
1180        data.mBoundingBox.Split(candidatePlane.mAxis, candidatePlane.mPosition, frontBox, backBox);
1181
1182        pFront = frontBox.GetVolume();
1183        pBack = pOverall - pFront;
1184               
1185
1186        // -- pvs rendering heuristics
1187        const int lowerPvsLimit = mViewCellsManager->GetMinPvsSize();
1188        const int upperPvsLimit = mViewCellsManager->GetMaxPvsSize();
1189
1190        // only render cost heuristics or combined with standard deviation
1191        const float penaltyOld = EvalPvsPenalty(totalPvs, lowerPvsLimit, upperPvsLimit);
1192    const float penaltyFront = EvalPvsPenalty(pvsFront, lowerPvsLimit, upperPvsLimit);
1193        const float penaltyBack = EvalPvsPenalty(pvsBack, lowerPvsLimit, upperPvsLimit);
1194                       
1195        const float oldRenderCost = pOverall * penaltyOld;
1196        const float newRenderCost = penaltyFront * pFront + penaltyBack * pBack;
1197
1198        //Debug << "decrease: " << oldRenderCost - newRenderCost << endl;
1199        const float renderCostDecrease = (oldRenderCost - newRenderCost) / mBox.GetVolume();
1200       
1201
1202#if 1
1203        // take render cost of node into account to avoid being stuck in a local minimum
1204        const float normalizedOldRenderCost = oldRenderCost / mBox.GetVolume();
1205        return 0.99f * renderCostDecrease + 0.01f * normalizedOldRenderCost;
1206#else
1207        return renderCostDecrease;
1208#endif
1209}
1210
1211
1212
1213float VspOspTree::EvalSplitCost(const VspOspTraversalData &data,
1214                                                                const AxisAlignedBox3 &box,
1215                                                                const int axis,
1216                                                                const float &position,                                                                           
1217                                                                float &pFront,
1218                                                                float &pBack) const
1219{
1220        float pvsTotal = 0;
1221        float pvsFront = 0;
1222        float pvsBack = 0;
1223       
1224        // create unique ids for pvs heuristics
1225        GenerateUniqueIdsForPvs();
1226
1227        const int pvsSize = data.mPvs;
1228
1229        RayInfoContainer::const_iterator rit, rit_end = data.mRays->end();
1230
1231        // this is the main ray classification loop!
1232        for(rit = data.mRays->begin(); rit != rit_end; ++ rit)
1233        {
1234                // determine the side of this ray with respect to the plane
1235                float t;
1236                const int side = (*rit).ComputeRayIntersection(axis, position, t);
1237       
1238                AddObjToPvs((*rit).mRay->mTerminationObject, side, pvsFront, pvsBack, pvsTotal);
1239                AddObjToPvs((*rit).mRay->mOriginObject, side, pvsFront, pvsBack, pvsTotal);
1240        }
1241
1242        //-- pvs heuristics
1243        float pOverall;
1244
1245        //-- compute heurstics
1246        //-- we take simplified computation for mid split
1247               
1248        pOverall = data.mProbability;
1249        pBack = pFront = pOverall * 0.5f;
1250       
1251       
1252#ifdef _DEBUG
1253        Debug << axis << " " << pvsSize << " " << pvsBack << " " << pvsFront << endl;
1254        Debug << pFront << " " << pBack << " " << pOverall << endl;
1255#endif
1256
1257       
1258        const float newCost = pvsBack * pBack + pvsFront * pFront;
1259        const float oldCost = (float)pvsSize * pOverall + Limits::Small;
1260
1261        return  (mCtDivCi + newCost) / oldCost;
1262}
1263
1264
1265void VspOspTree::AddObjToPvs(Intersectable *obj,
1266                                                         const int cf,
1267                                                         float &frontPvs,
1268                                                         float &backPvs,
1269                                                         float &totalPvs) const
1270{
1271        if (!obj)
1272                return;
1273
1274        const float renderCost = mViewCellsManager->EvalRenderCost(obj);
1275
1276        // new object
1277        if ((obj->mMailbox != sFrontId) &&
1278                (obj->mMailbox != sBackId) &&
1279                (obj->mMailbox != sFrontAndBackId))
1280        {
1281                totalPvs += renderCost;
1282        }
1283
1284        // TODO: does this really belong to no pvs?
1285        //if (cf == Ray::COINCIDENT) return;
1286
1287        // object belongs to both PVS
1288        if (cf >= 0)
1289        {
1290                if ((obj->mMailbox != sFrontId) &&
1291                        (obj->mMailbox != sFrontAndBackId))
1292                {
1293                        frontPvs += renderCost;
1294               
1295                        if (obj->mMailbox == sBackId)
1296                                obj->mMailbox = sFrontAndBackId;
1297                        else
1298                                obj->mMailbox = sFrontId;
1299                }
1300        }
1301
1302        if (cf <= 0)
1303        {
1304                if ((obj->mMailbox != sBackId) &&
1305                        (obj->mMailbox != sFrontAndBackId))
1306                {
1307                        backPvs += renderCost;
1308               
1309                        if (obj->mMailbox == sFrontId)
1310                                obj->mMailbox = sFrontAndBackId;
1311                        else
1312                                obj->mMailbox = sBackId;
1313                }
1314        }
1315}
1316
1317
1318void VspOspTree::CollectLeaves(vector<VspLeaf *> &leaves,
1319                                                           const bool onlyUnmailed,
1320                                                           const int maxPvsSize) const
1321{
1322        stack<VspNode *> nodeStack;
1323        nodeStack.push(mRoot);
1324
1325        while (!nodeStack.empty())
1326        {
1327                VspNode *node = nodeStack.top();
1328                nodeStack.pop();
1329               
1330                if (node->IsLeaf())
1331                {
1332                        // test if this leaf is in valid view space
1333                        VspLeaf *leaf = dynamic_cast<VspLeaf *>(node);
1334                        if (leaf->TreeValid() &&
1335                                (!onlyUnmailed || !leaf->Mailed()) &&
1336                                ((maxPvsSize < 0) || (leaf->GetViewCell()->GetPvs().GetSize() <= maxPvsSize)))
1337                        {
1338                                leaves.push_back(leaf);
1339                        }
1340                }
1341                else
1342                {
1343                        VspInterior *interior = dynamic_cast<VspInterior *>(node);
1344
1345                        nodeStack.push(interior->GetBack());
1346                        nodeStack.push(interior->GetFront());
1347                }
1348        }
1349}
1350
1351
1352AxisAlignedBox3 VspOspTree::GetBoundingBox() const
1353{
1354        return mBox;
1355}
1356
1357
1358VspNode *VspOspTree::GetRoot() const
1359{
1360        return mRoot;
1361}
1362
1363
1364void VspOspTree::EvaluateLeafStats(const VspOspTraversalData &data)
1365{
1366        // the node became a leaf -> evaluate stats for leafs
1367        VspLeaf *leaf = dynamic_cast<VspLeaf *>(data.mNode);
1368
1369
1370        if (data.mPvs > mVspStats.maxPvs)
1371        {
1372                mVspStats.maxPvs = data.mPvs;
1373        }
1374
1375        mVspStats.pvs += data.mPvs;
1376
1377        if (data.mDepth < mVspStats.minDepth)
1378        {
1379                mVspStats.minDepth = data.mDepth;
1380        }
1381       
1382        if (data.mDepth >= mTermMaxDepth)
1383        {
1384        ++ mVspStats.maxDepthNodes;
1385                //Debug << "new max depth: " << mVspStats.maxDepthNodes << endl;
1386        }
1387
1388        // accumulate rays to compute rays /  leaf
1389        mVspStats.accumRays += (int)data.mRays->size();
1390
1391        if (data.mPvs < mTermMinPvs)
1392                ++ mVspStats.minPvsNodes;
1393
1394        if ((int)data.mRays->size() < mTermMinRays)
1395                ++ mVspStats.minRaysNodes;
1396
1397        if (data.GetAvgRayContribution() > mTermMaxRayContribution)
1398                ++ mVspStats.maxRayContribNodes;
1399
1400        if (data.mProbability <= mTermMinProbability)
1401                ++ mVspStats.minProbabilityNodes;
1402       
1403        // accumulate depth to compute average depth
1404        mVspStats.accumDepth += data.mDepth;
1405
1406        ++ mCreatedViewCells;
1407
1408#ifdef _DEBUG
1409        Debug << "BSP stats: "
1410                  << "Depth: " << data.mDepth << " (max: " << mTermMaxDepth << "), "
1411                  << "PVS: " << data.mPvs << " (min: " << mTermMinPvs << "), "
1412          //              << "Area: " << data.mProbability << " (min: " << mTermMinProbability << "), "
1413                  << "#rays: " << (int)data.mRays->size() << " (max: " << mTermMinRays << "), "
1414                  << "#pvs: " << leaf->GetViewCell()->GetPvs().GetSize() << "=, "
1415                  << "#avg ray contrib (pvs): " << (float)data.mPvs / (float)data.mRays->size() << endl;
1416#endif
1417}
1418
1419
1420void VspOspTree::CollectViewCells(ViewCellContainer &viewCells, bool onlyValid) const
1421{
1422        ViewCell::NewMail();
1423        CollectViewCells(mRoot, onlyValid, viewCells, true);
1424}
1425
1426
1427void VspOspTree::CollapseViewCells()
1428{
1429// TODO
1430#if HAS_TO_BE_REDONE
1431        stack<VspNode *> nodeStack;
1432
1433        if (!mRoot)
1434                return;
1435
1436        nodeStack.push(mRoot);
1437       
1438        while (!nodeStack.empty())
1439        {
1440                VspNode *node = nodeStack.top();
1441                nodeStack.pop();
1442               
1443                if (node->IsLeaf())
1444        {
1445                        BspViewCell *viewCell = dynamic_cast<VspLeaf *>(node)->GetViewCell();
1446
1447                        if (!viewCell->GetValid())
1448                        {
1449                                BspViewCell *viewCell = dynamic_cast<VspLeaf *>(node)->GetViewCell();
1450       
1451                                ViewCellContainer leaves;
1452                                mViewCellsTree->CollectLeaves(viewCell, leaves);
1453
1454                                ViewCellContainer::const_iterator it, it_end = leaves.end();
1455
1456                                for (it = leaves.begin(); it != it_end; ++ it)
1457                                {
1458                                        VspLeaf *l = dynamic_cast<BspViewCell *>(*it)->mLeaf;
1459                                        l->SetViewCell(GetOrCreateOutOfBoundsCell());
1460                                        ++ mVspStats.invalidLeaves;
1461                                }
1462
1463                                // add to unbounded view cell
1464                                GetOrCreateOutOfBoundsCell()->GetPvs().AddPvs(viewCell->GetPvs());
1465                                DEL_PTR(viewCell);
1466                        }
1467                }
1468                else
1469                {
1470                        VspInterior *interior = dynamic_cast<VspInterior *>(node);
1471               
1472                        nodeStack.push(interior->GetFront());
1473                        nodeStack.push(interior->GetBack());
1474                }
1475        }
1476
1477        Debug << "invalid leaves: " << mVspStats.invalidLeaves << endl;
1478#endif
1479}
1480
1481
1482void VspOspTree::CollectRays(VssRayContainer &rays)
1483{
1484        vector<VspLeaf *> leaves;
1485
1486        vector<VspLeaf *>::const_iterator lit, lit_end = leaves.end();
1487
1488        for (lit = leaves.begin(); lit != lit_end; ++ lit)
1489        {
1490                VspLeaf *leaf = *lit;
1491                VssRayContainer::const_iterator rit, rit_end = leaf->mVssRays.end();
1492
1493                for (rit = leaf->mVssRays.begin(); rit != rit_end; ++ rit)
1494                        rays.push_back(*rit);
1495        }
1496}
1497
1498
1499void VspOspTree::ValidateTree()
1500{
1501        stack<VspNode *> nodeStack;
1502
1503        if (!mRoot)
1504                return;
1505
1506        nodeStack.push(mRoot);
1507       
1508        mVspStats.invalidLeaves = 0;
1509        while (!nodeStack.empty())
1510        {
1511                VspNode *node = nodeStack.top();
1512                nodeStack.pop();
1513               
1514                if (node->IsLeaf())
1515                {
1516                        VspLeaf *leaf = dynamic_cast<VspLeaf *>(node);
1517
1518                        if (!leaf->GetViewCell()->GetValid())
1519                                ++ mVspStats.invalidLeaves;
1520
1521                        // validity flags don't match => repair
1522                        if (leaf->GetViewCell()->GetValid() != leaf->TreeValid())
1523                        {
1524                                leaf->SetTreeValid(leaf->GetViewCell()->GetValid());
1525                                PropagateUpValidity(leaf);
1526                        }
1527                }
1528                else
1529                {
1530                        VspInterior *interior = dynamic_cast<VspInterior *>(node);
1531               
1532                        nodeStack.push(interior->GetFront());
1533                        nodeStack.push(interior->GetBack());
1534                }
1535        }
1536
1537        Debug << "invalid leaves: " << mVspStats.invalidLeaves << endl;
1538}
1539
1540
1541
1542void VspOspTree::CollectViewCells(VspNode *root,
1543                                                                  bool onlyValid,
1544                                                                  ViewCellContainer &viewCells,
1545                                                                  bool onlyUnmailed) const
1546{
1547        stack<VspNode *> nodeStack;
1548
1549        if (!root)
1550                return;
1551
1552        nodeStack.push(root);
1553       
1554        while (!nodeStack.empty())
1555        {
1556                VspNode *node = nodeStack.top();
1557                nodeStack.pop();
1558               
1559                if (node->IsLeaf())
1560                {
1561                        if (!onlyValid || node->TreeValid())
1562                        {
1563                                ViewCellLeaf *leafVc = dynamic_cast<VspLeaf *>(node)->GetViewCell();
1564
1565                                ViewCell *viewCell = mViewCellsTree->GetActiveViewCell(leafVc);
1566                                               
1567                                if (!onlyUnmailed || !viewCell->Mailed())
1568                                {
1569                                        viewCell->Mail();
1570                                        viewCells.push_back(viewCell);
1571                                }
1572                        }
1573                }
1574                else
1575                {
1576                        VspInterior *interior = dynamic_cast<VspInterior *>(node);
1577               
1578                        nodeStack.push(interior->GetFront());
1579                        nodeStack.push(interior->GetBack());
1580                }
1581        }
1582
1583}
1584
1585
1586int VspOspTree::FindNeighbors(VspLeaf *n,
1587                                                          vector<VspLeaf *> &neighbors,
1588                                                          const bool onlyUnmailed) const
1589{
1590        stack<VspNode *> nodeStack;
1591
1592        nodeStack.push(mRoot);
1593
1594        const AxisAlignedBox3 box = GetBBox(n);
1595
1596        while (!nodeStack.empty())
1597        {
1598                VspNode *node = nodeStack.top();
1599                nodeStack.pop();
1600
1601                if (node->IsLeaf())
1602                {
1603                        VspLeaf *leaf = dynamic_cast<VspLeaf *>(node);
1604
1605                        if (leaf != n && (!onlyUnmailed || !leaf->Mailed()))
1606                                neighbors.push_back(leaf);
1607                }
1608                else
1609                {
1610                        VspInterior *interior = dynamic_cast<VspInterior *>(node);
1611                       
1612                        if (interior->GetPosition() > box.Max(interior->GetAxis()))
1613                                nodeStack.push(interior->GetBack());
1614                        else
1615                        {
1616                                if (interior->GetPosition() < box.Min(interior->GetAxis()))
1617                                        nodeStack.push(interior->GetFront());
1618                                else
1619                                {
1620                                        // random decision
1621                                        nodeStack.push(interior->GetBack());
1622                                        nodeStack.push(interior->GetFront());
1623                                }
1624                        }
1625                }
1626        }
1627
1628        return (int)neighbors.size();
1629}
1630
1631
1632// Find random neighbor which was not mailed
1633VspLeaf *VspOspTree::GetRandomLeaf(const Plane3 &plane)
1634{
1635        stack<VspNode *> nodeStack;
1636        nodeStack.push(mRoot);
1637 
1638        int mask = rand();
1639 
1640        while (!nodeStack.empty())
1641        {
1642                VspNode *node = nodeStack.top();
1643               
1644                nodeStack.pop();
1645               
1646                if (node->IsLeaf())
1647                {
1648                        return dynamic_cast<VspLeaf *>(node);
1649                }
1650                else
1651                {
1652                        VspInterior *interior = dynamic_cast<VspInterior *>(node);
1653                        VspNode *next;
1654                       
1655                        if (GetBBox(interior->GetBack()).Side(plane) < 0)
1656                                next = interior->GetFront();
1657            else
1658                                if (GetBBox(interior->GetFront()).Side(plane) < 0)
1659                                        next = interior->GetBack();
1660                                else
1661                                {
1662                                        // random decision
1663                                        if (mask & 1)
1664                                                next = interior->GetBack();
1665                                        else
1666                                                next = interior->GetFront();
1667                                        mask = mask >> 1;
1668                                }
1669
1670                                nodeStack.push(next);
1671                }
1672        }
1673 
1674        return NULL;
1675}
1676
1677
1678VspLeaf *VspOspTree::GetRandomLeaf(const bool onlyUnmailed)
1679{
1680        stack<VspNode *> nodeStack;
1681
1682        nodeStack.push(mRoot);
1683
1684        int mask = rand();
1685
1686        while (!nodeStack.empty())
1687        {
1688                VspNode *node = nodeStack.top();
1689                nodeStack.pop();
1690
1691                if (node->IsLeaf())
1692                {
1693                        if ( (!onlyUnmailed || !node->Mailed()) )
1694                                return dynamic_cast<VspLeaf *>(node);
1695                }
1696                else
1697                {
1698                        VspInterior *interior = dynamic_cast<VspInterior *>(node);
1699
1700                        // random decision
1701                        if (mask & 1)
1702                                nodeStack.push(interior->GetBack());
1703                        else
1704                                nodeStack.push(interior->GetFront());
1705
1706                        mask = mask >> 1;
1707                }
1708        }
1709
1710        return NULL;
1711}
1712
1713
1714int VspOspTree::ComputePvsSize(const RayInfoContainer &rays) const
1715{
1716        int pvsSize = 0;
1717
1718        RayInfoContainer::const_iterator rit, rit_end = rays.end();
1719
1720        Intersectable::NewMail();
1721
1722        for (rit = rays.begin(); rit != rays.end(); ++ rit)
1723        {
1724                VssRay *ray = (*rit).mRay;
1725
1726                if (ray->mOriginObject)
1727                {
1728                        if (!ray->mOriginObject->Mailed())
1729                        {
1730                                ray->mOriginObject->Mail();
1731                                ++ pvsSize;
1732                        }
1733                }
1734                if (ray->mTerminationObject)
1735                {
1736                        if (!ray->mTerminationObject->Mailed())
1737                        {
1738                                ray->mTerminationObject->Mail();
1739                                ++ pvsSize;
1740                        }
1741                }
1742        }
1743
1744        return pvsSize;
1745}
1746
1747
1748float VspOspTree::GetEpsilon() const
1749{
1750        return mEpsilon;
1751}
1752
1753
1754int VspOspTree::CastLineSegment(const Vector3 &origin,
1755                                                                const Vector3 &termination,
1756                                                                ViewCellContainer &viewcells)
1757{
1758        int hits = 0;
1759
1760        float mint = 0.0f, maxt = 1.0f;
1761        const Vector3 dir = termination - origin;
1762
1763        stack<LineTraversalData> tStack;
1764
1765        Intersectable::NewMail();
1766        ViewCell::NewMail();
1767
1768        Vector3 entp = origin;
1769        Vector3 extp = termination;
1770
1771        VspNode *node = mRoot;
1772        VspNode *farChild;
1773
1774        float position;
1775        int axis;
1776
1777        while (1)
1778        {
1779                if (!node->IsLeaf())
1780                {
1781                        VspInterior *in = dynamic_cast<VspInterior *>(node);
1782                        position = in->GetPosition();
1783                        axis = in->GetAxis();
1784
1785                        if (entp[axis] <= position)
1786                        {
1787                                if (extp[axis] <= position)
1788                                {
1789                                        node = in->GetBack();
1790                                        // cases N1,N2,N3,P5,Z2,Z3
1791                                        continue;
1792                                } else
1793                                {
1794                                        // case N4
1795                                        node = in->GetBack();
1796                                        farChild = in->GetFront();
1797                                }
1798                        }
1799                        else
1800                        {
1801                                if (position <= extp[axis])
1802                                {
1803                                        node = in->GetFront();
1804                                        // cases P1,P2,P3,N5,Z1
1805                                        continue;
1806                                }
1807                                else
1808                                {
1809                                        node = in->GetFront();
1810                                        farChild = in->GetBack();
1811                                        // case P4
1812                                }
1813                        }
1814
1815                        // $$ modification 3.5.2004 - hints from Kamil Ghais
1816                        // case N4 or P4
1817                        const float tdist = (position - origin[axis]) / dir[axis];
1818                        tStack.push(LineTraversalData(farChild, extp, maxt)); //TODO
1819
1820                        extp = origin + dir * tdist;
1821                        maxt = tdist;
1822                }
1823                else
1824                {
1825                        // compute intersection with all objects in this leaf
1826                        VspLeaf *leaf = dynamic_cast<VspLeaf *>(node);
1827                        ViewCell *vc = leaf->GetViewCell();
1828
1829                        if (!vc->Mailed())
1830                        {
1831                                vc->Mail();
1832                                viewcells.push_back(vc);
1833                                ++ hits;
1834                        }
1835#if 0
1836                        leaf->mRays.push_back(RayInfo(new VssRay(origin, termination, NULL, NULL, 0)));
1837#endif
1838                        // get the next node from the stack
1839                        if (tStack.empty())
1840                                break;
1841
1842                        entp = extp;
1843                        mint = maxt;
1844                       
1845                        LineTraversalData &s  = tStack.top();
1846                        node = s.mNode;
1847                        extp = s.mExitPoint;
1848                        maxt = s.mMaxT;
1849                        tStack.pop();
1850                }
1851        }
1852
1853        return hits;
1854}
1855
1856
1857int VspOspTree::CastRay(Ray &ray)
1858{
1859        int hits = 0;
1860
1861        stack<LineTraversalData> tStack;
1862        const Vector3 dir = ray.GetDir();
1863
1864        float maxt, mint;
1865
1866        if (!mBox.GetRaySegment(ray, mint, maxt))
1867                return 0;
1868
1869        Intersectable::NewMail();
1870        ViewCell::NewMail();
1871
1872        Vector3 entp = ray.Extrap(mint);
1873        Vector3 extp = ray.Extrap(maxt);
1874
1875        const Vector3 origin = entp;
1876
1877        VspNode *node = mRoot;
1878        VspNode *farChild = NULL;
1879
1880        float position;
1881        int axis;
1882
1883        while (1)
1884        {
1885                if (!node->IsLeaf())
1886                {
1887                        VspInterior *in = dynamic_cast<VspInterior *>(node);
1888                        position = in->GetPosition();
1889                        axis = in->GetAxis();
1890
1891                        if (entp[axis] <= position)
1892                        {
1893                                if (extp[axis] <= position)
1894                                {
1895                                        node = in->GetBack();
1896                                        // cases N1,N2,N3,P5,Z2,Z3
1897                                        continue;
1898                                } else
1899                                {
1900                                        // case N4
1901                                        node = in->GetBack();
1902                                        farChild = in->GetFront();
1903                                }
1904                        }
1905                        else
1906                        {
1907                                if (position <= extp[axis])
1908                                {
1909                                        node = in->GetFront();
1910                                        // cases P1,P2,P3,N5,Z1
1911                                        continue;
1912                                }
1913                                else
1914                                {
1915                                        node = in->GetFront();
1916                                        farChild = in->GetBack();
1917                                        // case P4
1918                                }
1919                        }
1920
1921                        // $$ modification 3.5.2004 - hints from Kamil Ghais
1922                        // case N4 or P4
1923                        const float tdist = (position - origin[axis]) / dir[axis];
1924                        tStack.push(LineTraversalData(farChild, extp, maxt)); //TODO
1925                        extp = origin + dir * tdist;
1926                        maxt = tdist;
1927                }
1928                else
1929                {
1930                        // compute intersection with all objects in this leaf
1931                        VspLeaf *leaf = dynamic_cast<VspLeaf *>(node);
1932                        ViewCell *vc = leaf->GetViewCell();
1933
1934                        if (!vc->Mailed())
1935                        {
1936                                vc->Mail();
1937                                // todo: add view cells to ray
1938                                ++ hits;
1939                        }
1940
1941#if 0
1942                                leaf->mRays.push_back(RayInfo(new VssRay(origin, termination, NULL, NULL, 0)));
1943#endif
1944                        // get the next node from the stack
1945                        if (tStack.empty())
1946                                break;
1947
1948                        entp = extp;
1949                        mint = maxt;
1950                       
1951                        LineTraversalData &s  = tStack.top();
1952                        node = s.mNode;
1953                        extp = s.mExitPoint;
1954                        maxt = s.mMaxT;
1955                        tStack.pop();
1956                }
1957        }
1958
1959        return hits;
1960}
1961
1962
1963VspNode *VspOspTree::CollapseTree(VspNode *node, int &collapsed)
1964{
1965// TODO
1966#if HAS_TO_BE_REDONE
1967        if (node->IsLeaf())
1968                return node;
1969
1970        VspInterior *interior = dynamic_cast<VspInterior *>(node);
1971
1972        VspNode *front = CollapseTree(interior->GetFront(), collapsed);
1973        VspNode *back = CollapseTree(interior->GetBack(), collapsed);
1974
1975        if (front->IsLeaf() && back->IsLeaf())
1976        {
1977                VspLeaf *frontLeaf = dynamic_cast<VspLeaf *>(front);
1978                VspLeaf *backLeaf = dynamic_cast<VspLeaf *>(back);
1979
1980                //-- collapse tree
1981                if (frontLeaf->GetViewCell() == backLeaf->GetViewCell())
1982                {
1983                        BspViewCell *vc = frontLeaf->GetViewCell();
1984
1985                        VspLeaf *leaf = new VspLeaf(interior->GetParent(), vc);
1986                        leaf->SetTreeValid(frontLeaf->TreeValid());
1987
1988                        // replace a link from node's parent
1989                        if (leaf->GetParent())
1990                                leaf->GetParent()->ReplaceChildLink(node, leaf);
1991                        else
1992                                mRoot = leaf;
1993
1994                        ++ collapsed;
1995                        delete interior;
1996
1997                        return leaf;
1998                }
1999        }
2000#endif
2001        return node;
2002}
2003
2004
2005int VspOspTree::CollapseTree()
2006{
2007        int collapsed = 0;
2008        //TODO
2009#if HAS_TO_BE_REDONE
2010        (void)CollapseTree(mRoot, collapsed);
2011
2012        // revalidate leaves
2013        RepairViewCellsLeafLists();
2014#endif
2015        return collapsed;
2016}
2017
2018
2019void VspOspTree::RepairViewCellsLeafLists()
2020{
2021// TODO
2022#if HAS_TO_BE_REDONE
2023        // list not valid anymore => clear
2024        stack<VspNode *> nodeStack;
2025        nodeStack.push(mRoot);
2026
2027        ViewCell::NewMail();
2028
2029        while (!nodeStack.empty())
2030        {
2031                VspNode *node = nodeStack.top();
2032                nodeStack.pop();
2033
2034                if (node->IsLeaf())
2035                {
2036                        VspLeaf *leaf = dynamic_cast<VspLeaf *>(node);
2037
2038                        BspViewCell *viewCell = leaf->GetViewCell();
2039
2040                        if (!viewCell->Mailed())
2041                        {
2042                                viewCell->mLeaves.clear();
2043                                viewCell->Mail();
2044                        }
2045       
2046                        viewCell->mLeaves.push_back(leaf);
2047
2048                }
2049                else
2050                {
2051                        VspInterior *interior = dynamic_cast<VspInterior *>(node);
2052
2053                        nodeStack.push(interior->GetFront());
2054                        nodeStack.push(interior->GetBack());
2055                }
2056        }
2057// TODO
2058#endif
2059}
2060
2061
2062
2063ViewCell *VspOspTree::GetViewCell(const Vector3 &point, const bool active)
2064{
2065        if (mRoot == NULL)
2066                return NULL;
2067
2068        stack<VspNode *> nodeStack;
2069        nodeStack.push(mRoot);
2070 
2071        ViewCellLeaf *viewcell = NULL;
2072 
2073        while (!nodeStack.empty()) 
2074        {
2075                VspNode *node = nodeStack.top();
2076                nodeStack.pop();
2077       
2078                if (node->IsLeaf())
2079                {
2080                        viewcell = dynamic_cast<VspLeaf *>(node)->GetViewCell();
2081                        break;
2082                }
2083                else   
2084                {       
2085                        VspInterior *interior = dynamic_cast<VspInterior *>(node);
2086     
2087                        // random decision
2088                        if (interior->GetPosition() - point[interior->GetAxis()] < 0)
2089                                nodeStack.push(interior->GetBack());
2090                        else
2091                                nodeStack.push(interior->GetFront());
2092                }
2093        }
2094 
2095        if (active)
2096                return mViewCellsTree->GetActiveViewCell(viewcell);
2097        else
2098                return viewcell;
2099}
2100
2101
2102bool VspOspTree::ViewPointValid(const Vector3 &viewPoint) const
2103{
2104        VspNode *node = mRoot;
2105
2106        while (1)
2107        {
2108                // early exit
2109                if (node->TreeValid())
2110                        return true;
2111
2112                if (node->IsLeaf())
2113                        return false;
2114                       
2115                VspInterior *in = dynamic_cast<VspInterior *>(node);
2116                                       
2117                if (in->GetPosition() - viewPoint[in->GetAxis()] <= 0)
2118                {
2119                        node = in->GetBack();
2120                }
2121                else
2122                {
2123                        node = in->GetFront();
2124                }
2125        }
2126
2127        // should never come here
2128        return false;
2129}
2130
2131
2132void VspOspTree::PropagateUpValidity(VspNode *node)
2133{
2134        const bool isValid = node->TreeValid();
2135
2136        // propagative up invalid flag until only invalid nodes exist over this node
2137        if (!isValid)
2138        {
2139                while (!node->IsRoot() && node->GetParent()->TreeValid())
2140                {
2141                        node = node->GetParent();
2142                        node->SetTreeValid(false);
2143                }
2144        }
2145        else
2146        {
2147                // propagative up valid flag until one of the subtrees is invalid
2148                while (!node->IsRoot() && !node->TreeValid())
2149                {
2150            node = node->GetParent();
2151                        VspInterior *interior = dynamic_cast<VspInterior *>(node);
2152                       
2153                        // the parent is valid iff both leaves are valid
2154                        node->SetTreeValid(interior->GetBack()->TreeValid() &&
2155                                                           interior->GetFront()->TreeValid());
2156                }
2157        }
2158}
2159
2160#if ZIPPED_VIEWCELLS
2161bool VspOspTree::Export(ogzstream &stream)
2162#else
2163bool VspOspTree::Export(ofstream &stream)
2164#endif
2165{
2166        ExportNode(mRoot, stream);
2167
2168        return true;
2169}
2170
2171
2172#if ZIPPED_VIEWCELLS
2173void VspOspTree::ExportNode(VspNode *node, ogzstream &stream)
2174#else
2175void VspOspTree::ExportNode(VspNode *node, ofstream &stream)
2176#endif
2177{
2178        if (node->IsLeaf())
2179        {
2180                VspLeaf *leaf = dynamic_cast<VspLeaf *>(node);
2181                ViewCell *viewCell = mViewCellsTree->GetActiveViewCell(leaf->GetViewCell());
2182
2183                int id = -1;
2184                if (viewCell != mOutOfBoundsCell)
2185                        id = viewCell->GetId();
2186
2187                stream << "<Leaf viewCellId=\"" << id << "\" />" << endl;
2188        }
2189        else
2190        {       
2191                VspInterior *interior = dynamic_cast<VspInterior *>(node);
2192       
2193                AxisAlignedPlane plane = interior->GetPlane();
2194                stream << "<Interior plane=\"" << plane.mPosition << " "
2195                           << plane.mAxis << "\">" << endl;
2196
2197                ExportNode(interior->GetBack(), stream);
2198                ExportNode(interior->GetFront(), stream);
2199
2200                stream << "</Interior>" << endl;
2201        }
2202}
2203
2204
2205int VspOspTree::SplitRays(const AxisAlignedPlane &plane,
2206                                                  RayInfoContainer &rays,
2207                                                  RayInfoContainer &frontRays,
2208                                                  RayInfoContainer &backRays) const
2209{
2210        int splits = 0;
2211
2212        RayInfoContainer::const_iterator rit, rit_end = rays.end();
2213
2214        for (rit = rays.begin(); rit != rit_end; ++ rit)
2215        {
2216                RayInfo bRay = *rit;
2217               
2218                VssRay *ray = bRay.mRay;
2219                float t;
2220
2221                // get classification and receive new t
2222                //-- test if start point behind or in front of plane
2223                const int side = plane.ComputeRayIntersection(bRay, t);
2224
2225                if (side == 0)
2226                {
2227                        ++ splits;
2228
2229                        if (ray->HasPosDir(plane.mAxis))
2230                        {
2231                                backRays.push_back(RayInfo(ray, bRay.GetMinT(), t));
2232                                frontRays.push_back(RayInfo(ray, t, bRay.GetMaxT()));
2233                        }
2234                        else
2235                        {
2236                                frontRays.push_back(RayInfo(ray, bRay.GetMinT(), t));
2237                                backRays.push_back(RayInfo(ray, t, bRay.GetMaxT()));
2238                        }
2239                }
2240                else if (side == 1)
2241                {
2242                        frontRays.push_back(bRay);
2243                }
2244                else
2245                {
2246                        backRays.push_back(bRay);
2247                }
2248
2249        }
2250
2251        return splits;
2252}
2253
2254
2255AxisAlignedBox3 VspOspTree::GetBBox(VspNode *node) const
2256{
2257        if (!node->GetParent())
2258                return mBox;
2259
2260        if (!node->IsLeaf())
2261        {
2262                return (dynamic_cast<VspInterior *>(node))->GetBox();           
2263        }
2264
2265        VspInterior *parent = dynamic_cast<VspInterior *>(node->GetParent());
2266
2267        AxisAlignedBox3 box(parent->GetBox());
2268
2269        if (parent->GetFront() == node)
2270                box.SetMin(parent->GetAxis(), parent->GetPosition());
2271        else
2272                box.SetMax(parent->GetAxis(), parent->GetPosition());
2273
2274        return box;
2275}
2276
2277
2278}
Note: See TracBrowser for help on using the repository browser.