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

Revision 1123, 95.7 KB checked in by mattausch, 18 years ago (diff)
RevLine 
[1002]1#include <stack>
2#include <time.h>
3#include <iomanip>
4
[1010]5#include "ViewCell.h"
[1002]6#include "Plane3.h"
[1006]7#include "VspOspTree.h"
[1002]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"
[1020]18#include "KdTree.h"
[1002]19
[1020]20
[1002]21namespace GtpVisibilityPreprocessor {
22
23#define USE_FIXEDPOINT_T 0
24
[1123]25
[1002]26//-- static members
27
[1099]28VspTree *VspTree::VspSplitCandidate::sVspTree = NULL;
29OspTree *OspTree::OspSplitCandidate::sOspTree = NULL;
[1002]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;
[1077]39        else if (pvs > upper)
[1002]40                return (float)upper;
41
42        return (float)pvs;
43}
44
45
[1008]46int VspNode::sMailId = 1;
[1002]47
48
[1021]49
50void VspTreeStatistics::Print(ostream &app) const
51{
[1121]52        app << "=========== VspTree statistics ===============\n";
[1021]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
[1121]101        app << "========== END OF VspTree statistics ==========\n";
[1021]102}
103
104
[1008]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
[1011]181
[1008]182VspInterior::~VspInterior()
183{
184        DEL_PTR(mFront);
185        DEL_PTR(mBack);
186}
187
[1011]188
[1008]189bool VspInterior::IsLeaf() const
190{
191        return false;
192}
193
[1011]194
[1008]195VspNode *VspInterior::GetBack()
196{
197        return mBack;
198}
199
[1011]200
[1008]201VspNode *VspInterior::GetFront()
202{
203        return mFront;
204}
205
[1011]206
[1008]207AxisAlignedPlane VspInterior::GetPlane() const
208{
209        return mPlane;
210}
211
[1011]212
[1012]213float VspInterior::GetPosition() const
214{
215        return mPlane.mPosition;
216}
217
218
219int VspInterior::GetAxis() const
220{
221        return mPlane.mAxis;
222}
223
224
[1008]225void VspInterior::ReplaceChildLink(VspNode *oldChild, VspNode *newChild)
226{
227        if (mBack == oldChild)
228                mBack = newChild;
229        else
230                mFront = newChild;
231}
232
[1011]233
[1109]234void VspInterior::SetupChildLinks(VspNode *front, VspNode *back)
[1008]235{
[1109]236    mBack = back;
237    mFront = front;
[1008]238}
239
240
[1016]241AxisAlignedBox3 VspInterior::GetBoundingBox() const
[1012]242{
[1016]243        return mBoundingBox;
[1012]244}
245
246
[1016]247void VspInterior::SetBoundingBox(const AxisAlignedBox3 &box)
[1012]248{
[1016]249        mBoundingBox = box;
[1012]250}
251
252
253int VspInterior::Type() const
254{
255        return Interior;
256}
257
[1027]258
259
[1008]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);
[1074]273        CLEAR_CONTAINER(mVssRays);
[1008]274}
275
[1021]276
[1012]277int VspLeaf::Type() const
278{
279        return Leaf;
280}
[1008]281
[1012]282
[1008]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
[1002]317/******************************************************************************/
[1016]318/*                       class VspTree implementation                      */
[1002]319/******************************************************************************/
320
321
[1016]322VspTree::VspTree():
[1002]323mRoot(NULL),
324mOutOfBoundsCell(NULL),
325mStoreRays(false),
326mTimeStamp(1)
327{
328        bool randomize = false;
[1016]329        Environment::GetSingleton()->GetBoolValue("VspTree.Construction.randomize", randomize);
[1002]330        if (randomize)
331                Randomize(); // initialise random generator for heuristics
332
333        //-- termination criteria for autopartition
[1016]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);
[1023]339       
[1016]340        Environment::GetSingleton()->GetIntValue("VspTree.Termination.missTolerance", mTermMissTolerance);
341        Environment::GetSingleton()->GetIntValue("VspTree.Termination.maxViewCells", mMaxViewCells);
[1002]342
343        //-- max cost ratio for early tree termination
[1016]344        Environment::GetSingleton()->GetFloatValue("VspTree.Termination.maxCostRatio", mTermMaxCostRatio);
[1002]345
[1016]346        Environment::GetSingleton()->GetFloatValue("VspTree.Termination.minGlobalCostRatio", mTermMinGlobalCostRatio);
347        Environment::GetSingleton()->GetIntValue("VspTree.Termination.globalCostMissTolerance", mTermGlobalCostMissTolerance);
[1002]348
349        //-- factors for bsp tree split plane heuristics
[1016]350        Environment::GetSingleton()->GetFloatValue("VspTree.Termination.ct_div_ci", mCtDivCi);
[1002]351
352        //-- partition criteria
[1016]353        Environment::GetSingleton()->GetFloatValue("VspTree.Construction.epsilon", mEpsilon);
[1121]354        Environment::GetSingleton()->GetFloatValue("VspTree.Construction.renderCostDecreaseWeight", mRenderCostDecreaseWeight);
[1002]355
356        // if only the driving axis is used for axis aligned split
[1016]357        Environment::GetSingleton()->GetBoolValue("VspTree.splitUseOnlyDrivingAxis", mOnlyDrivingAxis);
[1002]358       
[1016]359        //Environment::GetSingleton()->GetFloatValue("VspTree.maxTotalMemory", mMaxTotalMemory);
360        Environment::GetSingleton()->GetFloatValue("VspTree.maxStaticMemory", mMaxMemory);
[1002]361
[1016]362        Environment::GetSingleton()->GetBoolValue("VspTree.useCostHeuristics", mUseCostHeuristics);
363        Environment::GetSingleton()->GetBoolValue("VspTree.simulateOctree", mCirculatingAxis);
[1002]364       
[1074]365        Environment::GetSingleton()->GetIntValue("VspTree.pvsCountMethod", mPvsCountMethod);
[1023]366
[1002]367        char subdivisionStatsLog[100];
[1016]368        Environment::GetSingleton()->GetStringValue("VspTree.subdivisionStats", subdivisionStatsLog);
[1002]369        mSubdivisionStats.open(subdivisionStatsLog);
370
[1016]371        Environment::GetSingleton()->GetFloatValue("VspTree.Construction.minBand", mMinBand);
372        Environment::GetSingleton()->GetFloatValue("VspTree.Construction.maxBand", mMaxBand);
[1006]373       
[1002]374
375        //-- debug output
376
[1076]377        Debug << "******* VSP options ******** " << endl;
[1106]378
[1002]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;
[1121]395        Debug << "render cost decrease weight: " << mRenderCostDecreaseWeight << endl;
396
[1002]397        Debug << "circulating axis: " << mCirculatingAxis << endl;
398        Debug << "minband: " << mMinBand << endl;
399        Debug << "maxband: " << mMaxBand << endl;
400
[1106]401        if (mPvsCountMethod == 0)
402                Debug << "pvs count method: per object" << endl;
403        else
404                Debug << "pvs count method: per kd node" << endl;
[1084]405
[1002]406        mSplitCandidates = new vector<SortableEntry>;
407
408        Debug << endl;
409}
410
411
[1016]412VspViewCell *VspTree::GetOutOfBoundsCell()
[1002]413{
414        return mOutOfBoundsCell;
415}
416
417
[1016]418VspViewCell *VspTree::GetOrCreateOutOfBoundsCell()
[1002]419{
420        if (!mOutOfBoundsCell)
421        {
[1008]422                mOutOfBoundsCell = new VspViewCell();
[1002]423                mOutOfBoundsCell->SetId(-1);
424                mOutOfBoundsCell->SetValid(false);
425        }
426
427        return mOutOfBoundsCell;
428}
429
430
[1016]431const VspTreeStatistics &VspTree::GetStatistics() const
[1002]432{
[1008]433        return mVspStats;
[1002]434}
435
436
[1016]437VspTree::~VspTree()
[1002]438{
439        DEL_PTR(mRoot);
440        DEL_PTR(mSplitCandidates);
441}
442
[1008]443
[1020]444void VspTree::PrepareConstruction(const VssRayContainer &sampleRays,
445                                                                  AxisAlignedBox3 *forcedBoundingBox)
[1002]446{
[1099]447        // store pointer to this tree
448        VspSplitCandidate::sVspTree = this;
449
[1008]450        mVspStats.nodes = 1;
[1027]451       
[1002]452        if (forcedBoundingBox)
[1027]453        {
[1020]454                mBoundingBox = *forcedBoundingBox;
[1027]455        }
456        else // compute vsp tree bounding box
[1002]457        {
[1027]458                mBoundingBox.Initialize();
[1002]459
[1027]460                VssRayContainer::const_iterator rit, rit_end = sampleRays.end();
461
462                //-- compute bounding box
463        for (rit = sampleRays.begin(); rit != rit_end; ++ rit)
[1002]464                {
[1027]465                        VssRay *ray = *rit;
466
467                        // compute bounding box of view space
[1020]468                        mBoundingBox.Include(ray->GetTermination());
469                        mBoundingBox.Include(ray->GetOrigin());
[1002]470                }
[1027]471        }
[1108]472
473        mTermMinProbability *= mBoundingBox.GetVolume();
474        mGlobalCostMisses = 0;
[1020]475}
[1002]476
477
[1020]478void VspTree::AddSubdivisionStats(const int viewCells,
479                                                                  const float renderCostDecr,
480                                                                  const float splitCandidateCost,
481                                                                  const float totalRenderCost,
[1027]482                                                                  const float avgRenderCost)
[1020]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;
[1002]490}
491
492
493// TODO: return memory usage in MB
[1016]494float VspTree::GetMemUsage() const
[1002]495{
496        return (float)
[1016]497                 (sizeof(VspTree) +
[1008]498                  mVspStats.Leaves() * sizeof(VspLeaf) +
[1010]499                  mCreatedViewCells * sizeof(VspViewCell) +
[1008]500                  mVspStats.pvs * sizeof(ObjectPvsData) +
501                  mVspStats.Interior() * sizeof(VspInterior) +
502                  mVspStats.accumRays * sizeof(RayInfo)) / (1024.0f * 1024.0f);
[1002]503}
504
505
[1016]506bool VspTree::LocalTerminationCriteriaMet(const VspTraversalData &data) const
[1002]507{
[1106]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               
[1002]516}
517
518
[1016]519bool VspTree::GlobalTerminationCriteriaMet(const VspTraversalData &data) const
[1002]520{
[1108]521        Debug << "cost misses " << mGlobalCostMisses << " " << mTermGlobalCostMissTolerance << endl;
522        Debug << "leaves " << mVspStats.Leaves() << " " <<  mMaxViewCells << endl;
523
[1106]524        return (
525//              mOutOfMemory ||
[1011]526                (mVspStats.Leaves() >= mMaxViewCells) ||
[1106]527        (mGlobalCostMisses >= mTermGlobalCostMissTolerance)
528                );             
[1002]529}
530
531
[1108]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();
[1109]555
[1108]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;
[1109]569        mVspStats.sampleContributions += (int)sampCon;
[1108]570
571        // finally evaluate statistics for this leaf
572        EvaluateLeafStats(tData);
573}
574
575
[1016]576VspNode *VspTree::Subdivide(SplitQueue &tQueue,
[1020]577                                                        VspSplitCandidate &splitCandidate,
578                                                        const bool globalCriteriaMet)
[1002]579{
[1016]580        VspTraversalData &tData = splitCandidate.mParentData;
[1002]581
[1008]582        VspNode *newNode = tData.mNode;
[1002]583
[1020]584        if (!LocalTerminationCriteriaMet(tData) && !globalCriteriaMet)
[1002]585        {       
[1016]586                VspTraversalData tFrontData;
587                VspTraversalData tBackData;
[1002]588
589                //-- continue subdivision
590               
591                // create new interior node and two leaf node
[1008]592                const AxisAlignedPlane splitPlane = splitCandidate.mSplitPlane;
[1006]593                newNode = SubdivideNode(splitPlane, tData, tFrontData, tBackData);
[1002]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                       
[1121]602       
[1002]603                if (1)
604                {
[1108]605                        //-- subdivision statistics
606
[1020]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();
[1002]613
614                        mTotalCost += costDecr;
615                        mTotalPvsSize += tFrontData.mPvs + tBackData.mPvs - tData.mPvs;
616
[1020]617                        AddSubdivisionStats(mVspStats.Leaves(),
618                                                                -costDecr,
619                                                                splitCandidate.GetPriority(),
620                                                                mTotalCost,
621                                                                (float)mTotalPvsSize / (float)mVspStats.Leaves());
[1002]622                }
623
[1121]624
[1108]625                //-- evaluate new split candidates for global greedy cost heuristics
[1099]626                VspSplitCandidate *frontCandidate = new VspSplitCandidate(tFrontData);
627                VspSplitCandidate *backCandidate = new VspSplitCandidate(tBackData);
[1002]628
[1099]629                EvalSplitCandidate(*frontCandidate);
630                EvalSplitCandidate(*backCandidate);
[1020]631
[1099]632                tQueue.Push(frontCandidate);
633                tQueue.Push(backCandidate);
[1121]634               
635                // delete old view cell
636                delete tData.mNode->GetViewCell();
[1002]637                // delete old leaf node
638                DEL_PTR(tData.mNode);
639        }
640
[1108]641        // subdivision terminated
[1002]642        if (newNode->IsLeaf())
643        {
[1108]644                //-- create new view cell
[1109]645                //CreateViewCell(tData);
[1020]646       
[1002]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();                     
[1108]654                                dynamic_cast<VspLeaf *>(newNode)->mVssRays.push_back((*it).mRay);
[1002]655                        }
656                }
657
[1108]658                // finally evaluate statistics for this leaf
[1020]659                EvaluateLeafStats(tData);
[1002]660        }
661
662        //-- cleanup
663        tData.Clear();
[1020]664       
[1002]665        return newNode;
666}
667
668
[1099]669void VspTree::EvalSplitCandidate(VspSplitCandidate &splitCandidate)
[1002]670{
[1020]671        float frontProb;
[1027]672        float backProb;
[1020]673       
[1099]674        VspLeaf *leaf = dynamic_cast<VspLeaf *>(splitCandidate.mParentData.mNode);
[1020]675       
[1002]676        // compute locally best split plane
[1008]677        const bool success =
[1099]678                SelectSplitPlane(splitCandidate.mParentData, splitCandidate.mSplitPlane, frontProb, backProb);
[1008]679
[1002]680        // compute global decrease in render cost
[1099]681        const float priority = EvalRenderCostDecrease(splitCandidate.mSplitPlane, splitCandidate.mParentData);
682        splitCandidate.SetPriority(priority);
683        splitCandidate.mMaxCostMisses = success ? splitCandidate.mParentData.mMaxCostMisses : splitCandidate.mParentData.mMaxCostMisses + 1;
[1074]684        //Debug << "p: " << tData.mNode << " depth: " << tData.mDepth << endl;
[1002]685}
686
687
[1016]688VspInterior *VspTree::SubdivideNode(const AxisAlignedPlane &splitPlane,
689                                                                        VspTraversalData &tData,
690                                                                        VspTraversalData &frontData,
691                                                                        VspTraversalData &backData)
[1002]692{
[1008]693        VspLeaf *leaf = dynamic_cast<VspLeaf *>(tData.mNode);
[1002]694       
695        //-- the front and back traversal data is filled with the new values
[1108]696
[1002]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
[1109]709       
[1074]710        //Debug << "f: " << frontData.mRays->size() << " b: " << backData.mRays->size() << "d: " << tData.mRays->size() << endl;
[1011]711        //-- compute pvs
[1002]712        frontData.mPvs = ComputePvsSize(*frontData.mRays);
713        backData.mPvs = ComputePvsSize(*backData.mRays);
714
715        // split front and back node geometry and compute area
[1006]716        tData.mBoundingBox.Split(splitPlane.mAxis, splitPlane.mPosition,
717                                                         frontData.mBoundingBox, backData.mBoundingBox);
[1002]718
719
[1006]720        frontData.mProbability = frontData.mBoundingBox.GetVolume();
721        backData.mProbability = tData.mProbability - frontData.mProbability;
[1002]722
723       
[1006]724    ///////////////////////////////////////////
[1002]725        // subdivide further
726
727        // store maximal and minimal depth
[1008]728        if (tData.mDepth > mVspStats.maxDepth)
[1002]729        {
[1008]730                Debug << "max depth increases to " << tData.mDepth << " at " << mVspStats.Leaves() << " leaves" << endl;
731                mVspStats.maxDepth = tData.mDepth;
[1002]732        }
733
[1109]734        // two more leaves
[1008]735        mVspStats.nodes += 2;
[1002]736   
[1008]737        VspInterior *interior = new VspInterior(splitPlane);
[1002]738
739#ifdef _DEBUG
740        Debug << interior << endl;
741#endif
742
743
744        //-- create front and back leaf
745
[1008]746        VspInterior *parent = leaf->GetParent();
[1002]747
748        // replace a link from node's parent
749        if (parent)
750        {
751                parent->ReplaceChildLink(leaf, interior);
752                interior->SetParent(parent);
[1121]753
754                // remove "parent" view cell from pvs of all objects (traverse trough rays)
755                RemoveParentViewCellReferences(tData.mNode->GetViewCell());
[1002]756        }
757        else // new root
758        {
759                mRoot = interior;
760        }
761
[1109]762        VspLeaf *frontLeaf = new VspLeaf(interior);
763        VspLeaf *backLeaf = new VspLeaf(interior);
764
[1002]765        // and setup child links
[1109]766        interior->SetupChildLinks(frontLeaf, backLeaf);
[1106]767       
[1016]768        // add bounding box
769        interior->SetBoundingBox(tData.mBoundingBox);
[1002]770
[1106]771        // set front and back leaf
[1109]772        frontData.mNode = frontLeaf;
773        backData.mNode = backLeaf;
[1106]774
[1109]775        // explicitely create front and back view cell
776        CreateViewCell(frontData);
777        CreateViewCell(backData);
778
[1121]779       
[1109]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
[1110]783        AddViewCellReferences(frontLeaf->GetViewCell());
784        AddViewCellReferences(backLeaf->GetViewCell());
[1109]785
[1002]786        interior->mTimeStamp = mTimeStamp ++;
787       
788        return interior;
789}
790
791
[1109]792/*void VspTree::RemoveParentViewCellReferences(VspTraversalData &parentData)
793{
794        RayInfoContainer::const_iterator it, it_end = parentData.mRays->end();
[1089]795
[1109]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)
[1089]805        {
[1109]806                VssRay *ray = (*it).mRay;
[1089]807
[1109]808                Intersectable *obj = ray->mTerminationObject;
809
810                if (obj && !obj->Mailed())
[1089]811                {
[1109]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                        }
[1089]828                }
[1109]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                }
[1089]851        }
[1109]852}*/
[1089]853
[1109]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)
[1089]863        {
[1109]864                Intersectable *object = (*oit).first;
865                // HACK: make sure that the view cell is removed from the pvs
866                const float high_contri = 9999999;
[1089]867
[1109]868                // remove reference count of view cells
[1110]869                object->mViewCellPvs.RemoveSample(parent, high_contri);
[1089]870        }
[1109]871}
[1089]872
[1109]873
874void VspTree::AddViewCellReferences(ViewCell *vc) const
875{
876        KdLeaf::NewMail();
877
878        // Add front view cell to the object pvsss
[1110]879        ObjectPvsMap::const_iterator oit, oit_end = vc->GetPvs().mEntries.end();
[1109]880
[1110]881        for (oit = vc->GetPvs().mEntries.begin(); oit != oit_end; ++ oit)
[1089]882        {
[1109]883                Intersectable *object = (*oit).first;
[1089]884
[1109]885                // increase reference count of view cells
[1110]886                object->mViewCellPvs.AddSample(vc, 1);
[1089]887        }
888}
889
890
[1016]891void VspTree::AddToPvs(VspLeaf *leaf,
[1089]892                                           const RayInfoContainer &rays,
893                                           float &sampleContributions,
894                                           int &contributingSamples)
[1002]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
[1109]912                Intersectable *obj = ray->mTerminationObject;
913
914                if (obj)
[1002]915                {
[1109]916                        if (vc->AddPvsSample(obj, ray->mPdf, contribution))
[1089]917                        {
[1002]918                                madeContrib = true;
[1089]919                        }
920
[1002]921                        sc += contribution;
922                }
[1109]923
[1110]924                obj = ray->mOriginObject;
[1109]925
926                if (obj)
[1002]927                {
[1109]928                        if (vc->AddPvsSample(obj, ray->mPdf, contribution))
[1089]929                        {
[1002]930                                madeContrib = true;
[1089]931                        }
932
[1002]933                        sc += contribution;
934                }
935
[1109]936
[1047]937                if (madeContrib)
938                        ++ contributingSamples;
[1002]939               
[1077]940                // store rays for visualization
941                if (0) leaf->mVssRays.push_back(new VssRay(*ray));
[1002]942        }
943}
944
945
[1047]946void VspTree::SortSplitCandidates(const RayInfoContainer &rays,
947                                                                  const int axis,
948                                                                  float minBand,
949                                                                  float maxBand)
[1002]950{
951        mSplitCandidates->clear();
952
953        int requestedSize = 2 * (int)(rays.size());
[1047]954
[1002]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
[1084]967        //-- insert all queries
[1002]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
[1072]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
[1074]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)
[1072]1062{       
[1074]1063        Intersectable::NewMail();
1064        KdNode::NewMail();
[1072]1065
[1074]1066        int pvsSize = 0;
[1072]1067
[1074]1068        RayInfoContainer::const_iterator ri, ri_end = rays.end();
[1072]1069
[1074]1070        //-- set all kd nodes as belonging to the front pvs
1071
1072        for (ri = rays.begin(); ri != ri_end; ++ ri)
[1072]1073        {
[1074]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                }
[1072]1119        }
[1074]1120
1121        return pvsSize;
[1072]1122}
1123
1124
[1074]1125void VspTree::RemoveContriFromPvs(KdLeaf *leaf, int &pvs) const
[1072]1126{
[1074]1127        // leaf falls out of right pvs
1128        if (-- leaf->mCounter == 0)
[1072]1129        {
[1074]1130                pvs -= ((int)leaf->mObjects.size() - (int)leaf->mMultipleObjects.size());
1131        }
[1072]1132
[1074]1133        //-- handle objects which are in several kd leaves separately
1134        ObjectContainer::const_iterator oit, oit_end = leaf->mMultipleObjects.end();
[1072]1135
[1074]1136        for (oit = leaf->mMultipleObjects.begin(); oit != oit_end; ++ oit)
1137        {
1138                Intersectable *object = *oit;
[1072]1139
[1074]1140                if (-- object->mCounter == 0)
1141                {
1142                        -- pvs;
1143                }
[1072]1144        }
1145}
1146
1147
[1074]1148void VspTree::AddContriToPvs(KdLeaf *leaf, int &pvs) const
[1072]1149{
[1121]1150        if (leaf->Mailed())
1151                return;
1152       
1153        leaf->Mail();
[1072]1154
[1121]1155        // add objects without those which are part of several kd leaves
1156        pvs += ((int)leaf->mObjects.size() - (int)leaf->mMultipleObjects.size());
[1072]1157
[1121]1158        //-- handle objects of several kd leaves separately
1159        ObjectContainer::const_iterator oit, oit_end = leaf->mMultipleObjects.end();
[1072]1160
[1121]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())
[1074]1167                {
[1121]1168                        object->Mail();
1169                        ++ pvs;
[1074]1170                }
1171        }                                       
1172}
[1072]1173
1174
[1074]1175void VspTree::EvalPvsIncr(const SortableEntry &ci,
1176                                                  int &pvsLeft,
1177                                                  int &pvsRight) const
1178{
1179        VssRay *ray = ci.ray;
[1072]1180
[1074]1181        Intersectable *oObject = ray->mOriginObject;
1182        Intersectable *tObject = ray->mTerminationObject;
[1072]1183
[1074]1184        if (oObject)
1185        {       
1186                if (mPvsCountMethod == PER_OBJECT)
[1072]1187                {
[1074]1188                        if (ci.type == SortableEntry::ERayMin)
[1072]1189                        {
[1074]1190                                if (!oObject->Mailed())
[1072]1191                                {
[1074]1192                                        oObject->Mail();
1193                                        ++ pvsLeft;
[1072]1194                                }
[1074]1195                        }
1196                        else if (ci.type == SortableEntry::ERayMax)
1197                        {
1198                                if (-- oObject->mCounter == 0)
1199                                        -- pvsRight;
1200                        }
[1072]1201                }
[1074]1202                else // per kd node
[1072]1203                {
[1074]1204                        set<KdLeaf *>::const_iterator kit, kit_end = oObject->mKdLeaves.end();
1205
1206                        for (kit = oObject->mKdLeaves.begin(); kit != kit_end; ++ kit)
[1072]1207                        {
[1074]1208                                KdLeaf *node = dynamic_cast<KdLeaf *>(*kit);
[1072]1209
[1074]1210                                // add contributions of the kd nodes
1211                                if (ci.type == SortableEntry::ERayMin)
[1072]1212                                {
[1074]1213                                        AddContriToPvs(node, pvsLeft);
[1072]1214                                }
[1074]1215                                else if (ci.type == SortableEntry::ERayMax)
[1072]1216                                {
[1074]1217                                        RemoveContriFromPvs(node, pvsRight);
[1072]1218                                }
[1074]1219                        }
[1072]1220                }
[1074]1221                       
[1072]1222
[1074]1223        }
1224       
1225        if (tObject)
1226        {       
1227                if (mPvsCountMethod == PER_OBJECT)
1228                {
1229                        if (ci.type == SortableEntry::ERayMin)
[1072]1230                        {
[1074]1231                                if (!tObject->Mailed())
1232                                {
1233                                        tObject->Mail();
1234                                        ++ pvsLeft;
1235                                }
[1072]1236                        }
[1074]1237                        else if (ci.type == SortableEntry::ERayMax)
[1072]1238                        {
[1074]1239                                if (-- tObject->mCounter == 0)
1240                                        -- pvsRight;
[1072]1241                        }
1242                }
[1074]1243                else // per kd node
[1072]1244                {
[1074]1245                        set<KdLeaf *>::const_iterator kit, kit_end = tObject->mKdLeaves.end();
[1077]1246                       
[1074]1247                        for (kit = tObject->mKdLeaves.begin(); kit != kit_end; ++ kit)
[1072]1248                        {
[1074]1249                                KdLeaf *node = dynamic_cast<KdLeaf *>(*kit);
[1072]1250
[1074]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                                }
[1072]1259                        }
1260                }
1261        }
1262}
1263
1264
[1074]1265float VspTree::EvalLocalCostHeuristics(const RayInfoContainer &rays,
1266                                                                           const AxisAlignedBox3 &box,
1267                                                                           int pvsSize,
1268                                                                           const int axis,
1269                                                                           float &position)
[1072]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
[1077]1281        // prepare the sweep
[1074]1282        // (note: returns pvs size, so there would be no need
1283        // to give pvs size as argument)
1284        pvsSize = PrepareHeuristics(rays);
1285
[1072]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
[1074]1304        float ratio = 99999999.0f;
[1072]1305        bool splitPlaneFound = false;
1306
1307        Intersectable::NewMail();
[1074]1308        KdLeaf::NewMail();
[1072]1309
1310
[1074]1311        vector<SortableEntry>::const_iterator ci, ci_end = mSplitCandidates->end();
[1072]1312
[1074]1313        //-- traverse through visibility events
1314
1315        for (ci = mSplitCandidates->begin(); ci != ci_end; ++ ci)
[1072]1316        {
[1074]1317                EvalPvsIncr(*ci, pvsl, pvsr);
1318
[1072]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
[1074]1324                        //Debug  << "pos=" << (*ci).value << "\t pvs=(" <<  pvsl << "," << pvsr << ")" << "\t cost= " << sum << endl;
[1072]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       
[1074]1339       
[1072]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        }
[1074]1360       
[1072]1361        //if (axis != 1)
[1074]1362        Debug << "axis=" << axis << " costRatio=" << ratio << " pos=" << position << " t=" << (position - minBox) / (maxBox - minBox)
1363              <<"\t pb=(" << pvsBack << ")\t pf=(" << pvsFront << ")" << endl;
[1072]1364
1365        return ratio;
[1077]1366}
[1072]1367
1368
[1027]1369float VspTree::SelectSplitPlane(const VspTraversalData &tData,
1370                                                                AxisAlignedPlane &plane,
1371                                                                float &pFront,
1372                                                                float &pBack)
[1002]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
[1006]1380        AxisAlignedBox3 box = tData.mBoundingBox;
[1002]1381               
1382        int sAxis = 0;
[1006]1383        int bestAxis = -1;
[1002]1384
1385        // if we use some kind of specialised fixed axis
1386    const bool useSpecialAxis =
[1023]1387                mOnlyDrivingAxis || mCirculatingAxis;
[1074]1388        //Debug << "data: " << tData.mBoundingBox << " pvs " << tData.mPvs << endl;
[1047]1389        if (mCirculatingAxis)
1390        {
1391                int parentAxis = 0;
1392                VspNode *parent = tData.mNode->GetParent();
[1002]1393
[1047]1394                if (parent)
1395                        parentAxis = dynamic_cast<VspInterior *>(parent)->GetAxis();
[1027]1396
1397                sAxis = (parentAxis + 1) % 3;
[1047]1398        }
[1002]1399        else if (mOnlyDrivingAxis)
[1047]1400        {
[1002]1401                sAxis = box.Size().DrivingAxis();
[1047]1402        }
[1106]1403       
[1047]1404        for (int axis = 0; axis < 3; ++ axis)
[1002]1405        {
1406                if (!useSpecialAxis || (axis == sAxis))
1407                {
1408                        if (mUseCostHeuristics)
1409                        {
[1106]1410                                //-- place split plane using heuristics
[1002]1411                                nCostRatio[axis] =
[1027]1412                                        EvalLocalCostHeuristics(*tData.mRays,
1413                                                                                        box,
[1002]1414                                                                                        tData.mPvs,
1415                                                                                        axis,
1416                                                                                        nPosition[axis]);                       
1417                        }
[1106]1418                        else
[1002]1419                        {
[1106]1420                                //-- split plane position is spatial median
[1002]1421                                nPosition[axis] = (box.Min()[axis] + box.Max()[axis]) * 0.5f;
[1027]1422                                nCostRatio[axis] = EvalLocalSplitCost(tData,
1423                                                                                                          box,
1424                                                                                                          axis,
1425                                                                                                          nPosition[axis],
1426                                                                                                          nProbFront[axis],
1427                                                                                                          nProbBack[axis]);
[1002]1428                        }
1429                                               
[1006]1430                        if (bestAxis == -1)
[1002]1431                        {
[1006]1432                                bestAxis = axis;
[1002]1433                        }
[1006]1434                        else if (nCostRatio[axis] < nCostRatio[bestAxis])
[1002]1435                        {
[1006]1436                                bestAxis = axis;
[1002]1437                        }
1438                }
1439        }
1440
[1047]1441
[1108]1442        //-- assign values of best split
[1047]1443       
[1006]1444        plane.mAxis = bestAxis;
[1108]1445        plane.mPosition = nPosition[bestAxis]; // split plane position
[1047]1446
[1002]1447        pFront = nProbFront[bestAxis];
1448        pBack = nProbBack[bestAxis];
1449
[1074]1450        //Debug << "val: " << nCostRatio[bestAxis] << " axis: " << bestAxis << endl;
[1002]1451        return nCostRatio[bestAxis];
1452}
1453
1454
[1016]1455float VspTree::EvalRenderCostDecrease(const AxisAlignedPlane &candidatePlane,
1456                                                                          const VspTraversalData &data) const
[1002]1457{
[1047]1458#if 0
1459        return (float)-data.mDepth;
[1027]1460#endif
[1002]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
[1089]1472        Intersectable::NewMail();
[1002]1473       
[1016]1474        RayInfoContainer::const_iterator rit, rit_end = data.mRays->end();
[1015]1475
[1016]1476        for (rit = data.mRays->begin(); rit != rit_end; ++ rit)
[1002]1477        {
[1015]1478                RayInfo rayInf = *rit;
[1002]1479
1480                float t;
1481                VssRay *ray = rayInf.mRay;
[1027]1482
[1008]1483                const int cf =
1484                        rayInf.ComputeRayIntersection(candidatePlane.mAxis,
1485                                                                                  candidatePlane.mPosition, t);
[1002]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
[1027]1492
[1011]1493        AxisAlignedBox3 frontBox;
1494        AxisAlignedBox3 backBox;
1495
1496        data.mBoundingBox.Split(candidatePlane.mAxis, candidatePlane.mPosition, frontBox, backBox);
1497
1498        pFront = frontBox.GetVolume();
[1006]1499        pBack = pOverall - pFront;
1500               
1501
[1020]1502        //-- pvs rendering heuristics
[1002]1503        const int lowerPvsLimit = mViewCellsManager->GetMinPvsSize();
1504        const int upperPvsLimit = mViewCellsManager->GetMaxPvsSize();
1505
[1020]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);
[1002]1510                       
1511        const float oldRenderCost = pOverall * penaltyOld;
1512        const float newRenderCost = penaltyFront * pFront + penaltyBack * pBack;
1513
1514        //Debug << "decrease: " << oldRenderCost - newRenderCost << endl;
[1020]1515        const float renderCostDecrease = (oldRenderCost - newRenderCost) / mBoundingBox.GetVolume();
[1002]1516       
[1020]1517        // take render cost of node into account
1518        // otherwise danger of being stuck in a local minimum!!
[1121]1519        const float factor = mRenderCostDecreaseWeight;
[1027]1520
[1020]1521        const float normalizedOldRenderCost = oldRenderCost / mBoundingBox.GetVolume();
[1027]1522        return factor * renderCostDecrease + (1.0f - factor) * normalizedOldRenderCost;
[1002]1523}
1524
1525
[1027]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
[1002]1532{
1533        float pvsTotal = 0;
1534        float pvsFront = 0;
1535        float pvsBack = 0;
1536       
1537        // create unique ids for pvs heuristics
[1089]1538        Intersectable::NewMail();
[1002]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
[1108]1555        //-- cost heuristics
[1002]1556        float pOverall;
[1027]1557       
[1002]1558        pOverall = data.mProbability;
[1108]1559
1560        // we take simplified computation for spatial mid split
[1006]1561        pBack = pFront = pOverall * 0.5f;
1562       
[1108]1563        const float newCost = pvsBack * pBack + pvsFront * pFront;
1564        const float oldCost = (float)pvsSize * pOverall + Limits::Small;
[1006]1565       
[1002]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
[1016]1575void VspTree::AddObjToPvs(Intersectable *obj,
[1027]1576                                                  const int cf,
1577                                                  float &frontPvs,
1578                                                  float &backPvs,
1579                                                  float &totalPvs) const
[1002]1580{
1581        if (!obj)
1582                return;
1583
[1076]1584        //const float renderCost = mViewCellsManager->EvalRenderCost(obj);
1585        const int renderCost = 1;
[1002]1586
[1089]1587        // object in no pvs => new
1588        if (!obj->Mailed() && !obj->Mailed(1) && !obj->Mailed(2))
[1002]1589        {
1590                totalPvs += renderCost;
1591        }
1592
1593        // TODO: does this really belong to no pvs?
1594        //if (cf == Ray::COINCIDENT) return;
1595
[1089]1596        if (cf >= 0) // front pvs
[1002]1597        {
[1089]1598                if (!obj->Mailed() && !obj->Mailed(2))
[1002]1599                {
1600                        frontPvs += renderCost;
1601               
[1089]1602                        // already in back pvs => in both pvss
1603                        if (obj->Mailed(1))
1604                                obj->Mail(2);
[1002]1605                        else
[1089]1606                                obj->Mail();
[1002]1607                }
1608        }
1609
[1089]1610        if (cf <= 0) // back pvs
[1002]1611        {
[1089]1612                if (!obj->Mailed(1) && !obj->Mailed(2))
[1002]1613                {
1614                        backPvs += renderCost;
1615               
[1089]1616                        // already in front pvs => in both pvss
1617                        if (obj->Mailed())
1618                                obj->Mail(2);
[1002]1619                        else
[1089]1620                                obj->Mail(1);
[1002]1621                }
1622        }
1623}
1624
1625
[1016]1626void VspTree::CollectLeaves(vector<VspLeaf *> &leaves,
[1002]1627                                                           const bool onlyUnmailed,
1628                                                           const int maxPvsSize) const
1629{
[1008]1630        stack<VspNode *> nodeStack;
[1002]1631        nodeStack.push(mRoot);
1632
1633        while (!nodeStack.empty())
1634        {
[1008]1635                VspNode *node = nodeStack.top();
[1002]1636                nodeStack.pop();
1637               
1638                if (node->IsLeaf())
1639                {
1640                        // test if this leaf is in valid view space
[1008]1641                        VspLeaf *leaf = dynamic_cast<VspLeaf *>(node);
[1002]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                {
[1008]1651                        VspInterior *interior = dynamic_cast<VspInterior *>(node);
[1002]1652
1653                        nodeStack.push(interior->GetBack());
1654                        nodeStack.push(interior->GetFront());
1655                }
1656        }
1657}
1658
1659
[1016]1660AxisAlignedBox3 VspTree::GetBoundingBox() const
[1002]1661{
[1020]1662        return mBoundingBox;
[1002]1663}
1664
1665
[1016]1666VspNode *VspTree::GetRoot() const
[1002]1667{
1668        return mRoot;
1669}
1670
1671
[1016]1672void VspTree::EvaluateLeafStats(const VspTraversalData &data)
[1002]1673{
1674        // the node became a leaf -> evaluate stats for leafs
[1008]1675        VspLeaf *leaf = dynamic_cast<VspLeaf *>(data.mNode);
[1002]1676
1677
[1008]1678        if (data.mPvs > mVspStats.maxPvs)
[1002]1679        {
[1008]1680                mVspStats.maxPvs = data.mPvs;
[1002]1681        }
1682
[1008]1683        mVspStats.pvs += data.mPvs;
[1002]1684
[1008]1685        if (data.mDepth < mVspStats.minDepth)
[1002]1686        {
[1008]1687                mVspStats.minDepth = data.mDepth;
[1002]1688        }
1689       
1690        if (data.mDepth >= mTermMaxDepth)
1691        {
[1008]1692        ++ mVspStats.maxDepthNodes;
1693                //Debug << "new max depth: " << mVspStats.maxDepthNodes << endl;
[1002]1694        }
1695
1696        // accumulate rays to compute rays /  leaf
[1008]1697        mVspStats.accumRays += (int)data.mRays->size();
[1002]1698
1699        if (data.mPvs < mTermMinPvs)
[1008]1700                ++ mVspStats.minPvsNodes;
[1002]1701
1702        if ((int)data.mRays->size() < mTermMinRays)
[1008]1703                ++ mVspStats.minRaysNodes;
[1002]1704
1705        if (data.GetAvgRayContribution() > mTermMaxRayContribution)
[1008]1706                ++ mVspStats.maxRayContribNodes;
[1002]1707
1708        if (data.mProbability <= mTermMinProbability)
[1008]1709                ++ mVspStats.minProbabilityNodes;
[1002]1710       
1711        // accumulate depth to compute average depth
[1008]1712        mVspStats.accumDepth += data.mDepth;
[1002]1713
1714        ++ mCreatedViewCells;
1715
[1074]1716#ifdef _DEBUG
[1002]1717        Debug << "BSP stats: "
1718                  << "Depth: " << data.mDepth << " (max: " << mTermMaxDepth << "), "
1719                  << "PVS: " << data.mPvs << " (min: " << mTermMinPvs << "), "
1720                  << "#rays: " << (int)data.mRays->size() << " (max: " << mTermMinRays << "), "
[1027]1721                  << "#pvs: " << leaf->GetViewCell()->GetPvs().GetSize() << "), "
[1002]1722                  << "#avg ray contrib (pvs): " << (float)data.mPvs / (float)data.mRays->size() << endl;
[1074]1723#endif
[1002]1724}
1725
1726
[1016]1727void VspTree::CollectViewCells(ViewCellContainer &viewCells, bool onlyValid) const
[1002]1728{
1729        ViewCell::NewMail();
1730        CollectViewCells(mRoot, onlyValid, viewCells, true);
1731}
1732
1733
[1016]1734void VspTree::CollapseViewCells()
[1002]1735{
[1111]1736// TODO matt
[1002]1737#if HAS_TO_BE_REDONE
[1008]1738        stack<VspNode *> nodeStack;
[1002]1739
1740        if (!mRoot)
1741                return;
1742
1743        nodeStack.push(mRoot);
1744       
1745        while (!nodeStack.empty())
1746        {
[1008]1747                VspNode *node = nodeStack.top();
[1002]1748                nodeStack.pop();
1749               
1750                if (node->IsLeaf())
1751        {
[1008]1752                        BspViewCell *viewCell = dynamic_cast<VspLeaf *>(node)->GetViewCell();
[1002]1753
1754                        if (!viewCell->GetValid())
1755                        {
[1008]1756                                BspViewCell *viewCell = dynamic_cast<VspLeaf *>(node)->GetViewCell();
[1002]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                                {
[1008]1765                                        VspLeaf *l = dynamic_cast<BspViewCell *>(*it)->mLeaf;
[1002]1766                                        l->SetViewCell(GetOrCreateOutOfBoundsCell());
[1008]1767                                        ++ mVspStats.invalidLeaves;
[1002]1768                                }
1769
1770                                // add to unbounded view cell
1771                                GetOrCreateOutOfBoundsCell()->GetPvs().AddPvs(viewCell->GetPvs());
1772                                DEL_PTR(viewCell);
1773                        }
1774                }
1775                else
1776                {
[1008]1777                        VspInterior *interior = dynamic_cast<VspInterior *>(node);
[1002]1778               
1779                        nodeStack.push(interior->GetFront());
1780                        nodeStack.push(interior->GetBack());
1781                }
1782        }
1783
[1008]1784        Debug << "invalid leaves: " << mVspStats.invalidLeaves << endl;
[1002]1785#endif
1786}
1787
1788
[1016]1789void VspTree::CollectRays(VssRayContainer &rays)
[1002]1790{
[1008]1791        vector<VspLeaf *> leaves;
[1002]1792
[1008]1793        vector<VspLeaf *>::const_iterator lit, lit_end = leaves.end();
[1002]1794
1795        for (lit = leaves.begin(); lit != lit_end; ++ lit)
1796        {
[1008]1797                VspLeaf *leaf = *lit;
[1002]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
[1021]1806void VspTree::SetViewCellsManager(ViewCellsManager *vcm)
1807{
1808        mViewCellsManager = vcm;
1809}
1810
1811
[1016]1812void VspTree::ValidateTree()
[1002]1813{
[1020]1814        mVspStats.invalidLeaves = 0;
1815
[1008]1816        stack<VspNode *> nodeStack;
[1002]1817
1818        if (!mRoot)
1819                return;
1820
1821        nodeStack.push(mRoot);
[1020]1822
[1002]1823        while (!nodeStack.empty())
1824        {
[1008]1825                VspNode *node = nodeStack.top();
[1002]1826                nodeStack.pop();
1827               
1828                if (node->IsLeaf())
1829                {
[1008]1830                        VspLeaf *leaf = dynamic_cast<VspLeaf *>(node);
[1002]1831
1832                        if (!leaf->GetViewCell()->GetValid())
[1008]1833                                ++ mVspStats.invalidLeaves;
[1002]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                {
[1008]1844                        VspInterior *interior = dynamic_cast<VspInterior *>(node);
[1002]1845               
1846                        nodeStack.push(interior->GetFront());
1847                        nodeStack.push(interior->GetBack());
1848                }
1849        }
1850
[1008]1851        Debug << "invalid leaves: " << mVspStats.invalidLeaves << endl;
[1002]1852}
1853
1854
1855
[1016]1856void VspTree::CollectViewCells(VspNode *root,
[1002]1857                                                                  bool onlyValid,
1858                                                                  ViewCellContainer &viewCells,
1859                                                                  bool onlyUnmailed) const
1860{
[1008]1861        stack<VspNode *> nodeStack;
[1002]1862
1863        if (!root)
1864                return;
1865
1866        nodeStack.push(root);
1867       
1868        while (!nodeStack.empty())
1869        {
[1008]1870                VspNode *node = nodeStack.top();
[1002]1871                nodeStack.pop();
1872               
1873                if (node->IsLeaf())
1874                {
1875                        if (!onlyValid || node->TreeValid())
1876                        {
[1008]1877                                ViewCellLeaf *leafVc = dynamic_cast<VspLeaf *>(node)->GetViewCell();
[1002]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                {
[1008]1890                        VspInterior *interior = dynamic_cast<VspInterior *>(node);
[1002]1891               
1892                        nodeStack.push(interior->GetFront());
1893                        nodeStack.push(interior->GetBack());
1894                }
1895        }
1896}
1897
1898
[1016]1899int VspTree::FindNeighbors(VspLeaf *n,
[1020]1900                                                   vector<VspLeaf *> &neighbors,
1901                                                   const bool onlyUnmailed) const
[1002]1902{
[1012]1903        stack<VspNode *> nodeStack;
1904        nodeStack.push(mRoot);
[1002]1905
[1012]1906        const AxisAlignedBox3 box = GetBBox(n);
[1002]1907
1908        while (!nodeStack.empty())
1909        {
[1008]1910                VspNode *node = nodeStack.top();
[1002]1911                nodeStack.pop();
1912
1913                if (node->IsLeaf())
1914                {
[1012]1915                        VspLeaf *leaf = dynamic_cast<VspLeaf *>(node);
1916
1917                        if (leaf != n && (!onlyUnmailed || !leaf->Mailed()))
1918                                neighbors.push_back(leaf);
[1002]1919                }
1920                else
1921                {
[1008]1922                        VspInterior *interior = dynamic_cast<VspInterior *>(node);
[1012]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        }
[1002]1939
[1012]1940        return (int)neighbors.size();
1941}
[1002]1942
1943
[1012]1944// Find random neighbor which was not mailed
[1016]1945VspLeaf *VspTree::GetRandomLeaf(const Plane3 &plane)
[1012]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)
[1027]1968                        {
[1002]1969                                next = interior->GetFront();
[1027]1970                        }
[1012]1971            else
[1027]1972                        {
[1012]1973                                if (GetBBox(interior->GetFront()).Side(plane) < 0)
[1027]1974                                {
[1002]1975                                        next = interior->GetBack();
[1027]1976                                }
[1012]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                                }
[1027]1986                        }
1987                       
1988                        nodeStack.push(next);
[1002]1989                }
1990        }
[1012]1991 
[1002]1992        return NULL;
1993}
1994
1995
[1016]1996VspLeaf *VspTree::GetRandomLeaf(const bool onlyUnmailed)
[1002]1997{
[1008]1998        stack<VspNode *> nodeStack;
[1002]1999
2000        nodeStack.push(mRoot);
2001
2002        int mask = rand();
2003
2004        while (!nodeStack.empty())
2005        {
[1008]2006                VspNode *node = nodeStack.top();
[1002]2007                nodeStack.pop();
2008
2009                if (node->IsLeaf())
2010                {
2011                        if ( (!onlyUnmailed || !node->Mailed()) )
[1008]2012                                return dynamic_cast<VspLeaf *>(node);
[1002]2013                }
2014                else
2015                {
[1008]2016                        VspInterior *interior = dynamic_cast<VspInterior *>(node);
[1002]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
[1072]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
[1016]2069int VspTree::ComputePvsSize(const RayInfoContainer &rays) const
[1002]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
[1016]2103float VspTree::GetEpsilon() const
[1002]2104{
2105        return mEpsilon;
2106}
2107
2108
[1016]2109int VspTree::CastLineSegment(const Vector3 &origin,
[1027]2110                                                         const Vector3 &termination,
2111                             ViewCellContainer &viewcells)
[1002]2112{
2113        int hits = 0;
2114
2115        float mint = 0.0f, maxt = 1.0f;
[1012]2116        const Vector3 dir = termination - origin;
[1002]2117
[1012]2118        stack<LineTraversalData> tStack;
2119
[1002]2120        Intersectable::NewMail();
2121        ViewCell::NewMail();
2122
2123        Vector3 entp = origin;
2124        Vector3 extp = termination;
2125
[1008]2126        VspNode *node = mRoot;
[1012]2127        VspNode *farChild;
[1002]2128
[1012]2129        float position;
2130        int axis;
2131
[1002]2132        while (1)
2133        {
2134                if (!node->IsLeaf())
2135                {
[1008]2136                        VspInterior *in = dynamic_cast<VspInterior *>(node);
[1012]2137                        position = in->GetPosition();
2138                        axis = in->GetAxis();
[1002]2139
[1012]2140                        if (entp[axis] <= position)
[1002]2141                        {
[1012]2142                                if (extp[axis] <= position)
2143                                {
[1047]2144                                        node = in->GetBack();
[1012]2145                                        // cases N1,N2,N3,P5,Z2,Z3
[1002]2146                                        continue;
[1012]2147                                } else
2148                                {
2149                                        // case N4
[1047]2150                                        node = in->GetBack();
2151                                        farChild = in->GetFront();
[1012]2152                                }
2153                        }
2154                        else
[1002]2155                        {
[1012]2156                                if (position <= extp[axis])
2157                                {
[1047]2158                                        node = in->GetFront();
[1012]2159                                        // cases P1,P2,P3,N5,Z1
[1002]2160                                        continue;
[1012]2161                                }
2162                                else
2163                                {
[1047]2164                                        node = in->GetFront();
2165                                        farChild = in->GetBack();
[1012]2166                                        // case P4
2167                                }
[1002]2168                        }
2169
[1012]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
[1002]2174
[1012]2175                        extp = origin + dir * tdist;
2176                        maxt = tdist;
2177                }
[1002]2178                else
2179                {
[1012]2180                        // compute intersection with all objects in this leaf
[1008]2181                        VspLeaf *leaf = dynamic_cast<VspLeaf *>(node);
[1012]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;
[1002]2199                       
[1012]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
[1016]2212int VspTree::CastRay(Ray &ray)
[1012]2213{
2214        int hits = 0;
2215
2216        stack<LineTraversalData> tStack;
2217        const Vector3 dir = ray.GetDir();
2218
2219        float maxt, mint;
2220
[1020]2221        if (!mBoundingBox.GetRaySegment(ray, mint, maxt))
[1012]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;
[1020]2253                                }
2254                                else
[1012]2255                                {
2256                                        // case N4
2257                                        node = in->GetBack();
2258                                        farChild = in->GetFront();
2259                                }
2260                        }
[1002]2261                        else
[1012]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                        }
[1002]2276
[1012]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())
[1002]2291                        {
[1012]2292                                vc->Mail();
2293                                // todo: add view cells to ray
[1002]2294                                ++ hits;
2295                        }
[1012]2296#if 0
[1020]2297                        leaf->mRays.push_back(RayInfo(new VssRay(origin, termination, NULL, NULL, 0)));
[1012]2298#endif
2299                        // get the next node from the stack
[1002]2300                        if (tStack.empty())
2301                                break;
2302
2303                        entp = extp;
[1012]2304                        mint = maxt;
[1002]2305                       
[1012]2306                        LineTraversalData &s  = tStack.top();
[1002]2307                        node = s.mNode;
2308                        extp = s.mExitPoint;
[1012]2309                        maxt = s.mMaxT;
[1002]2310                        tStack.pop();
2311                }
2312        }
[1012]2313
[1002]2314        return hits;
2315}
2316
2317
[1016]2318ViewCell *VspTree::GetViewCell(const Vector3 &point, const bool active)
[1002]2319{
2320        if (mRoot == NULL)
2321                return NULL;
2322
[1008]2323        stack<VspNode *> nodeStack;
[1002]2324        nodeStack.push(mRoot);
2325 
2326        ViewCellLeaf *viewcell = NULL;
2327 
2328        while (!nodeStack.empty()) 
2329        {
[1008]2330                VspNode *node = nodeStack.top();
[1002]2331                nodeStack.pop();
2332       
2333                if (node->IsLeaf())
2334                {
[1008]2335                        viewcell = dynamic_cast<VspLeaf *>(node)->GetViewCell();
[1002]2336                        break;
2337                }
2338                else   
2339                {       
[1008]2340                        VspInterior *interior = dynamic_cast<VspInterior *>(node);
[1012]2341     
[1002]2342                        // random decision
[1012]2343                        if (interior->GetPosition() - point[interior->GetAxis()] < 0)
[1002]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
[1016]2357bool VspTree::ViewPointValid(const Vector3 &viewPoint) const
[1002]2358{
[1008]2359        VspNode *node = mRoot;
[1002]2360
2361        while (1)
2362        {
2363                // early exit
2364                if (node->TreeValid())
2365                        return true;
2366
2367                if (node->IsLeaf())
2368                        return false;
2369                       
[1008]2370                VspInterior *in = dynamic_cast<VspInterior *>(node);
[1002]2371                                       
[1012]2372                if (in->GetPosition() - viewPoint[in->GetAxis()] <= 0)
[1002]2373                {
2374                        node = in->GetBack();
2375                }
2376                else
2377                {
2378                        node = in->GetFront();
2379                }
2380        }
[1012]2381
[1002]2382        // should never come here
2383        return false;
2384}
2385
2386
[1016]2387void VspTree::PropagateUpValidity(VspNode *node)
[1002]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();
[1008]2406                        VspInterior *interior = dynamic_cast<VspInterior *>(node);
[1002]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
[1016]2416bool VspTree::Export(ogzstream &stream)
[1002]2417#else
[1016]2418bool VspTree::Export(ofstream &stream)
[1002]2419#endif
2420{
2421        ExportNode(mRoot, stream);
2422
2423        return true;
2424}
2425
[1006]2426
[1002]2427#if ZIPPED_VIEWCELLS
[1016]2428void VspTree::ExportNode(VspNode *node, ogzstream &stream)
[1002]2429#else
[1016]2430void VspTree::ExportNode(VspNode *node, ofstream &stream)
[1002]2431#endif
2432{
2433        if (node->IsLeaf())
2434        {
[1008]2435                VspLeaf *leaf = dynamic_cast<VspLeaf *>(node);
[1002]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
[1012]2445        {       
[1008]2446                VspInterior *interior = dynamic_cast<VspInterior *>(node);
[1002]2447       
[1012]2448                AxisAlignedPlane plane = interior->GetPlane();
2449                stream << "<Interior plane=\"" << plane.mPosition << " "
2450                           << plane.mAxis << "\">" << endl;
[1002]2451
2452                ExportNode(interior->GetBack(), stream);
2453                ExportNode(interior->GetFront(), stream);
[1012]2454
[1002]2455                stream << "</Interior>" << endl;
2456        }
2457}
2458
[1011]2459
[1016]2460int VspTree::SplitRays(const AxisAlignedPlane &plane,
[1020]2461                                           RayInfoContainer &rays,
2462                                           RayInfoContainer &frontRays,
2463                                           RayInfoContainer &backRays) const
[1011]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
[1047]2478                const int side = bRay.ComputeRayIntersection(plane.mAxis, plane.mPosition, t);
2479                       
2480
[1027]2481#if 1
[1011]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                }
[1027]2505#else
2506                if (side == 0)
2507                {
2508                        ++ splits;
[1011]2509
[1027]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
[1011]2531        }
2532
2533        return splits;
2534}
2535
[1012]2536
[1016]2537AxisAlignedBox3 VspTree::GetBBox(VspNode *node) const
[1012]2538{
2539        if (!node->GetParent())
[1020]2540                return mBoundingBox;
[1012]2541
2542        if (!node->IsLeaf())
2543        {
[1016]2544                return (dynamic_cast<VspInterior *>(node))->GetBoundingBox();           
[1012]2545        }
2546
2547        VspInterior *parent = dynamic_cast<VspInterior *>(node->GetParent());
2548
[1016]2549        AxisAlignedBox3 box(parent->GetBoundingBox());
[1012]2550
[1047]2551        if (parent->GetFront() == node)
2552      box.SetMin(parent->GetAxis(), parent->GetPosition());
2553    else
2554      box.SetMax(parent->GetAxis(), parent->GetPosition());
[1012]2555
2556        return box;
2557}
2558
2559
[1021]2560int VspTree::ComputeBoxIntersections(const AxisAlignedBox3 &box,
2561                                                                         ViewCellContainer &viewCells) const
2562{
[1027]2563        stack<VspNode *> nodeStack;
2564 
[1021]2565        ViewCell::NewMail();
2566
2567        while (!nodeStack.empty())
2568        {
[1027]2569                VspNode *node = nodeStack.top();
[1021]2570                nodeStack.pop();
2571
[1027]2572                const AxisAlignedBox3 bbox = GetBBox(node);
2573
2574                if (bbox.Includes(box))
[1021]2575                {
2576                        // node geometry is contained in box
2577                        CollectViewCells(node, true, viewCells, true);
[1027]2578                }
2579                else if (Overlap(bbox, box))
2580                {
[1021]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                        {
[1027]2593                                VspInterior *interior = dynamic_cast<VspInterior *>(node);
[1021]2594                       
[1027]2595                                VspNode *first = interior->GetFront();
2596                                VspNode *second = interior->GetBack();
[1021]2597           
[1027]2598                                nodeStack.push(first);
2599                                nodeStack.push(second);
[1021]2600                        }
[1027]2601                }       
2602                // default: cull
[1021]2603        }
2604
2605        return (int)viewCells.size();
2606}
2607
2608
2609
[1121]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
[1016]2642/*****************************************************************/
[1021]2643/*                class OspTree implementation                   */
[1016]2644/*****************************************************************/
2645
[1027]2646
[1022]2647OspTree::OspTree():
[1099]2648mRoot(NULL),
[1022]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
[1020]2655
[1022]2656        //-- termination criteria for autopartition
[1106]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);
[1023]2661       
[1106]2662        Environment::GetSingleton()->GetIntValue("OspTree.Termination.missTolerance", mTermMissTolerance);
2663       
[1022]2664        //-- max cost ratio for early tree termination
[1106]2665        Environment::GetSingleton()->GetFloatValue("OspTree.Termination.maxCostRatio", mTermMaxCostRatio);
[1022]2666
[1123]2667        Environment::GetSingleton()->GetFloatValue("OspTree.Termination.minGlobalCostRatio",
2668                mTermMinGlobalCostRatio);
2669        Environment::GetSingleton()->GetIntValue("OspTree.Termination.mTermGlobalCostMissTolerance);
[1022]2670
2671
2672        //-- factors for bsp tree split plane heuristics
[1106]2673        Environment::GetSingleton()->GetFloatValue("OspTree.Termination.ct_div_ci", mCtDivCi);
[1022]2674
2675        //-- partition criteria
[1106]2676        Environment::GetSingleton()->GetFloatValue("OspTree.Construction.epsilon", mEpsilon);
[1022]2677
2678        // if only the driving axis is used for axis aligned split
[1106]2679        Environment::GetSingleton()->GetBoolValue("OspTree.splitUseOnlyDrivingAxis", mOnlyDrivingAxis);
2680        Environment::GetSingleton()->GetFloatValue("OspTree.maxStaticMemory", mMaxMemory);
2681        Environment::GetSingleton()->GetBoolValue("OspTree.useCostHeuristics", mUseCostHeuristics);
[1023]2682
2683
[1022]2684        char subdivisionStatsLog[100];
[1106]2685        Environment::GetSingleton()->GetStringValue("OspTree.subdivisionStats", subdivisionStatsLog);
[1022]2686        mSubdivisionStats.open(subdivisionStatsLog);
2687
[1106]2688        Environment::GetSingleton()->GetFloatValue("OspTree.Construction.splitBorder", mSplitBorder);
[1121]2689        Environment::GetSingleton()->GetFloatValue("OspTree.Construction.renderCostDecreaseWeight", mRenderCostDecreaseWeight);
2690       
[1022]2691        //-- debug output
2692
[1076]2693        Debug << "******* OSP options ******** " << endl;
[1106]2694
[1022]2695    Debug << "max depth: " << mTermMaxDepth << endl;
2696        Debug << "min probabiliy: " << mTermMinProbability << endl;
[1106]2697        Debug << "min objects: " << mTermMinObjects << endl;
[1022]2698        Debug << "max cost ratio: " << mTermMaxCostRatio << endl;
2699        Debug << "miss tolerance: " << mTermMissTolerance << endl;
[1106]2700        Debug << "max leaves: " << mTermMaxLeaves << endl;
2701
[1022]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;
[1023]2710       
[1106]2711        Debug << "split borders: " << mSplitBorder << endl;
[1121]2712        Debug << "render cost decrease weight: " << mRenderCostDecreaseWeight << endl;
[1022]2713
[1106]2714
[1022]2715        mSplitCandidates = new vector<SortableEntry>;
2716
2717        Debug << endl;
2718}
2719
2720
2721
[1020]2722void OspTree::SplitObjects(const AxisAlignedPlane & splitPlane,
2723                                                   const ObjectContainer &objects,
[1106]2724                                                   ObjectContainer &front,
2725                                                   ObjectContainer &back)
[1016]2726{
[1020]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        }
[1106]2740
2741        mOspStats.objectRefs -= (int)objects.size();
2742        mOspStats.objectRefs += (int)back.size() + (int)front.size();
[1016]2743}
2744
2745
[1106]2746KdInterior *OspTree::SubdivideNode(
[1020]2747                                                                   const AxisAlignedPlane &splitPlane,
[1108]2748                                                                   const OspTraversalData &tData,
[1106]2749                                                                   OspTraversalData &frontData,
2750                                                                   OspTraversalData &backData
2751                                                                   )
[1016]2752{
[1108]2753        KdLeaf *leaf = tData.mNode;
2754       
2755        // two new leaves
[1106]2756    mOspStats.nodes += 2;
[1016]2757
[1020]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;
[1108]2766        node->mBox = tData.mBoundingBox;
[1020]2767
[1108]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
[1020]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               
[1121]2795                if (box.Max(axis) > position + mEpsilon)
[1020]2796                        ++ objectsFront;
2797   
[1121]2798                if (box.Min(axis) < position - mEpsilon)
[1020]2799                        ++ objectsBack;
2800        }
2801
2802        KdLeaf *back = new KdLeaf(node, objectsBack);
2803        KdLeaf *front = new KdLeaf(node, objectsFront);
2804
[1106]2805        /////////////
2806        //-- create front and back leaf
2807
2808        KdInterior *parent = leaf->mParent;
2809
[1020]2810        // replace a link from node's parent
[1106]2811        if (parent)
2812        {
2813                parent->ReplaceChildLink(leaf, node);
2814                node->mParent = parent;
2815        }
2816        else // new root
2817        {
2818                mRoot = node;
2819        }
[1020]2820
2821        // and setup child links
2822        node->SetupChildLinks(back, front);
2823
[1106]2824        SplitObjects(splitPlane, leaf->mObjects, front->mObjects, back->mObjects);
[1089]2825
[1122]2826        ProcessLeafObjects(leaf, front, back);
2827   
[1106]2828        backData.mNode = back;
2829        frontData.mNode = front;
2830
[1020]2831        //delete leaf;
2832        return node;
[1016]2833}
2834
2835
[1020]2836KdNode *OspTree::Subdivide(SplitQueue &tQueue,
2837                                                   OspSplitCandidate &splitCandidate,
2838                                                   const bool globalCriteriaMet)
[1108]2839{
[1020]2840        OspTraversalData &tData = splitCandidate.mParentData;
2841        KdNode *newNode = tData.mNode;
[1016]2842
[1020]2843        if (!LocalTerminationCriteriaMet(tData) && !globalCriteriaMet)
2844        {       
2845                OspTraversalData tFrontData;
2846                OspTraversalData tBackData;
[1108]2847
[1020]2848                //-- continue subdivision
2849               
2850                // create new interior node and two leaf node
2851                const AxisAlignedPlane splitPlane = splitCandidate.mSplitPlane;
[1027]2852               
[1108]2853                newNode = SubdivideNode(splitPlane,
2854                                                                tData,
[1106]2855                                                                tFrontData,
2856                                                                tBackData);
[1020]2857       
2858                const int maxCostMisses = splitCandidate.mMaxCostMisses;
[1016]2859
[1020]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
[1099]2865                OspSplitCandidate *frontCandidate = new OspSplitCandidate(tFrontData);
2866                OspSplitCandidate *backCandidate = new OspSplitCandidate(tBackData);
[1016]2867
[1099]2868                EvalSplitCandidate(*frontCandidate);
2869                EvalSplitCandidate(*backCandidate);
[1016]2870
[1099]2871                tQueue.Push(frontCandidate);
2872                tQueue.Push(backCandidate);
[1020]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        }
[1108]2884       
[1020]2885        //-- cleanup
2886        tData.Clear();
[1108]2887       
[1020]2888        return newNode;
2889}
[1016]2890
2891
[1099]2892void OspTree::EvalSplitCandidate(OspSplitCandidate &splitCandidate)
[1020]2893{
2894        float frontProb;
[1089]2895        float backProb;
[1020]2896       
2897        // compute locally best split plane
[1084]2898        const bool success =
[1099]2899                SelectSplitPlane(splitCandidate.mParentData, splitCandidate.mSplitPlane, frontProb, backProb);
[1020]2900
[1106]2901        const float priority = EvalRenderCostDecrease(splitCandidate.mSplitPlane, splitCandidate.mParentData);
[1020]2902        // compute global decrease in render cost
[1106]2903        splitCandidate.SetPriority(priority);
2904
2905        splitCandidate.mMaxCostMisses =
2906                success ? splitCandidate.mParentData.mMaxCostMisses : splitCandidate.mParentData.mMaxCostMisses + 1;
[1016]2907}
2908
2909
[1021]2910bool OspTree::LocalTerminationCriteriaMet(const OspTraversalData &data) const
2911{
2912        // matt: TODO
[1106]2913        return (
2914                //(data.mNode->mObjects.size() < mTermMinObjects) ||
2915                //(data.mProbability <= mTermMinProbability) ||
2916                (data.mDepth >= mTermMaxDepth)
2917                 );
[1021]2918}
[1020]2919
[1021]2920
2921bool OspTree::GlobalTerminationCriteriaMet(const OspTraversalData &data) const
2922{
2923        // matt: TODO
[1106]2924        return (
2925                (mOspStats.Leaves() >= mTermMaxLeaves)
2926                //mOutOfMemory ||
2927                //(mGlobalCostMisses >= mTermGlobalCostMissTolerance)
2928                );
[1021]2929}
2930
2931
2932void OspTree::EvaluateLeafStats(const OspTraversalData &data)
2933{
2934        // the node became a leaf -> evaluate stats for leafs
[1106]2935        KdLeaf *leaf = data.mNode;
[1021]2936
[1106]2937        if (data.mPvs > mOspStats.maxPvs)
[1021]2938        {
[1106]2939                mOspStats.maxPvs = data.mPvs;
[1021]2940        }
2941
[1106]2942        mOspStats.pvs += data.mPvs;
[1021]2943
[1106]2944        if (data.mDepth < mOspStats.minDepth)
[1021]2945        {
[1106]2946                mOspStats.minDepth = data.mDepth;
[1021]2947        }
2948       
2949        if (data.mDepth >= mTermMaxDepth)
2950        {
[1106]2951        ++ mOspStats.maxDepthNodes;
[1021]2952                //Debug << "new max depth: " << mVspStats.maxDepthNodes << endl;
2953        }
2954
[1106]2955//      if (data.mPvs < mTermMinPvs)
2956//              ++ mOspStats.minPvsNodes;
[1021]2957
2958        if (data.mProbability <= mTermMinProbability)
[1106]2959                ++ mOspStats.minProbabilityNodes;
[1021]2960       
2961        // accumulate depth to compute average depth
[1106]2962        mOspStats.accumDepth += data.mDepth;
[1021]2963
[1106]2964        ++ mCreatedLeaves;
[1021]2965
2966#ifdef _DEBUG
2967        Debug << "BSP stats: "
2968                  << "Depth: " << data.mDepth << " (max: " << mTermMaxDepth << "), "
[1106]2969                 // << "PVS: " << data.mPvs << " (min: " << mTermMinPvs << "), "
2970                  << "Prob: " << data.mProbability << " (min: " << mTermMinProbability << "), "
2971                  << "#pvs: " << data.mPvs << ")\n";
[1021]2972#endif
2973}
2974
2975
[1084]2976float OspTree::EvalLocalCostHeuristics(KdLeaf *node,
[1077]2977                                                                           const AxisAlignedBox3 &box,
2978                                                                           const int axis,
[1084]2979                                                                           float &position,
[1106]2980                                                                           int &objectsFront,
2981                                                                           int &objectsBack)
[1077]2982{
[1121]2983        // sort so we can use a sweep
[1084]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       
[1106]2990        int pvsSize = PrepareHeuristics(node->mObjects);
[1084]2991        int pvsl = 0, pvsr = pvsSize;
2992 
[1077]2993        const float minBox = box.Min(axis);
2994        const float maxBox = box.Max(axis);
[1072]2995
[1077]2996        const float sizeBox = maxBox - minBox;
2997       
[1084]2998                       
[1077]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;
[1084]3005 
3006        float minBand = minBox + mSplitBorder * (maxBox - minBox);
3007        float maxBand = minBox + (1.0f - mSplitBorder) * (maxBox - minBox);
[1077]3008
[1084]3009    float minSum = 1e20f;
[1077]3010
[1084]3011        int pvsBack = pvsl;
3012        int pvsFront = pvsr;
[1077]3013
[1084]3014        float sum = (float)pvsSize * sizeBox;
[1121]3015float dummy1 = 0;
3016float dummy2 = 0;
[1084]3017
[1077]3018        vector<SortableEntry>::const_iterator ci, ci_end = mSplitCandidates->end();
[1121]3019cout << "***********" << endl;
[1077]3020        //-- traverse through visibility events
3021        for (ci = mSplitCandidates->begin(); ci != ci_end; ++ ci)
3022        {
[1121]3023                int pvslold = pvsl;
[1077]3024                EvalPvsIncr(*ci, pvsl, pvsr);
3025
[1121]3026                cout << "incr: " << ci->mObject->mViewCellPvs.GetSize() << " obj id "
3027                         << ci->mObject->GetId() << endl;
3028
[1077]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
[1121]3034                        cout << "pos=" << (*ci).value << "\t pvs=(" <<  pvsl << ","
3035                                 << pvsr << ")" << "\t cost= " << sum << endl;
[1077]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       
[1121]3050        //-- compute cost
[1077]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)
[1121]3072        Debug << "axis=" << axis << " costRatio=" << ratio << " pos="
3073                  << position << " t=" << (position - minBox) / (maxBox - minBox)
3074              << "\t pb=(" << pvsBack << ")\t pf=(" << pvsFront << ")" << endl;
[1077]3075
3076        return ratio;
3077}
3078
[1121]3079
[1084]3080void OspTree::SortSplitCandidates(KdLeaf *node, const int axis)
[1077]3081{
[1084]3082        mSplitCandidates->clear();
[1077]3083
[1121]3084        int requestedSize = 2 * (int)node->mObjects.size();
[1077]3085       
[1084]3086        // creates a sorted split candidates array
3087        if (mSplitCandidates->capacity() > 500000 &&
3088                requestedSize < (int)(mSplitCandidates->capacity()/10))
[1077]3089        {
[1084]3090                delete mSplitCandidates;
3091                mSplitCandidates = new vector<SortableEntry>;
[1077]3092        }
3093
[1084]3094        mSplitCandidates->reserve(requestedSize);
[1077]3095       
[1084]3096        ObjectContainer::const_iterator mi, mi_end = node->mObjects.end();
[1077]3097
[1084]3098        // insert all queries
3099        for(mi = node->mObjects.begin(); mi != mi_end; ++ mi)
[1077]3100        {
[1084]3101                AxisAlignedBox3 box = (*mi)->GetBox();
[1077]3102
[1084]3103                mSplitCandidates->push_back(SortableEntry(SortableEntry::BOX_MIN,
3104                                        box.Min(axis), *mi));
[1077]3105
3106
[1084]3107                mSplitCandidates->push_back(SortableEntry(SortableEntry::BOX_MAX,
3108                                        box.Max(axis), *mi));
[1077]3109        }
3110
[1084]3111        stable_sort(mSplitCandidates->begin(), mSplitCandidates->end());
[1077]3112}
3113
3114
3115int OspTree::PrepareHeuristics(Intersectable *object)
3116{
[1123]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();
[1077]3122}
3123
3124
[1121]3125const OspTreeStatistics &OspTree::GetStatistics() const
3126{
3127        return mOspStats;
3128}
3129
3130
[1084]3131int OspTree::PrepareHeuristics(const ObjectContainer &objects)
[1077]3132{       
3133        Intersectable::NewMail();
[1084]3134        ViewCell::NewMail();
[1077]3135
3136        int pvsSize = 0;
3137
[1084]3138        ObjectContainer::const_iterator oit, oit_end = objects.end();
[1077]3139
[1084]3140        //-- set all pvs as belonging to the front pvs
3141        for (oit = objects.begin(); oit != oit_end; ++ oit)
[1077]3142        {
[1084]3143                Intersectable *obj = *oit;
3144                pvsSize += PrepareHeuristics(obj);             
[1077]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        }
[1121]3167
3168        cout << "pvs left: " << pvsLeft << " pvs right " << pvsRight << endl;
[1077]3169}
3170
3171
[1106]3172void OspTree::SetViewCellsManager(ViewCellsManager *vcm)
3173{
3174        mViewCellsManager = vcm;
3175}
3176
3177
3178AxisAlignedBox3 OspTree::GetBoundingBox() const
3179{
3180        return mBoundingBox;
3181}
3182
3183
[1077]3184void OspTree::RemoveContriFromPvs(Intersectable *object, int &pvs) const
3185{
[1123]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();
[1077]3191}
3192
3193
3194void OspTree::AddContriToPvs(Intersectable *object, int &pvs) const
3195{
[1123]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();
[1077]3201}
3202
3203
[1089]3204float OspTree::SelectSplitPlane(const OspTraversalData &tData,
[1084]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
[1106]3225        // -- evaluate split cost for all three axis
[1084]3226        for (int axis = 0; axis < 3; ++ axis)
3227        {
3228                if (!mOnlyDrivingAxis || (axis == sAxis))
3229                {
[1106]3230                        if (1 || mUseCostHeuristics)
3231                        {
3232                                //-- place split plane using heuristics
3233                                int objectsFront, objectsBack;
[1084]3234
3235                                nCostRatio[axis] =
[1106]3236                                        EvalLocalCostHeuristics(tData.mNode,
3237                                                                           tData.mBoundingBox,
3238                                                                           axis,
3239                                                                           nPosition[axis],
3240                                                                           objectsFront,
3241                                                                           objectsBack);                       
[1084]3242                        }
[1106]3243                        /*      else
[1084]3244                        {
[1106]3245                                //-- split plane position is spatial median
3246
[1084]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]);
[1106]3255                        }*/
[1084]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
[1106]3269
[1084]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
[1089]3282float OspTree::EvalViewCellPvsIncr(Intersectable *object) const
3283{
[1123]3284        if (mUseSah)
3285                return 1;
3286        else
3287                return (float)object->mViewCellPvs.GetSize();
[1089]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
[1108]3303        float pOverall = data.mBoundingBox.GetVolume();//data.mProbability;
[1089]3304        float pFront = 0;
3305        float pBack = 0;
3306
3307        Intersectable::NewMail();
[1106]3308        KdLeaf::NewMail();
[1089]3309        ViewCell::NewMail();
3310
[1106]3311        KdLeaf *leaf = data.mNode;
[1089]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();
[1108]3318               
3319                const float pvsIncr = EvalViewCellPvsIncr(obj);
3320                totalPvs += pvsIncr;
[1121]3321                cout << "totalpvs " << totalPvs << " incr " << pvsIncr << endl;
[1089]3322
3323                if (box.Max(candidatePlane.mAxis) > candidatePlane.mPosition)
[1108]3324                {
3325                        pvsFront += pvsIncr;
3326                }
[1089]3327                if (box.Min(candidatePlane.mAxis) > candidatePlane.mPosition)
[1108]3328                {
3329                        pvsBack += pvsIncr;
3330                }
[1089]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
[1108]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);
[1089]3351                       
[1108]3352        Debug << "total pvs: " << totalPvs << " p " << pOverall << endl;
[1089]3353        const float oldRenderCost = pOverall * penaltyOld;
3354        const float newRenderCost = penaltyFront * pFront + penaltyBack * pBack;
3355
[1121]3356        Debug << "old: " << oldRenderCost << " new: " << newRenderCost << " decrease: "
3357                << oldRenderCost - newRenderCost << endl;
[1089]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!!
[1121]3362        const float factor = mRenderCostDecreaseWeight;
[1089]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{
[1099]3372        // store pointer to this tree
3373        OspSplitCandidate::sOspTree = this;
3374
[1089]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        }
[1108]3397
3398        mTermMinProbability *= mBoundingBox.GetVolume();
3399        mGlobalCostMisses = 0;
[1089]3400}
3401
3402
[1122]3403void OspTree::ProcessLeafObjects(KdLeaf *parent, KdLeaf *front, KdLeaf *back) const
[1089]3404{
[1122]3405        if (parent)
[1089]3406        {
[1122]3407                // remove the parents from the set
3408                ObjectContainer::const_iterator oit, oit_end = parent->mObjects.end();
[1089]3409
[1122]3410                for (oit = parent->mObjects.begin(); oit != oit_end; ++ oit)
[1089]3411                {
[1122]3412                        Intersectable *object = *oit;
3413
[1089]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                }
[1122]3420        }
[1089]3421
[1122]3422        //Intersectable::NewMail();
[1089]3423
[1122]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                }
[1089]3434        }
[1122]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       
[1089]3478}
3479
3480
[1106]3481void OspTree::CollectLeaves(vector<KdLeaf *> &leaves) const
3482{
3483        stack<KdNode *> nodeStack;
3484        nodeStack.push(mRoot);
[1089]3485
[1106]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}
[1021]3503
[1072]3504
[1106]3505AxisAlignedBox3 OspTree::GetBBox(KdNode *node) const
3506{
3507        if (!node->mParent)
3508                return mBoundingBox;
[1072]3509
[1106]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
[1121]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
[1106]3596/********************************************************************/
3597/*               class HierarchyManager implementation              */
3598/********************************************************************/
3599
3600
3601
[1027]3602HierarchyManager::HierarchyManager(VspTree &vspTree, OspTree &ospTree):
3603mVspTree(vspTree), mOspTree(ospTree)
[1016]3604{
[1020]3605}
3606
3607
3608SplitCandidate *HierarchyManager::NextSplitCandidate()
3609{
[1106]3610        SplitCandidate *splitCandidate = mTQueue.Top();
3611
[1121]3612        cout << "next candidate: " << splitCandidate->Type()
[1106]3613                << ", priority: " << splitCandidate->GetPriority() << endl;
3614
[1099]3615        mTQueue.Pop();
[1020]3616
3617        return splitCandidate;
3618}
3619
3620
[1121]3621void HierarchyManager::ProcessRays(const VssRayContainer &sampleRays,
3622                                                                   RayInfoContainer &rays,
3623                                                                   ObjectContainer &sampledObjects)
[1020]3624{
[1121]3625        VssRayContainer::const_iterator rit, rit_end = sampleRays.end();
[1020]3626
3627        long startTime = GetTime();
3628
3629        cout << "storing rays ... ";
3630
3631        Intersectable::NewMail();
3632
[1121]3633        //-- store rays and objects
[1016]3634        for (rit = sampleRays.begin(); rit != rit_end; ++ rit)
3635        {
3636                VssRay *ray = *rit;
3637                float minT, maxT;
[1121]3638                static Ray hray;
[1016]3639
3640                hray.Init(*ray);
[1027]3641               
[1016]3642                // TODO: not very efficient to implictly cast between rays types
[1027]3643                if (mVspTree.GetBoundingBox().GetRaySegment(hray, minT, maxT))
[1016]3644                {
3645                        float len = ray->Length();
3646
3647                        if (!len)
3648                                len = Limits::Small;
3649
[1020]3650                        rays.push_back(RayInfo(ray, minT / len, maxT / len));
[1016]3651                }
[1121]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                }
[1016]3665        }
[1020]3666
[1027]3667        cout << "finished in " << TimeDiff(startTime, GetTime()) * 1e-3 << " secs" << endl;
[1121]3668}
[1020]3669
[1089]3670
[1121]3671void HierarchyManager::PrepareVsp(const VssRayContainer &sampleRays,
3672                                                                  AxisAlignedBox3 *forcedViewSpace,
3673                                                                  RayInfoContainer &rays)
3674{
3675        // prepare bounding box
3676        mVspTree.PrepareConstruction(sampleRays, forcedViewSpace);
[1089]3677
[1121]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;
[1089]3686        // -- prepare view space partition
3687
3688        // add first candidate for view space partition
[1109]3689        VspLeaf *leaf = new VspLeaf();
3690        mVspTree.mRoot = leaf;
[1027]3691        const float prop = mVspTree.mBoundingBox.GetVolume();
3692
[1089]3693        // first vsp traversal data
[1109]3694        VspTree::VspTraversalData vData(leaf,
[1027]3695                                                                        0,
3696                                                                        &rays,
[1089]3697                                                                        pvsSize,       
[1027]3698                                                                        prop,
3699                                                                        mVspTree.mBoundingBox);
3700
[1089]3701
[1121]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
[1027]3718        // compute first split candidate
[1099]3719        VspTree::VspSplitCandidate *splitCandidate = new VspTree::VspSplitCandidate(vData);
3720    mVspTree.EvalSplitCandidate(*splitCandidate);
[1027]3721
[1099]3722        mTQueue.Push(splitCandidate);
[1121]3723}
[1027]3724
[1089]3725
[1121]3726void HierarchyManager::PrepareOsp(const ObjectContainer &objects,                                                                                 
3727                                                                  AxisAlignedBox3 *forcedViewSpace)
3728{
[1089]3729        mOspTree.PrepareConstruction(objects, forcedViewSpace);
3730
3731        // add first candidate for view space partition
[1109]3732        KdLeaf *kdleaf = new KdLeaf(NULL, 0);
3733        kdleaf->mObjects = objects;
[1089]3734
[1109]3735        mOspTree.mRoot = kdleaf;
[1089]3736       
[1121]3737        // TODO matt: this is something different than pvs size.
3738        const int pvsSize = mOspTree.ComputePvsSize(objects);
3739        const float prop = mOspTree.mBoundingBox.GetVolume();
[1089]3740
[1121]3741        cout <<  "here40 pvs size: " << pvsSize << endl;
[1089]3742        // first osp traversal data
[1109]3743        OspTree::OspTraversalData oData(kdleaf,
[1089]3744                                                                        0,
3745                                                                        pvsSize,
3746                                                                        prop,
3747                                                                        mOspTree.mBoundingBox);
3748
3749               
[1122]3750        mOspTree.ProcessLeafObjects(NULL, kdleaf, NULL);
[1089]3751
3752        // compute first split candidate
[1121]3753        OspTree::OspSplitCandidate *oSplitCandidate =
3754                new OspTree::OspSplitCandidate(oData);
[1099]3755    mOspTree.EvalSplitCandidate(*oSplitCandidate);
[1089]3756
[1121]3757        mTQueue.Push(oSplitCandidate);
[1016]3758}
3759
3760
[1121]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
[1020]3771bool HierarchyManager::GlobalTerminationCriteriaMet(SplitCandidate *candidate) const
[1016]3772{
[1108]3773        // TODO matt: make virtual by creating superclasss for traveral data
3774        return candidate->GlobalTerminationCriteriaMet();
[1016]3775}
3776
3777
[1101]3778void HierarchyManager::Construct2(const VssRayContainer &sampleRays,
[1106]3779                                                                  const ObjectContainer &objects,
3780                                                                  AxisAlignedBox3 *forcedViewSpace)
[1101]3781{
3782        RayInfoContainer *rays = new RayInfoContainer();
[1123]3783
3784        // makes no sense otherwise because only one kd cell available
3785        // during view space partition
3786        mVspTree.mPvsCountMethod = VspTree::PER_OBJECT;
3787
[1121]3788        PrepareVsp(sampleRays, forcedViewSpace, *rays);
[1101]3789
3790        long startTime = GetTime();
[1121]3791        cout << "starting vsp contruction ... " << endl;
[1101]3792        int i = 0;
3793
3794        while (!FinishedConstruction())
3795        {
3796                SplitCandidate *splitCandidate = NextSplitCandidate();
3797           
3798                const bool globalTerminationCriteriaMet =
3799                        GlobalTerminationCriteriaMet(splitCandidate);
3800
[1121]3801                cout << "vsp nodes: " << i ++ << endl;
[1101]3802
3803                // cost ratio of cost decrease / totalCost
3804                const float costRatio = splitCandidate->GetPriority() / mTotalCost;
[1121]3805                //cout << "cost ratio: " << costRatio << endl;
[1101]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
[1121]3816       
[1101]3817                VspNode *r = mVspTree.Subdivide(mTQueue, *sc, globalTerminationCriteriaMet);
[1121]3818
[1101]3819                DEL_PTR(splitCandidate);
3820        }
3821
3822        cout << "finished in " << TimeDiff(startTime, GetTime())*1e-3 << " secs" << endl;
3823        mVspTree.mVspStats.Stop();
[1106]3824       
[1121]3825        startTime = GetTime();
3826        cout << "starting osp contruction ... " << endl;
[1101]3827
3828        //-- object space partition
[1121]3829    PrepareOsp(objects, forcedViewSpace);
[1101]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;
[1108]3841               
[1101]3842                //Debug << "cost ratio: " << costRatio << endl;
[1108]3843                //cout << "kd nodes: " << i ++ << endl;
[1101]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        }
[1106]3858       
[1101]3859        cout << "finished in " << TimeDiff(startTime, GetTime())*1e-3 << " secs" << endl;
3860}
3861
3862
[1020]3863void HierarchyManager::Construct(const VssRayContainer &sampleRays,
3864                                                                 const ObjectContainer &objects,
3865                                                                 AxisAlignedBox3 *forcedViewSpace)
[1016]3866{
[1020]3867        RayInfoContainer *rays = new RayInfoContainer();
3868       
3869        // prepare vsp and osp trees for traversal
[1027]3870        PrepareConstruction(sampleRays, objects, forcedViewSpace, *rays);
[1016]3871
[1027]3872        mVspTree.mVspStats.Reset();
3873        mVspTree.mVspStats.Start();
[1016]3874
[1020]3875        cout << "Constructing view space / object space tree ... \n";
[1016]3876        const long startTime = GetTime();       
[1074]3877       
3878        int i = 0;
[1016]3879        while (!FinishedConstruction())
3880        {
[1121]3881                mCurrentCandidate = NextSplitCandidate();
3882           
[1020]3883                const bool globalTerminationCriteriaMet =
[1121]3884                        GlobalTerminationCriteriaMet(mCurrentCandidate);
[1020]3885
[1074]3886                cout << "view cells: " << i ++ << endl;
[1027]3887
[1016]3888                // cost ratio of cost decrease / totalCost
[1121]3889                const float costRatio = mCurrentCandidate->GetPriority() / mTotalCost;
[1016]3890                //Debug << "cost ratio: " << costRatio << endl;
3891
3892                if (costRatio < mTermMinGlobalCostRatio)
3893                        ++ mGlobalCostMisses;
[1027]3894
[1020]3895                //-- subdivide leaf node
[1101]3896
3897                // we have either a object space or view space split
[1121]3898                if (mCurrentCandidate->Type() == SplitCandidate::VIEW_SPACE)
[1016]3899                {
3900                        VspTree::VspSplitCandidate *sc =
[1121]3901                                dynamic_cast<VspTree::VspSplitCandidate *>(mCurrentCandidate);
[1016]3902
[1027]3903                        VspNode *r = mVspTree.Subdivide(mTQueue, *sc, globalTerminationCriteriaMet);
[1016]3904                }
[1020]3905                else // object space split
[1089]3906                {                       
3907                        OspTree::OspSplitCandidate *sc =
[1121]3908                                dynamic_cast<OspTree::OspSplitCandidate *>(mCurrentCandidate);
[1089]3909
3910                        KdNode *r = mOspTree.Subdivide(mTQueue, *sc, globalTerminationCriteriaMet);
[1016]3911                }
[1020]3912
[1121]3913                // reevaluate affected candidates
3914                RepairQueue();
3915
3916                DEL_PTR(mCurrentCandidate);
[1016]3917        }
3918
[1074]3919        cout << "finished in " << TimeDiff(startTime, GetTime())*1e-3 << " secs" << endl;
[1016]3920
[1027]3921        mVspTree.mVspStats.Stop();
[1016]3922}
3923
[1020]3924
3925bool HierarchyManager::FinishedConstruction()
3926{
[1099]3927        return mTQueue.Empty();
[1020]3928}
3929
3930
[1121]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()
[1093]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
[1099]3967        // *) search is not efficient => store queue position with each
[1093]3968        // split candidate
3969
[1121]3970        // collect list of "dirty" candidates
3971        vector<SplitCandidate *> dirtyList;
3972        CollectDirtyCandidates(dirtyList);
3973       
3974        //-- reevaluate the dirty list
[1099]3975        vector<SplitCandidate *>::const_iterator sit, sit_end = dirtyList.end();
[1093]3976
[1099]3977        for (sit = dirtyList.begin(); sit != sit_end; ++ sit)
[1093]3978        {
[1099]3979                SplitCandidate* sc = *sit;
[1093]3980
[1099]3981                mTQueue.Erase(sc);
[1093]3982
[1099]3983                // reevaluate
3984                sc->EvalPriority();
[1093]3985
[1099]3986                // reinsert
3987                mTQueue.Push(sc);
3988        }
[1093]3989}
3990
[1002]3991}
Note: See TracBrowser for help on using the repository browser.