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

#include "Plane3.h"
#include "VspBspTree.h"
#include "Mesh.h"
#include "common.h"
#include "ViewCell.h"
#include "Environment.h"
#include "Polygon3.h"
#include "Ray.h"
#include "AxisAlignedBox3.h"
#include <stack>
#include <time.h>
#include <iomanip>
#include "Exporter.h"
#include "Plane3.h"
#include "ViewCellBsp.h"

//-- static members
/** Evaluates split plane classification with respect to the plane's 
        contribution for a minimum number of ray splits.
*/
const float VspBspTree::sLeastRaySplitsTable[] = {0, 0, 1, 1, 0};
/** Evaluates split plane classification with respect to the plane's 
        contribution for balanced rays.
*/
const float VspBspTree::sBalancedRaysTable[] = {1, -1, 0, 0, 0};


int VspBspTree::sFrontId = 0; 
int VspBspTree::sBackId = 0;
int VspBspTree::sFrontAndBackId = 0;


/****************************************************************/
/*                  class VspBspTree implementation             */
/****************************************************************/

VspBspTree::VspBspTree(): 
mRoot(NULL),
mPvsUseArea(true)
{
        mRootCell = new BspViewCell();

        Randomize(); // initialise random generator for heuristics

        //-- termination criteria for autopartition
        environment->GetIntValue("VspBspTree.Termination.maxDepth", mTermMaxDepth);
        environment->GetIntValue("VspBspTree.Termination.minPvs", mTermMinPvs);
        environment->GetIntValue("VspBspTree.Termination.minRays", mTermMinRays);
        environment->GetFloatValue("VspBspTree.Termination.minArea", mTermMinArea);     
        environment->GetFloatValue("VspBspTree.Termination.maxRayContribution", mTermMaxRayContribution);
        environment->GetFloatValue("VspBspTree.Termination.minAccRayLenght", mTermMinAccRayLength);

        //-- factors for bsp tree split plane heuristics
        environment->GetFloatValue("VspBspTree.Factor.balancedRays", mBalancedRaysFactor);
        environment->GetFloatValue("VspBspTree.Factor.pvs", mPvsFactor);
        environment->GetFloatValue("VspBspTree.Termination.ct_div_ci", mCtDivCi);

        //-- termination criteria for axis aligned split
        environment->GetFloatValue("VspBspTree.Termination.AxisAligned.ct_div_ci", mAaCtDivCi);
        environment->GetFloatValue("VspBspTree.Termination.AxisAligned.maxCostRatio", mMaxCostRatio);
        
        environment->GetIntValue("VspBspTree.Termination.AxisAligned.minRays", 
                                                         mTermMinRaysForAxisAligned);
        
        //-- partition criteria
        environment->GetIntValue("VspBspTree.maxPolyCandidates", mMaxPolyCandidates);
        environment->GetIntValue("VspBspTree.maxRayCandidates", mMaxRayCandidates);
        environment->GetIntValue("VspBspTree.splitPlaneStrategy", mSplitPlaneStrategy);
        environment->GetFloatValue("VspBspTree.AxisAligned.splitBorder", mSplitBorder);

        environment->GetFloatValue("VspBspTree.Construction.epsilon", mEpsilon);
        environment->GetIntValue("VspBspTree.maxTests", mMaxTests);

    Debug << "BSP max depth: " << mTermMaxDepth << endl;
        Debug << "BSP min PVS: " << mTermMinPvs << endl;
        Debug << "BSP min area: " << mTermMinArea << endl;
        Debug << "BSP min rays: " << mTermMinRays << endl;
        Debug << "BSP max polygon candidates: " << mMaxPolyCandidates << endl;
        Debug << "BSP max plane candidates: " << mMaxRayCandidates << endl;

        Debug << "Split plane strategy: ";
        if (mSplitPlaneStrategy & RANDOM_POLYGON)
                Debug << "random polygon ";
        if (mSplitPlaneStrategy & AXIS_ALIGNED)
                Debug << "axis aligned ";
        if (mSplitPlaneStrategy & LEAST_RAY_SPLITS)
                Debug << "least ray splits ";
        if (mSplitPlaneStrategy & BALANCED_RAYS)
                Debug << "balanced rays ";
        if (mSplitPlaneStrategy & PVS)
                Debug << "pvs";

        Debug << endl;
}


const BspTreeStatistics &VspBspTree::GetStatistics() const 
{
        return mStat;
}


VspBspTree::~VspBspTree()
{
        DEL_PTR(mRoot);
        DEL_PTR(mRootCell);
}

int VspBspTree::AddMeshToPolygons(Mesh *mesh, 
                                                                  PolygonContainer &polys, 
                                                                  MeshInstance *parent)
{
        FaceContainer::const_iterator fi;
        
        // copy the face data to polygons
        for (fi = mesh->mFaces.begin(); fi != mesh->mFaces.end(); ++ fi)
        {
                Polygon3 *poly = new Polygon3((*fi), mesh);
                
                if (poly->Valid(mEpsilon))
                {
                        poly->mParent = parent; // set parent intersectable
                        polys.push_back(poly);
                }
                else
                        DEL_PTR(poly);
        }
        return (int)mesh->mFaces.size();
}

int VspBspTree::AddToPolygonSoup(const ViewCellContainer &viewCells, 
                                                          PolygonContainer &polys, 
                                                          int maxObjects)
{
        int limit = (maxObjects > 0) ? 
                Min((int)viewCells.size(), maxObjects) : (int)viewCells.size();
  
        int polysSize = 0;

        for (int i = 0; i < limit; ++ i)
        {
                if (viewCells[i]->GetMesh()) // copy the mesh data to polygons
                {
                        mBox.Include(viewCells[i]->GetBox()); // add to BSP tree aabb
                        polysSize += 
                                AddMeshToPolygons(viewCells[i]->GetMesh(), 
                                                                  polys, 
                                                                  viewCells[i]);
                }
        }

        return polysSize;
}

int VspBspTree::AddToPolygonSoup(const ObjectContainer &objects, 
                                                                 PolygonContainer &polys, 
                                                                 int maxObjects)
{
        int limit = (maxObjects > 0) ? 
                Min((int)objects.size(), maxObjects) : (int)objects.size();
  
        for (int i = 0; i < limit; ++i)
        {
                Intersectable *object = objects[i];//*it;
                Mesh *mesh = NULL;

                switch (object->Type()) // extract the meshes
                {
                case Intersectable::MESH_INSTANCE:
                        mesh = dynamic_cast<MeshInstance *>(object)->GetMesh();
                        break;
                case Intersectable::VIEW_CELL:
                        mesh = dynamic_cast<ViewCell *>(object)->GetMesh();
                        break;
                        // TODO: handle transformed mesh instances
                default:
                        Debug << "intersectable type not supported" << endl;
                        break;
                }
                
        if (mesh) // copy the mesh data to polygons
                {
                        mBox.Include(object->GetBox()); // add to BSP tree aabb
                        AddMeshToPolygons(mesh, polys, mRootCell);
                }
        }

        return (int)polys.size();
}

void VspBspTree::Construct(const VssRayContainer &sampleRays)
{
    mStat.nodes = 1;
        mBox.Initialize();      // initialise BSP tree bounding box
        
        PolygonContainer polys;
        RayInfoContainer *rays = new RayInfoContainer();

        VssRayContainer::const_iterator rit, rit_end = sampleRays.end();

        long startTime = GetTime();

        Debug << "**** Extracting polygons from rays ****\n";

        Intersectable::NewMail();

        //-- extract polygons intersected by the rays
        for (rit = sampleRays.begin(); rit != rit_end; ++ rit)
        {
                VssRay *ray = *rit;
        
                if (ray->mTerminationObject && !ray->mTerminationObject->Mailed())
                {
                        ray->mTerminationObject->Mail();
                        MeshInstance *obj = dynamic_cast<MeshInstance *>(ray->mTerminationObject);
                        AddMeshToPolygons(obj->GetMesh(), polys, obj);
                }

                if (ray->mOriginObject && !ray->mOriginObject->Mailed())
                {
                        ray->mOriginObject->Mail();
                        MeshInstance *obj = dynamic_cast<MeshInstance *>(ray->mOriginObject);
                        AddMeshToPolygons(obj->GetMesh(), polys, obj);
                }
        }

        // compute bounding box
        Polygon3::IncludeInBox(polys, mBox);

        //-- store rays
        for (rit = sampleRays.begin(); rit != rit_end; ++ rit)
        {
                VssRay *ray = *rit;
                
                float minT, maxT;

                // TODO: not very efficient to implictly cast between rays types ...
                if (mBox.GetRaySegment(*ray, minT, maxT))
                {
                        float len = ray->Length();
                        
                        if (!len) 
                                len = Limits::Small;
                        
                        rays->push_back(RayInfo(ray, minT / len, maxT / len));
                }
        }

        mStat.polys = (int)polys.size();

        Debug << "**** Finished polygon extraction ****" << endl;
        Debug << (int)polys.size() << " polys extracted from " << (int)sampleRays.size() << " rays" << endl;
        Debug << "extraction time: " << TimeDiff(startTime, GetTime())*1e-3 << "s" << endl;

        Construct(polys, rays);

        // clean up polygons
        CLEAR_CONTAINER(polys);
}

void VspBspTree::Construct(const PolygonContainer &polys, RayInfoContainer *rays)
{
        std::stack<VspBspTraversalData> tStack;

        mRoot = new BspLeaf();

        // constrruct root node geometry
        BspNodeGeometry *geom = new BspNodeGeometry();
        ConstructGeometry(mRoot, *geom);

        VspBspTraversalData tData(mRoot, 
                                                          new PolygonContainer(polys), 
                                                          0,
                                                          rays, 
                              ComputePvsSize(*rays), 
                                                          geom->GetArea(), 
                                                          geom);

        tStack.push(tData);

        mStat.Start();
        cout << "Contructing vsp bsp tree ... ";

        long startTime = GetTime();
        
        while (!tStack.empty()) 
        {
                tData = tStack.top();

            tStack.pop();

                // subdivide leaf node
                BspNode *r = Subdivide(tStack, tData);

                if (r == mRoot)
                        Debug << "BSP tree construction time spent at root: " 
                                  << TimeDiff(startTime, GetTime())*1e-3 << "s" << endl;
        }

        cout << "finished\n";

        mStat.Stop();
}

bool VspBspTree::TerminationCriteriaMet(const VspBspTraversalData &data) const
{
        return 
                (((int)data.mRays->size() <= mTermMinRays) ||
                 (data.mPvs <= mTermMinPvs) ||
                 (data.mArea <= mTermMinArea) ||
                // (data.GetAvgRayContribution() >= mTermMaxRayContribution) ||
                 (data.mDepth >= mTermMaxDepth));
}

BspNode *VspBspTree::Subdivide(VspBspTraversalStack &tStack, 
                                                           VspBspTraversalData &tData)
{
        //-- terminate traversal  
        if (TerminationCriteriaMet(tData))              
        {
                BspLeaf *leaf = dynamic_cast<BspLeaf *>(tData.mNode);
        
                BspViewCell *viewCell = new BspViewCell();
                
                leaf->SetViewCell(viewCell);
                viewCell->mBspLeaves.push_back(leaf);

                //-- add pvs
                if (viewCell != mRootCell)
                {
                        int conSamp = 0, sampCon = 0;
                        AddToPvs(leaf, *tData.mRays, conSamp, sampCon);
                        
                        mStat.contributingSamples += conSamp;
                        mStat.sampleContributions += sampCon;
                }

                EvaluateLeafStats(tData);
        
                //-- clean up
                
                DEL_PTR(tData.mPolygons);
                DEL_PTR(tData.mRays);
                DEL_PTR(tData.mGeometry);

                return leaf;
        }

        //-- continue subdivision
        PolygonContainer coincident;
        
        VspBspTraversalData tFrontData(new PolygonContainer(), 
                                                                   tData.mDepth + 1,
                                                                   new RayInfoContainer(), 
                                                                   new BspNodeGeometry());

        VspBspTraversalData tBackData(new PolygonContainer(), 
                                                                  tData.mDepth + 1,
                                                                  new RayInfoContainer(), 
                                                                  new BspNodeGeometry());

        // create new interior node and two leaf nodes
        BspInterior *interior = 
                SubdivideNode(tData, tFrontData, tBackData, coincident);

#ifdef _DEBUG   
//      if (frontPolys->empty() && backPolys->empty() && (coincident.size() > 2))
//      {       for (PolygonContainer::iterator it = coincident.begin(); it != coincident.end(); ++it)
//                      Debug << (*it) << " " << (*it)->GetArea() << " " << (*it)->mParent << endl ;
//              Debug << endl;}
#endif

        // push the children on the stack
        tStack.push(tFrontData);
        tStack.push(tBackData);

        // cleanup
        DEL_PTR(tData.mNode);   

        DEL_PTR(tData.mPolygons);
        DEL_PTR(tData.mRays);
        DEL_PTR(tData.mGeometry);

        return interior;
}

BspInterior *VspBspTree::SubdivideNode(VspBspTraversalData &tData,
                                                                           VspBspTraversalData &frontData,
                                                                           VspBspTraversalData &backData,
                                                                           PolygonContainer &coincident)
{
        mStat.nodes += 2;
        
        BspLeaf *leaf = dynamic_cast<BspLeaf *>(tData.mNode);
        
        // select subdivision plane
        BspInterior *interior = new BspInterior(SelectPlane(leaf, tData)); 

#ifdef _DEBUG
        Debug << interior << endl;
#endif
        
        // subdivide rays into front and back rays
        SplitRays(interior->GetPlane(), 
                          *tData.mRays, 
                          *frontData.mRays, 
                          *backData.mRays);
        
        // subdivide polygons with plane
        mStat.splits += SplitPolygons(interior->GetPlane(),
                                                                  *tData.mPolygons, 
                                          *frontData.mPolygons, 
                                                                  *backData.mPolygons,
                                                                  coincident);

    // compute pvs
        frontData.mPvs = ComputePvsSize(*frontData.mRays);
        backData.mPvs = ComputePvsSize(*backData.mRays);

        // split geometry and compute area
        if (1)
        {
                tData.mGeometry->SplitGeometry(*frontData.mGeometry, 
                                                                           *backData.mGeometry, 
                                                                           interior->GetPlane(),
                                                                           mBox,
                                                                           mEpsilon);
        
                
                frontData.mArea = frontData.mGeometry->GetArea();
                backData.mArea = backData.mGeometry->GetArea();
        }

        // compute accumulated ray length
        //frontData.mAccRayLength = AccumulatedRayLength(*frontData.mRays);
        //backData.mAccRayLength = AccumulatedRayLength(*backData.mRays);

        //-- create front and back leaf

        BspInterior *parent = leaf->GetParent();

        // replace a link from node's parent
        if (!leaf->IsRoot())
        {
                parent->ReplaceChildLink(leaf, interior);
                interior->SetParent(parent);
        }
        else // new root
        {
                mRoot = interior;
        }

        // and setup child links
        interior->SetupChildLinks(new BspLeaf(interior), new BspLeaf(interior));
        
        frontData.mNode = interior->GetFront();
        backData.mNode = interior->GetBack();
        
        //DEL_PTR(leaf);
        return interior;
}

void VspBspTree::AddToPvs(BspLeaf *leaf,
                                                  const RayInfoContainer &rays, 
                                                  int &sampleContributions,
                                                  int &contributingSamples)
{
        sampleContributions = 0;
        contributingSamples = 0;

    RayInfoContainer::const_iterator it, it_end = rays.end();

        ViewCell *vc = leaf->GetViewCell();

        // add contributions from samples to the PVS
        for (it = rays.begin(); it != it_end; ++ it)
        {
                int contribution = 0;
                VssRay *ray = (*it).mRay;
                        
                if (ray->mTerminationObject)
                        contribution += vc->GetPvs().AddSample(ray->mTerminationObject);
                
                if (ray->mOriginObject)
                        contribution += vc->GetPvs().AddSample(ray->mOriginObject);

                if (contribution)
                {
                        sampleContributions += contribution;
                        ++ contributingSamples;
                }
                //leaf->mVssRays.push_back(ray);
        }
}

void VspBspTree::SortSplitCandidates(const PolygonContainer &polys, 
                                                                         const int axis, 
                                                                         vector<SortableEntry> &splitCandidates) const
{
  splitCandidates.clear();
  
  int requestedSize = 2 * (int)polys.size();
  // creates a sorted split candidates array  
  splitCandidates.reserve(requestedSize);
  
  PolygonContainer::const_iterator it, it_end = polys.end();

  AxisAlignedBox3 box;

  // insert all queries 
  for(it = polys.begin(); it != it_end; ++ it)
  {
          box.Initialize();
          box.Include(*(*it));
          
          splitCandidates.push_back(SortableEntry(SortableEntry::POLY_MIN, box.Min(axis), *it));
      splitCandidates.push_back(SortableEntry(SortableEntry::POLY_MAX, box.Max(axis), *it));
  }
  
  stable_sort(splitCandidates.begin(), splitCandidates.end());
}


float VspBspTree::BestCostRatio(const PolygonContainer &polys,
                                                                const AxisAlignedBox3 &box,
                                                                const int axis,
                                                                float &position,
                                                                int &objectsBack,
                                int &objectsFront) const
{
        vector<SortableEntry> splitCandidates;

        SortSplitCandidates(polys, axis, splitCandidates);
        
        // go through the lists, count the number of objects left and right
        // and evaluate the following cost funcion:
        // C = ct_div_ci  + (ol + or)/queries
        
        int objectsLeft = 0, objectsRight = (int)polys.size();
        
        float minBox = box.Min(axis);
        float maxBox = box.Max(axis);
        float boxArea = box.SurfaceArea();
  
        float minBand = minBox + mSplitBorder * (maxBox - minBox);
        float maxBand = minBox + (1.0f - mSplitBorder) * (maxBox - minBox);
        
        float minSum = 1e20f;
        vector<SortableEntry>::const_iterator ci, ci_end = splitCandidates.end();

        for(ci = splitCandidates.begin(); ci != ci_end; ++ ci) 
        {
                switch ((*ci).type) 
                {
                        case SortableEntry::POLY_MIN:
                                ++ objectsLeft;
                                break;
                        case SortableEntry::POLY_MAX:
                            -- objectsRight;
                                break;
                        default:
                                break;
                }
                
                if ((*ci).value > minBand && (*ci).value < maxBand) 
                {
                        AxisAlignedBox3 lbox = box;
                        AxisAlignedBox3 rbox = box;
                        lbox.SetMax(axis, (*ci).value);
                        rbox.SetMin(axis, (*ci).value);

                        float sum = objectsLeft * lbox.SurfaceArea() + 
                                                objectsRight * rbox.SurfaceArea();
      
                        if (sum < minSum) 
                        {
                                minSum = sum;
                                position = (*ci).value;

                                objectsBack = objectsLeft;
                                objectsFront = objectsRight;
                        }
                }
        }
  
        float oldCost = (float)polys.size();
        float newCost = mAaCtDivCi + minSum / boxArea;
        float ratio = newCost / oldCost;


#if 0
  Debug << "====================" << endl;
  Debug << "costRatio=" << ratio << " pos=" << position<<" t=" << (position - minBox)/(maxBox - minBox)
      << "\t o=(" << objectsBack << "," << objectsFront << ")" << endl;
#endif
  return ratio;
}

bool VspBspTree::SelectAxisAlignedPlane(Plane3 &plane,
                                        const PolygonContainer &polys) const
{
        AxisAlignedBox3 box;
        box.Initialize();
        
        // create bounding box of region
        Polygon3::IncludeInBox(polys, box);
        
        int objectsBack = 0, objectsFront = 0;
        int axis = 0;
        float costRatio = MAX_FLOAT;
        Vector3 position;

        //-- area subdivision
        for (int i = 0; i < 3; ++ i) 
        {
                float p = 0;
                float r = BestCostRatio(polys, box, i, p, objectsBack, objectsFront);
                
                if (r < costRatio)
                {
                        costRatio = r;
                        axis = i;
                        position = p;
                }
        }
        
        if (costRatio >= mMaxCostRatio)
                return false;

        Vector3 norm(0,0,0); norm[axis] = 1.0f;
        plane = Plane3(norm, position);

        return true;
}

Plane3 VspBspTree::SelectPlane(BspLeaf *leaf, VspBspTraversalData &data)
{
        if ((mSplitPlaneStrategy & AXIS_ALIGNED) &&
                ((int)data.mRays->size() > mTermMinRaysForAxisAligned))
        {
                Plane3 plane;
                if (SelectAxisAlignedPlane(plane, *data.mPolygons)) // TODO: should be rays
                        return plane;
        }

        // simplest strategy: just take next polygon
        if (mSplitPlaneStrategy & RANDOM_POLYGON)
        {
        if (!data.mPolygons->empty())
                {
                        const int randIdx = Random((int)data.mPolygons->size());
                        Polygon3 *nextPoly = (*data.mPolygons)[randIdx];
                        return nextPoly->GetSupportingPlane();
                }
                else
                {
                        //-- choose plane on midpoint of a ray
                        const int candidateIdx = Random((int)data.mRays->size());
                                                                        
                        const Vector3 minPt = (*data.mRays)[candidateIdx].ExtrapOrigin();
                        const Vector3 maxPt = (*data.mRays)[candidateIdx].ExtrapTermination();

                        const Vector3 pt = (maxPt + minPt) * 0.5;

                        const Vector3 normal = (*data.mRays)[candidateIdx].mRay->GetDir();
                        
                        return Plane3(normal, pt);
                }

                return Plane3();
        }

        // use heuristics to find appropriate plane
        return SelectPlaneHeuristics(leaf, data); 
}

Plane3 VspBspTree::SelectPlaneHeuristics(BspLeaf *leaf, VspBspTraversalData &data)
{
        float lowestCost = MAX_FLOAT;
        Plane3 bestPlane;
        Plane3 plane;

        int limit = Min((int)data.mPolygons->size(), mMaxPolyCandidates);
        
        int candidateIdx = limit;
        
        for (int i = 0; i < limit; ++ i)
        {
                candidateIdx = GetNextCandidateIdx(candidateIdx, *data.mPolygons);
                
                Polygon3 *poly = (*data.mPolygons)[candidateIdx];

                // evaluate current candidate
                const float candidateCost = 
                        SplitPlaneCost(poly->GetSupportingPlane(), data);

                if (candidateCost < lowestCost)
                {
                        bestPlane = poly->GetSupportingPlane();
                        lowestCost = candidateCost;
                }
        }
        
        //Debug << "lowest: " << lowestCost << endl;

        //-- choose candidate planes extracted from rays
        // we currently use different two methods chosen with
        // equal probability
        
        // take 3 ray endpoints, where two are minimum and one a maximum
        // point or the other way round
        for (int i = 0; i < mMaxRayCandidates / 2; ++ i)
        {
                candidateIdx = Random((int)data.mRays->size());
        
                RayInfo rayInf = (*data.mRays)[candidateIdx];

                const Vector3 minPt = rayInf.ExtrapOrigin();
                const Vector3 maxPt = rayInf.ExtrapTermination();

                const Vector3 pt = (maxPt + minPt) * 0.5;
                const Vector3 normal = Normalize(rayInf.mRay->GetDir());

                plane = Plane3(normal, pt);

                const float candidateCost = SplitPlaneCost(plane, data);

                if (candidateCost < lowestCost)
                {
                        bestPlane = plane;      
                        lowestCost = candidateCost;
                }
        }

        // take plane normal as plane normal and the midpoint of the ray.
        // PROBLEM: does not resemble any point where visibility is likely to change
        //Debug << "lowest2: " << lowestCost << endl;
        for (int i = 0; i < mMaxRayCandidates / 2; ++ i)
        {
                Vector3 pt[3];
                int idx[3];
                int cmaxT = 0;
                int cminT = 0;
                bool chooseMin = false;

                for (int j = 0; j < 3; j ++)
                {
                        idx[j] = Random((int)data.mRays->size() * 2);
                                
                        if (idx[j] >= (int)data.mRays->size())
                        {
                                idx[j] -= (int)data.mRays->size();
                                
                                chooseMin = (cminT < 2);
                        }
                        else
                                chooseMin = (cmaxT < 2);

                        RayInfo rayInf = (*data.mRays)[idx[j]];
                        pt[j] = chooseMin ? rayInf.ExtrapOrigin() : rayInf.ExtrapTermination();
                }       
                        
                plane = Plane3(pt[0], pt[1], pt[2]);

                const float candidateCost = SplitPlaneCost(plane, data);

                if (candidateCost < lowestCost)
                {
                        //Debug << "choose ray plane 2: " << candidateCost << endl;
                        bestPlane = plane;
                        
                        lowestCost = candidateCost;
                }
        }       

#ifdef _DEBUG
        Debug << "plane lowest cost: " << lowestCost << endl;
#endif
        return bestPlane;
}

int VspBspTree::GetNextCandidateIdx(int currentIdx, PolygonContainer &polys)
{
        const int candidateIdx = Random(currentIdx --);

        // swap candidates to avoid testing same plane 2 times
        std::swap(polys[currentIdx], polys[candidateIdx]);
        
        return currentIdx;
        //return Random((int)polys.size());
}


inline void VspBspTree::GenerateUniqueIdsForPvs()
{
        Intersectable::NewMail(); sBackId = ViewCell::sMailId;
        Intersectable::NewMail(); sFrontId = ViewCell::sMailId;
        Intersectable::NewMail(); sFrontAndBackId = ViewCell::sMailId;
}

float VspBspTree::SplitPlaneCost(const Plane3 &candidatePlane, 
                                                                 const VspBspTraversalData &data)
{
        float val = 0;

        float sumBalancedRays = 0;
        float sumRaySplits = 0;
        
        int frontPvs = 0;
        int backPvs = 0;

        // probability that view point lies in child
        float pOverall = 0;
        float pFront = 0;
        float pBack = 0;

        const bool pvsUseLen = false;

        if (mSplitPlaneStrategy & PVS)
        {
                // create unique ids for pvs heuristics
                GenerateUniqueIdsForPvs();
        
                if (mPvsUseArea) // use front and back cell areas to approximate volume
                {       
                        // construct child geometry with regard to the candidate split plane
                        BspNodeGeometry frontCell;
                        BspNodeGeometry backCell;
                
                        data.mGeometry->SplitGeometry(frontCell, 
                                                                                  backCell, 
                                                                                  candidatePlane,
                                                                                  mBox,
                                                                                  mEpsilon);
                
                        pFront = frontCell.GetArea();
                        pBack = backCell.GetArea();

                        pOverall = data.mArea;
                }
        }
                
        int limit;
        bool useRand;

        // choose test polyongs randomly if over threshold
        if ((int)data.mRays->size() > mMaxTests)
        {
                useRand = true;
                limit = mMaxTests;
        }
        else
        {
                useRand = false;
                limit = (int)data.mRays->size();
        }

        for (int i = 0; i < limit; ++ i)
        {
                const int testIdx = useRand ? Random(limit) : i;
                RayInfo rayInf = (*data.mRays)[testIdx];

                VssRay *ray = rayInf.mRay;
                const int cf = rayInf.ComputeRayIntersection(candidatePlane, ray->mT);

                if (mSplitPlaneStrategy & LEAST_RAY_SPLITS)
                {
                        sumBalancedRays += cf;
                }
                
                if (mSplitPlaneStrategy & BALANCED_RAYS)
                {
                        if (cf == 0)
                                ++ sumRaySplits;
                }

                if (mSplitPlaneStrategy & PVS)
                {
                        // in case the ray intersects an object
                        // assure that we only count the object 
                        // once for the front and once for the back side of the plane
                        
                        // add the termination object
                        AddObjToPvs(ray->mTerminationObject, cf, frontPvs, backPvs);
                        
                        // add the source object
                        AddObjToPvs(ray->mOriginObject, cf, frontPvs, backPvs);
                        
                        // use number or length of rays to approximate volume
                        if (!mPvsUseArea)
                        {
                                float len = 1;

                                if (pvsUseLen) // use length of rays
                                        len = rayInf.SqrSegmentLength();
                        
                                pOverall += len;

                                if (cf == 1)
                                        pFront += len;
                                if (cf == -1)
                                        pBack += len;
                                if (cf == 0)
                                {
                                        // use length of rays to approximate volume
                                        if (pvsUseLen) 
                                        {
                                                float newLen = len * 
                                                        (rayInf.GetMaxT() - rayInf.mRay->mT) / 
                                                        (rayInf.GetMaxT() - rayInf.GetMinT());
                
                                                if (candidatePlane.Side(rayInf.ExtrapOrigin()) <= 0)
                                                {
                                                        pBack += newLen;
                                                        pFront += len - newLen;
                                                }
                                                else
                                                {
                                                        pFront += newLen;
                                                        pBack += len - newLen;
                                                }
                                        }
                                        else
                                        {
                                                ++ pFront;
                                                ++ pBack;
                                        }
                                }
                        }
                }
        }

        const float raysSize = (float)data.mRays->size() + Limits::Small;

        if (mSplitPlaneStrategy & LEAST_RAY_SPLITS)
                val += mLeastRaySplitsFactor * sumRaySplits / raysSize;

        if (mSplitPlaneStrategy & BALANCED_RAYS)
                val += mBalancedRaysFactor * fabs(sumBalancedRays) /  raysSize;

        const float denom = pOverall * (float)data.mPvs * 2.0f + Limits::Small;

        if (mSplitPlaneStrategy & PVS)
        {
                val += mPvsFactor * (frontPvs * pFront + (backPvs * pBack)) / denom;

                // give penalty to unbalanced split
                if (0)
                if (((pFront * 0.2 + Limits::Small) > pBack) || 
                        (pFront < (pBack * 0.2 + Limits::Small)))
                        val += 0.5;
        }

#ifdef _DEBUG
//      Debug << "totalpvs: " << pvs << " ptotal: " << pOverall
//                << " frontpvs: " << frontPvs << " pFront: " << pFront 
//                << " backpvs: " << backPvs << " pBack: " << pBack << endl << endl;
#endif
        return val;
}

void VspBspTree::AddObjToPvs(Intersectable *obj,
                                                         const int cf,
                                                         int &frontPvs,
                                                         int &backPvs) const
{
        if (!obj)
                return;
        // TODO: does this really belong to no pvs?
        //if (cf == Ray::COINCIDENT) return;

        // object belongs to both PVS
        if (cf >= 0)
        {
                if ((obj->mMailbox != sFrontId) && 
                        (obj->mMailbox != sFrontAndBackId))
                {
                        ++ frontPvs;

                        if (obj->mMailbox == sBackId)
                                obj->mMailbox = sFrontAndBackId;        
                        else
                                obj->mMailbox = sFrontId;                                                               
                }
        }
        
        if (cf <= 0)
        {
                if ((obj->mMailbox != sBackId) &&
                        (obj->mMailbox != sFrontAndBackId))
                {
                        ++ backPvs;

                        if (obj->mMailbox == sFrontId)
                                obj->mMailbox = sFrontAndBackId; 
                        else
                                obj->mMailbox = sBackId;                                
                }
        }
}

void VspBspTree::CollectLeaves(vector<BspLeaf *> &leaves) const
{
        stack<BspNode *> nodeStack;
        nodeStack.push(mRoot);
 
        while (!nodeStack.empty()) 
        {
                BspNode *node = nodeStack.top();
    
                nodeStack.pop();
    
                if (node->IsLeaf()) 
                {
                        BspLeaf *leaf = (BspLeaf *)node;                
                        leaves.push_back(leaf);
                } 
                else 
                {
                        BspInterior *interior = dynamic_cast<BspInterior *>(node);

                        nodeStack.push(interior->GetBack());
                        nodeStack.push(interior->GetFront());
                }
        }
}

AxisAlignedBox3 VspBspTree::GetBoundingBox() const
{
        return mBox;
}

BspNode *VspBspTree::GetRoot() const
{
        return mRoot;
}

void VspBspTree::EvaluateLeafStats(const VspBspTraversalData &data)
{
        // the node became a leaf -> evaluate stats for leafs
        BspLeaf *leaf = dynamic_cast<BspLeaf *>(data.mNode);

        // store maximal and minimal depth
        if (data.mDepth > mStat.maxDepth)
                mStat.maxDepth = data.mDepth;

        if (data.mDepth < mStat.minDepth)
                mStat.minDepth = data.mDepth;

        // accumulate depth to compute average depth
        mStat.accumDepth += data.mDepth;
        

        if (data.mDepth >= mTermMaxDepth)
                ++ mStat.maxDepthNodes;

        if (data.mPvs < mTermMinPvs)
                ++ mStat.minPvsNodes;

        if ((int)data.mRays->size() < mTermMinRays)
                ++ mStat.minRaysNodes;

        if (data.GetAvgRayContribution() > mTermMaxRayContribution)
                ++ mStat.maxRayContribNodes;
        
        if (data.mGeometry->GetArea() <= mTermMinArea) 
                ++ mStat.minAreaNodes;

#ifdef _DEBUG
        Debug << "BSP stats: "
                  << "Depth: " << data.mDepth << " (max: " << mTermMaxDepth << "), "
                  << "PVS: " << data.mPvs << " (min: " << mTermMinPvs << "), "
                  << "Area: " << data.mArea << " (min: " << mTermMinArea << "), "
                  << "#rays: " << (int)data.mRays->size() << " (max: " << mTermMinRays << "), " 
                  << "#pvs: " << leaf->GetViewCell()->GetPvs().GetSize() << "=, "
                  << "#avg ray contrib (pvs): " << (float)data.mPvs / (float)data.mRays->size() << endl;
#endif
}

int VspBspTree::CastRay(Ray &ray)
{
        int hits = 0;
  
        stack<BspRayTraversalData> tStack;
  
        float maxt, mint;

        if (!mBox.GetRaySegment(ray, mint, maxt))
                return 0;

        Intersectable::NewMail();

        Vector3 entp = ray.Extrap(mint);
        Vector3 extp = ray.Extrap(maxt);
  
        BspNode *node = mRoot;
        BspNode *farChild = NULL;
        
        while (1)
        {
                if (!node->IsLeaf()) 
                {
                        BspInterior *in = dynamic_cast<BspInterior *>(node);

                        Plane3 splitPlane = in->GetPlane();
                        const int entSide = splitPlane.Side(entp);
                        const int extSide = splitPlane.Side(extp);

                        if (entSide < 0)
                        {
                                node = in->GetBack();

                                if(extSide <= 0) // plane does not split ray => no far child
                                        continue;
                                        
                                farChild = in->GetFront(); // plane splits ray

                        } else if (entSide > 0)
                        {
                                node = in->GetFront();

                                if (extSide >= 0) // plane does not split ray => no far child
                                        continue;

                                farChild = in->GetBack(); // plane splits ray                   
                        }
                        else // ray and plane are coincident
                        {
                                // WHAT TO DO IN THIS CASE ?
                                //break;
                                node = in->GetFront();
                                continue;
                        }

                        // push data for far child
                        tStack.push(BspRayTraversalData(farChild, extp, maxt));

                        // find intersection of ray segment with plane
                        float t;
                        extp = splitPlane.FindIntersection(ray.GetLoc(), extp, &t);
                        maxt *= t;
                        
                } else // reached leaf => intersection with view cell
                {
                        BspLeaf *leaf = dynamic_cast<BspLeaf *>(node);
      
                        if (!leaf->GetViewCell()->Mailed())
                        {
                                //ray.bspIntersections.push_back(Ray::VspBspIntersection(maxt, leaf));
                                leaf->GetViewCell()->Mail();
                                ++ hits;
                        }
                        
                        //-- fetch the next far child from the stack
                        if (tStack.empty())
                                break;
      
                        entp = extp;
                        mint = maxt; // NOTE: need this?

                        if (ray.GetType() == Ray::LINE_SEGMENT && mint > 1.0f)
                                break;

                        BspRayTraversalData &s = tStack.top();

                        node = s.mNode;
                        extp = s.mExitPoint;
                        maxt = s.mMaxT;

                        tStack.pop();
                }
        }

        return hits;
}

bool VspBspTree::Export(const string filename)
{
        Exporter *exporter = Exporter::GetExporter(filename);

        if (exporter) 
        {
                //exporter->ExportVspBspTree(*this);
                return true;
        }       

        return false;
}

void VspBspTree::CollectViewCells(ViewCellContainer &viewCells) const
{
        stack<BspNode *> nodeStack;
        nodeStack.push(mRoot);

        ViewCell::NewMail();

        while (!nodeStack.empty()) 
        {
                BspNode *node = nodeStack.top();
                nodeStack.pop();

                if (node->IsLeaf()) 
                {
                        ViewCell *viewCell = dynamic_cast<BspLeaf *>(node)->GetViewCell();

                        if (!viewCell->Mailed()) 
                        {
                                viewCell->Mail();
                                viewCells.push_back(viewCell);
                        }
                }
                else 
                {
                        BspInterior *interior = dynamic_cast<BspInterior *>(node);

                        nodeStack.push(interior->GetFront());
                        nodeStack.push(interior->GetBack());
                }
        }
}

void VspBspTree::EvaluateViewCellsStats(BspViewCellsStatistics &stat) const
{
        stat.Reset();

        stack<BspNode *> nodeStack;
        nodeStack.push(mRoot);

        ViewCell::NewMail();

        // exclude root cell
        mRootCell->Mail();

        while (!nodeStack.empty()) 
        {
                BspNode *node = nodeStack.top();
                nodeStack.pop();

                if (node->IsLeaf()) 
                {
                        ++ stat.bspLeaves;

                        BspViewCell *viewCell = dynamic_cast<BspLeaf *>(node)->GetViewCell();

                        if (!viewCell->Mailed()) 
                        {
                                viewCell->Mail();
                                
                                ++ stat.viewCells;
                                const int pvsSize = viewCell->GetPvs().GetSize();

                stat.pvs += pvsSize;

                                if (pvsSize < 1)
                                        ++ stat.emptyPvs;

                                if (pvsSize > stat.maxPvs)
                                        stat.maxPvs = pvsSize;

                                if (pvsSize < stat.minPvs)
                                        stat.minPvs = pvsSize;

                                if ((int)viewCell->mBspLeaves.size() > stat.maxBspLeaves)
                                        stat.maxBspLeaves = (int)viewCell->mBspLeaves.size();           
                        }
                }
                else 
                {
                        BspInterior *interior = dynamic_cast<BspInterior *>(node);

                        nodeStack.push(interior->GetFront());
                        nodeStack.push(interior->GetBack());
                }
        }
}

BspTreeStatistics &VspBspTree::GetStat()
{
        return mStat;
}

float VspBspTree::AccumulatedRayLength(const RayInfoContainer &rays) const
{
        float len = 0;

        RayInfoContainer::const_iterator it, it_end = rays.end();

        for (it = rays.begin(); it != it_end; ++ it)
                len += (*it).SegmentLength();

        return len;
}

int VspBspTree::SplitRays(const Plane3 &plane,
                                                  RayInfoContainer &rays, 
                                                  RayInfoContainer &frontRays, 
                                                  RayInfoContainer &backRays)
{
        int splits = 0;

        while (!rays.empty())
        {
                RayInfo bRay = rays.back();
                VssRay *ray = bRay.mRay;

                rays.pop_back();
        
                // get classification and receive new t
                const int cf = bRay.ComputeRayIntersection(plane, ray->mT);

                switch (cf)
                {
                case -1:
                        backRays.push_back(bRay);
                        break;
                case 1:
                        frontRays.push_back(bRay);
                        break;
                case 0: 
                        //-- split ray

                        // if start point behind plane
                        if (plane.Side(bRay.ExtrapOrigin()) <= 0)
                        {       
                                backRays.push_back(RayInfo(ray, bRay.GetMinT(), ray->mT));
                                frontRays.push_back(RayInfo(ray, ray->mT, bRay.GetMaxT()));
                        }
                        else
                        {
                                frontRays.push_back(RayInfo(ray, bRay.GetMinT(), ray->mT));
                                backRays.push_back(RayInfo(ray, ray->mT, bRay.GetMaxT()));
                        }
                        break;
                default:
                        Debug << "Should not come here 4" << endl;
                        break;
                }
        }

        return splits;
}

void VspBspTree::ExtractHalfSpaces(BspNode *n, vector<Plane3> &halfSpaces) const
{
        BspNode *lastNode;

        do
        {
                lastNode = n;

                // want to get planes defining geometry of this node => don't take
                // split plane of node itself
                n = n->GetParent();
                
                if (n)
                {
                        BspInterior *interior = dynamic_cast<BspInterior *>(n);
                        Plane3 halfSpace = dynamic_cast<BspInterior *>(interior)->GetPlane();

            if (interior->GetFront() != lastNode)
                                halfSpace.ReverseOrientation();

                        halfSpaces.push_back(halfSpace);
                }
        }
        while (n);
}

void VspBspTree::ConstructGeometry(BspNode *n, 
                                                                   BspNodeGeometry &cell) const
{
        PolygonContainer polys;
        ConstructGeometry(n, polys);
        cell.mPolys = polys;
}

void VspBspTree::ConstructGeometry(BspNode *n, 
                                                                   PolygonContainer &cell) const
{
        vector<Plane3> halfSpaces;
        ExtractHalfSpaces(n, halfSpaces);

        PolygonContainer candidatePolys;

        // bounded planes are added to the polygons (reverse polygons
        // as they have to be outfacing
        for (int i = 0; i < (int)halfSpaces.size(); ++ i)
        {
                Polygon3 *p = GetBoundingBox().CrossSection(halfSpaces[i]);
                
                if (p->Valid(mEpsilon))
                {
                        candidatePolys.push_back(p->CreateReversePolygon());
                        DEL_PTR(p);
                }
        }

        // add faces of bounding box (also could be faces of the cell)
        for (int i = 0; i < 6; ++ i)
        {
                VertexContainer vertices;
        
                for (int j = 0; j < 4; ++ j)
                        vertices.push_back(mBox.GetFace(i).mVertices[j]);

                candidatePolys.push_back(new Polygon3(vertices));
        }

        for (int i = 0; i < (int)candidatePolys.size(); ++ i)
        {
                // polygon is split by all other planes
                for (int j = 0; (j < (int)halfSpaces.size()) && candidatePolys[i]; ++ j)
                {
                        if (i == j) // polygon and plane are coincident
                                continue;

                        VertexContainer splitPts;
                        Polygon3 *frontPoly, *backPoly;

                        const int cf = 
                                candidatePolys[i]->ClassifyPlane(halfSpaces[j],
                                                                                                 mEpsilon);
                        
                        switch (cf)
                        {
                                case Polygon3::SPLIT:
                                        frontPoly = new Polygon3();
                                        backPoly = new Polygon3();

                                        candidatePolys[i]->Split(halfSpaces[j], 
                                                                                         *frontPoly, 
                                                                                         *backPoly,
                                                                                         mEpsilon);

                                        DEL_PTR(candidatePolys[i]);

                                        if (frontPoly->Valid(mEpsilon))
                                                candidatePolys[i] = frontPoly;
                                        else
                                                DEL_PTR(frontPoly);

                                        DEL_PTR(backPoly);
                                        break;
                                case Polygon3::BACK_SIDE:
                                        DEL_PTR(candidatePolys[i]);
                                        break;
                                // just take polygon as it is
                                case Polygon3::FRONT_SIDE:
                                case Polygon3::COINCIDENT:
                                default:
                                        break;
                        }
                }
                
                if (candidatePolys[i])
                        cell.push_back(candidatePolys[i]);
        }
}

void VspBspTree::ConstructGeometry(BspViewCell *vc, PolygonContainer &vcGeom) const
{
        vector<BspLeaf *> leaves = vc->mBspLeaves;

        vector<BspLeaf *>::const_iterator it, it_end = leaves.end();

        for (it = leaves.begin(); it != it_end; ++ it)
                ConstructGeometry(*it, vcGeom);
}

int VspBspTree::FindNeighbors(BspNode *n, vector<BspLeaf *> &neighbors, 
                                                   const bool onlyUnmailed) const
{
        PolygonContainer cell;

        ConstructGeometry(n, cell);

        stack<BspNode *> nodeStack;
        nodeStack.push(mRoot);
                
        // planes needed to verify that we found neighbor leaf.
        vector<Plane3> halfSpaces;
        ExtractHalfSpaces(n, halfSpaces);

        while (!nodeStack.empty()) 
        {
                BspNode *node = nodeStack.top();
                nodeStack.pop();

                if (node->IsLeaf())
                {
            if (node != n && (!onlyUnmailed || !node->Mailed()))
                        {
                                // test all planes of current node if candidate really
                                // is neighbour
                                PolygonContainer neighborCandidate;
                                ConstructGeometry(node, neighborCandidate);
                                
                                bool isAdjacent = true;
                                for (int i = 0; (i < halfSpaces.size()) && isAdjacent; ++ i)
                                {
                                        const int cf = 
                                                Polygon3::ClassifyPlane(neighborCandidate, 
                                                                                                halfSpaces[i],
                                                                                                mEpsilon);

                                        if (cf == Polygon3::BACK_SIDE)
                                                isAdjacent = false;
                                }

                                if (isAdjacent)
                                        neighbors.push_back(dynamic_cast<BspLeaf *>(node));

                                CLEAR_CONTAINER(neighborCandidate);
                        }
                }
                else 
                {
                        BspInterior *interior = dynamic_cast<BspInterior *>(node);
        
                        const int cf = Polygon3::ClassifyPlane(cell, 
                                                                                                   interior->GetPlane(), 
                                                                                                   mEpsilon);

                        if (cf == Polygon3::FRONT_SIDE)
                                nodeStack.push(interior->GetFront());
                        else
                                if (cf == Polygon3::BACK_SIDE)
                                        nodeStack.push(interior->GetBack());
                                else 
                                {
                                        // random decision
                                        nodeStack.push(interior->GetBack());
                                        nodeStack.push(interior->GetFront());
                                }
                }
        }
        
        CLEAR_CONTAINER(cell);
        return (int)neighbors.size();
}

BspLeaf *VspBspTree::GetRandomLeaf(const Plane3 &halfspace)
{
    stack<BspNode *> nodeStack;
        nodeStack.push(mRoot);
        
        int mask = rand();
  
        while (!nodeStack.empty()) 
        {
                BspNode *node = nodeStack.top();
                nodeStack.pop();
          
                if (node->IsLeaf()) 
                {
                        return dynamic_cast<BspLeaf *>(node);
                } 
                else 
                {
                        BspInterior *interior = dynamic_cast<BspInterior *>(node);
                        
                        BspNode *next;
        
                        PolygonContainer cell;

                        // todo: not very efficient: constructs full cell everytime
                        ConstructGeometry(interior, cell);

                        const int cf = Polygon3::ClassifyPlane(cell, halfspace, mEpsilon);

                        if (cf == Polygon3::BACK_SIDE)
                                next = interior->GetFront();
                        else
                                if (cf == Polygon3::FRONT_SIDE)
                                        next = interior->GetFront();
                        else 
                        {
                                // random decision
                                if (mask & 1)
                                        next = interior->GetBack();
                                else
                                        next = interior->GetFront();
                                mask = mask >> 1;
                        }

                        nodeStack.push(next);
                }
        }
        
        return NULL;
}

BspLeaf *VspBspTree::GetRandomLeaf(const bool onlyUnmailed)
{
        stack<BspNode *> nodeStack;
        
        nodeStack.push(mRoot);

        int mask = rand();
        
        while (!nodeStack.empty()) 
        {
                BspNode *node = nodeStack.top();
                nodeStack.pop();
                
                if (node->IsLeaf()) 
                {
                        if ( (!onlyUnmailed || !node->Mailed()) )
                                return dynamic_cast<BspLeaf *>(node);
                }
                else 
                {
                        BspInterior *interior = dynamic_cast<BspInterior *>(node);

                        // random decision
                        if (mask & 1)
                                nodeStack.push(interior->GetBack());
                        else
                                nodeStack.push(interior->GetFront());

                        mask = mask >> 1;
                }
        }
        
        return NULL;
}

int VspBspTree::ComputePvsSize(const RayInfoContainer &rays) const
{
        int pvsSize = 0;

        RayInfoContainer::const_iterator rit, rit_end = rays.end();

        Intersectable::NewMail();

        for (rit = rays.begin(); rit != rays.end(); ++ rit)
        {
                VssRay *ray = (*rit).mRay;
                
                if (ray->mOriginObject)
                {
                        if (!ray->mOriginObject->Mailed())
                        {
                                ray->mOriginObject->Mail();
                                ++ pvsSize;
                        }
                }
                if (ray->mTerminationObject)
                {
                        if (!ray->mTerminationObject->Mailed())
                        {
                                ray->mTerminationObject->Mail();
                                ++ pvsSize;
                        }
                }
        }

        return pvsSize;
}

float VspBspTree::GetEpsilon() const
{
        return mEpsilon;
}

BspViewCell *VspBspTree::GetRootCell() const
{
        return mRootCell;
}

int VspBspTree::SplitPolygons(const Plane3 &plane,
                                                          PolygonContainer &polys, 
                                                          PolygonContainer &frontPolys, 
                                                          PolygonContainer &backPolys, 
                                                          PolygonContainer &coincident) const
{
        int splits = 0;

        while (!polys.empty())
        {
                Polygon3 *poly = polys.back();
                polys.pop_back();

                // classify polygon
                const int cf = poly->ClassifyPlane(plane, mEpsilon);

                switch (cf)
                {
                        case Polygon3::COINCIDENT:
                                coincident.push_back(poly);
                                break;                  
                        case Polygon3::FRONT_SIDE:      
                                frontPolys.push_back(poly);
                                break;
                        case Polygon3::BACK_SIDE:
                                backPolys.push_back(poly);
                                break;
                        case Polygon3::SPLIT:
                                backPolys.push_back(poly);
                                frontPolys.push_back(poly);
                                ++ splits;
                                break;
                        default:
                Debug << "SHOULD NEVER COME HERE\n";
                                break;
                }
        }

        return splits;
}