Context Navigation

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

Revision 574, 86.9 KB checked in by mattausch, 18 years ago (diff)
finished function for view cell construction removed bsp rays from vspbspmanager

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

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

Download in other formats: