source: GTP/trunk/Lib/Vis/Preprocessing/src/FromPointVisibilityTree.cpp @ 1002

Revision 1002, 97.0 KB checked in by mattausch, 18 years ago (diff)

debug run: fixing memory holes

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