source: trunk/VUT/GtpVisibilityPreprocessor/src/VspBspTree.cpp @ 542

Revision 542, 74.9 KB checked in by mattausch, 18 years ago (diff)
Line 
1#include <stack>
2#include <time.h>
3#include <iomanip>
4
5#include "Plane3.h"
6#include "VspBspTree.h"
7#include "Mesh.h"
8#include "common.h"
9#include "ViewCell.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 "ViewCellBsp.h"
17#include "ViewCellsManager.h"
18#include "Beam.h"
19
20//-- static members
21
22/** Evaluates split plane classification with respect to the plane's
23        contribution for a minimum number of ray splits.
24*/
25const float VspBspTree::sLeastRaySplitsTable[] = {0, 0, 1, 1, 0};
26/** Evaluates split plane classification with respect to the plane's
27        contribution for balanced rays.
28*/
29const float VspBspTree::sBalancedRaysTable[] = {1, -1, 0, 0, 0};
30
31int BspMergeCandidate::sMaxPvsSize = 0;
32
33int VspBspTree::sFrontId = 0;
34int VspBspTree::sBackId = 0;
35int VspBspTree::sFrontAndBackId = 0;
36
37float BspMergeCandidate::sOverallCost = 0;
38
39/********************************************************************/
40/*                  class VspBspTree implementation                 */
41/********************************************************************/
42
43VspBspTree::VspBspTree():
44mRoot(NULL),
45mPvsUseArea(true),
46mCostNormalizer(Limits::Small),
47mViewCellsManager(NULL),
48mOutOfBoundsCell(NULL),
49mShowInvalidSpace(false),
50mStoreRays(false)
51{
52        bool randomize = false;
53        environment->GetBoolValue("VspBspTree.Construction.randomize", randomize);
54        if (randomize)
55                Randomize(); // initialise random generator for heuristics
56
57        //-- termination criteria for autopartition
58        environment->GetIntValue("VspBspTree.Termination.maxDepth", mTermMaxDepth);
59        environment->GetIntValue("VspBspTree.Termination.minPvs", mTermMinPvs);
60        environment->GetIntValue("VspBspTree.Termination.minRays", mTermMinRays);
61        environment->GetFloatValue("VspBspTree.Termination.minArea", mTermMinArea);
62        environment->GetFloatValue("VspBspTree.Termination.maxRayContribution", mTermMaxRayContribution);
63        environment->GetFloatValue("VspBspTree.Termination.minAccRayLenght", mTermMinAccRayLength);
64        environment->GetFloatValue("VspBspTree.Termination.maxCostRatio", mTermMaxCostRatio);
65        environment->GetIntValue("VspBspTree.Termination.missTolerance", mTermMissTolerance);
66
67        //-- factors for bsp tree split plane heuristics
68        environment->GetFloatValue("VspBspTree.Factor.balancedRays", mBalancedRaysFactor);
69        environment->GetFloatValue("VspBspTree.Factor.pvs", mPvsFactor);
70        environment->GetFloatValue("VspBspTree.Termination.ct_div_ci", mCtDivCi);
71
72        //-- max cost ratio for early tree termination
73        environment->GetFloatValue("VspBspTree.Termination.maxCostRatio", mTermMaxCostRatio);
74
75        //-- partition criteria
76        environment->GetIntValue("VspBspTree.maxPolyCandidates", mMaxPolyCandidates);
77        environment->GetIntValue("VspBspTree.maxRayCandidates", mMaxRayCandidates);
78        environment->GetIntValue("VspBspTree.splitPlaneStrategy", mSplitPlaneStrategy);
79
80        environment->GetFloatValue("VspBspTree.Construction.epsilon", mEpsilon);
81        environment->GetIntValue("VspBspTree.maxTests", mMaxTests);
82
83        // maximum and minimum number of view cells
84        environment->GetIntValue("VspBspTree.Termination.maxViewCells", mMaxViewCells);
85        environment->GetIntValue("VspBspTree.PostProcess.minViewCells", mMergeMinViewCells);
86        environment->GetFloatValue("VspBspTree.PostProcess.maxCostRatio", mMergeMaxCostRatio);
87
88        environment->GetBoolValue("VspBspTree.splitUseOnlyDrivingAxis", mOnlyDrivingAxis);
89        environment->GetBoolValue("VspBspTree.PostProcess.useRaysForMerge", mUseRaysForMerge);
90
91        environment->GetFloatValue("ViewCells.maxPvsRatio", mMaxPvsRatio);
92
93        //-- termination criteria for axis aligned split
94        environment->GetFloatValue("VspBspTree.Termination.AxisAligned.maxRayContribution",
95                mTermMaxRayContriForAxisAligned);
96        environment->GetIntValue("VspBspTree.Termination.AxisAligned.minRays",
97                mTermMinRaysForAxisAligned);
98
99        //environment->GetFloatValue("VspBspTree.maxTotalMemory", mMaxTotalMemory);
100        environment->GetFloatValue("VspBspTree.maxStaticMemory", mMaxMemory);
101
102        //-- debug output
103        Debug << "******* VSP BSP options ******** " << endl;
104    Debug << "max depth: " << mTermMaxDepth << endl;
105        Debug << "min PVS: " << mTermMinPvs << endl;
106        Debug << "min area: " << mTermMinArea << endl;
107        Debug << "min rays: " << mTermMinRays << endl;
108        Debug << "max ray contri: " << mTermMaxRayContribution << endl;
109        //Debug << "VSP BSP mininam accumulated ray lenght: ", mTermMinAccRayLength) << endl;
110        Debug << "max cost ratio: " << mTermMaxCostRatio << endl;
111        Debug << "miss tolerance: " << mTermMissTolerance << endl;
112        Debug << "max view cells: " << mMaxViewCells << endl;
113        Debug << "max polygon candidates: " << mMaxPolyCandidates << endl;
114        Debug << "max plane candidates: " << mMaxRayCandidates << endl;
115        Debug << "randomize: " << randomize << endl;
116
117        Debug << "Split plane strategy: ";
118        if (mSplitPlaneStrategy & RANDOM_POLYGON)
119        {
120                Debug << "random polygon ";
121        }
122        if (mSplitPlaneStrategy & AXIS_ALIGNED)
123        {
124                Debug << "axis aligned ";
125        }
126        if (mSplitPlaneStrategy & LEAST_RAY_SPLITS)
127        {
128                mCostNormalizer += mLeastRaySplitsFactor;
129                Debug << "least ray splits ";
130        }
131        if (mSplitPlaneStrategy & BALANCED_RAYS)
132        {
133                mCostNormalizer += mBalancedRaysFactor;
134                Debug << "balanced rays ";
135        }
136        if (mSplitPlaneStrategy & PVS)
137        {
138                mCostNormalizer += mPvsFactor;
139                Debug << "pvs";
140        }
141
142
143        mSplitCandidates = new vector<SortableEntry>;
144
145        Debug << endl;
146}
147
148BspViewCell *VspBspTree::GetOutOfBoundsCell()
149{
150        return mOutOfBoundsCell;
151}
152
153
154BspViewCell *VspBspTree::GetOrCreateOutOfBoundsCell()
155{
156        if (!mOutOfBoundsCell)
157        {
158                mOutOfBoundsCell = new BspViewCell();
159                mOutOfBoundsCell->SetId(-1);
160        }
161        return mOutOfBoundsCell;
162}
163
164
165const BspTreeStatistics &VspBspTree::GetStatistics() const
166{
167        return mStat;
168}
169
170
171VspBspTree::~VspBspTree()
172{
173        DEL_PTR(mRoot);
174        DEL_PTR(mSplitCandidates);
175}
176
177int VspBspTree::AddMeshToPolygons(Mesh *mesh,
178                                                                  PolygonContainer &polys,
179                                                                  MeshInstance *parent)
180{
181        FaceContainer::const_iterator fi;
182
183        // copy the face data to polygons
184        for (fi = mesh->mFaces.begin(); fi != mesh->mFaces.end(); ++ fi)
185        {
186                Polygon3 *poly = new Polygon3((*fi), mesh);
187
188                if (poly->Valid(mEpsilon))
189                {
190                        poly->mParent = parent; // set parent intersectable
191                        polys.push_back(poly);
192                }
193                else
194                        DEL_PTR(poly);
195        }
196        return (int)mesh->mFaces.size();
197}
198
199int VspBspTree::AddToPolygonSoup(const ViewCellContainer &viewCells,
200                                                          PolygonContainer &polys,
201                                                          int maxObjects)
202{
203        int limit = (maxObjects > 0) ?
204                Min((int)viewCells.size(), maxObjects) : (int)viewCells.size();
205
206        int polysSize = 0;
207
208        for (int i = 0; i < limit; ++ i)
209        {
210                if (viewCells[i]->GetMesh()) // copy the mesh data to polygons
211                {
212                        mBox.Include(viewCells[i]->GetBox()); // add to BSP tree aabb
213                        polysSize +=
214                                AddMeshToPolygons(viewCells[i]->GetMesh(),
215                                                                  polys,
216                                                                  viewCells[i]);
217                }
218        }
219
220        return polysSize;
221}
222
223int VspBspTree::AddToPolygonSoup(const ObjectContainer &objects,
224                                                                 PolygonContainer &polys,
225                                                                 int maxObjects)
226{
227        int limit = (maxObjects > 0) ?
228                Min((int)objects.size(), maxObjects) : (int)objects.size();
229
230        for (int i = 0; i < limit; ++i)
231        {
232                Intersectable *object = objects[i];//*it;
233                Mesh *mesh = NULL;
234
235                switch (object->Type()) // extract the meshes
236                {
237                case Intersectable::MESH_INSTANCE:
238                        mesh = dynamic_cast<MeshInstance *>(object)->GetMesh();
239                        break;
240                case Intersectable::VIEW_CELL:
241                        mesh = dynamic_cast<ViewCell *>(object)->GetMesh();
242                        break;
243                        // TODO: handle transformed mesh instances
244                default:
245                        Debug << "intersectable type not supported" << endl;
246                        break;
247                }
248
249        if (mesh) // copy the mesh data to polygons
250                {
251                        mBox.Include(object->GetBox()); // add to BSP tree aabb
252                        AddMeshToPolygons(mesh, polys, NULL);
253                }
254        }
255
256        return (int)polys.size();
257}
258
259void VspBspTree::Construct(const VssRayContainer &sampleRays,
260                                                   AxisAlignedBox3 *forcedBoundingBox)
261{
262    mStat.nodes = 1;
263        mBox.Initialize();      // initialise BSP tree bounding box
264
265        if (forcedBoundingBox)
266                mBox = *forcedBoundingBox;
267       
268        PolygonContainer polys;
269        RayInfoContainer *rays = new RayInfoContainer();
270
271        VssRayContainer::const_iterator rit, rit_end = sampleRays.end();
272
273        long startTime = GetTime();
274
275        cout << "Extracting polygons from rays ... ";
276
277        Intersectable::NewMail();
278
279        int numObj = 0;
280
281        //-- extract polygons intersected by the rays
282        for (rit = sampleRays.begin(); rit != rit_end; ++ rit)
283        {
284                VssRay *ray = *rit;
285
286                if ((mBox.IsInside(ray->mTermination) || !forcedBoundingBox) &&
287                        ray->mTerminationObject &&
288                        !ray->mTerminationObject->Mailed())
289                {
290                        ray->mTerminationObject->Mail();
291                        MeshInstance *obj = dynamic_cast<MeshInstance *>(ray->mTerminationObject);
292                        AddMeshToPolygons(obj->GetMesh(), polys, obj);
293                        ++ numObj;
294                        //-- compute bounding box
295                        if (!forcedBoundingBox)
296                                mBox.Include(ray->mTermination);
297                }
298
299                if ((mBox.IsInside(ray->mOrigin) || !forcedBoundingBox) &&
300                        ray->mOriginObject &&
301                        !ray->mOriginObject->Mailed())
302                {
303                        ray->mOriginObject->Mail();
304                        MeshInstance *obj = dynamic_cast<MeshInstance *>(ray->mOriginObject);
305                        AddMeshToPolygons(obj->GetMesh(), polys, obj);
306                        ++ numObj;
307
308                        //-- compute bounding box
309                        if (!forcedBoundingBox)
310                                mBox.Include(ray->mOrigin);
311                }
312        }
313
314        mMaxPvs = (int)(mMaxPvsRatio * (float)numObj);
315       
316        Debug << "maximal pvs (i.e., pvs still considered as valid) : " << mMaxPvs << endl;
317        //-- store rays
318        for (rit = sampleRays.begin(); rit != rit_end; ++ rit)
319        {
320                VssRay *ray = *rit;
321
322                float minT, maxT;
323
324                // TODO: not very efficient to implictly cast between rays types
325                if (mBox.GetRaySegment(*ray, minT, maxT))
326                {
327                        float len = ray->Length();
328
329                        if (!len)
330                                len = Limits::Small;
331
332                        rays->push_back(RayInfo(ray, minT / len, maxT / len));
333                }
334        }
335
336        mTermMinArea *= mBox.SurfaceArea(); // normalize
337        mStat.polys = (int)polys.size();
338
339        cout << "finished" << endl;
340
341        Debug << "\nPolygon extraction: " << (int)polys.size() << " polys extracted from "
342                  << (int)sampleRays.size() << " rays in "
343                  << TimeDiff(startTime, GetTime())*1e-3 << " secs" << endl << endl;
344
345        Construct(polys, rays);
346
347        // clean up polygons
348        CLEAR_CONTAINER(polys);
349}
350
351
352// return memory usage in MB
353float VspBspTree::GetMemUsage() const
354{
355        return
356                (sizeof(VspBspTree) +
357                 mStat.Leaves() * sizeof(BspLeaf) +
358                 mStat.Interior() * sizeof(BspInterior) +
359                 mStat.accumRays * sizeof(RayInfo)) / (1024.0f * 1024.0f);
360}
361
362
363
364void VspBspTree::Construct(const PolygonContainer &polys, RayInfoContainer *rays)
365{
366        VspBspTraversalStack tStack;
367
368        mRoot = new BspLeaf();
369
370        // constrruct root node geometry
371        BspNodeGeometry *geom = new BspNodeGeometry();
372        ConstructGeometry(mRoot, *geom);
373
374        VspBspTraversalData tData(mRoot,
375                                                          new PolygonContainer(polys),
376                                                          0,
377                                                          rays,
378                              ComputePvsSize(*rays),
379                                                          geom->GetArea(),
380                                                          geom);
381
382        tStack.push(tData);
383
384        mStat.Start();
385        cout << "Contructing vsp bsp tree ... \n";
386
387        long startTime = GetTime();     
388        long interTime = GetTime();     
389        mOutOfMemory = false;
390
391        int nleaves = 500;
392
393        while (!tStack.empty())
394        {
395                tData = tStack.top();
396            tStack.pop();               
397
398                if (0 && !mOutOfMemory)
399                {
400                        float mem = GetMemUsage();
401
402                        if (mem > mMaxMemory)
403                        {
404                                mOutOfMemory = true;
405                                Debug << "memory limit reached: " << mem << endl;
406                        }
407                }
408
409                        // subdivide leaf node
410                BspNode *r = Subdivide(tStack, tData);
411
412                if (r == mRoot)
413                        Debug << "VSP BSP tree construction time spent at root: "
414                                  << TimeDiff(startTime, GetTime())*1e-3 << "s" << endl;
415
416                if (mStat.Leaves() >= nleaves)
417                {
418                        nleaves += 500;
419                               
420                        cout << "leaves=" << mStat.Leaves() << endl;
421                        Debug << "needed "
422                          << TimeDiff(interTime, GetTime())*1e-3 << " secs to create 500 leaves" << endl;
423                        interTime = GetTime();
424                }
425        }
426
427        Debug << "Used Memory: " << GetMemUsage() << " MB" << endl << endl;
428        cout << "finished\n";
429
430        mStat.Stop();
431}
432
433
434bool VspBspTree::TerminationCriteriaMet(const VspBspTraversalData &data) const
435{
436        return
437                (((int)data.mRays->size() <= mTermMinRays) ||
438                 (data.mPvs <= mTermMinPvs)   ||
439                 (data.mArea <= mTermMinArea) ||
440                 (mStat.Leaves() >= mMaxViewCells) ||
441                 (data.GetAvgRayContribution() > mTermMaxRayContribution) ||
442                 (data.mDepth >= mTermMaxDepth));
443}
444
445
446BspNode *VspBspTree::Subdivide(VspBspTraversalStack &tStack,
447                                                           VspBspTraversalData &tData)
448{
449        BspNode *newNode = tData.mNode;
450
451        if (!mOutOfMemory && !TerminationCriteriaMet(tData))
452        {
453                PolygonContainer coincident;
454
455                VspBspTraversalData tFrontData;
456                VspBspTraversalData tBackData;
457
458                // create new interior node and two leaf nodes
459                // or return leaf as it is (if maxCostRatio missed)
460                newNode = SubdivideNode(tData, tFrontData, tBackData, coincident);
461
462                if (!newNode->IsLeaf()) //-- continue subdivision
463                {
464                        // push the children on the stack
465                        tStack.push(tFrontData);
466                        tStack.push(tBackData);
467
468                        // delete old leaf node
469                        DEL_PTR(tData.mNode);
470                }
471        }
472
473        //-- terminate traversal and create new view cell
474        if (newNode->IsLeaf())
475        {
476                BspLeaf *leaf = dynamic_cast<BspLeaf *>(newNode);
477                BspViewCell *viewCell;
478
479                if (!CheckValid(tData))
480                {
481                        leaf->SetTreeValid(false);
482                        PropagateUpValidity(leaf);
483                        // view cell for invalid view space
484                        viewCell = GetOrCreateOutOfBoundsCell();
485                        ++ mStat.invalidLeaves;
486                }
487                else
488                {       // create new view cell for this leaf                   
489                        viewCell = new BspViewCell();
490                }
491               
492                leaf->SetViewCell(viewCell);
493       
494                //-- update pvs
495                int conSamp = 0, sampCon = 0;
496                AddToPvs(leaf, *tData.mRays, conSamp, sampCon);
497
498                mStat.contributingSamples += conSamp;
499                mStat.sampleContributions += sampCon;
500
501                //-- store additional info
502                if (mStoreRays)
503                {
504                        RayInfoContainer::const_iterator it, it_end = tData.mRays->end();
505                        for (it = tData.mRays->begin(); it != it_end; ++ it)
506                                leaf->mVssRays.push_back((*it).mRay);
507                }
508
509                viewCell->mLeaves.push_back(leaf);
510                viewCell->SetArea(tData.mArea);
511                leaf->mArea = tData.mArea;
512
513                EvaluateLeafStats(tData);               
514        }
515
516        //-- cleanup
517        tData.Clear();
518
519        return newNode;
520}
521
522
523
524
525BspNode *VspBspTree::SubdivideNode(VspBspTraversalData &tData,
526                                                                   VspBspTraversalData &frontData,
527                                                                   VspBspTraversalData &backData,
528                                                                   PolygonContainer &coincident)
529{
530        BspLeaf *leaf = dynamic_cast<BspLeaf *>(tData.mNode);
531       
532        // select subdivision plane
533        Plane3 splitPlane;
534       
535        int maxCostMisses = tData.mMaxCostMisses;
536
537        const bool success =
538                SelectPlane(splitPlane, leaf, tData, frontData, backData);
539
540        if (!success)
541        {
542                ++ maxCostMisses;
543
544                if (maxCostMisses > mTermMissTolerance)
545                {
546                        // terminate branch because of max cost
547                        ++ mStat.maxCostNodes;
548            return leaf;
549                }
550        }
551
552        mStat.nodes += 2;
553
554        //-- subdivide further
555        BspInterior *interior = new BspInterior(splitPlane);
556
557#ifdef _DEBUG
558        Debug << interior << endl;
559#endif
560
561        //-- the front and back traversal data is filled with the new values
562        frontData.mDepth = tData.mDepth + 1;
563        frontData.mPolygons = new PolygonContainer();
564        frontData.mRays = new RayInfoContainer();
565       
566        backData.mDepth = tData.mDepth + 1;
567        backData.mPolygons = new PolygonContainer();
568        backData.mRays = new RayInfoContainer();
569       
570        // subdivide rays
571        SplitRays(interior->GetPlane(),
572                          *tData.mRays,
573                          *frontData.mRays,
574                          *backData.mRays);
575
576        // subdivide polygons
577        SplitPolygons(interior->GetPlane(),
578                                  *tData.mPolygons,
579                      *frontData.mPolygons,
580                                  *backData.mPolygons,
581                                  coincident);
582
583
584        // how often was max cost ratio missed in this branch?
585        frontData.mMaxCostMisses = maxCostMisses;
586        backData.mMaxCostMisses = maxCostMisses;
587
588        // compute pvs
589        frontData.mPvs = ComputePvsSize(*frontData.mRays);
590        backData.mPvs = ComputePvsSize(*backData.mRays);
591
592        // split front and back node geometry and compute area
593        if (mPvsUseArea)
594        {
595                // if geometry was not already computed
596                if (!frontData.mGeometry && !backData.mGeometry)
597                {
598                        frontData.mGeometry = new BspNodeGeometry();
599                        backData.mGeometry = new BspNodeGeometry();
600
601                        tData.mGeometry->SplitGeometry(*frontData.mGeometry,
602                                                                                   *backData.mGeometry,
603                                                                                   interior->GetPlane(),
604                                                                                   mBox,
605                                                                                   mEpsilon);
606               
607                        frontData.mArea = frontData.mGeometry->GetArea();
608                        backData.mArea = backData.mGeometry->GetArea();
609                }
610        }
611        else
612        {
613                frontData.mArea = (float)frontData.mRays->size();
614                backData.mArea = (float)backData.mRays->size();
615        }
616
617        //-- create front and back leaf
618
619        BspInterior *parent = leaf->GetParent();
620
621        // replace a link from node's parent
622        if (parent)
623        {
624                parent->ReplaceChildLink(leaf, interior);
625                interior->SetParent(parent);
626        }
627        else // new root
628        {
629                mRoot = interior;
630        }
631
632        // and setup child links
633        interior->SetupChildLinks(new BspLeaf(interior), new BspLeaf(interior));
634
635        frontData.mNode = interior->GetFront();
636        backData.mNode = interior->GetBack();
637
638        //DEL_PTR(leaf);
639        return interior;
640}
641
642
643void VspBspTree::AddToPvs(BspLeaf *leaf,
644                                                  const RayInfoContainer &rays,
645                                                  int &sampleContributions,
646                                                  int &contributingSamples)
647{
648        sampleContributions = 0;
649        contributingSamples = 0;
650
651    RayInfoContainer::const_iterator it, it_end = rays.end();
652
653        ViewCell *vc = leaf->GetViewCell();
654
655        // add contributions from samples to the PVS
656        for (it = rays.begin(); it != it_end; ++ it)
657        {
658                int contribution = 0;
659                VssRay *ray = (*it).mRay;
660
661                if (ray->mTerminationObject)
662                        contribution += vc->GetPvs().AddSample(ray->mTerminationObject);
663
664                if (ray->mOriginObject)
665                        contribution += vc->GetPvs().AddSample(ray->mOriginObject);
666
667                if (contribution)
668                {
669                        sampleContributions += contribution;
670                        ++ contributingSamples;
671                }
672                //leaf->mVssRays.push_back(ray);
673        }
674}
675
676void VspBspTree::SortSplitCandidates(const RayInfoContainer &rays, const int axis)
677{
678        mSplitCandidates->clear();
679
680        int requestedSize = 2 * (int)(rays.size());
681        // creates a sorted split candidates array
682        if (mSplitCandidates->capacity() > 500000 &&
683                requestedSize < (int)(mSplitCandidates->capacity()/10) )
684        {
685        delete mSplitCandidates;
686                mSplitCandidates = new vector<SortableEntry>;
687        }
688
689        mSplitCandidates->reserve(requestedSize);
690
691        // insert all queries
692        for(RayInfoContainer::const_iterator ri = rays.begin(); ri < rays.end(); ++ ri)
693        {
694                bool positive = (*ri).mRay->HasPosDir(axis);
695                mSplitCandidates->push_back(SortableEntry(positive ? SortableEntry::ERayMin : SortableEntry::ERayMax,
696                                                                                                  (*ri).ExtrapOrigin(axis), (*ri).mRay));
697
698                mSplitCandidates->push_back(SortableEntry(positive ? SortableEntry::ERayMax : SortableEntry::ERayMin,
699                                                                                                  (*ri).ExtrapTermination(axis), (*ri).mRay));
700        }
701
702        stable_sort(mSplitCandidates->begin(), mSplitCandidates->end());
703}
704
705float VspBspTree::BestCostRatioHeuristics(const RayInfoContainer &rays,
706                                                                                  const AxisAlignedBox3 &box,
707                                                                                  const int pvsSize,
708                                                                                  const int &axis,
709                                          float &position)
710{
711        int raysBack;
712        int raysFront;
713        int pvsBack;
714        int pvsFront;
715
716        SortSplitCandidates(rays, axis);
717
718        // go through the lists, count the number of objects left and right
719        // and evaluate the following cost funcion:
720        // C = ct_div_ci  + (ql*rl + qr*rr)/queries
721
722        int rl = 0, rr = (int)rays.size();
723        int pl = 0, pr = pvsSize;
724
725        float minBox = box.Min(axis);
726        float maxBox = box.Max(axis);
727        float sizeBox = maxBox - minBox;
728
729        float minBand = minBox + 0.1f*(maxBox - minBox);
730        float maxBand = minBox + 0.9f*(maxBox - minBox);
731
732        float sum = rr*sizeBox;
733        float minSum = 1e20f;
734
735        Intersectable::NewMail();
736
737        // set all object as belonging to the fron pvs
738        for(RayInfoContainer::const_iterator ri = rays.begin(); ri != rays.end(); ++ ri)
739        {
740                if ((*ri).mRay->IsActive())
741                {
742                        Intersectable *object = (*ri).mRay->mTerminationObject;
743
744                        if (object)
745                        {
746                                if (!object->Mailed())
747                                {
748                                        object->Mail();
749                                        object->mCounter = 1;
750                                }
751                                else
752                                        ++ object->mCounter;
753                        }
754                }
755        }
756
757        Intersectable::NewMail();
758
759        for (vector<SortableEntry>::const_iterator ci = mSplitCandidates->begin();
760                ci < mSplitCandidates->end(); ++ ci)
761        {
762                VssRay *ray;
763
764                switch ((*ci).type)
765                {
766                        case SortableEntry::ERayMin:
767                                {
768                                        ++ rl;
769                                        ray = (VssRay *) (*ci).ray;
770
771                                        Intersectable *object = ray->mTerminationObject;
772
773                                        if (object && !object->Mailed())
774                                        {
775                                                object->Mail();
776                                                ++ pl;
777                                        }
778                                        break;
779                                }
780                        case SortableEntry::ERayMax:
781                                {
782                                        -- rr;
783                                        ray = (VssRay *) (*ci).ray;
784
785                                        Intersectable *object = ray->mTerminationObject;
786
787                                        if (object)
788                                        {
789                                                if (-- object->mCounter == 0)
790                                                        -- pr;
791                                        }
792
793                                        break;
794                                }
795                }
796
797                // Note: sufficient to compare size of bounding boxes of front and back side?
798                if ((*ci).value > minBand && (*ci).value < maxBand)
799                {
800                        sum = pl*((*ci).value - minBox) + pr*(maxBox - (*ci).value);
801
802                        //  cout<<"pos="<<(*ci).value<<"\t q=("<<ql<<","<<qr<<")\t r=("<<rl<<","<<rr<<")"<<endl;
803                        // cout<<"cost= "<<sum<<endl;
804
805                        if (sum < minSum)
806                        {
807                                minSum = sum;
808                                position = (*ci).value;
809
810                                raysBack = rl;
811                                raysFront = rr;
812
813                                pvsBack = pl;
814                                pvsFront = pr;
815
816                        }
817                }
818        }
819
820        float oldCost = (float)pvsSize;
821        float newCost = mCtDivCi + minSum / sizeBox;
822        float ratio = newCost / oldCost;
823
824        //Debug << "costRatio=" << ratio << " pos=" << position << " t=" << (position - minBox) / (maxBox - minBox)
825        //     <<"\t q=(" << queriesBack << "," << queriesFront << ")\t r=(" << raysBack << "," << raysFront << ")" << endl;
826
827        return ratio;
828}
829
830
831float VspBspTree::SelectAxisAlignedPlane(Plane3 &plane,
832                                                                                 const VspBspTraversalData &tData,
833                                                                                 int &axis,
834                                                                                 BspNodeGeometry **frontGeom,
835                                                                                 BspNodeGeometry **backGeom,
836                                                                                 float &frontArea,
837                                                                                 float &backArea,
838                                                                                 bool useKdSplit)
839{
840        const bool useCostHeuristics = false;
841
842        //-- regular split
843        float nPosition[3];
844        float nCostRatio[3];
845        float nFrontArea[3];
846        float nBackArea[3];
847
848        BspNodeGeometry *nFrontGeom[3];
849        BspNodeGeometry *nBackGeom[3];
850
851        int bestAxis = -1;
852
853        // create bounding box of node extent
854        AxisAlignedBox3 box;
855        box.Initialize();
856        //TODO: for kd split geometry already is box
857        if (1 && mPvsUseArea)
858        {
859                PolygonContainer::const_iterator it, it_end = tData.mGeometry->mPolys.end();
860
861                for(it = tData.mGeometry->mPolys.begin(); it < it_end; ++ it)
862                        box.Include(*(*it));
863        }
864        else
865        {
866                RayInfoContainer::const_iterator ri, ri_end = tData.mRays->end();
867
868                for(ri = tData.mRays->begin(); ri < ri_end; ++ ri)
869                        box.Include((*ri).ExtrapTermination());
870        }
871
872        int sAxis = box.Size().DrivingAxis();
873
874        for (axis = 0; axis < 3; ++ axis)
875        {
876                nFrontGeom[axis] = new BspNodeGeometry();
877                nBackGeom[axis] = new BspNodeGeometry();
878
879                if (!mOnlyDrivingAxis || axis == sAxis)
880                {
881                        if (!useCostHeuristics)
882                        {
883                                nPosition[axis] = (box.Min()[axis] + box.Max()[axis]) * 0.5f;
884                                Vector3 normal(0,0,0); normal[axis] = 1;
885                               
886                                if (useKdSplit)
887                                {
888                                        nCostRatio[axis] = EvalAxisAlignedSplitCost(tData,
889                                                                                                                                box,
890                                                                                                                                axis,
891                                                                                                                                nPosition[axis]);
892                                       
893                                        Vector3 pos = box.Max(); pos[axis] = nPosition[axis];
894                                        AxisAlignedBox3 frontBox(box.Min(), pos);
895                                        frontBox.ExtractPolys(nFrontGeom[axis]->mPolys);
896                                       
897                                        pos = box.Min(); pos[axis] = nPosition[axis];
898                                        AxisAlignedBox3 backBox(pos, box.Max());
899                                        backBox.ExtractPolys(nBackGeom[axis]->mPolys);
900
901                                }
902                                else
903                                {
904                                        nCostRatio[axis] = SplitPlaneCost(Plane3(normal, nPosition[axis]),
905                                                                                                          tData, *nFrontGeom[axis], *nBackGeom[axis],
906                                                                                                          nFrontArea[axis], nBackArea[axis]);
907                                }
908                        }
909                        else
910                        {
911                                nCostRatio[axis] =
912                                        BestCostRatioHeuristics(*tData.mRays,
913                                                                                    box,
914                                                                                        tData.mPvs,
915                                                                                        axis,
916                                                                                        nPosition[axis]);
917                        }
918
919                        if (bestAxis == -1)
920                        {
921                                bestAxis = axis;
922                        }
923
924                        else if (nCostRatio[axis] < nCostRatio[bestAxis])
925                        {
926                                bestAxis = axis;
927                        }
928
929                }
930        }
931
932        //-- assign values
933        axis = bestAxis;
934        frontArea = nFrontArea[bestAxis];
935        backArea = nBackArea[bestAxis];
936
937        // assign best split nodes geometry
938        *frontGeom = nFrontGeom[bestAxis];
939        *backGeom = nBackGeom[bestAxis];
940
941        // and delete other geometry
942        delete nFrontGeom[(bestAxis + 1) % 3];
943        delete nBackGeom[(bestAxis + 2) % 3];
944
945        //-- split plane
946    Vector3 normal(0,0,0); normal[bestAxis] = 1;
947        plane = Plane3(normal, nPosition[bestAxis]);
948
949        return nCostRatio[bestAxis];
950}
951
952
953float VspBspTree::EvalAxisAlignedSplitCost(const VspBspTraversalData &data,
954                                                                                   const AxisAlignedBox3 &box,
955                                                                                   const int axis,
956                                                                                   const float &position) const
957{
958        int pvsFront = 0;
959        int pvsBack = 0;
960        int pvsTotal = 0;
961
962        // create unique ids for pvs heuristics
963        GenerateUniqueIdsForPvs();
964
965        const int pvsSize = data.mPvs;
966
967        // this is the main ray classification loop!
968        for(RayInfoContainer::iterator ri = data.mRays->begin();
969                ri != data.mRays->end(); ++ ri)
970        {
971                if (!(*ri).mRay->IsActive())
972                        continue;
973
974                // determine the side of this ray with respect to the plane
975                float t;
976                int side = (*ri).ComputeRayIntersection(axis, position, t);
977       
978                AddObjToPvs((*ri).mRay->mTerminationObject, side, pvsFront, pvsBack, pvsTotal);
979                AddObjToPvs((*ri).mRay->mOriginObject, side, pvsFront, pvsBack, pvsTotal);
980        }
981
982        //-- only one of these options should be one
983        //-- pvs + probability heuristics
984        float pBack, pFront, pOverall;
985
986        if (1)
987        {
988                // box length substitute for probability
989                const float minBox = box.Min(axis);
990                const float maxBox = box.Max(axis);
991
992                pBack = position - minBox;
993                pFront = maxBox - position;
994                pOverall = maxBox - minBox;
995        }
996        else //-- area substitute for probability
997        {
998                pOverall = box.SurfaceArea();
999
1000                const bool useMidSplit = true;
1001                //const bool useMidSplit = false;       
1002                                       
1003                //-- simplified computation for mid split
1004                const int axis2 = (axis + 1) % 3;
1005                const int axis3 = (axis + 2) % 3;
1006
1007                const float faceArea =
1008                        (box.Max(axis2) - box.Min(axis2)) *
1009                        (box.Max(axis3) - box.Min(axis3));
1010
1011                pBack = pFront = pOverall * 0.5f + faceArea;
1012        }
1013
1014        //Debug << axis << " " << pvsSize << " " << pvsBack << " " << pvsFront << endl;
1015        //Debug << pFront << " " << pBack << " " << pOverall << endl;
1016
1017        // float sum = raysBack*(position - minBox) + raysFront*(maxBox - position);
1018        const float newCost = pvsBack * pBack + pvsFront * pFront;
1019        //  float oldCost = leaf->mRays.size();
1020        const float oldCost = (float)pvsSize * pOverall;
1021
1022        return  (mCtDivCi + newCost) / oldCost;
1023}
1024
1025
1026
1027bool VspBspTree::SelectPlane(Plane3 &bestPlane,
1028                                                         BspLeaf *leaf,
1029                                                         VspBspTraversalData &data,
1030                                                         VspBspTraversalData &frontData,
1031                                                         VspBspTraversalData &backData)
1032{
1033        // simplest strategy: just take next polygon
1034        if (mSplitPlaneStrategy & RANDOM_POLYGON)
1035        {
1036        if (!data.mPolygons->empty())
1037                {
1038                        const int randIdx =
1039                                (int)RandomValue(0, (Real)((int)data.mPolygons->size() - 1));
1040                        Polygon3 *nextPoly = (*data.mPolygons)[randIdx];
1041
1042                        bestPlane = nextPoly->GetSupportingPlane();
1043                        return true;
1044                }
1045        }
1046
1047        //-- use heuristics to find appropriate plane
1048
1049        // intermediate plane
1050        Plane3 plane;
1051        float lowestCost = MAX_FLOAT;
1052       
1053        // decides if the first few splits should be only axisAligned
1054        const bool onlyAxisAligned  =
1055                (mSplitPlaneStrategy & AXIS_ALIGNED) &&
1056                ((int)data.mRays->size() > mTermMinRaysForAxisAligned) &&
1057                ((int)data.GetAvgRayContribution() < mTermMaxRayContriForAxisAligned);
1058       
1059        const int limit = onlyAxisAligned ? 0 :
1060                Min((int)data.mPolygons->size(), mMaxPolyCandidates);
1061
1062        float candidateCost;
1063
1064        int maxIdx = (int)data.mPolygons->size();
1065
1066        for (int i = 0; i < limit; ++ i)
1067        {
1068                //-- assure that no index is taken twice
1069                const int candidateIdx = (int)RandomValue(0, (Real)(-- maxIdx));
1070                Polygon3 *poly = (*data.mPolygons)[candidateIdx];
1071
1072                // swap candidate to the end to avoid testing same plane
1073                std::swap((*data.mPolygons)[maxIdx], (*data.mPolygons)[candidateIdx]);
1074                //Polygon3 *poly = (*data.mPolygons)[(int)RandomValue(0, (int)polys.size() - 1)];
1075
1076                // evaluate current candidate
1077                BspNodeGeometry fGeom, bGeom;
1078                float fArea, bArea;
1079                plane = poly->GetSupportingPlane();
1080                candidateCost = SplitPlaneCost(plane, data, fGeom, bGeom, fArea, bArea);
1081               
1082                if (candidateCost < lowestCost)
1083                {
1084                        bestPlane = plane;
1085                        lowestCost = candidateCost;
1086                }
1087        }
1088
1089        //-- evaluate axis aligned splits
1090        int axis;
1091        BspNodeGeometry *fGeom, *bGeom;
1092        float fArea, bArea;
1093
1094        candidateCost = SelectAxisAlignedPlane(plane, data, axis,
1095                                                                                   &fGeom, &bGeom,
1096                                                                                   fArea, bArea,
1097                                                                                   onlyAxisAligned);     
1098
1099        if (candidateCost < lowestCost)
1100        {       
1101                bestPlane = plane;
1102                lowestCost = candidateCost;
1103
1104                // assign already computed values
1105                // we can do this because we always save the
1106                // computed values from the axis aligned splits
1107                frontData.mGeometry = fGeom;
1108                backData.mGeometry = bGeom;
1109                frontData.mArea = fArea;
1110                backData.mArea = bArea;
1111
1112                //! error also computed if cost ratio is missed
1113                ++ mStat.splits[axis];
1114        }
1115        else
1116        {
1117                DEL_PTR(fGeom);
1118                DEL_PTR(bGeom);
1119        }
1120
1121#ifdef _DEBUG
1122        Debug << "plane lowest cost: " << lowestCost << endl;
1123#endif
1124
1125        // cost ratio miss
1126        if (lowestCost > mTermMaxCostRatio)
1127                return false;
1128
1129        return true;
1130}
1131
1132
1133Plane3 VspBspTree::ChooseCandidatePlane(const RayInfoContainer &rays) const
1134{
1135        const int candidateIdx = (int)RandomValue(0, (Real)((int)rays.size() - 1));
1136
1137        const Vector3 minPt = rays[candidateIdx].ExtrapOrigin();
1138        const Vector3 maxPt = rays[candidateIdx].ExtrapTermination();
1139
1140        const Vector3 pt = (maxPt + minPt) * 0.5;
1141        const Vector3 normal = Normalize(rays[candidateIdx].mRay->GetDir());
1142
1143        return Plane3(normal, pt);
1144}
1145
1146
1147Plane3 VspBspTree::ChooseCandidatePlane2(const RayInfoContainer &rays) const
1148{
1149        Vector3 pt[3];
1150
1151        int idx[3];
1152        int cmaxT = 0;
1153        int cminT = 0;
1154        bool chooseMin = false;
1155
1156        for (int j = 0; j < 3; ++ j)
1157        {
1158                idx[j] = (int)RandomValue(0, (Real)((int)rays.size() * 2 - 1));
1159
1160                if (idx[j] >= (int)rays.size())
1161                {
1162                        idx[j] -= (int)rays.size();
1163
1164                        chooseMin = (cminT < 2);
1165                }
1166                else
1167                        chooseMin = (cmaxT < 2);
1168
1169                RayInfo rayInf = rays[idx[j]];
1170                pt[j] = chooseMin ? rayInf.ExtrapOrigin() : rayInf.ExtrapTermination();
1171        }
1172
1173        return Plane3(pt[0], pt[1], pt[2]);
1174}
1175
1176Plane3 VspBspTree::ChooseCandidatePlane3(const RayInfoContainer &rays) const
1177{
1178        Vector3 pt[3];
1179
1180        int idx1 = (int)RandomValue(0, (Real)((int)rays.size() - 1));
1181        int idx2 = (int)RandomValue(0, (Real)((int)rays.size() - 1));
1182
1183        // check if rays different
1184        if (idx1 == idx2)
1185                idx2 = (idx2 + 1) % (int)rays.size();
1186
1187        const RayInfo ray1 = rays[idx1];
1188        const RayInfo ray2 = rays[idx2];
1189
1190        // normal vector of the plane parallel to both lines
1191        const Vector3 norm = Normalize(CrossProd(ray1.mRay->GetDir(), ray2.mRay->GetDir()));
1192
1193        // vector from line 1 to line 2
1194        const Vector3 vd = ray2.ExtrapOrigin() - ray1.ExtrapOrigin();
1195
1196        // project vector on normal to get distance
1197        const float dist = DotProd(vd, norm);
1198
1199        // point on plane lies halfway between the two planes
1200        const Vector3 planePt = ray1.ExtrapOrigin() + norm * dist * 0.5;
1201
1202        return Plane3(norm, planePt);
1203}
1204
1205
1206inline void VspBspTree::GenerateUniqueIdsForPvs()
1207{
1208        Intersectable::NewMail(); sBackId = ViewCell::sMailId;
1209        Intersectable::NewMail(); sFrontId = ViewCell::sMailId;
1210        Intersectable::NewMail(); sFrontAndBackId = ViewCell::sMailId;
1211}
1212
1213
1214float VspBspTree::SplitPlaneCost(const Plane3 &candidatePlane,
1215                                                                 const VspBspTraversalData &data,
1216                                                                 BspNodeGeometry &geomFront,
1217                                                                 BspNodeGeometry &geomBack,
1218                                                                 float &areaFront,
1219                                                                 float &areaBack) const
1220{
1221        float cost = 0;
1222
1223        float sumBalancedRays = 0;
1224        float sumRaySplits = 0;
1225
1226        int pvsFront = 0;
1227        int pvsBack = 0;
1228
1229        // probability that view point lies in back / front node
1230        float pOverall = 0;
1231        float pFront = 0;
1232        float pBack = 0;
1233
1234        int raysFront = 0;
1235        int raysBack = 0;
1236        int totalPvs = 0;
1237
1238        int limit;
1239        bool useRand;
1240
1241        // choose test rays randomly if too much
1242        if ((int)data.mRays->size() > mMaxTests)
1243        {
1244                useRand = true;
1245                limit = mMaxTests;
1246        }
1247        else
1248        {
1249                useRand = false;
1250                limit = (int)data.mRays->size();
1251        }
1252       
1253        for (int i = 0; i < limit; ++ i)
1254        {
1255                const int testIdx = useRand ?
1256                        (int)RandomValue(0, (Real)((int)data.mRays->size() - 1)) : i;
1257                RayInfo rayInf = (*data.mRays)[testIdx];
1258
1259                float t;
1260                VssRay *ray = rayInf.mRay;
1261                const int cf = rayInf.ComputeRayIntersection(candidatePlane, t);
1262
1263        if (cf >= 0)
1264                        ++ raysFront;
1265                if (cf <= 0)
1266                        ++ raysBack;
1267
1268                if (mSplitPlaneStrategy & LEAST_RAY_SPLITS)
1269                {
1270                        sumBalancedRays += cf;
1271                }
1272
1273                if (mSplitPlaneStrategy & BALANCED_RAYS)
1274                {
1275                        if (cf == 0)
1276                                ++ sumRaySplits;
1277                }
1278
1279                if (mSplitPlaneStrategy & PVS)
1280                {
1281                        // find front and back pvs for origing and termination object
1282                        AddObjToPvs(ray->mTerminationObject, cf, pvsFront, pvsBack, totalPvs);
1283                        AddObjToPvs(ray->mOriginObject, cf, pvsFront, pvsBack, totalPvs);
1284                }
1285        }
1286
1287        const float raysSize = (float)data.mRays->size() + Limits::Small;
1288        if (mSplitPlaneStrategy & PVS)
1289        {
1290                // create unique ids for pvs heuristics
1291                GenerateUniqueIdsForPvs();
1292
1293                if (mPvsUseArea) // use front and back cell areas to approximate volume
1294                {
1295                        // construct child geometry with regard to the candidate split plane
1296                        data.mGeometry->SplitGeometry(geomFront,
1297                                                                                  geomBack,
1298                                                                                  candidatePlane,
1299                                                                                  mBox,
1300                                                                                  mEpsilon);
1301#if 0
1302                        AxisAlignedBox3 fbox;
1303                        AxisAlignedBox3 bbox;
1304                        AxisAlignedBox3 box;
1305                        fbox.Initialize();
1306                        bbox.Initialize();
1307                        box.Initialize();
1308
1309                        geomFront.IncludeInBox(fbox);
1310                        geomBack.IncludeInBox(bbox);
1311                        data.mGeometry->IncludeInBox(box);
1312
1313                        areaFront = fbox.GetVolume();
1314                        areaBack = bbox.GetVolume();
1315
1316                        pOverall = box.GetVolume();
1317#else
1318                        areaFront = geomFront.GetArea();
1319                        areaBack = geomBack.GetArea();
1320
1321                        pOverall = data.mArea;
1322#endif
1323                        pBack = areaBack;
1324                        pFront = areaFront;
1325                }
1326                else // use number of rays to approximate volume
1327                {
1328                        pOverall = (float)data.mRays->size();
1329                        pFront = (float)raysFront;
1330                        pBack = (float)raysBack;
1331                }
1332        }
1333
1334
1335        if (mSplitPlaneStrategy & LEAST_RAY_SPLITS)
1336                cost += mLeastRaySplitsFactor * sumRaySplits / raysSize;
1337
1338        if (mSplitPlaneStrategy & BALANCED_RAYS)
1339                cost += mBalancedRaysFactor * fabs(sumBalancedRays) / raysSize;
1340
1341        // pvs criterium
1342        if (mSplitPlaneStrategy & PVS)
1343        {
1344                const float oldCost = pOverall * (float)totalPvs + Limits::Small;
1345                cost += mPvsFactor * (pvsFront * pFront + pvsBack * pBack) / oldCost;
1346        }
1347
1348#ifdef _DEBUG
1349        Debug << "totalpvs: " << data.mPvs << " ptotal: " << pOverall
1350                  << " frontpvs: " << pvsFront << " pFront: " << pFront
1351                  << " backpvs: " << pvsBack << " pBack: " << pBack << endl << endl;
1352#endif
1353
1354        // normalize cost by sum of linear factors
1355        return cost / (float)mCostNormalizer;
1356}
1357
1358
1359void VspBspTree::AddObjToPvs(Intersectable *obj,
1360                                                         const int cf,
1361                                                         int &frontPvs,
1362                                                         int &backPvs,
1363                                                         int &totalPvs) const
1364{
1365        if (!obj)
1366                return;
1367       
1368        if ((obj->mMailbox != sFrontId) &&
1369                (obj->mMailbox != sBackId) &&
1370                (obj->mMailbox != sFrontAndBackId))
1371        {
1372                ++ totalPvs;
1373        }
1374
1375        // TODO: does this really belong to no pvs?
1376        //if (cf == Ray::COINCIDENT) return;
1377
1378        // object belongs to both PVS
1379        if (cf >= 0)
1380        {
1381                if ((obj->mMailbox != sFrontId) &&
1382                        (obj->mMailbox != sFrontAndBackId))
1383                {
1384                        ++ frontPvs;
1385               
1386                        if (obj->mMailbox == sBackId)
1387                                obj->mMailbox = sFrontAndBackId;
1388                        else
1389                                obj->mMailbox = sFrontId;
1390                }
1391        }
1392
1393        if (cf <= 0)
1394        {
1395                if ((obj->mMailbox != sBackId) &&
1396                        (obj->mMailbox != sFrontAndBackId))
1397                {
1398                        ++ backPvs;
1399               
1400                        if (obj->mMailbox == sFrontId)
1401                                obj->mMailbox = sFrontAndBackId;
1402                        else
1403                                obj->mMailbox = sBackId;
1404                }
1405        }
1406}
1407
1408
1409void VspBspTree::CollectLeaves(vector<BspLeaf *> &leaves,
1410                                                           const bool onlyUnmailed,
1411                                                           const int maxPvsSize) const
1412{
1413        stack<BspNode *> nodeStack;
1414        nodeStack.push(mRoot);
1415
1416        while (!nodeStack.empty())
1417        {
1418                BspNode *node = nodeStack.top();
1419                nodeStack.pop();
1420               
1421                if (node->IsLeaf())
1422                {
1423                        // test if this leaf is in valid view space
1424                        BspLeaf *leaf = dynamic_cast<BspLeaf *>(node);
1425                        if (leaf->TreeValid() &&
1426                                (!onlyUnmailed || !leaf->Mailed()) &&
1427                                ((maxPvsSize < 0) || (leaf->GetViewCell()->GetPvs().GetSize() <= maxPvsSize)))
1428                        {
1429                                leaves.push_back(leaf);
1430                        }
1431                }
1432                else
1433                {
1434                        BspInterior *interior = dynamic_cast<BspInterior *>(node);
1435
1436                        nodeStack.push(interior->GetBack());
1437                        nodeStack.push(interior->GetFront());
1438                }
1439        }
1440}
1441
1442
1443AxisAlignedBox3 VspBspTree::GetBoundingBox() const
1444{
1445        return mBox;
1446}
1447
1448
1449BspNode *VspBspTree::GetRoot() const
1450{
1451        return mRoot;
1452}
1453
1454
1455void VspBspTree::EvaluateLeafStats(const VspBspTraversalData &data)
1456{
1457        // the node became a leaf -> evaluate stats for leafs
1458        BspLeaf *leaf = dynamic_cast<BspLeaf *>(data.mNode);
1459
1460        // store maximal and minimal depth
1461        if (data.mDepth > mStat.maxDepth)
1462                mStat.maxDepth = data.mDepth;
1463
1464        if (data.mPvs > mStat.maxPvs)
1465                mStat.maxPvs = data.mPvs;
1466       
1467        if (data.mDepth < mStat.minDepth)
1468                mStat.minDepth = data.mDepth;
1469
1470        if (data.mDepth >= mTermMaxDepth)
1471                ++ mStat.maxDepthNodes;
1472        // accumulate rays to compute rays /  leaf
1473        mStat.accumRays += (int)data.mRays->size();
1474
1475        if (data.mPvs < mTermMinPvs)
1476                ++ mStat.minPvsNodes;
1477
1478        if ((int)data.mRays->size() < mTermMinRays)
1479                ++ mStat.minRaysNodes;
1480
1481        if (data.GetAvgRayContribution() > mTermMaxRayContribution)
1482                ++ mStat.maxRayContribNodes;
1483
1484        if (data.mArea <= mTermMinArea)
1485                ++ mStat.minAreaNodes;
1486       
1487        // accumulate depth to compute average depth
1488        mStat.accumDepth += data.mDepth;
1489
1490#ifdef _DEBUG
1491        Debug << "BSP stats: "
1492                  << "Depth: " << data.mDepth << " (max: " << mTermMaxDepth << "), "
1493                  << "PVS: " << data.mPvs << " (min: " << mTermMinPvs << "), "
1494                  << "Area: " << data.mArea << " (min: " << mTermMinArea << "), "
1495                  << "#rays: " << (int)data.mRays->size() << " (max: " << mTermMinRays << "), "
1496                  << "#pvs: " << leaf->GetViewCell()->GetPvs().GetSize() << "=, "
1497                  << "#avg ray contrib (pvs): " << (float)data.mPvs / (float)data.mRays->size() << endl;
1498#endif
1499}
1500
1501int VspBspTree::CastRay(Ray &ray)
1502{
1503        int hits = 0;
1504
1505        stack<BspRayTraversalData> tStack;
1506
1507        float maxt, mint;
1508
1509        if (!mBox.GetRaySegment(ray, mint, maxt))
1510                return 0;
1511
1512        Intersectable::NewMail();
1513
1514        Vector3 entp = ray.Extrap(mint);
1515        Vector3 extp = ray.Extrap(maxt);
1516
1517        BspNode *node = mRoot;
1518        BspNode *farChild = NULL;
1519
1520        while (1)
1521        {
1522                if (!node->IsLeaf())
1523                {
1524                        BspInterior *in = dynamic_cast<BspInterior *>(node);
1525
1526                        Plane3 splitPlane = in->GetPlane();
1527                        const int entSide = splitPlane.Side(entp);
1528                        const int extSide = splitPlane.Side(extp);
1529
1530                        if (entSide < 0)
1531                        {
1532                                node = in->GetBack();
1533
1534                                if(extSide <= 0) // plane does not split ray => no far child
1535                                        continue;
1536
1537                                farChild = in->GetFront(); // plane splits ray
1538
1539                        } else if (entSide > 0)
1540                        {
1541                                node = in->GetFront();
1542
1543                                if (extSide >= 0) // plane does not split ray => no far child
1544                                        continue;
1545
1546                                farChild = in->GetBack(); // plane splits ray
1547                        }
1548                        else // ray and plane are coincident
1549                        {
1550                                // WHAT TO DO IN THIS CASE ?
1551                                //break;
1552                                node = in->GetFront();
1553                                continue;
1554                        }
1555
1556                        // push data for far child
1557                        tStack.push(BspRayTraversalData(farChild, extp, maxt));
1558
1559                        // find intersection of ray segment with plane
1560                        float t;
1561                        extp = splitPlane.FindIntersection(ray.GetLoc(), extp, &t);
1562                        maxt *= t;
1563
1564                } else // reached leaf => intersection with view cell
1565                {
1566                        BspLeaf *leaf = dynamic_cast<BspLeaf *>(node);
1567
1568                        if (!leaf->GetViewCell()->Mailed())
1569                        {
1570                                //ray.bspIntersections.push_back(Ray::VspBspIntersection(maxt, leaf));
1571                                leaf->GetViewCell()->Mail();
1572                                ++ hits;
1573                        }
1574
1575                        //-- fetch the next far child from the stack
1576                        if (tStack.empty())
1577                                break;
1578
1579                        entp = extp;
1580                        mint = maxt; // NOTE: need this?
1581
1582                        if (ray.GetType() == Ray::LINE_SEGMENT && mint > 1.0f)
1583                                break;
1584
1585                        BspRayTraversalData &s = tStack.top();
1586
1587                        node = s.mNode;
1588                        extp = s.mExitPoint;
1589                        maxt = s.mMaxT;
1590
1591                        tStack.pop();
1592                }
1593        }
1594
1595        return hits;
1596}
1597
1598
1599void VspBspTree::CollectViewCells(ViewCellContainer &viewCells) const
1600{
1601        ViewCell::NewMail();
1602        CollectViewCells(mRoot, viewCells, true);
1603}
1604
1605
1606
1607void VspBspTree::CheckValidy()
1608{
1609        stack<BspNode *> nodeStack;
1610
1611        if (!mRoot)
1612                return;
1613
1614        nodeStack.push(mRoot);
1615       
1616        while (!nodeStack.empty())
1617        {
1618                BspNode *node = nodeStack.top();
1619                nodeStack.pop();
1620               
1621                if (node->IsLeaf())
1622                {
1623                        BspLeaf *leaf = dynamic_cast<BspLeaf *>(node);
1624
1625                        if (node->TreeValid())
1626                        {
1627                                BspViewCell *viewCell = dynamic_cast<BspLeaf *>(node)->GetViewCell();
1628                       
1629                                if (viewCell->GetPvs().GetSize() > mMaxPvs)
1630                                {
1631                                        while (!viewCell->mLeaves.empty())
1632                                        {
1633                                                BspLeaf *l = viewCell->mLeaves.back();
1634                                                l->SetTreeValid(false);
1635                                                PropagateUpValidity(l);
1636
1637                                                l->SetViewCell(GetOrCreateOutOfBoundsCell());
1638                                                viewCell->mLeaves.pop_back();
1639                                        }
1640
1641                                        mOutOfBoundsCell->GetPvs().AddPvs(viewCell->GetPvs());
1642
1643                                        DEL_PTR(viewCell);
1644                                }
1645                        }
1646                }
1647                else
1648                {
1649                        BspInterior *interior = dynamic_cast<BspInterior *>(node);
1650               
1651                        nodeStack.push(interior->GetFront());
1652                        nodeStack.push(interior->GetBack());
1653                }
1654        }
1655}
1656
1657void VspBspTree::CollectViewCells(BspNode *root,
1658                                                                  ViewCellContainer &viewCells,
1659                                                                  bool onlyUnmailed) const
1660{
1661        stack<BspNode *> nodeStack;
1662
1663        if (!root)
1664                return;
1665
1666        nodeStack.push(root);
1667       
1668        while (!nodeStack.empty())
1669        {
1670                BspNode *node = nodeStack.top();
1671                nodeStack.pop();
1672               
1673                if (node->IsLeaf())
1674                {
1675                        if (!mShowInvalidSpace && node->TreeValid())
1676                        {
1677                                ViewCell *viewCell = dynamic_cast<BspLeaf *>(node)->GetViewCell();
1678                       
1679                                if (!onlyUnmailed || !viewCell->Mailed())
1680                                {
1681                                        viewCell->Mail();
1682                                        viewCells.push_back(viewCell);
1683                                }
1684                        }
1685                }
1686                else
1687                {
1688                        BspInterior *interior = dynamic_cast<BspInterior *>(node);
1689               
1690                        nodeStack.push(interior->GetFront());
1691                        nodeStack.push(interior->GetBack());
1692                }
1693        }
1694}
1695
1696
1697float VspBspTree::AccumulatedRayLength(const RayInfoContainer &rays) const
1698{
1699        float len = 0;
1700
1701        RayInfoContainer::const_iterator it, it_end = rays.end();
1702
1703        for (it = rays.begin(); it != it_end; ++ it)
1704                len += (*it).SegmentLength();
1705
1706        return len;
1707}
1708
1709
1710int VspBspTree::SplitRays(const Plane3 &plane,
1711                                                  RayInfoContainer &rays,
1712                                                  RayInfoContainer &frontRays,
1713                                                  RayInfoContainer &backRays)
1714{
1715        int splits = 0;
1716
1717        while (!rays.empty())
1718        {
1719                RayInfo bRay = rays.back();
1720                rays.pop_back();
1721
1722                VssRay *ray = bRay.mRay;
1723                float t;
1724
1725                // get classification and receive new t
1726                const int cf = bRay.ComputeRayIntersection(plane, t);
1727
1728                switch (cf)
1729                {
1730                case -1:
1731                        backRays.push_back(bRay);
1732                        break;
1733                case 1:
1734                        frontRays.push_back(bRay);
1735                        break;
1736                case 0:
1737                        {
1738                                //-- split ray
1739                                //-- test if start point behind or in front of plane
1740                                const int side = plane.Side(bRay.ExtrapOrigin());
1741
1742                                if (side <= 0)
1743                                {
1744                                        backRays.push_back(RayInfo(ray, bRay.GetMinT(), t));
1745                                        frontRays.push_back(RayInfo(ray, t, bRay.GetMaxT()));
1746                                }
1747                                else
1748                                {
1749                                        frontRays.push_back(RayInfo(ray, bRay.GetMinT(), t));
1750                                        backRays.push_back(RayInfo(ray, t, bRay.GetMaxT()));
1751                                }
1752                        }
1753                        break;
1754                default:
1755                        Debug << "Should not come here" << endl;
1756                        break;
1757                }
1758        }
1759
1760        return splits;
1761}
1762
1763
1764void VspBspTree::ExtractHalfSpaces(BspNode *n, vector<Plane3> &halfSpaces) const
1765{
1766        BspNode *lastNode;
1767
1768        do
1769        {
1770                lastNode = n;
1771
1772                // want to get planes defining geometry of this node => don't take
1773                // split plane of node itself
1774                n = n->GetParent();
1775
1776                if (n)
1777                {
1778                        BspInterior *interior = dynamic_cast<BspInterior *>(n);
1779                        Plane3 halfSpace = dynamic_cast<BspInterior *>(interior)->GetPlane();
1780
1781            if (interior->GetFront() != lastNode)
1782                                halfSpace.ReverseOrientation();
1783
1784                        halfSpaces.push_back(halfSpace);
1785                }
1786        }
1787        while (n);
1788}
1789
1790
1791void VspBspTree::ConstructGeometry(BspNode *n,
1792                                                                   BspNodeGeometry &geom) const
1793{
1794        vector<Plane3> halfSpaces;
1795        ExtractHalfSpaces(n, halfSpaces);
1796
1797        PolygonContainer candidatePolys;
1798
1799        // bounded planes are added to the polygons (reverse polygons
1800        // as they have to be outfacing
1801        for (int i = 0; i < (int)halfSpaces.size(); ++ i)
1802        {
1803                Polygon3 *p = GetBoundingBox().CrossSection(halfSpaces[i]);
1804
1805                if (p->Valid(mEpsilon))
1806                {
1807                        candidatePolys.push_back(p->CreateReversePolygon());
1808                        DEL_PTR(p);
1809                }
1810        }
1811
1812        // add faces of bounding box (also could be faces of the cell)
1813        for (int i = 0; i < 6; ++ i)
1814        {
1815                VertexContainer vertices;
1816
1817                for (int j = 0; j < 4; ++ j)
1818                        vertices.push_back(mBox.GetFace(i).mVertices[j]);
1819
1820                candidatePolys.push_back(new Polygon3(vertices));
1821        }
1822
1823        for (int i = 0; i < (int)candidatePolys.size(); ++ i)
1824        {
1825                // polygon is split by all other planes
1826                for (int j = 0; (j < (int)halfSpaces.size()) && candidatePolys[i]; ++ j)
1827                {
1828                        if (i == j) // polygon and plane are coincident
1829                                continue;
1830
1831                        VertexContainer splitPts;
1832                        Polygon3 *frontPoly, *backPoly;
1833
1834                        const int cf =
1835                                candidatePolys[i]->ClassifyPlane(halfSpaces[j],
1836                                                                                                 mEpsilon);
1837
1838                        switch (cf)
1839                        {
1840                                case Polygon3::SPLIT:
1841                                        frontPoly = new Polygon3();
1842                                        backPoly = new Polygon3();
1843
1844                                        candidatePolys[i]->Split(halfSpaces[j],
1845                                                                                         *frontPoly,
1846                                                                                         *backPoly,
1847                                                                                         mEpsilon);
1848
1849                                        DEL_PTR(candidatePolys[i]);
1850
1851                                        if (frontPoly->Valid(mEpsilon))
1852                                                candidatePolys[i] = frontPoly;
1853                                        else
1854                                                DEL_PTR(frontPoly);
1855
1856                                        DEL_PTR(backPoly);
1857                                        break;
1858                                case Polygon3::BACK_SIDE:
1859                                        DEL_PTR(candidatePolys[i]);
1860                                        break;
1861                                // just take polygon as it is
1862                                case Polygon3::FRONT_SIDE:
1863                                case Polygon3::COINCIDENT:
1864                                default:
1865                                        break;
1866                        }
1867                }
1868
1869                if (candidatePolys[i])
1870                        geom.mPolys.push_back(candidatePolys[i]);
1871        }
1872}
1873
1874
1875void VspBspTree::ConstructGeometry(BspViewCell *vc,
1876                                                                   BspNodeGeometry &vcGeom) const
1877{
1878        vector<BspLeaf *> leaves = vc->mLeaves;
1879        vector<BspLeaf *>::const_iterator it, it_end = leaves.end();
1880
1881        for (it = leaves.begin(); it != it_end; ++ it)
1882                ConstructGeometry(*it, vcGeom);
1883}
1884
1885
1886int VspBspTree::FindNeighbors(BspNode *n, vector<BspLeaf *> &neighbors,
1887                                                   const bool onlyUnmailed) const
1888{
1889        BspNodeGeometry geom;
1890        ConstructGeometry(n, geom);
1891
1892        stack<BspNode *> nodeStack;
1893        nodeStack.push(mRoot);
1894
1895        // split planes from the root to this node
1896        // needed to verify that we found neighbor leaf
1897        vector<Plane3> halfSpaces;
1898        ExtractHalfSpaces(n, halfSpaces);
1899
1900        while (!nodeStack.empty())
1901        {
1902                BspNode *node = nodeStack.top();
1903                nodeStack.pop();
1904
1905                if (node->IsLeaf())
1906
1907                {       // test if this leaf is in valid view space
1908            if (node->TreeValid() &&
1909                                node != n &&
1910                                (!onlyUnmailed || !node->Mailed()))
1911                        {
1912                                // test all planes of current node if candidate really
1913                                // is neighbour
1914                                BspNodeGeometry neighborCandidate;
1915                                ConstructGeometry(node, neighborCandidate);
1916
1917                                bool isAdjacent = true;
1918                                for (int i = 0; (i < halfSpaces.size()) && isAdjacent; ++ i)
1919                                {
1920                                        const int cf =
1921                                                Polygon3::ClassifyPlane(neighborCandidate.mPolys,
1922                                                                                                halfSpaces[i],
1923                                                                                                mEpsilon);
1924
1925                                        if (cf == Polygon3::BACK_SIDE)
1926                                                isAdjacent = false;
1927                                }
1928                                // neighbor was found
1929                                if (isAdjacent)
1930                                        neighbors.push_back(dynamic_cast<BspLeaf *>(node));
1931                        }
1932                }
1933                else
1934                {
1935                        BspInterior *interior = dynamic_cast<BspInterior *>(node);
1936
1937                        const int cf = Polygon3::ClassifyPlane(geom.mPolys,
1938                                                                                                   interior->GetPlane(),
1939                                                                                                   mEpsilon);
1940
1941                        if (cf == Polygon3::FRONT_SIDE)
1942                                nodeStack.push(interior->GetFront());
1943                        else
1944                                if (cf == Polygon3::BACK_SIDE)
1945                                        nodeStack.push(interior->GetBack());
1946                                else
1947                                {
1948                                        // random decision
1949                                        nodeStack.push(interior->GetBack());
1950                                        nodeStack.push(interior->GetFront());
1951                                }
1952                }
1953        }
1954
1955        return (int)neighbors.size();
1956}
1957
1958
1959BspLeaf *VspBspTree::GetRandomLeaf(const Plane3 &halfspace)
1960{
1961    stack<BspNode *> nodeStack;
1962        nodeStack.push(mRoot);
1963
1964        int mask = rand();
1965
1966        while (!nodeStack.empty())
1967        {
1968                BspNode *node = nodeStack.top();
1969                nodeStack.pop();
1970
1971                if (node->IsLeaf())
1972                {
1973                        return dynamic_cast<BspLeaf *>(node);
1974                }
1975                else
1976                {
1977                        BspInterior *interior = dynamic_cast<BspInterior *>(node);
1978                        BspNode *next;
1979                        BspNodeGeometry geom;
1980
1981                        // todo: not very efficient: constructs full cell everytime
1982                        ConstructGeometry(interior, geom);
1983
1984                        const int cf =
1985                                Polygon3::ClassifyPlane(geom.mPolys, halfspace, mEpsilon);
1986
1987                        if (cf == Polygon3::BACK_SIDE)
1988                                next = interior->GetFront();
1989                        else
1990                                if (cf == Polygon3::FRONT_SIDE)
1991                                        next = interior->GetFront();
1992                        else
1993                        {
1994                                // random decision
1995                                if (mask & 1)
1996                                        next = interior->GetBack();
1997                                else
1998                                        next = interior->GetFront();
1999                                mask = mask >> 1;
2000                        }
2001
2002                        nodeStack.push(next);
2003                }
2004        }
2005
2006        return NULL;
2007}
2008
2009BspLeaf *VspBspTree::GetRandomLeaf(const bool onlyUnmailed)
2010{
2011        stack<BspNode *> nodeStack;
2012
2013        nodeStack.push(mRoot);
2014
2015        int mask = rand();
2016
2017        while (!nodeStack.empty())
2018        {
2019                BspNode *node = nodeStack.top();
2020                nodeStack.pop();
2021
2022                if (node->IsLeaf())
2023                {
2024                        if ( (!onlyUnmailed || !node->Mailed()) )
2025                                return dynamic_cast<BspLeaf *>(node);
2026                }
2027                else
2028                {
2029                        BspInterior *interior = dynamic_cast<BspInterior *>(node);
2030
2031                        // random decision
2032                        if (mask & 1)
2033                                nodeStack.push(interior->GetBack());
2034                        else
2035                                nodeStack.push(interior->GetFront());
2036
2037                        mask = mask >> 1;
2038                }
2039        }
2040
2041        return NULL;
2042}
2043
2044int VspBspTree::ComputePvsSize(const RayInfoContainer &rays) const
2045{
2046        int pvsSize = 0;
2047
2048        RayInfoContainer::const_iterator rit, rit_end = rays.end();
2049
2050        Intersectable::NewMail();
2051
2052        for (rit = rays.begin(); rit != rays.end(); ++ rit)
2053        {
2054                VssRay *ray = (*rit).mRay;
2055
2056                if (ray->mOriginObject)
2057                {
2058                        if (!ray->mOriginObject->Mailed())
2059                        {
2060                                ray->mOriginObject->Mail();
2061                                ++ pvsSize;
2062                        }
2063                }
2064                if (ray->mTerminationObject)
2065                {
2066                        if (!ray->mTerminationObject->Mailed())
2067                        {
2068                                ray->mTerminationObject->Mail();
2069                                ++ pvsSize;
2070                        }
2071                }
2072        }
2073
2074        return pvsSize;
2075}
2076
2077float VspBspTree::GetEpsilon() const
2078{
2079        return mEpsilon;
2080}
2081
2082
2083int VspBspTree::SplitPolygons(const Plane3 &plane,
2084                                                          PolygonContainer &polys,
2085                                                          PolygonContainer &frontPolys,
2086                                                          PolygonContainer &backPolys,
2087                                                          PolygonContainer &coincident) const
2088{
2089        int splits = 0;
2090
2091        while (!polys.empty())
2092        {
2093                Polygon3 *poly = polys.back();
2094                polys.pop_back();
2095
2096                // classify polygon
2097                const int cf = poly->ClassifyPlane(plane, mEpsilon);
2098
2099                switch (cf)
2100                {
2101                        case Polygon3::COINCIDENT:
2102                                coincident.push_back(poly);
2103                                break;
2104                        case Polygon3::FRONT_SIDE:
2105                                frontPolys.push_back(poly);
2106                                break;
2107                        case Polygon3::BACK_SIDE:
2108                                backPolys.push_back(poly);
2109                                break;
2110                        case Polygon3::SPLIT:
2111                                backPolys.push_back(poly);
2112                                frontPolys.push_back(poly);
2113                                ++ splits;
2114                                break;
2115                        default:
2116                Debug << "SHOULD NEVER COME HERE\n";
2117                                break;
2118                }
2119        }
2120
2121        return splits;
2122}
2123
2124
2125int VspBspTree::CastLineSegment(const Vector3 &origin,
2126                                                                const Vector3 &termination,
2127                                                                vector<ViewCell *> &viewcells)
2128{
2129        int hits = 0;
2130        stack<BspRayTraversalData> tStack;
2131
2132        float mint = 0.0f, maxt = 1.0f;
2133
2134        Intersectable::NewMail();
2135
2136        Vector3 entp = origin;
2137        Vector3 extp = termination;
2138
2139        BspNode *node = mRoot;
2140        BspNode *farChild = NULL;
2141
2142        float t;
2143        while (1)
2144        {
2145                if (!node->IsLeaf())
2146                {
2147                        BspInterior *in = dynamic_cast<BspInterior *>(node);
2148
2149                        Plane3 splitPlane = in->GetPlane();
2150                       
2151                        const int entSide = splitPlane.Side(entp);
2152                        const int extSide = splitPlane.Side(extp);
2153
2154                        if (entSide < 0)
2155                        {
2156                                node = in->GetBack();
2157                                // plane does not split ray => no far child
2158                                if (extSide <= 0)
2159                                        continue;
2160
2161                                farChild = in->GetFront(); // plane splits ray
2162                        }
2163                        else if (entSide > 0)
2164                        {
2165                                node = in->GetFront();
2166
2167                                if (extSide >= 0) // plane does not split ray => no far child
2168                                        continue;
2169
2170                                farChild = in->GetBack(); // plane splits ray
2171                        }
2172                        else // ray end point on plane
2173                        {       // NOTE: what to do if ray is coincident with plane?
2174                                if (extSide < 0)
2175                                        node = in->GetBack();
2176                                else
2177                                        node = in->GetFront();
2178                                                               
2179                                continue; // no far child
2180                        }
2181
2182                        // push data for far child
2183                        tStack.push(BspRayTraversalData(farChild, extp));
2184
2185                        // find intersection of ray segment with plane
2186                        extp = splitPlane.FindIntersection(origin, extp, &t);
2187                }
2188                else
2189                {
2190                        // reached leaf => intersection with view cell
2191                        BspLeaf *leaf = dynamic_cast<BspLeaf *>(node);
2192
2193                        if (!leaf->GetViewCell()->Mailed())
2194                        {
2195                                viewcells.push_back(leaf->GetViewCell());
2196                                leaf->GetViewCell()->Mail();
2197                                ++ hits;
2198                        }
2199
2200                        //-- fetch the next far child from the stack
2201                        if (tStack.empty())
2202                                break;
2203
2204                        entp = extp;
2205                       
2206                        BspRayTraversalData &s = tStack.top();
2207
2208                        node = s.mNode;
2209                        extp = s.mExitPoint;
2210
2211                        tStack.pop();
2212                }
2213        }
2214
2215        return hits;
2216}
2217
2218int VspBspTree::TreeDistance(BspNode *n1, BspNode *n2) const
2219{
2220        std::deque<BspNode *> path1;
2221        BspNode *p1 = n1;
2222
2223        // create path from node 1 to root
2224        while (p1)
2225        {
2226                if (p1 == n2) // second node on path
2227                        return (int)path1.size();
2228
2229                path1.push_front(p1);
2230                p1 = p1->GetParent();
2231        }
2232
2233        int depth = n2->GetDepth();
2234        int d = depth;
2235
2236        BspNode *p2 = n2;
2237
2238        // compare with same depth
2239        while (1)
2240        {
2241                if ((d < (int)path1.size()) && (p2 == path1[d]))
2242                        return (depth - d) + ((int)path1.size() - 1 - d);
2243
2244                -- d;
2245                p2 = p2->GetParent();
2246        }
2247
2248        return 0; // never come here
2249}
2250
2251BspNode *VspBspTree::CollapseTree(BspNode *node, int &collapsed)
2252{
2253        if (node->IsLeaf())
2254                return node;
2255
2256        BspInterior *interior = dynamic_cast<BspInterior *>(node);
2257
2258        BspNode *front = CollapseTree(interior->GetFront(), collapsed);
2259        BspNode *back = CollapseTree(interior->GetBack(), collapsed);
2260
2261        if (front->IsLeaf() && back->IsLeaf())
2262        {
2263                BspLeaf *frontLeaf = dynamic_cast<BspLeaf *>(front);
2264                BspLeaf *backLeaf = dynamic_cast<BspLeaf *>(back);
2265
2266                //-- collapse tree
2267                if (frontLeaf->GetViewCell() == backLeaf->GetViewCell())
2268                {
2269                        BspViewCell *vc = frontLeaf->GetViewCell();
2270
2271                        BspLeaf *leaf = new BspLeaf(interior->GetParent(), vc);
2272                        leaf->SetTreeValid(frontLeaf->TreeValid());
2273
2274                        // replace a link from node's parent
2275                        if (leaf->GetParent())
2276                                leaf->GetParent()->ReplaceChildLink(node, leaf);
2277                        else
2278                                mRoot = leaf;
2279
2280                        ++ collapsed;
2281                        delete interior;
2282
2283                        return leaf;
2284                }
2285        }
2286
2287        return node;
2288}
2289
2290
2291int VspBspTree::CollapseTree()
2292{
2293        int collapsed = 0;
2294       
2295        (void)CollapseTree(mRoot, collapsed);
2296
2297        // revalidate leaves
2298        RepairViewCellsLeafLists();
2299
2300        return collapsed;
2301}
2302
2303
2304void VspBspTree::RepairViewCellsLeafLists()
2305{
2306        // list not valid anymore => clear
2307        stack<BspNode *> nodeStack;
2308        nodeStack.push(mRoot);
2309
2310        ViewCell::NewMail();
2311
2312        while (!nodeStack.empty())
2313        {
2314                BspNode *node = nodeStack.top();
2315                nodeStack.pop();
2316
2317                if (node->IsLeaf())
2318                {
2319                        BspLeaf *leaf = dynamic_cast<BspLeaf *>(node);
2320
2321                        BspViewCell *viewCell = leaf->GetViewCell();
2322
2323                        if (!viewCell->Mailed())
2324                        {
2325                                viewCell->mLeaves.clear();
2326                                viewCell->Mail();
2327                        }
2328
2329                        viewCell->mLeaves.push_back(leaf);
2330                }
2331                else
2332                {
2333                        BspInterior *interior = dynamic_cast<BspInterior *>(node);
2334
2335                        nodeStack.push(interior->GetFront());
2336                        nodeStack.push(interior->GetBack());
2337                }
2338        }
2339}
2340
2341
2342typedef pair<BspNode *, BspNodeGeometry *> bspNodePair;
2343
2344
2345int VspBspTree::CastBeam(Beam &beam)
2346{
2347    stack<bspNodePair> nodeStack;
2348        BspNodeGeometry *rgeom = new BspNodeGeometry();
2349        ConstructGeometry(mRoot, *rgeom);
2350
2351        nodeStack.push(bspNodePair(mRoot, rgeom));
2352 
2353        ViewCell::NewMail();
2354
2355        while (!nodeStack.empty())
2356        {
2357                BspNode *node = nodeStack.top().first;
2358                BspNodeGeometry *geom = nodeStack.top().second;
2359                nodeStack.pop();
2360               
2361                AxisAlignedBox3 box;
2362                box.Initialize();
2363                geom->IncludeInBox(box);
2364
2365                const int side = beam.ComputeIntersection(box);
2366               
2367                switch (side)
2368                {
2369                case -1:
2370                        CollectViewCells(node, beam.mViewCells, true);
2371                        break;
2372                case 0:
2373                       
2374                        if (node->IsLeaf())
2375                        {
2376                                BspLeaf *leaf = dynamic_cast<BspLeaf *>(node);
2377                       
2378                                if (!leaf->GetViewCell()->Mailed() && leaf->TreeValid())
2379                                {
2380                                        leaf->GetViewCell()->Mail();
2381                                        beam.mViewCells.push_back(leaf->GetViewCell());
2382                                }
2383                        }
2384                        else
2385                        {
2386                                BspInterior *interior = dynamic_cast<BspInterior *>(node);
2387                       
2388                                BspNode *first = interior->GetFront();
2389                                BspNode *second = interior->GetBack();
2390           
2391                                BspNodeGeometry *firstGeom = new BspNodeGeometry();
2392                                BspNodeGeometry *secondGeom = new BspNodeGeometry();
2393
2394                                geom->SplitGeometry(*firstGeom,
2395                                                                        *secondGeom,
2396                                                                        interior->GetPlane(),
2397                                                                        mBox,
2398                                                                        mEpsilon);
2399
2400                                // decide on the order of the nodes
2401                                if (DotProd(beam.mPlanes[0].mNormal,
2402                                        interior->GetPlane().mNormal) > 0)
2403                                {
2404                                        swap(first, second);
2405                                        swap(firstGeom, secondGeom);
2406                                }
2407
2408                                nodeStack.push(bspNodePair(first, firstGeom));
2409                                nodeStack.push(bspNodePair(second, secondGeom));
2410                        }
2411                       
2412                        break;
2413                default:
2414                        // default: cull
2415                        break;
2416                }
2417               
2418                DEL_PTR(geom);
2419               
2420        }
2421
2422        return (int)beam.mViewCells.size();
2423}
2424
2425
2426bool VspBspTree::MergeViewCells(BspLeaf *l1, BspLeaf *l2) const
2427{
2428        //-- change pointer to view cells of all leaves associated
2429        //-- with the previous view cells
2430        BspViewCell *fVc = l1->GetViewCell();
2431        BspViewCell *bVc = l2->GetViewCell();
2432
2433        BspViewCell *vc = dynamic_cast<BspViewCell *>(
2434                mViewCellsManager->MergeViewCells(*fVc, *bVc));
2435
2436        // if merge was unsuccessful
2437        if (!vc) return false;
2438
2439        // set new size of view cell
2440        vc->SetArea(fVc->GetArea() + bVc->GetArea());
2441
2442        vector<BspLeaf *> fLeaves = fVc->mLeaves;
2443        vector<BspLeaf *> bLeaves = bVc->mLeaves;
2444
2445        vector<BspLeaf *>::const_iterator it;
2446
2447        //-- change view cells of all the other leaves the view cell belongs to
2448        for (it = fLeaves.begin(); it != fLeaves.end(); ++ it)
2449        {
2450                (*it)->SetViewCell(vc);
2451                vc->mLeaves.push_back(*it);
2452        }
2453
2454        for (it = bLeaves.begin(); it != bLeaves.end(); ++ it)
2455        {
2456                (*it)->SetViewCell(vc);
2457                vc->mLeaves.push_back(*it);
2458        }
2459
2460        vc->Mail();
2461
2462        // clean up old view cells
2463        DEL_PTR(fVc);
2464        DEL_PTR(bVc);
2465
2466        return true;
2467}
2468
2469
2470void VspBspTree::SetViewCellsManager(ViewCellsManager *vcm)
2471{
2472        mViewCellsManager = vcm;
2473}
2474
2475
2476int VspBspTree::CollectMergeCandidates(const vector<BspLeaf *> leaves)
2477{
2478        BspLeaf::NewMail();
2479       
2480        vector<BspLeaf *>::const_iterator it, it_end = leaves.end();
2481
2482        int candidates = 0;
2483
2484        // find merge candidates and push them into queue
2485        for (it = leaves.begin(); it != it_end; ++ it)
2486        {
2487                BspLeaf *leaf = *it;
2488               
2489                /// create leaf pvs (needed for post processing
2490                leaf->mPvs = new ObjectPvs(leaf->GetViewCell()->GetPvs());
2491
2492                BspMergeCandidate::sOverallCost +=
2493                        leaf->mArea * leaf->mPvs->GetSize();
2494
2495                // the same leaves must not be part of two merge candidates
2496                leaf->Mail();
2497                vector<BspLeaf *> neighbors;
2498                FindNeighbors(leaf, neighbors, true);
2499
2500                vector<BspLeaf *>::const_iterator nit, nit_end = neighbors.end();
2501
2502                // TODO: test if at least one ray goes from one leaf to the other
2503                for (nit = neighbors.begin(); nit != nit_end; ++ nit)
2504                {
2505                        if ((*nit)->GetViewCell() != leaf->GetViewCell())
2506                        {
2507                                mMergeQueue.push(BspMergeCandidate(leaf, *nit));
2508                                ++ candidates;
2509                        }
2510                }
2511        }
2512
2513        Debug << "found " << candidates << " new merge candidates" << endl;
2514
2515        return (int)leaves.size();
2516}
2517
2518
2519int VspBspTree::CollectMergeCandidates(const VssRayContainer &rays)
2520{
2521        vector<BspRay *> bspRays;
2522
2523        long startTime = GetTime();
2524        ConstructBspRays(bspRays, rays);
2525        Debug << (int)bspRays.size() << " bsp rays constructed in "
2526                  << TimeDiff(startTime, GetTime()) * 1e-3f << " secs" << endl;
2527
2528        map<BspLeaf *, vector<BspLeaf*> > neighborMap;
2529        vector<BspIntersection>::const_iterator iit;
2530
2531        int numLeaves = 0;
2532       
2533        BspLeaf::NewMail();
2534
2535        for (int i = 0; i < (int)bspRays.size(); ++ i)
2536        { 
2537                BspRay *ray = bspRays[i];
2538         
2539                // traverse leaves stored in the rays and compare and
2540                // merge consecutive leaves (i.e., the neighbors in the tree)
2541                if (ray->intersections.size() < 2)
2542                        continue;
2543         
2544                iit = ray->intersections.begin();
2545                BspLeaf *leaf = (*(iit ++)).mLeaf;
2546               
2547                // create leaf pvs (needed for post processing)
2548                if (!leaf->mPvs)
2549                {
2550                        leaf->mPvs =
2551                                new ObjectPvs(leaf->GetViewCell()->GetPvs());
2552
2553                        BspMergeCandidate::sOverallCost +=
2554                                leaf->mArea * leaf->mPvs->GetSize();
2555                       
2556                        ++ numLeaves;
2557                }
2558               
2559                // traverse intersections
2560                // consecutive leaves are neighbors => add them to queue
2561                for (; iit != ray->intersections.end(); ++ iit)
2562                {
2563                        // next pair
2564                        BspLeaf *prevLeaf = leaf;
2565            leaf = (*iit).mLeaf;
2566
2567                        // view space not valid or same view cell
2568                        if (!leaf->TreeValid() || !prevLeaf->TreeValid() ||
2569                                (leaf->GetViewCell() == prevLeaf->GetViewCell()))
2570                                continue;
2571
2572            // create leaf pvs (needed for post processing)
2573                        if (!leaf->mPvs)
2574                        {
2575                                leaf->mPvs =
2576                                        new ObjectPvs(leaf->GetViewCell()->GetPvs());
2577                               
2578                                BspMergeCandidate::sOverallCost +=
2579                                        leaf->mArea * leaf->mPvs->GetSize();
2580
2581                                ++ numLeaves;
2582                        }
2583               
2584                        vector<BspLeaf *> &neighbors = neighborMap[leaf];
2585                       
2586                        bool found = false;
2587
2588                        // both leaves inserted in queue already =>
2589                        // look if double pair already exists
2590                        if (leaf->Mailed() && prevLeaf->Mailed())
2591                        {
2592                                vector<BspLeaf *>::const_iterator it, it_end = neighbors.end();
2593                               
2594                for (it = neighbors.begin(); !found && (it != it_end); ++ it)
2595                                        if (*it == prevLeaf)
2596                                                found = true; // already in queue
2597                        }
2598                       
2599                        if (!found)
2600                        {
2601                                // this pair is not in map already
2602                                // => insert into the neighbor map and the queue
2603                                neighbors.push_back(prevLeaf);
2604                                neighborMap[prevLeaf].push_back(leaf);
2605
2606                                leaf->Mail();
2607                                prevLeaf->Mail();
2608
2609                                mMergeQueue.push(BspMergeCandidate(leaf, prevLeaf));
2610                        }
2611        }
2612        }
2613        Debug << "neighbormap size: " << (int)neighborMap.size() << endl;
2614        Debug << "mergequeue: " << (int)mMergeQueue.size() << endl;
2615        Debug << "leaves in queue: " << numLeaves << endl;
2616        Debug << "overall cost: " << BspMergeCandidate::sOverallCost << endl;
2617
2618        CLEAR_CONTAINER(bspRays);
2619
2620        //-- collect the leaves which haven't been found by ray casting
2621        if (0)
2622        {
2623                vector<BspLeaf *> leaves;
2624                CollectLeaves(leaves, true, mMaxPvs);
2625                Debug << "found " << (int)leaves.size() << " new leaves" << endl << endl;
2626                CollectMergeCandidates(leaves);
2627        }
2628
2629        return numLeaves;
2630}
2631
2632
2633void VspBspTree::ConstructBspRays(vector<BspRay *> &bspRays,
2634                                                                  const VssRayContainer &rays)
2635{
2636        VssRayContainer::const_iterator it, it_end = rays.end();
2637
2638        for (it = rays.begin(); it != rays.end(); ++ it)
2639        {
2640                VssRay *vssRay = *it;
2641                BspRay *ray = new BspRay(vssRay);
2642
2643                ViewCellContainer viewCells;
2644
2645                Ray hray(*vssRay);
2646                float tmin = 0, tmax = 1.0;
2647                // matt TODO: remove this!!
2648                //hray.Init(ray.GetOrigin(), ray.GetDir(), Ray::LINE_SEGMENT);
2649                if (!mBox.GetRaySegment(hray, tmin, tmax) || (tmin > tmax))
2650                        continue;
2651
2652                Vector3 origin = hray.Extrap(tmin);
2653                Vector3 termination = hray.Extrap(tmax);
2654       
2655                // cast line segment to get intersections with bsp leaves
2656                CastLineSegment(origin, termination, viewCells);
2657
2658                ViewCellContainer::const_iterator vit, vit_end = viewCells.end();
2659                for (vit = viewCells.begin(); vit != vit_end; ++ vit)
2660                {
2661                        BspViewCell *vc = dynamic_cast<BspViewCell *>(*vit);
2662                        vector<BspLeaf *>::const_iterator it, it_end = vc->mLeaves.end();
2663                        //NOTE: not sorted!
2664                        for (it = vc->mLeaves.begin(); it != it_end; ++ it)
2665                                ray->intersections.push_back(BspIntersection(0, *it));
2666                }
2667
2668                bspRays.push_back(ray);
2669        }
2670}
2671
2672
2673int VspBspTree::MergeViewCells(const VssRayContainer &rays)
2674{
2675        BspMergeCandidate::sMaxPvsSize = mMaxPvs;
2676
2677        MergeStatistics mergeStats;
2678        mergeStats.Start();
2679        // TODO: REMOVE LATER for performance!
2680        const bool showMergeStats = true;
2681        //BspMergeCandidate::sOverallCost = mBox.SurfaceArea() * mStat.maxPvs;
2682        long startTime = GetTime();
2683
2684        if (mUseRaysForMerge)
2685                mergeStats.nodes = CollectMergeCandidates(rays);
2686        else
2687        {
2688                vector<BspLeaf *> leaves;
2689                CollectLeaves(leaves);
2690                mergeStats.nodes = CollectMergeCandidates(leaves);
2691        }
2692
2693        mergeStats.collectTime = TimeDiff(startTime, GetTime());
2694        mergeStats.candidates = (int)mMergeQueue.size();
2695        startTime = GetTime();
2696
2697        int nViewCells = /*mergeStats.nodes*/ mStat.Leaves() - mStat.invalidLeaves;
2698
2699        //-- use priority queue to merge leaf pairs
2700        while (!mMergeQueue.empty() && (nViewCells > mMergeMinViewCells) &&
2701                   (mMergeQueue.top().GetMergeCost() <
2702                    mMergeMaxCostRatio * BspMergeCandidate::sOverallCost))
2703        {
2704                /*Debug << "abs mergecost: " << mMergeQueue.top().GetMergeCost() << " rel mergecost: "
2705                          << mMergeQueue.top().GetMergeCost() / BspMergeCandidate::sOverallCost
2706                          << " max ratio: " << mMergeMaxCostRatio << endl;*/
2707                BspMergeCandidate mc = mMergeQueue.top();
2708                mMergeQueue.pop();
2709
2710                // both view cells equal!
2711                if (mc.GetLeaf1()->GetViewCell() == mc.GetLeaf2()->GetViewCell())
2712                        continue;
2713
2714                if (mc.Valid())
2715                {
2716                        ViewCell::NewMail();
2717       
2718                        MergeViewCells(mc.GetLeaf1(), mc.GetLeaf2());
2719                        -- nViewCells;
2720                       
2721                        ++ mergeStats.merged;
2722                       
2723                        // increase absolute merge cost
2724                        BspMergeCandidate::sOverallCost += mc.GetMergeCost();
2725
2726                        if (showMergeStats)
2727                        {
2728                                if (mc.GetLeaf1()->IsSibling(mc.GetLeaf2()))
2729                                        ++ mergeStats.siblings;
2730
2731                                const int dist =
2732                                        TreeDistance(mc.GetLeaf1(), mc.GetLeaf2());
2733                                if (dist > mergeStats.maxTreeDist)
2734                                        mergeStats.maxTreeDist = dist;
2735                                mergeStats.accTreeDist += dist;
2736                        }
2737                }
2738                // merge candidate not valid, because one of the leaves was already
2739                // merged with another one => validate and reinsert into queue
2740                else
2741                {
2742                        mc.SetValid();
2743                        mMergeQueue.push(mc);
2744                }
2745        }
2746
2747        mergeStats.mergeTime = TimeDiff(startTime, GetTime());
2748        mergeStats.Stop();
2749
2750        if (showMergeStats)
2751                Debug << mergeStats << endl << endl;
2752       
2753        //TODO: should return sample contributions?
2754        return mergeStats.merged;
2755}
2756
2757
2758ViewCell *VspBspTree::GetViewCell(const Vector3 &point)
2759{
2760  if (mRoot == NULL)
2761        return NULL;
2762 
2763  stack<BspNode *> nodeStack;
2764  nodeStack.push(mRoot);
2765 
2766  ViewCell *viewcell = NULL;
2767 
2768  while (!nodeStack.empty())  {
2769        BspNode *node = nodeStack.top();
2770        nodeStack.pop();
2771       
2772        if (node->IsLeaf()) {
2773          viewcell = dynamic_cast<BspLeaf *>(node)->GetViewCell();
2774          break;
2775        } else {
2776         
2777          BspInterior *interior = dynamic_cast<BspInterior *>(node);
2778               
2779          // random decision
2780          if (interior->GetPlane().Side(point) < 0)
2781                nodeStack.push(interior->GetBack());
2782          else
2783                nodeStack.push(interior->GetFront());
2784        }
2785  }
2786 
2787  return viewcell;
2788}
2789
2790
2791int VspBspTree::RefineViewCells(const VssRayContainer &rays)
2792{
2793        int shuffled = 0;
2794
2795        Debug << "refining " << (int)mMergeQueue.size() << " candidates " << endl;
2796        BspLeaf::NewMail();
2797
2798        // Use priority queue of remaining leaf pairs
2799        // These candidates either share the same view cells or
2800        // are border leaves which share a boundary.
2801        // We test if they can be shuffled, i.e.,
2802        // either one leaf is made part of one view cell or the other
2803        // leaf is made part of the other view cell. It is tested if the
2804        // remaining view cells are "better" than the old ones.
2805        while (!mMergeQueue.empty())
2806        {
2807                BspMergeCandidate mc = mMergeQueue.top();
2808                mMergeQueue.pop();
2809
2810                // both view cells equal or already shuffled
2811                if ((mc.GetLeaf1()->GetViewCell() == mc.GetLeaf2()->GetViewCell()) ||
2812                        (mc.GetLeaf1()->Mailed()) || (mc.GetLeaf2()->Mailed()))
2813                        continue;
2814               
2815                // candidate for shuffling
2816                const bool wasShuffled =
2817                        ShuffleLeaves(mc.GetLeaf1(), mc.GetLeaf2());
2818               
2819                //-- stats
2820                if (wasShuffled)
2821                        ++ shuffled;
2822        }
2823
2824        return shuffled;
2825}
2826
2827
2828inline int AddedPvsSize(ObjectPvs pvs1, const ObjectPvs &pvs2)
2829{
2830        return pvs1.AddPvs(pvs2);
2831}
2832
2833/*inline int SubtractedPvsSize(BspViewCell *vc, BspLeaf *l, const ObjectPvs &pvs2)
2834{
2835        ObjectPvs pvs;
2836        vector<BspLeaf *>::const_iterator it, it_end = vc->mLeaves.end();
2837        for (it = vc->mLeaves.begin(); it != vc->mLeaves.end(); ++ it)
2838                if (*it != l)
2839                        pvs.AddPvs(*(*it)->mPvs);
2840        return pvs.GetSize();
2841}*/
2842
2843inline int SubtractedPvsSize(ObjectPvs pvs1, const ObjectPvs &pvs2)
2844{
2845        return pvs1.SubtractPvs(pvs2);
2846}
2847
2848
2849float GetShuffledVcCost(BspLeaf *leaf, BspViewCell *vc1, BspViewCell *vc2)
2850{
2851        //const int pvs1 = SubtractedPvsSize(vc1, leaf, *leaf->mPvs);
2852        const int pvs1 = SubtractedPvsSize(vc1->GetPvs(), *leaf->mPvs);
2853        const int pvs2 = AddedPvsSize(vc2->GetPvs(), *leaf->mPvs);
2854
2855        const float area1 = vc1->GetArea() - leaf->mArea;
2856        const float area2 = vc2->GetArea() + leaf->mArea;
2857
2858        const float cost1 = pvs1 * area1;
2859        const float cost2 = pvs2 * area2;
2860
2861        return cost1 + cost2;
2862}
2863
2864
2865void VspBspTree::ShuffleLeaf(BspLeaf *leaf,
2866                                                         BspViewCell *vc1,
2867                                                         BspViewCell *vc2) const
2868{
2869        // compute new pvs and area
2870        vc1->GetPvs().SubtractPvs(*leaf->mPvs);
2871        vc2->GetPvs().AddPvs(*leaf->mPvs);
2872       
2873        vc1->SetArea(vc1->GetArea() - leaf->mArea);
2874        vc2->SetArea(vc2->GetArea() + leaf->mArea);
2875
2876        /// add to second view cell
2877        vc2->mLeaves.push_back(leaf);
2878
2879        // erase leaf from old view cell
2880        vector<BspLeaf *>::iterator it = vc1->mLeaves.begin();
2881
2882        for (; *it != leaf; ++ it);
2883        vc1->mLeaves.erase(it);
2884
2885        /*vc1->GetPvs().mEntries.clear();
2886        for (; it != vc1->mLeaves.end(); ++ it)
2887        {
2888                if (*it == leaf)
2889                        vc1->mLeaves.erase(it);
2890                else
2891                        vc1->GetPvs().AddPvs(*(*it)->mPvs);
2892        }*/
2893
2894        leaf->SetViewCell(vc2); // finally change view cell
2895}
2896
2897
2898bool VspBspTree::ShuffleLeaves(BspLeaf *leaf1, BspLeaf *leaf2) const
2899{
2900        BspViewCell *vc1 = leaf1->GetViewCell();
2901        BspViewCell *vc2 = leaf2->GetViewCell();
2902
2903        const float cost1 = vc1->GetPvs().GetSize() * vc1->GetArea();
2904        const float cost2 = vc2->GetPvs().GetSize() * vc2->GetArea();
2905
2906        const float oldCost = cost1 + cost2;
2907       
2908        float shuffledCost1 = Limits::Infinity;
2909        float shuffledCost2 = Limits::Infinity;
2910
2911        // the view cell should not be empty after the shuffle
2912        if (vc1->mLeaves.size() > 1)
2913                shuffledCost1 = GetShuffledVcCost(leaf1, vc1, vc2);
2914        if (vc2->mLeaves.size() > 1)
2915                shuffledCost2 = GetShuffledVcCost(leaf2, vc2, vc1);
2916
2917        // shuffling unsuccessful
2918        if ((oldCost <= shuffledCost1) && (oldCost <= shuffledCost2))
2919                return false;
2920       
2921        if (shuffledCost1 < shuffledCost2)
2922        {
2923                ShuffleLeaf(leaf1, vc1, vc2);
2924                leaf1->Mail();
2925        }
2926        else
2927        {
2928                ShuffleLeaf(leaf2, vc2, vc1);
2929                leaf2->Mail();
2930        }
2931
2932        return true;
2933}
2934
2935bool VspBspTree::ViewPointValid(const Vector3 &viewPoint) const
2936{
2937        BspNode *node = mRoot;
2938
2939        while (1)
2940        {
2941                // early exit
2942                if (node->TreeValid())
2943                        return true;
2944
2945                if (node->IsLeaf())
2946                        return false;
2947                       
2948                BspInterior *in = dynamic_cast<BspInterior *>(node);
2949                                       
2950                if (in->GetPlane().Side(viewPoint) <= 0)
2951                {
2952                        node = in->GetBack();
2953                }
2954                else
2955                {
2956                        node = in->GetFront();
2957                }
2958        }
2959
2960        // should never come here
2961        return false;
2962}
2963
2964
2965bool VspBspTree::CheckValid(const VspBspTraversalData &data) const
2966{
2967        return data.mPvs <= mMaxPvs;
2968}
2969
2970
2971void VspBspTree::PropagateUpValidity(BspNode *node)
2972{
2973        while (!node->IsRoot() && node->GetParent()->TreeValid())
2974        {
2975                node = node->GetParent();
2976                node->SetTreeValid(false);
2977        }
2978}
2979
2980
2981bool VspBspTree::Export(ofstream &stream)
2982{
2983        ExportNode(mRoot, stream);
2984
2985        return true;
2986}
2987
2988
2989void VspBspTree::ExportNode(BspNode *node, ofstream &stream)
2990{
2991        if (node->IsLeaf())
2992        {
2993                BspLeaf *leaf = dynamic_cast<BspLeaf *>(node);
2994                       
2995                int id = -1;
2996                if (leaf->GetViewCell() != mOutOfBoundsCell)
2997                        id = leaf->GetViewCell()->GetId();
2998
2999                stream << "<Leaf viewCellId=\"" << id << "\" />" << endl;
3000        }
3001        else
3002        {
3003                BspInterior *interior = dynamic_cast<BspInterior *>(node);
3004       
3005                Plane3 plane = interior->GetPlane();
3006                stream << "<Interior plane=\"" << plane.mNormal.x << " "
3007                           << plane.mNormal.y << " " << plane.mNormal.z << " "
3008                           << plane.mD << "\">" << endl;
3009
3010                ExportNode(interior->GetBack(), stream);
3011                ExportNode(interior->GetFront(), stream);
3012
3013                stream << "</Interior>" << endl;
3014        }
3015}
3016
3017/************************************************************************/
3018/*                BspMergeCandidate implementation                      */
3019/************************************************************************/
3020
3021
3022BspMergeCandidate::BspMergeCandidate(BspLeaf *l1, BspLeaf *l2):
3023mMergeCost(0),
3024mLeaf1(l1),
3025mLeaf2(l2),
3026mLeaf1Id(l1->GetViewCell()->mMailbox),
3027mLeaf2Id(l2->GetViewCell()->mMailbox)
3028{
3029        EvalMergeCost();
3030}
3031
3032float BspMergeCandidate::GetCost(ViewCell *vc) const
3033{
3034        return vc->GetPvs().GetSize() * vc->GetArea();
3035}
3036
3037float BspMergeCandidate::GetLeaf1Cost() const
3038{
3039        BspViewCell *vc = mLeaf1->GetViewCell();
3040        return GetCost(vc);
3041}
3042
3043float BspMergeCandidate::GetLeaf2Cost() const
3044{
3045        BspViewCell *vc = mLeaf2->GetViewCell();
3046        return GetCost(vc);
3047}
3048
3049void BspMergeCandidate::EvalMergeCost()
3050{
3051        //-- compute pvs difference
3052        BspViewCell *vc1 = mLeaf1->GetViewCell();
3053        BspViewCell *vc2 = mLeaf2->GetViewCell();
3054
3055        const int diff1 = vc1->GetPvs().Diff(vc2->GetPvs());
3056        const int newPvs = diff1 + vc1->GetPvs().GetSize();
3057
3058        //-- compute ratio of old cost
3059        //   (i.e., added size of left and right view cell times pvs size)
3060        //   to new rendering cost (i.e, size of merged view cell times pvs size)
3061        const float oldCost = GetLeaf1Cost() + GetLeaf2Cost();
3062
3063        const float newCost =
3064                (float)newPvs * (vc1->GetArea() + vc2->GetArea());
3065
3066        if (newPvs > sMaxPvsSize) // strong penalty if pvs size too large
3067        {
3068                mMergeCost = 1e15f;
3069        }
3070        else
3071        {
3072                mMergeCost = newCost - oldCost;
3073        }
3074}
3075
3076
3077void BspMergeCandidate::SetLeaf1(BspLeaf *l)
3078{
3079        mLeaf1 = l;
3080}
3081
3082
3083void BspMergeCandidate::SetLeaf2(BspLeaf *l)
3084{
3085        mLeaf2 = l;
3086}
3087
3088
3089BspLeaf *BspMergeCandidate::GetLeaf1()
3090{
3091        return mLeaf1;
3092}
3093
3094
3095BspLeaf *BspMergeCandidate::GetLeaf2()
3096{
3097        return mLeaf2;
3098}
3099
3100
3101bool BspMergeCandidate::Valid() const
3102{
3103        return
3104                (mLeaf1->GetViewCell()->mMailbox == mLeaf1Id) &&
3105                (mLeaf2->GetViewCell()->mMailbox == mLeaf2Id);
3106}
3107
3108
3109float BspMergeCandidate::GetMergeCost() const
3110{
3111        return mMergeCost;
3112}
3113
3114
3115void BspMergeCandidate::SetValid()
3116{
3117        mLeaf1Id = mLeaf1->GetViewCell()->mMailbox;
3118        mLeaf2Id = mLeaf2->GetViewCell()->mMailbox;
3119
3120        EvalMergeCost();
3121}
3122
3123
3124/************************************************************************/
3125/*                  MergeStatistics implementation                      */
3126/************************************************************************/
3127
3128
3129void MergeStatistics::Print(ostream &app) const
3130{
3131        app << "===== Merge statistics ===============\n";
3132
3133        app << setprecision(4);
3134
3135        app << "#N_CTIME ( Overall time [s] )\n" << Time() << " \n";
3136
3137        app << "#N_CCTIME ( Collect candidates time [s] )\n" << collectTime * 1e-3f << " \n";
3138
3139        app << "#N_MTIME ( Merge time [s] )\n" << mergeTime * 1e-3f << " \n";
3140
3141        app << "#N_NODES ( Number of nodes before merge )\n" << nodes << "\n";
3142
3143        app << "#N_CANDIDATES ( Number of merge candidates )\n" << candidates << "\n";
3144
3145        app << "#N_MERGEDSIBLINGS ( Number of merged siblings )\n" << siblings << "\n";
3146        app << "#N_MERGEDNODES ( Number of merged nodes )\n" << merged << "\n";
3147
3148        app << "#MAX_TREEDIST ( Maximal distance in tree of merged leaves )\n" << maxTreeDist << "\n";
3149
3150        app << "#AVG_TREEDIST ( Average distance in tree of merged leaves )\n" << AvgTreeDist() << "\n";
3151
3152        app << "===== END OF BspTree statistics ==========\n";
3153}
Note: See TracBrowser for help on using the repository browser.