Context Navigation

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

Revision 607, 72.5 KB checked in by mattausch, 19 years ago (diff)
added method for associating spatial hierarchy leaf with view cell

Rev	Line
[478]	1	#include <stack>
	2	#include <time.h>
	3	#include <iomanip>
	4
[463]	5	#include "Plane3.h"
	6	#include "VspBspTree.h"
	7	#include "Mesh.h"
	8	#include "common.h"
	9	#include "ViewCell.h"
	10	#include "Environment.h"
	11	#include "Polygon3.h"
	12	#include "Ray.h"
	13	#include "AxisAlignedBox3.h"
	14	#include "Exporter.h"
	15	#include "Plane3.h"
	16	#include "ViewCellBsp.h"
[478]	17	#include "ViewCellsManager.h"
[532]	18	#include "Beam.h"
[463]	19
[590]	20
[463]	21	//-- static members
[508]	22
[482]	23	/** Evaluates split plane classification with respect to the plane's
[463]	24	contribution for a minimum number of ray splits.
	25	*/
	26	const float VspBspTree::sLeastRaySplitsTable[] = {0, 0, 1, 1, 0};
[482]	27	/** Evaluates split plane classification with respect to the plane's
[580]	28	contribution for balanced rays.
[463]	29	*/
	30	const float VspBspTree::sBalancedRaysTable[] = {1, -1, 0, 0, 0};
	31
	32
[482]	33	int VspBspTree::sFrontId = 0;
[463]	34	int VspBspTree::sBackId = 0;
	35	int VspBspTree::sFrontAndBackId = 0;
	36
	37
[557]	38
[578]	39
[579]	40	// pvs penalty can be different from pvs size
	41	inline float EvalPvsPenalty(const int pvs,
	42	const int lower,
	43	const int upper)
	44	{
	45	// clamp to minmax values
	46	if (pvs < lower)
	47	return (float)lower;
	48	if (pvs > upper)
	49	return (float)upper;
[578]	50
[579]	51	return (float)pvs;
	52	}
[463]	53
[557]	54
[579]	55
	56
[601]	57	/******************************************************************************/
	58	/* class VspBspTree implementation */
	59	/******************************************************************************/
[579]	60
	61
[482]	62	VspBspTree::VspBspTree():
[463]	63	mRoot(NULL),
[547]	64	mUseAreaForPvs(false),
[478]	65	mCostNormalizer(Limits::Small),
	66	mViewCellsManager(NULL),
[497]	67	mOutOfBoundsCell(NULL),
[580]	68	mStoreRays(false),
[605]	69	mRenderCostWeight(0.5),
	70	mUseRandomAxis(false)
[463]	71	{
[486]	72	bool randomize = false;
	73	environment->GetBoolValue("VspBspTree.Construction.randomize", randomize);
	74	if (randomize)
	75	Randomize(); // initialise random generator for heuristics
[463]	76
	77	//-- termination criteria for autopartition
	78	environment->GetIntValue("VspBspTree.Termination.maxDepth", mTermMaxDepth);
	79	environment->GetIntValue("VspBspTree.Termination.minPvs", mTermMinPvs);
	80	environment->GetIntValue("VspBspTree.Termination.minRays", mTermMinRays);
[547]	81	environment->GetFloatValue("VspBspTree.Termination.minProbability", mTermMinProbability);
[463]	82	environment->GetFloatValue("VspBspTree.Termination.maxRayContribution", mTermMaxRayContribution);
	83	environment->GetFloatValue("VspBspTree.Termination.minAccRayLenght", mTermMinAccRayLength);
[472]	84	environment->GetFloatValue("VspBspTree.Termination.maxCostRatio", mTermMaxCostRatio);
	85	environment->GetIntValue("VspBspTree.Termination.missTolerance", mTermMissTolerance);
[564]	86	environment->GetIntValue("VspBspTree.Termination.maxViewCells", mMaxViewCells);
	87	//-- max cost ratio for early tree termination
	88	environment->GetFloatValue("VspBspTree.Termination.maxCostRatio", mTermMaxCostRatio);
[600]	89	mTermMinPolygons = 25;
[463]	90	//-- factors for bsp tree split plane heuristics
	91	environment->GetFloatValue("VspBspTree.Factor.balancedRays", mBalancedRaysFactor);
	92	environment->GetFloatValue("VspBspTree.Factor.pvs", mPvsFactor);
	93	environment->GetFloatValue("VspBspTree.Termination.ct_div_ci", mCtDivCi);
	94
[482]	95
[463]	96	//-- partition criteria
	97	environment->GetIntValue("VspBspTree.maxPolyCandidates", mMaxPolyCandidates);
	98	environment->GetIntValue("VspBspTree.maxRayCandidates", mMaxRayCandidates);
	99	environment->GetIntValue("VspBspTree.splitPlaneStrategy", mSplitPlaneStrategy);
	100
	101	environment->GetFloatValue("VspBspTree.Construction.epsilon", mEpsilon);
	102	environment->GetIntValue("VspBspTree.maxTests", mMaxTests);
	103
[564]	104	// if only the driving axis is used for axis aligned split
[551]	105	environment->GetBoolValue("VspBspTree.splitUseOnlyDrivingAxis", mOnlyDrivingAxis);
[580]	106
[508]	107	//-- termination criteria for axis aligned split
	108	environment->GetFloatValue("VspBspTree.Termination.AxisAligned.maxRayContribution",
	109	mTermMaxRayContriForAxisAligned);
	110	environment->GetIntValue("VspBspTree.Termination.AxisAligned.minRays",
	111	mTermMinRaysForAxisAligned);
[487]	112
[508]	113	//environment->GetFloatValue("VspBspTree.maxTotalMemory", mMaxTotalMemory);
	114	environment->GetFloatValue("VspBspTree.maxStaticMemory", mMaxMemory);
	115
[580]	116	environment->GetFloatValue("VspBspTree.Construction.renderCostWeight", mRenderCostWeight);
[564]	117
[607]	118	environment->GetBoolValue("VspBspTree.usePolygonSplitIfAvailable", mUsePolygonSplitIfAvailable);
[547]	119
[600]	120	mSubdivisionStats.open("subdivisionStats.log");
[580]	121
	122
[478]	123	//-- debug output
[580]	124
[473]	125	Debug << "***** VSP BSP options ****** " << endl;
	126	Debug << "max depth: " << mTermMaxDepth << endl;
	127	Debug << "min PVS: " << mTermMinPvs << endl;
[547]	128	Debug << "min probabiliy: " << mTermMinProbability << endl;
[473]	129	Debug << "min rays: " << mTermMinRays << endl;
	130	Debug << "max ray contri: " << mTermMaxRayContribution << endl;
	131	//Debug << "VSP BSP mininam accumulated ray lenght: ", mTermMinAccRayLength) << endl;
	132	Debug << "max cost ratio: " << mTermMaxCostRatio << endl;
	133	Debug << "miss tolerance: " << mTermMissTolerance << endl;
	134	Debug << "max view cells: " << mMaxViewCells << endl;
	135	Debug << "max polygon candidates: " << mMaxPolyCandidates << endl;
	136	Debug << "max plane candidates: " << mMaxRayCandidates << endl;
[486]	137	Debug << "randomize: " << randomize << endl;
[582]	138
[551]	139	Debug << "using area for pvs: " << mUseAreaForPvs << endl;
[580]	140	Debug << "render cost weight: " << mRenderCostWeight << endl;
[463]	141	Debug << "Split plane strategy: ";
[564]	142
[463]	143	if (mSplitPlaneStrategy & RANDOM_POLYGON)
[474]	144	{
[463]	145	Debug << "random polygon ";
[474]	146	}
[463]	147	if (mSplitPlaneStrategy & AXIS_ALIGNED)
[472]	148	{
[463]	149	Debug << "axis aligned ";
[472]	150	}
[463]	151	if (mSplitPlaneStrategy & LEAST_RAY_SPLITS)
[472]	152	{
[474]	153	mCostNormalizer += mLeastRaySplitsFactor;
[463]	154	Debug << "least ray splits ";
[472]	155	}
[463]	156	if (mSplitPlaneStrategy & BALANCED_RAYS)
[472]	157	{
[474]	158	mCostNormalizer += mBalancedRaysFactor;
[463]	159	Debug << "balanced rays ";
[472]	160	}
[463]	161	if (mSplitPlaneStrategy & PVS)
[472]	162	{
[474]	163	mCostNormalizer += mPvsFactor;
[463]	164	Debug << "pvs";
[472]	165	}
[482]	166
[489]	167
[480]	168	mSplitCandidates = new vector<SortableEntry>;
[463]	169
	170	Debug << endl;
	171	}
	172
[580]	173
[508]	174	BspViewCell *VspBspTree::GetOutOfBoundsCell()
	175	{
	176	return mOutOfBoundsCell;
	177	}
[463]	178
[508]	179
[489]	180	BspViewCell *VspBspTree::GetOrCreateOutOfBoundsCell()
	181	{
	182	if (!mOutOfBoundsCell)
[508]	183	{
[489]	184	mOutOfBoundsCell = new BspViewCell();
[508]	185	mOutOfBoundsCell->SetId(-1);
[547]	186	mOutOfBoundsCell->SetValid(false);
[508]	187	}
[547]	188
[489]	189	return mOutOfBoundsCell;
	190	}
	191
	192
[482]	193	const BspTreeStatistics &VspBspTree::GetStatistics() const
[463]	194	{
[574]	195	return mBspStats;
[463]	196	}
	197
	198
	199	VspBspTree::~VspBspTree()
	200	{
	201	DEL_PTR(mRoot);
[480]	202	DEL_PTR(mSplitCandidates);
[463]	203	}
	204
[579]	205
[482]	206	int VspBspTree::AddMeshToPolygons(Mesh *mesh,
	207	PolygonContainer &polys,
[463]	208	MeshInstance *parent)
	209	{
	210	FaceContainer::const_iterator fi;
[482]	211
[463]	212	// copy the face data to polygons
	213	for (fi = mesh->mFaces.begin(); fi != mesh->mFaces.end(); ++ fi)
	214	{
	215	Polygon3 poly = new Polygon3((fi), mesh);
[482]	216
[463]	217	if (poly->Valid(mEpsilon))
	218	{
	219	poly->mParent = parent; // set parent intersectable
	220	polys.push_back(poly);
	221	}
	222	else
	223	DEL_PTR(poly);
	224	}
	225	return (int)mesh->mFaces.size();
	226	}
	227
[580]	228
[482]	229	int VspBspTree::AddToPolygonSoup(const ViewCellContainer &viewCells,
	230	PolygonContainer &polys,
[463]	231	int maxObjects)
	232	{
[482]	233	int limit = (maxObjects > 0) ?
[463]	234	Min((int)viewCells.size(), maxObjects) : (int)viewCells.size();
[482]	235
[463]	236	int polysSize = 0;
	237
	238	for (int i = 0; i < limit; ++ i)
	239	{
	240	if (viewCells[i]->GetMesh()) // copy the mesh data to polygons
	241	{
	242	mBox.Include(viewCells[i]->GetBox()); // add to BSP tree aabb
[482]	243	polysSize +=
	244	AddMeshToPolygons(viewCells[i]->GetMesh(),
	245	polys,
[463]	246	viewCells[i]);
	247	}
	248	}
	249
	250	return polysSize;
	251	}
	252
[580]	253
[482]	254	int VspBspTree::AddToPolygonSoup(const ObjectContainer &objects,
	255	PolygonContainer &polys,
[463]	256	int maxObjects)
	257	{
[482]	258	int limit = (maxObjects > 0) ?
[463]	259	Min((int)objects.size(), maxObjects) : (int)objects.size();
[482]	260
[463]	261	for (int i = 0; i < limit; ++i)
	262	{
	263	Intersectable object = objects[i];//it;
	264	Mesh *mesh = NULL;
	265
	266	switch (object->Type()) // extract the meshes
	267	{
	268	case Intersectable::MESH_INSTANCE:
	269	mesh = dynamic_cast<MeshInstance *>(object)->GetMesh();
	270	break;
	271	case Intersectable::VIEW_CELL:
	272	mesh = dynamic_cast<ViewCell *>(object)->GetMesh();
	273	break;
	274	// TODO: handle transformed mesh instances
	275	default:
	276	Debug << "intersectable type not supported" << endl;
	277	break;
	278	}
[482]	279
[463]	280	if (mesh) // copy the mesh data to polygons
	281	{
	282	mBox.Include(object->GetBox()); // add to BSP tree aabb
[485]	283	AddMeshToPolygons(mesh, polys, NULL);
[463]	284	}
	285	}
	286
	287	return (int)polys.size();
	288	}
	289
[580]	290
[483]	291	void VspBspTree::Construct(const VssRayContainer &sampleRays,
	292	AxisAlignedBox3 *forcedBoundingBox)
[463]	293	{
[574]	294	mBspStats.nodes = 1;
[463]	295	mBox.Initialize(); // initialise BSP tree bounding box
[482]	296
[484]	297	if (forcedBoundingBox)
	298	mBox = *forcedBoundingBox;
	299
[463]	300	PolygonContainer polys;
	301	RayInfoContainer *rays = new RayInfoContainer();
	302
	303	VssRayContainer::const_iterator rit, rit_end = sampleRays.end();
	304
	305	long startTime = GetTime();
	306
[480]	307	cout << "Extracting polygons from rays ... ";
[463]	308
	309	Intersectable::NewMail();
	310
[511]	311	int numObj = 0;
[542]	312
[463]	313	//-- extract polygons intersected by the rays
	314	for (rit = sampleRays.begin(); rit != rit_end; ++ rit)
	315	{
	316	VssRay ray = rit;
[482]	317
[484]	318	if ((mBox.IsInside(ray->mTermination) \|\| !forcedBoundingBox) &&
	319	ray->mTerminationObject &&
	320	!ray->mTerminationObject->Mailed())
[463]	321	{
	322	ray->mTerminationObject->Mail();
	323	MeshInstance obj = dynamic_cast<MeshInstance >(ray->mTerminationObject);
	324	AddMeshToPolygons(obj->GetMesh(), polys, obj);
[511]	325	++ numObj;
[484]	326	//-- compute bounding box
	327	if (!forcedBoundingBox)
	328	mBox.Include(ray->mTermination);
[463]	329	}
	330
[484]	331	if ((mBox.IsInside(ray->mOrigin) \|\| !forcedBoundingBox) &&
	332	ray->mOriginObject &&
	333	!ray->mOriginObject->Mailed())
[463]	334	{
	335	ray->mOriginObject->Mail();
	336	MeshInstance obj = dynamic_cast<MeshInstance >(ray->mOriginObject);
	337	AddMeshToPolygons(obj->GetMesh(), polys, obj);
[511]	338	++ numObj;
	339
[484]	340	//-- compute bounding box
	341	if (!forcedBoundingBox)
	342	mBox.Include(ray->mOrigin);
[463]	343	}
	344	}
[535]	345
[547]	346	Debug << "maximal pvs (i.e., pvs still considered as valid) : " << mViewCellsManager->GetMaxPvsSize() << endl;
[580]	347
[463]	348	//-- store rays
	349	for (rit = sampleRays.begin(); rit != rit_end; ++ rit)
	350	{
	351	VssRay ray = rit;
[482]	352
[463]	353	float minT, maxT;
	354
[564]	355	static Ray hray;
	356	hray.Init(*ray);
[562]	357
[483]	358	// TODO: not very efficient to implictly cast between rays types
[564]	359	if (mBox.GetRaySegment(hray, minT, maxT))
[463]	360	{
	361	float len = ray->Length();
[482]	362
	363	if (!len)
[463]	364	len = Limits::Small;
[482]	365
[463]	366	rays->push_back(RayInfo(ray, minT / len, maxT / len));
	367	}
	368	}
	369
[587]	370	// normalize
[547]	371	if (mUseAreaForPvs)
[587]	372	mTermMinProbability *= mBox.SurfaceArea();
[547]	373	else
	374	mTermMinProbability *= mBox.GetVolume();
	375
[574]	376	mBspStats.polys = (int)polys.size();
[463]	377
[475]	378	cout << "finished" << endl;
[463]	379
[482]	380	Debug << "\nPolygon extraction: " << (int)polys.size() << " polys extracted from "
[475]	381	<< (int)sampleRays.size() << " rays in "
	382	<< TimeDiff(startTime, GetTime())*1e-3 << " secs" << endl << endl;
	383
[463]	384	Construct(polys, rays);
	385
	386	// clean up polygons
	387	CLEAR_CONTAINER(polys);
	388	}
	389
[508]	390
	391	// return memory usage in MB
[600]	392	float VspBspTree::GetMemUsage(/const VspBspTraversalQueue &tstack/) const
[508]	393	{
	394	return
	395	(sizeof(VspBspTree) +
[587]	396	(float)mBspStats.Leaves() * sizeof(BspLeaf) +
	397	// the nodes in the stack is the minimal additional number of leaves
	398	//(float)tstack.size() * sizeof(BspLeaf) +
[574]	399	mBspStats.Interior() * sizeof(BspInterior) +
	400	mBspStats.accumRays * sizeof(RayInfo)) / (1024.0f * 1024.0f);
[508]	401	}
	402
	403
	404
[463]	405	void VspBspTree::Construct(const PolygonContainer &polys, RayInfoContainer *rays)
	406	{
[600]	407	VspBspTraversalQueue tQueue;
[463]	408
	409	mRoot = new BspLeaf();
	410
	411	// constrruct root node geometry
	412	BspNodeGeometry *geom = new BspNodeGeometry();
	413	ConstructGeometry(mRoot, *geom);
	414
[547]	415	const float prop = mUseAreaForPvs ? geom->GetArea() : geom->GetVolume();
	416
[482]	417	VspBspTraversalData tData(mRoot,
	418	new PolygonContainer(polys),
[463]	419	0,
[482]	420	rays,
	421	ComputePvsSize(*rays),
[547]	422	prop,
[463]	423	geom);
	424
[578]	425	// first node is kd node, i.e. an axis aligned box
	426	if (1)
	427	tData.mIsKdNode = true;
	428	else
	429	tData.mIsKdNode = false;
[562]	430
[600]	431	tQueue.push(tData);
[463]	432
[605]	433	mTotalCost = tData.mPvs * tData.mProbability / mBox.GetVolume();
[607]	434	mTotalPvsSize = tData.mPvs;
	435
[601]	436	mSubdivisionStats
[605]	437	<< "#ViewCells\n1\n" << endl
[601]	438	<< "#RenderCostDecrease\n0\n" << endl
[607]	439	<< "#TotalRenderCost\n" << mTotalCost << endl
	440	<< "#AvgRenderCost\n" << mTotalPvsSize << endl;
[601]	441
[600]	442	Debug << "total cost: " << mTotalCost << endl;
	443
[601]	444
[574]	445	mBspStats.Start();
[542]	446	cout << "Contructing vsp bsp tree ... \n";
[463]	447
[542]	448	long startTime = GetTime();
[587]	449
	450	// used for intermediate time measurements and progress
[542]	451	long interTime = GetTime();
[547]	452
[542]	453	mOutOfMemory = false;
[482]	454
[542]	455	int nleaves = 500;
	456
[600]	457	while (!tQueue.empty())
[463]	458	{
[600]	459	tData = tQueue.top();
	460	tQueue.pop();
[463]	461
[508]	462	if (0 && !mOutOfMemory)
	463	{
	464	float mem = GetMemUsage();
[478]	465
[508]	466	if (mem > mMaxMemory)
	467	{
	468	mOutOfMemory = true;
	469	Debug << "memory limit reached: " << mem << endl;
	470	}
	471	}
	472
[587]	473	// subdivide leaf node
[600]	474	BspNode *r = Subdivide(tQueue, tData);
[463]	475
	476	if (r == mRoot)
[482]	477	Debug << "VSP BSP tree construction time spent at root: "
[542]	478	<< TimeDiff(startTime, GetTime())*1e-3 << "s" << endl;
	479
[574]	480	if (mBspStats.Leaves() >= nleaves)
[542]	481	{
	482	nleaves += 500;
	483
[574]	484	cout << "leaves=" << mBspStats.Leaves() << endl;
[542]	485	Debug << "needed "
[587]	486	<< TimeDiff(interTime, GetTime())*1e-3
	487	<< " secs to create 500 leaves" << endl;
[542]	488	interTime = GetTime();
	489	}
[463]	490	}
	491
[542]	492	Debug << "Used Memory: " << GetMemUsage() << " MB" << endl << endl;
[463]	493	cout << "finished\n";
	494
[574]	495	mBspStats.Stop();
[463]	496	}
	497
[508]	498
[478]	499	bool VspBspTree::TerminationCriteriaMet(const VspBspTraversalData &data) const
[463]	500	{
[482]	501	return
[463]	502	(((int)data.mRays->size() <= mTermMinRays) \|\|
[473]	503	(data.mPvs <= mTermMinPvs) \|\|
[547]	504	(data.mProbability <= mTermMinProbability) \|\|
[574]	505	(mBspStats.Leaves() >= mMaxViewCells) \|\|
[600]	506	#if 0
	507	(((int)data.mPolygons->size() <= mTermMinPolygons) && !data.mPolygons->empty())\|\|
	508	#endif
[535]	509	(data.GetAvgRayContribution() > mTermMaxRayContribution) \|\|
[463]	510	(data.mDepth >= mTermMaxDepth));
	511	}
	512
[508]	513
[600]	514	BspNode *VspBspTree::Subdivide(VspBspTraversalQueue &tQueue,
[463]	515	VspBspTraversalData &tData)
	516	{
[473]	517	BspNode *newNode = tData.mNode;
	518
[542]	519	if (!mOutOfMemory && !TerminationCriteriaMet(tData))
[473]	520	{
	521	PolygonContainer coincident;
[482]	522
[473]	523	VspBspTraversalData tFrontData;
	524	VspBspTraversalData tBackData;
	525
	526	// create new interior node and two leaf nodes
	527	// or return leaf as it is (if maxCostRatio missed)
	528	newNode = SubdivideNode(tData, tFrontData, tBackData, coincident);
	529
	530	if (!newNode->IsLeaf()) //-- continue subdivision
	531	{
[600]	532	if (1)
	533	{
[605]	534	float cFront = tFrontData.mPvs * tFrontData.mProbability;
	535	float cBack = tBackData.mPvs * tBackData.mProbability;
	536	float cData = tData.mPvs * tData.mProbability;;
	537
[607]	538	float costDecr =
[605]	539	(cFront + cBack - cData) / mBox.GetVolume();
[607]	540
	541	int pvsFront = tFrontData.mPvs * tFrontData.mProbability;
	542	int pvsBack = tBackData.mPvs * tBackData.mProbability;
	543	int pvsData = tData.mPvs * tData.mProbability;;
	544
[600]	545	mTotalCost += costDecr;
[607]	546	mTotalPvsSize += tFrontData.mPvs + tBackData.mPvs - tData.mPvs;
[600]	547
	548	mSubdivisionStats
[605]	549	<< "#ViewCells\n" << mBspStats.Leaves() << endl
[600]	550	<< "#RenderCostDecrease\n" << -costDecr << endl
[607]	551	<< "#TotalRenderCost\n" << mTotalCost << endl
	552	<< "#AvgRenderCost\n" << mTotalPvsSize / mBspStats.Leaves() << endl;
[600]	553	}
	554
[473]	555	// push the children on the stack
[600]	556	tQueue.push(tFrontData);
	557	tQueue.push(tBackData);
[473]	558
	559	// delete old leaf node
[482]	560	DEL_PTR(tData.mNode);
[473]	561	}
	562	}
[482]	563
[478]	564	//-- terminate traversal and create new view cell
[473]	565	if (newNode->IsLeaf())
[463]	566	{
[473]	567	BspLeaf leaf = dynamic_cast<BspLeaf >(newNode);
[547]	568	BspViewCell *viewCell = new BspViewCell();
[489]	569
[463]	570	leaf->SetViewCell(viewCell);
[487]	571
	572	//-- update pvs
[556]	573	int conSamp = 0;
	574	float sampCon = 0.0f;
	575	AddToPvs(leaf, *tData.mRays, sampCon, conSamp);
[487]	576
[574]	577	mBspStats.contributingSamples += conSamp;
	578	mBspStats.sampleContributions +=(int) sampCon;
[487]	579
	580	//-- store additional info
[478]	581	if (mStoreRays)
	582	{
	583	RayInfoContainer::const_iterator it, it_end = tData.mRays->end();
	584	for (it = tData.mRays->begin(); it != it_end; ++ it)
	585	leaf->mVssRays.push_back((*it).mRay);
	586	}
[564]	587	// should I check here?
[562]	588	if (0 && !mViewCellsManager->CheckValidity(viewCell, 0, mViewCellsManager->GetMaxPvsSize()))
[547]	589	{
	590	viewCell->SetValid(false);
	591	leaf->SetTreeValid(false);
	592	PropagateUpValidity(leaf);
[463]	593
[574]	594	++ mBspStats.invalidLeaves;
[547]	595	}
	596
[580]	597	viewCell->mLeaf = leaf;
[547]	598
	599	if (mUseAreaForPvs)
	600	viewCell->SetArea(tData.mProbability);
	601	else
	602	viewCell->SetVolume(tData.mProbability);
	603
	604	leaf->mProbability = tData.mProbability;
	605
[487]	606	EvaluateLeafStats(tData);
[463]	607	}
[482]	608
[473]	609	//-- cleanup
[478]	610	tData.Clear();
[463]	611
[472]	612	return newNode;
[463]	613	}
	614
[487]	615
[472]	616	BspNode *VspBspTree::SubdivideNode(VspBspTraversalData &tData,
	617	VspBspTraversalData &frontData,
	618	VspBspTraversalData &backData,
	619	PolygonContainer &coincident)
[463]	620	{
	621	BspLeaf leaf = dynamic_cast<BspLeaf >(tData.mNode);
[508]	622
[463]	623	// select subdivision plane
[472]	624	Plane3 splitPlane;
[542]	625
[508]	626	int maxCostMisses = tData.mMaxCostMisses;
[542]	627
	628	const bool success =
	629	SelectPlane(splitPlane, leaf, tData, frontData, backData);
	630
[508]	631	if (!success)
[472]	632	{
	633	++ maxCostMisses;
[463]	634
[480]	635	if (maxCostMisses > mTermMissTolerance)
[473]	636	{
	637	// terminate branch because of max cost
[574]	638	++ mBspStats.maxCostNodes;
[472]	639	return leaf;
[474]	640	}
[472]	641	}
	642
[574]	643	mBspStats.nodes += 2;
[473]	644
	645	//-- subdivide further
[482]	646	BspInterior *interior = new BspInterior(splitPlane);
[472]	647
[463]	648	#ifdef _DEBUG
	649	Debug << interior << endl;
	650	#endif
[482]	651
[473]	652	//-- the front and back traversal data is filled with the new values
	653	frontData.mDepth = tData.mDepth + 1;
[508]	654	frontData.mPolygons = new PolygonContainer();
[473]	655	frontData.mRays = new RayInfoContainer();
[508]	656
[473]	657	backData.mDepth = tData.mDepth + 1;
[508]	658	backData.mPolygons = new PolygonContainer();
[473]	659	backData.mRays = new RayInfoContainer();
[508]	660
[473]	661	// subdivide rays
[482]	662	SplitRays(interior->GetPlane(),
	663	*tData.mRays,
	664	*frontData.mRays,
[463]	665	*backData.mRays);
[482]	666
[473]	667	// subdivide polygons
[491]	668	SplitPolygons(interior->GetPlane(),
	669	*tData.mPolygons,
	670	*frontData.mPolygons,
	671	*backData.mPolygons,
	672	coincident);
[463]	673
[482]	674
[473]	675	// how often was max cost ratio missed in this branch?
[472]	676	frontData.mMaxCostMisses = maxCostMisses;
	677	backData.mMaxCostMisses = maxCostMisses;
	678
	679	// compute pvs
[463]	680	frontData.mPvs = ComputePvsSize(*frontData.mRays);
	681	backData.mPvs = ComputePvsSize(*backData.mRays);
	682
[508]	683	// split front and back node geometry and compute area
[547]	684
	685	// if geometry was not already computed
[602]	686	if (!frontData.mGeometry && !backData.mGeometry)
[463]	687	{
[547]	688	frontData.mGeometry = new BspNodeGeometry();
	689	backData.mGeometry = new BspNodeGeometry();
[482]	690
[547]	691	tData.mGeometry->SplitGeometry(*frontData.mGeometry,
	692	*backData.mGeometry,
	693	interior->GetPlane(),
	694	mBox,
	695	mEpsilon);
[508]	696
[547]	697	if (mUseAreaForPvs)
	698	{
	699	frontData.mProbability = frontData.mGeometry->GetArea();
	700	backData.mProbability = backData.mGeometry->GetArea();
[508]	701	}
[547]	702	else
	703	{
	704	frontData.mProbability = frontData.mGeometry->GetVolume();
[602]	705	backData.mProbability = tData.mProbability - frontData.mProbability;
[547]	706	}
[463]	707	}
[547]	708
[463]	709
	710	//-- create front and back leaf
	711
	712	BspInterior *parent = leaf->GetParent();
	713
	714	// replace a link from node's parent
[487]	715	if (parent)
[463]	716	{
	717	parent->ReplaceChildLink(leaf, interior);
	718	interior->SetParent(parent);
	719	}
	720	else // new root
	721	{
	722	mRoot = interior;
	723	}
	724
	725	// and setup child links
	726	interior->SetupChildLinks(new BspLeaf(interior), new BspLeaf(interior));
[482]	727
[463]	728	frontData.mNode = interior->GetFront();
	729	backData.mNode = interior->GetBack();
[473]	730
[463]	731	//DEL_PTR(leaf);
	732	return interior;
	733	}
	734
[508]	735
[463]	736	void VspBspTree::AddToPvs(BspLeaf *leaf,
[482]	737	const RayInfoContainer &rays,
[556]	738	float &sampleContributions,
[463]	739	int &contributingSamples)
	740	{
[556]	741	sampleContributions = 0;
	742	contributingSamples = 0;
	743
	744	RayInfoContainer::const_iterator it, it_end = rays.end();
	745
	746	ViewCell *vc = leaf->GetViewCell();
	747
	748	// add contributions from samples to the PVS
	749	for (it = rays.begin(); it != it_end; ++ it)
[463]	750	{
[556]	751	float sc = 0.0f;
	752	VssRay ray = (it).mRay;
	753	bool madeContrib = false;
	754	float contribution;
	755	if (ray->mTerminationObject) {
	756	if (vc->GetPvs().AddSample(ray->mTerminationObject, ray->mPdf, contribution))
	757	madeContrib = true;
	758	sc += contribution;
	759	}
	760
	761	if (ray->mOriginObject) {
	762	if (vc->GetPvs().AddSample(ray->mOriginObject, ray->mPdf, contribution))
	763	madeContrib = true;
	764	sc += contribution;
	765	}
	766
	767	sampleContributions += sc;
	768	if (madeContrib)
	769	++ contributingSamples;
	770
	771	//leaf->mVssRays.push_back(ray);
[463]	772	}
	773	}
	774
[580]	775
[480]	776	void VspBspTree::SortSplitCandidates(const RayInfoContainer &rays, const int axis)
[463]	777	{
[480]	778	mSplitCandidates->clear();
[463]	779
[480]	780	int requestedSize = 2 * (int)(rays.size());
	781	// creates a sorted split candidates array
	782	if (mSplitCandidates->capacity() > 500000 &&
	783	requestedSize < (int)(mSplitCandidates->capacity()/10) )
	784	{
[482]	785	delete mSplitCandidates;
[480]	786	mSplitCandidates = new vector<SortableEntry>;
	787	}
[463]	788
[480]	789	mSplitCandidates->reserve(requestedSize);
	790
	791	// insert all queries
	792	for(RayInfoContainer::const_iterator ri = rays.begin(); ri < rays.end(); ++ ri)
[473]	793	{
[480]	794	bool positive = (*ri).mRay->HasPosDir(axis);
	795	mSplitCandidates->push_back(SortableEntry(positive ? SortableEntry::ERayMin : SortableEntry::ERayMax,
[482]	796	(ri).ExtrapOrigin(axis), (ri).mRay));
[480]	797
	798	mSplitCandidates->push_back(SortableEntry(positive ? SortableEntry::ERayMax : SortableEntry::ERayMin,
[482]	799	(ri).ExtrapTermination(axis), (ri).mRay));
[473]	800	}
[480]	801
	802	stable_sort(mSplitCandidates->begin(), mSplitCandidates->end());
[463]	803	}
	804
[580]	805
[480]	806	float VspBspTree::BestCostRatioHeuristics(const RayInfoContainer &rays,
	807	const AxisAlignedBox3 &box,
	808	const int pvsSize,
[481]	809	const int &axis,
[480]	810	float &position)
[463]	811	{
[480]	812	int raysBack;
	813	int raysFront;
	814	int pvsBack;
	815	int pvsFront;
[463]	816
[480]	817	SortSplitCandidates(rays, axis);
	818
[463]	819	// go through the lists, count the number of objects left and right
	820	// and evaluate the following cost funcion:
[480]	821	// C = ct_div_ci + (qlrl + qrrr)/queries
	822
	823	int rl = 0, rr = (int)rays.size();
	824	int pl = 0, pr = pvsSize;
	825
[463]	826	float minBox = box.Min(axis);
	827	float maxBox = box.Max(axis);
[480]	828	float sizeBox = maxBox - minBox;
	829
	830	float minBand = minBox + 0.1f*(maxBox - minBox);
	831	float maxBand = minBox + 0.9f*(maxBox - minBox);
	832
	833	float sum = rr*sizeBox;
[463]	834	float minSum = 1e20f;
	835
[480]	836	Intersectable::NewMail();
	837
	838	// set all object as belonging to the fron pvs
	839	for(RayInfoContainer::const_iterator ri = rays.begin(); ri != rays.end(); ++ ri)
[463]	840	{
[480]	841	if ((*ri).mRay->IsActive())
[463]	842	{
[480]	843	Intersectable object = (ri).mRay->mTerminationObject;
	844
	845	if (object)
[482]	846	{
[480]	847	if (!object->Mailed())
	848	{
	849	object->Mail();
	850	object->mCounter = 1;
	851	}
	852	else
	853	++ object->mCounter;
[482]	854	}
[463]	855	}
[480]	856	}
	857
	858	Intersectable::NewMail();
	859
	860	for (vector<SortableEntry>::const_iterator ci = mSplitCandidates->begin();
	861	ci < mSplitCandidates->end(); ++ ci)
	862	{
	863	VssRay *ray;
	864
	865	switch ((*ci).type)
[463]	866	{
[480]	867	case SortableEntry::ERayMin:
	868	{
	869	++ rl;
[482]	870	ray = (VssRay ) (ci).ray;
	871
[480]	872	Intersectable *object = ray->mTerminationObject;
[482]	873
[480]	874	if (object && !object->Mailed())
	875	{
	876	object->Mail();
	877	++ pl;
	878	}
	879	break;
	880	}
	881	case SortableEntry::ERayMax:
	882	{
	883	-- rr;
[482]	884	ray = (VssRay ) (ci).ray;
[463]	885
[480]	886	Intersectable *object = ray->mTerminationObject;
	887
	888	if (object)
	889	{
	890	if (-- object->mCounter == 0)
[482]	891	-- pr;
[480]	892	}
	893
	894	break;
	895	}
	896	}
	897
	898	// Note: sufficient to compare size of bounding boxes of front and back side?
	899	if ((ci).value > minBand && (ci).value < maxBand)
	900	{
	901	sum = pl((ci).value - minBox) + pr(maxBox - (ci).value);
	902
	903	// cout<<"pos="<<(*ci).value<<"\t q=("<<ql<<","<<qr<<")\t r=("<<rl<<","<<rr<<")"<<endl;
	904	// cout<<"cost= "<<sum<<endl;
	905
	906	if (sum < minSum)
[463]	907	{
	908	minSum = sum;
	909	position = (*ci).value;
	910
[480]	911	raysBack = rl;
	912	raysFront = rr;
	913
	914	pvsBack = pl;
	915	pvsFront = pr;
	916
[463]	917	}
	918	}
	919	}
	920
[480]	921	float oldCost = (float)pvsSize;
	922	float newCost = mCtDivCi + minSum / sizeBox;
	923	float ratio = newCost / oldCost;
[463]	924
[480]	925	//Debug << "costRatio=" << ratio << " pos=" << position << " t=" << (position - minBox) / (maxBox - minBox)
	926	// <<"\t q=(" << queriesBack << "," << queriesFront << ")\t r=(" << raysBack << "," << raysFront << ")" << endl;
	927
	928	return ratio;
[463]	929	}
	930
[480]	931
[482]	932	float VspBspTree::SelectAxisAlignedPlane(Plane3 &plane,
[491]	933	const VspBspTraversalData &tData,
[495]	934	int &axis,
[508]	935	BspNodeGeometry **frontGeom,
	936	BspNodeGeometry **backGeom,
[547]	937	float &pFront,
	938	float &pBack,
	939	const bool useKdSplit)
[463]	940	{
[508]	941	const bool useCostHeuristics = false;
	942
	943	//-- regular split
	944	float nPosition[3];
	945	float nCostRatio[3];
[547]	946	float nProbFront[3];
	947	float nProbBack[3];
[508]	948
	949	BspNodeGeometry *nFrontGeom[3];
	950	BspNodeGeometry *nBackGeom[3];
	951
	952	int bestAxis = -1;
	953
[545]	954	// create bounding box of node geometry
[463]	955	AxisAlignedBox3 box;
	956	box.Initialize();
[547]	957
[562]	958	//TODO: for kd split geometry already is box => only take minmax vertices
[551]	959	if (1)
[509]	960	{
	961	PolygonContainer::const_iterator it, it_end = tData.mGeometry->mPolys.end();
[480]	962
[509]	963	for(it = tData.mGeometry->mPolys.begin(); it < it_end; ++ it)
	964	box.Include((it));
	965	}
	966	else
	967	{
	968	RayInfoContainer::const_iterator ri, ri_end = tData.mRays->end();
[480]	969
[509]	970	for(ri = tData.mRays->begin(); ri < ri_end; ++ ri)
	971	box.Include((*ri).ExtrapTermination());
	972	}
[463]	973
[605]	974	const int sAxis = mUseRandomAxis ? Random(3) : box.Size().DrivingAxis();
[480]	975
[495]	976	for (axis = 0; axis < 3; ++ axis)
[480]	977	{
[508]	978	nFrontGeom[axis] = new BspNodeGeometry();
	979	nBackGeom[axis] = new BspNodeGeometry();
	980
[545]	981	if (!mOnlyDrivingAxis \|\| (axis == sAxis))
[463]	982	{
[480]	983	if (!useCostHeuristics)
	984	{
	985	nPosition[axis] = (box.Min()[axis] + box.Max()[axis]) * 0.5f;
[543]	986	Vector3 normal(0,0,0); normal[axis] = 1.0f;
[547]	987
	988	// allows faster split because we have axis aligned kd tree boxes
[542]	989	if (useKdSplit)
	990	{
	991	nCostRatio[axis] = EvalAxisAlignedSplitCost(tData,
	992	box,
	993	axis,
[547]	994	nPosition[axis],
	995	nProbFront[axis],
	996	nProbBack[axis]);
[542]	997
[543]	998	Vector3 pos;
[542]	999
[543]	1000	pos = box.Max(); pos[axis] = nPosition[axis];
[547]	1001	AxisAlignedBox3 bBox(box.Min(), pos);
	1002	bBox.ExtractPolys(nBackGeom[axis]->mPolys);
[543]	1003
[542]	1004	pos = box.Min(); pos[axis] = nPosition[axis];
[547]	1005	AxisAlignedBox3 fBox(pos, box.Max());
	1006	fBox.ExtractPolys(nFrontGeom[axis]->mPolys);
[542]	1007	}
	1008	else
	1009	{
[573]	1010	nCostRatio[axis] =
	1011	EvalSplitPlaneCost(Plane3(normal, nPosition[axis]),
	1012	tData, nFrontGeom[axis], nBackGeom[axis],
	1013	nProbFront[axis], nProbBack[axis]);
[542]	1014	}
[480]	1015	}
[481]	1016	else
	1017	{
[482]	1018	nCostRatio[axis] =
[481]	1019	BestCostRatioHeuristics(*tData.mRays,
	1020	box,
	1021	tData.mPvs,
	1022	axis,
	1023	nPosition[axis]);
	1024	}
	1025
[480]	1026	if (bestAxis == -1)
	1027	{
	1028	bestAxis = axis;
	1029	}
	1030
	1031	else if (nCostRatio[axis] < nCostRatio[bestAxis])
	1032	{
	1033	bestAxis = axis;
	1034	}
	1035
[463]	1036	}
	1037	}
	1038
[495]	1039	//-- assign values
	1040	axis = bestAxis;
[547]	1041	pFront = nProbFront[bestAxis];
	1042	pBack = nProbBack[bestAxis];
[495]	1043
[508]	1044	// assign best split nodes geometry
	1045	*frontGeom = nFrontGeom[bestAxis];
	1046	*backGeom = nBackGeom[bestAxis];
[542]	1047
[508]	1048	// and delete other geometry
	1049	delete nFrontGeom[(bestAxis + 1) % 3];
	1050	delete nBackGeom[(bestAxis + 2) % 3];
[495]	1051
[508]	1052	//-- split plane
	1053	Vector3 normal(0,0,0); normal[bestAxis] = 1;
[480]	1054	plane = Plane3(normal, nPosition[bestAxis]);
[508]	1055
[480]	1056	return nCostRatio[bestAxis];
[463]	1057	}
	1058
[480]	1059
[508]	1060	bool VspBspTree::SelectPlane(Plane3 &bestPlane,
	1061	BspLeaf *leaf,
	1062	VspBspTraversalData &data,
	1063	VspBspTraversalData &frontData,
	1064	VspBspTraversalData &backData)
[491]	1065	{
[508]	1066	// simplest strategy: just take next polygon
	1067	if (mSplitPlaneStrategy & RANDOM_POLYGON)
[491]	1068	{
[508]	1069	if (!data.mPolygons->empty())
	1070	{
	1071	const int randIdx =
	1072	(int)RandomValue(0, (Real)((int)data.mPolygons->size() - 1));
	1073	Polygon3 nextPoly = (data.mPolygons)[randIdx];
[491]	1074
[508]	1075	bestPlane = nextPoly->GetSupportingPlane();
	1076	return true;
	1077	}
[491]	1078	}
	1079
[508]	1080	//-- use heuristics to find appropriate plane
[491]	1081
[508]	1082	// intermediate plane
	1083	Plane3 plane;
	1084	float lowestCost = MAX_FLOAT;
[517]	1085
	1086	// decides if the first few splits should be only axisAligned
[508]	1087	const bool onlyAxisAligned =
	1088	(mSplitPlaneStrategy & AXIS_ALIGNED) &&
	1089	((int)data.mRays->size() > mTermMinRaysForAxisAligned) &&
	1090	((int)data.GetAvgRayContribution() < mTermMaxRayContriForAxisAligned);
	1091
	1092	const int limit = onlyAxisAligned ? 0 :
	1093	Min((int)data.mPolygons->size(), mMaxPolyCandidates);
[491]	1094
[508]	1095	float candidateCost;
[491]	1096
[508]	1097	int maxIdx = (int)data.mPolygons->size();
[491]	1098
[508]	1099	for (int i = 0; i < limit; ++ i)
[491]	1100	{
[562]	1101	// the already taken candidates are stored behind maxIdx
	1102	// => assure that no index is taken twice
[508]	1103	const int candidateIdx = (int)RandomValue(0, (Real)(-- maxIdx));
	1104	Polygon3 poly = (data.mPolygons)[candidateIdx];
[491]	1105
[508]	1106	// swap candidate to the end to avoid testing same plane
	1107	std::swap((data.mPolygons)[maxIdx], (data.mPolygons)[candidateIdx]);
	1108	//Polygon3 poly = (data.mPolygons)[(int)RandomValue(0, (int)polys.size() - 1)];
[491]	1109
[508]	1110	// evaluate current candidate
	1111	BspNodeGeometry fGeom, bGeom;
	1112	float fArea, bArea;
	1113	plane = poly->GetSupportingPlane();
[573]	1114	candidateCost = EvalSplitPlaneCost(plane, data, fGeom, bGeom, fArea, bArea);
[491]	1115
[508]	1116	if (candidateCost < lowestCost)
[491]	1117	{
[508]	1118	bestPlane = plane;
	1119	lowestCost = candidateCost;
[491]	1120	}
	1121	}
	1122
[607]	1123	// cost ratio miss
	1124	if (!mUsePolygonSplitIfAvailable \|\| data.mPolygons->empty())
	1125	{
	1126	frontData.mIsKdNode = backData.mIsKdNode = false;
	1127	if (lowestCost > mTermMaxCostRatio)
	1128	return false;
	1129
	1130	return true;
	1131	}
	1132
[508]	1133	//-- evaluate axis aligned splits
	1134	int axis;
	1135	BspNodeGeometry fGeom, bGeom;
[547]	1136	float pFront, pBack;
[491]	1137
[508]	1138	candidateCost = SelectAxisAlignedPlane(plane, data, axis,
	1139	&fGeom, &bGeom,
[547]	1140	pFront, pBack,
[562]	1141	data.mIsKdNode);
[491]	1142
[562]	1143	bool isAxisAlignedSplit = false;
	1144
[508]	1145	if (candidateCost < lowestCost)
	1146	{
	1147	bestPlane = plane;
	1148	lowestCost = candidateCost;
[491]	1149
[542]	1150	// assign already computed values
	1151	// we can do this because we always save the
[562]	1152	// computed values from the axis aligned splits
[545]	1153	frontData.mGeometry = fGeom;
[508]	1154	backData.mGeometry = bGeom;
[547]	1155
	1156	frontData.mProbability = pFront;
	1157	backData.mProbability = pBack;
	1158
[508]	1159	//! error also computed if cost ratio is missed
[574]	1160	++ mBspStats.splits[axis];
[562]	1161	isAxisAlignedSplit = true;
[508]	1162	}
	1163	else
[463]	1164	{
[508]	1165	DEL_PTR(fGeom);
	1166	DEL_PTR(bGeom);
[463]	1167	}
	1168
[562]	1169	frontData.mIsKdNode = backData.mIsKdNode = (data.mIsKdNode && isAxisAlignedSplit);
	1170
[508]	1171	#ifdef _DEBUG
	1172	Debug << "plane lowest cost: " << lowestCost << endl;
	1173	#endif
	1174
	1175	// cost ratio miss
	1176	if (lowestCost > mTermMaxCostRatio)
	1177	return false;
	1178
	1179	return true;
[463]	1180	}
	1181
[480]	1182
[473]	1183	Plane3 VspBspTree::ChooseCandidatePlane(const RayInfoContainer &rays) const
[482]	1184	{
[473]	1185	const int candidateIdx = (int)RandomValue(0, (Real)((int)rays.size() - 1));
[482]	1186
[473]	1187	const Vector3 minPt = rays[candidateIdx].ExtrapOrigin();
	1188	const Vector3 maxPt = rays[candidateIdx].ExtrapTermination();
	1189
	1190	const Vector3 pt = (maxPt + minPt) * 0.5;
	1191	const Vector3 normal = Normalize(rays[candidateIdx].mRay->GetDir());
	1192
	1193	return Plane3(normal, pt);
	1194	}
	1195
[480]	1196
[473]	1197	Plane3 VspBspTree::ChooseCandidatePlane2(const RayInfoContainer &rays) const
[482]	1198	{
[473]	1199	Vector3 pt[3];
[482]	1200
[473]	1201	int idx[3];
	1202	int cmaxT = 0;
	1203	int cminT = 0;
	1204	bool chooseMin = false;
	1205
	1206	for (int j = 0; j < 3; ++ j)
	1207	{
	1208	idx[j] = (int)RandomValue(0, (Real)((int)rays.size() * 2 - 1));
[482]	1209
[473]	1210	if (idx[j] >= (int)rays.size())
	1211	{
	1212	idx[j] -= (int)rays.size();
[482]	1213
[473]	1214	chooseMin = (cminT < 2);
	1215	}
	1216	else
	1217	chooseMin = (cmaxT < 2);
	1218
	1219	RayInfo rayInf = rays[idx[j]];
	1220	pt[j] = chooseMin ? rayInf.ExtrapOrigin() : rayInf.ExtrapTermination();
[482]	1221	}
[473]	1222
	1223	return Plane3(pt[0], pt[1], pt[2]);
	1224	}
	1225
[580]	1226
[473]	1227	Plane3 VspBspTree::ChooseCandidatePlane3(const RayInfoContainer &rays) const
[482]	1228	{
[473]	1229	Vector3 pt[3];
[482]	1230
[473]	1231	int idx1 = (int)RandomValue(0, (Real)((int)rays.size() - 1));
	1232	int idx2 = (int)RandomValue(0, (Real)((int)rays.size() - 1));
	1233
	1234	// check if rays different
	1235	if (idx1 == idx2)
	1236	idx2 = (idx2 + 1) % (int)rays.size();
	1237
	1238	const RayInfo ray1 = rays[idx1];
	1239	const RayInfo ray2 = rays[idx2];
	1240
	1241	// normal vector of the plane parallel to both lines
	1242	const Vector3 norm = Normalize(CrossProd(ray1.mRay->GetDir(), ray2.mRay->GetDir()));
	1243
	1244	// vector from line 1 to line 2
[479]	1245	const Vector3 vd = ray2.ExtrapOrigin() - ray1.ExtrapOrigin();
[482]	1246
[473]	1247	// project vector on normal to get distance
	1248	const float dist = DotProd(vd, norm);
	1249
	1250	// point on plane lies halfway between the two planes
	1251	const Vector3 planePt = ray1.ExtrapOrigin() + norm * dist * 0.5;
	1252
	1253	return Plane3(norm, planePt);
	1254	}
	1255
[495]	1256
[463]	1257	inline void VspBspTree::GenerateUniqueIdsForPvs()
	1258	{
[580]	1259	Intersectable::NewMail(); sBackId = Intersectable::sMailId;
	1260	Intersectable::NewMail(); sFrontId = Intersectable::sMailId;
	1261	Intersectable::NewMail(); sFrontAndBackId = Intersectable::sMailId;
[463]	1262	}
	1263
[495]	1264
[573]	1265	float VspBspTree::EvalSplitPlaneCost(const Plane3 &candidatePlane,
	1266	const VspBspTraversalData &data,
	1267	BspNodeGeometry &geomFront,
	1268	BspNodeGeometry &geomBack,
	1269	float &pFront,
	1270	float &pBack) const
[463]	1271	{
[472]	1272	float cost = 0;
[463]	1273
	1274	float sumBalancedRays = 0;
	1275	float sumRaySplits = 0;
[482]	1276
[508]	1277	int pvsFront = 0;
	1278	int pvsBack = 0;
[463]	1279
[508]	1280	// probability that view point lies in back / front node
[463]	1281	float pOverall = 0;
[547]	1282	pFront = 0;
	1283	pBack = 0;
[463]	1284
[508]	1285	int raysFront = 0;
	1286	int raysBack = 0;
	1287	int totalPvs = 0;
	1288
[463]	1289	int limit;
	1290	bool useRand;
	1291
[491]	1292	// choose test rays randomly if too much
[463]	1293	if ((int)data.mRays->size() > mMaxTests)
	1294	{
	1295	useRand = true;
	1296	limit = mMaxTests;
	1297	}
	1298	else
	1299	{
	1300	useRand = false;
	1301	limit = (int)data.mRays->size();
	1302	}
[508]	1303
[463]	1304	for (int i = 0; i < limit; ++ i)
	1305	{
[508]	1306	const int testIdx = useRand ?
	1307	(int)RandomValue(0, (Real)((int)data.mRays->size() - 1)) : i;
[463]	1308	RayInfo rayInf = (*data.mRays)[testIdx];
	1309
	1310	float t;
[508]	1311	VssRay *ray = rayInf.mRay;
[463]	1312	const int cf = rayInf.ComputeRayIntersection(candidatePlane, t);
	1313
[562]	1314	if (0)
	1315	{
	1316	if (cf >= 0)
	1317	++ raysFront;
	1318	if (cf <= 0)
	1319	++ raysBack;
	1320	}
[508]	1321
[463]	1322	if (mSplitPlaneStrategy & LEAST_RAY_SPLITS)
	1323	{
	1324	sumBalancedRays += cf;
	1325	}
[482]	1326
[463]	1327	if (mSplitPlaneStrategy & BALANCED_RAYS)
	1328	{
	1329	if (cf == 0)
	1330	++ sumRaySplits;
	1331	}
	1332
	1333	if (mSplitPlaneStrategy & PVS)
	1334	{
[508]	1335	// find front and back pvs for origing and termination object
	1336	AddObjToPvs(ray->mTerminationObject, cf, pvsFront, pvsBack, totalPvs);
	1337	AddObjToPvs(ray->mOriginObject, cf, pvsFront, pvsBack, totalPvs);
[463]	1338	}
	1339	}
	1340
	1341	const float raysSize = (float)data.mRays->size() + Limits::Small;
[547]	1342
[542]	1343	if (mSplitPlaneStrategy & PVS)
	1344	{
	1345	// create unique ids for pvs heuristics
	1346	GenerateUniqueIdsForPvs();
[482]	1347
[547]	1348	// construct child geometry with regard to the candidate split plane
	1349	data.mGeometry->SplitGeometry(geomFront,
	1350	geomBack,
	1351	candidatePlane,
	1352	mBox,
	1353	mEpsilon);
	1354
	1355	pOverall = data.mProbability;
	1356
	1357	if (!mUseAreaForPvs) // use front and back cell areas to approximate volume
[542]	1358	{
[545]	1359	pFront = geomFront.GetVolume();
[603]	1360	pBack = pOverall - pFront;
[542]	1361	}
[547]	1362	else
[542]	1363	{
[547]	1364	pFront = geomFront.GetArea();
	1365	pBack = geomBack.GetArea();
[542]	1366	}
	1367	}
	1368
	1369
[463]	1370	if (mSplitPlaneStrategy & LEAST_RAY_SPLITS)
[472]	1371	cost += mLeastRaySplitsFactor * sumRaySplits / raysSize;
[463]	1372
	1373	if (mSplitPlaneStrategy & BALANCED_RAYS)
[495]	1374	cost += mBalancedRaysFactor * fabs(sumBalancedRays) / raysSize;
[482]	1375
[580]	1376	// -- pvs rendering heuristics
[463]	1377	if (mSplitPlaneStrategy & PVS)
	1378	{
[580]	1379	const int lowerPvsLimit = mViewCellsManager->GetMinPvsSize();
	1380	const int upperPvsLimit = mViewCellsManager->GetMaxPvsSize();
	1381
	1382	// only render cost heuristics or combined with standard deviation
[607]	1383	const float penaltyOld = EvalPvsPenalty(totalPvs, lowerPvsLimit, upperPvsLimit);
	1384	const float penaltyFront = EvalPvsPenalty(pvsFront, lowerPvsLimit, upperPvsLimit);
	1385	const float penaltyBack = EvalPvsPenalty(pvsBack, lowerPvsLimit, upperPvsLimit);
[579]	1386
[607]	1387	const float oldRenderCost = pOverall * (float)penaltyOld + Limits::Small;
	1388	const float newRenderCost = penaltyFront * pFront + penaltyBack * pBack;
[579]	1389
[607]	1390	float oldCost, newCost;
[580]	1391
[607]	1392	// only render cost
	1393	if (1)
	1394	{
	1395	oldCost = oldRenderCost;
	1396	newCost = newRenderCost;
	1397	}
	1398	else // also considering standard deviation
	1399	{
	1400	// standard deviation is difference of back and front pvs
	1401	const float expectedCost = 0.5f * (penaltyFront + penaltyBack);
[580]	1402
[607]	1403	const float newDeviation = 0.5f *
	1404	fabs(penaltyFront - expectedCost) + fabs(penaltyBack - expectedCost);
[580]	1405
[607]	1406	const float oldDeviation = penaltyOld + Limits::Small;
[580]	1407
[607]	1408	newCost = mRenderCostWeight * newRenderCost + (1.0f - mRenderCostWeight) * newDeviation;
	1409	oldCost = mRenderCostWeight * oldRenderCost + (1.0f - mRenderCostWeight) * oldDeviation;
	1410	}
[580]	1411
[607]	1412	cost += mPvsFactor * newCost / oldCost;
	1413
[463]	1414	}
	1415
	1416	#ifdef _DEBUG
[474]	1417	Debug << "totalpvs: " << data.mPvs << " ptotal: " << pOverall
[508]	1418	<< " frontpvs: " << pvsFront << " pFront: " << pFront
	1419	<< " backpvs: " << pvsBack << " pBack: " << pBack << endl << endl;
[600]	1420	Debug << "cost: " << cost << endl;
[463]	1421	#endif
[482]	1422
[600]	1423
[474]	1424	// normalize cost by sum of linear factors
[562]	1425	if(0)
	1426	return cost / (float)mCostNormalizer;
	1427	else
	1428	return cost;
[463]	1429	}
	1430
[508]	1431
[580]	1432	float VspBspTree::EvalAxisAlignedSplitCost(const VspBspTraversalData &data,
	1433	const AxisAlignedBox3 &box,
	1434	const int axis,
	1435	const float &position,
	1436	float &pFront,
	1437	float &pBack) const
	1438	{
	1439	int pvsTotal = 0;
	1440	int pvsFront = 0;
	1441	int pvsBack = 0;
	1442
	1443	// create unique ids for pvs heuristics
	1444	GenerateUniqueIdsForPvs();
	1445
	1446	const int pvsSize = data.mPvs;
	1447
	1448	RayInfoContainer::const_iterator rit, rit_end = data.mRays->end();
	1449
	1450	// this is the main ray classification loop!
	1451	for(rit = data.mRays->begin(); rit != rit_end; ++ rit)
	1452	{
	1453	//if (!(*rit).mRay->IsActive()) continue;
	1454
	1455	// determine the side of this ray with respect to the plane
	1456	float t;
	1457	const int side = (*rit).ComputeRayIntersection(axis, position, t);
	1458
	1459	AddObjToPvs((*rit).mRay->mTerminationObject, side, pvsFront, pvsBack, pvsTotal);
	1460	AddObjToPvs((*rit).mRay->mOriginObject, side, pvsFront, pvsBack, pvsTotal);
	1461	}
	1462
	1463	//-- pvs heuristics
	1464	float pOverall;
	1465
	1466	//-- compute heurstics
	1467	// we take simplified computation for mid split
	1468
	1469	pOverall = data.mProbability;
	1470
	1471	if (!mUseAreaForPvs)
	1472	{ // volume
	1473	pBack = pFront = pOverall * 0.5f;
	1474
	1475	#if 0
	1476	// box length substitute for probability
	1477	const float minBox = box.Min(axis);
	1478	const float maxBox = box.Max(axis);
	1479
	1480	pBack = position - minBox;
	1481	pFront = maxBox - position;
	1482	pOverall = maxBox - minBox;
	1483	#endif
	1484	}
	1485	else //-- area substitute for probability
	1486	{
	1487	const int axis2 = (axis + 1) % 3;
	1488	const int axis3 = (axis + 2) % 3;
	1489
	1490	const float faceArea =
	1491	(box.Max(axis2) - box.Min(axis2)) *
	1492	(box.Max(axis3) - box.Min(axis3));
	1493
	1494	pBack = pFront = pOverall * 0.5f + faceArea;
	1495	}
	1496
	1497	#ifdef _DEBUG
	1498	Debug << axis << " " << pvsSize << " " << pvsBack << " " << pvsFront << endl;
	1499	Debug << pFront << " " << pBack << " " << pOverall << endl;
	1500	#endif
	1501
	1502
	1503	const float newCost = pvsBack * pBack + pvsFront * pFront;
	1504	const float oldCost = (float)pvsSize * pOverall + Limits::Small;
	1505
	1506	return (mCtDivCi + newCost) / oldCost;
	1507	}
	1508
	1509
[463]	1510	void VspBspTree::AddObjToPvs(Intersectable *obj,
	1511	const int cf,
	1512	int &frontPvs,
[508]	1513	int &backPvs,
	1514	int &totalPvs) const
[463]	1515	{
	1516	if (!obj)
	1517	return;
[508]	1518
	1519	if ((obj->mMailbox != sFrontId) &&
	1520	(obj->mMailbox != sBackId) &&
	1521	(obj->mMailbox != sFrontAndBackId))
	1522	{
	1523	++ totalPvs;
	1524	}
	1525
[463]	1526	// TODO: does this really belong to no pvs?
	1527	//if (cf == Ray::COINCIDENT) return;
	1528
	1529	// object belongs to both PVS
	1530	if (cf >= 0)
	1531	{
[482]	1532	if ((obj->mMailbox != sFrontId) &&
[463]	1533	(obj->mMailbox != sFrontAndBackId))
	1534	{
	1535	++ frontPvs;
[508]	1536
[463]	1537	if (obj->mMailbox == sBackId)
[482]	1538	obj->mMailbox = sFrontAndBackId;
[463]	1539	else
[482]	1540	obj->mMailbox = sFrontId;
[463]	1541	}
	1542	}
[482]	1543
[463]	1544	if (cf <= 0)
	1545	{
	1546	if ((obj->mMailbox != sBackId) &&
	1547	(obj->mMailbox != sFrontAndBackId))
	1548	{
	1549	++ backPvs;
[508]	1550
[463]	1551	if (obj->mMailbox == sFrontId)
[482]	1552	obj->mMailbox = sFrontAndBackId;
[463]	1553	else
[482]	1554	obj->mMailbox = sBackId;
[463]	1555	}
	1556	}
	1557	}
	1558
[491]	1559
[503]	1560	void VspBspTree::CollectLeaves(vector<BspLeaf *> &leaves,
	1561	const bool onlyUnmailed,
	1562	const int maxPvsSize) const
[463]	1563	{
	1564	stack<BspNode *> nodeStack;
	1565	nodeStack.push(mRoot);
[482]	1566
	1567	while (!nodeStack.empty())
[463]	1568	{
	1569	BspNode *node = nodeStack.top();
	1570	nodeStack.pop();
[489]	1571
[482]	1572	if (node->IsLeaf())
[463]	1573	{
[490]	1574	// test if this leaf is in valid view space
[503]	1575	BspLeaf leaf = dynamic_cast<BspLeaf >(node);
	1576	if (leaf->TreeValid() &&
[508]	1577	(!onlyUnmailed \|\| !leaf->Mailed()) &&
[503]	1578	((maxPvsSize < 0) \|\| (leaf->GetViewCell()->GetPvs().GetSize() <= maxPvsSize)))
[490]	1579	{
	1580	leaves.push_back(leaf);
	1581	}
[482]	1582	}
	1583	else
[463]	1584	{
	1585	BspInterior interior = dynamic_cast<BspInterior >(node);
	1586
	1587	nodeStack.push(interior->GetBack());
	1588	nodeStack.push(interior->GetFront());
	1589	}
	1590	}
	1591	}
	1592
[489]	1593
[463]	1594	AxisAlignedBox3 VspBspTree::GetBoundingBox() const
	1595	{
	1596	return mBox;
	1597	}
	1598
[489]	1599
[463]	1600	BspNode *VspBspTree::GetRoot() const
	1601	{
	1602	return mRoot;
	1603	}
	1604
[489]	1605
[463]	1606	void VspBspTree::EvaluateLeafStats(const VspBspTraversalData &data)
	1607	{
	1608	// the node became a leaf -> evaluate stats for leafs
	1609	BspLeaf leaf = dynamic_cast<BspLeaf >(data.mNode);
	1610
	1611	// store maximal and minimal depth
[574]	1612	if (data.mDepth > mBspStats.maxDepth)
	1613	mBspStats.maxDepth = data.mDepth;
[463]	1614
[574]	1615	if (data.mPvs > mBspStats.maxPvs)
	1616	mBspStats.maxPvs = data.mPvs;
[508]	1617
[574]	1618	if (data.mDepth < mBspStats.minDepth)
	1619	mBspStats.minDepth = data.mDepth;
[508]	1620
[463]	1621	if (data.mDepth >= mTermMaxDepth)
[574]	1622	++ mBspStats.maxDepthNodes;
[508]	1623	// accumulate rays to compute rays / leaf
[574]	1624	mBspStats.accumRays += (int)data.mRays->size();
[463]	1625
[437]	1626	if (data.mPvs < mTermMinPvs)
[574]	1627	++ mBspStats.minPvsNodes;
[437]	1628
	1629	if ((int)data.mRays->size() < mTermMinRays)
[574]	1630	++ mBspStats.minRaysNodes;
[437]	1631
	1632	if (data.GetAvgRayContribution() > mTermMaxRayContribution)
[574]	1633	++ mBspStats.maxRayContribNodes;
[482]	1634
[547]	1635	if (data.mProbability <= mTermMinProbability)
[574]	1636	++ mBspStats.minProbabilityNodes;
[508]	1637
[474]	1638	// accumulate depth to compute average depth
[574]	1639	mBspStats.accumDepth += data.mDepth;
[463]	1640
	1641	#ifdef _DEBUG
	1642	Debug << "BSP stats: "
	1643	<< "Depth: " << data.mDepth << " (max: " << mTermMaxDepth << "), "
	1644	<< "PVS: " << data.mPvs << " (min: " << mTermMinPvs << "), "
[587]	1645	// << "Area: " << data.mProbability << " (min: " << mTermMinProbability << "), "
[482]	1646	<< "#rays: " << (int)data.mRays->size() << " (max: " << mTermMinRays << "), "
[463]	1647	<< "#pvs: " << leaf->GetViewCell()->GetPvs().GetSize() << "=, "
	1648	<< "#avg ray contrib (pvs): " << (float)data.mPvs / (float)data.mRays->size() << endl;
	1649	#endif
	1650	}
	1651
	1652	int VspBspTree::CastRay(Ray &ray)
	1653	{
	1654	int hits = 0;
[482]	1655
[600]	1656	stack<BspRayTraversalData> tQueue;
[482]	1657
[463]	1658	float maxt, mint;
	1659
	1660	if (!mBox.GetRaySegment(ray, mint, maxt))
	1661	return 0;
	1662
	1663	Intersectable::NewMail();
[600]	1664	ViewCell::NewMail();
[463]	1665	Vector3 entp = ray.Extrap(mint);
	1666	Vector3 extp = ray.Extrap(maxt);
[482]	1667
[463]	1668	BspNode *node = mRoot;
	1669	BspNode *farChild = NULL;
[482]	1670
[463]	1671	while (1)
	1672	{
[482]	1673	if (!node->IsLeaf())
[463]	1674	{
	1675	BspInterior in = dynamic_cast<BspInterior >(node);
	1676
	1677	Plane3 splitPlane = in->GetPlane();
	1678	const int entSide = splitPlane.Side(entp);
	1679	const int extSide = splitPlane.Side(extp);
	1680
	1681	if (entSide < 0)
	1682	{
	1683	node = in->GetBack();
	1684
	1685	if(extSide <= 0) // plane does not split ray => no far child
	1686	continue;
[482]	1687
[463]	1688	farChild = in->GetFront(); // plane splits ray
	1689
	1690	} else if (entSide > 0)
	1691	{
	1692	node = in->GetFront();
	1693
	1694	if (extSide >= 0) // plane does not split ray => no far child
	1695	continue;
	1696
[482]	1697	farChild = in->GetBack(); // plane splits ray
[463]	1698	}
	1699	else // ray and plane are coincident
	1700	{
	1701	// WHAT TO DO IN THIS CASE ?
	1702	//break;
	1703	node = in->GetFront();
	1704	continue;
	1705	}
	1706
	1707	// push data for far child
[600]	1708	tQueue.push(BspRayTraversalData(farChild, extp, maxt));
[463]	1709
	1710	// find intersection of ray segment with plane
	1711	float t;
	1712	extp = splitPlane.FindIntersection(ray.GetLoc(), extp, &t);
	1713	maxt *= t;
[482]	1714
[463]	1715	} else // reached leaf => intersection with view cell
	1716	{
	1717	BspLeaf leaf = dynamic_cast<BspLeaf >(node);
[482]	1718
[463]	1719	if (!leaf->GetViewCell()->Mailed())
	1720	{
	1721	//ray.bspIntersections.push_back(Ray::VspBspIntersection(maxt, leaf));
	1722	leaf->GetViewCell()->Mail();
	1723	++ hits;
	1724	}
[482]	1725
[463]	1726	//-- fetch the next far child from the stack
[600]	1727	if (tQueue.empty())
[463]	1728	break;
[482]	1729
[463]	1730	entp = extp;
	1731	mint = maxt; // NOTE: need this?
	1732
	1733	if (ray.GetType() == Ray::LINE_SEGMENT && mint > 1.0f)
	1734	break;
	1735
[600]	1736	BspRayTraversalData &s = tQueue.top();
[463]	1737
	1738	node = s.mNode;
	1739	extp = s.mExitPoint;
	1740	maxt = s.mMaxT;
	1741
[600]	1742	tQueue.pop();
[463]	1743	}
	1744	}
	1745
	1746	return hits;
	1747	}
	1748
[532]	1749
[547]	1750	void VspBspTree::CollectViewCells(ViewCellContainer &viewCells, bool onlyValid) const
[463]	1751	{
[532]	1752	ViewCell::NewMail();
[551]	1753
	1754	CollectViewCells(mRoot, onlyValid, viewCells, true);
[532]	1755	}
	1756
	1757
[574]	1758	void VspBspTree::CollapseViewCells()
[542]	1759	{
[590]	1760	// TODO
	1761	#if VC_HISTORY
[542]	1762	stack<BspNode *> nodeStack;
	1763
	1764	if (!mRoot)
	1765	return;
	1766
	1767	nodeStack.push(mRoot);
	1768
	1769	while (!nodeStack.empty())
	1770	{
	1771	BspNode *node = nodeStack.top();
	1772	nodeStack.pop();
	1773
	1774	if (node->IsLeaf())
[574]	1775	{
	1776	BspViewCell viewCell = dynamic_cast<BspLeaf >(node)->GetViewCell();
[542]	1777
[574]	1778	if (!viewCell->GetValid())
[542]	1779	{
	1780	BspViewCell viewCell = dynamic_cast<BspLeaf >(node)->GetViewCell();
[580]	1781
	1782	ViewCellContainer leaves;
[590]	1783	mViewCellsTree->CollectLeaves(viewCell, leaves);
[580]	1784
	1785	ViewCellContainer::const_iterator it, it_end = leaves.end();
	1786
	1787	for (it = leaves.begin(); it != it_end; ++ it)
[542]	1788	{
[580]	1789	BspLeaf l = dynamic_cast<BspViewCell >(*it)->mLeaf;
[574]	1790	l->SetViewCell(GetOrCreateOutOfBoundsCell());
	1791	++ mBspStats.invalidLeaves;
	1792	}
[542]	1793
[574]	1794	// add to unbounded view cell
	1795	GetOrCreateOutOfBoundsCell()->GetPvs().AddPvs(viewCell->GetPvs());
	1796	DEL_PTR(viewCell);
	1797	}
	1798	}
	1799	else
	1800	{
	1801	BspInterior interior = dynamic_cast<BspInterior >(node);
	1802
	1803	nodeStack.push(interior->GetFront());
	1804	nodeStack.push(interior->GetBack());
	1805	}
	1806	}
[542]	1807
[574]	1808	Debug << "invalid leaves: " << mBspStats.invalidLeaves << endl;
[590]	1809	#endif
[574]	1810	}
	1811
	1812
	1813	void VspBspTree::ValidateTree()
	1814	{
	1815	stack<BspNode *> nodeStack;
	1816
	1817	if (!mRoot)
	1818	return;
	1819
	1820	nodeStack.push(mRoot);
	1821
	1822	mBspStats.invalidLeaves = 0;
	1823	while (!nodeStack.empty())
	1824	{
	1825	BspNode *node = nodeStack.top();
	1826	nodeStack.pop();
	1827
	1828	if (node->IsLeaf())
	1829	{
	1830	BspLeaf leaf = dynamic_cast<BspLeaf >(node);
	1831
	1832	if (!leaf->GetViewCell()->GetValid())
	1833	++ mBspStats.invalidLeaves;
	1834
	1835	// validity flags don't match => repair
	1836	if (leaf->GetViewCell()->GetValid() != leaf->TreeValid())
	1837	{
	1838	leaf->SetTreeValid(leaf->GetViewCell()->GetValid());
	1839	PropagateUpValidity(leaf);
[542]	1840	}
	1841	}
	1842	else
	1843	{
	1844	BspInterior interior = dynamic_cast<BspInterior >(node);
	1845
	1846	nodeStack.push(interior->GetFront());
	1847	nodeStack.push(interior->GetBack());
	1848	}
	1849	}
[562]	1850
[574]	1851	Debug << "invalid leaves: " << mBspStats.invalidLeaves << endl;
[542]	1852	}
	1853
[547]	1854
[574]	1855
[547]	1856	void VspBspTree::CollectViewCells(BspNode *root,
	1857	bool onlyValid,
[532]	1858	ViewCellContainer &viewCells,
	1859	bool onlyUnmailed) const
	1860	{
[498]	1861	stack<BspNode *> nodeStack;
[463]	1862
[538]	1863	if (!root)
[508]	1864	return;
[463]	1865
[538]	1866	nodeStack.push(root);
[498]	1867
	1868	while (!nodeStack.empty())
	1869	{
	1870	BspNode *node = nodeStack.top();
	1871	nodeStack.pop();
	1872
	1873	if (node->IsLeaf())
	1874	{
[564]	1875	if (!onlyValid \|\| node->TreeValid())
[498]	1876	{
[590]	1877	ViewCell *viewCell =
	1878	mViewCellsTree->GetActiveViewCell(dynamic_cast<BspLeaf *>(node)->GetViewCell());
	1879
[532]	1880	if (!onlyUnmailed \|\| !viewCell->Mailed())
[498]	1881	{
	1882	viewCell->Mail();
	1883	viewCells.push_back(viewCell);
	1884	}
	1885	}
	1886	}
	1887	else
	1888	{
	1889	BspInterior interior = dynamic_cast<BspInterior >(node);
[508]	1890
[498]	1891	nodeStack.push(interior->GetFront());
	1892	nodeStack.push(interior->GetBack());
	1893	}
	1894	}
[590]	1895
[463]	1896	}
	1897
	1898
	1899	float VspBspTree::AccumulatedRayLength(const RayInfoContainer &rays) const
	1900	{
	1901	float len = 0;
	1902
	1903	RayInfoContainer::const_iterator it, it_end = rays.end();
	1904
	1905	for (it = rays.begin(); it != it_end; ++ it)
	1906	len += (*it).SegmentLength();
	1907
	1908	return len;
	1909	}
	1910
[479]	1911
[463]	1912	int VspBspTree::SplitRays(const Plane3 &plane,
[482]	1913	RayInfoContainer &rays,
	1914	RayInfoContainer &frontRays,
[463]	1915	RayInfoContainer &backRays)
	1916	{
	1917	int splits = 0;
	1918
[574]	1919	RayInfoContainer::const_iterator it, it_end = rays.end();
	1920
	1921	for (it = rays.begin(); it != it_end; ++ it)
[463]	1922	{
[574]	1923	RayInfo bRay = *it;
	1924
[463]	1925	VssRay *ray = bRay.mRay;
[473]	1926	float t;
[463]	1927
[485]	1928	// get classification and receive new t
[463]	1929	const int cf = bRay.ComputeRayIntersection(plane, t);
[482]	1930
[463]	1931	switch (cf)
	1932	{
	1933	case -1:
	1934	backRays.push_back(bRay);
	1935	break;
	1936	case 1:
	1937	frontRays.push_back(bRay);
	1938	break;
[482]	1939	case 0:
	1940	{
[485]	1941	//-- split ray
	1942	//-- test if start point behind or in front of plane
	1943	const int side = plane.Side(bRay.ExtrapOrigin());
	1944
	1945	if (side <= 0)
	1946	{
	1947	backRays.push_back(RayInfo(ray, bRay.GetMinT(), t));
	1948	frontRays.push_back(RayInfo(ray, t, bRay.GetMaxT()));
	1949	}
	1950	else
	1951	{
	1952	frontRays.push_back(RayInfo(ray, bRay.GetMinT(), t));
	1953	backRays.push_back(RayInfo(ray, t, bRay.GetMaxT()));
	1954	}
[463]	1955	}
	1956	break;
	1957	default:
[485]	1958	Debug << "Should not come here" << endl;
[463]	1959	break;
	1960	}
	1961	}
	1962
	1963	return splits;
	1964	}
	1965
[479]	1966
[463]	1967	void VspBspTree::ExtractHalfSpaces(BspNode *n, vector<Plane3> &halfSpaces) const
	1968	{
	1969	BspNode *lastNode;
	1970
	1971	do
	1972	{
	1973	lastNode = n;
	1974
	1975	// want to get planes defining geometry of this node => don't take
	1976	// split plane of node itself
	1977	n = n->GetParent();
[482]	1978
[463]	1979	if (n)
	1980	{
	1981	BspInterior interior = dynamic_cast<BspInterior >(n);
	1982	Plane3 halfSpace = dynamic_cast<BspInterior *>(interior)->GetPlane();
	1983
	1984	if (interior->GetFront() != lastNode)
	1985	halfSpace.ReverseOrientation();
	1986
	1987	halfSpaces.push_back(halfSpace);
	1988	}
	1989	}
	1990	while (n);
	1991	}
	1992
[485]	1993
[482]	1994	void VspBspTree::ConstructGeometry(BspNode *n,
[503]	1995	BspNodeGeometry &geom) const
[463]	1996	{
[437]	1997	vector<Plane3> halfSpaces;
	1998	ExtractHalfSpaces(n, halfSpaces);
	1999
	2000	PolygonContainer candidatePolys;
	2001
	2002	// bounded planes are added to the polygons (reverse polygons
	2003	// as they have to be outfacing
	2004	for (int i = 0; i < (int)halfSpaces.size(); ++ i)
	2005	{
	2006	Polygon3 *p = GetBoundingBox().CrossSection(halfSpaces[i]);
[482]	2007
[448]	2008	if (p->Valid(mEpsilon))
[437]	2009	{
	2010	candidatePolys.push_back(p->CreateReversePolygon());
	2011	DEL_PTR(p);
	2012	}
	2013	}
	2014
	2015	// add faces of bounding box (also could be faces of the cell)
	2016	for (int i = 0; i < 6; ++ i)
	2017	{
	2018	VertexContainer vertices;
[482]	2019
[437]	2020	for (int j = 0; j < 4; ++ j)
	2021	vertices.push_back(mBox.GetFace(i).mVertices[j]);
	2022
	2023	candidatePolys.push_back(new Polygon3(vertices));
	2024	}
	2025
	2026	for (int i = 0; i < (int)candidatePolys.size(); ++ i)
	2027	{
	2028	// polygon is split by all other planes
	2029	for (int j = 0; (j < (int)halfSpaces.size()) && candidatePolys[i]; ++ j)
	2030	{
	2031	if (i == j) // polygon and plane are coincident
	2032	continue;
	2033
	2034	VertexContainer splitPts;
	2035	Polygon3 frontPoly, backPoly;
	2036
[482]	2037	const int cf =
[448]	2038	candidatePolys[i]->ClassifyPlane(halfSpaces[j],
	2039	mEpsilon);
[482]	2040
[437]	2041	switch (cf)
	2042	{
	2043	case Polygon3::SPLIT:
	2044	frontPoly = new Polygon3();
	2045	backPoly = new Polygon3();
	2046
[482]	2047	candidatePolys[i]->Split(halfSpaces[j],
	2048	*frontPoly,
[448]	2049	*backPoly,
	2050	mEpsilon);
[437]	2051
	2052	DEL_PTR(candidatePolys[i]);
	2053
[448]	2054	if (frontPoly->Valid(mEpsilon))
[437]	2055	candidatePolys[i] = frontPoly;
	2056	else
	2057	DEL_PTR(frontPoly);
	2058
	2059	DEL_PTR(backPoly);
	2060	break;
	2061	case Polygon3::BACK_SIDE:
	2062	DEL_PTR(candidatePolys[i]);
	2063	break;
	2064	// just take polygon as it is
	2065	case Polygon3::FRONT_SIDE:
	2066	case Polygon3::COINCIDENT:
	2067	default:
	2068	break;
	2069	}
	2070	}
[482]	2071
[437]	2072	if (candidatePolys[i])
[503]	2073	geom.mPolys.push_back(candidatePolys[i]);
[437]	2074	}
[463]	2075	}
	2076
[485]	2077
[582]	2078	void VspBspTree::ConstructGeometry(ViewCell *vc,
[503]	2079	BspNodeGeometry &vcGeom) const
[589]	2080	{
[580]	2081	ViewCellContainer leaves;
[590]	2082
	2083	mViewCellsTree->CollectLeaves(vc, leaves);
[463]	2084
[580]	2085	ViewCellContainer::const_iterator it, it_end = leaves.end();
	2086
[463]	2087	for (it = leaves.begin(); it != it_end; ++ it)
[580]	2088	{
	2089	BspLeaf l = dynamic_cast<BspViewCell >(*it)->mLeaf;
[590]	2090
[580]	2091	ConstructGeometry(l, vcGeom);
	2092	}
[463]	2093	}
	2094
[485]	2095
[551]	2096	typedef pair<BspNode , BspNodeGeometry > bspNodePair;
	2097
[557]	2098
[482]	2099	int VspBspTree::FindNeighbors(BspNode n, vector<BspLeaf > &neighbors,
[562]	2100	const bool onlyUnmailed) const
[463]	2101	{
[551]	2102	stack<bspNodePair> nodeStack;
	2103
	2104	BspNodeGeometry nodeGeom;
	2105	ConstructGeometry(n, nodeGeom);
	2106
[500]	2107	// split planes from the root to this node
	2108	// needed to verify that we found neighbor leaf
[557]	2109	// TODO: really needed?
[463]	2110	vector<Plane3> halfSpaces;
	2111	ExtractHalfSpaces(n, halfSpaces);
	2112
[551]	2113
	2114	BspNodeGeometry *rgeom = new BspNodeGeometry();
	2115	ConstructGeometry(mRoot, *rgeom);
	2116
	2117	nodeStack.push(bspNodePair(mRoot, rgeom));
	2118
[482]	2119	while (!nodeStack.empty())
[463]	2120	{
[551]	2121	BspNode *node = nodeStack.top().first;
	2122	BspNodeGeometry *geom = nodeStack.top().second;
[562]	2123
[463]	2124	nodeStack.pop();
	2125
[557]	2126	if (node->IsLeaf())
[562]	2127	{
[557]	2128	// test if this leaf is in valid view space
	2129	if (node->TreeValid() &&
	2130	(node != n) &&
	2131	(!onlyUnmailed \|\| !node->Mailed()))
	2132	{
	2133	bool isAdjacent = true;
[551]	2134
[570]	2135	if (1)
[557]	2136	{
[562]	2137	// test all planes of current node if still adjacent
	2138	for (int i = 0; (i < halfSpaces.size()) && isAdjacent; ++ i)
	2139	{
	2140	const int cf =
	2141	Polygon3::ClassifyPlane(geom->mPolys,
	2142	halfSpaces[i],
	2143	mEpsilon);
[482]	2144
[562]	2145	if (cf == Polygon3::BACK_SIDE)
	2146	{
	2147	isAdjacent = false;
	2148	}
[557]	2149	}
	2150	}
[562]	2151	else
[557]	2152	{
[562]	2153	// TODO: why is this wrong??
	2154	// test all planes of current node if still adjacent
	2155	for (int i = 0; (i < (int)nodeGeom.mPolys.size()) && isAdjacent; ++ i)
	2156	{
	2157	Polygon3 *poly = nodeGeom.mPolys[i];
[555]	2158
[562]	2159	const int cf =
	2160	Polygon3::ClassifyPlane(geom->mPolys,
	2161	poly->GetSupportingPlane(),
	2162	mEpsilon);
[557]	2163
[562]	2164	if (cf == Polygon3::BACK_SIDE)
	2165	{
	2166	isAdjacent = false;
	2167	}
[557]	2168	}
[570]	2169	}
[557]	2170	// neighbor was found
	2171	if (isAdjacent)
[562]	2172	{
[551]	2173	neighbors.push_back(dynamic_cast<BspLeaf *>(node));
[562]	2174	}
[463]	2175	}
[562]	2176	}
	2177	else
	2178	{
	2179	BspInterior interior = dynamic_cast<BspInterior >(node);
[482]	2180
[562]	2181	const int cf = Polygon3::ClassifyPlane(nodeGeom.mPolys,
	2182	interior->GetPlane(),
	2183	mEpsilon);
[551]	2184
[562]	2185	BspNode *front = interior->GetFront();
	2186	BspNode *back = interior->GetBack();
[551]	2187
[562]	2188	BspNodeGeometry *fGeom = new BspNodeGeometry();
	2189	BspNodeGeometry *bGeom = new BspNodeGeometry();
[463]	2190
[562]	2191	geom->SplitGeometry(*fGeom,
	2192	*bGeom,
	2193	interior->GetPlane(),
	2194	mBox,
	2195	mEpsilon);
[551]	2196
[562]	2197	if (cf == Polygon3::FRONT_SIDE)
	2198	{
	2199	nodeStack.push(bspNodePair(interior->GetFront(), fGeom));
	2200	DEL_PTR(bGeom);
	2201	}
	2202	else
	2203	{
	2204	if (cf == Polygon3::BACK_SIDE)
[551]	2205	{
[562]	2206	nodeStack.push(bspNodePair(interior->GetBack(), bGeom));
	2207	DEL_PTR(fGeom);
[551]	2208	}
[482]	2209	else
[562]	2210	{ // random decision
	2211	nodeStack.push(bspNodePair(front, fGeom));
	2212	nodeStack.push(bspNodePair(back, bGeom));
[463]	2213	}
[551]	2214	}
[463]	2215	}
[562]	2216
[551]	2217	DEL_PTR(geom);
[463]	2218	}
[482]	2219
[463]	2220	return (int)neighbors.size();
	2221	}
	2222
[489]	2223
[600]	2224
	2225	int VspBspTree::FindApproximateNeighbors(BspNode n, vector<BspLeaf > &neighbors,
	2226	const bool onlyUnmailed) const
	2227	{
	2228	stack<bspNodePair> nodeStack;
	2229
	2230	BspNodeGeometry nodeGeom;
	2231	ConstructGeometry(n, nodeGeom);
	2232
	2233	// split planes from the root to this node
	2234	// needed to verify that we found neighbor leaf
	2235	// TODO: really needed?
	2236	vector<Plane3> halfSpaces;
	2237	ExtractHalfSpaces(n, halfSpaces);
	2238
	2239
	2240	BspNodeGeometry *rgeom = new BspNodeGeometry();
	2241	ConstructGeometry(mRoot, *rgeom);
	2242
	2243	nodeStack.push(bspNodePair(mRoot, rgeom));
	2244
	2245	while (!nodeStack.empty())
	2246	{
	2247	BspNode *node = nodeStack.top().first;
	2248	BspNodeGeometry *geom = nodeStack.top().second;
	2249
	2250	nodeStack.pop();
	2251
	2252	if (node->IsLeaf())
	2253	{
	2254	// test if this leaf is in valid view space
	2255	if (node->TreeValid() &&
	2256	(node != n) &&
	2257	(!onlyUnmailed \|\| !node->Mailed()))
	2258	{
	2259	bool isAdjacent = true;
	2260
	2261	// neighbor was found
	2262	if (isAdjacent)
	2263	{
	2264	neighbors.push_back(dynamic_cast<BspLeaf *>(node));
	2265	}
	2266	}
	2267	}
	2268	else
	2269	{
	2270	BspInterior interior = dynamic_cast<BspInterior >(node);
	2271
	2272	const int cf = Polygon3::ClassifyPlane(nodeGeom.mPolys,
	2273	interior->GetPlane(),
	2274	mEpsilon);
	2275
	2276	BspNode *front = interior->GetFront();
	2277	BspNode *back = interior->GetBack();
	2278
	2279	BspNodeGeometry *fGeom = new BspNodeGeometry();
	2280	BspNodeGeometry *bGeom = new BspNodeGeometry();
	2281
	2282	geom->SplitGeometry(*fGeom,
	2283	*bGeom,
	2284	interior->GetPlane(),
	2285	mBox,
	2286	mEpsilon);
	2287
	2288	if (cf == Polygon3::FRONT_SIDE)
	2289	{
	2290	nodeStack.push(bspNodePair(interior->GetFront(), fGeom));
	2291	DEL_PTR(bGeom);
	2292	}
	2293	else
	2294	{
	2295	if (cf == Polygon3::BACK_SIDE)
	2296	{
	2297	nodeStack.push(bspNodePair(interior->GetBack(), bGeom));
	2298	DEL_PTR(fGeom);
	2299	}
	2300	else
	2301	{ // random decision
	2302	nodeStack.push(bspNodePair(front, fGeom));
	2303	nodeStack.push(bspNodePair(back, bGeom));
	2304	}
	2305	}
	2306	}
	2307
	2308	DEL_PTR(geom);
	2309	}
	2310
	2311	return (int)neighbors.size();
	2312	}
	2313
	2314
	2315
[463]	2316	BspLeaf *VspBspTree::GetRandomLeaf(const Plane3 &halfspace)
	2317	{
	2318	stack<BspNode *> nodeStack;
	2319	nodeStack.push(mRoot);
[482]	2320
[463]	2321	int mask = rand();
[482]	2322
	2323	while (!nodeStack.empty())
[463]	2324	{
	2325	BspNode *node = nodeStack.top();
	2326	nodeStack.pop();
[482]	2327
	2328	if (node->IsLeaf())
[463]	2329	{
	2330	return dynamic_cast<BspLeaf *>(node);
[482]	2331	}
	2332	else
[463]	2333	{
	2334	BspInterior interior = dynamic_cast<BspInterior >(node);
	2335	BspNode *next;
[503]	2336	BspNodeGeometry geom;
[482]	2337
[463]	2338	// todo: not very efficient: constructs full cell everytime
[498]	2339	ConstructGeometry(interior, geom);
[463]	2340
[503]	2341	const int cf =
	2342	Polygon3::ClassifyPlane(geom.mPolys, halfspace, mEpsilon);
[463]	2343
	2344	if (cf == Polygon3::BACK_SIDE)
	2345	next = interior->GetFront();
	2346	else
	2347	if (cf == Polygon3::FRONT_SIDE)
	2348	next = interior->GetFront();
[482]	2349	else
[463]	2350	{
	2351	// random decision
	2352	if (mask & 1)
	2353	next = interior->GetBack();
	2354	else
	2355	next = interior->GetFront();
	2356	mask = mask >> 1;
	2357	}
	2358
	2359	nodeStack.push(next);
	2360	}
	2361	}
[482]	2362
[463]	2363	return NULL;
	2364	}
	2365
	2366	BspLeaf *VspBspTree::GetRandomLeaf(const bool onlyUnmailed)
	2367	{
	2368	stack<BspNode *> nodeStack;
[482]	2369
[463]	2370	nodeStack.push(mRoot);
	2371
	2372	int mask = rand();
[482]	2373
	2374	while (!nodeStack.empty())
[463]	2375	{
	2376	BspNode *node = nodeStack.top();
	2377	nodeStack.pop();
[482]	2378
	2379	if (node->IsLeaf())
[463]	2380	{
	2381	if ( (!onlyUnmailed \|\| !node->Mailed()) )
	2382	return dynamic_cast<BspLeaf *>(node);
	2383	}
[482]	2384	else
[463]	2385	{
	2386	BspInterior interior = dynamic_cast<BspInterior >(node);
	2387
	2388	// random decision
	2389	if (mask & 1)
	2390	nodeStack.push(interior->GetBack());
	2391	else
	2392	nodeStack.push(interior->GetFront());
	2393
	2394	mask = mask >> 1;
	2395	}
	2396	}
[482]	2397
[463]	2398	return NULL;
	2399	}
	2400
	2401	int VspBspTree::ComputePvsSize(const RayInfoContainer &rays) const
	2402	{
	2403	int pvsSize = 0;
	2404
	2405	RayInfoContainer::const_iterator rit, rit_end = rays.end();
	2406
	2407	Intersectable::NewMail();
	2408
	2409	for (rit = rays.begin(); rit != rays.end(); ++ rit)
	2410	{
	2411	VssRay ray = (rit).mRay;
[482]	2412
[463]	2413	if (ray->mOriginObject)
	2414	{
	2415	if (!ray->mOriginObject->Mailed())
	2416	{
	2417	ray->mOriginObject->Mail();
	2418	++ pvsSize;
	2419	}
	2420	}
	2421	if (ray->mTerminationObject)
	2422	{
	2423	if (!ray->mTerminationObject->Mailed())
	2424	{
	2425	ray->mTerminationObject->Mail();
	2426	++ pvsSize;
	2427	}
	2428	}
	2429	}
	2430
	2431	return pvsSize;
	2432	}
	2433
	2434	float VspBspTree::GetEpsilon() const
	2435	{
	2436	return mEpsilon;
	2437	}
	2438
	2439
	2440	int VspBspTree::SplitPolygons(const Plane3 &plane,
[482]	2441	PolygonContainer &polys,
	2442	PolygonContainer &frontPolys,
	2443	PolygonContainer &backPolys,
[463]	2444	PolygonContainer &coincident) const
	2445	{
	2446	int splits = 0;
	2447
[574]	2448	PolygonContainer::const_iterator it, it_end = polys.end();
	2449
	2450	for (it = polys.begin(); it != polys.end(); ++ it)
[463]	2451	{
[574]	2452	Polygon3 poly = it;
[463]	2453
	2454	// classify polygon
	2455	const int cf = poly->ClassifyPlane(plane, mEpsilon);
	2456
	2457	switch (cf)
	2458	{
	2459	case Polygon3::COINCIDENT:
	2460	coincident.push_back(poly);
[482]	2461	break;
	2462	case Polygon3::FRONT_SIDE:
[463]	2463	frontPolys.push_back(poly);
	2464	break;
	2465	case Polygon3::BACK_SIDE:
	2466	backPolys.push_back(poly);
	2467	break;
	2468	case Polygon3::SPLIT:
	2469	backPolys.push_back(poly);
	2470	frontPolys.push_back(poly);
	2471	++ splits;
	2472	break;
	2473	default:
	2474	Debug << "SHOULD NEVER COME HERE\n";
	2475	break;
	2476	}
	2477	}
	2478
	2479	return splits;
	2480	}
[466]	2481
	2482
[469]	2483	int VspBspTree::CastLineSegment(const Vector3 &origin,
	2484	const Vector3 &termination,
	2485	vector<ViewCell *> &viewcells)
[466]	2486	{
[469]	2487	int hits = 0;
[600]	2488	stack<BspRayTraversalData> tQueue;
[482]	2489
[469]	2490	float mint = 0.0f, maxt = 1.0f;
[482]	2491
[469]	2492	Intersectable::NewMail();
[600]	2493	ViewCell::NewMail();
[482]	2494
[469]	2495	Vector3 entp = origin;
	2496	Vector3 extp = termination;
[482]	2497
[469]	2498	BspNode *node = mRoot;
	2499	BspNode *farChild = NULL;
[482]	2500
[485]	2501	float t;
[482]	2502	while (1)
[469]	2503	{
[482]	2504	if (!node->IsLeaf())
[469]	2505	{
	2506	BspInterior in = dynamic_cast<BspInterior >(node);
[482]	2507
[469]	2508	Plane3 splitPlane = in->GetPlane();
[485]	2509
[469]	2510	const int entSide = splitPlane.Side(entp);
	2511	const int extSide = splitPlane.Side(extp);
[482]	2512
[485]	2513	if (entSide < 0)
[469]	2514	{
[576]	2515	node = in->GetBack();
	2516	// plane does not split ray => no far child
	2517	if (extSide <= 0)
	2518	continue;
	2519
	2520	farChild = in->GetFront(); // plane splits ray
[485]	2521	}
	2522	else if (entSide > 0)
[469]	2523	{
[576]	2524	node = in->GetFront();
[482]	2525
[576]	2526	if (extSide >= 0) // plane does not split ray => no far child
	2527	continue;
[482]	2528
	2529	farChild = in->GetBack(); // plane splits ray
[469]	2530	}
[485]	2531	else // ray end point on plane
	2532	{ // NOTE: what to do if ray is coincident with plane?
	2533	if (extSide < 0)
	2534	node = in->GetBack();
[487]	2535	else
[485]	2536	node = in->GetFront();
[487]	2537
[485]	2538	continue; // no far child
[469]	2539	}
[482]	2540
[469]	2541	// push data for far child
[600]	2542	tQueue.push(BspRayTraversalData(farChild, extp));
[482]	2543
[469]	2544	// find intersection of ray segment with plane
	2545	extp = splitPlane.FindIntersection(origin, extp, &t);
[485]	2546	}
	2547	else
[469]	2548	{
	2549	// reached leaf => intersection with view cell
	2550	BspLeaf leaf = dynamic_cast<BspLeaf >(node);
[482]	2551
[590]	2552	ViewCell *viewCell = mViewCellsTree->GetActiveViewCell(leaf->GetViewCell());
	2553	if (!viewCell->Mailed())
[469]	2554	{
[590]	2555	viewcells.push_back(viewCell);
	2556	viewCell->Mail();
[469]	2557	++ hits;
	2558	}
[482]	2559
[469]	2560	//-- fetch the next far child from the stack
[600]	2561	if (tQueue.empty())
[469]	2562	break;
[482]	2563
[469]	2564	entp = extp;
[485]	2565
[600]	2566	BspRayTraversalData &s = tQueue.top();
[482]	2567
[469]	2568	node = s.mNode;
	2569	extp = s.mExitPoint;
[482]	2570
[600]	2571	tQueue.pop();
[469]	2572	}
[466]	2573	}
[487]	2574
[469]	2575	return hits;
[466]	2576	}
[478]	2577
[576]	2578
	2579
	2580
[485]	2581	int VspBspTree::TreeDistance(BspNode n1, BspNode n2) const
[482]	2582	{
	2583	std::deque<BspNode *> path1;
	2584	BspNode *p1 = n1;
[479]	2585
[482]	2586	// create path from node 1 to root
	2587	while (p1)
	2588	{
	2589	if (p1 == n2) // second node on path
	2590	return (int)path1.size();
	2591
	2592	path1.push_front(p1);
	2593	p1 = p1->GetParent();
	2594	}
	2595
	2596	int depth = n2->GetDepth();
	2597	int d = depth;
	2598
	2599	BspNode *p2 = n2;
	2600
	2601	// compare with same depth
	2602	while (1)
	2603	{
	2604	if ((d < (int)path1.size()) && (p2 == path1[d]))
	2605	return (depth - d) + ((int)path1.size() - 1 - d);
	2606
	2607	-- d;
	2608	p2 = p2->GetParent();
	2609	}
	2610
	2611	return 0; // never come here
	2612	}
	2613
[580]	2614
[501]	2615	BspNode VspBspTree::CollapseTree(BspNode node, int &collapsed)
[479]	2616	{
[590]	2617	// TODO
	2618	#if VC_HISTORY
[495]	2619	if (node->IsLeaf())
[479]	2620	return node;
	2621
[492]	2622	BspInterior interior = dynamic_cast<BspInterior >(node);
	2623
[501]	2624	BspNode *front = CollapseTree(interior->GetFront(), collapsed);
	2625	BspNode *back = CollapseTree(interior->GetBack(), collapsed);
[492]	2626
[479]	2627	if (front->IsLeaf() && back->IsLeaf())
	2628	{
	2629	BspLeaf frontLeaf = dynamic_cast<BspLeaf >(front);
	2630	BspLeaf backLeaf = dynamic_cast<BspLeaf >(back);
	2631
	2632	//-- collapse tree
	2633	if (frontLeaf->GetViewCell() == backLeaf->GetViewCell())
	2634	{
	2635	BspViewCell *vc = frontLeaf->GetViewCell();
	2636
	2637	BspLeaf *leaf = new BspLeaf(interior->GetParent(), vc);
[489]	2638	leaf->SetTreeValid(frontLeaf->TreeValid());
[482]	2639
[479]	2640	// replace a link from node's parent
	2641	if (leaf->GetParent())
	2642	leaf->GetParent()->ReplaceChildLink(node, leaf);
[517]	2643	else
	2644	mRoot = leaf;
	2645
[501]	2646	++ collapsed;
[479]	2647	delete interior;
	2648
	2649	return leaf;
	2650	}
	2651	}
[590]	2652	#endif
[495]	2653	return node;
	2654	}
	2655
	2656
[501]	2657	int VspBspTree::CollapseTree()
[495]	2658	{
[501]	2659	int collapsed = 0;
[580]	2660	//TODO
	2661	#if VC_HISTORY
[501]	2662	(void)CollapseTree(mRoot, collapsed);
[517]	2663
[485]	2664	// revalidate leaves
[517]	2665	RepairViewCellsLeafLists();
[580]	2666	#endif
[501]	2667	return collapsed;
[479]	2668	}
	2669
	2670
[517]	2671	void VspBspTree::RepairViewCellsLeafLists()
[492]	2672	{
[590]	2673	// TODO
	2674	#if VC_HISTORY
[479]	2675	// list not valid anymore => clear
[492]	2676	stack<BspNode *> nodeStack;
	2677	nodeStack.push(mRoot);
	2678
	2679	ViewCell::NewMail();
	2680
	2681	while (!nodeStack.empty())
	2682	{
	2683	BspNode *node = nodeStack.top();
	2684	nodeStack.pop();
	2685
	2686	if (node->IsLeaf())
	2687	{
	2688	BspLeaf leaf = dynamic_cast<BspLeaf >(node);
	2689
	2690	BspViewCell *viewCell = leaf->GetViewCell();
[590]	2691
[492]	2692	if (!viewCell->Mailed())
	2693	{
	2694	viewCell->mLeaves.clear();
	2695	viewCell->Mail();
	2696	}
[580]	2697
[492]	2698	viewCell->mLeaves.push_back(leaf);
[590]	2699
[492]	2700	}
	2701	else
	2702	{
	2703	BspInterior interior = dynamic_cast<BspInterior >(node);
	2704
	2705	nodeStack.push(interior->GetFront());
	2706	nodeStack.push(interior->GetBack());
	2707	}
[479]	2708	}
[590]	2709	// TODO
	2710	#endif
[479]	2711	}
	2712
	2713
[532]	2714	typedef pair<BspNode , BspNodeGeometry > bspNodePair;
	2715
[535]	2716
[532]	2717	int VspBspTree::CastBeam(Beam &beam)
	2718	{
	2719	stack<bspNodePair> nodeStack;
	2720	BspNodeGeometry *rgeom = new BspNodeGeometry();
	2721	ConstructGeometry(mRoot, *rgeom);
	2722
	2723	nodeStack.push(bspNodePair(mRoot, rgeom));
	2724
	2725	ViewCell::NewMail();
	2726
	2727	while (!nodeStack.empty())
	2728	{
	2729	BspNode *node = nodeStack.top().first;
	2730	BspNodeGeometry *geom = nodeStack.top().second;
	2731	nodeStack.pop();
	2732
	2733	AxisAlignedBox3 box;
[535]	2734	box.Initialize();
	2735	geom->IncludeInBox(box);
[532]	2736
[535]	2737	const int side = beam.ComputeIntersection(box);
[532]	2738
	2739	switch (side)
	2740	{
	2741	case -1:
[547]	2742	CollectViewCells(node, true, beam.mViewCells, true);
[532]	2743	break;
	2744	case 0:
[535]	2745
[532]	2746	if (node->IsLeaf())
	2747	{
[535]	2748	BspLeaf leaf = dynamic_cast<BspLeaf >(node);
	2749
[532]	2750	if (!leaf->GetViewCell()->Mailed() && leaf->TreeValid())
[535]	2751	{
	2752	leaf->GetViewCell()->Mail();
[532]	2753	beam.mViewCells.push_back(leaf->GetViewCell());
[535]	2754	}
[532]	2755	}
	2756	else
	2757	{
	2758	BspInterior interior = dynamic_cast<BspInterior >(node);
[535]	2759
[538]	2760	BspNode *first = interior->GetFront();
	2761	BspNode *second = interior->GetBack();
[535]	2762
	2763	BspNodeGeometry *firstGeom = new BspNodeGeometry();
	2764	BspNodeGeometry *secondGeom = new BspNodeGeometry();
	2765
[538]	2766	geom->SplitGeometry(*firstGeom,
	2767	*secondGeom,
	2768	interior->GetPlane(),
	2769	mBox,
	2770	mEpsilon);
[535]	2771
[532]	2772	// decide on the order of the nodes
	2773	if (DotProd(beam.mPlanes[0].mNormal,
	2774	interior->GetPlane().mNormal) > 0)
	2775	{
	2776	swap(first, second);
[535]	2777	swap(firstGeom, secondGeom);
[532]	2778	}
	2779
[535]	2780	nodeStack.push(bspNodePair(first, firstGeom));
	2781	nodeStack.push(bspNodePair(second, secondGeom));
[532]	2782	}
[535]	2783
[532]	2784	break;
[538]	2785	default:
[532]	2786	// default: cull
[538]	2787	break;
[532]	2788	}
[538]	2789
[532]	2790	DEL_PTR(geom);
[535]	2791
[532]	2792	}
	2793
[538]	2794	return (int)beam.mViewCells.size();
[532]	2795	}
	2796
	2797
[485]	2798	void VspBspTree::SetViewCellsManager(ViewCellsManager *vcm)
[478]	2799	{
[485]	2800	mViewCellsManager = vcm;
	2801	}
	2802
	2803
[580]	2804	int VspBspTree::CollectMergeCandidates(const vector<BspLeaf *> leaves,
	2805	vector<MergeCandidate> &candidates)
[485]	2806	{
[478]	2807	BspLeaf::NewMail();
[508]	2808
[478]	2809	vector<BspLeaf *>::const_iterator it, it_end = leaves.end();
	2810
[580]	2811	int numCandidates = 0;
[509]	2812
[478]	2813	// find merge candidates and push them into queue
	2814	for (it = leaves.begin(); it != it_end; ++ it)
	2815	{
	2816	BspLeaf leaf = it;
[485]	2817
	2818	// the same leaves must not be part of two merge candidates
	2819	leaf->Mail();
	2820	vector<BspLeaf *> neighbors;
	2821	FindNeighbors(leaf, neighbors, true);
[600]	2822	//FindApproximateNeighbors(leaf, neighbors, true);
[485]	2823	vector<BspLeaf *>::const_iterator nit, nit_end = neighbors.end();
	2824
	2825	// TODO: test if at least one ray goes from one leaf to the other
	2826	for (nit = neighbors.begin(); nit != nit_end; ++ nit)
[508]	2827	{
	2828	if ((*nit)->GetViewCell() != leaf->GetViewCell())
[509]	2829	{
[580]	2830	MergeCandidate mc(leaf->GetViewCell(), (*nit)->GetViewCell());
	2831	candidates.push_back(mc);
[564]	2832
[580]	2833	++ numCandidates;
	2834	if ((numCandidates % 1000) == 0)
[566]	2835	{
[580]	2836	cout << "collected " << numCandidates << " merge candidates" << endl;
[566]	2837	}
[509]	2838	}
[485]	2839	}
	2840	}
	2841
[580]	2842	Debug << "merge queue: " << (int)candidates.size() << endl;
	2843	Debug << "leaves in queue: " << numCandidates << endl;
	2844
[508]	2845
[485]	2846	return (int)leaves.size();
	2847	}
	2848
	2849
[580]	2850	int VspBspTree::CollectMergeCandidates(const VssRayContainer &rays,
	2851	vector<MergeCandidate> &candidates)
[485]	2852	{
[547]	2853	ViewCell::NewMail();
[503]	2854	long startTime = GetTime();
[574]	2855
[485]	2856	map<BspLeaf , vector<BspLeaf> > neighborMap;
[574]	2857	ViewCellContainer::const_iterator iit;
[485]	2858
[503]	2859	int numLeaves = 0;
[485]	2860
	2861	BspLeaf::NewMail();
	2862
[574]	2863	for (int i = 0; i < (int)rays.size(); ++ i)
[485]	2864	{
[574]	2865	VssRay *ray = rays[i];
[547]	2866
[485]	2867	// traverse leaves stored in the rays and compare and
	2868	// merge consecutive leaves (i.e., the neighbors in the tree)
[574]	2869	if (ray->mViewCells.size() < 2)
[485]	2870	continue;
[580]	2871	//TODO viewcellhierarchy
[574]	2872	iit = ray->mViewCells.begin();
	2873	BspViewCell bspVc = dynamic_cast<BspViewCell >(*(iit ++));
[580]	2874	BspLeaf *leaf = bspVc->mLeaf;
[485]	2875
	2876	// traverse intersections
[489]	2877	// consecutive leaves are neighbors => add them to queue
[574]	2878	for (; iit != ray->mViewCells.end(); ++ iit)
[485]	2879	{
[489]	2880	// next pair
	2881	BspLeaf *prevLeaf = leaf;
[574]	2882	bspVc = dynamic_cast<BspViewCell >(iit);
[580]	2883	leaf = bspVc->mLeaf;
[489]	2884
[508]	2885	// view space not valid or same view cell
	2886	if (!leaf->TreeValid() \|\| !prevLeaf->TreeValid() \|\|
	2887	(leaf->GetViewCell() == prevLeaf->GetViewCell()))
[489]	2888	continue;
	2889
[580]	2890	vector<BspLeaf *> &neighbors = neighborMap[leaf];
[485]	2891
	2892	bool found = false;
[478]	2893
[485]	2894	// both leaves inserted in queue already =>
	2895	// look if double pair already exists
	2896	if (leaf->Mailed() && prevLeaf->Mailed())
[478]	2897	{
[485]	2898	vector<BspLeaf *>::const_iterator it, it_end = neighbors.end();
	2899
	2900	for (it = neighbors.begin(); !found && (it != it_end); ++ it)
	2901	if (*it == prevLeaf)
	2902	found = true; // already in queue
	2903	}
[547]	2904
[485]	2905	if (!found)
	2906	{
[564]	2907	// this pair is not in map yet
[485]	2908	// => insert into the neighbor map and the queue
	2909	neighbors.push_back(prevLeaf);
	2910	neighborMap[prevLeaf].push_back(leaf);
[478]	2911
[485]	2912	leaf->Mail();
	2913	prevLeaf->Mail();
[547]	2914
[580]	2915	MergeCandidate mc(leaf->GetViewCell(), prevLeaf->GetViewCell());
	2916
	2917	candidates.push_back(mc);
[564]	2918
[580]	2919	if (((int)candidates.size() % 1000) == 0)
[564]	2920	{
[580]	2921	cout << "collected " << (int)candidates.size() << " merge candidates" << endl;
[564]	2922	}
[478]	2923	}
[485]	2924	}
	2925	}
[564]	2926
[485]	2927	Debug << "neighbormap size: " << (int)neighborMap.size() << endl;
[580]	2928	Debug << "merge queue: " << (int)candidates.size() << endl;
[503]	2929	Debug << "leaves in queue: " << numLeaves << endl;
[485]	2930
[580]	2931
[503]	2932	//-- collect the leaves which haven't been found by ray casting
[542]	2933	if (0)
	2934	{
[551]	2935	cout << "finding additional merge candidates using geometry" << endl;
[542]	2936	vector<BspLeaf *> leaves;
[547]	2937	CollectLeaves(leaves, true);
[542]	2938	Debug << "found " << (int)leaves.size() << " new leaves" << endl << endl;
[580]	2939	CollectMergeCandidates(leaves, candidates);
[542]	2940	}
[503]	2941
	2942	return numLeaves;
[485]	2943	}
	2944
	2945
	2946
	2947
[508]	2948	ViewCell *VspBspTree::GetViewCell(const Vector3 &point)
[492]	2949	{
	2950	if (mRoot == NULL)
	2951	return NULL;
	2952
	2953	stack<BspNode *> nodeStack;
	2954	nodeStack.push(mRoot);
	2955
	2956	ViewCell *viewcell = NULL;
	2957
	2958	while (!nodeStack.empty()) {
	2959	BspNode *node = nodeStack.top();
	2960	nodeStack.pop();
	2961
	2962	if (node->IsLeaf()) {
	2963	viewcell = dynamic_cast<BspLeaf *>(node)->GetViewCell();
	2964	break;
	2965	} else {
	2966
	2967	BspInterior interior = dynamic_cast<BspInterior >(node);
	2968
	2969	// random decision
	2970	if (interior->GetPlane().Side(point) < 0)
	2971	nodeStack.push(interior->GetBack());
	2972	else
	2973	nodeStack.push(interior->GetFront());
	2974	}
	2975	}
	2976
	2977	return viewcell;
	2978	}
	2979
	2980
[487]	2981	bool VspBspTree::ViewPointValid(const Vector3 &viewPoint) const
	2982	{
	2983	BspNode *node = mRoot;
[485]	2984
[487]	2985	while (1)
	2986	{
	2987	// early exit
	2988	if (node->TreeValid())
	2989	return true;
	2990
	2991	if (node->IsLeaf())
	2992	return false;
	2993
	2994	BspInterior in = dynamic_cast<BspInterior >(node);
[490]	2995
	2996	if (in->GetPlane().Side(viewPoint) <= 0)
[487]	2997	{
	2998	node = in->GetBack();
	2999	}
	3000	else
	3001	{
	3002	node = in->GetFront();
	3003	}
	3004	}
	3005
	3006	// should never come here
	3007	return false;
	3008	}
	3009
	3010
	3011	void VspBspTree::PropagateUpValidity(BspNode *node)
	3012	{
[574]	3013	const bool isValid = node->TreeValid();
	3014
	3015	// propagative up invalid flag until only invalid nodes exist over this node
	3016	if (!isValid)
[487]	3017	{
[574]	3018	while (!node->IsRoot() && node->GetParent()->TreeValid())
	3019	{
	3020	node = node->GetParent();
	3021	node->SetTreeValid(false);
	3022	}
[487]	3023	}
[574]	3024	else
	3025	{
	3026	// propagative up valid flag until one of the subtrees is invalid
	3027	while (!node->IsRoot() && !node->TreeValid())
	3028	{
	3029	node = node->GetParent();
	3030	BspInterior interior = dynamic_cast<BspInterior >(node);
	3031
	3032	// the parent is valid iff both leaves are valid
	3033	node->SetTreeValid(interior->GetBack()->TreeValid() &&
	3034	interior->GetFront()->TreeValid());
	3035	}
	3036	}
[487]	3037	}
	3038
[508]	3039
	3040	bool VspBspTree::Export(ofstream &stream)
[503]	3041	{
[508]	3042	ExportNode(mRoot, stream);
[487]	3043
[503]	3044	return true;
	3045	}
	3046
	3047
[508]	3048	void VspBspTree::ExportNode(BspNode *node, ofstream &stream)
	3049	{
	3050	if (node->IsLeaf())
[503]	3051	{
[508]	3052	BspLeaf leaf = dynamic_cast<BspLeaf >(node);
[590]	3053	ViewCell *viewCell = mViewCellsTree->GetActiveViewCell(leaf->GetViewCell());
	3054
[508]	3055	int id = -1;
[590]	3056	if (viewCell != mOutOfBoundsCell)
	3057	id = viewCell->GetId();
[503]	3058
[508]	3059	stream << "<Leaf viewCellId=\"" << id << "\" />" << endl;
[503]	3060	}
[508]	3061	else
[503]	3062	{
[508]	3063	BspInterior interior = dynamic_cast<BspInterior >(node);
	3064
	3065	Plane3 plane = interior->GetPlane();
	3066	stream << "<Interior plane=\"" << plane.mNormal.x << " "
	3067	<< plane.mNormal.y << " " << plane.mNormal.z << " "
	3068	<< plane.mD << "\">" << endl;
[503]	3069
[508]	3070	ExportNode(interior->GetBack(), stream);
	3071	ExportNode(interior->GetFront(), stream);
[503]	3072
[508]	3073	stream << "</Interior>" << endl;
[503]	3074	}
	3075	}

Note: See TracBrowser for help on using the repository browser.

Download in other formats: