source: trunk/VUT/GtpVisibilityPreprocessor/src/ViewCellBsp.cpp @ 396

Revision 396, 59.4 KB checked in by mattausch, 19 years ago (diff)

fixed bug in bsp geometry construction

Line 
1#include "Plane3.h"
2#include "ViewCellBsp.h"
3#include "Mesh.h"
4#include "common.h"
5#include "ViewCell.h"
6#include "Environment.h"
7#include "Polygon3.h"
8#include "Ray.h"
9#include "AxisAlignedBox3.h"
10#include <stack>
11#include <time.h>
12#include <iomanip>
13#include "Exporter.h"
14#include "Plane3.h"
15
16int BspTree::sTermMaxPolygons = 10;
17int BspTree::sTermMinPvs = 20;
18int BspTree::sTermMaxDepth = 20;
19int BspTree::sMaxPolyCandidates = 10;
20int BspTree::sMaxRayCandidates = 10;
21int BspTree::sSplitPlaneStrategy = BALANCED_POLYS;
22int BspTree::sConstructionMethod = FROM_INPUT_VIEW_CELLS;
23int BspTree::sTermMaxPolysForAxisAligned = 50;
24int BspTree::sTermMaxObjectsForAxisAligned = 50;
25int BspTree::sTermMaxRaysForAxisAligned = -1;
26int BspTree::sTermMaxRays = -1;
27int BspTree::sMinPvsDif = 10;
28int BspTree::sMinPvs = 10;
29int BspTree::sMaxPvs = 500;
30
31float BspTree::sCt_div_ci = 1.0f;
32float BspTree::sSplitBorder = 1.0f;
33float BspTree::sMaxCostRatio = 0.9f;
34
35//-- factors for bsp tree split plane heuristics
36
37float BspTree::sLeastSplitsFactor = 1.0f;
38float BspTree::sBalancedPolysFactor = 1.0f;
39float BspTree::sBalancedViewCellsFactor = 1.0f;
40
41// NOTE:  very important criterium for 2.5d scenes
42float BspTree::sVerticalSplitsFactor = 1.0f;
43float BspTree::sLargestPolyAreaFactor = 1.0f;
44float BspTree::sBlockedRaysFactor = 1.0f;
45float BspTree::sLeastRaySplitsFactor = 1.0f;
46float BspTree::sBalancedRaysFactor = 1.0f;
47float BspTree::sPvsFactor = 1.0f;
48
49int BspNode::mailID = 1;
50
51/** Evaluates split plane classification with respect to the plane's
52        contribution for a minimum number splits in the tree.
53*/
54float BspTree::sLeastPolySplitsTable[] = {0, 0, 1, 0};
55/** Evaluates split plane classification with respect to the plane's
56        contribution for a balanced tree.
57*/
58float BspTree::sBalancedPolysTable[] = {1, -1, 0, 0};
59
60/** Evaluates split plane classification with respect to the plane's
61        contribution for a minimum number of ray splits.
62*/
63float BspTree::sLeastRaySplitsTable[] = {0, 0, 1, 1, 0};
64/** Evaluates split plane classification with respect to the plane's
65        contribution for balanced rays.
66*/
67float BspTree::sBalancedRaysTable[] = {1, -1, 0, 0, 0};
68
69/****************************************************************/
70/*                  class BspNode implementation                */
71/****************************************************************/
72
73BspNode::BspNode():
74mParent(NULL)
75{}
76
77BspNode::BspNode(BspInterior *parent):
78mParent(parent)
79{}
80
81
82bool BspNode::IsRoot() const
83{
84        return mParent == NULL;
85}
86
87BspInterior *BspNode::GetParent()
88{
89        return mParent;
90}
91
92void BspNode::SetParent(BspInterior *parent)
93{
94        mParent = parent;
95}
96
97/****************************************************************/
98/*              class BspInterior implementation                */
99/****************************************************************/
100
101
102BspInterior::BspInterior(const Plane3 &plane):
103mPlane(plane), mFront(NULL), mBack(NULL)
104{}
105
106bool BspInterior::IsLeaf() const
107{
108        return false;
109}
110
111BspNode *BspInterior::GetBack()
112{
113        return mBack;
114}
115
116BspNode *BspInterior::GetFront()
117{
118        return mFront;
119}
120
121Plane3 *BspInterior::GetPlane()
122{
123        return &mPlane;
124}
125
126void BspInterior::ReplaceChildLink(BspNode *oldChild, BspNode *newChild)
127{
128        if (mBack == oldChild)
129                mBack = newChild;
130        else
131                mFront = newChild;
132}
133
134void BspInterior::SetupChildLinks(BspNode *b, BspNode *f)
135{
136    mBack = b;
137    mFront = f;
138}
139
140int BspInterior::SplitPolygons(PolygonContainer &polys,
141                                                           PolygonContainer &frontPolys,
142                                                           PolygonContainer &backPolys,
143                                                           PolygonContainer &coincident)
144{
145        Polygon3 *splitPoly = NULL;
146
147        int splits = 0;
148
149#ifdef _Debug
150        Debug << "splitting polygons of node " << this << " with plane " << mPlane << endl;
151#endif
152        while (!polys.empty())
153        {
154                Polygon3 *poly = polys.back();
155                polys.pop_back();
156
157                //Debug << "New polygon with plane: " << poly->GetSupportingPlane() << "\n";
158
159                // classify polygon
160                const int cf = poly->ClassifyPlane(mPlane);
161
162                Polygon3 *front_piece = NULL;
163                Polygon3 *back_piece = NULL;
164       
165                VertexContainer splitVertices;
166
167                switch (cf)
168                {
169                        case Polygon3::COINCIDENT:
170                                coincident.push_back(poly);
171                                break;                 
172                        case Polygon3::FRONT_SIDE:     
173                                frontPolys.push_back(poly);
174                                break;
175                        case Polygon3::BACK_SIDE:
176                                backPolys.push_back(poly);
177                                break;
178                        case Polygon3::SPLIT:
179                                front_piece = new Polygon3(poly->mParent);
180                                back_piece = new Polygon3(poly->mParent);
181
182                                //-- split polygon into front and back part
183                                poly->Split(mPlane,
184                                                        *front_piece,
185                                                        *back_piece,
186                                                        splitVertices);
187                                       
188                                ++ splits; // increase number of splits
189
190                                //-- inherit rays from parent polygon for blocked ray criterium
191                                poly->InheritRays(*front_piece, *back_piece);
192                                //Debug << "p: " << poly->mPiercingRays.size() << " f: " << front_piece->mPiercingRays.size() << " b: " << back_piece->mPiercingRays.size() << endl;
193                               
194                                // check if polygons still valid
195                                if (front_piece->Valid())
196                                        frontPolys.push_back(front_piece);
197                                else
198                                        DEL_PTR(front_piece);
199                               
200                                if (back_piece->Valid())
201                                        backPolys.push_back(back_piece);
202                                else                           
203                                        DEL_PTR(back_piece);
204                               
205#ifdef _DEBUG
206                                Debug << "split " << *poly << endl << *front_piece << endl << *back_piece << endl;
207#endif
208                                DEL_PTR(poly);                 
209                                break;
210                        default:
211                Debug << "SHOULD NEVER COME HERE\n";
212                                break;
213                }
214        }
215
216        return splits;
217}
218
219/****************************************************************/
220/*                  class BspLeaf implementation                */
221/****************************************************************/
222BspLeaf::BspLeaf(): mViewCell(NULL)
223{
224}
225
226BspLeaf::BspLeaf(BspViewCell *viewCell):
227mViewCell(viewCell)
228{
229}
230
231BspLeaf::BspLeaf(BspInterior *parent):
232BspNode(parent), mViewCell(NULL)
233{}
234
235BspLeaf::BspLeaf(BspInterior *parent, BspViewCell *viewCell):
236BspNode(parent), mViewCell(viewCell)
237{
238}
239
240BspViewCell *BspLeaf::GetViewCell() const
241{
242        return mViewCell;
243}
244
245void BspLeaf::SetViewCell(BspViewCell *viewCell)
246{
247        mViewCell = viewCell;
248}
249
250bool BspLeaf::IsLeaf() const
251{
252        return true;
253}
254
255void BspLeaf::AddToPvs(const BoundedRayContainer &rays,
256                                           int &sampleContributions,
257                                           int &contributingSamples)
258{
259        sampleContributions = 0;
260        contributingSamples = 0;
261
262    BoundedRayContainer::const_iterator it, it_end = rays.end();
263
264        // add contributions from samples to the PVS
265        for (it = rays.begin(); it != it_end; ++ it)
266        {
267                int contribution = 0;
268                Ray *ray = (*it)->mRay;
269                       
270                if (!ray->intersections.empty())
271                {       
272            //contribution += mViewCell->GetPvs().AddSample(ray->intersections[0].mObject);
273                }
274
275                //contribution += mViewCell->GetPvs().AddSample(ray->sourceObject.mObject);
276
277                if (contribution > 0)
278                {
279                        sampleContributions += contribution;
280                        ++ contributingSamples;
281                }
282
283                // warning: not ordered
284                ray->bspIntersections.push_back(Ray::BspIntersection((*it)->mMinT, this));
285        }
286}
287
288
289/****************************************************************/
290/*                  class BspTree implementation                */
291/****************************************************************/
292
293BspTree::BspTree(BspViewCell *viewCell):
294mRootCell(viewCell), mRoot(NULL), mGenerateViewCells(false),
295mStorePiercingRays(true)
296{
297        Randomize(); // initialise random generator for heuristics
298}
299 
300const BspTreeStatistics &BspTree::GetStatistics() const
301{
302        return mStat;
303}
304
305void BspViewCellsStatistics::Print(ostream &app) const
306{
307        app << "===== BspViewCells statistics ===============\n";
308
309        app << setprecision(4);
310
311        app << "#N_OVERALLPVS ( objects in PVS )\n" << pvs << endl;
312
313        app << "#N_PMAXPVS ( largest PVS )\n" << maxPvs << endl;
314
315        app << "#N_PMINPVS ( smallest PVS )\n" << minPvs << endl;
316
317        app << "#N_PAVGPVS ( average PVS )\n" << AvgPvs() << endl;
318
319        app << "#N_PEMPTYPVS ( view cells with PVS smaller 2 )\n" << emptyPvs << endl;
320
321        app << "#N_VIEWCELLS ( number of view cells)\n" << viewCells << endl;
322
323        app << "#N_AVGBSPLEAVES (average number of BSP leaves per view cell )\n" << AvgBspLeaves() << endl;
324
325        app << "#N_MAXBSPLEAVES ( maximal number of BSP leaves per view cell )\n" << maxBspLeaves << endl;
326       
327        app << "===== END OF BspViewCells statistics ==========\n";
328}
329
330void BspTreeStatistics::Print(ostream &app) const
331{
332        app << "===== BspTree statistics ===============\n";
333
334        app << setprecision(4);
335
336        app << "#N_NODES ( Number of nodes )\n" << nodes << "\n";
337
338        app << "#N_INTERIORS ( Number of interior nodes )\n" << Interior() << "\n";
339
340        app << "#N_LEAVES ( Number of leaves )\n" << Leaves() << "\n";
341
342        app << "#N_SPLITS ( Number of splits )\n" << splits << "\n";
343
344        app << "#N_PMAXDEPTHLEAVES ( Percentage of leaves at maxdepth )\n"<<
345        maxDepthNodes * 100 / (double)Leaves() << endl;
346
347        app << "#N_PMAXDEPTH ( Maximal reached depth )\n" << maxDepth << endl;
348
349        app << "#N_PMINDEPTH ( Minimal reached depth )\n" << minDepth << endl;
350
351        app << "#AVGDEPTH ( average depth )\n" << AvgDepth() << endl;
352
353        app << "#N_INPUT_POLYGONS (number of input polygons )\n" <<     polys << endl;
354
355        //app << "#N_PVS: " << pvs << endl;
356
357        app << "#N_ROUTPUT_INPUT_POLYGONS ( ratio polygons after subdivision / input polygons )\n" <<
358                 (polys + splits) / (double)polys << endl;
359       
360        app << "===== END OF BspTree statistics ==========\n";
361}
362
363
364BspTree::~BspTree()
365{
366        std::stack<BspNode *> tStack;
367
368        tStack.push(mRoot);
369
370        while (!tStack.empty())
371        {
372                BspNode *node = tStack.top();
373
374            tStack.pop();
375       
376                if (!node->IsLeaf())
377                {
378                        BspInterior *interior = dynamic_cast<BspInterior *>(node);
379
380                        // push the children on the stack (there are always two children)
381                        tStack.push(interior->GetBack());
382                        tStack.push(interior->GetFront());
383                }
384
385                DEL_PTR(node);
386        }
387}
388
389void BspTree::InsertViewCell(ViewCell *viewCell)
390{
391        PolygonContainer *polys = new PolygonContainer();
392
393        // extract polygons that guide the split process
394        mStat.polys += AddMeshToPolygons(viewCell->GetMesh(), *polys, viewCell);
395        mBox.Include(viewCell->GetBox()); // add to BSP aabb
396
397        InsertPolygons(polys);
398}
399
400void BspTree::InsertPolygons(PolygonContainer *polys)
401{       
402        std::stack<BspTraversalData> tStack;
403       
404        // traverse existing tree or create new tree
405    if (!mRoot)
406                mRoot = new BspLeaf();
407
408        tStack.push(BspTraversalData(mRoot, polys, 0, mRootCell, new BoundedRayContainer(), 0,
409                                                                 mBox.SurfaceArea(), new BspNodeGeometry()));
410
411        while (!tStack.empty())
412        {
413                // filter polygons donw the tree
414                BspTraversalData tData = tStack.top();
415            tStack.pop();
416                       
417                if (!tData.mNode->IsLeaf())
418                {
419                        BspInterior *interior = dynamic_cast<BspInterior *>(tData.mNode);
420
421                        //-- filter view cell polygons down the tree until a leaf is reached
422                        if (!tData.mPolygons->empty())
423                        {
424                                PolygonContainer *frontPolys = new PolygonContainer();
425                                PolygonContainer *backPolys = new PolygonContainer();
426                                PolygonContainer coincident;
427
428                                int splits = 0;
429               
430                                // split viewcell polygons with respect to split plane
431                                splits += interior->SplitPolygons(*tData.mPolygons,
432                                                                                                  *frontPolys,
433                                                                                                  *backPolys,
434                                                                                                  coincident);
435                               
436                                // extract view cells associated with the split polygons
437                                ViewCell *frontViewCell = mRootCell;
438                                ViewCell *backViewCell = mRootCell;
439                       
440                                BspTraversalData frontData(interior->GetFront(),
441                                                                                   frontPolys,
442                                                                                   tData.mDepth + 1,
443                                                                                   mRootCell,   
444                                                                                   tData.mRays,
445                                                                                   tData.mPvs,
446                                                                                   mBox.SurfaceArea(),
447                                                                                   new BspNodeGeometry());
448
449                                BspTraversalData backData(interior->GetBack(),
450                                                                                  backPolys,
451                                                                                  tData.mDepth + 1,
452                                                                                  mRootCell,   
453                                                                                  tData.mRays,
454                                                                                  tData.mPvs,
455                                                                                  mBox.SurfaceArea(),
456                                                                                  new BspNodeGeometry());
457
458                                if (!mGenerateViewCells)
459                                {
460                                        ExtractViewCells(frontData,
461                                                                         backData,
462                                                                         coincident,
463                                                                         interior->mPlane);
464                                }
465
466                                // don't need coincident polygons anymore
467                                CLEAR_CONTAINER(coincident);
468
469                                mStat.splits += splits;
470
471                                // push the children on the stack
472                                tStack.push(frontData);
473                                tStack.push(backData);
474                        }
475
476                        // cleanup
477                        DEL_PTR(tData.mPolygons);
478                }
479                else
480                {
481                        // reached leaf => subdivide current viewcell
482                        BspNode *subRoot = Subdivide(tStack, tData);
483                }
484        }
485}
486
487int BspTree::AddMeshToPolygons(Mesh *mesh, PolygonContainer &polys, MeshInstance *parent)
488{
489        FaceContainer::const_iterator fi;
490       
491        // copy the face data to polygons
492        for (fi = mesh->mFaces.begin(); fi != mesh->mFaces.end(); ++ fi)
493        {
494                Polygon3 *poly = new Polygon3((*fi), mesh);
495               
496                if (poly->Valid())
497                {
498                        poly->mParent = parent; // set parent intersectable
499                        polys.push_back(poly);
500                }
501                else
502                        DEL_PTR(poly);
503        }
504        return (int)mesh->mFaces.size();
505}
506
507int BspTree::AddToPolygonSoup(const ViewCellContainer &viewCells,
508                                                          PolygonContainer &polys,
509                                                          int maxObjects)
510{
511        int limit = (maxObjects > 0) ?
512                Min((int)viewCells.size(), maxObjects) : (int)viewCells.size();
513 
514        int polysSize = 0;
515
516        for (int i = 0; i < limit; ++ i)
517        {
518                if (viewCells[i]->GetMesh()) // copy the mesh data to polygons
519                {
520                        mBox.Include(viewCells[i]->GetBox()); // add to BSP tree aabb
521                        polysSize += AddMeshToPolygons(viewCells[i]->GetMesh(), polys, viewCells[i]);
522                }
523        }
524
525        return polysSize;
526}
527
528int BspTree::AddToPolygonSoup(const ObjectContainer &objects, PolygonContainer &polys, int maxObjects)
529{
530        int limit = (maxObjects > 0) ? Min((int)objects.size(), maxObjects) : (int)objects.size();
531 
532        for (int i = 0; i < limit; ++i)
533        {
534                Intersectable *object = objects[i];//*it;
535                Mesh *mesh = NULL;
536
537                switch (object->Type()) // extract the meshes
538                {
539                case Intersectable::MESH_INSTANCE:
540                        mesh = dynamic_cast<MeshInstance *>(object)->GetMesh();
541                        break;
542                case Intersectable::VIEW_CELL:
543                        mesh = dynamic_cast<ViewCell *>(object)->GetMesh();
544                        break;
545                        // TODO: handle transformed mesh instances
546                default:
547                        Debug << "intersectable type not supported" << endl;
548                        break;
549                }
550               
551        if (mesh) // copy the mesh data to polygons
552                {
553                        mBox.Include(object->GetBox()); // add to BSP tree aabb
554                        AddMeshToPolygons(mesh, polys, mRootCell);
555                }
556        }
557
558        return (int)polys.size();
559}
560
561void BspTree::Construct(const ViewCellContainer &viewCells)
562{
563        mStat.nodes = 1;
564        mBox.Initialize();      // initialise bsp tree bounding box
565
566        // copy view cell meshes into one big polygon soup
567        PolygonContainer *polys = new PolygonContainer();
568        mStat.polys = AddToPolygonSoup(viewCells, *polys);
569
570        // construct tree from the view cell polygons
571        Construct(polys, new BoundedRayContainer());
572}
573
574
575void BspTree::Construct(const ObjectContainer &objects)
576{
577        mStat.nodes = 1;
578        mBox.Initialize();      // initialise bsp tree bounding box
579       
580        PolygonContainer *polys = new PolygonContainer();
581       
582        // copy mesh instance polygons into one big polygon soup
583        mStat.polys = AddToPolygonSoup(objects, *polys);
584
585        // construct tree from polygon soup
586        Construct(polys, new BoundedRayContainer());
587}
588
589void BspTree::Construct(const RayContainer &sampleRays)
590{
591        mStat.nodes = 1;
592        mBox.Initialize();      // initialise BSP tree bounding box
593       
594        PolygonContainer *polys = new PolygonContainer();
595        BoundedRayContainer *rays = new BoundedRayContainer();
596
597        RayContainer::const_iterator rit, rit_end = sampleRays.end();
598
599        long startTime = GetTime();
600
601        Debug << "**** Extracting polygons from rays ****\n";
602
603        std::map<Face *, Polygon3 *> facePolyMap;
604
605        //-- extract polygons intersected by the rays
606        for (rit = sampleRays.begin(); rit != rit_end; ++ rit)
607        {
608                Ray *ray = *rit;
609       
610                // get ray-face intersection. Store polygon representing the rays together
611                // with rays intersecting the face.
612                if (!ray->intersections.empty())
613                {
614                        MeshInstance *obj = dynamic_cast<MeshInstance *>(ray->intersections[0].mObject);
615                        Face *face = obj->GetMesh()->mFaces[ray->intersections[0].mFace];
616
617                        std::map<Face *, Polygon3 *>::iterator it = facePolyMap.find(face);
618
619                        if (it != facePolyMap.end())
620                        {
621                                //store rays if needed for heuristics
622                                if (sSplitPlaneStrategy & BLOCKED_RAYS)
623                                        (*it).second->mPiercingRays.push_back(ray);
624                        }
625                        else
626                        {       //store rays if needed for heuristics
627                                Polygon3 *poly = new Polygon3(face, obj->GetMesh());
628                                poly->mParent = obj;
629                                polys->push_back(poly);
630
631                                if (sSplitPlaneStrategy & BLOCKED_RAYS)
632                                        poly->mPiercingRays.push_back(ray);
633
634                                facePolyMap[face] = poly;
635                        }
636                }
637        }
638       
639        facePolyMap.clear();
640
641        // compue bounding box
642        Polygon3::IncludeInBox(*polys, mBox);
643
644        //-- store rays
645        for (rit = sampleRays.begin(); rit != rit_end; ++ rit)
646        {
647                Ray *ray = *rit;
648        ray->SetId(-1); // reset id
649
650                float minT, maxT;
651                if (BoundRay(*ray, minT, maxT))
652                        rays->push_back(new BoundedRay(ray, minT, maxT));
653        }
654
655        mStat.polys = (int)polys->size();
656
657        Debug << "**** Finished polygon extraction ****" << endl;
658        Debug << (int)polys->size() << " polys extracted from " << (int)sampleRays.size() << " rays" << endl;
659        Debug << "extraction time: " << TimeDiff(startTime, GetTime())*1e-3 << "s" << endl;
660
661        Construct(polys, rays);
662}
663
664void BspTree::Construct(PolygonContainer *polys, BoundedRayContainer *rays)
665{
666        std::stack<BspTraversalData> tStack;
667
668        mRoot = new BspLeaf();
669
670        BspNodeGeometry *cell = new BspNodeGeometry();
671        ConstructGeometry(mRoot, *cell);
672
673        BspTraversalData tData(mRoot, polys, 0, mRootCell, rays,
674                                                   ComputePvsSize(*rays), cell->GetArea(), cell);
675
676        tStack.push(tData);
677
678        long startTime = GetTime();
679        cout << "**** Contructing bsp tree ****\n";
680
681        while (!tStack.empty())
682        {
683                tData = tStack.top();
684            tStack.pop();
685
686                // subdivide leaf node
687                BspNode *subRoot = Subdivide(tStack, tData);
688        }
689
690        cout << "**** Finished tree construction ****\n";
691        Debug << "BSP tree contruction time: " << TimeDiff(startTime, GetTime())*1e-3 << "s" << endl;
692}
693
694
695BspNode *BspTree::Subdivide(BspTraversalStack &tStack, BspTraversalData &tData)
696{
697        //-- terminate traversal 
698        if (((int)tData.mPolygons->size() <= sTermMaxPolygons) ||
699                ((int)tData.mRays->size() <= sTermMaxRays) ||
700                (tData.mPvs <= sTermMinPvs) ||
701                 (tData.mDepth >= sTermMaxDepth))
702               
703        {
704                BspLeaf *leaf = dynamic_cast<BspLeaf *>(tData.mNode);
705       
706                BspViewCell *viewCell;
707
708                // generate new view cell for each leaf
709                if (mGenerateViewCells)
710                {
711                        viewCell = dynamic_cast<BspViewCell *>(ViewCell::Generate());
712                }
713                else
714                {
715                        // add view cell to leaf
716                        viewCell = dynamic_cast<BspViewCell *>(tData.mViewCell);
717                }
718
719                leaf->SetViewCell(viewCell);
720                viewCell->mBspLeaves.push_back(leaf);
721
722                //-- add pvs
723                if (viewCell != mRootCell)
724                {
725                        int conSamp = 0, sampCon = 0;
726                        leaf->AddToPvs(*tData.mRays, conSamp, sampCon);
727                       
728                        mStat.contributingSamples += conSamp;
729                        mStat.sampleContributions += sampCon;
730                }
731
732                EvaluateLeafStats(tData);
733               
734                //-- clean up
735               
736                // discard polygons
737                CLEAR_CONTAINER(*tData.mPolygons);
738                // discard rays
739                CLEAR_CONTAINER(*tData.mRays);
740
741                DEL_PTR(tData.mPolygons);
742                DEL_PTR(tData.mRays);
743                DEL_PTR(tData.mCell);
744               
745                return leaf;
746        }
747
748        //-- continue subdivision
749
750        PolygonContainer coincident;
751       
752        PolygonContainer *frontPolys = new PolygonContainer();
753        PolygonContainer *backPolys = new PolygonContainer();
754
755        BoundedRayContainer *frontRays = new BoundedRayContainer();
756        BoundedRayContainer *backRays = new BoundedRayContainer();
757       
758        BspTraversalData tFrontData(NULL, new PolygonContainer(), tData.mDepth + 1, mRootCell,
759                                                                new BoundedRayContainer(), 0, NULL, 0);
760        BspTraversalData tBackData(NULL, new PolygonContainer(), tData.mDepth + 1, mRootCell,
761                                                           new BoundedRayContainer(), 0, NULL, 0);
762
763        // create new interior node and two leaf nodes
764        BspInterior *interior = SubdivideNode(tData,
765                                                                                  tFrontData,
766                                                                                  tBackData,
767                                                                                  coincident);
768
769#ifdef _DEBUG   
770        if (frontPolys->empty() && backPolys->empty() && (coincident.size() > 2))
771        {       for (PolygonContainer::iterator it = coincident.begin(); it != coincident.end(); ++it)
772                        Debug << (*it) << " " << (*it)->GetArea() << " " << (*it)->mParent << endl ;
773                Debug << endl;}
774#endif
775
776        // extract view cells from coincident polygons according to plane normal
777    // only if front or back polygons are empty
778        if (!mGenerateViewCells)
779        {
780                ExtractViewCells(tFrontData,
781                                                 tBackData,
782                                                 coincident,
783                                                 interior->mPlane);
784                                                 
785        }
786
787        // don't need coincident polygons anymory
788        CLEAR_CONTAINER(coincident);
789
790        // push the children on the stack
791        tStack.push(tFrontData);
792        tStack.push(tBackData);
793
794        // cleanup
795        DEL_PTR(tData.mNode);
796        DEL_PTR(tData.mPolygons);
797        DEL_PTR(tData.mRays);
798        DEL_PTR(tData.mCell);           
799
800        return interior;
801}
802
803void BspTree::ExtractViewCells(BspTraversalData &frontData,
804                                                           BspTraversalData &backData,
805                                                           const PolygonContainer &coincident,
806                                                           const Plane3 splitPlane) const
807{
808        // if not empty, tree is further subdivided => don't have to find view cell
809        bool foundFront = !frontData.mPolygons->empty();
810        bool foundBack = !frontData.mPolygons->empty();
811
812        PolygonContainer::const_iterator it =
813                coincident.begin(), it_end = coincident.end();
814
815        //-- find first view cells in front and back leafs
816        for (; !(foundFront && foundBack) && (it != it_end); ++ it)
817        {
818                if (DotProd((*it)->GetNormal(), splitPlane.mNormal) > 0)
819                {
820                        backData.mViewCell = dynamic_cast<ViewCell *>((*it)->mParent);
821                        foundBack = true;
822                }
823                else
824                {
825                        frontData.mViewCell = dynamic_cast<ViewCell *>((*it)->mParent);
826                        foundFront = true;
827                }
828        }
829}
830
831BspInterior *BspTree::SubdivideNode(BspTraversalData &tData,
832                                                                        BspTraversalData &frontData,
833                                                                        BspTraversalData &backData,
834                                                                        PolygonContainer &coincident)
835{
836        mStat.nodes += 2;
837       
838        BspLeaf *leaf = dynamic_cast<BspLeaf *>(tData.mNode);
839        // select subdivision plane
840        BspInterior *interior =
841                new BspInterior(SelectPlane(leaf, tData, frontData, backData));
842
843#ifdef _DEBUG
844        Debug << interior << endl;
845#endif
846       
847        // subdivide rays into front and back rays
848        SplitRays(interior->mPlane, *tData.mRays, *frontData.mRays, *backData.mRays);
849       
850        // subdivide polygons with plane
851        mStat.splits += interior->SplitPolygons(*tData.mPolygons,
852                                                    *frontData.mPolygons,
853                                                                                        *backData.mPolygons,
854                                                                                        coincident);
855
856        BspInterior *parent = leaf->GetParent();
857
858        // replace a link from node's parent
859        if (!leaf->IsRoot())
860        {
861                parent->ReplaceChildLink(leaf, interior);
862                interior->SetParent(parent);
863        }
864        else // new root
865        {
866                mRoot = interior;
867        }
868
869        // and setup child links
870        interior->SetupChildLinks(new BspLeaf(interior), new BspLeaf(interior));
871       
872        frontData.mNode = interior->mFront;
873        backData.mNode = interior->mBack;
874       
875        return interior;
876}
877
878void BspTree::SortSplitCandidates(const PolygonContainer &polys,
879                                                                  const int axis,
880                                                                  vector<SortableEntry> &splitCandidates) const
881{
882  splitCandidates.clear();
883 
884  int requestedSize = 2 * (int)polys.size();
885  // creates a sorted split candidates array 
886  splitCandidates.reserve(requestedSize);
887 
888  PolygonContainer::const_iterator it, it_end = polys.end();
889
890  AxisAlignedBox3 box;
891
892  // insert all queries
893  for(it = polys.begin(); it != it_end; ++ it)
894  {
895          box.Initialize();
896          box.Include(*(*it));
897         
898          splitCandidates.push_back(SortableEntry(SortableEntry::POLY_MIN, box.Min(axis), *it));
899      splitCandidates.push_back(SortableEntry(SortableEntry::POLY_MAX, box.Max(axis), *it));
900  }
901 
902  stable_sort(splitCandidates.begin(), splitCandidates.end());
903}
904
905
906float BspTree::BestCostRatio(const PolygonContainer &polys,
907                                                         const AxisAlignedBox3 &box,
908                                                         const int axis,
909                                                         float &position,
910                                                         int &objectsBack,
911                                                         int &objectsFront) const
912{
913        vector<SortableEntry> splitCandidates;
914
915        SortSplitCandidates(polys, axis, splitCandidates);
916       
917        // go through the lists, count the number of objects left and right
918        // and evaluate the following cost funcion:
919        // C = ct_div_ci  + (ol + or)/queries
920       
921        int objectsLeft = 0, objectsRight = (int)polys.size();
922       
923        float minBox = box.Min(axis);
924        float maxBox = box.Max(axis);
925        float boxArea = box.SurfaceArea();
926 
927        float minBand = minBox + sSplitBorder * (maxBox - minBox);
928        float maxBand = minBox + (1.0f - sSplitBorder) * (maxBox - minBox);
929       
930        float minSum = 1e20f;
931        vector<SortableEntry>::const_iterator ci, ci_end = splitCandidates.end();
932
933        for(ci = splitCandidates.begin(); ci != ci_end; ++ ci)
934        {
935                switch ((*ci).type)
936                {
937                        case SortableEntry::POLY_MIN:
938                                ++ objectsLeft;
939                                break;
940                        case SortableEntry::POLY_MAX:
941                            -- objectsRight;
942                                break;
943                        default:
944                                break;
945                }
946               
947                if ((*ci).value > minBand && (*ci).value < maxBand)
948                {
949                        AxisAlignedBox3 lbox = box;
950                        AxisAlignedBox3 rbox = box;
951                        lbox.SetMax(axis, (*ci).value);
952                        rbox.SetMin(axis, (*ci).value);
953
954                        float sum = objectsLeft * lbox.SurfaceArea() +
955                                                objectsRight * rbox.SurfaceArea();
956     
957                        if (sum < minSum)
958                        {
959                                minSum = sum;
960                                position = (*ci).value;
961
962                                objectsBack = objectsLeft;
963                                objectsFront = objectsRight;
964                        }
965                }
966        }
967 
968        float oldCost = (float)polys.size();
969        float newCost = sCt_div_ci + minSum / boxArea;
970        float ratio = newCost / oldCost;
971
972
973#if 0
974  Debug << "====================" << endl;
975  Debug << "costRatio=" << ratio << " pos=" << position<<" t=" << (position - minBox)/(maxBox - minBox)
976      << "\t o=(" << objectsBack << "," << objectsFront << ")" << endl;
977#endif
978  return ratio;
979}
980
981bool BspTree::SelectAxisAlignedPlane(Plane3 &plane,
982                                                                         const PolygonContainer &polys) const
983{
984        AxisAlignedBox3 box;
985        box.Initialize();
986       
987        // create bounding box of region
988        Polygon3::IncludeInBox(polys, box);
989       
990        int objectsBack = 0, objectsFront = 0;
991        int axis = 0;
992        float costRatio = MAX_FLOAT;
993        Vector3 position;
994
995        //-- area subdivision
996        for (int i = 0; i < 3; ++ i)
997        {
998                float p = 0;
999                float r = BestCostRatio(polys, box, i, p, objectsBack, objectsFront);
1000               
1001                if (r < costRatio)
1002                {
1003                        costRatio = r;
1004                        axis = i;
1005                        position = p;
1006                }
1007        }
1008       
1009        if (costRatio >= sMaxCostRatio)
1010                return false;
1011
1012        Vector3 norm(0,0,0); norm[axis] = 1.0f;
1013        plane = Plane3(norm, position);
1014
1015        return true;
1016}
1017
1018Plane3 BspTree::SelectPlane(BspLeaf *leaf,
1019                                                        BspTraversalData &data,
1020                                                        BspTraversalData &frontData,
1021                                                        BspTraversalData &backData)
1022{
1023        if (data.mPolygons->empty() && data.mRays->empty())
1024        {
1025                Debug << "Warning: No autopartition polygon candidate available\n";
1026       
1027                // return axis aligned split
1028                AxisAlignedBox3 box;
1029                box.Initialize();
1030       
1031                // create bounding box of region
1032                Polygon3::IncludeInBox(*data.mPolygons, box);
1033
1034                const int axis = box.Size().DrivingAxis();
1035                const Vector3 position = (box.Min()[axis] + box.Max()[axis])*0.5f;
1036
1037                Vector3 norm(0,0,0); norm[axis] = 1.0f;
1038                return Plane3(norm, position);
1039        }
1040       
1041        if ((sSplitPlaneStrategy & AXIS_ALIGNED) &&
1042                ((int)data.mPolygons->size() > sTermMaxPolysForAxisAligned) &&
1043                ((int)data.mRays->size() > sTermMaxRaysForAxisAligned) &&
1044                ((sTermMaxObjectsForAxisAligned < 0) ||
1045                  (Polygon3::ParentObjectsSize(*data.mPolygons) > sTermMaxObjectsForAxisAligned)))
1046        {
1047                Plane3 plane;
1048                if (SelectAxisAlignedPlane(plane, *data.mPolygons))
1049                        return plane;
1050        }
1051
1052        // simplest strategy: just take next polygon
1053        if (sSplitPlaneStrategy & RANDOM_POLYGON)
1054        {
1055                Polygon3 *nextPoly = (*data.mPolygons)[Random((int)data.mPolygons->size())];
1056                return nextPoly->GetSupportingPlane();
1057        }
1058
1059        // use heuristics to find appropriate plane
1060        return SelectPlaneHeuristics(leaf, data, frontData, backData);
1061}
1062
1063Plane3 BspTree::SelectPlaneHeuristics(BspLeaf *leaf,
1064                                                                          BspTraversalData &data,
1065                                                                          BspTraversalData &frontData,
1066                                                                          BspTraversalData &backData)
1067{
1068        float lowestCost = MAX_FLOAT;
1069        Plane3 bestPlane;
1070       
1071        int limit = sMaxPolyCandidates > 0 ?
1072                Min((int)data.mPolygons->size(), sMaxPolyCandidates) : (int)data.mPolygons->size();
1073       
1074        int candidateIdx = limit;
1075
1076        for (int i = 0; i < limit; ++ i)
1077        {
1078                candidateIdx = GetNextCandidateIdx(candidateIdx, *data.mPolygons);
1079               
1080                Polygon3 *poly = (*data.mPolygons)[candidateIdx];
1081                // evaluate current candidate
1082                const float candidateCost =
1083                        SplitPlaneCost(poly->GetSupportingPlane(), data, frontData, backData);
1084
1085                Debug << "cost: " << candidateCost << ", lowest: " << lowestCost << endl;
1086
1087                if (candidateCost < lowestCost)
1088                {
1089                        bestPlane = poly->GetSupportingPlane();
1090                        lowestCost = candidateCost;
1091                }
1092        }
1093       
1094        //limit = maxRaysCandidates > 0 ? Min((int)rays.size(), maxRayCandidates) : (int)rays.size();   
1095        const BoundedRayContainer *rays = data.mRays;
1096
1097        for (int i = 0; i < sMaxRayCandidates; ++ i)
1098        {
1099                Plane3 plane;
1100
1101                if (1)
1102                {
1103                        Vector3 pt[3];
1104                        int idx[3];
1105                       
1106                        for (int j = 0; j < 3; j ++)
1107                        {
1108                                idx[j] = Random((int)rays->size() * 2);
1109                               
1110                                BoundedRay *bRay = (*rays)[idx[j]];
1111                                Ray *ray = bRay->mRay;
1112
1113                                if (idx[j] < (int)rays->size())
1114                                        pt[j] = ray->Extrap(bRay->mMinT);
1115                                else
1116                                {
1117                                        idx[j] -= (int)rays->size();
1118                                        pt[j] = ray->Extrap(bRay->mMaxT);
1119                                }
1120                        }       
1121               
1122                        plane = Plane3(pt[0], pt[1], pt[2]);
1123                }
1124                else
1125                {
1126                        candidateIdx = Random((int)rays->size());
1127                        BoundedRay *bRay = (*rays)[candidateIdx];
1128
1129                        Ray *ray = bRay->mRay;
1130                                               
1131                        const Vector3 minPt = ray->Extrap(bRay->mMinT);
1132                        const Vector3 maxPt = ray->Extrap(bRay->mMaxT);
1133
1134                        const Vector3 pt = (maxPt + minPt) * 0.5;
1135
1136                        const Vector3 normal = ray->GetDir();
1137                       
1138                        plane = Plane3(normal, pt);
1139                }
1140               
1141                const float candidateCost =
1142                        SplitPlaneCost(plane, data, frontData, backData);
1143
1144                if (candidateCost < lowestCost)
1145                {
1146                        Debug << "choose ray plane: " << lowestCost << endl;
1147                        bestPlane = plane;
1148                       
1149                        lowestCost = candidateCost;
1150                }
1151                else
1152                        Debug << "ray cost: " << candidateCost << ", lowest: " << lowestCost << endl;
1153        }
1154
1155        //Debug << "Plane lowest cost: " << lowestCost << endl;
1156        return bestPlane;
1157}
1158
1159int BspTree::GetNextCandidateIdx(int currentIdx, PolygonContainer &polys)
1160{
1161        const int candidateIdx = Random(currentIdx --);
1162
1163        // swap candidates to avoid testing same plane 2 times
1164        std::swap(polys[currentIdx], polys[candidateIdx]);
1165       
1166        return currentIdx;
1167        //return Random((int)polys.size());
1168}
1169
1170float BspTree::SplitPlaneCost(const Plane3 &candidatePlane,
1171                                                          const PolygonContainer &polys) const
1172{
1173        float val = 0;
1174
1175        float sumBalancedPolys = 0;
1176        float sumSplits = 0;
1177        float sumPolyArea = 0;
1178        float sumBalancedViewCells = 0;
1179        float sumBlockedRays = 0;
1180        float totalBlockedRays = 0;
1181        //float totalArea = 0;
1182        int totalViewCells = 0;
1183
1184        // need three unique ids for each type of view cell
1185        // for balanced view cells criterium
1186        ViewCell::NewMail();
1187        const int backId = ViewCell::sMailId;
1188        ViewCell::NewMail();
1189        const int frontId = ViewCell::sMailId;
1190        ViewCell::NewMail();
1191        const int frontAndBackId = ViewCell::sMailId;
1192
1193        PolygonContainer::const_iterator it, it_end = polys.end();
1194
1195        for (it = polys.begin(); it != it_end; ++ it)
1196        {
1197                const int classification = (*it)->ClassifyPlane(candidatePlane);
1198
1199                if (sSplitPlaneStrategy & BALANCED_POLYS)
1200                        sumBalancedPolys += sBalancedPolysTable[classification];
1201               
1202                if (sSplitPlaneStrategy & LEAST_SPLITS)
1203                        sumSplits += sLeastPolySplitsTable[classification];
1204
1205                if (sSplitPlaneStrategy & LARGEST_POLY_AREA)
1206                {
1207                        if (classification == Polygon3::COINCIDENT)
1208                                sumPolyArea += (*it)->GetArea();
1209                        //totalArea += area;
1210                }
1211               
1212                if (sSplitPlaneStrategy & BLOCKED_RAYS)
1213                {
1214                        const float blockedRays = (float)(*it)->mPiercingRays.size();
1215               
1216                        if (classification == Polygon3::COINCIDENT)
1217                                sumBlockedRays += blockedRays;
1218                       
1219                        totalBlockedRays += blockedRays;
1220                }
1221
1222                // assign view cells to back or front according to classificaion
1223                if (sSplitPlaneStrategy & BALANCED_VIEW_CELLS)
1224                {
1225                        MeshInstance *viewCell = (*it)->mParent;
1226               
1227                        // assure that we only count a view cell
1228                        // once for the front and once for the back side of the plane
1229                        if (classification == Polygon3::FRONT_SIDE)
1230                        {
1231                                if ((viewCell->mMailbox != frontId) &&
1232                                        (viewCell->mMailbox != frontAndBackId))
1233                                {
1234                                        sumBalancedViewCells += 1.0;
1235
1236                                        if (viewCell->mMailbox != backId)
1237                                                viewCell->mMailbox = frontId;
1238                                        else
1239                                                viewCell->mMailbox = frontAndBackId;
1240                                       
1241                                        ++ totalViewCells;
1242                                }
1243                        }
1244                        else if (classification == Polygon3::BACK_SIDE)
1245                        {
1246                                if ((viewCell->mMailbox != backId) &&
1247                                    (viewCell->mMailbox != frontAndBackId))
1248                                {
1249                                        sumBalancedViewCells -= 1.0;
1250
1251                                        if (viewCell->mMailbox != frontId)
1252                                                viewCell->mMailbox = backId;
1253                                        else
1254                                                viewCell->mMailbox = frontAndBackId;
1255
1256                                        ++ totalViewCells;
1257                                }
1258                        }
1259                }
1260        }
1261
1262        // all values should be approx. between 0 and 1 so they can be combined
1263        // and scaled with the factors according to their importance
1264        if (sSplitPlaneStrategy & BALANCED_POLYS)
1265                val += sBalancedPolysFactor * fabs(sumBalancedPolys) / (float)polys.size();
1266       
1267        if (sSplitPlaneStrategy & LEAST_SPLITS)   
1268                val += sLeastSplitsFactor * sumSplits / (float)polys.size();
1269
1270        if (sSplitPlaneStrategy & LARGEST_POLY_AREA)
1271                // HACK: polys.size should be total area so scaling is between 0 and 1
1272                val += sLargestPolyAreaFactor * (float)polys.size() / sumPolyArea;
1273
1274        if (sSplitPlaneStrategy & BLOCKED_RAYS)
1275                if (totalBlockedRays != 0)
1276                        val += sBlockedRaysFactor * (totalBlockedRays - sumBlockedRays) / totalBlockedRays;
1277
1278        if (sSplitPlaneStrategy & BALANCED_VIEW_CELLS)
1279                if (totalViewCells != 0)
1280                        val += sBalancedViewCellsFactor * fabs(sumBalancedViewCells) / (float)totalViewCells;
1281       
1282        return val;
1283}
1284
1285bool BspTree::BoundRay(const Ray &ray, float &minT, float &maxT) const
1286{
1287        maxT = 1e6;
1288        minT = 0;
1289       
1290        // test with tree bounding box
1291        if (!mBox.GetMinMaxT(ray, &minT, &maxT))
1292                return false;
1293
1294        if (minT < 0) // start ray from origin
1295                minT = 0;
1296
1297        // bound ray or line segment
1298        if ((ray.GetType() == Ray::LOCAL_RAY) &&
1299            !ray.intersections.empty() &&
1300                (ray.intersections[0].mT <= maxT))
1301        {
1302                maxT = ray.intersections[0].mT;
1303        }
1304
1305        return true;
1306}
1307
1308float BspTree::SplitPlaneCost(const Plane3 &candidatePlane,
1309                                                          const BoundedRayContainer &rays,
1310                                                          const int pvs,
1311                                                          const float area,
1312                                                          const BspNodeGeometry &cell,
1313                                                          BspTraversalData &frontData,
1314                                                          BspTraversalData &backData) const
1315{
1316        float val = 0;
1317
1318        float sumBalancedRays = 0;
1319        float sumRaySplits = 0;
1320       
1321        int backId = 0;
1322        int frontId = 0;
1323        int frontAndBackId = 0;
1324
1325        frontData.mPvs = 0;
1326        backData.mPvs = 0;
1327        frontData.mArea = 0;
1328        backData.mArea = 0;
1329
1330        // probability that view point lies in child
1331        float pOverall = 0;
1332        float pFront = 0;
1333        float pBack = 0;
1334
1335        if (sSplitPlaneStrategy & PVS)
1336        {
1337                // create three unique ids for pvs heuristics
1338                Intersectable::NewMail(); backId = ViewCell::sMailId;
1339                Intersectable::NewMail(); frontId = ViewCell::sMailId;
1340                Intersectable::NewMail(); frontAndBackId = ViewCell::sMailId;
1341
1342                if (1) // use front and back cell areas to approximate volume
1343                {       
1344                        // construct child geometry with regard to the candidate split plane
1345                        frontData.mCell = new BspNodeGeometry();
1346                        backData.mCell = new BspNodeGeometry();
1347
1348                        cell.SplitGeometry(*frontData.mCell, *backData.mCell, *this, candidatePlane);
1349               
1350                        pFront = frontData.mArea = frontData.mCell->GetArea();
1351                        pBack = backData.mArea = backData.mCell->GetArea();
1352
1353                        pOverall = area;
1354                }
1355                else
1356                {
1357                        pOverall = (float)rays.size();
1358                }
1359        }
1360                       
1361        BoundedRayContainer::const_iterator rit, rit_end = rays.end();
1362
1363        for (rit = rays.begin(); rit != rays.end(); ++ rit)
1364        {
1365                Ray *ray = (*rit)->mRay;
1366                const float minT = (*rit)->mMinT;
1367                const float maxT = (*rit)->mMaxT;
1368
1369                const int cf =
1370                        ray->ClassifyPlane(candidatePlane, minT, maxT);
1371
1372                if (sSplitPlaneStrategy & LEAST_RAY_SPLITS)
1373                {
1374                        sumBalancedRays += sBalancedRaysTable[cf];
1375                }
1376               
1377                if (sSplitPlaneStrategy & BALANCED_RAYS)
1378                {
1379                        sumRaySplits += sLeastRaySplitsTable[cf];
1380                }
1381
1382                if (sSplitPlaneStrategy & PVS)
1383                {
1384                        if (!ray->intersections.empty())
1385                        {
1386                                // in case the ray intersects an objcrs
1387                                // assure that we only count a object
1388                                // once for the front and once for the back side of the plane
1389                                AddToPvs(*ray->intersections[0].mObject, frontData.mPvs, backData.mPvs,
1390                                                 cf, frontId, backId, frontAndBackId);
1391                        }
1392
1393                        // the source object in the origin of the ray
1394                        AddToPvs(*ray->sourceObject.mObject, frontData.mPvs, backData.mPvs,
1395                                         cf, frontId, backId, frontAndBackId);
1396
1397                        // use #rays to approximate volume
1398                        if (0)
1399                        {
1400                                if ((cf == Ray::BACK) || (cf == Ray::FRONT_BACK) || (cf == Ray::BACK_FRONT))
1401                                        ++ pFront;
1402                                if ((cf == Ray::FRONT) || (cf == Ray::FRONT_BACK) || (cf == Ray::BACK_FRONT))
1403                                        ++ pBack;
1404                        }
1405                }
1406        }
1407
1408        if ((sSplitPlaneStrategy & LEAST_RAY_SPLITS) && !rays.empty())
1409                        val += sLeastRaySplitsFactor * sumRaySplits / (float)rays.size();
1410
1411        if ((sSplitPlaneStrategy & BALANCED_RAYS) && !rays.empty())
1412                        val += sBalancedRaysFactor * fabs(sumBalancedRays) / (float)rays.size();
1413
1414        if ((sSplitPlaneStrategy & PVS) && area && pvs)
1415                val += sPvsFactor * (frontData.mPvs * pFront + (backData.mPvs * pBack)) /
1416                           (pOverall * (float)pvs * 2);
1417
1418        Debug << "totalpvs: " << pvs << " ptotal: " << pOverall
1419                  << " frontpvs: " << frontData.mPvs << " pFront: " << pFront
1420                  << " backpvs: " << backData.mPvs << " pBack: " << pBack
1421                  << " val " << val << endl;
1422
1423        return val;
1424}
1425
1426void BspTree::AddToPvs(Intersectable &obj,
1427                                           int &frontPvs,
1428                                           int &backPvs,
1429                                           const int cf,
1430                                           const int frontId,
1431                                           const int backId,
1432                                           const int frontAndBackId) const
1433{
1434        if (cf == Ray::COINCIDENT) //TODO: really belongs to no pvs?
1435                return;
1436
1437        if (cf == Ray::FRONT)
1438        {
1439                if ((obj.mMailbox != frontId) &&
1440                        (obj.mMailbox != frontAndBackId))
1441                {
1442                        ++ frontPvs;
1443
1444                        if (obj.mMailbox != backId)
1445                                obj.mMailbox = frontId;
1446                        else
1447                                obj.mMailbox = frontAndBackId;                                 
1448                }
1449        }
1450        else if (cf == Ray::BACK)
1451        {
1452                if ((obj.mMailbox != backId) &&
1453                        (obj.mMailbox != frontAndBackId))
1454                {
1455                        ++ backPvs;
1456
1457                        if (obj.mMailbox != frontId)
1458                                obj.mMailbox = backId;
1459                        else
1460                                obj.mMailbox = frontAndBackId;
1461                }
1462        }
1463        // object belongs to both pvs
1464        else if ((cf == Ray::FRONT_BACK) || (cf == Ray::BACK_FRONT))
1465        {
1466                if (obj.mMailbox !=  frontAndBackId)
1467                {
1468                        if (obj.mMailbox != frontId)
1469                                ++ frontPvs;
1470                        if (obj.mMailbox != backId)
1471                                ++ backPvs;
1472               
1473                        obj.mMailbox = frontAndBackId;
1474                }
1475        }
1476}
1477
1478float BspTree::SplitPlaneCost(const Plane3 &candidatePlane,
1479                                                          BspTraversalData &data,
1480                                                          BspTraversalData &frontData,
1481                                                          BspTraversalData &backData) const
1482{
1483        float val = 0;
1484
1485        if (sSplitPlaneStrategy & VERTICAL_AXIS)
1486        {
1487                Vector3 tinyAxis(0,0,0); tinyAxis[mBox.Size().TinyAxis()] = 1.0f;
1488                // we put a penalty on the dot product between the "tiny" vertical axis
1489                // and the split plane axis
1490                val += sVerticalSplitsFactor *
1491                           fabs(DotProd(candidatePlane.mNormal, tinyAxis));
1492        }
1493
1494        // the following criteria loop over all polygons to find the cost value
1495        if ((sSplitPlaneStrategy & BALANCED_POLYS)      ||
1496                (sSplitPlaneStrategy & LEAST_SPLITS)        ||
1497                (sSplitPlaneStrategy & LARGEST_POLY_AREA)   ||
1498                (sSplitPlaneStrategy & BALANCED_VIEW_CELLS) ||
1499                (sSplitPlaneStrategy & BLOCKED_RAYS))
1500        {
1501                val += SplitPlaneCost(candidatePlane, *data.mPolygons);
1502        }
1503
1504        // the following criteria loop over all rays to find the cost value
1505        if ((sSplitPlaneStrategy & BALANCED_RAYS)      ||
1506                (sSplitPlaneStrategy & LEAST_RAY_SPLITS)   ||
1507                (sSplitPlaneStrategy & PVS))
1508        {
1509                val += SplitPlaneCost(candidatePlane, *data.mRays, data.mPvs,
1510                                                          data.mArea, *data.mCell, frontData, backData);
1511        }
1512
1513        // return linear combination of the sums
1514        return val;
1515}
1516
1517void BspTree::ParseEnvironment()
1518{
1519        //-- parse bsp cell tree construction method
1520        char constructionMethodStr[60];
1521       
1522        environment->GetStringValue("BspTree.Construction.input", constructionMethodStr);
1523
1524        sConstructionMethod = FROM_INPUT_VIEW_CELLS;
1525       
1526        if (strcmp(constructionMethodStr, "fromViewCells") == 0)
1527                sConstructionMethod = FROM_INPUT_VIEW_CELLS;
1528        else if (strcmp(constructionMethodStr, "fromSceneGeometry") == 0)
1529                sConstructionMethod = FROM_SCENE_GEOMETRY;
1530        else if (strcmp(constructionMethodStr, "fromRays") == 0)
1531                sConstructionMethod = FROM_RAYS;
1532        else
1533        {
1534                cerr << "Wrong construction method " << constructionMethodStr << endl;
1535                exit(1);
1536    }
1537
1538        Debug << "Construction method: " << constructionMethodStr << endl;
1539
1540        //-- termination criteria for autopartition
1541        environment->GetIntValue("BspTree.Termination.maxDepth", sTermMaxDepth);
1542        environment->GetIntValue("BspTree.Termination.minPvs", sTermMinPvs);
1543        environment->GetIntValue("BspTree.Termination.maxPolygons", sTermMaxPolygons);
1544        environment->GetIntValue("BspTree.Termination.maxRays", sTermMaxRays);
1545
1546        //-- termination criteria for axis aligned split
1547        environment->GetFloatValue("BspTree.Termination.AxisAligned.ct_div_ci", sCt_div_ci);
1548        environment->GetFloatValue("BspTree.Termination.AxisAligned.maxCostRatio", sMaxCostRatio);
1549        environment->GetIntValue("BspTree.Termination.AxisAligned.maxPolys",
1550                                                         sTermMaxPolysForAxisAligned);
1551        environment->GetIntValue("BspTree.Termination.AxisAligned.maxRays",
1552                                                         sTermMaxRaysForAxisAligned);
1553        environment->GetIntValue("BspTree.Termination.AxisAligned.maxObjects",
1554                                                         sTermMaxObjectsForAxisAligned);
1555        //-- partition criteria
1556        environment->GetIntValue("BspTree.maxPolyCandidates", sMaxPolyCandidates);
1557        environment->GetIntValue("BspTree.maxRayCandidates", sMaxRayCandidates);
1558        environment->GetIntValue("BspTree.splitPlaneStrategy", sSplitPlaneStrategy);
1559        environment->GetFloatValue("BspTree.AxisAligned.splitBorder", sSplitBorder);
1560
1561        environment->GetFloatValue("BspTree.Construction.sideTolerance", Vector3::sDistTolerance);
1562        Vector3::sDistToleranceSqrt = Vector3::sDistTolerance * Vector3::sDistTolerance;
1563
1564        // post processing stuff
1565        environment->GetIntValue("ViewCells.PostProcessing.minPvsDif", sMinPvsDif);
1566        environment->GetIntValue("ViewCells.PostProcessing.minPvs", sMinPvs);
1567        environment->GetIntValue("ViewCells.PostProcessing.maxPvs", sMaxPvs);
1568
1569    Debug << "BSP max depth: " << sTermMaxDepth << endl;
1570        Debug << "BSP min PVS: " << sTermMinPvs << endl;
1571        Debug << "BSP max polys: " << sTermMaxPolygons << endl;
1572        Debug << "BSP max rays: " << sTermMaxRays << endl;
1573        Debug << "BSP max polygon candidates: " << sMaxPolyCandidates << endl;
1574        Debug << "BSP max plane candidates: " << sMaxRayCandidates << endl;
1575
1576        Debug << "Split plane strategy: ";
1577        if (sSplitPlaneStrategy & RANDOM_POLYGON)
1578                Debug << "random polygon ";
1579        if (sSplitPlaneStrategy & AXIS_ALIGNED)
1580                Debug << "axis aligned ";
1581        if (sSplitPlaneStrategy & LEAST_SPLITS)     
1582                Debug << "least splits ";
1583        if (sSplitPlaneStrategy & BALANCED_POLYS)
1584                Debug << "balanced polygons ";
1585        if (sSplitPlaneStrategy & BALANCED_VIEW_CELLS)
1586                Debug << "balanced view cells ";
1587        if (sSplitPlaneStrategy & LARGEST_POLY_AREA)
1588                Debug << "largest polygon area ";
1589        if (sSplitPlaneStrategy & VERTICAL_AXIS)
1590                Debug << "vertical axis ";
1591        if (sSplitPlaneStrategy & BLOCKED_RAYS)
1592                Debug << "blocked rays ";
1593        if (sSplitPlaneStrategy & LEAST_RAY_SPLITS)
1594                Debug << "least ray splits ";
1595        if (sSplitPlaneStrategy & BALANCED_RAYS)
1596                Debug << "balanced rays ";
1597        if (sSplitPlaneStrategy & PVS)
1598                Debug << "pvs";
1599        Debug << endl;
1600}
1601
1602void BspTree::CollectLeaves(vector<BspLeaf *> &leaves) const
1603{
1604        stack<BspNode *> nodeStack;
1605        nodeStack.push(mRoot);
1606 
1607        while (!nodeStack.empty())
1608        {
1609                BspNode *node = nodeStack.top();
1610   
1611                nodeStack.pop();
1612   
1613                if (node->IsLeaf())
1614                {
1615                        BspLeaf *leaf = (BspLeaf *)node;               
1616                        leaves.push_back(leaf);
1617                }
1618                else
1619                {
1620                        BspInterior *interior = dynamic_cast<BspInterior *>(node);
1621
1622                        nodeStack.push(interior->GetBack());
1623                        nodeStack.push(interior->GetFront());
1624                }
1625        }
1626}
1627
1628AxisAlignedBox3 BspTree::GetBoundingBox() const
1629{
1630        return mBox;
1631}
1632
1633BspNode *BspTree::GetRoot() const
1634{
1635        return mRoot;
1636}
1637
1638void BspTree::EvaluateLeafStats(const BspTraversalData &data)
1639{
1640        // the node became a leaf -> evaluate stats for leafs
1641        BspLeaf *leaf = dynamic_cast<BspLeaf *>(data.mNode);
1642
1643        if (data.mDepth >= sTermMaxDepth)
1644                ++ mStat.maxDepthNodes;
1645       
1646        // store maximal and minimal depth
1647        if (data.mDepth > mStat.maxDepth)
1648                mStat.maxDepth = data.mDepth;
1649
1650        if (data.mDepth < mStat.minDepth)
1651                mStat.minDepth = data.mDepth;
1652
1653        // store minimal and maximal pvs
1654        /*if (data.mPvs > mStat.pvs)
1655                mStat.pvs = data.mPvs;
1656       
1657        if (data.mPvs < mStat.pvs)
1658                mStat.pvs = data.mPvs;*/
1659
1660        // accumulate depth to compute average
1661        mStat.accumDepth += data.mDepth;
1662       
1663//#ifdef _DEBUG
1664        Debug << "BSP stats: "
1665                  << "Depth: " << data.mDepth << " (max: " << sTermMaxDepth << "), "
1666                   << "PVS: " << data.mPvs << " (max: " << sTermMinPvs << "), "
1667                  << "#polygons: " << (int)data.mPolygons->size() << " (max: " << sTermMaxPolygons << "), "
1668                  << "#rays: " << (int)data.mRays->size() << " (max: " << sTermMaxRays << "), "
1669                  << "#pvs: " << leaf->GetViewCell()->GetPvs().GetSize() << endl;
1670//#endif
1671}
1672
1673int BspTree::CastRay(Ray &ray)
1674{
1675        int hits = 0;
1676 
1677        stack<BspRayTraversalData> tStack;
1678 
1679        float maxt, mint;
1680
1681        if (!BoundRay(ray, mint, maxt))
1682                return 0;
1683
1684        Intersectable::NewMail();
1685
1686        Vector3 entp = ray.Extrap(mint);
1687        Vector3 extp = ray.Extrap(maxt);
1688 
1689        BspNode *node = mRoot;
1690        BspNode *farChild = NULL;
1691       
1692        while (1)
1693        {
1694                if (!node->IsLeaf())
1695                {
1696                        BspInterior *in = (BspInterior *) node;
1697                       
1698                        Plane3 *splitPlane = in->GetPlane();
1699
1700                        int entSide = splitPlane->Side(entp);
1701                        int extSide = splitPlane->Side(extp);
1702
1703                        Vector3 intersection;
1704
1705                        if (entSide < 0)
1706                        {
1707                                node = in->GetBack();
1708
1709                                if(extSide <= 0) // plane does not split ray => no far child
1710                                        continue;
1711                                       
1712                                farChild = in->GetFront(); // plane splits ray
1713
1714                        } else if (entSide > 0)
1715                        {
1716                                node = in->GetFront();
1717
1718                                if (extSide >= 0) // plane does not split ray => no far child
1719                                        continue;
1720
1721                                farChild = in->GetBack(); // plane splits ray                   
1722                        }
1723                        else // ray and plane are coincident
1724                        // WHAT TO DO IN THIS CASE ?
1725                        {
1726                                break;
1727                                //node = in->GetFront();
1728                                //continue;
1729                        }
1730
1731                        // push data for far child
1732                        tStack.push(BspRayTraversalData(farChild, extp, maxt));
1733
1734                        // find intersection of ray segment with plane
1735                        float t;
1736                        extp = splitPlane->FindIntersection(ray.GetLoc(), extp, &t);
1737                        maxt *= t;
1738                       
1739                } else // reached leaf => intersection with view cell
1740                {
1741                        BspLeaf *leaf = dynamic_cast<BspLeaf *>(node);
1742     
1743                        if (!leaf->mViewCell->Mailed())
1744                        {
1745                                ray.bspIntersections.push_back(Ray::BspIntersection(maxt, leaf));
1746                                leaf->mViewCell->Mail();
1747                                ++ hits;
1748                        }
1749                       
1750                        //-- fetch the next far child from the stack
1751                        if (tStack.empty())
1752                                break;
1753     
1754                        entp = extp;
1755                        mint = maxt; // NOTE: need this?
1756
1757                        if (ray.GetType() == Ray::LINE_SEGMENT && mint > 1.0f)
1758                                break;
1759
1760                        BspRayTraversalData &s = tStack.top();
1761
1762                        node = s.mNode;
1763                        extp = s.mExitPoint;
1764                        maxt = s.mMaxT;
1765
1766                        tStack.pop();
1767                }
1768        }
1769
1770        return hits;
1771}
1772
1773bool BspTree::Export(const string filename)
1774{
1775        Exporter *exporter = Exporter::GetExporter(filename);
1776
1777        if (exporter)
1778        {
1779                exporter->ExportBspTree(*this);
1780                return true;
1781        }       
1782
1783        return false;
1784}
1785
1786void BspTree::CollectViewCells(ViewCellContainer &viewCells) const
1787{
1788        stack<BspNode *> nodeStack;
1789        nodeStack.push(mRoot);
1790
1791        ViewCell::NewMail();
1792
1793        while (!nodeStack.empty())
1794        {
1795                BspNode *node = nodeStack.top();
1796                nodeStack.pop();
1797
1798                if (node->IsLeaf())
1799                {
1800                        ViewCell *viewCell = dynamic_cast<BspLeaf *>(node)->mViewCell;
1801
1802                        if (!viewCell->Mailed())
1803                        {
1804                                viewCell->Mail();
1805                                viewCells.push_back(viewCell);
1806                        }
1807                }
1808                else
1809                {
1810                        BspInterior *interior = dynamic_cast<BspInterior *>(node);
1811
1812                        nodeStack.push(interior->mFront);
1813                        nodeStack.push(interior->mBack);
1814                }
1815        }
1816}
1817
1818void BspTree::EvaluateViewCellsStats(BspViewCellsStatistics &stat) const
1819{
1820        stat.Reset();
1821
1822        stack<BspNode *> nodeStack;
1823        nodeStack.push(mRoot);
1824
1825        ViewCell::NewMail();
1826
1827        // exclude root cell
1828        mRootCell->Mail();
1829
1830        while (!nodeStack.empty())
1831        {
1832                BspNode *node = nodeStack.top();
1833                nodeStack.pop();
1834
1835                if (node->IsLeaf())
1836                {
1837                        ++ stat.bspLeaves;
1838
1839                        BspViewCell *viewCell = dynamic_cast<BspLeaf *>(node)->mViewCell;
1840
1841                        if (!viewCell->Mailed())
1842                        {
1843                                viewCell->Mail();
1844                               
1845                                ++ stat.viewCells;
1846                                int pvsSize = viewCell->GetPvs().GetSize();
1847
1848                stat.pvs += pvsSize;
1849
1850                                if (pvsSize < 2)
1851                                        ++ stat.emptyPvs;
1852
1853                                if (pvsSize > stat.maxPvs)
1854                                        stat.maxPvs = pvsSize;
1855
1856                                if (pvsSize < stat.minPvs)
1857                                        stat.minPvs = pvsSize;
1858
1859                                if ((int)viewCell->mBspLeaves.size() > stat.maxBspLeaves)
1860                                        stat.maxBspLeaves = (int)viewCell->mBspLeaves.size();           
1861                        }
1862                }
1863                else
1864                {
1865                        BspInterior *interior = dynamic_cast<BspInterior *>(node);
1866
1867                        nodeStack.push(interior->mFront);
1868                        nodeStack.push(interior->mBack);
1869                }
1870        }
1871}
1872
1873bool BspTree::MergeViewCells(BspLeaf *front, BspLeaf *back) const
1874{
1875        BspViewCell *viewCell =
1876                dynamic_cast<BspViewCell *>(ViewCell::Merge(*front->mViewCell, *back->mViewCell));
1877       
1878        if (!viewCell)
1879                return false;
1880
1881        // change view cells of all leaves associated with the
1882        // previous view cells
1883
1884        BspViewCell *fVc = front->mViewCell;
1885        BspViewCell *bVc = back->mViewCell;
1886
1887        vector<BspLeaf *> fLeaves = fVc->mBspLeaves;
1888        vector<BspLeaf *> bLeaves = bVc->mBspLeaves;
1889
1890        vector<BspLeaf *>::const_iterator it;
1891       
1892        for (it = fLeaves.begin(); it != fLeaves.end(); ++ it)
1893        {
1894                (*it)->SetViewCell(viewCell);
1895                viewCell->mBspLeaves.push_back(*it);
1896        }
1897        for (it = bLeaves.begin(); it != bLeaves.end(); ++ it)
1898        {
1899                (*it)->SetViewCell(viewCell);
1900                viewCell->mBspLeaves.push_back(*it);
1901        }
1902       
1903        DEL_PTR(fVc);
1904        DEL_PTR(bVc);
1905
1906        return true;
1907}
1908
1909bool BspTree::ShouldMerge(BspLeaf *front, BspLeaf *back) const
1910{
1911        ViewCell *fvc = front->mViewCell;
1912        ViewCell *bvc = back->mViewCell;
1913
1914        if ((fvc == mRootCell) || (bvc == mRootCell) || (fvc == bvc))
1915                return false;
1916
1917        int fdiff = fvc->GetPvs().Diff(bvc->GetPvs());
1918
1919        if (fvc->GetPvs().GetSize() + fdiff < sMaxPvs)
1920        {
1921                if ((fvc->GetPvs().GetSize() < sMinPvs) ||     
1922                        (bvc->GetPvs().GetSize() < sMinPvs) ||
1923                        ((fdiff < sMinPvsDif) && (bvc->GetPvs().Diff(fvc->GetPvs()) < sMinPvsDif)))
1924                {
1925                        return true;
1926                }
1927        }
1928       
1929        return false;
1930}
1931
1932void BspTree::SetGenerateViewCells(int generateViewCells)
1933{
1934        mGenerateViewCells = generateViewCells;
1935}
1936
1937BspTreeStatistics &BspTree::GetStat()
1938{
1939        return mStat;
1940}
1941
1942int BspTree::SplitRays(const Plane3 &plane,
1943                                           BoundedRayContainer &rays,
1944                                           BoundedRayContainer &frontRays,
1945                                           BoundedRayContainer &backRays)
1946{
1947        int splits = 0;
1948
1949        while (!rays.empty())
1950        {
1951                BoundedRay *bRay = rays.back();
1952                Ray *ray = bRay->mRay;
1953                rays.pop_back();
1954               
1955                const int cf =
1956                        ray->ClassifyPlane(plane, bRay->mMinT, bRay->mMaxT);
1957               
1958                ray->SetId(cf);
1959
1960                switch (cf)
1961                {
1962                case Ray::COINCIDENT:
1963                        //frontRays.push_back(bRay);
1964                        DEL_PTR(bRay);
1965                        break;
1966                case Ray::BACK:
1967                        backRays.push_back(bRay);
1968                        break;
1969                case Ray::FRONT:
1970                        frontRays.push_back(bRay);
1971                        break;
1972                case Ray::FRONT_BACK:
1973                        {
1974                                // find intersection of ray segment with plane
1975                                const Vector3 extp = ray->Extrap(bRay->mMaxT);
1976                                const float t = plane.FindT(ray->GetLoc(), extp);
1977                               
1978                                const float newT = t * bRay->mMaxT;
1979
1980                                frontRays.push_back(new BoundedRay(ray, bRay->mMinT, newT));
1981                                backRays.push_back(new BoundedRay(ray, newT, bRay->mMaxT));
1982
1983                                DEL_PTR(bRay);
1984
1985                                ++ splits;
1986                        }
1987                        break;
1988                case Ray::BACK_FRONT:
1989                        {
1990                                // find intersection of ray segment with plane
1991                                const Vector3 extp = ray->Extrap(bRay->mMaxT);
1992                                const float t = plane.FindT(ray->GetLoc(), extp);
1993                                const float newT = t * bRay->mMaxT;
1994
1995                                backRays.push_back(new BoundedRay(ray, bRay->mMinT, newT));
1996                                frontRays.push_back(new BoundedRay(ray, newT, bRay->mMaxT));
1997                               
1998                                DEL_PTR(bRay);
1999
2000                                ++ splits;
2001                        }
2002                        break;
2003                default:
2004                        Debug << "Should not come here" << endl;
2005                        break;
2006                }
2007        }
2008
2009        return splits;
2010}
2011
2012void BspTree::ExtractHalfSpaces(BspNode *n, vector<Plane3> &halfSpaces) const
2013{
2014        BspNode *lastNode;
2015        do
2016        {
2017                lastNode = n;
2018
2019                // want to get planes defining geometry of this node => don't take
2020                // split plane of node itself
2021                n = n->GetParent();
2022               
2023                if (n)
2024                {
2025                        BspInterior *interior = dynamic_cast<BspInterior *>(n);
2026                        Plane3 halfSpace = *dynamic_cast<BspInterior *>(interior)->GetPlane();
2027
2028            if (interior->mFront != lastNode)
2029                                halfSpace.ReverseOrientation();
2030
2031                        halfSpaces.push_back(halfSpace);
2032                }
2033        }
2034        while (n);
2035}
2036
2037void BspTree::ConstructGeometry(BspNode *n, BspNodeGeometry &cell) const
2038{
2039        PolygonContainer polys;
2040        ConstructGeometry(n, polys);
2041        cell.mPolys = polys;
2042}
2043
2044void BspTree::ConstructGeometry(BspViewCell *vc, PolygonContainer &cell) const
2045{
2046        vector<BspLeaf *> leaves = vc->mBspLeaves;
2047
2048        vector<BspLeaf *>::const_iterator it, it_end = leaves.end();
2049
2050        for (it = leaves.begin(); it != it_end; ++ it)
2051                ConstructGeometry(*it, cell);
2052}
2053
2054
2055void BspTree::ConstructGeometry(BspNode *n, PolygonContainer &cell) const
2056{
2057        vector<Plane3> halfSpaces;
2058        ExtractHalfSpaces(n, halfSpaces);
2059
2060        PolygonContainer candidatePolys;
2061
2062        // bounded planes are added to the polygons (reverse polygons
2063        // as they have to be outfacing
2064        for (int i = 0; i < (int)halfSpaces.size(); ++ i)
2065        {
2066                Polygon3 *p = GetBoundingBox().CrossSection(halfSpaces[i]);
2067               
2068                if (p->Valid())
2069                {
2070                        candidatePolys.push_back(p->CreateReversePolygon());
2071                        DEL_PTR(p);
2072                }
2073        }
2074
2075        // add faces of bounding box (also could be faces of the cell)
2076        for (int i = 0; i < 6; ++ i)
2077        {
2078                VertexContainer vertices;
2079       
2080                for (int j = 0; j < 4; ++ j)
2081                        vertices.push_back(mBox.GetFace(i).mVertices[j]);
2082
2083                candidatePolys.push_back(new Polygon3(vertices));
2084        }
2085
2086        for (int i = 0; i < (int)candidatePolys.size(); ++ i)
2087        {
2088                // polygon is split by all other planes
2089                for (int j = 0; (j < (int)halfSpaces.size()) && candidatePolys[i]; ++ j)
2090                {
2091                        if (i == j) // polygon and plane are coincident
2092                                continue;
2093
2094                        VertexContainer splitPts;
2095                        Polygon3 *frontPoly, *backPoly;
2096
2097                        const int cf = candidatePolys[i]->ClassifyPlane(halfSpaces[j]);
2098                       
2099                        switch (cf)
2100                        {
2101                                case Polygon3::SPLIT:
2102                                        frontPoly = new Polygon3();
2103                                        backPoly = new Polygon3();
2104
2105                                        candidatePolys[i]->Split(halfSpaces[j], *frontPoly,
2106                                                                                         *backPoly, splitPts);
2107
2108                                        DEL_PTR(candidatePolys[i]);
2109
2110                                        if (frontPoly->Valid())
2111                                                candidatePolys[i] = frontPoly;
2112                                        else
2113                                                DEL_PTR(frontPoly);
2114
2115                                        DEL_PTR(backPoly);
2116                                        break;
2117                                case Polygon3::BACK_SIDE:
2118                                        DEL_PTR(candidatePolys[i]);
2119                                        break;
2120                                // just take polygon as it is
2121                                case Polygon3::FRONT_SIDE:
2122                                case Polygon3::COINCIDENT:
2123                                default:
2124                                        break;
2125                        }
2126                }
2127               
2128                if (candidatePolys[i])
2129                        cell.push_back(candidatePolys[i]);
2130        }
2131}
2132
2133
2134int BspTree::FindNeighbors(BspNode *n, vector<BspLeaf *> &neighbors,
2135                                                   const bool onlyUnmailed) const
2136{
2137        PolygonContainer cell;
2138
2139        ConstructGeometry(n, cell);
2140
2141        stack<BspNode *> nodeStack;
2142        nodeStack.push(mRoot);
2143               
2144        // planes needed to verify that we found neighbor leaf.
2145        vector<Plane3> halfSpaces;
2146        ExtractHalfSpaces(n, halfSpaces);
2147
2148        while (!nodeStack.empty())
2149        {
2150                BspNode *node = nodeStack.top();
2151                nodeStack.pop();
2152
2153                if (node->IsLeaf())
2154                {
2155            if (node != n && (!onlyUnmailed || !node->Mailed()))
2156                        {
2157                                // test all planes of current node if candidate really
2158                                // is neighbour
2159                                PolygonContainer neighborCandidate;
2160                                ConstructGeometry(node, neighborCandidate);
2161                               
2162                                bool isAdjacent = true;
2163                                for (int i = 0; (i < halfSpaces.size()) && isAdjacent; ++ i)
2164                                {
2165                                        const int cf =
2166                                                Polygon3::ClassifyPlane(neighborCandidate, halfSpaces[i]);
2167
2168                                        if (cf == Polygon3::BACK_SIDE)
2169                                                isAdjacent = false;
2170                                }
2171
2172                                if (isAdjacent)
2173                                        neighbors.push_back(dynamic_cast<BspLeaf *>(node));
2174
2175                                CLEAR_CONTAINER(neighborCandidate);
2176                        }
2177                }
2178                else
2179                {
2180                        BspInterior *interior = dynamic_cast<BspInterior *>(node);
2181       
2182                        const int cf = Polygon3::ClassifyPlane(cell, interior->mPlane);
2183
2184                        if (cf == Polygon3::FRONT_SIDE)
2185                                nodeStack.push(interior->mFront);
2186                        else
2187                                if (cf == Polygon3::BACK_SIDE)
2188                                        nodeStack.push(interior->mBack);
2189                                else
2190                                {
2191                                        // random decision
2192                                        nodeStack.push(interior->mBack);
2193                                        nodeStack.push(interior->mFront);
2194                                }
2195                }
2196        }
2197       
2198        CLEAR_CONTAINER(cell);
2199        return (int)neighbors.size();
2200}
2201
2202BspLeaf *BspTree::GetRandomLeaf(const Plane3 &halfspace)
2203{
2204    stack<BspNode *> nodeStack;
2205        nodeStack.push(mRoot);
2206       
2207        int mask = rand();
2208 
2209        while (!nodeStack.empty())
2210        {
2211                BspNode *node = nodeStack.top();
2212                nodeStack.pop();
2213         
2214                if (node->IsLeaf())
2215                {
2216                        return dynamic_cast<BspLeaf *>(node);
2217                }
2218                else
2219                {
2220                        BspInterior *interior = dynamic_cast<BspInterior *>(node);
2221                       
2222                        BspNode *next;
2223       
2224                        PolygonContainer cell;
2225
2226                        // todo: not very efficient: constructs full cell everytime
2227                        ConstructGeometry(interior, cell);
2228
2229                        const int cf = Polygon3::ClassifyPlane(cell, halfspace);
2230
2231                        if (cf == Polygon3::BACK_SIDE)
2232                                next = interior->mFront;
2233                        else
2234                                if (cf == Polygon3::FRONT_SIDE)
2235                                        next = interior->mFront;
2236                        else
2237                        {
2238                                // random decision
2239                                if (mask & 1)
2240                                        next = interior->mBack;
2241                                else
2242                                        next = interior->mFront;
2243                                mask = mask >> 1;
2244                        }
2245
2246                        nodeStack.push(next);
2247                }
2248        }
2249       
2250        return NULL;
2251}
2252
2253BspLeaf *BspTree::GetRandomLeaf(const bool onlyUnmailed)
2254{
2255        stack<BspNode *> nodeStack;
2256       
2257        nodeStack.push(mRoot);
2258
2259        int mask = rand();
2260       
2261        while (!nodeStack.empty())
2262        {
2263                BspNode *node = nodeStack.top();
2264                nodeStack.pop();
2265               
2266                if (node->IsLeaf())
2267                {
2268                        if ( (!onlyUnmailed || !node->Mailed()) )
2269                                return dynamic_cast<BspLeaf *>(node);
2270                }
2271                else
2272                {
2273                        BspInterior *interior = dynamic_cast<BspInterior *>(node);
2274
2275                        // random decision
2276                        if (mask & 1)
2277                                nodeStack.push(interior->mBack);
2278                        else
2279                                nodeStack.push(interior->mFront);
2280
2281                        mask = mask >> 1;
2282                }
2283        }
2284       
2285        return NULL;
2286}
2287
2288int BspTree::ComputePvsSize(const BoundedRayContainer &rays) const
2289{
2290        int pvsSize = 0;
2291
2292        BoundedRayContainer::const_iterator rit, rit_end = rays.end();
2293
2294        Intersectable::NewMail();
2295
2296        for (rit = rays.begin(); rit != rays.end(); ++ rit)
2297        {
2298                Ray *ray = (*rit)->mRay;
2299               
2300                if (!ray->intersections.empty())
2301                {
2302                        if (!ray->intersections[0].mObject->Mailed())
2303                        {
2304                                ray->intersections[0].mObject->Mail();
2305                                ++ pvsSize;
2306                        }
2307                }
2308                if (!ray->sourceObject.mObject->Mailed())
2309                {
2310                        ray->sourceObject.mObject->Mail();
2311                        ++ pvsSize;
2312                }
2313        }
2314
2315        return pvsSize;
2316}
2317
2318/*************************************************************
2319 *            BspNodeGeometry Implementation                 *
2320 *************************************************************/
2321
2322BspNodeGeometry::~BspNodeGeometry()
2323{
2324        CLEAR_CONTAINER(mPolys);
2325}
2326
2327float BspNodeGeometry::GetArea() const
2328{
2329        return Polygon3::GetArea(mPolys);
2330}
2331
2332void BspNodeGeometry::SplitGeometry(BspNodeGeometry &front,
2333                                                                        BspNodeGeometry &back,
2334                                                                        const BspTree &tree,                                             
2335                                                                        const Plane3 &splitPlane) const
2336{       
2337        // get cross section of new polygon
2338        Polygon3 *planePoly = tree.GetBoundingBox().CrossSection(splitPlane);
2339
2340        planePoly = SplitPolygon(planePoly, tree);
2341
2342        //-- plane poly splits all other cell polygons
2343        for (int i = 0; i < (int)mPolys.size(); ++ i)
2344        {
2345                const int cf = mPolys[i]->ClassifyPlane(splitPlane);
2346                       
2347                // split new polygon with all previous planes
2348                switch (cf)
2349                {
2350                        case Polygon3::SPLIT:
2351                                {
2352                                        Polygon3 *poly = new Polygon3(mPolys[i]->mVertices);
2353
2354                                        Polygon3 *frontPoly = new Polygon3();
2355                                        Polygon3 *backPoly = new Polygon3();
2356                               
2357                                        VertexContainer splitPts;
2358                                               
2359                                        poly->Split(splitPlane, *frontPoly, *backPoly, splitPts);
2360
2361                                        DEL_PTR(poly);
2362
2363                                        if (frontPoly->Valid())
2364                                                front.mPolys.push_back(frontPoly);
2365                                        if (backPoly->Valid())
2366                                                back.mPolys.push_back(backPoly);
2367                                }
2368                               
2369                                break;
2370                        case Polygon3::BACK_SIDE:
2371                                back.mPolys.push_back(new Polygon3(mPolys[i]->mVertices));                     
2372                                break;
2373                        case Polygon3::FRONT_SIDE:
2374                                front.mPolys.push_back(new Polygon3(mPolys[i]->mVertices));     
2375                                break;
2376                        case Polygon3::COINCIDENT:
2377                                //front.mPolys.push_back(CreateReversePolygon(mPolys[i]));
2378                                //back.mPolys.push_back(new Polygon3(mPolys[i]->mVertices));
2379                                break;
2380                        default:
2381                                break;
2382                }
2383        }
2384
2385        //-- finally add the new polygon to the child cells
2386        if (planePoly)
2387        {
2388                // add polygon with normal pointing into positive half space to back cell
2389                back.mPolys.push_back(planePoly);
2390                // add polygon with reverse orientation to front cell
2391                front.mPolys.push_back(planePoly->CreateReversePolygon());
2392        }
2393
2394        //Debug << "returning new geometry " << mPolys.size() << " CHILD: " << childCell->mPolys.size() << endl;
2395        //Debug << "old area " << this->GetArea() << " new: " << childCell->GetArea() << endl;
2396}
2397
2398Polygon3 *BspNodeGeometry::SplitPolygon(Polygon3 *planePoly,
2399                                                                                const BspTree &tree) const
2400{
2401        // polygon is split by all other planes
2402        for (int i = 0; (i < (int)mPolys.size()) && planePoly; ++ i)
2403        {
2404                Plane3 plane = mPolys[i]->GetSupportingPlane();
2405
2406                const int cf =
2407                        planePoly->ClassifyPlane(plane);
2408                       
2409                // split new polygon with all previous planes
2410                switch (cf)
2411                {
2412                        case Polygon3::SPLIT:
2413                                {
2414                                        VertexContainer splitPts;
2415                               
2416                                        Polygon3 *frontPoly = new Polygon3();
2417                                        Polygon3 *backPoly = new Polygon3();
2418
2419                                        planePoly->Split(plane, *frontPoly, *backPoly, splitPts);
2420                                        DEL_PTR(planePoly);
2421
2422                                        if (backPoly->Valid())
2423                                                planePoly = backPoly;
2424                                        else
2425                                                DEL_PTR(backPoly);
2426                                }
2427                                break;
2428                        case Polygon3::FRONT_SIDE:
2429                                DEL_PTR(planePoly);
2430                break;
2431                        // polygon is taken as it is
2432                        case Polygon3::BACK_SIDE:
2433                        case Polygon3::COINCIDENT:
2434                        default:
2435                                break;
2436                }
2437        }
2438
2439        return planePoly;
2440}
Note: See TracBrowser for help on using the repository browser.