source: GTP/trunk/Lib/Vis/Preprocessing/src/ViewCell.cpp @ 1027

#include <time.h>
#include <iomanip>
#include <stack>


#include "ViewCell.h"
#include "Mesh.h"
#include "Intersectable.h"
#include "MeshKdTree.h"
#include "Triangle3.h"
#include "common.h"
#include "Environment.h"
#include "ViewCellsManager.h"
#include "Exporter.h"


namespace GtpVisibilityPreprocessor {



template <typename T> class myless
{
public:
        bool operator() (T v1, T v2) const
        {
                return (v1->GetMergeCost() < v2->GetMergeCost());
        }
};


typedef priority_queue<ViewCell *, vector<ViewCell *>, myless<vector<ViewCell *>::value_type> > TraversalQueue;

int ViewCell::sMailId = 21843194198;
int ViewCell::sReservedMailboxes = 1;


//int upperPvsLimit = 120;
//int lowerPvsLimit = 5;

float MergeCandidate::sRenderCostWeight = 0;


// pvs penalty can be different from pvs size
inline static float EvalPvsPenalty(const int pvs, 
                                                                   const int lower, 
                                                                   const int upper)
{
        // clamp to minmax values
        if (pvs < lower)
                return (float)lower;
        if (pvs > upper)
                return (float)upper;

        return (float)pvs;
}

/// Counts differences between pvss.
inline int CountDiffPvs(ViewCell *vc)
{
        int count = 0;

        ObjectPvsMap::const_iterator it, it_end = vc->GetPvs().mEntries.end();
        for (it = vc->GetPvs().mEntries.begin(); it != it_end; ++ it)
        {
                if (!(*it).first->Mailed())
                {
                        (*it).first->Mail();
                        ++ count;
                }
        }

        return count;
}

/// Fast computation of merged pvs size
static int ComputeMergedPvsSize(const ObjectPvs &pvs1, const ObjectPvs &pvs2)
{
        // add first pvs
        int pvs = pvs1.GetSize();

        ObjectPvsMap::const_iterator it, it_end =  pvs1.mEntries.end();

        Intersectable::NewMail();

        // mail all objects in first pvs
        for (it = pvs1.mEntries.begin(); it != it_end; ++ it)
        {
                (*it).first->Mail();
        }

        it_end = pvs2.mEntries.end();

        // look if they are in second pvs
        for (it = pvs2.mEntries.begin(); it != it_end; ++ it)
        {
                Intersectable *obj = (*it).first;
                if (!obj->Mailed())
                        ++ pvs;
        }

        return pvs;
}


ViewCell::ViewCell(): 
MeshInstance(NULL), 
mPiercingRays(0),
mArea(-1),
mVolume(-1),
mValid(true),
mParent(NULL),
mMergeCost(0),
mPvsSize(0),
mPvsSizeValid(false)
{
}

ViewCell::ViewCell(Mesh *mesh): 
MeshInstance(mesh), 
mPiercingRays(0),
mArea(-1),
mVolume(-1),
mValid(true),
mParent(NULL),
mMergeCost(0),
mPvsSize(0),
mPvsSizeValid(false)
{
}


const ObjectPvs &ViewCell::GetPvs() const
{
        return mPvs;
}


ObjectPvs &ViewCell::GetPvs()
{
        return mPvs;
}


void ViewCell::SetPvs(const ObjectPvs &pvs)
{
        mPvs = pvs;
}


int ViewCell::Type() const 
{ 
        return VIEW_CELL; 
}


float ViewCell::GetVolume() const
{
        return mVolume;
}


void ViewCell::SetVolume(float volume) 
{
        mVolume = volume;
}


void ViewCell::SetMesh(Mesh *mesh)
{
        mMesh = mesh;
}


float ViewCell::GetArea() const
{
        return mArea;
}


void ViewCell::SetArea(float area) 
{
        mArea = area;
}


void ViewCell::SetColor(const RgbColor &color) 
{
        mColor = color;
}


RgbColor ViewCell::GetColor() const 
{
        return mColor;
}


void ViewCell::SetValid(const bool valid)
{
        mValid = valid;
}


bool ViewCell::GetValid() const
{
        return mValid;
}


void ViewCell::SetParent(ViewCellInterior *parent)
{
        mParent = parent;
}


bool ViewCell::IsRoot() const
{
        return !mParent;
}


ViewCellInterior *ViewCell::GetParent() const
{
        return mParent;
}


void ViewCell::SetMergeCost(const float mergeCost)
{
        mMergeCost = mergeCost;
}


float ViewCell::GetRenderCost() const
{
        //return (float)mPvs.GetSize() * GetVolume();
        return (float)mPvsSize * GetVolume();
}


float ViewCell::GetMergeCost() const
{
        return mMergeCost;
}


bool ViewCell::AddPvsSample(Intersectable *sample,
                                                        const float pdf, 
                                                        float &contribution)
{
        const bool result = mPvs.AddSample(sample, pdf, contribution);

        // update pvs size
        mPvsSize = mPvs.GetSize();
        mPvsSizeValid = true;

        return result;
}




/************************************************************************/
/*                class ViewCellInterior implementation                 */
/************************************************************************/


ViewCellInterior::ViewCellInterior()
{
}


ViewCellInterior::~ViewCellInterior()
{
        ViewCellContainer::const_iterator it, it_end = mChildren.end();

        for (it = mChildren.begin(); it != it_end; ++ it)
                delete (*it);
}


ViewCellInterior::ViewCellInterior(Mesh *mesh): 
ViewCell(mesh)
{
}


bool ViewCellInterior::IsLeaf() const
{
        return false;
}


void ViewCellInterior::SetupChildLink(ViewCell *vc)
{
    mChildren.push_back(vc);
    vc->SetParent(this);
}


void ViewCellInterior::RemoveChildLink(ViewCell *l)
{
        // erase leaf from old view cell
        ViewCellContainer::iterator it = mChildren.begin();

        for (; (*it) != l; ++ it);
        if (it == mChildren.end())
                Debug << "error" << endl;
        else
                mChildren.erase(it);
}


/************************************************************************/
/*                class ViewCellsStatistics implementation              */
/************************************************************************/




void ViewCellsStatistics::Print(ostream &app) const
{
        app << "=========== View Cells Statistics ===============\n";

        app << setprecision(4);

        //app << "#N_CTIME  ( Construction time [s] )\n" << Time() << " \n";

        app << "#N_OVERALLPVS ( objects in PVS )\n" << pvsSize << endl;

        app << "#N_PMAXPVS ( largest PVS )\n" << maxPvs << endl;

        app << "#N_PMINPVS ( smallest PVS )\n" << minPvs << endl;

        app << "#N_PAVGPVS ( average PVS )\n" << AvgPvs() << endl;

        app << "#N_PEMPTYPVS ( view cells with empty PVS )\n" << emptyPvs << endl;

        app << "#N_VIEWCELLS ( number of view cells)\n" << viewCells << endl;

        app << "#N_AVGLEAVES (average number of leaves per view cell )\n" << AvgLeaves() << endl;

        app << "#N_MAXLEAVES ( maximal number of leaves per view cell )\n" << maxLeaves << endl;
        
        app << "#N_INVALID ( number of invalid view cells )\n" << invalid << endl;

        app << "========== End of View Cells Statistics ==========\n";
}


/*************************************************************************/
/*                    class ViewCellsTree implementation                 */
/*************************************************************************/


ViewCellsTree::ViewCellsTree(ViewCellsManager *vcm):
mRoot(NULL),
mUseAreaForPvs(false),
mViewCellsManager(vcm),
#if 0
mViewCellsStorage(PVS_IN_INTERIORS)
#else
mViewCellsStorage(PVS_IN_LEAVES)
#endif
{
        Environment::GetSingleton()->GetBoolValue("ViewCells.Visualization.exportMergedViewCells", mExportMergedViewCells);
        Environment::GetSingleton()->GetFloatValue("ViewCells.maxStaticMemory", mMaxMemory);

        //-- merge options
        Environment::GetSingleton()->GetFloatValue("ViewCells.PostProcess.renderCostWeight", mRenderCostWeight);
        Environment::GetSingleton()->GetIntValue("ViewCells.PostProcess.minViewCells", mMergeMinViewCells);
        Environment::GetSingleton()->GetFloatValue("ViewCells.PostProcess.maxCostRatio", mMergeMaxCostRatio);
        Environment::GetSingleton()->GetBoolValue("ViewCells.PostProcess.refine", mRefineViewCells);    

        Environment::GetSingleton()->GetIntValue("ViewCells.PostProcess.maxMergesPerPass", mMaxMergesPerPass);
        Environment::GetSingleton()->GetFloatValue("ViewCells.PostProcess.avgCostMaxDeviation", mAvgCostMaxDeviation);

        Debug << "============= view cell tree options ================\n";
        Debug << "minimum view cells: " << mMergeMinViewCells << endl;
        Debug << "max cost ratio: " << mMergeMaxCostRatio << endl;
        Debug << "max memory: " << mMaxMemory << endl;
        Debug << "refining view cells: " << mRefineViewCells << endl;
        Debug << "=========== end view cell tree options ===============\n";

        MergeCandidate::sRenderCostWeight = mRenderCostWeight;
}


// return memory usage in MB
float ViewCellsTree::GetMemUsage() const
{
        // TODO
        return 0;
                /*(sizeof(ViewCellsTree) +
                 mBspStats.Leaves() * sizeof(BspLeaf) +
                 mBspStats.Interior() * sizeof(BspInterior) +
                 mBspStats.accumRays * sizeof(RayInfo)) / (1024.0f * 1024.0f);*/
}


int ViewCellsTree::GetNumInitialViewCells(ViewCell *vc) const
{
        int vcSize = 0;

        stack<ViewCell *> tstack;

        tstack.push(vc);

        while (!tstack.empty())
        {
                ViewCell *vc = tstack.top();
                tstack.pop();

                if (vc->IsLeaf())
                {
                        ++ vcSize;
                }
                else
                {
                        ViewCellInterior *interior = dynamic_cast<ViewCellInterior *>(vc);
                        ViewCellContainer::const_iterator it, it_end = interior->mChildren.end();
                        for (it = interior->mChildren.begin(); it != it_end; ++ it)
                                tstack.push(*it);
                        
                }
        }

        return vcSize;
}


void ViewCellsTree::CollectLeaves(ViewCell *vc, ViewCellContainer &leaves) const
{
        stack<ViewCell *> tstack;

        tstack.push(vc);

        while (!tstack.empty())
        {
                ViewCell *vc = tstack.top();
                tstack.pop();

                if (vc->IsLeaf())
                {
                        leaves.push_back(vc);
                }
                else
                {
                        ViewCellInterior *interior = dynamic_cast<ViewCellInterior *>(vc);
                        ViewCellContainer::const_iterator it, it_end = interior->mChildren.end();

                        for (it = interior->mChildren.begin(); it != it_end; ++ it)
                        {
                                tstack.push(*it);
                        }
                        
                }
        }
}


ViewCellsTree::~ViewCellsTree()
{
        DEL_PTR(mRoot);
}


int ViewCellsTree::ConstructMergeTree(const VssRayContainer &rays, 
                                                                          const ObjectContainer &objects)
{
        mNumActiveViewCells = (int)mViewCellsManager->GetViewCells().size();

        float variance = 0;
        int totalPvs = 0;
        float totalRenderCost = 0;

        //-- compute statistics values of initial view cells
        mViewCellsManager->EvaluateRenderStatistics(totalRenderCost,
                                                                                                mExpectedCost,
                                                                                                mDeviation,
                                                                                                variance,
                                                                                                totalPvs,
                                                                                                mAvgRenderCost);


        //-- fill merge queue
        vector<MergeCandidate> candidates;

        mViewCellsManager->CollectMergeCandidates(rays, candidates);

        while(!candidates.empty())
        {
                MergeCandidate mc = candidates.back();
                candidates.pop_back();
                EvalMergeCost(mc);
                mMergeQueue.push(mc);
        }

        Debug << "************************* merge ***********************************" << endl;  
        Debug << "deviation: " << mDeviation << endl;
        Debug << "avg render cost: " << mAvgRenderCost << endl;
        Debug << "expected cost: " << mExpectedCost << endl;


        ViewCellsManager::PvsStatistics pvsStats;
        mViewCellsManager->GetPvsStatistics(pvsStats);

        //static float expectedValue = pvsStats.avgPvs;
        
        //-- the current view cells are kept in this container
        //-- we start with the current view cells from the view cell manager. 
        //-- The active view cells will change with subsequent merges
        
        // todo: should rather take initial view cells
    ViewCellContainer &activeViewCells = mViewCellsManager->GetViewCells();
        
        
        ViewCell::NewMail();

        MergeStatistics mergeStats;
        mergeStats.Start();
        
        long startTime = GetTime();

        mergeStats.collectTime = TimeDiff(startTime, GetTime());
        mergeStats.candidates = (int)mMergeQueue.size();
        startTime = GetTime();

        // frequency stats are updated
        const int statsOut = 500;
        
        // passes are needed for statistics, because we don't want to record
        // every merge
        int pass = 0;
        int mergedPerPass = 0;
        float realExpectedCost = mExpectedCost;
        float realAvgRenderCost = mAvgRenderCost;
        int realNumActiveViewCells = mNumActiveViewCells;
        
        // maximal ratio of old expected render cost to expected render
        // when the the render queue has to be reset.
        int numMergedViewCells = 0;
                

        cout << "actual merge starts now ... " << endl;

        //-- use priority queue to merge leaf pairs

        while (!mMergeQueue.empty())
        {
                //-- reset merge queue if the ratio of current expected cost / real expected cost
                //   too small or after a given number of merges
                if ((mergedPerPass > mMaxMergesPerPass) ||
                        (mAvgCostMaxDeviation > mAvgRenderCost / realAvgRenderCost))
                {
                        Debug << "************ reset queue *****************\n"
                                  << "ratios: " << mAvgCostMaxDeviation 
                                  << " real avg render cost " << realAvgRenderCost << " average render cost " << mAvgRenderCost
                                  << " merged per pass : " << mergedPerPass << " of maximal " << mMaxMergesPerPass << endl;

                        Debug << "Values before reset: "  
                                  << " erc: " << mExpectedCost 
                                  << " avgrc: " << mAvgRenderCost
                                  << " dev: " << mDeviation << endl;
        
                        // adjust render cost
                        ++ pass;

                        mergedPerPass = 0;
                        mExpectedCost = realExpectedCost;
                        mAvgRenderCost = realAvgRenderCost;
                        mNumActiveViewCells = realNumActiveViewCells;
                        
                        const int numMergedViewCells = UpdateActiveViewCells(activeViewCells);
                
                        
                        //-- resets / refines the view cells
                        //-- priorities are recomputed
                        //-- the candidates are put back into merge queue
                        if (mRefineViewCells)
                                RefineViewCells(rays, objects);
                        else
                                ResetMergeQueue();

                        Debug << "Values after reset: "  
                                  << " erc: " << mExpectedCost 
                                  << " avg: " << mAvgRenderCost 
                                  << " dev: " << mDeviation << endl;

                        if (mExportMergedViewCells)
                        {
                                ExportMergedViewCells(activeViewCells, objects, numMergedViewCells);
                        }
                }


                MergeCandidate mc = mMergeQueue.top();
                mMergeQueue.pop();
        
                // both view cells equal because of previous merges
                // NOTE: do I really still need this? probably cannot happen!!
                if (mc.mLeftViewCell == mc.mRightViewCell)
                        continue;

                if (mc.IsValid())
                {
                        ViewCell::NewMail();

                        //-- update statistical values
                        -- realNumActiveViewCells;
                        ++ mergeStats.merged;
                        ++ mergedPerPass;

                        const float renderCostIncr = mc.GetRenderCost();
                        const float mergeCostIncr = mc.GetMergeCost();

                        totalRenderCost += renderCostIncr;
                        mDeviation += mc.GetDeviationIncr();

                                                
                        //-- merge the view cells of leaf1 and leaf2
                        int pvsDiff;
                        ViewCellInterior *mergedVc = 
                                MergeViewCells(mc.mLeftViewCell, mc.mRightViewCell, pvsDiff);
                        

                        // total render cost and deviation has changed
                        // real expected cost will be larger than expected cost used for the
                        // cost heuristics, but cannot recompute costs on each increase of the 
                        // expected cost
                        totalPvs += pvsDiff;
                        realExpectedCost = totalRenderCost / (float)realNumActiveViewCells;
                        realAvgRenderCost = (float)totalPvs / (float)realNumActiveViewCells;
        
                        // set merge cost to this node for priority traversal
                        mergedVc->SetMergeCost(totalRenderCost);
                        
                        // HACK
                        //mergedVc->SetMergeCost(1.0f / (float)realNumActiveViewCells);

                        // check if "siblings (back and front node of the same parent)
                        if (0 && mViewCellsManager->EqualToSpatialNode(mergedVc))
                                ++ mergeStats.siblings;
                        // set the coŽst for rendering a view cell
                        mergedVc->SetCost(realExpectedCost);

                        if ((mergeStats.merged % statsOut) == 0)
                                cout << "merged " << mergeStats.merged << " view cells" << endl;

                }
                else
                { 
                        // merge candidate not valid, because one of the leaves was already
                        // merged with another one => validate and reinsert into queue
                        if (ValidateMergeCandidate(mc))
                        {
                                EvalMergeCost(mc);
                                mMergeQueue.push(mc);
                        }
                }
                
        }

        // adjust stats and reset queue one final time
        mExpectedCost = realExpectedCost;
        mAvgRenderCost = realAvgRenderCost;
        mNumActiveViewCells = realNumActiveViewCells;

        UpdateActiveViewCells(activeViewCells);

        // refine view cells and reset costs
        if (mRefineViewCells)
                RefineViewCells(rays, objects);
        else
                ResetMergeQueue();


        // create a root node if the merge was not done till root level,
        // else take the single node as new root
        if ((int)activeViewCells.size() > 1)
        {
                Debug << "creating root of view cell hierarchy for " 
                          << (int)activeViewCells.size() << " view cells" << endl;
                
                ViewCellInterior *root = mViewCellsManager->MergeViewCells(activeViewCells);
        
                ViewCellContainer::const_iterator it, it_end = root->mChildren.end();

                for (it = root->mChildren.begin(); it != it_end; ++ it)
                        (*it)->SetParent(root);
        
                root->SetMergeCost(totalRenderCost);
                // $$JB keep this 0 temporarilly
                root->SetCost(0.0f);

                mRoot = root;
        }
        // normal case
        else if (!activeViewCells.empty())
        {
                Debug << "setting root of the merge history" << endl;
                mRoot = activeViewCells[0];
        }
        else
        {
                Debug << "big error, root is NULL" << endl;
        }
        
        //-- empty merge queue just in case
        while (!mMergeQueue.empty())
        {
                mMergeQueue.pop();
        }

        // test if voluje is equal
        Debug << "volume of the root view cell: " << mRoot->GetVolume() << " " << mViewCellsManager->GetViewSpaceBox().GetVolume() << endl;

        //hack!!
        //mRoot->GetPvs().Clear();
        
        // TODO: delete because makes no sense here
        mergeStats.expectedRenderCost = realExpectedCost;
        mergeStats.deviation = mDeviation;

        // we want to optimize this heuristics
        mergeStats.heuristics = 
                mDeviation * (1.0f - mRenderCostWeight) + 
                mExpectedCost * mRenderCostWeight;

        mergeStats.mergeTime = TimeDiff(startTime, GetTime());
        mergeStats.Stop();
        Debug << mergeStats << endl << endl;

        // assign colors for the view cells so that at least one is always consistent
        AssignRandomColors();

        //TODO: should return sample contributions?
        return mergeStats.merged;
}


ViewCell *ViewCellsTree::GetRoot() const
{
        return mRoot;
}


void ViewCellsTree::ResetMergeQueue()
{
        cout << "reset merge queue ... ";
        
        vector<MergeCandidate> buf;
        buf.reserve(mMergeQueue.size());
                        
        
        // store merge candidates in intermediate buffer
        while (!mMergeQueue.empty())
        {
                MergeCandidate mc = mMergeQueue.top();
                mMergeQueue.pop();
                
                // recalculate cost
                if (ValidateMergeCandidate(mc))
                {
                        EvalMergeCost(mc);
                        buf.push_back(mc);                              
                }
        }

        vector<MergeCandidate>::const_iterator bit, bit_end = buf.end();

        // reinsert back into queue
        for (bit = buf.begin(); bit != bit_end; ++ bit)
        {      
                mMergeQueue.push(*bit); 
        }

        cout << "finished" << endl;
}


float ViewCellsTree::ComputeMergedPvsCost(const ObjectPvs &pvs1, 
                                                                                  const ObjectPvs &pvs2) const
{
        // computes render cost of merge
        float renderCost = 0;

        // compute new pvs size
        ObjectPvsMap::const_iterator it, it_end =  pvs1.mEntries.end();

        Intersectable::NewMail();

        // first mail all objects in first pvs
        for (it = pvs1.mEntries.begin(); it != it_end; ++ it)
        {
                Intersectable *obj = (*it).first;

                obj->Mail();
                renderCost += mViewCellsManager->EvalRenderCost(obj);
        }

        it_end = pvs2.mEntries.end();

        
        for (it = pvs2.mEntries.begin(); it != it_end; ++ it)
        {
                Intersectable *obj = (*it).first;

                // test if object already considered    
                if (!obj->Mailed())
                {
                        renderCost += mViewCellsManager->EvalRenderCost(obj);
                }
        }

        return renderCost;
}


int ViewCellsTree::UpdateActiveViewCells(ViewCellContainer &viewCells)
{
        int numMergedViewCells = 0;

        Debug << "updating active vc: " << (int)viewCells.size() << endl;
        // find all already merged view cells and remove them from view cells
                
        // sort out all view cells which are not active anymore, i.e., they
        // were already part of a merge
        int i = 0;

        ViewCell::NewMail();

        while (1)
        {
                // remove all merged view cells from end of the vector
                while (!viewCells.empty() && (viewCells.back()->GetParent()))
                {
                        viewCells.pop_back();
                }

                // all merged view cells have been found
                if (i >= (int)viewCells.size()) 
                        break;

                // already merged this view cell, put it to end of vector
                if (viewCells[i]->GetParent())
                        swap(viewCells[i], viewCells.back());
                
                // mail view cell that it has not been merged
                viewCells[i]->Mail();

                // increase loop counter
                ++ i;
        }


        // add new view cells to container only if they don't have been 
        // merged in the mean time
        ViewCellContainer::const_iterator ait, ait_end = mMergedViewCells.end();

        for (ait = mMergedViewCells.begin(); ait != ait_end; ++ ait)
        {
                ViewCell *vc = mMergedViewCells.back();
                if (!vc->GetParent() && !vc->Mailed())
                {
                        vc->Mail();
                        viewCells.push_back(vc);
                        ++ numMergedViewCells;
                }
        }

        // dispose old merged view cells
        mMergedViewCells.clear();

        // update standard deviation
        ViewCellContainer::const_iterator vit, vit_end = viewCells.end();
        
        mDeviation = 0;

        for (vit = viewCells.begin(); vit != vit_end; ++ vit)
        {
                const int lower = mViewCellsManager->GetMinPvsSize();
                const int upper = mViewCellsManager->GetMaxPvsSize();

                const float penalty = EvalPvsPenalty((*vit)->GetPvs().GetSize(), lower, upper);
                
                mDeviation += fabs(mAvgRenderCost - penalty);
        }

        mDeviation /= (float)viewCells.size();
        
        return numMergedViewCells;
}


void ViewCellsTree::ExportMergedViewCells(ViewCellContainer &viewCells, 
                                                                                  const ObjectContainer &objects,
                                                                                  const int numMergedViewCells)
{
        

        char s[64];

        sprintf(s, "merged_viewcells%07d.x3d", (int)viewCells.size());
        Exporter *exporter = Exporter::GetExporter(s);

        if (exporter)
        {
                cout << "exporting " << (int)viewCells.size() << " merged view cells ... ";
                exporter->ExportGeometry(objects);
                //Debug << "vc size " << (int)viewCells.size() << " merge queue size: " << (int)mMergeQueue.size() << endl;
                ViewCellContainer::const_iterator it, it_end = viewCells.end();

                int i = 0;
                for (it = viewCells.begin(); it != it_end; ++ it)
                {
                        Material m;
                        // assign special material to new view cells
                        // new view cells are on the back of container
                        if (i ++ >= (viewCells.size() - numMergedViewCells))
                        {
                                //m = RandomMaterial();
                                m.mDiffuseColor.r = RandomValue(0.5f, 1.0f);
                                m.mDiffuseColor.g = RandomValue(0.5f, 1.0f);
                                m.mDiffuseColor.b = RandomValue(0.5f, 1.0f);
                        }
                        else
                        {
                                float col = RandomValue(0.1f, 0.4f);
                                m.mDiffuseColor.r = col;
                                m.mDiffuseColor.g = col;
                                m.mDiffuseColor.b = col;
                        }

                        exporter->SetForcedMaterial(m);
                        mViewCellsManager->ExportViewCellGeometry(exporter, *it);
                }
                delete exporter;
                cout << "finished" << endl;
        }
}


// TODO: should be done in view cells manager
ViewCellInterior *ViewCellsTree::MergeViewCells(ViewCell *l, 
                                                                                                ViewCell *r, 
                                                                                                int &pvsDiff) //const
{
        ViewCellInterior *vc = mViewCellsManager->MergeViewCells(l, r);

        // if merge was unsuccessful
        if (!vc) return NULL;

        // set to the new parent view cell
        l->SetParent(vc);
        r->SetParent(vc);

        // set new size of view cell
        if (mUseAreaForPvs)
        {
                vc->SetArea(l->GetArea() + l->GetArea());
        }
        else
        {
                vc->SetVolume(r->GetVolume() + l->GetVolume());
        }

        
        // important so other merge candidates sharing this view cell
        // are notified that the merge cost must be updated!!
        vc->Mail();

        const int pvs1 = l->GetPvs().GetSize();
        const int pvs2 = r->GetPvs().GetSize();


        // new view cells are stored in this vector
        mMergedViewCells.push_back(vc);

        pvsDiff = vc->GetPvs().GetSize() - pvs1 - pvs2;



        //Ždon't store intermediate pvs
        if (mViewCellsStorage == PVS_IN_LEAVES)
        {
                l->mPvsSize = l->GetPvs().GetSize(); 
                l->mPvsSizeValid = true;
                
                if (!l->IsLeaf())
                        l->GetPvs().Clear();
                
                r->mPvsSize = r->GetPvs().GetSize(); 
                r->mPvsSizeValid = true;
                
                if (!r->IsLeaf())
                        r->GetPvs().Clear();
        
}


        return vc;
}


int ViewCellsTree::RefineViewCells(const VssRayContainer &rays, 
                                                                   const ObjectContainer &objects)
{
        Debug << "refining " << (int)mMergeQueue.size() << " candidates " << endl;

        // intermediate buffer for shuffled view cells
        vector<MergeCandidate> buf;
        buf.reserve(mMergeQueue.size());
                        
        // Use priority queue of remaining leaf pairs 
        // The candidates either share the same view cells or
        // are border leaves which share a boundary.
        // We test if they can be shuffled, i.e.,
        // either one leaf is made part of one view cell or the other
        // leaf is made part of the other view cell. It is tested if the
        // remaining view cells are "better" than the old ones.
        
        const int numPasses = 3;
        int pass = 0;
        int passShuffled = 0;
        int shuffled = 0;
        int shuffledViewCells = 0;

        ViewCell::NewMail();
        
        while (!mMergeQueue.empty())
        {
                MergeCandidate mc = mMergeQueue.top();
                mMergeQueue.pop();

                // both view cells equal or already shuffled
                if ((mc.GetLeftViewCell() == mc.GetRightViewCell()) ||
                        mc.GetLeftViewCell()->IsLeaf() || mc.GetRightViewCell()->IsLeaf())
                {                       
                        continue;
                }

                // candidate for shuffling
                const bool wasShuffled = ShuffleLeaves(mc);
                
                // shuffled or put into other queue for further refine
                if (wasShuffled)
                {
                        ++ passShuffled;

                        if (!mc.GetLeftViewCell()->Mailed())
                        {
                                mc.GetLeftViewCell()->Mail();
                                ++ shuffledViewCells;
                        }
                        if (!mc.GetRightViewCell()->Mailed())
                        {
                                mc.GetRightViewCell()->Mail();
                                ++ shuffledViewCells;
                        }
                }

                // put back into intermediate vector
                buf.push_back(mc);
        }


        //-- in the end, the candidates must be in the mergequeue again
        //   with the correct cost

        cout << "reset merge queue ... ";
        
        
        vector<MergeCandidate>::const_iterator bit, bit_end = buf.end();
        
        for (bit = buf.begin(); bit != bit_end; ++ bit)
        {    
                MergeCandidate mc = *bit;
                // recalculate cost
                if (ValidateMergeCandidate(mc))
                {
                        EvalMergeCost(mc);
                        mMergeQueue.push(mc);   
                }
        }

        cout << "finished" << endl;

        return shuffledViewCells;
}


inline int AddedPvsSize(ObjectPvs pvs1, const ObjectPvs &pvs2)
{
        return pvs1.AddPvs(pvs2);
}


// recomputes pvs size minus pvs of leaf l
#if 0
inline int SubtractedPvsSize(BspViewCell *vc, BspLeaf *l, const ObjectPvs &pvs2)
{
        ObjectPvs pvs;
        vector<BspLeaf *>::const_iterator it, it_end = vc->mLeaves.end();
        for (it = vc->mLeaves.begin(); it != vc->mLeaves.end(); ++ it)
                if (*it != l)
                        pvs.AddPvs(*(*it)->mPvs);
        return pvs.GetSize();
}
#endif


// computes pvs1 minus pvs2
inline int SubtractedPvsSize(ObjectPvs pvs1, const ObjectPvs &pvs2)
{
        return pvs1.SubtractPvs(pvs2);
}


float ViewCellsTree::EvalShuffleCost(ViewCell *leaf, 
                                                                         ViewCellInterior *vc1, 
                                                                         ViewCellInterior *vc2) const
{
        //const int pvs1 = SubtractedPvsSize(vc1, leaf, *leaf->mPvs);
        const int pvs1 = SubtractedPvsSize(vc1->GetPvs(), leaf->GetPvs());
        const int pvs2 = AddedPvsSize(vc2->GetPvs(), leaf->GetPvs());

        const int lowerPvsLimit = mViewCellsManager->GetMinPvsSize();
        const int upperPvsLimit = mViewCellsManager->GetMaxPvsSize();

        const float pvsPenalty1 = 
                EvalPvsPenalty(pvs1, lowerPvsLimit, upperPvsLimit);

        const float pvsPenalty2 =
                EvalPvsPenalty(pvs2, lowerPvsLimit, upperPvsLimit);


        // don't shuffle leaves with pvs > max
        if (0 && (pvs1 + pvs2 > mViewCellsManager->GetMaxPvsSize()))
        {
                return 1e20f;
        }

        float p1, p2;

    if (mUseAreaForPvs)
        {
                p1 = vc1->GetArea() - leaf->GetArea();
                p2 = vc2->GetArea() + leaf->GetArea();
        }
        else
        {
                p1 = vc1->GetVolume() - leaf->GetVolume();
                p2 = vc2->GetVolume() + leaf->GetVolume();
        }

        const float renderCost1 = pvsPenalty1 * p1;
        const float renderCost2 = pvsPenalty2 * p2;

        float dev1, dev2;

        if (1)
        {
                dev1 = fabs(mAvgRenderCost - pvsPenalty1);
                dev2 = fabs(mAvgRenderCost - pvsPenalty2);
        }
        else
        {
                dev1 = fabs(mExpectedCost - renderCost1);
                dev2 = fabs(mExpectedCost - renderCost2);
        }
        
        return mRenderCostWeight * (renderCost1 + renderCost2) +
                  (1.0f - mRenderCostWeight) * (dev1 + dev2) / (float)mNumActiveViewCells;
}


void ViewCellsTree::ShuffleLeaf(ViewCell *leaf,
                                                                ViewCellInterior *vc1, 
                                                                ViewCellInterior *vc2) const
{
        // compute new pvs and area
        // TODO change
        vc1->GetPvs().SubtractPvs(leaf->GetPvs());
        vc2->GetPvs().AddPvs(leaf->GetPvs());
        
        if (mUseAreaForPvs)
        {
                vc1->SetArea(vc1->GetArea() - leaf->GetArea());
                vc2->SetArea(vc2->GetArea() + leaf->GetArea());
        }
        else
        {
                vc1->SetVolume(vc1->GetVolume() - leaf->GetVolume());
                vc2->SetVolume(vc2->GetVolume() + leaf->GetVolume());
        }

        
        ViewCellInterior *p = dynamic_cast<ViewCellInterior *>(leaf->GetParent());

        p->RemoveChildLink(leaf);
        vc2->SetupChildLink(leaf);
}


bool ViewCellsTree::ShuffleLeaves(MergeCandidate &mc) const
{
        float cost1, cost2;

        ViewCellInterior *vc1 = dynamic_cast<ViewCellInterior *>(mc.GetLeftViewCell());
        ViewCellInterior *vc2 = dynamic_cast<ViewCellInterior *>(mc.GetRightViewCell());

        ViewCell *leaf1 = mc.GetInitialLeftViewCell();
        ViewCell *leaf2 = mc.GetInitialRightViewCell();

        //-- first test if shuffling would decrease cost
        cost1 = GetCostHeuristics(vc1);
        cost2 = GetCostHeuristics(vc2);

        const float oldCost = cost1 + cost2;
        
        float shuffledCost1 = Limits::Infinity;
        float shuffledCost2 = Limits::Infinity;

        if (leaf1)
                shuffledCost1 = EvalShuffleCost(leaf1, vc1, vc2);       
        if (leaf2)
                shuffledCost2 = EvalShuffleCost(leaf2, vc2, vc1);

        // if cost of shuffle is less than old cost => shuffle
        if ((oldCost <= shuffledCost1) && (oldCost <= shuffledCost2))
                return false;
        
        if (shuffledCost1 < shuffledCost2)
        {
                if (leaf1)
                        ShuffleLeaf(leaf1, vc1, vc2);
                mc.mInitialLeftViewCell = NULL;
        }
        else
        {
                if (leaf2)
                        ShuffleLeaf(leaf2, vc2, vc1);
                mc.mInitialRightViewCell = NULL;
        }

        return true;
}


float ViewCellsTree::GetVariance(ViewCell *vc) const
{
        const int upper = mViewCellsManager->GetMaxPvsSize();
        const int lower = mViewCellsManager->GetMinPvsSize();

        if (1)
        {
                const float penalty = 
                        EvalPvsPenalty(vc->GetPvs().GetSize(), lower, upper);
                return (mAvgRenderCost - penalty) * (mAvgRenderCost - penalty) / 
                        (float)mNumActiveViewCells;
        }

    const float leafCost = GetRenderCost(vc);
        return (mExpectedCost - leafCost) * (mExpectedCost - leafCost);
}


float ViewCellsTree::GetDeviation(ViewCell *vc) const
{
        const int upper = mViewCellsManager->GetMaxPvsSize();
        const int lower = mViewCellsManager->GetMinPvsSize();

        if (1)
        {
                const float penalty = EvalPvsPenalty(vc->GetPvs().GetSize(), lower, upper);
                return fabs(mAvgRenderCost - penalty) / (float)mNumActiveViewCells;
        }

    const float renderCost = GetRenderCost(vc);
        return fabs(mExpectedCost - renderCost);
}



float ViewCellsTree::GetRenderCost(ViewCell *vc) const
{
        if (mUseAreaForPvs)
                return vc->GetPvs().GetSize() * vc->GetArea();

        return vc->GetPvs().GetSize() * vc->GetVolume();
}


float ViewCellsTree::GetCostHeuristics(ViewCell *vc) const
{
        return GetRenderCost(vc) * mRenderCostWeight + 
                   GetDeviation(vc) * (1.0f - mRenderCostWeight);
}


bool ViewCellsTree::ValidateMergeCandidate(MergeCandidate &mc) const
{
        // propagate up so we have only the active view cells
        while (mc.mLeftViewCell->mParent)
        {
                mc.mLeftViewCell = mc.mLeftViewCell->mParent;
        }

        while (mc.mRightViewCell->mParent)
        {
                mc.mRightViewCell = mc.mRightViewCell->mParent;
        }

        // this view cell was already merged
        //return mc.mLeftViewCell && (mc.mLeftViewCell != mc.mRightViewCell);
        return mc.mLeftViewCell != mc.mRightViewCell;
}


void ViewCellsTree::EvalMergeCost(MergeCandidate &mc) const
{
        //-- compute pvs difference
        int newPvs;
        if (1) // not valid if not using const cost per object!!
                newPvs = ComputeMergedPvsSize(mc.mLeftViewCell->GetPvs(), mc.mRightViewCell->GetPvs());
        else
                newPvs = (int)ComputeMergedPvsCost(mc.mLeftViewCell->GetPvs(), mc.mRightViewCell->GetPvs());

        const float newPenalty = EvalPvsPenalty(newPvs,
                                                                                        mViewCellsManager->GetMinPvsSize(),
                                                                                        mViewCellsManager->GetMaxPvsSize());

        ViewCell *vc1 = mc.mLeftViewCell;
        ViewCell *vc2 = mc.mRightViewCell;

        //-- compute ratio of old cost
        //   (i.e., added size of left and right view cell times pvs size)
        //   to new rendering cost (i.e, size of merged view cell times pvs size)
        const float oldCost = GetRenderCost(vc1) + GetRenderCost(vc2);

    const float newCost = mUseAreaForPvs ? 
                (float)newPenalty * (vc1->GetArea() + vc2->GetArea()) :
                (float)newPenalty * (vc1->GetVolume() + vc2->GetVolume());


        // strong penalty if pvs size too large
        if (0 && (newPvs > mViewCellsManager->GetMaxPvsSize()))
        {
                mc.mRenderCost = 1e20f;
        }
        else
        {
                mc.mRenderCost = (newCost - oldCost) / 
                        mViewCellsManager->GetViewSpaceBox().GetVolume();
        }       
        

        //-- merge cost also takes deviation into account
        float newDev, oldDev;

        if (1)
                newDev = fabs(mAvgRenderCost - newPenalty) / (float)mNumActiveViewCells;
        else
                newDev = fabs(mExpectedCost - newCost) / (float)mNumActiveViewCells;
        
        oldDev = GetDeviation(vc1) + GetDeviation(vc2);

        // compute deviation increase
        mc.mDeviationIncr = newDev - oldDev;
        
        //Debug << "render cost: " << mc.mRenderCost * mRenderCostWeight << endl;
        //Debug << "standard deviation: " << mc.mDeviationIncr * mRenderCostWeight << endl;
}


void ViewCellsTree::SetViewCellsStorage(int stype)
{
        if (mViewCellsStorage == stype)
                return;

        // TODO
        switch (stype)
        {
        case COMPRESSED:
                CompressViewCellsPvs(mRoot);
                break;
        default:
                break;
        }

        mViewCellsStorage = stype;
}


void ViewCellsTree::CompressViewCellsPvs(ViewCell *root)
{
        if (!root->IsLeaf())
        {
                ViewCellInterior *interior = dynamic_cast<ViewCellInterior *>(root);

        ViewCellContainer::const_iterator it, it_end = interior->mChildren.end();
                
                // compress child sets first
                for (it = interior->mChildren.begin(); it != it_end; ++ it)
                {
                        CompressViewCellsPvs(*it);
                }

                // compress root node
                PullUpVisibility(interior);
        }
}


void ViewCellsTree::ExportStats(const string &mergeStats)
{
        TraversalQueue tqueue;

        tqueue.push(mRoot);
        int numViewCells = 1;
        
        const AxisAlignedBox3 box = mViewCellsManager->GetViewSpaceBox();
        const float vol = box.GetVolume();

        const int rootPvs = GetPvsSize(mRoot);

        Debug << "vsb volume: " << vol << endl;
        Debug << "root volume: " << mRoot->GetVolume() << endl;
        Debug << "root pvs: " << rootPvs << endl;

        int totalPvs;
        float totalRenderCost, avgRenderCost, expectedCost;

        float deviation = 0;
        totalPvs = rootPvs;
        totalRenderCost = avgRenderCost = expectedCost = (float)rootPvs;

        ofstream stats;
        stats.open(mergeStats.c_str());

        //-- first view cell
        stats 
                << "#Pass\n" << 0 << endl
                //<< "#Merged\n" << mergeStats.merged << endl 
                << "#ViewCells\n" << numViewCells << endl 
        << "#RenderCostDecrease\n" << 0 << endl // TODO
                << "#TotalRenderCost\n" << totalRenderCost << endl
                << "#CurrentPvs\n" << rootPvs << endl
                << "#ExpectedCost\n" << expectedCost << endl
                << "#AvgRenderCost\n" << avgRenderCost << endl
                << "#Deviation\n" << deviation << endl
                << "#TotalPvs\n" << totalPvs << endl
                << "#PvsSizeDecrease\n" << 0 << endl // TODO
                << "#Volume\n" << mRoot->GetVolume() << endl
                << endl;


        while (!tqueue.empty())
        {
                ViewCell *vc = tqueue.top();
                tqueue.pop();

                if (!vc->IsLeaf()) 
                {       
                        ViewCellInterior *interior = dynamic_cast<ViewCellInterior *>(vc);
                        
                        const int parentPvs = GetPvsSize(interior);
                        const float parentCost = (float)parentPvs * interior->GetVolume();
                        float childCost = 0;
                        int childPvs = 0;

                        -- numViewCells;

                        ViewCellContainer::const_iterator it, it_end = interior->mChildren.end();

                        for (it = interior->mChildren.begin(); it != it_end; ++ it)
                        {
                                int pvsSize = GetPvsSize(*it);
                                childCost += (float) pvsSize * (*it)->GetVolume();
                                childPvs += pvsSize;

                                tqueue.push(*it);
                                ++ numViewCells;
                        }

                
                        const float costDecr = (parentCost - childCost) / vol;

                        totalRenderCost -= costDecr;
                        totalPvs += childPvs - parentPvs;
                        
                        expectedCost = totalRenderCost / (float)numViewCells;
                        avgRenderCost = (float)totalPvs / (float)numViewCells;

                        stats 
                                << "#Pass\n" << 0 << endl
                                //<< "#Merged\n" << mergeStats.merged << endl 
                                << "#ViewCells\n" << numViewCells << endl 
                                << "#RenderCostDecrease\n" << costDecr << endl // TODO
                                << "#TotalRenderCost\n" << totalRenderCost << endl
                                << "#CurrentPvs\n" << parentPvs << endl
                                << "#ExpectedCost\n" << expectedCost << endl
                                << "#AvgRenderCost\n" << avgRenderCost << endl
                                << "#Deviation\n" << deviation << endl
                                << "#TotalPvs\n" << totalPvs << endl
                                << "#PvsSizeDecrease\n" << childPvs - parentPvs << endl // TODO
                                << "#Volume\n" << vc->GetVolume() << endl
                                << endl;

                }
        }

        stats.close();
}


void ViewCellsTree::SetRoot(ViewCell *root)
{
        mRoot = root;
}


void ViewCellsTree::CollectBestViewCellSet(ViewCellContainer &viewCells, 
                                                                                   const int numViewCells)
{
        TraversalQueue tqueue;
        tqueue.push(mRoot);
        
        while (!tqueue.empty())
        {
                ViewCell *vc = tqueue.top();
                tqueue.pop();

                // save the view cells if it is a leaf or if enough view cells have already been traversed
                // because of the priority queue, this will be the optimal set of v
                if (vc->IsLeaf() || ((viewCells.size() + tqueue.size() + 1) >= numViewCells))
                {
                        viewCells.push_back(vc);
                }
                else 
                {       
                        ViewCellInterior *interior = dynamic_cast<ViewCellInterior *>(vc);

                        ViewCellContainer::const_iterator it, it_end = interior->mChildren.end();

                        for (it = interior->mChildren.begin(); it != it_end; ++ it)
                        {
                                tqueue.push(*it);
                        }
                }
        }
}
        

void ViewCellsTree::PullUpVisibility(ViewCellInterior *interior)
{
        Intersectable::NewMail((int)interior->mChildren.size());

        ViewCellContainer::const_iterator cit, cit_end = interior->mChildren.end();

        ObjectPvsMap::const_iterator oit;

        // mail all objects in the leaf sets
        // we are interested in the objects which are present in all leaves
        // => count how often an object is part of a child set
        for (cit = interior->mChildren.begin(); cit != cit_end; ++ cit)
        {
                ViewCell *vc = *cit;

                ObjectPvsMap::const_iterator oit_end = vc->GetPvs().mEntries.end();

                for (oit = vc->GetPvs().mEntries.begin(); oit != oit_end; ++ oit)
                {
                        Intersectable *obj = (*oit).first;
                        if ((cit == interior->mChildren.begin()) && !obj->Mailed())
                                obj->Mail();
                        
                        int incm = obj->IncMail();
                }
        }

        interior->GetPvs().mEntries.clear();
        
        
        // only the objects which are present in all leaf pvs 
        // should remain in the parent pvs
        // these are the objects which have been mailed in all children
        for (cit = interior->mChildren.begin(); cit != cit_end; ++ cit)
        {
                ViewCell *vc = *cit;

                ObjectPvsMap::const_iterator oit_end = vc->GetPvs().mEntries.end();

                for (oit = vc->GetPvs().mEntries.begin(); oit != oit_end; ++ oit)
                {               
                        if ((*oit).first->Mailed((int)interior->mChildren.size()))
                        {       
                                interior->GetPvs().AddSample((*oit).first, (*oit).second.mSumPdf);
                        }
                }
        }



        // delete all the objects from the leaf sets which were moved to parent pvs
        ObjectPvsMap::const_iterator oit_end = interior->GetPvs().mEntries.end();

        for (oit = interior->GetPvs().mEntries.begin(); oit != oit_end; ++ oit)
        {
                for (cit = interior->mChildren.begin(); cit != cit_end; ++ cit)
                {
                        if (!(*cit)->GetPvs().RemoveSample((*oit).first, Limits::Infinity))
                                Debug << "should not come here!" << endl;
                }
        }

        /*int dummy = interior->GetPvs().GetSize();

        for (cit = interior->mChildren.begin(); cit != cit_end; ++ cit)
        {
                dummy += (*cit)->GetPvs().GetSize();
        }*/

}

// TODO
void ViewCellsTree::GetPvs(ViewCell *vc, ObjectPvs &pvs) const
{
        // pvs is stored in each cell
        if (mViewCellsStorage == PVS_IN_INTERIORS)
        {
                pvs = vc->GetPvs();
                return;
        }

        Intersectable::NewMail();


        int pvsSize = 0;
        ViewCell *root = vc;
        
        // add pvs from this view cell to root
        while (root->GetParent())
        {
                root = root->GetParent();
                pvs.AddPvs(root->GetPvs());
        }

        // add pvs from leaves
        stack<ViewCell *> tstack;
        tstack.push(vc);

        while (!tstack.empty())
        {
                vc = tstack.top();
                tstack.pop();

                // add new pvs
                pvs.AddPvs(vc->GetPvs());

                if (!vc->IsLeaf())
                {
                        ViewCellInterior *interior = dynamic_cast<ViewCellInterior *>(vc);

                        ViewCellContainer::const_iterator it, it_end = interior->mChildren.end();

                        for (it = interior->mChildren.begin(); it != it_end; ++ it)
                        {
                                tstack.push(*it);
                        }               
                }
        }
}


int ViewCellsTree::GetPvsSize(ViewCell *vc) const
{
        int pvsSize = 0;

        Intersectable::NewMail();

        //////////////////////////
        // for interiors, pvs can be stored using different methods

        switch (mViewCellsStorage)
        {
        case PVS_IN_LEAVES: //-- store pvs only in leaves
                {                       
                        // pvs is always stored directly in leaves
                        if (vc->IsLeaf())
                        {
                                pvsSize = vc->GetPvs().GetSize();
                                break;
                        }
        
                        // the stored pvs size is the valid pvs size => just return scalar
                        if (vc->mPvsSizeValid)
                        {
                                pvsSize = vc->mPvsSize;
                                break;
                        }
        
                        //-- if no valid pvs size stored as a scalar => compute new pvs size
                        ViewCellContainer leaves;
                        CollectLeaves(vc, leaves);

                        ViewCellContainer::const_iterator it, it_end = leaves.end();

                        Intersectable::NewMail();

                        // sum different intersectables
                        for (it = leaves.begin(); it != it_end; ++ it)
                        {
                                ObjectPvsMap::iterator oit, oit_end = (*it)->GetPvs().mEntries.end();

                                // mail all from first pvs
                                for (oit = (*it)->GetPvs().mEntries.begin(); oit != oit_end; ++ oit)
                                {
                                        Intersectable *intersect = (*oit).first;

                                        if (!intersect->Mailed())
                                        {
                                                ++ pvsSize;
                                                intersect->Mail();                                      
                                        }
                                }
                        }

                        break;
                }
        case COMPRESSED:
                {
                        ////////////////////////
                        //-- compressed pvs

                        if (vc->mPvsSizeValid)
                        {
                                return vc->mPvsSize;
                        }

                        // if no pvs size stored, compute new one
                        int pvsSize = 0;
                        ViewCell *root = vc;
        
                        // also add pvs from this view cell to root
                        while (root->GetParent())
                        {
                                root = root->GetParent();
                                pvsSize += CountDiffPvs(root);
                        }

                        stack<ViewCell *> tstack;
                        tstack.push(vc);

                        while (!tstack.empty())
                        {
                                vc = tstack.top();
                                tstack.pop();
        
                                pvsSize += CountDiffPvs(vc);

                                if (!vc->IsLeaf())
                                {
                                        ViewCellInterior *interior = dynamic_cast<ViewCellInterior *>(vc);

                                        ViewCellContainer::const_iterator it, it_end = interior->mChildren.end();

                                        for (it = interior->mChildren.begin(); it != it_end; ++ it)
                                        {
                                                tstack.push(*it);
                                        }               
                                }
                        }
                        break;
                }
        case PVS_IN_INTERIORS:
        default:
                Debug << "in interiors: " << vc->mPvsSize << " $$ " << vc->GetPvs().GetSize() << endl;
                pvsSize = vc->GetPvs().GetSize();               
        }

        return pvsSize;  
}


float ViewCellsTree::GetMemoryCost(ViewCell *vc) const
{
        const float entrySize = 
                sizeof(PvsData<Intersectable *>) + sizeof(Intersectable *);

        return (float)GetNumPvsEntries(vc) * entrySize;
}


int ViewCellsTree::GetNumPvsEntries(ViewCell *vc) const
{
        int pvsSize = 0;
        // only count leaves for uncompressed method for fairness
        if ((mViewCellsStorage == PVS_IN_INTERIORS) || vc->IsLeaf())
        {
                pvsSize = vc->GetPvs().GetSize();
        }

        if (!vc->IsLeaf())
        {
                ViewCellInterior *interior = dynamic_cast<ViewCellInterior *>(vc);

                ViewCellContainer::const_iterator it, it_end = interior->mChildren.end();

                for (it = interior->mChildren.begin(); it != it_end; ++ it)
                {
                        pvsSize += GetNumPvsEntries(*it);
                }
        }

        return pvsSize;         
}


int ViewCellsTree::ViewCellsStorage() const
{
        return mViewCellsStorage;
}


ViewCell *ViewCellsTree::GetActiveViewCell(ViewCellLeaf *vc) const
{
        return vc->GetActiveViewCell();
}


void ViewCellsTree::PropagatePvs(ViewCell *vc)
{       
        ViewCell *viewCell = vc;

        // propagate pvs up
        while (viewCell->GetParent())
        {
                viewCell->GetParent()->GetPvs().Merge(vc->GetPvs());
                viewCell = viewCell->GetParent();
        }

        if (vc->IsLeaf())
                return;

        // propagate pvs to the leaves
        stack<ViewCell *> tstack;
        tstack.push(vc);

        while (!tstack.empty())
        {
                ViewCell *viewCell = tstack.top();
                tstack.pop();

                if (viewCell != vc)
                {
                        viewCell->GetPvs().Merge(vc->GetPvs());
                }

                if (!viewCell->IsLeaf())
                {
                        ViewCellInterior *interior = dynamic_cast<ViewCellInterior *>(viewCell);

                        ViewCellContainer::const_iterator it, it_end = interior->mChildren.end();

                        for (it = interior->mChildren.begin(); it != it_end; ++ it)
                        {
                                tstack.push(*it);
                        }
                }
        }
}


void ViewCellsTree::AssignRandomColors()
{
        TraversalQueue tqueue;
        
        tqueue.push(mRoot);
        
        mRoot->SetColor(RandomColor(0.3f, 1.0f));
        
        while (!tqueue.empty())
        {
                ViewCell *vc = tqueue.top();
                tqueue.pop();

                // save the view cells if it is a leaf or if enough view cells have already been traversed
                // because of the priority queue, this will be the optimal set of v
                if (!vc->IsLeaf()) 
                {       
                        ViewCellInterior *interior = dynamic_cast<ViewCellInterior *>(vc);
                  
                        ViewCellContainer::const_iterator it, it_end = interior->mChildren.end();
                  
                        float maxProbability = -1.0f;
                  
                        ViewCell *maxViewCell = NULL;
                        for (it = interior->mChildren.begin(); it != it_end; ++ it) 
                        {
                                ViewCell *v = *it;
                          
                                // set random color
                                v->SetColor(RandomColor(0.3f, 1.0f));
                          
                                if (v->GetVolume() > maxProbability) 
                                {
                                        maxProbability = v->GetVolume();
                                        maxViewCell = v;
                                }

                                if (maxViewCell)
                                {
                                        maxViewCell->SetColor(vc->GetColor());
                                }
                                
                                tqueue.push(v);
                        }
                }        
        }
}


/** Get costs resulting from each merge step. */
void
ViewCellsTree::GetCostFunction(vector<float> &costFunction)
{
  TraversalQueue tqueue;
  tqueue.push(mRoot);
  while (!tqueue.empty()) {
        ViewCell *vc = tqueue.top();
        tqueue.pop();
        // save the view cells if it is a leaf or if enough view cells have already been traversed
        // because of the priority queue, this will be the optimal set of v
        if (!vc->IsLeaf()) {    
          ViewCellInterior *interior = dynamic_cast<ViewCellInterior *>(vc);
          costFunction.push_back(interior->GetMergeCost());
          
          ViewCellContainer::const_iterator it, it_end = interior->mChildren.end();
          
          for (it = interior->mChildren.begin(); it != it_end; ++ it) {
                tqueue.push(*it);
          }
          
        }
  }
}


void  ViewCellsTree::UpdateViewCellsStats(ViewCell *vc, ViewCellsStatistics &vcStat)
{
        ++ vcStat.viewCells;
                
        const int pvsSize = GetPvsSize(vc);

        vcStat.pvsSize += pvsSize;

        if (pvsSize == 0)
                ++ vcStat.emptyPvs;

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

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

        if (!vc->GetValid())
                ++ vcStat.invalid;
}

#if ZIPPED_VIEWCELLS
bool ViewCellsTree::Export(ogzstream &stream, const bool exportPvs)
#else
bool ViewCellsTree::Export(ofstream &stream, const bool exportPvs)
#endif
{
        // export recursivly all view cells from the root
        ExportViewCell(mRoot, stream, exportPvs);

        return true;
}


void ViewCellsTree::CreateUniqueViewCellsIds()
{
        stack<ViewCell *> tstack;

        int currentId = 0;

        tstack.push(mRoot);

        while (!tstack.empty())
        {
                ViewCell *vc = tstack.top();
                tstack.pop();

                vc->SetId(currentId ++);

                if (!vc->IsLeaf())
                {
                        ViewCellInterior *interior = dynamic_cast<ViewCellInterior *>(vc);
                        ViewCellContainer::const_iterator it, it_end = interior->mChildren.end();
                        for (it = interior->mChildren.begin(); it != it_end; ++ it)
                        {
                                tstack.push(*it);
                        }
                }
        }
}

#if ZIPPED_VIEWCELLS
void ViewCellsTree::ExportPvs(ViewCell *viewCell, ogzstream &stream)
#else
void ViewCellsTree::ExportPvs(ViewCell *viewCell, ofstream &stream)
#endif
{
        ObjectPvsMap::iterator it, it_end = viewCell->GetPvs().mEntries.end();

        for (it = viewCell->GetPvs().mEntries.begin(); it != it_end; ++ it)
        {
                stream << (*it).first->GetId() << " ";
        }
}

#if ZIPPED_VIEWCELLS
void ViewCellsTree::ExportViewCell(ViewCell *viewCell, ogzstream &stream, const bool exportPvs)
#else
void ViewCellsTree::ExportViewCell(ViewCell *viewCell, ofstream &stream, const bool exportPvs)
#endif
{
        if (viewCell->IsLeaf())
        {
                stream << "<Leaf ";
                stream << "id=\"" << viewCell->GetId() << "\" ";
                stream << "active=\"" << dynamic_cast<ViewCellLeaf *>(viewCell)->GetActiveViewCell()->GetId() << "\" ";
                stream << "mergecost=\"" << viewCell->GetMergeCost() << "\" ";
                stream << "pvs=\"";
                
                //-- export pvs
                if (exportPvs)
                {
                        ExportPvs(viewCell, stream);
                }

                stream << "\" />" << endl;
                //stream << " </Leaf>" << endl;
        }
        else
        {
                ViewCellInterior *interior = dynamic_cast<ViewCellInterior *>(viewCell);
        
                stream << "<Interior ";
                stream << "id=\"" << viewCell->GetId() << "\" ";
                stream << "mergecost=\"" << viewCell->GetMergeCost() << "\" ";
                stream << "pvs=\"";

                //-- NOTE: do not export pvss for interior view cells because
                // they can be completely reconstructed from the leaf pvss
                if (0)
                        ExportPvs(viewCell, stream);
                
                stream << "\" >" << endl;

                //-- recursivly export child view cells
                ViewCellContainer::const_iterator it, it_end = interior->mChildren.end();

                for (it = interior->mChildren.begin(); it != it_end; ++ it)
                {
                        ExportViewCell(*it, stream, exportPvs);
                }

                stream << "</Interior>" << endl;
        }
}


void ViewCellsTree::ResetPvs()
{
        stack<ViewCell *> tstack;

        tstack.push(mRoot);

        while (!tstack.empty())
        {
                ViewCell *vc = tstack.top();
                tstack.pop();

                vc->GetPvs().mEntries.clear();
                
                if (!vc->IsLeaf())
                {
                        ViewCellInterior *interior = dynamic_cast<ViewCellInterior *>(vc);
                        ViewCellContainer::const_iterator it, it_end = interior->mChildren.end();

                        for (it = interior->mChildren.begin(); it != it_end; ++ it)
                        {
                                tstack.push(*it);
                        }
                }
        }
}


void ViewCellsTree::SetActiveSetToLeaves()
{
        // todo
}

/**************************************************************************/
/*                     MergeCandidate implementation                      */
/**************************************************************************/


MergeCandidate::MergeCandidate(ViewCell *l, ViewCell *r):
mRenderCost(0),
mDeviationIncr(0),
mLeftViewCell(l),
mRightViewCell(r),
mInitialLeftViewCell(l),
mInitialRightViewCell(r)
{
        //EvalMergeCost();
}


void MergeCandidate::SetRightViewCell(ViewCell *v)
{
        mRightViewCell = v;
}


void MergeCandidate::SetLeftViewCell(ViewCell *v)
{
        mLeftViewCell = v;
}


ViewCell *MergeCandidate::GetRightViewCell() const
{
        return mRightViewCell;
}


ViewCell *MergeCandidate::GetLeftViewCell() const
{
        return mLeftViewCell;
}


ViewCell *MergeCandidate::GetInitialRightViewCell() const
{
        return mInitialRightViewCell;
}


ViewCell *MergeCandidate::GetInitialLeftViewCell() const
{
        return mInitialLeftViewCell;
}


bool MergeCandidate::IsValid() const
{
        // if one has a parent, it was already merged
        return !(mLeftViewCell->GetParent() || mRightViewCell->GetParent());
}


float MergeCandidate::GetRenderCost() const
{
        return mRenderCost;
}


float MergeCandidate::GetDeviationIncr() const
{
        return mDeviationIncr;
}


float MergeCandidate::GetMergeCost() const
{
        return mRenderCost * sRenderCostWeight + 
                   mDeviationIncr * (1.0f - sRenderCostWeight);
}






/************************************************************************/
/*                    MergeStatistics implementation                    */
/************************************************************************/


void MergeStatistics::Print(ostream &app) const
{
        app << "===== Merge statistics ===============\n";

        app << setprecision(4);

        app << "#N_CTIME ( Overall time [s] )\n" << Time() << " \n";

        app << "#N_CCTIME ( Collect candidates time [s] )\n" << collectTime * 1e-3f << " \n";

        app << "#N_MTIME ( Merge time [s] )\n" << mergeTime * 1e-3f << " \n";

        app << "#N_NODES ( Number of nodes before merge )\n" << nodes << "\n";

        app << "#N_CANDIDATES ( Number of merge candidates )\n" << candidates << "\n";

        app << "#N_MERGEDSIBLINGS ( Number of merged siblings )\n" << siblings << "\n";

        app << "#OVERALLCOST ( overall merge cost )\n" << overallCost << "\n";

        app << "#N_MERGEDNODES ( Number of merged nodes )\n" << merged << "\n";

        app << "#MAX_TREEDIST ( Maximal distance in tree of merged leaves )\n" << maxTreeDist << "\n";

        app << "#AVG_TREEDIST ( Average distance in tree of merged leaves )\n" << AvgTreeDist() << "\n";

        app << "#EXPECTEDCOST ( expected render cost )\n" << expectedRenderCost << "\n";

        app << "#DEVIATION ( deviation )\n" << deviation << "\n";

        app << "#HEURISTICS ( heuristics )\n" << heuristics << "\n";
        

        app << "===== END OF BspTree statistics ==========\n";
}

}