Context Navigation

source: GTP/trunk/Lib/Vis/Preprocessing/src/VspBspTree.cpp @ 1563

Revision 1563, 98.9 KB checked in by mattausch, 18 years ago (diff)
fixed bug with view space box

Line
1	#include <stack>
2	#include <time.h>
3	#include <iomanip>
4
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"
17	#include "ViewCellsManager.h"
18	#include "Beam.h"
19
20
21
22	namespace GtpVisibilityPreprocessor {
23
24
25	#define USE_FIXEDPOINT_T 0
26	#define COUNT_ORIGIN_OBJECTS 1
27
28
29	//-- static members
30
31	int VspBspTree::sFrontId = 0;
32	int VspBspTree::sBackId = 0;
33	int VspBspTree::sFrontAndBackId = 0;
34
35
36
37	typedef pair<BspNode , BspNodeGeometry > bspNodePair;
38
39
40	// pvs penalty can be different from pvs size
41	inline static 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;
50
51	return (float)pvs;
52	}
53
54
55
56
57	/******************************************************************************/
58	/* class VspBspTree implementation */
59	/******************************************************************************/
60
61
62	VspBspTree::VspBspTree():
63	mRoot(NULL),
64	mUseAreaForPvs(false),
65	mCostNormalizer(Limits::Small),
66	mViewCellsManager(NULL),
67	mOutOfBoundsCell(NULL),
68	mStoreRays(false),
69	mRenderCostWeight(0.5),
70	mUseRandomAxis(false),
71	mTimeStamp(1)
72	{
73	bool randomize = false;
74	Environment::GetSingleton()->GetBoolValue("VspBspTree.Construction.randomize", randomize);
75	if (randomize) Randomize(); // initialise random generator for heuristics
76
77	//////////////////
78	//-- termination criteria for autopartition
79
80	Environment::GetSingleton()->GetIntValue("VspBspTree.Termination.maxDepth", mTermMaxDepth);
81	Environment::GetSingleton()->GetIntValue("VspBspTree.Termination.minPvs", mTermMinPvs);
82	Environment::GetSingleton()->GetIntValue("VspBspTree.Termination.minRays", mTermMinRays);
83	Environment::GetSingleton()->GetFloatValue("VspBspTree.Termination.minProbability", mTermMinProbability);
84	Environment::GetSingleton()->GetFloatValue("VspBspTree.Termination.maxRayContribution", mTermMaxRayContribution);
85	Environment::GetSingleton()->GetFloatValue("VspBspTree.Termination.minAccRayLenght", mTermMinAccRayLength);
86
87	Environment::GetSingleton()->GetIntValue("VspBspTree.Termination.missTolerance", mTermMissTolerance);
88	Environment::GetSingleton()->GetIntValue("VspBspTree.Termination.maxViewCells", mMaxViewCells);
89
90	////////////////////////
91	//-- cost ratios for early tree termination
92	Environment::GetSingleton()->GetFloatValue("VspBspTree.Termination.maxCostRatio", mTermMaxCostRatio);
93	Environment::GetSingleton()->GetFloatValue("VspBspTree.Termination.minGlobalCostRatio", mTermMinGlobalCostRatio);
94	Environment::GetSingleton()->GetIntValue("VspBspTree.Termination.globalCostMissTolerance", mTermGlobalCostMissTolerance);
95
96	///////////
97	//-- factors for bsp tree split plane heuristics
98
99	Environment::GetSingleton()->GetFloatValue("VspBspTree.Factor.pvs", mPvsFactor);
100	Environment::GetSingleton()->GetFloatValue("VspBspTree.Termination.ct_div_ci", mCtDivCi);
101
102	//////////
103	//-- partition criteria
104
105	Environment::GetSingleton()->GetIntValue("VspBspTree.maxPolyCandidates", mMaxPolyCandidates);
106	Environment::GetSingleton()->GetIntValue("VspBspTree.maxRayCandidates", mMaxRayCandidates);
107	Environment::GetSingleton()->GetIntValue("VspBspTree.splitPlaneStrategy", mSplitPlaneStrategy);
108
109	Environment::GetSingleton()->GetFloatValue("VspBspTree.Construction.epsilon", mEpsilon);
110	Environment::GetSingleton()->GetIntValue("VspBspTree.maxTests", mMaxTests);
111
112	// if only the driving axis is used for axis aligned split
113	Environment::GetSingleton()->GetBoolValue("VspBspTree.splitUseOnlyDrivingAxis", mOnlyDrivingAxis);
114
115	//////////////////////
116	//-- termination criteria for axis aligned split
117	Environment::GetSingleton()->GetFloatValue("VspBspTree.Termination.AxisAligned.maxRayContribution",
118	mTermMaxRayContriForAxisAligned);
119	Environment::GetSingleton()->GetIntValue("VspBspTree.Termination.AxisAligned.minRays",
120	mTermMinRaysForAxisAligned);
121
122	Environment::GetSingleton()->GetFloatValue("VspBspTree.maxStaticMemory", mMaxMemory);
123
124	Environment::GetSingleton()->GetFloatValue("VspBspTree.Construction.renderCostWeight", mRenderCostWeight);
125	Environment::GetSingleton()->GetFloatValue("VspBspTree.Construction.renderCostDecreaseWeight", mRenderCostDecreaseWeight);
126	Environment::GetSingleton()->GetBoolValue("VspBspTree.usePolygonSplitIfAvailable", mUsePolygonSplitIfAvailable);
127
128	Environment::GetSingleton()->GetBoolValue("VspBspTree.useCostHeuristics", mUseCostHeuristics);
129	Environment::GetSingleton()->GetBoolValue("VspBspTree.useSplitCostQueue", mUseSplitCostQueue);
130	Environment::GetSingleton()->GetBoolValue("VspBspTree.simulateOctree", mCirculatingAxis);
131	Environment::GetSingleton()->GetBoolValue("VspBspTree.useRandomAxis", mUseRandomAxis);
132	Environment::GetSingleton()->GetIntValue("VspBspTree.nodePriorityQueueType", mNodePriorityQueueType);
133
134
135	char subdivisionStatsLog[100];
136	Environment::GetSingleton()->GetStringValue("VspBspTree.subdivisionStats", subdivisionStatsLog);
137	mSubdivisionStats.open(subdivisionStatsLog);
138
139	Environment::GetSingleton()->GetFloatValue("VspBspTree.Construction.minBand", mMinBand);
140	Environment::GetSingleton()->GetFloatValue("VspBspTree.Construction.maxBand", mMaxBand);
141	Environment::GetSingleton()->GetBoolValue("VspBspTree.Construction.useDrivingAxisForMaxCost", mUseDrivingAxisIfMaxCostViolated);
142
143	/////////
144	//-- debug output
145
146	Debug << "***** VSP BSP options ****** " << endl;
147	Debug << "max depth: " << mTermMaxDepth << endl;
148	Debug << "min PVS: " << mTermMinPvs << endl;
149	Debug << "min probabiliy: " << mTermMinProbability << endl;
150	Debug << "min rays: " << mTermMinRays << endl;
151	Debug << "max ray contri: " << mTermMaxRayContribution << endl;
152	Debug << "max cost ratio: " << mTermMaxCostRatio << endl;
153	Debug << "miss tolerance: " << mTermMissTolerance << endl;
154	Debug << "max view cells: " << mMaxViewCells << endl;
155	Debug << "max polygon candidates: " << mMaxPolyCandidates << endl;
156	Debug << "randomize: " << randomize << endl;
157
158	Debug << "using area for pvs: " << mUseAreaForPvs << endl;
159	Debug << "render cost weight: " << mRenderCostWeight << endl;
160	Debug << "min global cost ratio: " << mTermMinGlobalCostRatio << endl;
161	Debug << "global cost miss tolerance: " << mTermGlobalCostMissTolerance << endl;
162	Debug << "only driving axis: " << mOnlyDrivingAxis << endl;
163	Debug << "max memory: " << mMaxMemory << endl;
164	Debug << "use poly split if available: " << mUsePolygonSplitIfAvailable << endl;
165	Debug << "use cost heuristics: " << mUseCostHeuristics << endl;
166	Debug << "use split cost queue: " << mUseSplitCostQueue << endl;
167	Debug << "subdivision stats log: " << subdivisionStatsLog << endl;
168	Debug << "use random axis: " << mUseRandomAxis << endl;
169	Debug << "priority queue type: " << mNodePriorityQueueType << endl;
170	Debug << "circulating axis: " << mCirculatingAxis << endl;
171	Debug << "minband: " << mMinBand << endl;
172	Debug << "maxband: " << mMaxBand << endl;
173	Debug << "use driving axis for max cost: " << mUseDrivingAxisIfMaxCostViolated << endl;
174	Debug << "render cost decrease weight: " << mRenderCostDecreaseWeight << endl;
175
176	Debug << "Split plane strategy: ";
177
178	if (mSplitPlaneStrategy & RANDOM_POLYGON)
179	{
180	Debug << "random polygon ";
181	}
182	if (mSplitPlaneStrategy & AXIS_ALIGNED)
183	{
184	Debug << "axis aligned ";
185	}
186	if (mSplitPlaneStrategy & LEAST_RAY_SPLITS)
187	{
188	mCostNormalizer += mLeastRaySplitsFactor;
189	Debug << "least ray splits ";
190	}
191	if (mSplitPlaneStrategy & BALANCED_RAYS)
192	{
193	mCostNormalizer += mBalancedRaysFactor;
194	Debug << "balanced rays ";
195	}
196	if (mSplitPlaneStrategy & PVS)
197	{
198	mCostNormalizer += mPvsFactor;
199	Debug << "pvs";
200	}
201
202	Debug << endl;
203
204	mLocalSubdivisionCandidates = new vector<SortableEntry>;
205	}
206
207
208	BspViewCell *VspBspTree::GetOutOfBoundsCell()
209	{
210	return mOutOfBoundsCell;
211	}
212
213
214	BspViewCell *VspBspTree::GetOrCreateOutOfBoundsCell()
215	{
216	if (!mOutOfBoundsCell)
217	{
218	mOutOfBoundsCell = new BspViewCell();
219	mOutOfBoundsCell->SetId(OUT_OF_BOUNDS_ID);
220	mOutOfBoundsCell->SetValid(false);
221	}
222
223	return mOutOfBoundsCell;
224	}
225
226
227	const BspTreeStatistics &VspBspTree::GetStatistics() const
228	{
229	return mBspStats;
230	}
231
232
233	VspBspTree::~VspBspTree()
234	{
235	DEL_PTR(mRoot);
236	DEL_PTR(mLocalSubdivisionCandidates);
237	}
238
239
240	int VspBspTree::AddMeshToPolygons(Mesh *mesh,
241	PolygonContainer &polys,
242	MeshInstance *parent)
243	{
244	FaceContainer::const_iterator fi;
245
246	// copy the face data to polygons
247	for (fi = mesh->mFaces.begin(); fi != mesh->mFaces.end(); ++ fi)
248	{
249	Polygon3 poly = new Polygon3((fi), mesh);
250
251	if (poly->Valid(mEpsilon))
252	{
253	poly->mParent = parent; // set parent intersectable
254	polys.push_back(poly);
255	}
256	else
257	{
258	DEL_PTR(poly);
259	}
260	}
261	return (int)mesh->mFaces.size();
262	}
263
264
265	int VspBspTree::AddToPolygonSoup(const ViewCellContainer &viewCells,
266	PolygonContainer &polys,
267	int maxObjects)
268	{
269	int limit = (maxObjects > 0) ?
270	Min((int)viewCells.size(), maxObjects) : (int)viewCells.size();
271
272	int polysSize = 0;
273
274	for (int i = 0; i < limit; ++ i)
275	{
276	Mesh *mesh = viewCells[i]->GetMesh();
277	if (mesh)
278	{ // // copy the mesh into polygons and add to BSP tree aabb
279	mBoundingBox.Include(viewCells[i]->GetBox());
280	polysSize += AddMeshToPolygons(mesh,
281	polys,
282	viewCells[i]);
283	}
284	}
285
286	return polysSize;
287	}
288
289
290	int VspBspTree::AddToPolygonSoup(const ObjectContainer &objects,
291	PolygonContainer &polys,
292	int maxObjects)
293	{
294	int limit = (maxObjects > 0) ?
295	Min((int)objects.size(), maxObjects) : (int)objects.size();
296
297	for (int i = 0; i < limit; ++i)
298	{
299	Intersectable object = objects[i];//it;
300	Mesh *mesh = NULL;
301
302	switch (object->Type()) // extract the meshes
303	{
304	case Intersectable::MESH_INSTANCE:
305	mesh = dynamic_cast<MeshInstance *>(object)->GetMesh();
306	break;
307	case Intersectable::VIEW_CELL:
308	mesh = dynamic_cast<ViewCell *>(object)->GetMesh();
309	break;
310	case Intersectable::TRANSFORMED_MESH_INSTANCE:
311	{
312	TransformedMeshInstance mi = dynamic_cast<TransformedMeshInstance >(object);
313
314	if (!mi->GetMesh())
315	break;
316	mesh = new Mesh();
317	mi->GetTransformedMesh(*mesh);
318
319	break;
320	}
321	default:
322	Debug << "intersectable type not supported" << endl;
323	break;
324	}
325
326	if (mesh) // copy the mesh data to polygons
327	{
328	mBoundingBox.Include(object->GetBox()); // add to BSP tree aabb
329	AddMeshToPolygons(mesh, polys, NULL);
330
331	// cleanup
332	if (object->Type() == Intersectable::TRANSFORMED_MESH_INSTANCE)
333	DEL_PTR(mesh);
334	}
335	}
336
337	return (int)polys.size();
338	}
339
340
341	void VspBspTree::ComputeBoundingBox(const VssRayContainer &sampleRays,
342	AxisAlignedBox3 *forcedBoundingBox)
343	{
344	if (forcedBoundingBox)
345	{
346	mBoundingBox = *forcedBoundingBox;
347	}
348	else // compute vsp tree bounding box
349	{
350	mBoundingBox.Initialize();
351
352	VssRayContainer::const_iterator rit, rit_end = sampleRays.end();
353
354	//-- compute bounding box
355	for (rit = sampleRays.begin(); rit != rit_end; ++ rit)
356	{
357	VssRay ray = rit;
358
359	// compute bounding box of view space
360	mBoundingBox.Include(ray->GetTermination());
361	mBoundingBox.Include(ray->GetOrigin());
362	}
363	}
364	}
365
366
367	void VspBspTree::Construct(const VssRayContainer &sampleRays,
368	AxisAlignedBox3 *forcedBoundingBox)
369	{
370	// Compute the bounding box from the rays
371	ComputeBoundingBox(sampleRays, forcedBoundingBox);
372
373	PolygonContainer polys;
374	RayInfoContainer *rays = new RayInfoContainer();
375
376	//-- extract polygons from rays if there are polygon candidates
377	if (mMaxPolyCandidates)
378	{
379	int numObj = 0;
380
381	Intersectable::NewMail();
382
383	cout << "Extracting polygons from rays ... ";
384
385	long startTime = GetTime();
386
387	VssRayContainer::const_iterator rit, rit_end = sampleRays.end();
388
389	//-- extract polygons intersected by the rays
390	for (rit = sampleRays.begin(); rit != rit_end; ++ rit)
391	{
392	VssRay ray = rit;
393	Intersectable *obj = ray->mTerminationObject;
394
395	if ((mBoundingBox.IsInside(ray->mTermination) \|\| !forcedBoundingBox) &&
396	obj && !obj->Mailed())
397	{
398	obj->Mail();
399
400	// transformed mesh instance and mesh instance handle mesh differently
401	if (obj->Type() == Intersectable::TRANSFORMED_MESH_INSTANCE)
402	{
403	Mesh mesh;
404
405	TransformedMeshInstance *tmobj =
406	dynamic_cast<TransformedMeshInstance *>(obj);
407
408	tmobj->GetTransformedMesh(mesh);
409	AddMeshToPolygons(&mesh, polys, tmobj);
410	}
411	else // MeshInstance
412	{
413	MeshInstance mobj = dynamic_cast<MeshInstance >(obj);
414	AddMeshToPolygons(mobj->GetMesh(), polys, mobj);
415	}
416
417	++ numObj;
418
419	//-- compute bounding box
420	if (!forcedBoundingBox)
421	mBoundingBox.Include(ray->mTermination);
422	}
423
424	if ((mBoundingBox.IsInside(ray->mOrigin) \|\| !forcedBoundingBox) &&
425	ray->mOriginObject &&
426	!ray->mOriginObject->Mailed())
427	{
428	ray->mOriginObject->Mail();
429	MeshInstance obj = dynamic_cast<MeshInstance >(ray->mOriginObject);
430	AddMeshToPolygons(obj->GetMesh(), polys, obj);
431
432	++ numObj;
433	}
434	}
435
436	// throw out unnecessary polygons
437	PreprocessPolygons(polys);
438
439	cout << "finished" << endl;
440
441	Debug << "\n" << (int)polys.size() << " polys extracted from "
442	<< (int)sampleRays.size() << " rays in "
443	<< TimeDiff(startTime, GetTime())*1e-3 << " secs" << endl << endl;
444	}
445
446	Debug << "maximal pvs (i.e., pvs still considered as valid): "
447	<< mViewCellsManager->GetMaxPvsSize() << endl;
448
449	VssRayContainer::const_iterator rit, rit_end = sampleRays.end();
450
451	//-- store rays
452	for (rit = sampleRays.begin(); rit != rit_end; ++ rit)
453	{
454	VssRay ray = rit;
455
456	float minT, maxT;
457
458	static Ray hray;
459	hray.Init(*ray);
460
461	// TODO: not very efficient to implictly cast between rays types
462	if (mBoundingBox.GetRaySegment(hray, minT, maxT))
463	{
464	float len = ray->Length();
465
466	if (!len)
467	len = Limits::Small;
468
469	rays->push_back(RayInfo(ray, minT / len, maxT / len));
470	}
471	}
472
473	// normalize
474	if (mUseAreaForPvs)
475	mTermMinProbability *= mBoundingBox.SurfaceArea();
476	else
477	mTermMinProbability *= mBoundingBox.GetVolume();
478
479
480	mBspStats.nodes = 1;
481	mBspStats.polys = (int)polys.size();
482	mBspStats.mGlobalCostMisses = 0;
483
484
485	// use split cost priority queue
486	if (mUseSplitCostQueue)
487	{
488	ConstructWithSplitQueue(polys, rays);
489	}
490	else
491	{
492	Construct(polys, rays);
493	}
494
495	// clean up polygons
496	CLEAR_CONTAINER(polys);
497	}
498
499
500	// TODO: return memory usage in MB
501	float VspBspTree::GetMemUsage() const
502	{
503	return (float)
504	(sizeof(VspBspTree) +
505	mBspStats.Leaves() * sizeof(BspLeaf) +
506	mCreatedViewCells * sizeof(BspViewCell) +
507	mBspStats.pvs * sizeof(PvsData) +
508	mBspStats.Interior() * sizeof(BspInterior) +
509	mBspStats.accumRays * sizeof(RayInfo)) / (1024.0f * 1024.0f);
510	}
511
512
513	void VspBspTree::Construct(const PolygonContainer &polys, RayInfoContainer *rays)
514	{
515	VspBspTraversalQueue tQueue;
516
517	/// create new vsp tree
518	mRoot = new BspLeaf();
519
520	// constrruct root node geometry
521	BspNodeGeometry *geom = new BspNodeGeometry();
522	ConstructGeometry(mRoot, *geom);
523
524	/// we use the overall probability as normalizer
525	/// either the overall area or the volume
526	const float prop = mUseAreaForPvs ? geom->GetArea() : geom->GetVolume();
527
528	/// first traversal data
529	VspBspTraversalData tData(mRoot,
530	new PolygonContainer(polys),
531	0,
532	rays,
533	ComputePvsSize(*rays),
534	prop,
535	geom);
536
537	// evaluate the priority of this traversal data
538	EvalPriority(tData);
539
540	// first node is kd node, i.e. an axis aligned box
541	if (1)
542	tData.mIsKdNode = true;
543	else
544	tData.mIsKdNode = false;
545
546	tQueue.push(tData);
547
548
549	mTotalCost = tData.mPvs * tData.mProbability / mBoundingBox.GetVolume();
550	mTotalPvsSize = tData.mPvs;
551
552	// first subdivison statistics
553	AddSubdivisionStats(1, 0, 0, mTotalCost, (float)mTotalPvsSize);
554
555	mBspStats.Start();
556	cout << "Constructing vsp bsp tree ... \n";
557
558	const long startTime = GetTime();
559	// used for intermediate time measurements and progress
560	long interTime = GetTime();
561
562	int nLeaves = 500;
563	int nViewCells = 500;
564
565	mOutOfMemory = false;
566	mCreatedViewCells = 0;
567
568	while (!tQueue.empty())
569	{
570	tData = tQueue.top();
571	tQueue.pop();
572
573	if (0 && !mOutOfMemory)
574	{
575	float mem = GetMemUsage();
576
577	if (mem > mMaxMemory)
578	{
579	mOutOfMemory = true;
580	Debug << "memory limit reached: " << mem << endl;
581	}
582	}
583
584	// subdivide leaf node
585	const BspNode *r = Subdivide(tQueue, tData);
586
587	if (r == mRoot)
588	Debug << "VSP BSP tree construction time spent at root: "
589	<< TimeDiff(startTime, GetTime())*1e-3 << "s" << endl;
590
591	if (mBspStats.Leaves() >= nLeaves)
592	{
593	nLeaves += 500;
594
595	cout << "leaves=" << mBspStats.Leaves() << endl;
596	Debug << "needed "
597	<< TimeDiff(interTime, GetTime())*1e-3
598	<< " secs to create 500 view cells" << endl;
599	interTime = GetTime();
600	}
601
602	if (mCreatedViewCells >= nViewCells)
603	{
604	nViewCells += 500;
605
606	cout << "generated " << mCreatedViewCells << " viewcells" << endl;
607	}
608	}
609
610	Debug << "Used Memory: " << GetMemUsage() << " MB" << endl << endl;
611	cout << "finished in " << TimeDiff(startTime, GetTime())*1e-3 << "secs" << endl;
612
613	mBspStats.Stop();
614	}
615
616
617
618	void VspBspTree::ConstructWithSplitQueue(const PolygonContainer &polys,
619	RayInfoContainer *rays)
620	{
621	VspBspSplitQueue tQueue;
622
623	mRoot = new BspLeaf();
624
625	// constrruct root node geometry
626	BspNodeGeometry *geom = new BspNodeGeometry();
627	ConstructGeometry(mRoot, *geom);
628
629	const float prop = mUseAreaForPvs ? geom->GetArea() : geom->GetVolume();
630
631	VspBspTraversalData tData(mRoot,
632	new PolygonContainer(polys),
633	0,
634	rays,
635	ComputePvsSize(*rays),
636	prop,
637	geom);
638
639
640	// first node is kd node, i.e. an axis aligned box
641	if (1)
642	tData.mIsKdNode = true;
643	else
644	tData.mIsKdNode = false;
645
646	// compute first split candidate
647	VspBspSubdivisionCandidate splitCandidate;
648	splitCandidate.mParentData = tData;
649
650	EvalSubdivisionCandidate(splitCandidate);
651
652	tQueue.push(splitCandidate);
653
654	mTotalCost = tData.mPvs * tData.mProbability / mBoundingBox.GetVolume();
655	mTotalPvsSize = tData.mPvs;
656
657	// first subdivison statistics
658	AddSubdivisionStats(1, 0, 0, mTotalCost, (float)mTotalPvsSize);
659
660	mBspStats.Start();
661	cout << "Constructing vsp bsp tree ... \n";
662
663	long startTime = GetTime();
664	int nLeaves = 500;
665	int nViewCells = 500;
666
667	// used for intermediate time measurements and progress
668	long interTime = GetTime();
669
670	mOutOfMemory = false;
671
672	mCreatedViewCells = 0;
673
674	while (!tQueue.empty())
675	{
676	splitCandidate = tQueue.top();
677	tQueue.pop();
678
679	// cost ratio of cost decrease / totalCost
680	float costRatio = splitCandidate.mRenderCostDecr / mTotalCost;
681
682	//Debug << "cost ratio: " << costRatio << endl;
683	if (costRatio < mTermMinGlobalCostRatio)
684	{
685	++ mBspStats.mGlobalCostMisses;
686	}
687
688	if (0 && !mOutOfMemory)
689	{
690	float mem = GetMemUsage();
691	if (mem > mMaxMemory)
692	{
693	mOutOfMemory = true;
694	Debug << "memory limit reached: " << mem << endl;
695	}
696	}
697
698	// subdivide leaf node
699	BspNode *r = Subdivide(tQueue, splitCandidate);
700
701	if (r == mRoot)
702	{
703	Debug << "VSP BSP tree construction time spent at root: "
704	<< TimeDiff(startTime, GetTime())*1e-3 << "s" << endl;
705	}
706
707	if (mBspStats.Leaves() >= nLeaves)
708	{
709	nLeaves += 500;
710
711	cout << "leaves=" << mBspStats.Leaves() << endl;
712	Debug << "needed "
713	<< TimeDiff(interTime, GetTime())*1e-3
714	<< " secs to create 500 view cells" << endl;
715	interTime = GetTime();
716	}
717
718	if (mCreatedViewCells == nViewCells)
719	{
720	nViewCells += 500;
721	cout << "generated " << mCreatedViewCells << " viewcells" << endl;
722	}
723	}
724
725	Debug << "Used Memory: " << GetMemUsage() << " MB" << endl << endl;
726	cout << "finished\n";
727
728	mBspStats.Stop();
729	}
730
731
732	bool VspBspTree::LocalTerminationCriteriaMet(const VspBspTraversalData &data) const
733	{
734	return
735	(((int)data.mRays->size() <= mTermMinRays) \|\|
736	(data.mPvs <= mTermMinPvs) \|\|
737	(data.mProbability <= mTermMinProbability) \|\|
738	(data.GetAvgRayContribution() > mTermMaxRayContribution) \|\|
739	(data.mDepth >= mTermMaxDepth));
740	}
741
742
743	void VspBspTree::AddSubdivisionStats(const int viewCells,
744	const float renderCostDecr,
745	const float splitCandidateCost,
746	const float totalRenderCost,
747	const float avgRenderCost)
748	{
749	mSubdivisionStats
750	<< "#ViewCells\n" << viewCells << endl
751	<< "#RenderCostDecrease\n" << renderCostDecr << endl
752	<< "#SubdivisionCandidateCost\n" << splitCandidateCost << endl
753	<< "#TotalRenderCost\n" << totalRenderCost << endl
754	<< "#AvgRenderCost\n" << avgRenderCost << endl;
755	}
756
757
758	bool VspBspTree::GlobalTerminationCriteriaMet(const VspBspTraversalData &data) const
759	{
760	return
761	(0
762	\|\| mOutOfMemory
763	\|\| (mBspStats.Leaves() >= mMaxViewCells)
764	\|\| (mBspStats.mGlobalCostMisses >= mTermGlobalCostMissTolerance)
765	);
766	}
767
768
769	BspNode *VspBspTree::Subdivide(VspBspTraversalQueue &tQueue,
770	VspBspTraversalData &tData)
771	{
772	BspNode *newNode = tData.mNode;
773
774	if (!LocalTerminationCriteriaMet(tData) && !GlobalTerminationCriteriaMet(tData))
775	{
776	PolygonContainer coincident;
777
778	VspBspTraversalData tFrontData;
779	VspBspTraversalData tBackData;
780
781	// create new interior node and two leaf nodes
782	// or return leaf as it is (if maxCostRatio missed)
783	int splitAxis;
784	bool splitFurther = true;
785	int maxCostMisses = tData.mMaxCostMisses;
786
787	Plane3 splitPlane;
788	BspLeaf leaf = dynamic_cast<BspLeaf >(tData.mNode);
789
790	// choose next split plane
791	if (!SelectPlane(splitPlane, leaf, tData, tFrontData, tBackData, splitAxis))
792	{
793	++ maxCostMisses;
794
795	if (maxCostMisses > mTermMissTolerance)
796	{
797	// terminate branch because of max cost
798	++ mBspStats.maxCostNodes;
799	splitFurther = false;
800	}
801	}
802
803	// if this a valid split => subdivide this node further
804
805	if (splitFurther)
806	{
807	newNode = SubdivideNode(splitPlane, tData, tFrontData, tBackData, coincident);
808
809	if (splitAxis < 3)
810	++ mBspStats.splits[splitAxis];
811	else
812	++ mBspStats.polySplits;
813
814	// if it was a kd node (i.e., a box) and the split axis is axis aligned, it is still a kd node
815	tFrontData.mIsKdNode = tBackData.mIsKdNode = (tData.mIsKdNode && (splitAxis < 3));
816
817	tFrontData.mAxis = tBackData.mAxis = splitAxis;
818
819	// how often was max cost ratio missed in this branch?
820	tFrontData.mMaxCostMisses = maxCostMisses;
821	tBackData.mMaxCostMisses = maxCostMisses;
822
823	EvalPriority(tFrontData);
824	EvalPriority(tBackData);
825
826	// evaluate subdivision stats
827	if (1)
828	EvalSubdivisionStats(tData, tFrontData, tBackData);
829
830
831	// push the children on the stack
832	tQueue.push(tFrontData);
833	tQueue.push(tBackData);
834
835	// delete old leaf node
836	DEL_PTR(tData.mNode);
837	}
838	}
839
840	//-- terminate traversal and create new view cell
841	if (newNode->IsLeaf())
842	{
843	BspLeaf leaf = dynamic_cast<BspLeaf >(newNode);
844
845	BspViewCell *viewCell = new BspViewCell();
846	leaf->SetViewCell(viewCell);
847
848	//-- update pvs
849	int conSamp = 0;
850	float sampCon = 0.0f;
851	AddToPvs(leaf, *tData.mRays, sampCon, conSamp);
852
853	// update scalar pvs size lookup
854	ObjectPvs &pvs = viewCell->GetPvs();
855	mViewCellsManager->UpdateScalarPvsSize(viewCell, pvs.CountObjectsInPvs(), pvs.GetSize());
856
857
858	mBspStats.contributingSamples += conSamp;
859	mBspStats.sampleContributions += (int)sampCon;
860
861	//-- store additional info
862	if (mStoreRays)
863	{
864	RayInfoContainer::const_iterator it, it_end = tData.mRays->end();
865	for (it = tData.mRays->begin(); it != it_end; ++ it)
866	{
867	(*it).mRay->Ref();
868	leaf->mVssRays.push_back((*it).mRay);
869	}
870	}
871
872	// should I check here?
873	if (0 && !mViewCellsManager->CheckValidity(viewCell, 0,
874	mViewCellsManager->GetMaxPvsSize()))
875	{
876	viewCell->SetValid(false);
877	leaf->SetTreeValid(false);
878	PropagateUpValidity(leaf);
879
880	++ mBspStats.invalidLeaves;
881	}
882
883	viewCell->mLeaves.push_back(leaf);
884
885	if (mUseAreaForPvs)
886	viewCell->SetArea(tData.mProbability);
887	else
888	viewCell->SetVolume(tData.mProbability);
889
890	leaf->mProbability = tData.mProbability;
891
892	// finally evaluate stats until this leaf
893	if (0)
894	EvaluateLeafStats(tData);
895	}
896
897	//-- cleanup
898	tData.Clear();
899
900	return newNode;
901	}
902
903
904	// subdivide node using a split plane queue
905	BspNode *VspBspTree::Subdivide(VspBspSplitQueue &tQueue,
906	VspBspSubdivisionCandidate &splitCandidate)
907	{
908	VspBspTraversalData &tData = splitCandidate.mParentData;
909
910	BspNode *newNode = tData.mNode;
911
912	if (!LocalTerminationCriteriaMet(tData) && !GlobalTerminationCriteriaMet(tData))
913	{
914	PolygonContainer coincident;
915
916	VspBspTraversalData tFrontData;
917	VspBspTraversalData tBackData;
918
919	//-- continue subdivision
920
921	// create new interior node and two leaf node
922	const Plane3 splitPlane = splitCandidate.mSplitPlane;
923
924	newNode = SubdivideNode(splitPlane, tData, tFrontData, tBackData, coincident);
925
926	const int splitAxis = splitCandidate.mSplitAxis;
927	const int maxCostMisses = splitCandidate.mMaxCostMisses;
928
929	if (splitAxis < 3)
930	++ mBspStats.splits[splitAxis];
931	else
932	++ mBspStats.polySplits;
933
934	tFrontData.mIsKdNode = tBackData.mIsKdNode = (tData.mIsKdNode && (splitAxis < 3));
935	tFrontData.mAxis = tBackData.mAxis = splitAxis;
936
937	// how often was max cost ratio missed in this branch?
938	tFrontData.mMaxCostMisses = maxCostMisses;
939	tBackData.mMaxCostMisses = maxCostMisses;
940
941	// statistics
942	if (1)
943	{
944	float cFront = (float)tFrontData.mPvs * tFrontData.mProbability;
945	float cBack = (float)tBackData.mPvs * tBackData.mProbability;
946	float cData = (float)tData.mPvs * tData.mProbability;
947
948
949	float costDecr =
950	(cFront + cBack - cData) / mBoundingBox.GetVolume();
951
952	mTotalCost += costDecr;
953	mTotalPvsSize += tFrontData.mPvs + tBackData.mPvs - tData.mPvs;
954
955	AddSubdivisionStats(mBspStats.Leaves(),
956	-costDecr,
957	splitCandidate.GetPriority(),
958	mTotalCost,
959	(float)mTotalPvsSize / (float)mBspStats.Leaves());
960	}
961
962
963	//-- push the new split candidates on the stack
964	VspBspSubdivisionCandidate frontCandidate;
965	frontCandidate.mParentData = tFrontData;
966
967	VspBspSubdivisionCandidate backCandidate;
968	backCandidate.mParentData = tBackData;
969
970	EvalSubdivisionCandidate(frontCandidate);
971	EvalSubdivisionCandidate(backCandidate);
972
973	tQueue.push(frontCandidate);
974	tQueue.push(backCandidate);
975
976	// delete old leaf node
977	DEL_PTR(tData.mNode);
978	}
979
980
981	//////////////////
982	//-- terminate traversal and create new view cell
983	if (newNode->IsLeaf())
984	{
985	BspLeaf leaf = dynamic_cast<BspLeaf >(newNode);
986
987	BspViewCell *viewCell = new BspViewCell();
988	leaf->SetViewCell(viewCell);
989
990	//-- update pvs
991	int conSamp = 0;
992	float sampCon = 0.0f;
993	AddToPvs(leaf, *tData.mRays, sampCon, conSamp);
994
995	// update scalar pvs size value
996	ObjectPvs &pvs = viewCell->GetPvs();
997	mViewCellsManager->UpdateScalarPvsSize(viewCell, pvs.CountObjectsInPvs(), pvs.GetSize());
998
999	mBspStats.contributingSamples += conSamp;
1000	mBspStats.sampleContributions +=(int) sampCon;
1001
1002	viewCell->mLeaves.push_back(leaf);
1003
1004	///////////
1005	//-- store additional info
1006	if (mStoreRays)
1007	{
1008	RayInfoContainer::const_iterator it, it_end = tData.mRays->end();
1009	for (it = tData.mRays->begin(); it != it_end; ++ it)
1010	{
1011	(*it).mRay->Ref();
1012	leaf->mVssRays.push_back((*it).mRay);
1013	}
1014	}
1015
1016	if (mUseAreaForPvs)
1017	viewCell->SetArea(tData.mProbability);
1018	else
1019	viewCell->SetVolume(tData.mProbability);
1020
1021	leaf->mProbability = tData.mProbability;
1022
1023	// finally evaluate stats until this leaf
1024	if (0)
1025	EvaluateLeafStats(tData);
1026	}
1027
1028	//-- cleanup
1029	tData.Clear();
1030
1031	return newNode;
1032	}
1033
1034
1035	void VspBspTree::EvalPriority(VspBspTraversalData &tData) const
1036	{
1037	switch (mNodePriorityQueueType)
1038	{
1039	case BREATH_FIRST:
1040	tData.mPriority = (float)-tData.mDepth;
1041	break;
1042	case DEPTH_FIRST:
1043	tData.mPriority = (float)tData.mDepth;
1044	break;
1045	default:
1046	tData.mPriority = tData.mPvs * tData.mProbability;
1047	//Debug << "priority: " << tData.mPriority << endl;
1048	break;
1049	}
1050	}
1051
1052
1053	void VspBspTree::EvalSubdivisionCandidate(VspBspSubdivisionCandidate &splitCandidate)
1054	{
1055	VspBspTraversalData frontData;
1056	VspBspTraversalData backData;
1057
1058	BspLeaf leaf = dynamic_cast<BspLeaf >(splitCandidate.mParentData.mNode);
1059
1060	// compute locally best split plane
1061	const bool costRatioViolated =
1062	SelectPlane(splitCandidate.mSplitPlane,
1063	leaf,
1064	splitCandidate.mParentData,
1065	frontData,
1066	backData,
1067	splitCandidate.mSplitAxis);
1068
1069	// max cost threshold violated?
1070	splitCandidate.mMaxCostMisses = costRatioViolated ?
1071	splitCandidate.mParentData.mMaxCostMisses :
1072	splitCandidate.mParentData.mMaxCostMisses + 1;
1073
1074	float oldRenderCost;
1075
1076	// compute global decrease in render cost
1077	const float renderCostDecr = EvalRenderCostDecrease(splitCandidate.mSplitPlane,
1078	splitCandidate.mParentData,
1079	oldRenderCost);
1080
1081	splitCandidate.mRenderCostDecr = renderCostDecr;
1082
1083	// TODO: geometry could be reused
1084	delete frontData.mGeometry;
1085	delete backData.mGeometry;
1086
1087	// set priority for queue
1088	#if 0
1089	const float priority = (float)-data.mDepth;
1090	#else
1091
1092	// take render cost of node into account
1093	// otherwise danger of being stuck in a local minimum!!
1094	const float factor = mRenderCostDecreaseWeight;
1095	const float priority = factor * renderCostDecr + (1.0f - factor) * oldRenderCost;
1096	#endif
1097
1098	splitCandidate.mPriority = priority;
1099	}
1100
1101
1102	void VspBspTree::EvalSubdivisionStats(const VspBspTraversalData &tData,
1103	const VspBspTraversalData &tFrontData,
1104	const VspBspTraversalData &tBackData)
1105	{
1106	const float cFront = (float)tFrontData.mPvs * tFrontData.mProbability;
1107	const float cBack = (float)tBackData.mPvs * tBackData.mProbability;
1108	const float cData = (float)tData.mPvs * tData.mProbability;
1109
1110	const float costDecr =
1111	(cFront + cBack - cData) / mBoundingBox.GetVolume();
1112
1113	mTotalCost += costDecr;
1114	mTotalPvsSize += tFrontData.mPvs + tBackData.mPvs - tData.mPvs;
1115
1116	AddSubdivisionStats(mBspStats.Leaves(),
1117	-costDecr,
1118	0,
1119	mTotalCost,
1120	(float)mTotalPvsSize / (float)mBspStats.Leaves());
1121	}
1122
1123
1124	BspInterior *VspBspTree::SubdivideNode(const Plane3 &splitPlane,
1125	VspBspTraversalData &tData,
1126	VspBspTraversalData &frontData,
1127	VspBspTraversalData &backData,
1128	PolygonContainer &coincident)
1129	{
1130	BspLeaf leaf = dynamic_cast<BspLeaf >(tData.mNode);
1131
1132	//-- the front and back traversal data is filled with the new values
1133	frontData.mDepth = tData.mDepth + 1;
1134	frontData.mPolygons = new PolygonContainer();
1135	frontData.mRays = new RayInfoContainer();
1136
1137	backData.mDepth = tData.mDepth + 1;
1138	backData.mPolygons = new PolygonContainer();
1139	backData.mRays = new RayInfoContainer();
1140
1141
1142	//-- subdivide rays
1143	SplitRays(splitPlane,
1144	*tData.mRays,
1145	*frontData.mRays,
1146	*backData.mRays);
1147
1148
1149	// compute pvs
1150	frontData.mPvs = ComputePvsSize(*frontData.mRays);
1151	backData.mPvs = ComputePvsSize(*backData.mRays);
1152
1153	// split front and back node geometry and compute area
1154
1155	// if geometry was not already computed
1156	if (!frontData.mGeometry && !backData.mGeometry)
1157	{
1158	frontData.mGeometry = new BspNodeGeometry();
1159	backData.mGeometry = new BspNodeGeometry();
1160
1161	tData.mGeometry->SplitGeometry(*frontData.mGeometry,
1162	*backData.mGeometry,
1163	splitPlane,
1164	mBoundingBox,
1165	//0.0f);
1166	mEpsilon);
1167
1168	if (mUseAreaForPvs)
1169	{
1170	frontData.mProbability = frontData.mGeometry->GetArea();
1171	backData.mProbability = backData.mGeometry->GetArea();
1172	}
1173	else
1174	{
1175	frontData.mProbability = frontData.mGeometry->GetVolume();
1176	backData.mProbability = tData.mProbability - frontData.mProbability;
1177
1178	// should never come here: wrong volume !!!
1179	if (0)
1180	{
1181	if (frontData.mProbability < -0.00001)
1182	Debug << "fatal error f: " << frontData.mProbability << endl;
1183	if (backData.mProbability < -0.00001)
1184	Debug << "fatal error b: " << backData.mProbability << endl;
1185
1186	// clamp because of precision issues
1187	if (frontData.mProbability < 0) frontData.mProbability = 0;
1188	if (backData.mProbability < 0) backData.mProbability = 0;
1189	}
1190	}
1191	}
1192
1193
1194	// subdivide polygons
1195	SplitPolygons(splitPlane,
1196	*tData.mPolygons,
1197	*frontData.mPolygons,
1198	*backData.mPolygons,
1199	coincident);
1200
1201
1202
1203	///////////////////////////////////////
1204	// subdivide further
1205
1206	// store maximal and minimal depth
1207	if (tData.mDepth > mBspStats.maxDepth)
1208	{
1209	Debug << "max depth increases to " << tData.mDepth << " at " << mBspStats.Leaves() << " leaves" << endl;
1210	mBspStats.maxDepth = tData.mDepth;
1211	}
1212
1213	mBspStats.nodes += 2;
1214
1215
1216	BspInterior *interior = new BspInterior(splitPlane);
1217
1218	#ifdef _DEBUG
1219	Debug << interior << endl;
1220	#endif
1221
1222
1223	//-- create front and back leaf
1224
1225	BspInterior *parent = leaf->GetParent();
1226
1227	// replace a link from node's parent
1228	if (parent)
1229	{
1230	parent->ReplaceChildLink(leaf, interior);
1231	interior->SetParent(parent);
1232	}
1233	else // new root
1234	{
1235	mRoot = interior;
1236	}
1237
1238	// and setup child links
1239	interior->SetupChildLinks(new BspLeaf(interior), new BspLeaf(interior));
1240
1241	frontData.mNode = interior->GetFront();
1242	backData.mNode = interior->GetBack();
1243
1244	interior->mTimeStamp = mTimeStamp ++;
1245
1246
1247	//DEL_PTR(leaf);
1248	return interior;
1249	}
1250
1251
1252	void VspBspTree::AddToPvs(BspLeaf *leaf,
1253	const RayInfoContainer &rays,
1254	float &sampleContributions,
1255	int &contributingSamples)
1256	{
1257	sampleContributions = 0;
1258	contributingSamples = 0;
1259
1260	RayInfoContainer::const_iterator it, it_end = rays.end();
1261
1262	ViewCellLeaf *vc = leaf->GetViewCell();
1263
1264	// add contributions from samples to the PVS
1265	for (it = rays.begin(); it != it_end; ++ it)
1266	{
1267	float sc = 0.0f;
1268	VssRay ray = (it).mRay;
1269
1270	bool madeContrib = false;
1271	float contribution;
1272
1273	if (ray->mTerminationObject)
1274	{
1275	if (vc->AddPvsSample(ray->mTerminationObject, ray->mPdf, contribution))
1276	madeContrib = true;
1277	sc += contribution;
1278	}
1279
1280	// only count termination objects?
1281	if (COUNT_ORIGIN_OBJECTS && ray->mOriginObject)
1282	{
1283	if (vc->AddPvsSample(ray->mOriginObject, ray->mPdf, contribution))
1284	madeContrib = true;
1285
1286	sc += contribution;
1287	}
1288
1289	sampleContributions += sc;
1290
1291	if (madeContrib)
1292	++ contributingSamples;
1293	}
1294	}
1295
1296
1297	void VspBspTree::SortSubdivisionCandidates(const RayInfoContainer &rays,
1298	const int axis,
1299	float minBand,
1300	float maxBand)
1301	{
1302	mLocalSubdivisionCandidates->clear();
1303
1304	int requestedSize = 2 * (int)(rays.size());
1305	// creates a sorted split candidates array
1306	if (mLocalSubdivisionCandidates->capacity() > 500000 &&
1307	requestedSize < (int)(mLocalSubdivisionCandidates->capacity() / 10) )
1308	{
1309	delete mLocalSubdivisionCandidates;
1310	mLocalSubdivisionCandidates = new vector<SortableEntry>;
1311	}
1312
1313	mLocalSubdivisionCandidates->reserve(requestedSize);
1314
1315	if (0)
1316	{ // float values => don't compare with exact values
1317	minBand += Limits::Small;
1318	maxBand -= Limits::Small;
1319	}
1320
1321	// insert all queries
1322	for (RayInfoContainer::const_iterator ri = rays.begin(); ri < rays.end(); ++ ri)
1323	{
1324	const bool positive = (*ri).mRay->HasPosDir(axis);
1325	float pos = (*ri).ExtrapOrigin(axis);
1326
1327	// clamp to min / max band
1328	if (0) ClipValue(pos, minBand, maxBand);
1329
1330	mLocalSubdivisionCandidates->push_back(SortableEntry(positive ? SortableEntry::ERayMin : SortableEntry::ERayMax,
1331	pos, (*ri).mRay));
1332
1333	pos = (*ri).ExtrapTermination(axis);
1334
1335	// clamp to min / max band
1336	if (0) ClipValue(pos, minBand, maxBand);
1337
1338	mLocalSubdivisionCandidates->push_back(SortableEntry(positive ? SortableEntry::ERayMax : SortableEntry::ERayMin,
1339	pos, (*ri).mRay));
1340	}
1341
1342	stable_sort(mLocalSubdivisionCandidates->begin(), mLocalSubdivisionCandidates->end());
1343	}
1344
1345
1346	float VspBspTree::BestCostRatioHeuristics(const RayInfoContainer &rays,
1347	const AxisAlignedBox3 &box,
1348	const int pvsSize,
1349	const int axis,
1350	float &position)
1351	{
1352	RayInfoContainer usedRays;
1353
1354	if (mMaxTests < rays.size())
1355	{
1356	GetRayInfoSets(rays, mMaxTests, usedRays);
1357	}
1358	else
1359	{
1360	usedRays = rays;
1361	}
1362
1363	const float minBox = box.Min(axis);
1364	const float maxBox = box.Max(axis);
1365
1366	const float sizeBox = maxBox - minBox;
1367
1368	const float minBand = minBox + mMinBand * sizeBox;
1369	const float maxBand = minBox + mMaxBand * sizeBox;
1370
1371	SortSubdivisionCandidates(usedRays, axis, minBand, maxBand);
1372
1373	// go through the lists, count the number of objects left and right
1374	// and evaluate the following cost funcion:
1375	// C = ct_div_ci + (qlrl + qrrr)/queries
1376
1377	int pvsl = 0;
1378	int pvsr = pvsSize;
1379
1380	int pvsBack = pvsl;
1381	int pvsFront = pvsr;
1382
1383	float sum = (float)pvsSize * sizeBox;
1384	float minSum = 1e20f;
1385
1386
1387	// if no border can be found, take mid split
1388	position = minBox + 0.5f * sizeBox;
1389
1390	// the relative cost ratio
1391	float ratio = 99999999.0f;
1392	bool splitPlaneFound = false;
1393
1394	Intersectable::NewMail();
1395
1396	RayInfoContainer::const_iterator ri, ri_end = usedRays.end();
1397
1398	// set all object as belonging to the front pvs
1399	for(ri = usedRays.begin(); ri != ri_end; ++ ri)
1400	{
1401	Intersectable oObject = (ri).mRay->mOriginObject;
1402	Intersectable tObject = (ri).mRay->mTerminationObject;
1403
1404	if (COUNT_ORIGIN_OBJECTS && oObject)
1405	{
1406	if (!oObject->Mailed())
1407	{
1408	oObject->Mail();
1409	oObject->mCounter = 1;
1410	}
1411	else
1412	{
1413	++ oObject->mCounter;
1414	}
1415	}
1416
1417	if (tObject)
1418	{
1419	if (!tObject->Mailed())
1420	{
1421	tObject->Mail();
1422	tObject->mCounter = 1;
1423	}
1424	else
1425	{
1426	++ tObject->mCounter;
1427	}
1428	}
1429	}
1430
1431	Intersectable::NewMail();
1432
1433	vector<SortableEntry>::const_iterator ci, ci_end = mLocalSubdivisionCandidates->end();
1434
1435	for (ci = mLocalSubdivisionCandidates->begin(); ci != ci_end; ++ ci)
1436	{
1437	VssRay *ray;
1438	ray = (*ci).ray;
1439
1440	Intersectable *oObject = ray->mOriginObject;
1441	Intersectable *tObject = ray->mTerminationObject;
1442
1443
1444	switch ((*ci).type)
1445	{
1446	case SortableEntry::ERayMin:
1447	{
1448	if (COUNT_ORIGIN_OBJECTS && oObject && !oObject->Mailed())
1449	{
1450	oObject->Mail();
1451	++ pvsl;
1452	}
1453
1454	if (tObject && !tObject->Mailed())
1455	{
1456	tObject->Mail();
1457	++ pvsl;
1458	}
1459
1460	break;
1461	}
1462	case SortableEntry::ERayMax:
1463	{
1464	if (COUNT_ORIGIN_OBJECTS && oObject)
1465	{
1466	if (-- oObject->mCounter == 0)
1467	-- pvsr;
1468	}
1469
1470	if (tObject)
1471	{
1472	if (-- tObject->mCounter == 0)
1473	-- pvsr;
1474	}
1475
1476	break;
1477	}
1478	}
1479
1480
1481	// Note: we compare size of bounding boxes of front and back side because
1482	// of efficiency reasons (otherwise a new geometry would have to be computed
1483	// in each step and incremential evaluation would be difficult.
1484	// but then errors happen if the geometry is not an axis aligned box
1485	// (i.e., if a geometry aligned split was taken before)
1486	// question: is it sufficient to make this approximation?
1487	if (((ci).value >= minBand) && ((ci).value <= maxBand))
1488	{
1489	sum = pvsl * ((ci).value - minBox) + pvsr (maxBox - (*ci).value);
1490
1491	float currentPos;
1492
1493	// HACK: current positition is BETWEEN visibility events
1494	if (0 && ((ci + 1) != ci_end))
1495	{
1496	currentPos = ((ci).value + ((ci + 1)).value) * 0.5f;
1497	}
1498	else
1499	currentPos = (*ci).value;
1500
1501	//Debug << "pos=" << (*ci).value << "\t pvs=(" << pvsl << "," << pvsr << ")" << endl;
1502	//Debug << "cost= " << sum << endl;
1503
1504	if (sum < minSum)
1505	{
1506	splitPlaneFound = true;
1507
1508	minSum = sum;
1509	position = currentPos;
1510
1511	pvsBack = pvsl;
1512	pvsFront = pvsr;
1513	}
1514	}
1515	}
1516
1517	// -- compute cost
1518	const int lowerPvsLimit = mViewCellsManager->GetMinPvsSize();
1519	const int upperPvsLimit = mViewCellsManager->GetMaxPvsSize();
1520
1521	const float pOverall = sizeBox;
1522
1523	const float pBack = position - minBox;
1524	const float pFront = maxBox - position;
1525
1526	const float penaltyOld = EvalPvsPenalty(pvsSize, lowerPvsLimit, upperPvsLimit);
1527	const float penaltyFront = EvalPvsPenalty(pvsFront, lowerPvsLimit, upperPvsLimit);
1528	const float penaltyBack = EvalPvsPenalty(pvsBack, lowerPvsLimit, upperPvsLimit);
1529
1530	const float oldRenderCost = penaltyOld * pOverall;
1531	const float newRenderCost = penaltyFront * pFront + penaltyBack * pBack;
1532
1533	if (splitPlaneFound)
1534	{
1535	ratio = mPvsFactor * newRenderCost / (oldRenderCost + Limits::Small);
1536	}
1537	//if (axis != 1)
1538	//Debug << "axis=" << axis << " costRatio=" << ratio << " pos=" << position << " t=" << (position - minBox) / (maxBox - minBox)
1539	// <<"\t pb=(" << pvsBack << ")\t pf=(" << pvsFront << ")" << endl;
1540
1541	return ratio;
1542	}
1543
1544
1545	float VspBspTree::SelectAxisAlignedPlane(Plane3 &plane,
1546	const VspBspTraversalData &tData,
1547	int &axis,
1548	BspNodeGeometry **frontGeom,
1549	BspNodeGeometry **backGeom,
1550	float &pFront,
1551	float &pBack,
1552	const bool isKdNode)
1553	{
1554	float nPosition[3];
1555	float nCostRatio[3];
1556	float nProbFront[3];
1557	float nProbBack[3];
1558
1559	BspNodeGeometry *nFrontGeom[3];
1560	BspNodeGeometry *nBackGeom[3];
1561
1562	// set to NULL, so I can find out which gemetry was stored
1563	for (int i = 0; i < 3; ++ i)
1564	{
1565	nFrontGeom[i] = NULL;
1566	nBackGeom[i] = NULL;
1567	}
1568
1569	// create bounding box of node geometry
1570	AxisAlignedBox3 box;
1571
1572	//TODO: for kd split geometry already is box => only take minmax vertices
1573	if (1)
1574	{ // get bounding box from geometry
1575	tData.mGeometry->GetBoundingBox(box);
1576	}
1577	else
1578	{
1579	box.Initialize();
1580	RayInfoContainer::const_iterator ri, ri_end = tData.mRays->end();
1581
1582	for(ri = tData.mRays->begin(); ri < ri_end; ++ ri)
1583	box.Include((*ri).ExtrapTermination());
1584	}
1585
1586
1587	int sAxis = 0;
1588	int bestAxis;
1589
1590	// if max cost ratio is exceeded, take split along longest axis instead
1591	const float maxCostRatioForArbitraryAxis = 0.9f;
1592
1593	if (mUseDrivingAxisIfMaxCostViolated)
1594	bestAxis = box.Size().DrivingAxis();
1595	else
1596	bestAxis = -1;
1597
1598	#if 0
1599	// maximum cost ratio for axis to be valid:
1600	// if exceeded, spatial mid split is used instead
1601	const maxCostRatioForHeur = 0.99f;
1602	#endif
1603
1604	// if we use some kind of specialised fixed axis
1605	const bool useSpecialAxis =
1606	mOnlyDrivingAxis \|\| mUseRandomAxis \|\| mCirculatingAxis;
1607
1608	if (mUseRandomAxis)
1609	sAxis = Random(3);
1610	else if (mCirculatingAxis)
1611	sAxis = (tData.mAxis + 1) % 3;
1612	else if (mOnlyDrivingAxis)
1613	sAxis = box.Size().DrivingAxis();
1614
1615
1616	//Debug << "use special axis: " << useSpecialAxis << endl;
1617	//Debug << "axis: " << sAxis << " drivingaxis: " << box.Size().DrivingAxis();
1618
1619	for (axis = 0; axis < 3 ; ++ axis)
1620	{
1621	if (!useSpecialAxis \|\| (axis == sAxis))
1622	{
1623	if (mUseCostHeuristics)
1624	{
1625	//-- place split plane using heuristics
1626	nCostRatio[axis] =
1627	BestCostRatioHeuristics(*tData.mRays,
1628	box,
1629	tData.mPvs,
1630	axis,
1631	nPosition[axis]);
1632	}
1633	else
1634	{ //-- split plane position is spatial median
1635
1636	nPosition[axis] = (box.Min()[axis] + box.Max()[axis]) * 0.5f;
1637	Vector3 normal(0,0,0); normal[axis] = 1.0f;
1638
1639	// allows faster split because we have axis aligned kd tree boxes
1640	if (isKdNode)
1641	{
1642	nCostRatio[axis] = EvalAxisAlignedSplitCost(tData,
1643	box,
1644	axis,
1645	nPosition[axis],
1646	nProbFront[axis],
1647	nProbBack[axis]);
1648
1649	Vector3 pos;
1650
1651	// create back geometry from box
1652	// NOTE: the geometry is returned from the function so we
1653	// don't have to recompute it when possible
1654	pos = box.Max(); pos[axis] = nPosition[axis];
1655	AxisAlignedBox3 bBox(box.Min(), pos);
1656	PolygonContainer fPolys;
1657	bBox.ExtractPolys(fPolys);
1658
1659	nBackGeom[axis] = new BspNodeGeometry(fPolys);
1660
1661	//-- create front geometry from box
1662	pos = box.Min(); pos[axis] = nPosition[axis];
1663	AxisAlignedBox3 fBox(pos, box.Max());
1664
1665	PolygonContainer bPolys;
1666	fBox.ExtractPolys(bPolys);
1667	nFrontGeom[axis] = new BspNodeGeometry(bPolys);
1668	}
1669	else
1670	{
1671	nFrontGeom[axis] = new BspNodeGeometry();
1672	nBackGeom[axis] = new BspNodeGeometry();
1673
1674	nCostRatio[axis] =
1675	EvalSplitPlaneCost(Plane3(normal, nPosition[axis]),
1676	tData, nFrontGeom[axis], nBackGeom[axis],
1677	nProbFront[axis], nProbBack[axis]);
1678	}
1679	}
1680
1681
1682	if (mUseDrivingAxisIfMaxCostViolated)
1683	{
1684	// we take longest axis split if cost ratio exceeds threshold
1685	if (nCostRatio[axis] < min(maxCostRatioForArbitraryAxis, nCostRatio[bestAxis]))
1686	{
1687	bestAxis = axis;
1688	}
1689	/*else if (nCostRatio[axis] < nCostRatio[bestAxis])
1690	{
1691	Debug << "taking split along longest axis (" << bestAxis << ") instead of (" << axis << ")" << endl;
1692	}*/
1693
1694	}
1695	else
1696	{
1697	if (bestAxis == -1)
1698	{
1699	bestAxis = axis;
1700	}
1701	else if (nCostRatio[axis] < nCostRatio[bestAxis])
1702	{
1703	bestAxis = axis;
1704	}
1705	}
1706	}
1707	}
1708
1709	//-- assign values
1710
1711	axis = bestAxis;
1712	pFront = nProbFront[bestAxis];
1713	pBack = nProbBack[bestAxis];
1714
1715	// assign best split nodes geometry
1716	*frontGeom = nFrontGeom[bestAxis];
1717	*backGeom = nBackGeom[bestAxis];
1718
1719	// and delete other geometry
1720	DEL_PTR(nFrontGeom[(bestAxis + 1) % 3]);
1721	DEL_PTR(nBackGeom[(bestAxis + 2) % 3]);
1722
1723	//-- split plane
1724	Vector3 normal(0,0,0); normal[bestAxis] = 1;
1725	plane = Plane3(normal, nPosition[bestAxis]);
1726
1727	//Debug << "best axis: " << bestAxis << " pos " << nPosition[bestAxis] << endl;
1728
1729	return nCostRatio[bestAxis];
1730	}
1731
1732
1733	bool VspBspTree::SelectPlane(Plane3 &bestPlane,
1734	BspLeaf *leaf,
1735	VspBspTraversalData &data,
1736	VspBspTraversalData &frontData,
1737	VspBspTraversalData &backData,
1738	int &splitAxis)
1739	{
1740	// HACK matt: subdivide regularily to certain depth
1741	if (data.mDepth < 0) // question matt: why depth < 0 ?
1742	{
1743	cout << "depth: " << data.mDepth << endl;
1744
1745	// return axis aligned split
1746	AxisAlignedBox3 box;
1747	box.Initialize();
1748
1749	// create bounding box of region
1750	data.mGeometry->GetBoundingBox(box);
1751
1752	const int axis = box.Size().DrivingAxis();
1753	const Vector3 position = (box.Min()[axis] + box.Max()[axis]) * 0.5f;
1754
1755	Vector3 norm(0,0,0); norm[axis] = 1.0f;
1756	bestPlane = Plane3(norm, position);
1757	splitAxis = axis;
1758
1759	return true;
1760	}
1761
1762	// simplest strategy: just take next polygon
1763	if (mSplitPlaneStrategy & RANDOM_POLYGON)
1764	{
1765	if (!data.mPolygons->empty())
1766	{
1767	const int randIdx =
1768	(int)RandomValue(0, (Real)((int)data.mPolygons->size() - 1));
1769	Polygon3 nextPoly = (data.mPolygons)[randIdx];
1770
1771	bestPlane = nextPoly->GetSupportingPlane();
1772	return true;
1773	}
1774	}
1775
1776	//-- use heuristics to find appropriate plane
1777
1778	// intermediate plane
1779	Plane3 plane;
1780	float lowestCost = MAX_FLOAT;
1781
1782	// decides if the first few splits should be only axisAligned
1783	const bool onlyAxisAligned =
1784	(mSplitPlaneStrategy & AXIS_ALIGNED) &&
1785	((int)data.mRays->size() > mTermMinRaysForAxisAligned) &&
1786	((int)data.GetAvgRayContribution() < mTermMaxRayContriForAxisAligned);
1787
1788	const int limit = onlyAxisAligned ? 0 :
1789	Min((int)data.mPolygons->size(), mMaxPolyCandidates);
1790
1791	float candidateCost;
1792
1793	int maxIdx = (int)data.mPolygons->size();
1794
1795	for (int i = 0; i < limit; ++ i)
1796	{
1797	// the already taken candidates are stored behind maxIdx
1798	// => assure that no index is taken twice
1799	const int candidateIdx = (int)RandomValue(0, (Real)(-- maxIdx));
1800	Polygon3 poly = (data.mPolygons)[candidateIdx];
1801
1802	// swap candidate to the end to avoid testing same plane
1803	std::swap((data.mPolygons)[maxIdx], (data.mPolygons)[candidateIdx]);
1804	//Polygon3 poly = (data.mPolygons)[(int)RandomValue(0, (int)polys.size() - 1)];
1805
1806	// evaluate current candidate
1807	BspNodeGeometry fGeom, bGeom;
1808	float fArea, bArea;
1809	plane = poly->GetSupportingPlane();
1810	candidateCost = EvalSplitPlaneCost(plane, data, fGeom, bGeom, fArea, bArea);
1811
1812	if (candidateCost < lowestCost)
1813	{
1814	bestPlane = plane;
1815	lowestCost = candidateCost;
1816	}
1817	}
1818
1819
1820	//-- evaluate axis aligned splits
1821
1822	int axis;
1823	BspNodeGeometry fGeom, bGeom;
1824	float pFront, pBack;
1825
1826	candidateCost = 99999999.0f;
1827
1828	// as a variant, we take axis aligned split only if there is
1829	// more polygon available to guide the split
1830	if (!mUsePolygonSplitIfAvailable \|\| data.mPolygons->empty())
1831	{
1832	candidateCost = SelectAxisAlignedPlane(plane,
1833	data,
1834	axis,
1835	&fGeom,
1836	&bGeom,
1837	pFront,
1838	pBack,
1839	data.mIsKdNode);
1840	}
1841
1842	splitAxis = 3;
1843
1844	if (candidateCost < lowestCost)
1845	{
1846	bestPlane = plane;
1847	lowestCost = candidateCost;
1848	splitAxis = axis;
1849
1850	// assign already computed values
1851	// we can do this because we always save the
1852	// computed values from the axis aligned splits
1853
1854	if (fGeom && bGeom)
1855	{
1856	frontData.mGeometry = fGeom;
1857	backData.mGeometry = bGeom;
1858
1859	frontData.mProbability = pFront;
1860	backData.mProbability = pBack;
1861	}
1862	}
1863	else
1864	{
1865	DEL_PTR(fGeom);
1866	DEL_PTR(bGeom);
1867	}
1868
1869	#ifdef _DEBUG
1870	Debug << "plane lowest cost: " << lowestCost << endl;
1871	#endif
1872
1873	// exeeded relative max cost ratio
1874	if (lowestCost > mTermMaxCostRatio)
1875	{
1876	return false;
1877	}
1878
1879	return true;
1880	}
1881
1882
1883	Plane3 VspBspTree::ChooseCandidatePlane(const RayInfoContainer &rays) const
1884	{
1885	const int candidateIdx = (int)RandomValue(0, (Real)((int)rays.size() - 1));
1886
1887	const Vector3 minPt = rays[candidateIdx].ExtrapOrigin();
1888	const Vector3 maxPt = rays[candidateIdx].ExtrapTermination();
1889
1890	const Vector3 pt = (maxPt + minPt) * 0.5;
1891	const Vector3 normal = Normalize(rays[candidateIdx].mRay->GetDir());
1892
1893	return Plane3(normal, pt);
1894	}
1895
1896
1897	Plane3 VspBspTree::ChooseCandidatePlane2(const RayInfoContainer &rays) const
1898	{
1899	Vector3 pt[3];
1900
1901	int idx[3];
1902	int cmaxT = 0;
1903	int cminT = 0;
1904	bool chooseMin = false;
1905
1906	for (int j = 0; j < 3; ++ j)
1907	{
1908	idx[j] = (int)RandomValue(0, (Real)((int)rays.size() * 2 - 1));
1909
1910	if (idx[j] >= (int)rays.size())
1911	{
1912	idx[j] -= (int)rays.size();
1913
1914	chooseMin = (cminT < 2);
1915	}
1916	else
1917	chooseMin = (cmaxT < 2);
1918
1919	RayInfo rayInf = rays[idx[j]];
1920	pt[j] = chooseMin ? rayInf.ExtrapOrigin() : rayInf.ExtrapTermination();
1921	}
1922
1923	return Plane3(pt[0], pt[1], pt[2]);
1924	}
1925
1926
1927	Plane3 VspBspTree::ChooseCandidatePlane3(const RayInfoContainer &rays) const
1928	{
1929	Vector3 pt[3];
1930
1931	int idx1 = (int)RandomValue(0, (Real)((int)rays.size() - 1));
1932	int idx2 = (int)RandomValue(0, (Real)((int)rays.size() - 1));
1933
1934	// check if rays different
1935	if (idx1 == idx2)
1936	idx2 = (idx2 + 1) % (int)rays.size();
1937
1938	const RayInfo ray1 = rays[idx1];
1939	const RayInfo ray2 = rays[idx2];
1940
1941	// normal vector of the plane parallel to both lines
1942	const Vector3 norm = Normalize(CrossProd(ray1.mRay->GetDir(), ray2.mRay->GetDir()));
1943
1944	// vector from line 1 to line 2
1945	const Vector3 vd = ray2.ExtrapOrigin() - ray1.ExtrapOrigin();
1946
1947	// project vector on normal to get distance
1948	const float dist = DotProd(vd, norm);
1949
1950	// point on plane lies halfway between the two planes
1951	const Vector3 planePt = ray1.ExtrapOrigin() + norm * dist * 0.5;
1952
1953	return Plane3(norm, planePt);
1954	}
1955
1956
1957	inline void VspBspTree::GenerateUniqueIdsForPvs()
1958	{
1959	Intersectable::NewMail(); sBackId = Intersectable::sMailId;
1960	Intersectable::NewMail(); sFrontId = Intersectable::sMailId;
1961	Intersectable::NewMail(); sFrontAndBackId = Intersectable::sMailId;
1962	}
1963
1964
1965	float VspBspTree::EvalRenderCostDecrease(const Plane3 &candidatePlane,
1966	const VspBspTraversalData &data,
1967	float &normalizedOldRenderCost) const
1968	{
1969	float pvsFront = 0;
1970	float pvsBack = 0;
1971	float totalPvs = 0;
1972
1973	// probability that view point lies in back / front node
1974	float pOverall = data.mProbability;
1975	float pFront = 0;
1976	float pBack = 0;
1977
1978
1979	// create unique ids for pvs heuristics
1980	GenerateUniqueIdsForPvs();
1981
1982	for (int i = 0; i < data.mRays->size(); ++ i)
1983	{
1984	RayInfo rayInf = (*data.mRays)[i];
1985
1986	float t;
1987	VssRay *ray = rayInf.mRay;
1988	const int cf = rayInf.ComputeRayIntersection(candidatePlane, t);
1989
1990	// find front and back pvs for origing and termination object
1991	AddObjToPvs(ray->mTerminationObject, cf, pvsFront, pvsBack, totalPvs);
1992
1993	if (COUNT_ORIGIN_OBJECTS)
1994	AddObjToPvs(ray->mOriginObject, cf, pvsFront, pvsBack, totalPvs);
1995	}
1996
1997
1998	BspNodeGeometry geomFront;
1999	BspNodeGeometry geomBack;
2000
2001	// construct child geometry with regard to the candidate split plane
2002	data.mGeometry->SplitGeometry(geomFront,
2003	geomBack,
2004	candidatePlane,
2005	mBoundingBox,
2006	//0.0f);
2007	mEpsilon);
2008
2009	if (!mUseAreaForPvs) // use front and back cell areas to approximate volume
2010	{
2011	pFront = geomFront.GetVolume();
2012	pBack = pOverall - pFront;
2013
2014	// something is wrong with the volume
2015	if (0 && ((pFront < 0.0) \|\| (pBack < 0.0)))
2016	{
2017	Debug << "ERROR in volume:\n"
2018	<< "volume f :" << pFront << " b: " << pBack << " p: " << pOverall
2019	<< ", real volume f: " << pFront << " b: " << geomBack.GetVolume()
2020	<< ", #polygons f: " << geomFront.Size() << " b: " << geomBack.Size() << " p: " << data.mGeometry->Size() << endl;
2021	}
2022	}
2023	else
2024	{
2025	pFront = geomFront.GetArea();
2026	pBack = geomBack.GetArea();
2027	}
2028
2029
2030	// -- pvs rendering heuristics
2031
2032	// upper and lower bounds
2033	const int lowerPvsLimit = mViewCellsManager->GetMinPvsSize();
2034	const int upperPvsLimit = mViewCellsManager->GetMaxPvsSize();
2035
2036	const float penaltyOld = EvalPvsPenalty((int)totalPvs, lowerPvsLimit, upperPvsLimit);
2037	const float penaltyFront = EvalPvsPenalty((int)pvsFront, lowerPvsLimit, upperPvsLimit);
2038	const float penaltyBack = EvalPvsPenalty((int)pvsBack, lowerPvsLimit, upperPvsLimit);
2039
2040	const float oldRenderCost = pOverall * penaltyOld;
2041	const float newRenderCost = penaltyFront * pFront + penaltyBack * pBack;
2042
2043	const float renderCostDecrease = (oldRenderCost - newRenderCost) / mBoundingBox.GetVolume();
2044
2045	// take render cost of node into account to avoid being stuck in a local minimum
2046	normalizedOldRenderCost = oldRenderCost / mBoundingBox.GetVolume();
2047
2048	return renderCostDecrease;
2049	}
2050
2051
2052	float VspBspTree::EvalSplitPlaneCost(const Plane3 &candidatePlane,
2053	const VspBspTraversalData &data,
2054	BspNodeGeometry &geomFront,
2055	BspNodeGeometry &geomBack,
2056	float &pFront,
2057	float &pBack) const
2058	{
2059	float totalPvs = 0;
2060	float pvsFront = 0;
2061	float pvsBack = 0;
2062
2063	// overall probability is used as normalizer
2064	float pOverall = 0;
2065
2066	// probability that view point lies in back / front node
2067	pFront = 0;
2068	pBack = 0;
2069
2070	int numTests; // the number of tests
2071
2072	// if random samples shold be taken instead of testing all the rays
2073	bool useRand;
2074
2075	if ((int)data.mRays->size() > mMaxTests)
2076	{
2077	useRand = true;
2078	numTests = mMaxTests;
2079	}
2080	else
2081	{
2082	useRand = false;
2083	numTests = (int)data.mRays->size();
2084	}
2085
2086	// create unique ids for pvs heuristics
2087	GenerateUniqueIdsForPvs();
2088
2089	for (int i = 0; i < numTests; ++ i)
2090	{
2091	const int testIdx = useRand ?
2092	(int)RandomValue(0, (Real)((int)data.mRays->size() - 1)) : i;
2093	RayInfo rayInf = (*data.mRays)[testIdx];
2094
2095	float t;
2096	VssRay *ray = rayInf.mRay;
2097	const int cf = rayInf.ComputeRayIntersection(candidatePlane, t);
2098
2099	// find front and back pvs for origing and termination object
2100	AddObjToPvs(ray->mTerminationObject, cf, pvsFront, pvsBack, totalPvs);
2101
2102	if (COUNT_ORIGIN_OBJECTS)
2103	AddObjToPvs(ray->mOriginObject, cf, pvsFront, pvsBack, totalPvs);
2104	}
2105
2106	// construct child geometry with regard to the candidate split plane
2107	bool splitSuccessFull = data.mGeometry->SplitGeometry(geomFront,
2108	geomBack,
2109	candidatePlane,
2110	mBoundingBox,
2111	//0.0f);
2112	mEpsilon);
2113
2114	pOverall = data.mProbability;
2115
2116	if (!mUseAreaForPvs) // use front and back cell areas to approximate volume
2117	{
2118	pFront = geomFront.GetVolume();
2119	pBack = pOverall - pFront;
2120
2121	// HACK: precision issues possible for unbalanced split => don't take this split!
2122	if (1 &&
2123	(!splitSuccessFull \|\| (pFront <= 0) \|\| (pBack <= 0) \|\|
2124	!geomFront.Valid() \|\| !geomBack.Valid()))
2125	{
2126	//Debug << "error f: " << pFront << " b: " << pBack << endl;
2127
2128	// high penalty for degenerated / wrong split
2129	return 99999.9f;
2130	}
2131	}
2132	else
2133	{
2134	pFront = geomFront.GetArea();
2135	pBack = geomBack.GetArea();
2136	}
2137
2138
2139	// -- pvs rendering heuristics
2140	const int lowerPvsLimit = mViewCellsManager->GetMinPvsSize();
2141	const int upperPvsLimit = mViewCellsManager->GetMaxPvsSize();
2142
2143	// only render cost heuristics or combined with standard deviation
2144	const float penaltyOld = EvalPvsPenalty((int)totalPvs, lowerPvsLimit, upperPvsLimit);
2145	const float penaltyFront = EvalPvsPenalty((int)pvsFront, lowerPvsLimit, upperPvsLimit);
2146	const float penaltyBack = EvalPvsPenalty((int)pvsBack, lowerPvsLimit, upperPvsLimit);
2147
2148	const float oldRenderCost = pOverall * penaltyOld;
2149	const float newRenderCost = penaltyFront * pFront + penaltyBack * pBack;
2150
2151	float oldCost, newCost;
2152
2153	// only render cost
2154	if (1)
2155	{
2156	oldCost = oldRenderCost;
2157	newCost = newRenderCost;
2158	}
2159	else // also considering standard deviation
2160	{
2161	// standard deviation is difference of back and front pvs
2162	const float expectedCost = 0.5f * (penaltyFront + penaltyBack);
2163
2164	const float newDeviation = 0.5f *
2165	fabs(penaltyFront - expectedCost) + fabs(penaltyBack - expectedCost);
2166
2167	const float oldDeviation = penaltyOld;
2168
2169	newCost = mRenderCostWeight * newRenderCost + (1.0f - mRenderCostWeight) * newDeviation;
2170	oldCost = mRenderCostWeight * oldRenderCost + (1.0f - mRenderCostWeight) * oldDeviation;
2171	}
2172
2173	const float cost = mPvsFactor * newCost / (oldCost + Limits::Small);
2174
2175
2176	#ifdef _DEBUG
2177	Debug << "totalpvs: " << data.mPvs << " ptotal: " << pOverall
2178	<< " frontpvs: " << pvsFront << " pFront: " << pFront
2179	<< " backpvs: " << pvsBack << " pBack: " << pBack << endl << endl;
2180	Debug << "cost: " << cost << endl;
2181	#endif
2182
2183	return cost;
2184	}
2185
2186
2187	int VspBspTree::ComputeBoxIntersections(const AxisAlignedBox3 &box,
2188	ViewCellContainer &viewCells) const
2189	{
2190	stack<bspNodePair> nodeStack;
2191	BspNodeGeometry *rgeom = new BspNodeGeometry();
2192
2193	ConstructGeometry(mRoot, *rgeom);
2194
2195	nodeStack.push(bspNodePair(mRoot, rgeom));
2196
2197	ViewCell::NewMail();
2198
2199	while (!nodeStack.empty())
2200	{
2201	BspNode *node = nodeStack.top().first;
2202	BspNodeGeometry *geom = nodeStack.top().second;
2203	nodeStack.pop();
2204
2205	const int side = geom->ComputeIntersection(box);
2206
2207	switch (side)
2208	{
2209	case -1:
2210	// node geometry is contained in box
2211	CollectViewCells(node, true, viewCells, true);
2212	break;
2213
2214	case 0:
2215	if (node->IsLeaf())
2216	{
2217	BspLeaf leaf = dynamic_cast<BspLeaf >(node);
2218
2219	if (!leaf->GetViewCell()->Mailed() && leaf->TreeValid())
2220	{
2221	leaf->GetViewCell()->Mail();
2222	viewCells.push_back(leaf->GetViewCell());
2223	}
2224	}
2225	else
2226	{
2227	BspInterior interior = dynamic_cast<BspInterior >(node);
2228
2229	BspNode *first = interior->GetFront();
2230	BspNode *second = interior->GetBack();
2231
2232	BspNodeGeometry *firstGeom = new BspNodeGeometry();
2233	BspNodeGeometry *secondGeom = new BspNodeGeometry();
2234
2235	geom->SplitGeometry(*firstGeom,
2236	*secondGeom,
2237	interior->GetPlane(),
2238	mBoundingBox,
2239	//0.0000001f);
2240	mEpsilon);
2241
2242	nodeStack.push(bspNodePair(first, firstGeom));
2243	nodeStack.push(bspNodePair(second, secondGeom));
2244	}
2245
2246	break;
2247	default:
2248	// default: cull
2249	break;
2250	}
2251
2252	DEL_PTR(geom);
2253
2254	}
2255
2256	return (int)viewCells.size();
2257	}
2258
2259
2260	float VspBspTree::EvalAxisAlignedSplitCost(const VspBspTraversalData &data,
2261	const AxisAlignedBox3 &box,
2262	const int axis,
2263	const float &position,
2264	float &pFront,
2265	float &pBack) const
2266	{
2267	float pvsTotal = 0;
2268	float pvsFront = 0;
2269	float pvsBack = 0;
2270
2271	// create unique ids for pvs heuristics
2272	GenerateUniqueIdsForPvs();
2273
2274	const int pvsSize = data.mPvs;
2275
2276	RayInfoContainer::const_iterator rit, rit_end = data.mRays->end();
2277
2278	// this is the main ray classification loop!
2279	for(rit = data.mRays->begin(); rit != rit_end; ++ rit)
2280	{
2281	// determine the side of this ray with respect to the plane
2282	float t;
2283	const int side = (*rit).ComputeRayIntersection(axis, position, t);
2284
2285	AddObjToPvs((*rit).mRay->mTerminationObject, side, pvsFront, pvsBack, pvsTotal);
2286
2287	if (COUNT_ORIGIN_OBJECTS)
2288	AddObjToPvs((*rit).mRay->mOriginObject, side, pvsFront, pvsBack, pvsTotal);
2289	}
2290
2291
2292	//-- pvs heuristics
2293
2294	float pOverall = data.mProbability;
2295
2296	// note: we use a simplified computation assuming that we always do a
2297	// spatial mid split
2298
2299	if (!mUseAreaForPvs)
2300	{
2301	// volume
2302	pBack = pFront = pOverall * 0.5f;
2303	#if 0
2304	// box length substitute for probability
2305	const float minBox = box.Min(axis);
2306	const float maxBox = box.Max(axis);
2307
2308	pBack = position - minBox;
2309	pFront = maxBox - position;
2310	pOverall = maxBox - minBox;
2311	#endif
2312	}
2313	else //-- area substitute for probability
2314	{
2315	const int axis2 = (axis + 1) % 3;
2316	const int axis3 = (axis + 2) % 3;
2317
2318	const float faceArea =
2319	(box.Max(axis2) - box.Min(axis2)) *
2320	(box.Max(axis3) - box.Min(axis3));
2321
2322	pBack = pFront = pOverall * 0.5f + faceArea;
2323	}
2324
2325	#ifdef _DEBUG
2326	Debug << "axis: " << axis << " " << pvsSize << " " << pvsBack << " " << pvsFront << endl;
2327	Debug << "p: " << pFront << " " << pBack << " " << pOverall << endl;
2328	#endif
2329
2330
2331	const float newCost = pvsBack * pBack + pvsFront * pFront;
2332	const float oldCost = (float)pvsSize * pOverall + Limits::Small;
2333
2334	return (mCtDivCi + newCost) / oldCost;
2335	}
2336
2337
2338	inline void VspBspTree::AddObjToPvs(Intersectable *obj,
2339	const int cf,
2340	float &frontPvs,
2341	float &backPvs,
2342	float &totalPvs) const
2343	{
2344	if (!obj)
2345	return;
2346	#if 0
2347	const float renderCost = mViewCellsManager->EvalRenderCost(obj);
2348	#else
2349	const int renderCost = 1;
2350	#endif
2351	// new object
2352	if ((obj->mMailbox != sFrontId) &&
2353	(obj->mMailbox != sBackId) &&
2354	(obj->mMailbox != sFrontAndBackId))
2355	{
2356	totalPvs += renderCost;
2357	}
2358
2359	// TODO: does this really belong to no pvs?
2360	//if (cf == Ray::COINCIDENT) return;
2361
2362	// object belongs to both PVS
2363	if (cf >= 0)
2364	{
2365	if ((obj->mMailbox != sFrontId) &&
2366	(obj->mMailbox != sFrontAndBackId))
2367	{
2368	frontPvs += renderCost;
2369
2370	if (obj->mMailbox == sBackId)
2371	obj->mMailbox = sFrontAndBackId;
2372	else
2373	obj->mMailbox = sFrontId;
2374	}
2375	}
2376
2377	if (cf <= 0)
2378	{
2379	if ((obj->mMailbox != sBackId) &&
2380	(obj->mMailbox != sFrontAndBackId))
2381	{
2382	backPvs += renderCost;
2383
2384	if (obj->mMailbox == sFrontId)
2385	obj->mMailbox = sFrontAndBackId;
2386	else
2387	obj->mMailbox = sBackId;
2388	}
2389	}
2390	}
2391
2392
2393	void VspBspTree::CollectLeaves(vector<BspLeaf *> &leaves,
2394	const bool onlyUnmailed,
2395	const int maxPvsSize) const
2396	{
2397	stack<BspNode *> nodeStack;
2398	nodeStack.push(mRoot);
2399
2400	while (!nodeStack.empty())
2401	{
2402	BspNode *node = nodeStack.top();
2403	nodeStack.pop();
2404
2405	if (node->IsLeaf())
2406	{
2407	// test if this leaf is in valid view space
2408	BspLeaf leaf = dynamic_cast<BspLeaf >(node);
2409	if (leaf->TreeValid() &&
2410	(!onlyUnmailed \|\| !leaf->Mailed()) &&
2411	((maxPvsSize < 0) \|\| (leaf->GetViewCell()->GetPvs().CountObjectsInPvs() <= maxPvsSize)))
2412	{
2413	leaves.push_back(leaf);
2414	}
2415	}
2416	else
2417	{
2418	BspInterior interior = dynamic_cast<BspInterior >(node);
2419
2420	nodeStack.push(interior->GetBack());
2421	nodeStack.push(interior->GetFront());
2422	}
2423	}
2424	}
2425
2426
2427	AxisAlignedBox3 VspBspTree::GetBoundingBox() const
2428	{
2429	return mBoundingBox;
2430	}
2431
2432
2433	BspNode *VspBspTree::GetRoot() const
2434	{
2435	return mRoot;
2436	}
2437
2438
2439	void VspBspTree::EvaluateLeafStats(const VspBspTraversalData &data)
2440	{
2441	// the node became a leaf -> evaluate stats for leafs
2442	BspLeaf leaf = dynamic_cast<BspLeaf >(data.mNode);
2443
2444
2445	if (data.mPvs > mBspStats.maxPvs)
2446	{
2447	mBspStats.maxPvs = data.mPvs;
2448	}
2449
2450	mBspStats.pvs += data.mPvs;
2451
2452	if (data.mDepth < mBspStats.minDepth)
2453	{
2454	mBspStats.minDepth = data.mDepth;
2455	}
2456
2457	if (data.mDepth >= mTermMaxDepth)
2458	{
2459	++ mBspStats.maxDepthNodes;
2460	//Debug << "new max depth: " << mBspStats.maxDepthNodes << endl;
2461	}
2462
2463	// accumulate rays to compute rays / leaf
2464	mBspStats.accumRays += (int)data.mRays->size();
2465
2466	if (data.mPvs < mTermMinPvs)
2467	++ mBspStats.minPvsNodes;
2468
2469	if ((int)data.mRays->size() < mTermMinRays)
2470	++ mBspStats.minRaysNodes;
2471
2472	if (data.GetAvgRayContribution() > mTermMaxRayContribution)
2473	++ mBspStats.maxRayContribNodes;
2474
2475	if (data.mProbability <= mTermMinProbability)
2476	++ mBspStats.minProbabilityNodes;
2477
2478	// accumulate depth to compute average depth
2479	mBspStats.accumDepth += data.mDepth;
2480
2481	++ mCreatedViewCells;
2482
2483	#ifdef _DEBUG
2484	Debug << "BSP stats: "
2485	<< "Depth: " << data.mDepth << " (max: " << mTermMaxDepth << "), "
2486	<< "PVS: " << data.mPvs << " (min: " << mTermMinPvs << "), "
2487	<< "#rays: " << (int)data.mRays->size() << " (max: " << mTermMinRays << "), "
2488	<< "#pvs: " << leaf->GetViewCell()->GetPvs().CountObjectsInPvs() << "), "
2489	<< "#avg ray contrib (pvs): " << (float)data.mPvs / (float)data.mRays->size() << endl;
2490	#endif
2491	}
2492
2493
2494	int VspBspTree::CastRay(Ray &ray)
2495	{
2496	int hits = 0;
2497
2498	stack<BspRayTraversalData> tQueue;
2499
2500	float maxt, mint;
2501
2502	if (!mBoundingBox.GetRaySegment(ray, mint, maxt))
2503	return 0;
2504
2505	Intersectable::NewMail();
2506	ViewCell::NewMail();
2507
2508	Vector3 entp = ray.Extrap(mint);
2509	Vector3 extp = ray.Extrap(maxt);
2510
2511	BspNode *node = mRoot;
2512	BspNode *farChild = NULL;
2513
2514	while (1)
2515	{
2516	if (!node->IsLeaf())
2517	{
2518	BspInterior in = dynamic_cast<BspInterior >(node);
2519
2520	Plane3 splitPlane = in->GetPlane();
2521	const int entSide = splitPlane.Side(entp);
2522	const int extSide = splitPlane.Side(extp);
2523
2524	if (entSide < 0)
2525	{
2526	node = in->GetBack();
2527
2528	if(extSide <= 0) // plane does not split ray => no far child
2529	continue;
2530
2531	farChild = in->GetFront(); // plane splits ray
2532
2533	} else if (entSide > 0)
2534	{
2535	node = in->GetFront();
2536
2537	if (extSide >= 0) // plane does not split ray => no far child
2538	continue;
2539
2540	farChild = in->GetBack(); // plane splits ray
2541	}
2542	else // ray and plane are coincident
2543	{
2544	// matt: WHAT TO DO IN THIS CASE ?
2545	//break;
2546	node = in->GetFront();
2547	continue;
2548	}
2549
2550	// push data for far child
2551	tQueue.push(BspRayTraversalData(farChild, extp, maxt));
2552
2553	// find intersection of ray segment with plane
2554	float t;
2555	extp = splitPlane.FindIntersection(ray.GetLoc(), extp, &t);
2556	maxt *= t;
2557	}
2558	else // reached leaf => intersection with view cell
2559	{
2560	BspLeaf leaf = dynamic_cast<BspLeaf >(node);
2561
2562	if (!leaf->GetViewCell()->Mailed())
2563	{
2564	//ray.bspIntersections.push_back(Ray::VspBspIntersection(maxt, leaf));
2565	leaf->GetViewCell()->Mail();
2566	++ hits;
2567	}
2568
2569	//-- fetch the next far child from the stack
2570	if (tQueue.empty())
2571	break;
2572
2573	entp = extp;
2574	mint = maxt; // NOTE: need this?
2575
2576	if (ray.GetType() == Ray::LINE_SEGMENT && mint > 1.0f)
2577	break;
2578
2579	BspRayTraversalData &s = tQueue.top();
2580
2581	node = s.mNode;
2582	extp = s.mExitPoint;
2583	maxt = s.mMaxT;
2584
2585	tQueue.pop();
2586	}
2587	}
2588
2589	return hits;
2590	}
2591
2592
2593	void VspBspTree::CollectViewCells(ViewCellContainer &viewCells,
2594	bool onlyValid) const
2595	{
2596	ViewCell::NewMail();
2597	CollectViewCells(mRoot, onlyValid, viewCells, true);
2598	}
2599
2600
2601	void VspBspTree::CollectViewCells(BspNode *root,
2602	bool onlyValid,
2603	ViewCellContainer &viewCells,
2604	bool onlyUnmailed) const
2605	{
2606	stack<BspNode *> nodeStack;
2607
2608	if (!root)
2609	return;
2610
2611	nodeStack.push(root);
2612
2613	while (!nodeStack.empty())
2614	{
2615	BspNode *node = nodeStack.top();
2616	nodeStack.pop();
2617
2618	if (node->IsLeaf())
2619	{
2620	if (!onlyValid \|\| node->TreeValid())
2621	{
2622	ViewCellLeaf leafVc = dynamic_cast<BspLeaf >(node)->GetViewCell();
2623
2624	ViewCell *viewCell = mViewCellsTree->GetActiveViewCell(leafVc);
2625
2626	if (!onlyUnmailed \|\| !viewCell->Mailed())
2627	{
2628	viewCell->Mail();
2629	viewCells.push_back(viewCell);
2630	}
2631	}
2632	}
2633	else
2634	{
2635	BspInterior interior = dynamic_cast<BspInterior >(node);
2636
2637	nodeStack.push(interior->GetFront());
2638	nodeStack.push(interior->GetBack());
2639	}
2640	}
2641	}
2642
2643
2644	void VspBspTree::CollapseViewCells()
2645	{
2646	// TODO
2647	#if HAS_TO_BE_REDONE
2648	stack<BspNode *> nodeStack;
2649
2650	if (!mRoot)
2651	return;
2652
2653	nodeStack.push(mRoot);
2654
2655	while (!nodeStack.empty())
2656	{
2657	BspNode *node = nodeStack.top();
2658	nodeStack.pop();
2659
2660	if (node->IsLeaf())
2661	{
2662	BspViewCell viewCell = dynamic_cast<BspLeaf >(node)->GetViewCell();
2663
2664	if (!viewCell->GetValid())
2665	{
2666	BspViewCell viewCell = dynamic_cast<BspLeaf >(node)->GetViewCell();
2667
2668	ViewCellContainer leaves;
2669	mViewCellsTree->CollectLeaves(viewCell, leaves);
2670
2671	ViewCellContainer::const_iterator it, it_end = leaves.end();
2672
2673	for (it = leaves.begin(); it != it_end; ++ it)
2674	{
2675	BspLeaf l = dynamic_cast<BspViewCell >(*it)->mLeaf;
2676	l->SetViewCell(GetOrCreateOutOfBoundsCell());
2677	++ mBspStats.invalidLeaves;
2678	}
2679
2680	// add to unbounded view cell
2681	GetOrCreateOutOfBoundsCell()->GetPvs().AddPvs(viewCell->GetPvs());
2682	DEL_PTR(viewCell);
2683	}
2684	}
2685	else
2686	{
2687	BspInterior interior = dynamic_cast<BspInterior >(node);
2688
2689	nodeStack.push(interior->GetFront());
2690	nodeStack.push(interior->GetBack());
2691	}
2692	}
2693
2694	Debug << "invalid leaves: " << mBspStats.invalidLeaves << endl;
2695	#endif
2696	}
2697
2698
2699	void VspBspTree::CollectRays(VssRayContainer &rays)
2700	{
2701	vector<BspLeaf *> leaves;
2702
2703	vector<BspLeaf *>::const_iterator lit, lit_end = leaves.end();
2704
2705	for (lit = leaves.begin(); lit != lit_end; ++ lit)
2706	{
2707	BspLeaf leaf = lit;
2708	VssRayContainer::const_iterator rit, rit_end = leaf->mVssRays.end();
2709
2710	for (rit = leaf->mVssRays.begin(); rit != rit_end; ++ rit)
2711	rays.push_back(*rit);
2712	}
2713	}
2714
2715
2716	void VspBspTree::ValidateTree()
2717	{
2718	stack<BspNode *> nodeStack;
2719
2720	if (!mRoot)
2721	return;
2722
2723	nodeStack.push(mRoot);
2724
2725	mBspStats.invalidLeaves = 0;
2726	while (!nodeStack.empty())
2727	{
2728	BspNode *node = nodeStack.top();
2729	nodeStack.pop();
2730
2731	if (node->IsLeaf())
2732	{
2733	BspLeaf leaf = dynamic_cast<BspLeaf >(node);
2734
2735	if (!leaf->GetViewCell()->GetValid())
2736	++ mBspStats.invalidLeaves;
2737
2738	// validity flags don't match => repair
2739	if (leaf->GetViewCell()->GetValid() != leaf->TreeValid())
2740	{
2741	leaf->SetTreeValid(leaf->GetViewCell()->GetValid());
2742	PropagateUpValidity(leaf);
2743	}
2744	}
2745	else
2746	{
2747	BspInterior interior = dynamic_cast<BspInterior >(node);
2748
2749	nodeStack.push(interior->GetFront());
2750	nodeStack.push(interior->GetBack());
2751	}
2752	}
2753
2754	Debug << "invalid leaves: " << mBspStats.invalidLeaves << endl;
2755	}
2756
2757
2758	void VspBspTree::PreprocessPolygons(PolygonContainer &polys)
2759	{
2760	// preprocess: throw out polygons coincident to the view space box (not needed)
2761	PolygonContainer boxPolys;
2762
2763	mBoundingBox.ExtractPolys(boxPolys);
2764	vector<Plane3> boxPlanes;
2765
2766	PolygonContainer::iterator pit, pit_end = boxPolys.end();
2767
2768	// extract planes of box
2769	// TODO: can be done more elegantly than first extracting polygons
2770	// and take their planes
2771	for (pit = boxPolys.begin(); pit != pit_end; ++ pit)
2772	{
2773	boxPlanes.push_back((*pit)->GetSupportingPlane());
2774	}
2775
2776	pit_end = polys.end();
2777
2778	for (pit = polys.begin(); pit != pit_end; ++ pit)
2779	{
2780	vector<Plane3>::const_iterator bit, bit_end = boxPlanes.end();
2781
2782	for (bit = boxPlanes.begin(); (bit != bit_end) && (*pit); ++ bit)
2783	{
2784	const int cf = (pit)->ClassifyPlane(bit, mEpsilon);
2785
2786	if (cf == Polygon3::COINCIDENT)
2787	{
2788	DEL_PTR(*pit);
2789	//Debug << "coincident!!" << endl;
2790	}
2791	}
2792	}
2793
2794	// remove deleted entries after swapping them to end of vector
2795	for (int i = 0; i < (int)polys.size(); ++ i)
2796	{
2797	while (!polys[i] && (i < (int)polys.size()))
2798	{
2799	swap(polys[i], polys.back());
2800	polys.pop_back();
2801	}
2802	}
2803	}
2804
2805
2806	float VspBspTree::AccumulatedRayLength(const RayInfoContainer &rays) const
2807	{
2808	float len = 0;
2809
2810	RayInfoContainer::const_iterator it, it_end = rays.end();
2811
2812	for (it = rays.begin(); it != it_end; ++ it)
2813	len += (*it).SegmentLength();
2814
2815	return len;
2816	}
2817
2818
2819	int VspBspTree::SplitRays(const Plane3 &plane,
2820	RayInfoContainer &rays,
2821	RayInfoContainer &frontRays,
2822	RayInfoContainer &backRays) const
2823	{
2824	int splits = 0;
2825
2826	RayInfoContainer::const_iterator it, it_end = rays.end();
2827
2828	for (it = rays.begin(); it != it_end; ++ it)
2829	{
2830	RayInfo bRay = *it;
2831
2832	VssRay *ray = bRay.mRay;
2833	float t;
2834
2835	// get classification and receive new t
2836	const int cf = bRay.ComputeRayIntersection(plane, t);
2837
2838	switch (cf)
2839	{
2840	case -1:
2841	backRays.push_back(bRay);
2842	break;
2843	case 1:
2844	frontRays.push_back(bRay);
2845	break;
2846	case 0:
2847	{
2848	//-- split ray
2849	// test if start point behind or in front of plane
2850	const int side = plane.Side(bRay.ExtrapOrigin());
2851
2852	++ splits;
2853
2854	if (side <= 0)
2855	{
2856	backRays.push_back(RayInfo(ray, bRay.GetMinT(), t));
2857	frontRays.push_back(RayInfo(ray, t, bRay.GetMaxT()));
2858	}
2859	else
2860	{
2861	frontRays.push_back(RayInfo(ray, bRay.GetMinT(), t));
2862	backRays.push_back(RayInfo(ray, t, bRay.GetMaxT()));
2863	}
2864	}
2865	break;
2866	default:
2867	Debug << "Should not come here" << endl;
2868	break;
2869	}
2870	}
2871
2872	return splits;
2873	}
2874
2875
2876	void VspBspTree::ExtractHalfSpaces(BspNode *n, vector<Plane3> &halfSpaces) const
2877	{
2878	BspNode *lastNode;
2879
2880	do
2881	{
2882	lastNode = n;
2883
2884	// want to get planes defining geometry of this node => don't take
2885	// split plane of node itself
2886	n = n->GetParent();
2887
2888	if (n)
2889	{
2890	BspInterior interior = dynamic_cast<BspInterior >(n);
2891	Plane3 halfSpace = dynamic_cast<BspInterior *>(interior)->GetPlane();
2892
2893	if (interior->GetBack() != lastNode)
2894	halfSpace.ReverseOrientation();
2895
2896	halfSpaces.push_back(halfSpace);
2897	}
2898	}
2899	while (n);
2900	}
2901
2902
2903	void VspBspTree::ConstructGeometry(BspNode *n,
2904	BspNodeGeometry &geom) const
2905	{
2906	vector<Plane3> halfSpaces;
2907	ExtractHalfSpaces(n, halfSpaces);
2908
2909	PolygonContainer candidatePolys;
2910	vector<Plane3> candidatePlanes;
2911
2912	vector<Plane3>::const_iterator pit, pit_end = halfSpaces.end();
2913
2914	// bounded planes are added to the polygons
2915	for (pit = halfSpaces.begin(); pit != pit_end; ++ pit)
2916	{
2917	Polygon3 p = GetBoundingBox().CrossSection(pit);
2918
2919	if (p->Valid(mEpsilon))
2920	{
2921	candidatePolys.push_back(p);
2922	candidatePlanes.push_back(*pit);
2923	}
2924	}
2925
2926	// add faces of bounding box (also could be faces of the cell)
2927	for (int i = 0; i < 6; ++ i)
2928	{
2929	VertexContainer vertices;
2930
2931	for (int j = 0; j < 4; ++ j)
2932	{
2933	vertices.push_back(mBoundingBox.GetFace(i).mVertices[j]);
2934	}
2935
2936	Polygon3 *poly = new Polygon3(vertices);
2937
2938	candidatePolys.push_back(poly);
2939	candidatePlanes.push_back(poly->GetSupportingPlane());
2940	}
2941
2942	for (int i = 0; i < (int)candidatePolys.size(); ++ i)
2943	{
2944	// polygon is split by all other planes
2945	for (int j = 0; (j < (int)halfSpaces.size()) && candidatePolys[i]; ++ j)
2946	{
2947	if (i == j) // polygon and plane are coincident
2948	continue;
2949
2950	VertexContainer splitPts;
2951	Polygon3 frontPoly, backPoly;
2952
2953	const int cf =
2954	candidatePolys[i]->ClassifyPlane(halfSpaces[j],
2955	mEpsilon);
2956
2957	switch (cf)
2958	{
2959	case Polygon3::SPLIT:
2960	frontPoly = new Polygon3();
2961	backPoly = new Polygon3();
2962
2963	candidatePolys[i]->Split(halfSpaces[j],
2964	*frontPoly,
2965	*backPoly,
2966	mEpsilon);
2967
2968	DEL_PTR(candidatePolys[i]);
2969
2970	if (backPoly->Valid(mEpsilon))
2971	{
2972	candidatePolys[i] = backPoly;
2973	}
2974	else
2975	{
2976	DEL_PTR(backPoly);
2977	}
2978
2979	// outside, don't need this
2980	DEL_PTR(frontPoly);
2981	break;
2982
2983	// polygon outside of halfspace
2984	case Polygon3::FRONT_SIDE:
2985	DEL_PTR(candidatePolys[i]);
2986	break;
2987
2988	// just take polygon as it is
2989	case Polygon3::BACK_SIDE:
2990	case Polygon3::COINCIDENT:
2991	default:
2992	break;
2993	}
2994	}
2995
2996	if (candidatePolys[i])
2997	{
2998	geom.Add(candidatePolys[i], candidatePlanes[i]);
2999	}
3000	}
3001	}
3002
3003
3004	bool VspBspTree::IsOutOfBounds(ViewCell *vc) const
3005	{
3006	return vc->GetId() == OUT_OF_BOUNDS_ID;
3007	}
3008
3009
3010	void VspBspTree::SetViewCellsTree(ViewCellsTree *viewCellsTree)
3011	{
3012	mViewCellsTree = viewCellsTree;
3013	}
3014
3015
3016	void VspBspTree::ConstructGeometry(ViewCell *vc,
3017	BspNodeGeometry &vcGeom) const
3018	{
3019	// if false, cannot construct geometry for interior leaf
3020	if (!mViewCellsTree)
3021	return;
3022
3023	ViewCellContainer leaves;
3024	mViewCellsTree->CollectLeaves(vc, leaves);
3025
3026	ViewCellContainer::const_iterator it, it_end = leaves.end();
3027
3028	for (it = leaves.begin(); it != it_end; ++ it)
3029	{
3030	if (IsOutOfBounds(*it))
3031	continue;
3032
3033	BspViewCell bspVc = dynamic_cast<BspViewCell >(*it);
3034	vector<BspLeaf *>::const_iterator bit, bit_end = bspVc->mLeaves.end();
3035
3036	for (bit = bspVc->mLeaves.begin(); bit != bit_end; ++ bit)
3037	{
3038	BspLeaf l = bit;
3039	ConstructGeometry(l, vcGeom);
3040	}
3041	}
3042	}
3043
3044
3045	int VspBspTree::FindNeighbors(BspNode n, vector<BspLeaf > &neighbors,
3046	const bool onlyUnmailed) const
3047	{
3048	stack<bspNodePair> nodeStack;
3049
3050	BspNodeGeometry nodeGeom;
3051	ConstructGeometry(n, nodeGeom);
3052	// const float eps = 0.5f;
3053	const float eps = 0.01f;
3054	// split planes from the root to this node
3055	// needed to verify that we found neighbor leaf
3056	// TODO: really needed?
3057	vector<Plane3> halfSpaces;
3058	ExtractHalfSpaces(n, halfSpaces);
3059
3060
3061	BspNodeGeometry *rgeom = new BspNodeGeometry();
3062	ConstructGeometry(mRoot, *rgeom);
3063
3064	nodeStack.push(bspNodePair(mRoot, rgeom));
3065
3066	while (!nodeStack.empty())
3067	{
3068	BspNode *node = nodeStack.top().first;
3069	BspNodeGeometry *geom = nodeStack.top().second;
3070
3071	nodeStack.pop();
3072
3073	if (node->IsLeaf())
3074	{
3075	// test if this leaf is in valid view space
3076	if (node->TreeValid() &&
3077	(node != n) &&
3078	(!onlyUnmailed \|\| !node->Mailed()))
3079	{
3080	bool isAdjacent = true;
3081
3082	if (1)
3083	{
3084	// test all planes of current node if still adjacent
3085	for (int i = 0; (i < halfSpaces.size()) && isAdjacent; ++ i)
3086	{
3087	const int cf =
3088	Polygon3::ClassifyPlane(geom->GetPolys(),
3089	halfSpaces[i],
3090	eps);
3091
3092	if (cf == Polygon3::FRONT_SIDE)
3093	{
3094	isAdjacent = false;
3095	}
3096	}
3097	}
3098	else
3099	{
3100	// TODO: why is this wrong??
3101	// test all planes of current node if still adjacent
3102	for (int i = 0; (i < nodeGeom.Size()) && isAdjacent; ++ i)
3103	{
3104	Polygon3 *poly = nodeGeom.GetPolys()[i];
3105
3106	const int cf =
3107	Polygon3::ClassifyPlane(geom->GetPolys(),
3108	poly->GetSupportingPlane(),
3109	eps);
3110
3111	if (cf == Polygon3::FRONT_SIDE)
3112	{
3113	isAdjacent = false;
3114	}
3115	}
3116	}
3117	// neighbor was found
3118	if (isAdjacent)
3119	{
3120	neighbors.push_back(dynamic_cast<BspLeaf *>(node));
3121	}
3122	}
3123	}
3124	else
3125	{
3126	BspInterior interior = dynamic_cast<BspInterior >(node);
3127
3128	const int cf = Polygon3::ClassifyPlane(nodeGeom.GetPolys(),
3129	interior->GetPlane(),
3130	eps);
3131
3132	BspNode *front = interior->GetFront();
3133	BspNode *back = interior->GetBack();
3134
3135	BspNodeGeometry *fGeom = new BspNodeGeometry();
3136	BspNodeGeometry *bGeom = new BspNodeGeometry();
3137
3138	geom->SplitGeometry(*fGeom,
3139	*bGeom,
3140	interior->GetPlane(),
3141	mBoundingBox,
3142	//0.0000001f);
3143	eps);
3144
3145	if (cf == Polygon3::BACK_SIDE)
3146	{
3147	nodeStack.push(bspNodePair(interior->GetBack(), bGeom));
3148	DEL_PTR(fGeom);
3149	}
3150	else
3151	{
3152	if (cf == Polygon3::FRONT_SIDE)
3153	{
3154	nodeStack.push(bspNodePair(interior->GetFront(), fGeom));
3155	DEL_PTR(bGeom);
3156	}
3157	else
3158	{ // random decision
3159	nodeStack.push(bspNodePair(front, fGeom));
3160	nodeStack.push(bspNodePair(back, bGeom));
3161	}
3162	}
3163	}
3164
3165	DEL_PTR(geom);
3166	}
3167
3168	return (int)neighbors.size();
3169	}
3170
3171
3172
3173	int VspBspTree::FindApproximateNeighbors(BspNode *n,
3174	vector<BspLeaf *> &neighbors,
3175	const bool onlyUnmailed) const
3176	{
3177	stack<bspNodePair> nodeStack;
3178
3179	BspNodeGeometry nodeGeom;
3180	ConstructGeometry(n, nodeGeom);
3181
3182	float eps = 0.01f;
3183	// split planes from the root to this node
3184	// needed to verify that we found neighbor leaf
3185	// TODO: really needed?
3186	vector<Plane3> halfSpaces;
3187	ExtractHalfSpaces(n, halfSpaces);
3188
3189
3190	BspNodeGeometry *rgeom = new BspNodeGeometry();
3191	ConstructGeometry(mRoot, *rgeom);
3192
3193	nodeStack.push(bspNodePair(mRoot, rgeom));
3194
3195	while (!nodeStack.empty())
3196	{
3197	BspNode *node = nodeStack.top().first;
3198	BspNodeGeometry *geom = nodeStack.top().second;
3199
3200	nodeStack.pop();
3201
3202	if (node->IsLeaf())
3203	{
3204	// test if this leaf is in valid view space
3205	if (node->TreeValid() &&
3206	(node != n) &&
3207	(!onlyUnmailed \|\| !node->Mailed()))
3208	{
3209	bool isAdjacent = true;
3210
3211	// neighbor was found
3212	if (isAdjacent)
3213	{
3214	neighbors.push_back(dynamic_cast<BspLeaf *>(node));
3215	}
3216	}
3217	}
3218	else
3219	{
3220	BspInterior interior = dynamic_cast<BspInterior >(node);
3221
3222	const int cf = Polygon3::ClassifyPlane(nodeGeom.GetPolys(),
3223	interior->GetPlane(),
3224	eps);
3225
3226	BspNode *front = interior->GetFront();
3227	BspNode *back = interior->GetBack();
3228
3229	BspNodeGeometry *fGeom = new BspNodeGeometry();
3230	BspNodeGeometry *bGeom = new BspNodeGeometry();
3231
3232	geom->SplitGeometry(*fGeom,
3233	*bGeom,
3234	interior->GetPlane(),
3235	mBoundingBox,
3236	//0.0000001f);
3237	eps);
3238
3239	if (cf == Polygon3::BACK_SIDE)
3240	{
3241	nodeStack.push(bspNodePair(interior->GetBack(), bGeom));
3242	DEL_PTR(fGeom);
3243	}
3244	else
3245	{
3246	if (cf == Polygon3::FRONT_SIDE)
3247	{
3248	nodeStack.push(bspNodePair(interior->GetFront(), fGeom));
3249	DEL_PTR(bGeom);
3250	}
3251	else
3252	{ // random decision
3253	nodeStack.push(bspNodePair(front, fGeom));
3254	nodeStack.push(bspNodePair(back, bGeom));
3255	}
3256	}
3257	}
3258
3259	DEL_PTR(geom);
3260	}
3261
3262	return (int)neighbors.size();
3263	}
3264
3265
3266
3267	BspLeaf *VspBspTree::GetRandomLeaf(const Plane3 &halfspace)
3268	{
3269	stack<BspNode *> nodeStack;
3270	nodeStack.push(mRoot);
3271
3272	int mask = rand();
3273
3274	while (!nodeStack.empty())
3275	{
3276	BspNode *node = nodeStack.top();
3277	nodeStack.pop();
3278
3279	if (node->IsLeaf())
3280	{
3281	return dynamic_cast<BspLeaf *>(node);
3282	}
3283	else
3284	{
3285	BspInterior interior = dynamic_cast<BspInterior >(node);
3286	BspNode *next;
3287	BspNodeGeometry geom;
3288
3289	// todo: not very efficient: constructs full cell everytime
3290	ConstructGeometry(interior, geom);
3291
3292	const int cf =
3293	Polygon3::ClassifyPlane(geom.GetPolys(), halfspace, mEpsilon);
3294
3295	if (cf == Polygon3::BACK_SIDE)
3296	next = interior->GetFront();
3297	else
3298	if (cf == Polygon3::FRONT_SIDE)
3299	next = interior->GetFront();
3300	else
3301	{
3302	// random decision
3303	if (mask & 1)
3304	next = interior->GetBack();
3305	else
3306	next = interior->GetFront();
3307	mask = mask >> 1;
3308	}
3309
3310	nodeStack.push(next);
3311	}
3312	}
3313
3314	return NULL;
3315	}
3316
3317
3318	BspLeaf *VspBspTree::GetRandomLeaf(const bool onlyUnmailed)
3319	{
3320	stack<BspNode *> nodeStack;
3321
3322	nodeStack.push(mRoot);
3323
3324	int mask = rand();
3325
3326	while (!nodeStack.empty())
3327	{
3328	BspNode *node = nodeStack.top();
3329	nodeStack.pop();
3330
3331	if (node->IsLeaf())
3332	{
3333	if ( (!onlyUnmailed \|\| !node->Mailed()) )
3334	return dynamic_cast<BspLeaf *>(node);
3335	}
3336	else
3337	{
3338	BspInterior interior = dynamic_cast<BspInterior >(node);
3339
3340	// random decision
3341	if (mask & 1)
3342	nodeStack.push(interior->GetBack());
3343	else
3344	nodeStack.push(interior->GetFront());
3345
3346	mask = mask >> 1;
3347	}
3348	}
3349
3350	return NULL;
3351	}
3352
3353
3354	int VspBspTree::ComputePvsSize(const RayInfoContainer &rays) const
3355	{
3356	int pvsSize = 0;
3357
3358	RayInfoContainer::const_iterator rit, rit_end = rays.end();
3359
3360	Intersectable::NewMail();
3361
3362	for (rit = rays.begin(); rit != rays.end(); ++ rit)
3363	{
3364	VssRay ray = (rit).mRay;
3365
3366	if (COUNT_ORIGIN_OBJECTS && ray->mOriginObject)
3367	{
3368	if (!ray->mOriginObject->Mailed())
3369	{
3370	ray->mOriginObject->Mail();
3371	++ pvsSize;
3372	}
3373	}
3374
3375	if (ray->mTerminationObject)
3376	{
3377	if (!ray->mTerminationObject->Mailed())
3378	{
3379	ray->mTerminationObject->Mail();
3380	++ pvsSize;
3381	}
3382	}
3383	}
3384
3385	return pvsSize;
3386	}
3387
3388
3389	float VspBspTree::GetEpsilon() const
3390	{
3391	return mEpsilon;
3392	}
3393
3394
3395	int VspBspTree::SplitPolygons(const Plane3 &plane,
3396	PolygonContainer &polys,
3397	PolygonContainer &frontPolys,
3398	PolygonContainer &backPolys,
3399	PolygonContainer &coincident) const
3400	{
3401	int splits = 0;
3402
3403	PolygonContainer::const_iterator it, it_end = polys.end();
3404
3405	for (it = polys.begin(); it != polys.end(); ++ it)
3406	{
3407	Polygon3 poly = it;
3408
3409	// classify polygon
3410	const int cf = poly->ClassifyPlane(plane, mEpsilon);
3411
3412	switch (cf)
3413	{
3414	case Polygon3::COINCIDENT:
3415	coincident.push_back(poly);
3416	break;
3417	case Polygon3::FRONT_SIDE:
3418	frontPolys.push_back(poly);
3419	break;
3420	case Polygon3::BACK_SIDE:
3421	backPolys.push_back(poly);
3422	break;
3423	case Polygon3::SPLIT:
3424	backPolys.push_back(poly);
3425	frontPolys.push_back(poly);
3426	++ splits;
3427	break;
3428	default:
3429	Debug << "SHOULD NEVER COME HERE\n";
3430	break;
3431	}
3432	}
3433
3434	return splits;
3435	}
3436
3437
3438	int VspBspTree::CastLineSegment(const Vector3 &origin,
3439	const Vector3 &termination,
3440	ViewCellContainer &viewcells)
3441	{
3442	int hits = 0;
3443	stack<BspRayTraversalData> tStack;
3444
3445	float mint = 0.0f, maxt = 1.0f;
3446
3447	//ViewCell::NewMail();
3448
3449	Vector3 entp = origin;
3450	Vector3 extp = termination;
3451
3452	BspNode *node = mRoot;
3453	BspNode *farChild = NULL;
3454
3455	float t;
3456	const float thresh = 1e-6f; // matt: change this to adjustable value
3457
3458	while (1)
3459	{
3460	if (!node->IsLeaf())
3461	{
3462	BspInterior in = dynamic_cast<BspInterior >(node);
3463
3464	Plane3 splitPlane = in->GetPlane();
3465
3466	const int entSide = splitPlane.Side(entp, thresh);
3467	const int extSide = splitPlane.Side(extp, thresh);
3468
3469	if (entSide < 0)
3470	{
3471	node = in->GetBack();
3472
3473	// plane does not split ray => no far child
3474	if (extSide <= 0)
3475	continue;
3476
3477	farChild = in->GetFront(); // plane splits ray
3478	}
3479	else if (entSide > 0)
3480	{
3481	node = in->GetFront();
3482
3483	if (extSide >= 0) // plane does not split ray => no far child
3484	continue;
3485
3486	farChild = in->GetBack(); // plane splits ray
3487	}
3488	else // one of the ray end points is on the plane
3489	{
3490	// NOTE: what to do if ray is coincident with plane?
3491	if (extSide < 0)
3492	node = in->GetBack();
3493	else //if (extSide > 0)
3494	node = in->GetFront();
3495	//else break; // coincident => count no intersections
3496
3497	continue; // no far child
3498	}
3499
3500	// push data for far child
3501	tStack.push(BspRayTraversalData(farChild, extp));
3502
3503	// find intersection of ray segment with plane
3504	extp = splitPlane.FindIntersection(origin, extp, &t);
3505	}
3506	else
3507	{
3508	// reached leaf => intersection with view cell
3509	BspLeaf leaf = dynamic_cast<BspLeaf >(node);
3510	ViewCell *viewCell;
3511
3512	// question: always contribute to leaf or to currently active view cell?
3513	if (0)
3514	viewCell = mViewCellsTree->GetActiveViewCell(leaf->GetViewCell());
3515	else
3516	viewCell = leaf->GetViewCell();
3517
3518	if (!viewCell->Mailed())
3519	{
3520	viewcells.push_back(viewCell);
3521	viewCell->Mail();
3522	++ hits;
3523	}
3524
3525	//-- fetch the next far child from the stack
3526	if (tStack.empty())
3527	break;
3528
3529	entp = extp;
3530
3531	const BspRayTraversalData &s = tStack.top();
3532
3533	node = s.mNode;
3534	extp = s.mExitPoint;
3535
3536	tStack.pop();
3537	}
3538	}
3539
3540	return hits;
3541	}
3542
3543
3544
3545
3546	int VspBspTree::TreeDistance(BspNode n1, BspNode n2) const
3547	{
3548	std::deque<BspNode *> path1;
3549	BspNode *p1 = n1;
3550
3551	// create path from node 1 to root
3552	while (p1)
3553	{
3554	if (p1 == n2) // second node on path
3555	return (int)path1.size();
3556
3557	path1.push_front(p1);
3558	p1 = p1->GetParent();
3559	}
3560
3561	int depth = n2->GetDepth();
3562	int d = depth;
3563
3564	BspNode *p2 = n2;
3565
3566	// compare with same depth
3567	while (1)
3568	{
3569	if ((d < (int)path1.size()) && (p2 == path1[d]))
3570	return (depth - d) + ((int)path1.size() - 1 - d);
3571
3572	-- d;
3573	p2 = p2->GetParent();
3574	}
3575
3576	return 0; // never come here
3577	}
3578
3579
3580	BspNode VspBspTree::CollapseTree(BspNode node, int &collapsed)
3581	{
3582	// TODO
3583	#if HAS_TO_BE_REDONE
3584	if (node->IsLeaf())
3585	return node;
3586
3587	BspInterior interior = dynamic_cast<BspInterior >(node);
3588
3589	BspNode *front = CollapseTree(interior->GetFront(), collapsed);
3590	BspNode *back = CollapseTree(interior->GetBack(), collapsed);
3591
3592	if (front->IsLeaf() && back->IsLeaf())
3593	{
3594	BspLeaf frontLeaf = dynamic_cast<BspLeaf >(front);
3595	BspLeaf backLeaf = dynamic_cast<BspLeaf >(back);
3596
3597	//-- collapse tree
3598	if (frontLeaf->GetViewCell() == backLeaf->GetViewCell())
3599	{
3600	BspViewCell *vc = frontLeaf->GetViewCell();
3601
3602	BspLeaf *leaf = new BspLeaf(interior->GetParent(), vc);
3603	leaf->SetTreeValid(frontLeaf->TreeValid());
3604
3605	// replace a link from node's parent
3606	if (leaf->GetParent())
3607	leaf->GetParent()->ReplaceChildLink(node, leaf);
3608	else
3609	mRoot = leaf;
3610
3611	++ collapsed;
3612	delete interior;
3613
3614	return leaf;
3615	}
3616	}
3617	#endif
3618	return node;
3619	}
3620
3621
3622	int VspBspTree::CollapseTree()
3623	{
3624	int collapsed = 0;
3625	//TODO
3626	#if HAS_TO_BE_REDONE
3627	(void)CollapseTree(mRoot, collapsed);
3628
3629	// revalidate leaves
3630	RepairViewCellsLeafLists();
3631	#endif
3632	return collapsed;
3633	}
3634
3635
3636	void VspBspTree::RepairViewCellsLeafLists()
3637	{
3638	// TODO
3639	#if HAS_TO_BE_REDONE
3640	// list not valid anymore => clear
3641	stack<BspNode *> nodeStack;
3642	nodeStack.push(mRoot);
3643
3644	ViewCell::NewMail();
3645
3646	while (!nodeStack.empty())
3647	{
3648	BspNode *node = nodeStack.top();
3649	nodeStack.pop();
3650
3651	if (node->IsLeaf())
3652	{
3653	BspLeaf leaf = dynamic_cast<BspLeaf >(node);
3654
3655	BspViewCell *viewCell = leaf->GetViewCell();
3656
3657	if (!viewCell->Mailed())
3658	{
3659	viewCell->mLeaves.clear();
3660	viewCell->Mail();
3661	}
3662
3663	viewCell->mLeaves.push_back(leaf);
3664
3665	}
3666	else
3667	{
3668	BspInterior interior = dynamic_cast<BspInterior >(node);
3669
3670	nodeStack.push(interior->GetFront());
3671	nodeStack.push(interior->GetBack());
3672	}
3673	}
3674	// TODO
3675	#endif
3676	}
3677
3678
3679	int VspBspTree::CastBeam(Beam &beam)
3680	{
3681	stack<bspNodePair> nodeStack;
3682	BspNodeGeometry *rgeom = new BspNodeGeometry();
3683	ConstructGeometry(mRoot, *rgeom);
3684
3685	nodeStack.push(bspNodePair(mRoot, rgeom));
3686
3687	ViewCell::NewMail();
3688
3689	while (!nodeStack.empty())
3690	{
3691	BspNode *node = nodeStack.top().first;
3692	BspNodeGeometry *geom = nodeStack.top().second;
3693	nodeStack.pop();
3694
3695	AxisAlignedBox3 box;
3696	geom->GetBoundingBox(box);
3697
3698	const int side = beam.ComputeIntersection(box);
3699
3700	switch (side)
3701	{
3702	case -1:
3703	CollectViewCells(node, true, beam.mViewCells, true);
3704	break;
3705	case 0:
3706
3707	if (node->IsLeaf())
3708	{
3709	BspLeaf leaf = dynamic_cast<BspLeaf >(node);
3710
3711	if (!leaf->GetViewCell()->Mailed() && leaf->TreeValid())
3712	{
3713	leaf->GetViewCell()->Mail();
3714	beam.mViewCells.push_back(leaf->GetViewCell());
3715	}
3716	}
3717	else
3718	{
3719	BspInterior interior = dynamic_cast<BspInterior >(node);
3720
3721	BspNode *first = interior->GetFront();
3722	BspNode *second = interior->GetBack();
3723
3724	BspNodeGeometry *firstGeom = new BspNodeGeometry();
3725	BspNodeGeometry *secondGeom = new BspNodeGeometry();
3726
3727	geom->SplitGeometry(*firstGeom,
3728	*secondGeom,
3729	interior->GetPlane(),
3730	mBoundingBox,
3731	//0.0000001f);
3732	mEpsilon);
3733
3734	// decide on the order of the nodes
3735	if (DotProd(beam.mPlanes[0].mNormal,
3736	interior->GetPlane().mNormal) > 0)
3737	{
3738	swap(first, second);
3739	swap(firstGeom, secondGeom);
3740	}
3741
3742	nodeStack.push(bspNodePair(first, firstGeom));
3743	nodeStack.push(bspNodePair(second, secondGeom));
3744	}
3745
3746	break;
3747	default:
3748	// default: cull
3749	break;
3750	}
3751
3752	DEL_PTR(geom);
3753
3754	}
3755
3756	return (int)beam.mViewCells.size();
3757	}
3758
3759
3760	void VspBspTree::SetViewCellsManager(ViewCellsManager *vcm)
3761	{
3762	mViewCellsManager = vcm;
3763	}
3764
3765
3766	int VspBspTree::CollectMergeCandidates(const vector<BspLeaf *> leaves,
3767	vector<MergeCandidate> &candidates)
3768	{
3769	BspLeaf::NewMail();
3770
3771	vector<BspLeaf *>::const_iterator it, it_end = leaves.end();
3772
3773	int numCandidates = 0;
3774
3775	// find merge candidates and push them into queue
3776	for (it = leaves.begin(); it != it_end; ++ it)
3777	{
3778	BspLeaf leaf = it;
3779
3780	// the same leaves must not be part of two merge candidates
3781	leaf->Mail();
3782
3783	vector<BspLeaf *> neighbors;
3784
3785	// appoximate neighbor search has slightl relaxed constraints
3786	if (1)
3787	FindNeighbors(leaf, neighbors, true);
3788	else
3789	FindApproximateNeighbors(leaf, neighbors, true);
3790
3791	vector<BspLeaf *>::const_iterator nit, nit_end = neighbors.end();
3792
3793	// TODO: test if at least one ray goes from one leaf to the other
3794	for (nit = neighbors.begin(); nit != nit_end; ++ nit)
3795	{
3796	if ((*nit)->GetViewCell() != leaf->GetViewCell())
3797	{
3798	MergeCandidate mc(leaf->GetViewCell(), (*nit)->GetViewCell());
3799
3800	if (!leaf->GetViewCell()->GetPvs().Empty() \|\|
3801	!(*nit)->GetViewCell()->GetPvs().Empty() \|\|
3802	leaf->IsSibling(*nit))
3803	{
3804	candidates.push_back(mc);
3805	}
3806
3807	++ numCandidates;
3808	if ((numCandidates % 1000) == 0)
3809	{
3810	cout << "collected " << numCandidates << " merge candidates" << endl;
3811	}
3812	}
3813	}
3814	}
3815
3816	Debug << "merge queue: " << (int)candidates.size() << endl;
3817	Debug << "leaves in queue: " << numCandidates << endl;
3818
3819
3820	return (int)leaves.size();
3821	}
3822
3823
3824	int VspBspTree::CollectMergeCandidates(const VssRayContainer &rays,
3825	vector<MergeCandidate> &candidates)
3826	{
3827	ViewCell::NewMail();
3828	long startTime = GetTime();
3829
3830	map<BspLeaf , vector<BspLeaf> > neighborMap;
3831	ViewCellContainer::const_iterator iit;
3832
3833	int numLeaves = 0;
3834
3835	BspLeaf::NewMail();
3836
3837	for (int i = 0; i < (int)rays.size(); ++ i)
3838	{
3839	VssRay *ray = rays[i];
3840
3841	// traverse leaves stored in the rays and compare and
3842	// merge consecutive leaves (i.e., the neighbors in the tree)
3843	if (ray->mViewCells.size() < 2)
3844	continue;
3845
3846	iit = ray->mViewCells.begin();
3847	BspViewCell bspVc = dynamic_cast<BspViewCell >(*(iit ++));
3848	BspLeaf *leaf = bspVc->mLeaves[0];
3849
3850	// traverse intersections
3851	// consecutive leaves are neighbors => add them to queue
3852	for (; iit != ray->mViewCells.end(); ++ iit)
3853	{
3854	// next pair
3855	BspLeaf *prevLeaf = leaf;
3856	bspVc = dynamic_cast<BspViewCell >(iit);
3857	leaf = bspVc->mLeaves[0]; // exactly one leaf
3858
3859	// view space not valid or same view cell
3860	if (!leaf->TreeValid() \|\| !prevLeaf->TreeValid() \|\|
3861	(leaf->GetViewCell() == prevLeaf->GetViewCell()))
3862	continue;
3863
3864	vector<BspLeaf *> &neighbors = neighborMap[leaf];
3865
3866	bool found = false;
3867
3868	// both leaves inserted in queue already =>
3869	// look if double pair already exists
3870	if (leaf->Mailed() && prevLeaf->Mailed())
3871	{
3872	vector<BspLeaf *>::const_iterator it, it_end = neighbors.end();
3873
3874	for (it = neighbors.begin(); !found && (it != it_end); ++ it)
3875	if (*it == prevLeaf)
3876	found = true; // already in queue
3877	}
3878
3879	if (!found)
3880	{
3881	// this pair is not in map yet
3882	// => insert into the neighbor map and the queue
3883	neighbors.push_back(prevLeaf);
3884	neighborMap[prevLeaf].push_back(leaf);
3885
3886	leaf->Mail();
3887	prevLeaf->Mail();
3888
3889	MergeCandidate mc(leaf->GetViewCell(), prevLeaf->GetViewCell());
3890
3891	candidates.push_back(mc);
3892
3893	if (((int)candidates.size() % 1000) == 0)
3894	{
3895	cout << "collected " << (int)candidates.size() << " merge candidates" << endl;
3896	}
3897	}
3898	}
3899	}
3900
3901	Debug << "neighbormap size: " << (int)neighborMap.size() << endl;
3902	Debug << "merge queue: " << (int)candidates.size() << endl;
3903	Debug << "leaves in queue: " << numLeaves << endl;
3904
3905
3906	//-- collect the leaves which haven't been found by ray casting
3907	if (0)
3908	{
3909	cout << "finding additional merge candidates using geometry" << endl;
3910	vector<BspLeaf *> leaves;
3911	CollectLeaves(leaves, true);
3912	Debug << "found " << (int)leaves.size() << " new leaves" << endl << endl;
3913	CollectMergeCandidates(leaves, candidates);
3914	}
3915
3916	return numLeaves;
3917	}
3918
3919
3920
3921
3922	ViewCell *VspBspTree::GetViewCell(const Vector3 &point, const bool active)
3923	{
3924	if (mRoot == NULL)
3925	return NULL;
3926
3927	stack<BspNode *> nodeStack;
3928	nodeStack.push(mRoot);
3929
3930	ViewCellLeaf *viewcell = NULL;
3931
3932	while (!nodeStack.empty())
3933	{
3934	BspNode *node = nodeStack.top();
3935	nodeStack.pop();
3936
3937	if (node->IsLeaf())
3938	{
3939	viewcell = dynamic_cast<BspLeaf *>(node)->GetViewCell();
3940	break;
3941	}
3942	else
3943	{
3944	BspInterior interior = dynamic_cast<BspInterior >(node);
3945
3946	// random decision
3947	if (interior->GetPlane().Side(point) < 0)
3948	nodeStack.push(interior->GetBack());
3949	else
3950	nodeStack.push(interior->GetFront());
3951	}
3952	}
3953
3954	if (active)
3955	return mViewCellsTree->GetActiveViewCell(viewcell);
3956	else
3957	return viewcell;
3958	}
3959
3960
3961	bool VspBspTree::ViewPointValid(const Vector3 &viewPoint) const
3962	{
3963	BspNode *node = mRoot;
3964
3965	while (1)
3966	{
3967	// early exit
3968	if (node->TreeValid())
3969	return true;
3970
3971	if (node->IsLeaf())
3972	return false;
3973
3974	BspInterior in = dynamic_cast<BspInterior >(node);
3975
3976	if (in->GetPlane().Side(viewPoint) <= 0)
3977	{
3978	node = in->GetBack();
3979	}
3980	else
3981	{
3982	node = in->GetFront();
3983	}
3984	}
3985
3986	// should never come here
3987	return false;
3988	}
3989
3990
3991	void VspBspTree::PropagateUpValidity(BspNode *node)
3992	{
3993	const bool isValid = node->TreeValid();
3994
3995	// propagative up invalid flag until only invalid nodes exist over this node
3996	if (!isValid)
3997	{
3998	while (!node->IsRoot() && node->GetParent()->TreeValid())
3999	{
4000	node = node->GetParent();
4001	node->SetTreeValid(false);
4002	}
4003	}
4004	else
4005	{
4006	// propagative up valid flag until one of the subtrees is invalid
4007	while (!node->IsRoot() && !node->TreeValid())
4008	{
4009	node = node->GetParent();
4010	BspInterior interior = dynamic_cast<BspInterior >(node);
4011
4012	// the parent is valid iff both leaves are valid
4013	node->SetTreeValid(interior->GetBack()->TreeValid() &&
4014	interior->GetFront()->TreeValid());
4015	}
4016	}
4017	}
4018
4019
4020	bool VspBspTree::Export(OUT_STREAM &stream)
4021	{
4022	ExportNode(mRoot, stream);
4023	return true;
4024	}
4025
4026
4027	void VspBspTree::ExportNode(BspNode *node, OUT_STREAM &stream)
4028	{
4029	if (node->IsLeaf())
4030	{
4031	BspLeaf leaf = dynamic_cast<BspLeaf >(node);
4032	ViewCell *viewCell = mViewCellsTree->GetActiveViewCell(leaf->GetViewCell());
4033
4034	int id = -1;
4035	if (viewCell != mOutOfBoundsCell)
4036	id = viewCell->GetId();
4037
4038	stream << "<Leaf viewCellId=\"" << id << "\" />" << endl;
4039	}
4040	else
4041	{
4042	BspInterior interior = dynamic_cast<BspInterior >(node);
4043
4044	Plane3 plane = interior->GetPlane();
4045	stream << "<Interior plane=\"" << plane.mNormal.x << " "
4046	<< plane.mNormal.y << " " << plane.mNormal.z << " "
4047	<< plane.mD << "\">" << endl;
4048
4049	ExportNode(interior->GetBack(), stream);
4050	ExportNode(interior->GetFront(), stream);
4051
4052	stream << "</Interior>" << endl;
4053	}
4054	}
4055
4056	}

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

Download in other formats: