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


#include "ViewCell.h"
#include "Mesh.h"
#include "Intersectable.h"
#include "KdTree.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::CountKdPvs(const ViewCellLeaf *vc) const
{
	ObjectPvsMap::const_iterator oit, oit_end = vc->GetPvs().mEntries.end();

	int pvsSize = 0;

	Intersectable::NewMail();

	for (oit = vc->GetPvs().mEntries.begin(); oit != oit_end; ++ oit)
	{
		Intersectable *obj = (*oit).first;

		set<KdLeaf *>::const_iterator kit, kit_end = obj->mKdLeaves.end();

		for (kit = obj->mKdLeaves.begin(); kit != kit_end; ++ kit)
		{
			KdLeaf *leaf = *kit;

			pvsSize += ((int)leaf->mObjects.size() - (int)leaf->mMultipleObjects.size());

			//-- handle objects of several kd leaves separately
			ObjectContainer::const_iterator oit, oit_end = leaf->mMultipleObjects.end();

			for (oit = leaf->mMultipleObjects.begin(); oit != oit_end; ++ oit)
			{
				Intersectable *object = *oit;

				// object not previously in left pvs
				if (!object->Mailed())
				{
					object->Mail();
					++ pvsSize;
				}
			}		
		}
	}

	return pvsSize;
}


int ViewCellsTree::GetPvsSize(ViewCell *vc, const bool countKdPvs) 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())
			{
				if (countKdPvs)
					pvsSize = vc->GetPvs().GetSize();
				else
					pvsSize = CountKdPvs(dynamic_cast<ViewCellLeaf *>(vc));
				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";
}

}