Context Navigation

source: trunk/VUT/GtpVisibilityPreprocessor/src/ViewCellBsp.cpp @ 426

Revision 426, 63.1 KB checked in by mattausch, 19 years ago (diff)
fixed ray bug in vspkdtree added visualizations

Line
1	#include "Plane3.h"
2	#include "ViewCellBsp.h"
3	#include "Mesh.h"
4	#include "common.h"
5	#include "ViewCell.h"
6	#include "Environment.h"
7	#include "Polygon3.h"
8	#include "Ray.h"
9	#include "AxisAlignedBox3.h"
10	#include <stack>
11	#include <time.h>
12	#include <iomanip>
13	#include "Exporter.h"
14	#include "Plane3.h"
15
16	//-- static members
17
18	int BspNode::sMailId = 1;
19
20	/** Evaluates split plane classification with respect to the plane's
21	contribution for a minimum number splits in the tree.
22	*/
23	const float BspTree::sLeastPolySplitsTable[] = {0, 0, 1, 0};
24	/** Evaluates split plane classification with respect to the plane's
25	contribution for a balanced tree.
26	*/
27	const float BspTree::sBalancedPolysTable[] = {1, -1, 0, 0};
28
29	/** Evaluates split plane classification with respect to the plane's
30	contribution for a minimum number of ray splits.
31	*/
32	const float BspTree::sLeastRaySplitsTable[] = {0, 0, 1, 1, 0};
33	/** Evaluates split plane classification with respect to the plane's
34	contribution for balanced rays.
35	*/
36	const float BspTree::sBalancedRaysTable[] = {1, -1, 0, 0, 0};
37
38	int BspTree::sConstructionMethod = FROM_INPUT_VIEW_CELLS;
39
40
41	/****************************************************************/
42	/* class BspNode implementation */
43	/****************************************************************/
44
45	BspNode::BspNode():
46	mParent(NULL)
47	{}
48
49	BspNode::BspNode(BspInterior *parent):
50	mParent(parent)
51	{}
52
53
54	bool BspNode::IsRoot() const
55	{
56	return mParent == NULL;
57	}
58
59	BspInterior *BspNode::GetParent()
60	{
61	return mParent;
62	}
63
64	void BspNode::SetParent(BspInterior *parent)
65	{
66	mParent = parent;
67	}
68
69	/****************************************************************/
70	/* class BspInterior implementation */
71	/****************************************************************/
72
73
74	BspInterior::BspInterior(const Plane3 &plane):
75	mPlane(plane), mFront(NULL), mBack(NULL)
76	{}
77
78	BspInterior::~BspInterior()
79	{
80	DEL_PTR(mFront);
81	DEL_PTR(mBack);
82	}
83
84	bool BspInterior::IsLeaf() const
85	{
86	return false;
87	}
88
89	BspNode *BspInterior::GetBack()
90	{
91	return mBack;
92	}
93
94	BspNode *BspInterior::GetFront()
95	{
96	return mFront;
97	}
98
99	Plane3 *BspInterior::GetPlane()
100	{
101	return &mPlane;
102	}
103
104	void BspInterior::ReplaceChildLink(BspNode oldChild, BspNode newChild)
105	{
106	if (mBack == oldChild)
107	mBack = newChild;
108	else
109	mFront = newChild;
110	}
111
112	void BspInterior::SetupChildLinks(BspNode b, BspNode f)
113	{
114	mBack = b;
115	mFront = f;
116	}
117
118	int BspInterior::SplitPolygons(PolygonContainer &polys,
119	PolygonContainer &frontPolys,
120	PolygonContainer &backPolys,
121	PolygonContainer &coincident)
122	{
123	Polygon3 *splitPoly = NULL;
124
125	int splits = 0;
126
127	#ifdef _Debug
128	Debug << "splitting polygons of node " << this << " with plane " << mPlane << endl;
129	#endif
130	while (!polys.empty())
131	{
132	Polygon3 *poly = polys.back();
133	polys.pop_back();
134
135	//Debug << "New polygon with plane: " << poly->GetSupportingPlane() << "\n";
136
137	// classify polygon
138	const int cf = poly->ClassifyPlane(mPlane);
139
140	Polygon3 *front_piece = NULL;
141	Polygon3 *back_piece = NULL;
142
143	VertexContainer splitVertices;
144
145	switch (cf)
146	{
147	case Polygon3::COINCIDENT:
148	coincident.push_back(poly);
149	break;
150	case Polygon3::FRONT_SIDE:
151	frontPolys.push_back(poly);
152	break;
153	case Polygon3::BACK_SIDE:
154	backPolys.push_back(poly);
155	break;
156	case Polygon3::SPLIT:
157	front_piece = new Polygon3(poly->mParent);
158	back_piece = new Polygon3(poly->mParent);
159
160	//-- split polygon into front and back part
161	poly->Split(mPlane,
162	*front_piece,
163	*back_piece,
164	splitVertices);
165
166	++ splits; // increase number of splits
167
168	//-- inherit rays from parent polygon for blocked ray criterium
169	poly->InheritRays(front_piece, back_piece);
170	//Debug << "p: " << poly->mPiercingRays.size() << " f: " << front_piece->mPiercingRays.size() << " b: " << back_piece->mPiercingRays.size() << endl;
171
172	// check if polygons still valid
173	if (front_piece->Valid())
174	frontPolys.push_back(front_piece);
175	else
176	DEL_PTR(front_piece);
177
178	if (back_piece->Valid())
179	backPolys.push_back(back_piece);
180	else
181	DEL_PTR(back_piece);
182
183	#ifdef _DEBUG
184	Debug << "split " << poly << endl << front_piece << endl << *back_piece << endl;
185	#endif
186	DEL_PTR(poly);
187	break;
188	default:
189	Debug << "SHOULD NEVER COME HERE\n";
190	break;
191	}
192	}
193
194	return splits;
195	}
196
197	/****************************************************************/
198	/* class BspLeaf implementation */
199	/****************************************************************/
200	BspLeaf::BspLeaf(): mViewCell(NULL)
201	{
202	}
203
204	BspLeaf::BspLeaf(BspViewCell *viewCell):
205	mViewCell(viewCell)
206	{
207	}
208
209	BspLeaf::BspLeaf(BspInterior *parent):
210	BspNode(parent), mViewCell(NULL)
211	{}
212
213	BspLeaf::BspLeaf(BspInterior parent, BspViewCell viewCell):
214	BspNode(parent), mViewCell(viewCell)
215	{
216	}
217
218	BspViewCell *BspLeaf::GetViewCell() const
219	{
220	return mViewCell;
221	}
222
223	void BspLeaf::SetViewCell(BspViewCell *viewCell)
224	{
225	mViewCell = viewCell;
226	}
227
228	bool BspLeaf::IsLeaf() const
229	{
230	return true;
231	}
232
233	void BspLeaf::AddToPvs(const BoundedRayContainer &rays,
234	int &sampleContributions,
235	int &contributingSamples,
236	bool storeLeavesWithRays)
237	{
238	sampleContributions = 0;
239	contributingSamples = 0;
240
241	BoundedRayContainer::const_iterator it, it_end = rays.end();
242
243	// add contributions from samples to the PVS
244	for (it = rays.begin(); it != it_end; ++ it)
245	{
246	int contribution = 0;
247	Ray ray = (it)->mRay;
248
249	if (!ray->intersections.empty())
250	contribution += mViewCell->GetPvs().AddSample(ray->intersections[0].mObject);
251
252	if (ray->sourceObject.mObject)
253	contribution += mViewCell->GetPvs().AddSample(ray->sourceObject.mObject);
254
255	if (contribution > 0)
256	{
257	sampleContributions += contribution;
258	++ contributingSamples;
259	}
260
261	if (storeLeavesWithRays)
262	// warning: intersections not ordered
263	ray->bspIntersections.push_back(Ray::BspIntersection((*it)->mMinT, this));
264	}
265	}
266
267
268	/****************************************************************/
269	/* class BspTree implementation */
270	/****************************************************************/
271
272	BspTree::BspTree(BspViewCell *viewCell):
273	mRootCell(viewCell),
274	mRoot(NULL),
275	mGenerateViewCells(false),
276	mStoreLeavesWithRays(false),
277	mPvsUseArea(true)
278	{
279	Randomize(); // initialise random generator for heuristics
280
281	//-- termination criteria for autopartition
282	environment->GetIntValue("BspTree.Termination.maxDepth", mTermMaxDepth);
283	environment->GetIntValue("BspTree.Termination.minPvs", mTermMinPvs);
284	environment->GetIntValue("BspTree.Termination.minPolygons", mTermMinPolys);
285	environment->GetIntValue("BspTree.Termination.minRays", mTermMinRays);
286	environment->GetFloatValue("BspTree.Termination.minArea", mTermMinArea);
287	environment->GetFloatValue("BspTree.Termination.maxRayContribution", mTermMaxRayContribution);
288	environment->GetFloatValue("BspTree.Termination.minAccRayLenght", mTermMinAccRayLength);
289
290	//-- factors for bsp tree split plane heuristics
291	environment->GetFloatValue("BspTree.Factor.verticalSplits", mVerticalSplitsFactor);
292	environment->GetFloatValue("BspTree.Factor.largestPolyArea", mLargestPolyAreaFactor);
293	environment->GetFloatValue("BspTree.Factor.blockedRays", mBlockedRaysFactor);
294	environment->GetFloatValue("BspTree.Factor.leastRaySplits", mLeastRaySplitsFactor);
295	environment->GetFloatValue("BspTree.Factor.balancedRays", mBalancedRaysFactor);
296	environment->GetFloatValue("BspTree.Factor.pvsFactor", mPvsFactor);
297	environment->GetFloatValue("BspTree.Factor.leastSplits" , mLeastSplitsFactor);
298	environment->GetFloatValue("BspTree.Factor.balancedPolys", mBalancedPolysFactor);
299	environment->GetFloatValue("BspTree.Factor.balancedViewCells", mBalancedViewCellsFactor);
300	environment->GetFloatValue("BspTree.Termination.ct_div_ci", mCtDivCi);
301
302	//-- termination criteria for axis aligned split
303	environment->GetFloatValue("BspTree.Termination.AxisAligned.ct_div_ci", mAaCtDivCi);
304	environment->GetFloatValue("BspTree.Termination.AxisAligned.maxCostRatio", mMaxCostRatio);
305	environment->GetIntValue("BspTree.Termination.AxisAligned.minPolys",
306	mTermMinPolysForAxisAligned);
307	environment->GetIntValue("BspTree.Termination.AxisAligned.maxRays",
308	mTermMinRaysForAxisAligned);
309	environment->GetIntValue("BspTree.Termination.AxisAligned.maxObjects",
310	mTermMinObjectsForAxisAligned);
311	//-- partition criteria
312	environment->GetIntValue("BspTree.maxPolyCandidates", mMaxPolyCandidates);
313	environment->GetIntValue("BspTree.maxRayCandidates", mMaxRayCandidates);
314	environment->GetIntValue("BspTree.splitPlaneStrategy", mSplitPlaneStrategy);
315	environment->GetFloatValue("BspTree.AxisAligned.splitBorder", mSplitBorder);
316
317	environment->GetFloatValue("BspTree.Construction.sideTolerance", Vector3::sDistTolerance);
318	Vector3::sDistToleranceSqrt = Vector3::sDistTolerance * Vector3::sDistTolerance;
319
320	// post processing stuff
321	environment->GetIntValue("ViewCells.PostProcessing.minPvsDif", mMinPvsDif);
322	environment->GetIntValue("ViewCells.PostProcessing.minPvs", mMinPvs);
323	environment->GetIntValue("ViewCells.PostProcessing.maxPvs", mMaxPvs);
324
325	Debug << "BSP max depth: " << mTermMaxDepth << endl;
326	Debug << "BSP min PVS: " << mTermMinPvs << endl;
327	Debug << "BSP min area: " << mTermMinArea << endl;
328	Debug << "BSP max polys: " << mTermMinPolys << endl;
329	Debug << "BSP max rays: " << mTermMinRays << endl;
330	Debug << "BSP max polygon candidates: " << mMaxPolyCandidates << endl;
331	Debug << "BSP max plane candidates: " << mMaxRayCandidates << endl;
332
333	Debug << "Split plane strategy: ";
334	if (mSplitPlaneStrategy & RANDOM_POLYGON)
335	Debug << "random polygon ";
336	if (mSplitPlaneStrategy & AXIS_ALIGNED)
337	Debug << "axis aligned ";
338	if (mSplitPlaneStrategy & LEAST_SPLITS)
339	Debug << "least splits ";
340	if (mSplitPlaneStrategy & BALANCED_POLYS)
341	Debug << "balanced polygons ";
342	if (mSplitPlaneStrategy & BALANCED_VIEW_CELLS)
343	Debug << "balanced view cells ";
344	if (mSplitPlaneStrategy & LARGEST_POLY_AREA)
345	Debug << "largest polygon area ";
346	if (mSplitPlaneStrategy & VERTICAL_AXIS)
347	Debug << "vertical axis ";
348	if (mSplitPlaneStrategy & BLOCKED_RAYS)
349	Debug << "blocked rays ";
350	if (mSplitPlaneStrategy & LEAST_RAY_SPLITS)
351	Debug << "least ray splits ";
352	if (mSplitPlaneStrategy & BALANCED_RAYS)
353	Debug << "balanced rays ";
354	if (mSplitPlaneStrategy & PVS)
355	Debug << "pvs";
356
357	Debug << endl;
358	}
359
360	void BspTree::ParseEnvironment()
361	{
362	//-- parse bsp cell tree construction method
363	char constructionMethodStr[60];
364
365	environment->GetStringValue("BspTree.Construction.input", constructionMethodStr);
366
367	sConstructionMethod = FROM_INPUT_VIEW_CELLS;
368
369	if (strcmp(constructionMethodStr, "fromViewCells") == 0)
370	sConstructionMethod = FROM_INPUT_VIEW_CELLS;
371	else if (strcmp(constructionMethodStr, "fromSceneGeometry") == 0)
372	sConstructionMethod = FROM_SCENE_GEOMETRY;
373	else if (strcmp(constructionMethodStr, "fromSamples") == 0)
374	sConstructionMethod = FROM_SAMPLES;
375	else
376	{
377	cerr << "Wrong construction method " << constructionMethodStr << endl;
378	exit(1);
379	}
380
381	Debug << "Construction method: " << constructionMethodStr << endl;
382	}
383
384	const BspTreeStatistics &BspTree::GetStatistics() const
385	{
386	return mStat;
387	}
388
389	void BspTreeStatistics::Print(ostream &app) const
390	{
391	app << "===== BspTree statistics ===============\n";
392
393	app << setprecision(4);
394
395	app << "#N_CTIME ( Construction time [s] )\n" << Time() << " \n";
396
397	app << "#N_NODES ( Number of nodes )\n" << nodes << "\n";
398
399	app << "#N_INTERIORS ( Number of interior nodes )\n" << Interior() << "\n";
400
401	app << "#N_LEAVES ( Number of leaves )\n" << Leaves() << "\n";
402
403	app << "#N_SPLITS ( Number of splits )\n" << splits << "\n";
404
405	app << "#N_PMAXDEPTHLEAVES ( Percentage of leaves at maximum depth )\n"
406	<< maxDepthNodes * 100 / (double)Leaves() << endl;
407
408	app << "#N_PMAXDEPTHLEAVES ( Percentage of leaves at maximum depth )\n"
409	<< maxDepthNodes * 100 / (double)Leaves() << endl;
410
411	app << "#N_PMINPVSLEAVES ( Percentage of leaves with mininimal PVS )\n"
412	<< minPvsNodes * 100 / (double)Leaves() << endl;
413
414	app << "#N_PMINRAYSLEAVES ( Percentage of leaves with minimal number of rays)\n"
415	<< minRaysNodes * 100 / (double)Leaves() << endl;
416
417	app << "#N_PMINAREALEAVES ( Percentage of leaves with mininum area )\n"
418	<< minAreaNodes * 100 / (double)Leaves() << endl;
419
420	app << "#N_PMAXRAYCONTRIBLEAVES ( Percentage of leaves with maximal ray contribution )\n"
421	<< maxRayContribNodes * 100 / (double)Leaves() << endl;
422
423	app << "#N_PMAXDEPTH ( Maximal reached depth )\n" << maxDepth << endl;
424
425	app << "#N_PMINDEPTH ( Minimal reached depth )\n" << minDepth << endl;
426
427	app << "#AVGDEPTH ( average depth )\n" << AvgDepth() << endl;
428
429	app << "#N_INPUT_POLYGONS (number of input polygons )\n" << polys << endl;
430
431	//app << "#N_PVS: " << pvs << endl;
432
433	app << "#N_ROUTPUT_INPUT_POLYGONS ( ratio polygons after subdivision / input polygons )\n" <<
434	(polys + splits) / (double)polys << endl;
435
436	app << "===== END OF BspTree statistics ==========\n";
437	}
438
439
440	BspTree::~BspTree()
441	{
442	DEL_PTR(mRoot);
443	}
444
445	void BspTree::InsertViewCell(ViewCell *viewCell)
446	{
447	PolygonContainer *polys = new PolygonContainer();
448
449	// extract polygons that guide the split process
450	mStat.polys += AddMeshToPolygons(viewCell->GetMesh(), *polys, viewCell);
451	mBox.Include(viewCell->GetBox()); // add to BSP aabb
452
453	InsertPolygons(polys);
454	}
455
456	void BspTree::InsertPolygons(PolygonContainer *polys)
457	{
458	std::stack<BspTraversalData> tStack;
459
460	// traverse existing tree or create new tree
461	if (!mRoot)
462	mRoot = new BspLeaf();
463
464	tStack.push(BspTraversalData(mRoot, polys, 0, mRootCell, new BoundedRayContainer(), 0,
465	mBox.SurfaceArea(), new BspNodeGeometry()));
466
467	while (!tStack.empty())
468	{
469	// filter polygons donw the tree
470	BspTraversalData tData = tStack.top();
471	tStack.pop();
472
473	if (!tData.mNode->IsLeaf())
474	{
475	BspInterior interior = dynamic_cast<BspInterior >(tData.mNode);
476
477	//-- filter view cell polygons down the tree until a leaf is reached
478	if (!tData.mPolygons->empty())
479	{
480	PolygonContainer *frontPolys = new PolygonContainer();
481	PolygonContainer *backPolys = new PolygonContainer();
482	PolygonContainer coincident;
483
484	int splits = 0;
485
486	// split viewcell polygons with respect to split plane
487	splits += interior->SplitPolygons(*tData.mPolygons,
488	*frontPolys,
489	*backPolys,
490	coincident);
491
492	// extract view cells associated with the split polygons
493	ViewCell *frontViewCell = mRootCell;
494	ViewCell *backViewCell = mRootCell;
495
496	BspTraversalData frontData(interior->GetFront(),
497	frontPolys,
498	tData.mDepth + 1,
499	mRootCell,
500	tData.mRays,
501	tData.mPvs,
502	mBox.SurfaceArea(),
503	new BspNodeGeometry());
504
505	BspTraversalData backData(interior->GetBack(),
506	backPolys,
507	tData.mDepth + 1,
508	mRootCell,
509	tData.mRays,
510	tData.mPvs,
511	mBox.SurfaceArea(),
512	new BspNodeGeometry());
513
514	if (!mGenerateViewCells)
515	{
516	ExtractViewCells(frontData,
517	backData,
518	coincident,
519	interior->mPlane);
520	}
521
522	// don't need coincident polygons anymore
523	CLEAR_CONTAINER(coincident);
524
525	mStat.splits += splits;
526
527	// push the children on the stack
528	tStack.push(frontData);
529	tStack.push(backData);
530	}
531
532	// cleanup
533	DEL_PTR(tData.mPolygons);
534	}
535	else
536	{
537	// reached leaf => subdivide current viewcell
538	BspNode *subRoot = Subdivide(tStack, tData);
539	}
540	}
541	}
542
543	int BspTree::AddMeshToPolygons(Mesh mesh, PolygonContainer &polys, MeshInstance parent)
544	{
545	FaceContainer::const_iterator fi;
546
547	// copy the face data to polygons
548	for (fi = mesh->mFaces.begin(); fi != mesh->mFaces.end(); ++ fi)
549	{
550	Polygon3 poly = new Polygon3((fi), mesh);
551
552	if (poly->Valid())
553	{
554	poly->mParent = parent; // set parent intersectable
555	polys.push_back(poly);
556	}
557	else
558	DEL_PTR(poly);
559	}
560	return (int)mesh->mFaces.size();
561	}
562
563	int BspTree::AddToPolygonSoup(const ViewCellContainer &viewCells,
564	PolygonContainer &polys,
565	int maxObjects)
566	{
567	int limit = (maxObjects > 0) ?
568	Min((int)viewCells.size(), maxObjects) : (int)viewCells.size();
569
570	int polysSize = 0;
571
572	for (int i = 0; i < limit; ++ i)
573	{
574	if (viewCells[i]->GetMesh()) // copy the mesh data to polygons
575	{
576	mBox.Include(viewCells[i]->GetBox()); // add to BSP tree aabb
577	polysSize += AddMeshToPolygons(viewCells[i]->GetMesh(), polys, viewCells[i]);
578	}
579	}
580
581	return polysSize;
582	}
583
584	int BspTree::AddToPolygonSoup(const ObjectContainer &objects, PolygonContainer &polys, int maxObjects)
585	{
586	int limit = (maxObjects > 0) ? Min((int)objects.size(), maxObjects) : (int)objects.size();
587
588	for (int i = 0; i < limit; ++i)
589	{
590	Intersectable object = objects[i];//it;
591	Mesh *mesh = NULL;
592
593	switch (object->Type()) // extract the meshes
594	{
595	case Intersectable::MESH_INSTANCE:
596	mesh = dynamic_cast<MeshInstance *>(object)->GetMesh();
597	break;
598	case Intersectable::VIEW_CELL:
599	mesh = dynamic_cast<ViewCell *>(object)->GetMesh();
600	break;
601	// TODO: handle transformed mesh instances
602	default:
603	Debug << "intersectable type not supported" << endl;
604	break;
605	}
606
607	if (mesh) // copy the mesh data to polygons
608	{
609	mBox.Include(object->GetBox()); // add to BSP tree aabb
610	AddMeshToPolygons(mesh, polys, mRootCell);
611	}
612	}
613
614	return (int)polys.size();
615	}
616
617	void BspTree::Construct(const ViewCellContainer &viewCells)
618	{
619	mStat.nodes = 1;
620	mBox.Initialize(); // initialise bsp tree bounding box
621
622	// copy view cell meshes into one big polygon soup
623	PolygonContainer *polys = new PolygonContainer();
624	mStat.polys = AddToPolygonSoup(viewCells, *polys);
625
626	// construct tree from the view cell polygons
627	Construct(polys, new BoundedRayContainer());
628	}
629
630
631	void BspTree::Construct(const ObjectContainer &objects)
632	{
633	mStat.nodes = 1;
634	mBox.Initialize(); // initialise bsp tree bounding box
635
636	PolygonContainer *polys = new PolygonContainer();
637
638	// copy mesh instance polygons into one big polygon soup
639	mStat.polys = AddToPolygonSoup(objects, *polys);
640
641	// construct tree from polygon soup
642	Construct(polys, new BoundedRayContainer());
643	}
644
645	void BspTree::Construct(const RayContainer &sampleRays)
646	{
647	mStat.nodes = 1;
648	mBox.Initialize(); // initialise BSP tree bounding box
649
650	PolygonContainer *polys = new PolygonContainer();
651	BoundedRayContainer *rays = new BoundedRayContainer();
652
653	RayContainer::const_iterator rit, rit_end = sampleRays.end();
654
655	long startTime = GetTime();
656
657	Debug << "** Extracting polygons from rays **\n";
658
659	std::map<Face , Polygon3 > facePolyMap;
660
661	//-- extract polygons intersected by the rays
662	for (rit = sampleRays.begin(); rit != rit_end; ++ rit)
663	{
664	Ray ray = rit;
665
666	// get ray-face intersection. Store polygon representing the rays together
667	// with rays intersecting the face.
668	if (!ray->intersections.empty())
669	{
670	MeshInstance obj = dynamic_cast<MeshInstance >(ray->intersections[0].mObject);
671	Face *face = obj->GetMesh()->mFaces[ray->intersections[0].mFace];
672
673	std::map<Face , Polygon3 >::iterator it = facePolyMap.find(face);
674
675	if (it != facePolyMap.end())
676	{
677	//store rays if needed for heuristics
678	if (mSplitPlaneStrategy & BLOCKED_RAYS)
679	(*it).second->mPiercingRays.push_back(ray);
680	}
681	else
682	{ //store rays if needed for heuristics
683	Polygon3 *poly = new Polygon3(face, obj->GetMesh());
684	poly->mParent = obj;
685	polys->push_back(poly);
686
687	if (mSplitPlaneStrategy & BLOCKED_RAYS)
688	poly->mPiercingRays.push_back(ray);
689
690	facePolyMap[face] = poly;
691	}
692	}
693	}
694
695	facePolyMap.clear();
696
697	// compue bounding box
698	Polygon3::IncludeInBox(*polys, mBox);
699
700	//-- store rays
701	for (rit = sampleRays.begin(); rit != rit_end; ++ rit)
702	{
703	Ray ray = rit;
704	ray->SetId(-1); // reset id
705
706	float minT, maxT;
707	if (BoundRay(*ray, minT, maxT))
708	rays->push_back(new BoundedRay(ray, minT, maxT));
709	}
710
711	mStat.polys = (int)polys->size();
712
713	Debug << "** Finished polygon extraction **" << endl;
714	Debug << (int)polys->size() << " polys extracted from " << (int)sampleRays.size() << " rays" << endl;
715	Debug << "extraction time: " << TimeDiff(startTime, GetTime())*1e-3 << "s" << endl;
716
717	Construct(polys, rays);
718	}
719
720	void BspTree::Construct(PolygonContainer polys, BoundedRayContainer rays)
721	{
722	std::stack<BspTraversalData> tStack;
723
724	mRoot = new BspLeaf();
725
726	BspNodeGeometry *cell = new BspNodeGeometry();
727	ConstructGeometry(mRoot, *cell);
728
729	BspTraversalData tData(mRoot, polys, 0, mRootCell, rays,
730	ComputePvsSize(*rays), cell->GetArea(), cell);
731
732	tStack.push(tData);
733
734	mStat.Start();
735	cout << "** Contructing bsp tree **\n";
736
737	while (!tStack.empty())
738	{
739	tData = tStack.top();
740	tStack.pop();
741
742	// subdivide leaf node
743	BspNode *subRoot = Subdivide(tStack, tData);
744	}
745
746	mStat.Stop();
747	}
748
749	bool BspTree::TerminationCriteriaMet(const BspTraversalData &data) const
750	{
751	return
752	(((int)data.mPolygons->size() <= mTermMinPolys) \|\|
753	((int)data.mRays->size() <= mTermMinRays) \|\|
754	(data.mPvs <= mTermMinPvs) \|\|
755	(data.mArea <= mTermMinArea) \|\|
756	(data.mDepth >= mTermMaxDepth) \|\|
757	(data.GetAvgRayContribution() < mTermMaxRayContribution));
758	}
759
760	BspNode *BspTree::Subdivide(BspTraversalStack &tStack, BspTraversalData &tData)
761	{
762	//-- terminate traversal
763	if (TerminationCriteriaMet(tData))
764	{
765	BspLeaf leaf = dynamic_cast<BspLeaf >(tData.mNode);
766
767	BspViewCell *viewCell;
768
769	// generate new view cell for each leaf
770	if (mGenerateViewCells)
771	{
772	viewCell = dynamic_cast<BspViewCell *>(ViewCell::Generate());
773	}
774	else
775	{
776	// add view cell to leaf
777	viewCell = dynamic_cast<BspViewCell *>(tData.mViewCell);
778	}
779
780	leaf->SetViewCell(viewCell);
781	viewCell->mBspLeaves.push_back(leaf);
782
783	//-- add pvs
784	if (viewCell != mRootCell)
785	{
786	int conSamp = 0, sampCon = 0;
787	leaf->AddToPvs(*tData.mRays, conSamp, sampCon);
788
789	mStat.contributingSamples += conSamp;
790	mStat.sampleContributions += sampCon;
791	}
792
793	EvaluateLeafStats(tData);
794
795	//-- clean up
796
797	// discard polygons
798	CLEAR_CONTAINER(*tData.mPolygons);
799	// discard rays
800	CLEAR_CONTAINER(*tData.mRays);
801
802	DEL_PTR(tData.mPolygons);
803	DEL_PTR(tData.mRays);
804	DEL_PTR(tData.mGeometry);
805
806	return leaf;
807	}
808
809	//-- continue subdivision
810	PolygonContainer coincident;
811
812	PolygonContainer *frontPolys = new PolygonContainer();
813	PolygonContainer *backPolys = new PolygonContainer();
814
815	BoundedRayContainer *frontRays = new BoundedRayContainer();
816	BoundedRayContainer *backRays = new BoundedRayContainer();
817
818	BspTraversalData tFrontData(NULL, new PolygonContainer(), tData.mDepth + 1, mRootCell,
819	new BoundedRayContainer(), 0, 0, new BspNodeGeometry());
820	BspTraversalData tBackData(NULL, new PolygonContainer(), tData.mDepth + 1, mRootCell,
821	new BoundedRayContainer(), 0, 0, new BspNodeGeometry());
822
823	// create new interior node and two leaf nodes
824	BspInterior *interior = SubdivideNode(tData,
825	tFrontData,
826	tBackData,
827	coincident);
828
829	#ifdef _DEBUG
830	if (frontPolys->empty() && backPolys->empty() && (coincident.size() > 2))
831	{ for (PolygonContainer::iterator it = coincident.begin(); it != coincident.end(); ++it)
832	Debug << (it) << " " << (it)->GetArea() << " " << (*it)->mParent << endl ;
833	Debug << endl;}
834	#endif
835
836	// extract view cells from coincident polygons according to plane normal
837	// only if front or back polygons are empty
838	if (!mGenerateViewCells)
839	{
840	ExtractViewCells(tFrontData,
841	tBackData,
842	coincident,
843	interior->mPlane);
844	}
845
846	// don't need coincident polygons anymory
847	CLEAR_CONTAINER(coincident);
848
849	// push the children on the stack
850	tStack.push(tFrontData);
851	tStack.push(tBackData);
852
853	// cleanup
854	DEL_PTR(tData.mNode);
855	DEL_PTR(tData.mPolygons);
856	DEL_PTR(tData.mRays);
857	DEL_PTR(tData.mGeometry);
858
859	return interior;
860	}
861
862	void BspTree::ExtractViewCells(BspTraversalData &frontData,
863	BspTraversalData &backData,
864	const PolygonContainer &coincident,
865	const Plane3 splitPlane) const
866	{
867	// if not empty, tree is further subdivided => don't have to find view cell
868	bool foundFront = !frontData.mPolygons->empty();
869	bool foundBack = !frontData.mPolygons->empty();
870
871	PolygonContainer::const_iterator it =
872	coincident.begin(), it_end = coincident.end();
873
874	//-- find first view cells in front and back leafs
875	for (; !(foundFront && foundBack) && (it != it_end); ++ it)
876	{
877	if (DotProd((*it)->GetNormal(), splitPlane.mNormal) > 0)
878	{
879	backData.mViewCell = dynamic_cast<ViewCell >((it)->mParent);
880	foundBack = true;
881	}
882	else
883	{
884	frontData.mViewCell = dynamic_cast<ViewCell >((it)->mParent);
885	foundFront = true;
886	}
887	}
888	}
889
890	BspInterior *BspTree::SubdivideNode(BspTraversalData &tData,
891	BspTraversalData &frontData,
892	BspTraversalData &backData,
893	PolygonContainer &coincident)
894	{
895	mStat.nodes += 2;
896
897	BspLeaf leaf = dynamic_cast<BspLeaf >(tData.mNode);
898	// select subdivision plane
899	BspInterior *interior =
900	new BspInterior(SelectPlane(leaf, tData));
901
902	#ifdef _DEBUG
903	Debug << interior << endl;
904	#endif
905
906	// subdivide rays into front and back rays
907	SplitRays(interior->mPlane, tData.mRays, frontData.mRays, *backData.mRays);
908
909	// subdivide polygons with plane
910	mStat.splits += interior->SplitPolygons(*tData.mPolygons,
911	*frontData.mPolygons,
912	*backData.mPolygons,
913	coincident);
914
915	// compute pvs
916	frontData.mPvs = ComputePvsSize(*frontData.mRays);
917	backData.mPvs = ComputePvsSize(*backData.mRays);
918
919	// split geometry and compute area
920	if (1)
921	{
922	tData.mGeometry->SplitGeometry(*frontData.mGeometry,
923	*backData.mGeometry,
924	*this,
925	interior->mPlane);
926
927
928	frontData.mArea = frontData.mGeometry->GetArea();
929	backData.mArea = backData.mGeometry->GetArea();
930	}
931
932	// compute accumulated ray length
933	//frontData.mAccRayLength = AccumulatedRayLength(*frontData.mRays);
934	//backData.mAccRayLength = AccumulatedRayLength(*backData.mRays);
935
936	//-- create front and back leaf
937
938	BspInterior *parent = leaf->GetParent();
939
940	// replace a link from node's parent
941	if (!leaf->IsRoot())
942	{
943	parent->ReplaceChildLink(leaf, interior);
944	interior->SetParent(parent);
945	}
946	else // new root
947	{
948	mRoot = interior;
949	}
950
951	// and setup child links
952	interior->SetupChildLinks(new BspLeaf(interior), new BspLeaf(interior));
953
954	frontData.mNode = interior->mFront;
955	backData.mNode = interior->mBack;
956
957	//DEL_PTR(leaf);
958	return interior;
959	}
960
961	void BspTree::SortSplitCandidates(const PolygonContainer &polys,
962	const int axis,
963	vector<SortableEntry> &splitCandidates) const
964	{
965	splitCandidates.clear();
966
967	int requestedSize = 2 * (int)polys.size();
968	// creates a sorted split candidates array
969	splitCandidates.reserve(requestedSize);
970
971	PolygonContainer::const_iterator it, it_end = polys.end();
972
973	AxisAlignedBox3 box;
974
975	// insert all queries
976	for(it = polys.begin(); it != it_end; ++ it)
977	{
978	box.Initialize();
979	box.Include((it));
980
981	splitCandidates.push_back(SortableEntry(SortableEntry::POLY_MIN, box.Min(axis), *it));
982	splitCandidates.push_back(SortableEntry(SortableEntry::POLY_MAX, box.Max(axis), *it));
983	}
984
985	stable_sort(splitCandidates.begin(), splitCandidates.end());
986	}
987
988
989	float BspTree::BestCostRatio(const PolygonContainer &polys,
990	const AxisAlignedBox3 &box,
991	const int axis,
992	float &position,
993	int &objectsBack,
994	int &objectsFront) const
995	{
996	vector<SortableEntry> splitCandidates;
997
998	SortSplitCandidates(polys, axis, splitCandidates);
999
1000	// go through the lists, count the number of objects left and right
1001	// and evaluate the following cost funcion:
1002	// C = ct_div_ci + (ol + or)/queries
1003
1004	int objectsLeft = 0, objectsRight = (int)polys.size();
1005
1006	float minBox = box.Min(axis);
1007	float maxBox = box.Max(axis);
1008	float boxArea = box.SurfaceArea();
1009
1010	float minBand = minBox + mSplitBorder * (maxBox - minBox);
1011	float maxBand = minBox + (1.0f - mSplitBorder) * (maxBox - minBox);
1012
1013	float minSum = 1e20f;
1014	vector<SortableEntry>::const_iterator ci, ci_end = splitCandidates.end();
1015
1016	for(ci = splitCandidates.begin(); ci != ci_end; ++ ci)
1017	{
1018	switch ((*ci).type)
1019	{
1020	case SortableEntry::POLY_MIN:
1021	++ objectsLeft;
1022	break;
1023	case SortableEntry::POLY_MAX:
1024	-- objectsRight;
1025	break;
1026	default:
1027	break;
1028	}
1029
1030	if ((ci).value > minBand && (ci).value < maxBand)
1031	{
1032	AxisAlignedBox3 lbox = box;
1033	AxisAlignedBox3 rbox = box;
1034	lbox.SetMax(axis, (*ci).value);
1035	rbox.SetMin(axis, (*ci).value);
1036
1037	float sum = objectsLeft * lbox.SurfaceArea() +
1038	objectsRight * rbox.SurfaceArea();
1039
1040	if (sum < minSum)
1041	{
1042	minSum = sum;
1043	position = (*ci).value;
1044
1045	objectsBack = objectsLeft;
1046	objectsFront = objectsRight;
1047	}
1048	}
1049	}
1050
1051	float oldCost = (float)polys.size();
1052	float newCost = mAaCtDivCi + minSum / boxArea;
1053	float ratio = newCost / oldCost;
1054
1055
1056	#if 0
1057	Debug << "====================" << endl;
1058	Debug << "costRatio=" << ratio << " pos=" << position<<" t=" << (position - minBox)/(maxBox - minBox)
1059	<< "\t o=(" << objectsBack << "," << objectsFront << ")" << endl;
1060	#endif
1061	return ratio;
1062	}
1063
1064	bool BspTree::SelectAxisAlignedPlane(Plane3 &plane,
1065	const PolygonContainer &polys) const
1066	{
1067	AxisAlignedBox3 box;
1068	box.Initialize();
1069
1070	// create bounding box of region
1071	Polygon3::IncludeInBox(polys, box);
1072
1073	int objectsBack = 0, objectsFront = 0;
1074	int axis = 0;
1075	float costRatio = MAX_FLOAT;
1076	Vector3 position;
1077
1078	//-- area subdivision
1079	for (int i = 0; i < 3; ++ i)
1080	{
1081	float p = 0;
1082	float r = BestCostRatio(polys, box, i, p, objectsBack, objectsFront);
1083
1084	if (r < costRatio)
1085	{
1086	costRatio = r;
1087	axis = i;
1088	position = p;
1089	}
1090	}
1091
1092	if (costRatio >= mMaxCostRatio)
1093	return false;
1094
1095	Vector3 norm(0,0,0); norm[axis] = 1.0f;
1096	plane = Plane3(norm, position);
1097
1098	return true;
1099	}
1100
1101	Plane3 BspTree::SelectPlane(BspLeaf *leaf, BspTraversalData &data)
1102	{
1103	if (data.mPolygons->empty() && data.mRays->empty())
1104	{
1105	Debug << "Warning: No autopartition polygon candidate available\n";
1106
1107	// return axis aligned split
1108	AxisAlignedBox3 box;
1109	box.Initialize();
1110
1111	// create bounding box of region
1112	Polygon3::IncludeInBox(*data.mPolygons, box);
1113
1114	const int axis = box.Size().DrivingAxis();
1115	const Vector3 position = (box.Min()[axis] + box.Max()[axis])*0.5f;
1116
1117	Vector3 norm(0,0,0); norm[axis] = 1.0f;
1118	return Plane3(norm, position);
1119	}
1120
1121	if ((mSplitPlaneStrategy & AXIS_ALIGNED) &&
1122	((int)data.mPolygons->size() > mTermMinPolysForAxisAligned) &&
1123	((int)data.mRays->size() > mTermMinRaysForAxisAligned) &&
1124	((mTermMinObjectsForAxisAligned < 0) \|\|
1125	(Polygon3::ParentObjectsSize(*data.mPolygons) > mTermMinObjectsForAxisAligned)))
1126	{
1127	Plane3 plane;
1128	if (SelectAxisAlignedPlane(plane, *data.mPolygons))
1129	return plane;
1130	}
1131
1132	// simplest strategy: just take next polygon
1133	if (mSplitPlaneStrategy & RANDOM_POLYGON)
1134	{
1135	if (!data.mPolygons->empty())
1136	{
1137	Polygon3 nextPoly = (data.mPolygons)[Random((int)data.mPolygons->size())];
1138	return nextPoly->GetSupportingPlane();
1139	}
1140	else
1141	{
1142	const int candidateIdx = Random((int)data.mRays->size());
1143	BoundedRay bRay = (data.mRays)[candidateIdx];
1144
1145	Ray *ray = bRay->mRay;
1146
1147	const Vector3 minPt = ray->Extrap(bRay->mMinT);
1148	const Vector3 maxPt = ray->Extrap(bRay->mMaxT);
1149
1150	const Vector3 pt = (maxPt + minPt) * 0.5;
1151
1152	const Vector3 normal = ray->GetDir();
1153
1154	return Plane3(normal, pt);
1155	}
1156
1157	return Plane3();
1158	}
1159
1160	// use heuristics to find appropriate plane
1161	return SelectPlaneHeuristics(leaf, data);
1162	}
1163
1164	Plane3 BspTree::SelectPlaneHeuristics(BspLeaf *leaf, BspTraversalData &data)
1165	{
1166	float lowestCost = MAX_FLOAT;
1167	Plane3 bestPlane;
1168	Plane3 plane;
1169
1170	int limit = Min((int)data.mPolygons->size(), mMaxPolyCandidates);
1171
1172	int candidateIdx = limit;
1173
1174	for (int i = 0; i < limit; ++ i)
1175	{
1176	candidateIdx = GetNextCandidateIdx(candidateIdx, *data.mPolygons);
1177
1178	Polygon3 poly = (data.mPolygons)[candidateIdx];
1179	// evaluate current candidate
1180	const float candidateCost =
1181	SplitPlaneCost(poly->GetSupportingPlane(), data);
1182
1183	if (candidateCost < lowestCost)
1184	{
1185	bestPlane = poly->GetSupportingPlane();
1186	lowestCost = candidateCost;
1187	}
1188	}
1189
1190	//Debug << "lowest: " << lowestCost << endl;
1191
1192	//-- choose candidate planes extracted from rays
1193	// we currently use two methods
1194	// 1) take 3 ray endpoints, where two are minimum and one a maximum
1195	// point or the other way round
1196	// 2) take plane normal as plane normal and the midpoint of the ray.
1197	// PROBLEM: does not resemble any point where visibility is likely to change
1198	const BoundedRayContainer *rays = data.mRays;
1199
1200	for (int i = 0; i < mMaxRayCandidates / 2; ++ i)
1201	{
1202	candidateIdx = Random((int)rays->size());
1203	BoundedRay bRay = (rays)[candidateIdx];
1204
1205	Ray *ray = bRay->mRay;
1206
1207	const Vector3 minPt = ray->Extrap(bRay->mMinT);
1208	const Vector3 maxPt = ray->Extrap(bRay->mMaxT);
1209
1210	const Vector3 pt = (maxPt + minPt) * 0.5;
1211
1212	const Vector3 normal = ray->GetDir();
1213
1214	plane = Plane3(normal, pt);
1215
1216	const float candidateCost = SplitPlaneCost(plane, data);
1217
1218	if (candidateCost < lowestCost)
1219	{
1220	bestPlane = plane;
1221
1222	lowestCost = candidateCost;
1223	}
1224	}
1225
1226	//Debug << "lowest: " << lowestCost << endl;
1227	for (int i = 0; i < mMaxRayCandidates / 2; ++ i)
1228	{
1229	Vector3 pt[3];
1230	int idx[3];
1231	int cmaxT = 0;
1232	int cminT = 0;
1233	bool chooseMin = false;
1234
1235	for (int j = 0; j < 3; j ++)
1236	{
1237	idx[j] = Random((int)rays->size() * 2);
1238
1239	if (idx[j] >= (int)rays->size())
1240	{
1241	idx[j] -= (int)rays->size();
1242
1243	chooseMin = (cminT < 2);
1244	}
1245	else
1246	chooseMin = (cmaxT < 2);
1247
1248	BoundedRay bRay = (rays)[idx[j]];
1249	pt[j] = chooseMin ? bRay->mRay->Extrap(bRay->mMinT) : bRay->mRay->Extrap(bRay->mMaxT);
1250	}
1251
1252	plane = Plane3(pt[0], pt[1], pt[2]);
1253
1254	const float candidateCost = SplitPlaneCost(plane, data);
1255
1256	if (candidateCost < lowestCost)
1257	{
1258	//Debug << "choose ray plane 2: " << candidateCost << endl;
1259	bestPlane = plane;
1260
1261	lowestCost = candidateCost;
1262	}
1263	}
1264
1265	#ifdef _DEBUG
1266	Debug << "plane lowest cost: " << lowestCost << endl;
1267	#endif
1268	return bestPlane;
1269	}
1270
1271	int BspTree::GetNextCandidateIdx(int currentIdx, PolygonContainer &polys)
1272	{
1273	const int candidateIdx = Random(currentIdx --);
1274
1275	// swap candidates to avoid testing same plane 2 times
1276	std::swap(polys[currentIdx], polys[candidateIdx]);
1277
1278	return currentIdx;
1279	//return Random((int)polys.size());
1280	}
1281
1282	float BspTree::SplitPlaneCost(const Plane3 &candidatePlane,
1283	const PolygonContainer &polys) const
1284	{
1285	float val = 0;
1286
1287	float sumBalancedPolys = 0;
1288	float sumSplits = 0;
1289	float sumPolyArea = 0;
1290	float sumBalancedViewCells = 0;
1291	float sumBlockedRays = 0;
1292	float totalBlockedRays = 0;
1293	//float totalArea = 0;
1294	int totalViewCells = 0;
1295
1296	// need three unique ids for each type of view cell
1297	// for balanced view cells criterium
1298	ViewCell::NewMail();
1299	const int backId = ViewCell::sMailId;
1300	ViewCell::NewMail();
1301	const int frontId = ViewCell::sMailId;
1302	ViewCell::NewMail();
1303	const int frontAndBackId = ViewCell::sMailId;
1304
1305	PolygonContainer::const_iterator it, it_end = polys.end();
1306
1307	for (it = polys.begin(); it != it_end; ++ it)
1308	{
1309	const int classification = (*it)->ClassifyPlane(candidatePlane);
1310
1311	if (mSplitPlaneStrategy & BALANCED_POLYS)
1312	sumBalancedPolys += sBalancedPolysTable[classification];
1313
1314	if (mSplitPlaneStrategy & LEAST_SPLITS)
1315	sumSplits += sLeastPolySplitsTable[classification];
1316
1317	if (mSplitPlaneStrategy & LARGEST_POLY_AREA)
1318	{
1319	if (classification == Polygon3::COINCIDENT)
1320	sumPolyArea += (*it)->GetArea();
1321	//totalArea += area;
1322	}
1323
1324	if (mSplitPlaneStrategy & BLOCKED_RAYS)
1325	{
1326	const float blockedRays = (float)(*it)->mPiercingRays.size();
1327
1328	if (classification == Polygon3::COINCIDENT)
1329	sumBlockedRays += blockedRays;
1330
1331	totalBlockedRays += blockedRays;
1332	}
1333
1334	// assign view cells to back or front according to classificaion
1335	if (mSplitPlaneStrategy & BALANCED_VIEW_CELLS)
1336	{
1337	MeshInstance viewCell = (it)->mParent;
1338
1339	// assure that we only count a view cell
1340	// once for the front and once for the back side of the plane
1341	if (classification == Polygon3::FRONT_SIDE)
1342	{
1343	if ((viewCell->mMailbox != frontId) &&
1344	(viewCell->mMailbox != frontAndBackId))
1345	{
1346	sumBalancedViewCells += 1.0;
1347
1348	if (viewCell->mMailbox != backId)
1349	viewCell->mMailbox = frontId;
1350	else
1351	viewCell->mMailbox = frontAndBackId;
1352
1353	++ totalViewCells;
1354	}
1355	}
1356	else if (classification == Polygon3::BACK_SIDE)
1357	{
1358	if ((viewCell->mMailbox != backId) &&
1359	(viewCell->mMailbox != frontAndBackId))
1360	{
1361	sumBalancedViewCells -= 1.0;
1362
1363	if (viewCell->mMailbox != frontId)
1364	viewCell->mMailbox = backId;
1365	else
1366	viewCell->mMailbox = frontAndBackId;
1367
1368	++ totalViewCells;
1369	}
1370	}
1371	}
1372	}
1373
1374	const float polysSize = (float)polys.size() + Limits::Small;
1375
1376	// all values should be approx. between 0 and 1 so they can be combined
1377	// and scaled with the factors according to their importance
1378	if (mSplitPlaneStrategy & BALANCED_POLYS)
1379	val += mBalancedPolysFactor * fabs(sumBalancedPolys) / polysSize;
1380
1381	if (mSplitPlaneStrategy & LEAST_SPLITS)
1382	val += mLeastSplitsFactor * sumSplits / polysSize;
1383
1384	if (mSplitPlaneStrategy & LARGEST_POLY_AREA)
1385	// HACK: polys.size should be total area so scaling is between 0 and 1
1386	val += mLargestPolyAreaFactor * (float)polys.size() / sumPolyArea;
1387
1388	if (mSplitPlaneStrategy & BLOCKED_RAYS)
1389	if (totalBlockedRays != 0)
1390	val += mBlockedRaysFactor * (totalBlockedRays - sumBlockedRays) / totalBlockedRays;
1391
1392	if (mSplitPlaneStrategy & BALANCED_VIEW_CELLS)
1393	val += mBalancedViewCellsFactor * fabs(sumBalancedViewCells) /
1394	((float)totalViewCells + Limits::Small);
1395
1396	return val;
1397	}
1398
1399	bool BspTree::BoundRay(const Ray &ray, float &minT, float &maxT) const
1400	{
1401	maxT = 1e6;
1402	minT = 0;
1403
1404	// test with tree bounding box
1405	if (!mBox.GetMinMaxT(ray, &minT, &maxT))
1406	return false;
1407
1408	if (minT < 0) // start ray from origin
1409	minT = 0;
1410
1411	// bound ray or line segment
1412	if ((ray.GetType() == Ray::LOCAL_RAY) &&
1413	!ray.intersections.empty() &&
1414	(ray.intersections[0].mT <= maxT))
1415	{
1416	maxT = ray.intersections[0].mT;
1417	}
1418
1419	return true;
1420	}
1421
1422	float BspTree::SplitPlaneCost(const Plane3 &candidatePlane,
1423	const BoundedRayContainer &rays,
1424	const int pvs,
1425	const float area,
1426	const BspNodeGeometry &cell) const
1427	{
1428	float val = 0;
1429
1430	float sumBalancedRays = 0;
1431	float sumRaySplits = 0;
1432
1433	int backId = 0;
1434	int frontId = 0;
1435	int frontAndBackId = 0;
1436
1437	int frontPvs = 0;
1438	int backPvs = 0;
1439
1440	// probability that view point lies in child
1441	float pOverall = 0;
1442	float pFront = 0;
1443	float pBack = 0;
1444
1445	if (mSplitPlaneStrategy & PVS)
1446	{
1447	// create three unique ids for pvs heuristics
1448	Intersectable::NewMail(); backId = ViewCell::sMailId;
1449	Intersectable::NewMail(); frontId = ViewCell::sMailId;
1450	Intersectable::NewMail(); frontAndBackId = ViewCell::sMailId;
1451
1452	if (mPvsUseArea) // use front and back cell areas to approximate volume
1453	{
1454	// construct child geometry with regard to the candidate split plane
1455	BspNodeGeometry frontCell;
1456	BspNodeGeometry backCell;
1457
1458	cell.SplitGeometry(frontCell, backCell, *this, candidatePlane);
1459
1460	pFront = frontCell.GetArea();
1461	pBack = backCell.GetArea();
1462
1463	pOverall = area;
1464	}
1465	}
1466
1467	BoundedRayContainer::const_iterator rit, rit_end = rays.end();
1468
1469	for (rit = rays.begin(); rit != rays.end(); ++ rit)
1470	{
1471	Ray ray = (rit)->mRay;
1472	const float minT = (*rit)->mMinT;
1473	const float maxT = (*rit)->mMaxT;
1474
1475	Vector3 entP, extP;
1476
1477	const int cf =
1478	ray->ClassifyPlane(candidatePlane, minT, maxT, entP, extP);
1479
1480	if (mSplitPlaneStrategy & LEAST_RAY_SPLITS)
1481	{
1482	sumBalancedRays += sBalancedRaysTable[cf];
1483	}
1484
1485	if (mSplitPlaneStrategy & BALANCED_RAYS)
1486	{
1487	sumRaySplits += sLeastRaySplitsTable[cf];
1488	}
1489
1490	if (mSplitPlaneStrategy & PVS)
1491	{
1492	if (!ray->intersections.empty())
1493	{
1494	// in case the ray intersects an objcrs
1495	// assure that we only count a object
1496	// once for the front and once for the back side of the plane
1497	IncPvs(*ray->intersections[0].mObject, frontPvs, backPvs,
1498	cf, frontId, backId, frontAndBackId);
1499	}
1500
1501	// the source object in the origin of the ray
1502	if (ray->sourceObject.mObject)
1503	{
1504	IncPvs(*ray->sourceObject.mObject, frontPvs, backPvs,
1505	cf, frontId, backId, frontAndBackId);
1506	}
1507
1508	if (!mPvsUseArea) // use front and back cell areas to approximate volume
1509	{
1510	float len = Distance(entP, extP);
1511	pOverall += len;
1512
1513	// use length of rays to approximate volume
1514	switch (cf)
1515	{
1516	case Ray::COINCIDENT:
1517	pBack += len;
1518	pFront += len;
1519	break;
1520	case Ray::BACK:
1521	pBack += len;
1522	break;
1523	case Ray::FRONT:
1524	pFront += len;
1525	break;
1526	case Ray::FRONT_BACK:
1527	{
1528	// find intersection of ray segment with plane
1529	const Vector3 extp = ray->Extrap(maxT);
1530	const float t = candidatePlane.FindT(ray->GetLoc(), extp);
1531
1532	const float newT = t * maxT;
1533	float newLen = Distance(ray->Extrap(newT), extp);
1534
1535	pFront += len - newLen;
1536	pBack += newLen;
1537	}
1538	break;
1539	case Ray::BACK_FRONT:
1540	{
1541	// find intersection of ray segment with plane
1542	const Vector3 extp = ray->Extrap(maxT);
1543	const float t = candidatePlane.FindT(ray->GetLoc(), extp);
1544
1545	const float newT = t * maxT;
1546	float newLen = Distance(ray->Extrap(newT), extp);
1547
1548	pFront += len;
1549	pBack += len - newLen;
1550	}
1551	break;
1552	default:
1553	Debug << "Should not come here" << endl;
1554	break;
1555	}
1556	}
1557	}
1558	}
1559
1560	const float raysSize = (float)rays.size() + Limits::Small;
1561	if (mSplitPlaneStrategy & LEAST_RAY_SPLITS)
1562	val += mLeastRaySplitsFactor * sumRaySplits / raysSize;
1563
1564	if (mSplitPlaneStrategy & BALANCED_RAYS)
1565	val += mBalancedRaysFactor * fabs(sumBalancedRays) / raysSize;
1566
1567	float denom = pOverall * (float)pvs * 2.0f + Limits::Small;
1568	if ((mSplitPlaneStrategy & PVS) && area && pvs)
1569	{
1570	val += mPvsFactor * (frontPvs * pFront + (backPvs * pBack)) / denom;
1571
1572	// give penalty to unbalanced split
1573	if (0)
1574	if (((pFront * 0.2 + Limits::Small) > pBack) \|\| (pFront < (pBack * 0.2 + Limits::Small)))
1575	val += 0.5;
1576	}
1577
1578	#ifdef _DEBUG
1579	Debug << "totalpvs: " << pvs << " ptotal: " << pOverall
1580	<< " frontpvs: " << frontPvs << " pFront: " << pFront
1581	<< " backpvs: " << backPvs << " pBack: " << pBack << endl << endl;
1582	#endif
1583	return val;
1584	}
1585
1586	void BspTree::IncPvs(Intersectable &obj,
1587	int &frontPvs,
1588	int &backPvs,
1589	const int cf,
1590	const int frontId,
1591	const int backId,
1592	const int frontAndBackId) const
1593	{
1594	// TODO: does this really belong to no pvs?
1595	//if (cf == Ray::COINCIDENT) return;
1596
1597	if (cf == Ray::FRONT)
1598	{
1599	if ((obj.mMailbox != frontId) &&
1600	(obj.mMailbox != frontAndBackId))
1601	{
1602	++ frontPvs;
1603
1604	if (obj.mMailbox != backId)
1605	obj.mMailbox = frontId;
1606	else
1607	obj.mMailbox = frontAndBackId;
1608	}
1609	}
1610	else if (cf == Ray::BACK)
1611	{
1612	if ((obj.mMailbox != backId) &&
1613	(obj.mMailbox != frontAndBackId))
1614	{
1615	++ backPvs;
1616
1617	if (obj.mMailbox != frontId)
1618	obj.mMailbox = backId;
1619	else
1620	obj.mMailbox = frontAndBackId;
1621	}
1622	}
1623	// object belongs to both PVS
1624	else if ((cf == Ray::FRONT_BACK) \|\| (cf == Ray::BACK_FRONT) \|\|(cf == Ray::COINCIDENT))
1625	{
1626	if (obj.mMailbox != frontAndBackId)
1627	{
1628	if (obj.mMailbox != frontId)
1629	++ frontPvs;
1630	if (obj.mMailbox != backId)
1631	++ backPvs;
1632
1633	obj.mMailbox = frontAndBackId;
1634	}
1635	}
1636	}
1637
1638	float BspTree::SplitPlaneCost(const Plane3 &candidatePlane,
1639	BspTraversalData &data) const
1640	{
1641	float val = 0;
1642
1643	if (mSplitPlaneStrategy & VERTICAL_AXIS)
1644	{
1645	Vector3 tinyAxis(0,0,0); tinyAxis[mBox.Size().TinyAxis()] = 1.0f;
1646	// we put a penalty on the dot product between the "tiny" vertical axis
1647	// and the split plane axis
1648	val += mVerticalSplitsFactor *
1649	fabs(DotProd(candidatePlane.mNormal, tinyAxis));
1650	}
1651
1652	// the following criteria loop over all polygons to find the cost value
1653	if ((mSplitPlaneStrategy & BALANCED_POLYS) \|\|
1654	(mSplitPlaneStrategy & LEAST_SPLITS) \|\|
1655	(mSplitPlaneStrategy & LARGEST_POLY_AREA) \|\|
1656	(mSplitPlaneStrategy & BALANCED_VIEW_CELLS) \|\|
1657	(mSplitPlaneStrategy & BLOCKED_RAYS))
1658	{
1659	val += SplitPlaneCost(candidatePlane, *data.mPolygons);
1660	}
1661
1662	// the following criteria loop over all rays to find the cost value
1663	if ((mSplitPlaneStrategy & BALANCED_RAYS) \|\|
1664	(mSplitPlaneStrategy & LEAST_RAY_SPLITS) \|\|
1665	(mSplitPlaneStrategy & PVS))
1666	{
1667	val += SplitPlaneCost(candidatePlane, *data.mRays, data.mPvs,
1668	data.mArea, *data.mGeometry);
1669	}
1670
1671	// return linear combination of the sums
1672	return val;
1673	}
1674
1675	void BspTree::CollectLeaves(vector<BspLeaf *> &leaves) const
1676	{
1677	stack<BspNode *> nodeStack;
1678	nodeStack.push(mRoot);
1679
1680	while (!nodeStack.empty())
1681	{
1682	BspNode *node = nodeStack.top();
1683
1684	nodeStack.pop();
1685
1686	if (node->IsLeaf())
1687	{
1688	BspLeaf leaf = (BspLeaf )node;
1689	leaves.push_back(leaf);
1690	}
1691	else
1692	{
1693	BspInterior interior = dynamic_cast<BspInterior >(node);
1694
1695	nodeStack.push(interior->GetBack());
1696	nodeStack.push(interior->GetFront());
1697	}
1698	}
1699	}
1700
1701	AxisAlignedBox3 BspTree::GetBoundingBox() const
1702	{
1703	return mBox;
1704	}
1705
1706	BspNode *BspTree::GetRoot() const
1707	{
1708	return mRoot;
1709	}
1710
1711	void BspTree::EvaluateLeafStats(const BspTraversalData &data)
1712	{
1713	// the node became a leaf -> evaluate stats for leafs
1714	BspLeaf leaf = dynamic_cast<BspLeaf >(data.mNode);
1715
1716	// store maximal and minimal depth
1717	if (data.mDepth > mStat.maxDepth)
1718	mStat.maxDepth = data.mDepth;
1719
1720	if (data.mDepth < mStat.minDepth)
1721	mStat.minDepth = data.mDepth;
1722
1723	// accumulate depth to compute average depth
1724	mStat.accumDepth += data.mDepth;
1725
1726
1727	if (data.mDepth >= mTermMaxDepth)
1728	++ mStat.maxDepthNodes;
1729
1730	if (data.mPvs < mTermMinPvs)
1731	++ mStat.minPvsNodes;
1732
1733	if ((int)data.mRays->size() < mTermMinRays)
1734	++ mStat.minRaysNodes;
1735
1736	if (data.GetAvgRayContribution() > mTermMaxRayContribution)
1737	++ mStat.maxRayContribNodes;
1738
1739	if (data.mGeometry->GetArea() <= mTermMinArea)
1740	++ mStat.minAreaNodes;
1741
1742	#ifdef _DEBUG
1743	Debug << "BSP stats: "
1744	<< "Depth: " << data.mDepth << " (max: " << mTermMaxDepth << "), "
1745	<< "PVS: " << data.mPvs << " (min: " << mTermMinPvs << "), "
1746	<< "Area: " << data.mArea << " (min: " << mTermMinArea << "), "
1747	<< "#polygons: " << (int)data.mPolygons->size() << " (max: " << mTermMinPolys << "), "
1748	<< "#rays: " << (int)data.mRays->size() << " (max: " << mTermMinRays << "), "
1749	<< "#pvs: " << leaf->GetViewCell()->GetPvs().GetSize() << "=, "
1750	<< "#avg ray contrib (pvs): " << (float)data.mPvs / (float)data.mRays->size() << endl;
1751	#endif
1752	}
1753
1754	int BspTree::CastRay(Ray &ray)
1755	{
1756	int hits = 0;
1757
1758	stack<BspRayTraversalData> tStack;
1759
1760	float maxt, mint;
1761
1762	if (!BoundRay(ray, mint, maxt))
1763	return 0;
1764
1765	Intersectable::NewMail();
1766
1767	Vector3 entp = ray.Extrap(mint);
1768	Vector3 extp = ray.Extrap(maxt);
1769
1770	BspNode *node = mRoot;
1771	BspNode *farChild = NULL;
1772
1773	while (1)
1774	{
1775	if (!node->IsLeaf())
1776	{
1777	BspInterior in = (BspInterior ) node;
1778
1779	Plane3 *splitPlane = in->GetPlane();
1780
1781	int entSide = splitPlane->Side(entp);
1782	int extSide = splitPlane->Side(extp);
1783
1784	Vector3 intersection;
1785
1786	if (entSide < 0)
1787	{
1788	node = in->GetBack();
1789
1790	if(extSide <= 0) // plane does not split ray => no far child
1791	continue;
1792
1793	farChild = in->GetFront(); // plane splits ray
1794
1795	} else if (entSide > 0)
1796	{
1797	node = in->GetFront();
1798
1799	if (extSide >= 0) // plane does not split ray => no far child
1800	continue;
1801
1802	farChild = in->GetBack(); // plane splits ray
1803	}
1804	else // ray and plane are coincident
1805	{
1806	// WHAT TO DO IN THIS CASE ?
1807	//break;
1808	node = in->GetFront();
1809	continue;
1810	}
1811
1812	// push data for far child
1813	tStack.push(BspRayTraversalData(farChild, extp, maxt));
1814
1815	// find intersection of ray segment with plane
1816	float t;
1817	extp = splitPlane->FindIntersection(ray.GetLoc(), extp, &t);
1818	maxt *= t;
1819
1820	} else // reached leaf => intersection with view cell
1821	{
1822	BspLeaf leaf = dynamic_cast<BspLeaf >(node);
1823
1824	if (!leaf->mViewCell->Mailed())
1825	{
1826	ray.bspIntersections.push_back(Ray::BspIntersection(maxt, leaf));
1827	leaf->mViewCell->Mail();
1828	++ hits;
1829	}
1830
1831	//-- fetch the next far child from the stack
1832	if (tStack.empty())
1833	break;
1834
1835	entp = extp;
1836	mint = maxt; // NOTE: need this?
1837
1838	if (ray.GetType() == Ray::LINE_SEGMENT && mint > 1.0f)
1839	break;
1840
1841	BspRayTraversalData &s = tStack.top();
1842
1843	node = s.mNode;
1844	extp = s.mExitPoint;
1845	maxt = s.mMaxT;
1846
1847	tStack.pop();
1848	}
1849	}
1850
1851	return hits;
1852	}
1853
1854	bool BspTree::Export(const string filename)
1855	{
1856	Exporter *exporter = Exporter::GetExporter(filename);
1857
1858	if (exporter)
1859	{
1860	exporter->ExportBspTree(*this);
1861	return true;
1862	}
1863
1864	return false;
1865	}
1866
1867	void BspTree::CollectViewCells(ViewCellContainer &viewCells) const
1868	{
1869	stack<BspNode *> nodeStack;
1870	nodeStack.push(mRoot);
1871
1872	ViewCell::NewMail();
1873
1874	while (!nodeStack.empty())
1875	{
1876	BspNode *node = nodeStack.top();
1877	nodeStack.pop();
1878
1879	if (node->IsLeaf())
1880	{
1881	ViewCell viewCell = dynamic_cast<BspLeaf >(node)->mViewCell;
1882
1883	if (!viewCell->Mailed())
1884	{
1885	viewCell->Mail();
1886	viewCells.push_back(viewCell);
1887	}
1888	}
1889	else
1890	{
1891	BspInterior interior = dynamic_cast<BspInterior >(node);
1892
1893	nodeStack.push(interior->mFront);
1894	nodeStack.push(interior->mBack);
1895	}
1896	}
1897	}
1898
1899	void BspTree::EvaluateViewCellsStats(BspViewCellsStatistics &stat) const
1900	{
1901	stat.Reset();
1902
1903	stack<BspNode *> nodeStack;
1904	nodeStack.push(mRoot);
1905
1906	ViewCell::NewMail();
1907
1908	// exclude root cell
1909	mRootCell->Mail();
1910
1911	while (!nodeStack.empty())
1912	{
1913	BspNode *node = nodeStack.top();
1914	nodeStack.pop();
1915
1916	if (node->IsLeaf())
1917	{
1918	++ stat.bspLeaves;
1919
1920	BspViewCell viewCell = dynamic_cast<BspLeaf >(node)->mViewCell;
1921
1922	if (!viewCell->Mailed())
1923	{
1924	viewCell->Mail();
1925
1926	++ stat.viewCells;
1927	const int pvsSize = viewCell->GetPvs().GetSize();
1928
1929	stat.pvs += pvsSize;
1930
1931	if (pvsSize < 1)
1932	++ stat.emptyPvs;
1933
1934	if (pvsSize > stat.maxPvs)
1935	stat.maxPvs = pvsSize;
1936
1937	if (pvsSize < stat.minPvs)
1938	stat.minPvs = pvsSize;
1939
1940	if ((int)viewCell->mBspLeaves.size() > stat.maxBspLeaves)
1941	stat.maxBspLeaves = (int)viewCell->mBspLeaves.size();
1942	}
1943	}
1944	else
1945	{
1946	BspInterior interior = dynamic_cast<BspInterior >(node);
1947
1948	nodeStack.push(interior->mFront);
1949	nodeStack.push(interior->mBack);
1950	}
1951	}
1952	}
1953
1954	bool BspTree::MergeViewCells(BspLeaf front, BspLeaf back) const
1955	{
1956	BspViewCell *viewCell =
1957	dynamic_cast<BspViewCell >(ViewCell::Merge(front->mViewCell, *back->mViewCell));
1958
1959	if (!viewCell)
1960	return false;
1961
1962	// change view cells of all leaves associated with the
1963	// previous view cells
1964
1965	BspViewCell *fVc = front->mViewCell;
1966	BspViewCell *bVc = back->mViewCell;
1967
1968	vector<BspLeaf *> fLeaves = fVc->mBspLeaves;
1969	vector<BspLeaf *> bLeaves = bVc->mBspLeaves;
1970
1971	vector<BspLeaf *>::const_iterator it;
1972
1973	for (it = fLeaves.begin(); it != fLeaves.end(); ++ it)
1974	{
1975	(*it)->SetViewCell(viewCell);
1976	viewCell->mBspLeaves.push_back(*it);
1977	}
1978	for (it = bLeaves.begin(); it != bLeaves.end(); ++ it)
1979	{
1980	(*it)->SetViewCell(viewCell);
1981	viewCell->mBspLeaves.push_back(*it);
1982	}
1983
1984	DEL_PTR(fVc);
1985	DEL_PTR(bVc);
1986
1987	return true;
1988	}
1989
1990	bool BspTree::ShouldMerge(BspLeaf front, BspLeaf back) const
1991	{
1992	ViewCell *fvc = front->mViewCell;
1993	ViewCell *bvc = back->mViewCell;
1994
1995	if ((fvc == mRootCell) \|\| (bvc == mRootCell) \|\| (fvc == bvc))
1996	return false;
1997
1998	int fdiff = fvc->GetPvs().Diff(bvc->GetPvs());
1999
2000	if (fvc->GetPvs().GetSize() + fdiff < mMaxPvs)
2001	{
2002	if ((fvc->GetPvs().GetSize() < mMinPvs) \|\|
2003	(bvc->GetPvs().GetSize() < mMinPvs) \|\|
2004	((fdiff < mMinPvsDif) && (bvc->GetPvs().Diff(fvc->GetPvs()) < mMinPvsDif)))
2005	{
2006	return true;
2007	}
2008	}
2009
2010	return false;
2011	}
2012
2013	void BspTree::SetGenerateViewCells(int generateViewCells)
2014	{
2015	mGenerateViewCells = generateViewCells;
2016	}
2017
2018	BspTreeStatistics &BspTree::GetStat()
2019	{
2020	return mStat;
2021	}
2022
2023	float BspTree::AccumulatedRayLength(BoundedRayContainer &rays) const
2024	{
2025	float len = 0;
2026
2027	BoundedRayContainer::const_iterator it, it_end = rays.end();
2028
2029	for (it = rays.begin(); it != it_end; ++ it)
2030	{
2031	len += SqrDistance((it)->mRay->Extrap((it)->mMinT),
2032	(it)->mRay->Extrap((it)->mMaxT));
2033	}
2034
2035	return len;
2036	}
2037
2038	int BspTree::SplitRays(const Plane3 &plane,
2039	BoundedRayContainer &rays,
2040	BoundedRayContainer &frontRays,
2041	BoundedRayContainer &backRays)
2042	{
2043	int splits = 0;
2044
2045	while (!rays.empty())
2046	{
2047	BoundedRay *bRay = rays.back();
2048	Ray *ray = bRay->mRay;
2049	float minT = bRay->mMinT;
2050	float maxT = bRay->mMaxT;
2051
2052	rays.pop_back();
2053
2054	Vector3 entP, extP;
2055
2056	const int cf =
2057	ray->ClassifyPlane(plane, minT, maxT, entP, extP);
2058
2059	// set id to ray classification
2060	ray->SetId(cf);
2061
2062	switch (cf)
2063	{
2064	case Ray::COINCIDENT: // TODO: should really discard ray?
2065	frontRays.push_back(bRay);
2066	//DEL_PTR(bRay);
2067	break;
2068	case Ray::BACK:
2069	backRays.push_back(bRay);
2070	break;
2071	case Ray::FRONT:
2072	frontRays.push_back(bRay);
2073	break;
2074	case Ray::FRONT_BACK:
2075	{
2076	// find intersection of ray segment with plane
2077	const float t = plane.FindT(ray->GetLoc(), extP);
2078
2079	const float newT = t * maxT;
2080
2081	frontRays.push_back(new BoundedRay(ray, minT, newT));
2082	backRays.push_back(new BoundedRay(ray, newT, maxT));
2083
2084	DEL_PTR(bRay);
2085	}
2086	break;
2087	case Ray::BACK_FRONT:
2088	{
2089	// find intersection of ray segment with plane
2090	const float t = plane.FindT(ray->GetLoc(), extP);
2091	const float newT = t * bRay->mMaxT;
2092
2093	backRays.push_back(new BoundedRay(ray, bRay->mMinT, newT));
2094	frontRays.push_back(new BoundedRay(ray, newT, bRay->mMaxT));
2095	DEL_PTR(bRay);
2096
2097	++ splits;
2098	}
2099	break;
2100	default:
2101	Debug << "Should not come here" << endl;
2102	break;
2103	}
2104	}
2105
2106	return splits;
2107	}
2108
2109	void BspTree::ExtractHalfSpaces(BspNode *n, vector<Plane3> &halfSpaces) const
2110	{
2111	BspNode *lastNode;
2112	do
2113	{
2114	lastNode = n;
2115
2116	// want to get planes defining geometry of this node => don't take
2117	// split plane of node itself
2118	n = n->GetParent();
2119
2120	if (n)
2121	{
2122	BspInterior interior = dynamic_cast<BspInterior >(n);
2123	Plane3 halfSpace = dynamic_cast<BspInterior >(interior)->GetPlane();
2124
2125	if (interior->mFront != lastNode)
2126	halfSpace.ReverseOrientation();
2127
2128	halfSpaces.push_back(halfSpace);
2129	}
2130	}
2131	while (n);
2132	}
2133
2134	void BspTree::ConstructGeometry(BspNode *n, BspNodeGeometry &cell) const
2135	{
2136	PolygonContainer polys;
2137	ConstructGeometry(n, polys);
2138	cell.mPolys = polys;
2139	}
2140
2141	void BspTree::ConstructGeometry(BspViewCell *vc, PolygonContainer &cell) const
2142	{
2143	vector<BspLeaf *> leaves = vc->mBspLeaves;
2144
2145	vector<BspLeaf *>::const_iterator it, it_end = leaves.end();
2146
2147	for (it = leaves.begin(); it != it_end; ++ it)
2148	ConstructGeometry(*it, cell);
2149	}
2150
2151
2152	void BspTree::ConstructGeometry(BspNode *n, PolygonContainer &cell) const
2153	{
2154	vector<Plane3> halfSpaces;
2155	ExtractHalfSpaces(n, halfSpaces);
2156
2157	PolygonContainer candidatePolys;
2158
2159	// bounded planes are added to the polygons (reverse polygons
2160	// as they have to be outfacing
2161	for (int i = 0; i < (int)halfSpaces.size(); ++ i)
2162	{
2163	Polygon3 *p = GetBoundingBox().CrossSection(halfSpaces[i]);
2164
2165	if (p->Valid())
2166	{
2167	candidatePolys.push_back(p->CreateReversePolygon());
2168	DEL_PTR(p);
2169	}
2170	}
2171
2172	// add faces of bounding box (also could be faces of the cell)
2173	for (int i = 0; i < 6; ++ i)
2174	{
2175	VertexContainer vertices;
2176
2177	for (int j = 0; j < 4; ++ j)
2178	vertices.push_back(mBox.GetFace(i).mVertices[j]);
2179
2180	candidatePolys.push_back(new Polygon3(vertices));
2181	}
2182
2183	for (int i = 0; i < (int)candidatePolys.size(); ++ i)
2184	{
2185	// polygon is split by all other planes
2186	for (int j = 0; (j < (int)halfSpaces.size()) && candidatePolys[i]; ++ j)
2187	{
2188	if (i == j) // polygon and plane are coincident
2189	continue;
2190
2191	VertexContainer splitPts;
2192	Polygon3 frontPoly, backPoly;
2193
2194	const int cf = candidatePolys[i]->ClassifyPlane(halfSpaces[j]);
2195
2196	switch (cf)
2197	{
2198	case Polygon3::SPLIT:
2199	frontPoly = new Polygon3();
2200	backPoly = new Polygon3();
2201
2202	candidatePolys[i]->Split(halfSpaces[j], *frontPoly,
2203	*backPoly, splitPts);
2204
2205	DEL_PTR(candidatePolys[i]);
2206
2207	if (frontPoly->Valid())
2208	candidatePolys[i] = frontPoly;
2209	else
2210	DEL_PTR(frontPoly);
2211
2212	DEL_PTR(backPoly);
2213	break;
2214	case Polygon3::BACK_SIDE:
2215	DEL_PTR(candidatePolys[i]);
2216	break;
2217	// just take polygon as it is
2218	case Polygon3::FRONT_SIDE:
2219	case Polygon3::COINCIDENT:
2220	default:
2221	break;
2222	}
2223	}
2224
2225	if (candidatePolys[i])
2226	cell.push_back(candidatePolys[i]);
2227	}
2228	}
2229
2230
2231	int BspTree::FindNeighbors(BspNode n, vector<BspLeaf > &neighbors,
2232	const bool onlyUnmailed) const
2233	{
2234	PolygonContainer cell;
2235
2236	ConstructGeometry(n, cell);
2237
2238	stack<BspNode *> nodeStack;
2239	nodeStack.push(mRoot);
2240
2241	// planes needed to verify that we found neighbor leaf.
2242	vector<Plane3> halfSpaces;
2243	ExtractHalfSpaces(n, halfSpaces);
2244
2245	while (!nodeStack.empty())
2246	{
2247	BspNode *node = nodeStack.top();
2248	nodeStack.pop();
2249
2250	if (node->IsLeaf())
2251	{
2252	if (node != n && (!onlyUnmailed \|\| !node->Mailed()))
2253	{
2254	// test all planes of current node if candidate really
2255	// is neighbour
2256	PolygonContainer neighborCandidate;
2257	ConstructGeometry(node, neighborCandidate);
2258
2259	bool isAdjacent = true;
2260	for (int i = 0; (i < halfSpaces.size()) && isAdjacent; ++ i)
2261	{
2262	const int cf =
2263	Polygon3::ClassifyPlane(neighborCandidate, halfSpaces[i]);
2264
2265	if (cf == Polygon3::BACK_SIDE)
2266	isAdjacent = false;
2267	}
2268
2269	if (isAdjacent)
2270	neighbors.push_back(dynamic_cast<BspLeaf *>(node));
2271
2272	CLEAR_CONTAINER(neighborCandidate);
2273	}
2274	}
2275	else
2276	{
2277	BspInterior interior = dynamic_cast<BspInterior >(node);
2278
2279	const int cf = Polygon3::ClassifyPlane(cell, interior->mPlane);
2280
2281	if (cf == Polygon3::FRONT_SIDE)
2282	nodeStack.push(interior->mFront);
2283	else
2284	if (cf == Polygon3::BACK_SIDE)
2285	nodeStack.push(interior->mBack);
2286	else
2287	{
2288	// random decision
2289	nodeStack.push(interior->mBack);
2290	nodeStack.push(interior->mFront);
2291	}
2292	}
2293	}
2294
2295	CLEAR_CONTAINER(cell);
2296	return (int)neighbors.size();
2297	}
2298
2299	BspLeaf *BspTree::GetRandomLeaf(const Plane3 &halfspace)
2300	{
2301	stack<BspNode *> nodeStack;
2302	nodeStack.push(mRoot);
2303
2304	int mask = rand();
2305
2306	while (!nodeStack.empty())
2307	{
2308	BspNode *node = nodeStack.top();
2309	nodeStack.pop();
2310
2311	if (node->IsLeaf())
2312	{
2313	return dynamic_cast<BspLeaf *>(node);
2314	}
2315	else
2316	{
2317	BspInterior interior = dynamic_cast<BspInterior >(node);
2318
2319	BspNode *next;
2320
2321	PolygonContainer cell;
2322
2323	// todo: not very efficient: constructs full cell everytime
2324	ConstructGeometry(interior, cell);
2325
2326	const int cf = Polygon3::ClassifyPlane(cell, halfspace);
2327
2328	if (cf == Polygon3::BACK_SIDE)
2329	next = interior->mFront;
2330	else
2331	if (cf == Polygon3::FRONT_SIDE)
2332	next = interior->mFront;
2333	else
2334	{
2335	// random decision
2336	if (mask & 1)
2337	next = interior->mBack;
2338	else
2339	next = interior->mFront;
2340	mask = mask >> 1;
2341	}
2342
2343	nodeStack.push(next);
2344	}
2345	}
2346
2347	return NULL;
2348	}
2349
2350	BspLeaf *BspTree::GetRandomLeaf(const bool onlyUnmailed)
2351	{
2352	stack<BspNode *> nodeStack;
2353
2354	nodeStack.push(mRoot);
2355
2356	int mask = rand();
2357
2358	while (!nodeStack.empty())
2359	{
2360	BspNode *node = nodeStack.top();
2361	nodeStack.pop();
2362
2363	if (node->IsLeaf())
2364	{
2365	if ( (!onlyUnmailed \|\| !node->Mailed()) )
2366	return dynamic_cast<BspLeaf *>(node);
2367	}
2368	else
2369	{
2370	BspInterior interior = dynamic_cast<BspInterior >(node);
2371
2372	// random decision
2373	if (mask & 1)
2374	nodeStack.push(interior->mBack);
2375	else
2376	nodeStack.push(interior->mFront);
2377
2378	mask = mask >> 1;
2379	}
2380	}
2381
2382	return NULL;
2383	}
2384
2385	int BspTree::ComputePvsSize(const BoundedRayContainer &rays) const
2386	{
2387	int pvsSize = 0;
2388
2389	BoundedRayContainer::const_iterator rit, rit_end = rays.end();
2390
2391	Intersectable::NewMail();
2392
2393	for (rit = rays.begin(); rit != rays.end(); ++ rit)
2394	{
2395	Ray ray = (rit)->mRay;
2396
2397	if (!ray->intersections.empty())
2398	{
2399	if (!ray->intersections[0].mObject->Mailed())
2400	{
2401	ray->intersections[0].mObject->Mail();
2402	++ pvsSize;
2403	}
2404	}
2405	if (ray->sourceObject.mObject)
2406	{
2407	if (!ray->sourceObject.mObject->Mailed())
2408	{
2409	ray->sourceObject.mObject->Mail();
2410	++ pvsSize;
2411	}
2412	}
2413	}
2414
2415	return pvsSize;
2416	}
2417
2418	/*************************************************************
2419	* BspNodeGeometry Implementation *
2420	*************************************************************/
2421
2422	BspNodeGeometry::~BspNodeGeometry()
2423	{
2424	CLEAR_CONTAINER(mPolys);
2425	}
2426
2427	float BspNodeGeometry::GetArea() const
2428	{
2429	return Polygon3::GetArea(mPolys);
2430	}
2431
2432	void BspNodeGeometry::SplitGeometry(BspNodeGeometry &front,
2433	BspNodeGeometry &back,
2434	const BspTree &tree,
2435	const Plane3 &splitPlane) const
2436	{
2437	// get cross section of new polygon
2438	Polygon3 *planePoly = tree.GetBoundingBox().CrossSection(splitPlane);
2439
2440	planePoly = SplitPolygon(planePoly, tree);
2441
2442	//-- plane poly splits all other cell polygons
2443	for (int i = 0; i < (int)mPolys.size(); ++ i)
2444	{
2445	const int cf = mPolys[i]->ClassifyPlane(splitPlane, 0.00001f);
2446
2447	// split new polygon with all previous planes
2448	switch (cf)
2449	{
2450	case Polygon3::SPLIT:
2451	{
2452	Polygon3 *poly = new Polygon3(mPolys[i]->mVertices);
2453
2454	Polygon3 *frontPoly = new Polygon3();
2455	Polygon3 *backPoly = new Polygon3();
2456
2457	VertexContainer splitPts;
2458
2459	poly->Split(splitPlane, frontPoly, backPoly, splitPts);
2460
2461	DEL_PTR(poly);
2462
2463	if (frontPoly->Valid())
2464	front.mPolys.push_back(frontPoly);
2465	if (backPoly->Valid())
2466	back.mPolys.push_back(backPoly);
2467	}
2468
2469	break;
2470	case Polygon3::BACK_SIDE:
2471	back.mPolys.push_back(new Polygon3(mPolys[i]->mVertices));
2472	break;
2473	case Polygon3::FRONT_SIDE:
2474	front.mPolys.push_back(new Polygon3(mPolys[i]->mVertices));
2475	break;
2476	case Polygon3::COINCIDENT:
2477	//front.mPolys.push_back(CreateReversePolygon(mPolys[i]));
2478	back.mPolys.push_back(new Polygon3(mPolys[i]->mVertices));
2479	break;
2480	default:
2481	break;
2482	}
2483	}
2484
2485	//-- finally add the new polygon to the child cells
2486	if (planePoly)
2487	{
2488	// add polygon with normal pointing into positive half space to back cell
2489	back.mPolys.push_back(planePoly);
2490	// add polygon with reverse orientation to front cell
2491	front.mPolys.push_back(planePoly->CreateReversePolygon());
2492	}
2493
2494	//Debug << "returning new geometry " << mPolys.size() << " f: " << front.mPolys.size() << " b: " << back.mPolys.size() << endl;
2495	//Debug << "old area " << GetArea() << " f: " << front.GetArea() << " b: " << back.GetArea() << endl;
2496	}
2497
2498	Polygon3 BspNodeGeometry::SplitPolygon(Polygon3 planePoly,
2499	const BspTree &tree) const
2500	{
2501	// polygon is split by all other planes
2502	for (int i = 0; (i < (int)mPolys.size()) && planePoly; ++ i)
2503	{
2504	Plane3 plane = mPolys[i]->GetSupportingPlane();
2505
2506	const int cf =
2507	planePoly->ClassifyPlane(plane, 0.00001f);
2508
2509	// split new polygon with all previous planes
2510	switch (cf)
2511	{
2512	case Polygon3::SPLIT:
2513	{
2514	VertexContainer splitPts;
2515
2516	Polygon3 *frontPoly = new Polygon3();
2517	Polygon3 *backPoly = new Polygon3();
2518
2519	planePoly->Split(plane, frontPoly, backPoly, splitPts);
2520	DEL_PTR(planePoly);
2521
2522	if (backPoly->Valid())
2523	planePoly = backPoly;
2524	else
2525	DEL_PTR(backPoly);
2526	}
2527	break;
2528	case Polygon3::FRONT_SIDE:
2529	DEL_PTR(planePoly);
2530	break;
2531	// polygon is taken as it is
2532	case Polygon3::BACK_SIDE:
2533	case Polygon3::COINCIDENT:
2534	default:
2535	break;
2536	}
2537	}
2538
2539	return planePoly;
2540	}
2541
2542	void BspViewCellsStatistics::Print(ostream &app) const
2543	{
2544	app << "===== BspViewCells statistics ===============\n";
2545
2546	app << setprecision(4);
2547
2548	//app << "#N_CTIME ( Construction time [s] )\n" << Time() << " \n";
2549
2550	app << "#N_OVERALLPVS ( objects in PVS )\n" << pvs << endl;
2551
2552	app << "#N_PMAXPVS ( largest PVS )\n" << maxPvs << endl;
2553
2554	app << "#N_PMINPVS ( smallest PVS )\n" << minPvs << endl;
2555
2556	app << "#N_PAVGPVS ( average PVS )\n" << AvgPvs() << endl;
2557
2558	app << "#N_PEMPTYPVS ( view cells with PVS smaller 2 )\n" << emptyPvs << endl;
2559
2560	app << "#N_VIEWCELLS ( number of view cells)\n" << viewCells << endl;
2561
2562	app << "#N_AVGBSPLEAVES (average number of BSP leaves per view cell )\n" << AvgBspLeaves() << endl;
2563
2564	app << "#N_MAXBSPLEAVES ( maximal number of BSP leaves per view cell )\n" << maxBspLeaves << endl;
2565
2566	app << "===== END OF BspViewCells statistics ==========\n";
2567	}

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

Download in other formats: