Context Navigation

source: GTP/trunk/Lib/Vis/Preprocessing/src/BvHierarchy.cpp @ 1830

Revision 1830, 71.1 KB checked in by mattausch, 18 years ago (diff)

Line
1	#include <stack>
2	#include <time.h>
3	#include <iomanip>
4
5	#include "BvHierarchy.h"
6	#include "ViewCell.h"
7	#include "Plane3.h"
8	#include "Mesh.h"
9	#include "common.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 "ViewCellsManager.h"
17	#include "Beam.h"
18	#include "VspTree.h"
19	#include "HierarchyManager.h"
20
21
22	namespace GtpVisibilityPreprocessor {
23
24
25	#define PROBABILIY_IS_BV_VOLUME 1
26	#define USE_FIXEDPOINT_T 0
27	#define USE_VOLUMES_FOR_HEURISTICS 1
28
29	//int BvhNode::sMailId = 10000; //2147483647;
30	//int BvhNode::sReservedMailboxes = 1;
31
32	BvHierarchy *BvHierarchy::BvhSubdivisionCandidate::sBvHierarchy = NULL;
33
34
35	/// sorting operator
36	inline static bool ilt(Intersectable obj1, Intersectable obj2)
37	{
38	return obj1->mId < obj2->mId;
39	}
40
41
42	/// sorting operator
43	inline static bool smallerSize(Intersectable obj1, Intersectable obj2)
44	{
45	return obj1->GetBox().SurfaceArea() < obj2->GetBox().SurfaceArea();
46	}
47
48	/***************************************************************/
49	/* class BvhNode implementation */
50	/***************************************************************/
51
52	BvhNode::BvhNode():
53	mParent(NULL),
54	mTimeStamp(0),
55	mRenderCost(-1)
56
57	{
58
59	}
60
61	BvhNode::BvhNode(const AxisAlignedBox3 &bbox):
62	mParent(NULL),
63	mBoundingBox(bbox),
64	mTimeStamp(0),
65	mRenderCost(-1)
66	{
67	}
68
69
70	BvhNode::BvhNode(const AxisAlignedBox3 &bbox, BvhInterior *parent):
71	mBoundingBox(bbox),
72	mParent(parent),
73	mTimeStamp(0),
74	mRenderCost(-1)
75	{
76	}
77
78
79	bool BvhNode::IsRoot() const
80	{
81	return mParent == NULL;
82	}
83
84
85	BvhInterior *BvhNode::GetParent()
86	{
87	return mParent;
88	}
89
90
91	void BvhNode::SetParent(BvhInterior *parent)
92	{
93	mParent = parent;
94	}
95
96
97	int BvhNode::GetRandomEdgePoint(Vector3 &point,
98	Vector3 &normal)
99	{
100	// get random edge
101	const int idx = Random(12);
102	Vector3 a, b;
103	mBoundingBox.GetEdge(idx, &a, &b);
104
105	const float w = RandomValue(0.0f, 1.0f);
106
107	point = a * w + b * (1.0f - w);
108
109	// TODO
110	normal = Vector3(0);
111
112	return idx;
113	}
114
115
116
117	/******************************************************************/
118	/* class BvhInterior implementation */
119	/******************************************************************/
120
121
122	BvhLeaf::BvhLeaf(const AxisAlignedBox3 &bbox):
123	BvhNode(bbox),
124	mSubdivisionCandidate(NULL),
125	mGlList(0)
126	{
127	mActiveNode = this;
128	}
129
130
131	BvhLeaf::BvhLeaf(const AxisAlignedBox3 &bbox, BvhInterior *parent):
132	BvhNode(bbox, parent),
133	mGlList(0)
134
135	{
136	mActiveNode = this;
137	}
138
139
140	BvhLeaf::BvhLeaf(const AxisAlignedBox3 &bbox,
141	BvhInterior *parent,
142	const int numObjects):
143	BvhNode(bbox, parent),
144	mGlList(0)
145
146	{
147	mObjects.reserve(numObjects);
148	mActiveNode = this;
149	}
150
151
152	bool BvhLeaf::IsLeaf() const
153	{
154	return true;
155	}
156
157
158	BvhLeaf::~BvhLeaf()
159	{
160	}
161
162
163	void BvhLeaf::CollectObjects(ObjectContainer &objects)
164	{
165	ObjectContainer::const_iterator oit, oit_end = mObjects.end();
166	for (oit = mObjects.begin(); oit != oit_end; ++ oit)
167	{
168	objects.push_back(*oit);
169	}
170	}
171
172	/******************************************************************/
173	/* class BvhInterior implementation */
174	/******************************************************************/
175
176
177	BvhInterior::BvhInterior(const AxisAlignedBox3 &bbox):
178	BvhNode(bbox), mFront(NULL), mBack(NULL)
179	{
180	}
181
182
183	BvhInterior::BvhInterior(const AxisAlignedBox3 &bbox, BvhInterior *parent):
184	BvhNode(bbox, parent), mFront(NULL), mBack(NULL)
185	{
186	}
187
188
189	void BvhInterior::ReplaceChildLink(BvhNode oldChild, BvhNode newChild)
190	{
191	if (mBack == oldChild)
192	mBack = newChild;
193	else
194	mFront = newChild;
195	}
196
197
198	bool BvhInterior::IsLeaf() const
199	{
200	return false;
201	}
202
203
204	BvhInterior::~BvhInterior()
205	{
206	DEL_PTR(mFront);
207	DEL_PTR(mBack);
208	}
209
210
211	void BvhInterior::SetupChildLinks(BvhNode front, BvhNode back)
212	{
213	mBack = back;
214	mFront = front;
215	}
216
217
218	void BvhInterior::CollectObjects(ObjectContainer &objects)
219	{
220	mFront->CollectObjects(objects);
221	mBack->CollectObjects(objects);
222	}
223
224
225	/*******************************************************************/
226	/* class BvHierarchy implementation */
227	/*******************************************************************/
228
229
230	BvHierarchy::BvHierarchy():
231	mRoot(NULL),
232	mTimeStamp(1),
233	mIsInitialSubdivision(false)
234	{
235	ReadEnvironment();
236	mSubdivisionCandidates = new SortableEntryContainer;
237	// for (int i = 0; i < 4; ++ i)
238	// mSortedObjects[i] = NULL;
239	}
240
241
242	BvHierarchy::~BvHierarchy()
243	{
244	// delete the local subdivision candidates
245	DEL_PTR(mSubdivisionCandidates);
246
247	// delete the presorted objects
248	/*for (int i = 0; i < 4; ++ i)
249	{
250	DEL_PTR(mSortedObjects[i]);
251	}*/
252
253	// delete the tree
254	DEL_PTR(mRoot);
255	}
256
257
258	void BvHierarchy::ReadEnvironment()
259	{
260	bool randomize = false;
261	Environment::GetSingleton()->GetBoolValue("BvHierarchy.Construction.randomize", randomize);
262
263	// initialise random generator for heuristics
264	if (randomize)
265	Randomize();
266
267	//////////////////////////////
268	//-- termination criteria for autopartition
269
270	Environment::GetSingleton()->GetIntValue("BvHierarchy.Termination.maxDepth", mTermMaxDepth);
271	Environment::GetSingleton()->GetIntValue("BvHierarchy.Termination.maxLeaves", mTermMaxLeaves);
272	Environment::GetSingleton()->GetIntValue("BvHierarchy.Termination.minObjects", mTermMinObjects);
273	Environment::GetSingleton()->GetIntValue("BvHierarchy.Termination.minRays", mTermMinRays);
274	Environment::GetSingleton()->GetFloatValue(
275	"BvHierarchy.Termination.minProbability", mTermMinProbability);
276	Environment::GetSingleton()->GetIntValue("BvHierarchy.Termination.missTolerance", mTermMissTolerance);
277
278
279	//////////////////////////////
280	//-- max cost ratio for early tree termination
281
282	Environment::GetSingleton()->GetFloatValue("BvHierarchy.Termination.maxCostRatio", mTermMaxCostRatio);
283	Environment::GetSingleton()->GetFloatValue("BvHierarchy.Termination.minGlobalCostRatio",
284	mTermMinGlobalCostRatio);
285	Environment::GetSingleton()->GetIntValue("BvHierarchy.Termination.globalCostMissTolerance",
286	mTermGlobalCostMissTolerance);
287
288
289	//////////////////////////////
290	//-- factors for subdivision heuristics
291
292	// if only the driving axis is used for splits
293	Environment::GetSingleton()->GetBoolValue("BvHierarchy.splitUseOnlyDrivingAxis", mOnlyDrivingAxis);
294	Environment::GetSingleton()->GetFloatValue("BvHierarchy.maxStaticMemory", mMaxMemory);
295	Environment::GetSingleton()->GetBoolValue("BvHierarchy.useCostHeuristics", mUseCostHeuristics);
296	Environment::GetSingleton()->GetBoolValue("BvHierarchy.useSah", mUseSah);
297
298	char subdivisionStatsLog[100];
299	Environment::GetSingleton()->GetStringValue("BvHierarchy.subdivisionStats", subdivisionStatsLog);
300	mSubdivisionStats.open(subdivisionStatsLog);
301
302	Environment::GetSingleton()->GetFloatValue("BvHierarchy.Construction.renderCostDecreaseWeight", mRenderCostDecreaseWeight);
303	Environment::GetSingleton()->GetBoolValue("BvHierarchy.Construction.useGlobalSorting", mUseGlobalSorting);
304	Environment::GetSingleton()->GetIntValue("BvHierarchy.minRaysForVisibility", mMinRaysForVisibility);
305	Environment::GetSingleton()->GetIntValue("BvHierarchy.maxTests", mMaxTests);
306	Environment::GetSingleton()->GetBoolValue("BvHierarchy.Construction.useInitialSubdivision", mApplyInitialPartition);
307	Environment::GetSingleton()->GetIntValue("BvHierarchy.Construction.Initial.minObjects", mInitialMinObjects);
308	Environment::GetSingleton()->GetFloatValue("BvHierarchy.Construction.Initial.maxAreaRatio", mInitialMaxAreaRatio);
309	Environment::GetSingleton()->GetFloatValue("BvHierarchy.Construction.Initial.minArea", mInitialMinArea);
310
311	//mMemoryConst = (float)(sizeof(VspLeaf) + sizeof(VspViewCell));
312	//mMemoryConst = (float)sizeof(BvhLeaf);
313	mMemoryConst = 16;//(float)sizeof(ObjectContainer);
314
315	mUseBboxAreaForSah = true;
316
317	/////////////
318	//-- debug output
319
320	Debug << "***** Bvh hierarchy options ****** " << endl;
321	Debug << "max depth: " << mTermMaxDepth << endl;
322	Debug << "min probabiliy: " << mTermMinProbability<< endl;
323	Debug << "min objects: " << mTermMinObjects << endl;
324	Debug << "max cost ratio: " << mTermMaxCostRatio << endl;
325	Debug << "miss tolerance: " << mTermMissTolerance << endl;
326	Debug << "max leaves: " << mTermMaxLeaves << endl;
327	Debug << "randomize: " << randomize << endl;
328	Debug << "min global cost ratio: " << mTermMinGlobalCostRatio << endl;
329	Debug << "global cost miss tolerance: " << mTermGlobalCostMissTolerance << endl;
330	Debug << "only driving axis: " << mOnlyDrivingAxis << endl;
331	Debug << "max memory: " << mMaxMemory << endl;
332	Debug << "use cost heuristics: " << mUseCostHeuristics << endl;
333	Debug << "use surface area heuristics: " << mUseSah << endl;
334	Debug << "subdivision stats log: " << subdivisionStatsLog << endl;
335	Debug << "split borders: " << mSplitBorder << endl;
336	Debug << "render cost decrease weight: " << mRenderCostDecreaseWeight << endl;
337	Debug << "use global sort: " << mUseGlobalSorting << endl;
338	Debug << "minimal rays for visibility: " << mMinRaysForVisibility << endl;
339	Debug << "bvh mem const: " << mMemoryConst << endl;
340	Debug << "apply initial partition: " << mApplyInitialPartition << endl;
341	Debug << "min objects: " << mInitialMinObjects << endl;
342	Debug << "max area ratio: " << mInitialMaxAreaRatio << endl;
343	Debug << "min area: " << mInitialMinArea << endl;
344
345	Debug << endl;
346	}
347
348
349	void BvHierarchy::AssociateObjectsWithLeaf(BvhLeaf *leaf)
350	{
351	ObjectContainer::const_iterator oit, oit_end = leaf->mObjects.end();
352
353	for (oit = leaf->mObjects.begin(); oit != oit_end; ++ oit)
354	{
355	(*oit)->mBvhLeaf = leaf;
356	}
357	}
358
359
360	static int CountRays(const ObjectContainer &objects)
361	{
362	int nRays = 0;
363
364	ObjectContainer::const_iterator oit, oit_end = objects.end();
365
366	for (oit = objects.begin(); oit != oit_end; ++ oit)
367	{
368	nRays += (int)(*oit)->GetOrCreateRays()->size();
369	}
370
371	return nRays;
372	}
373
374
375	BvhInterior *BvHierarchy::SubdivideNode(const BvhSubdivisionCandidate &sc,
376	BvhTraversalData &frontData,
377	BvhTraversalData &backData)
378	{
379	const BvhTraversalData &tData = sc.mParentData;
380	BvhLeaf *leaf = tData.mNode;
381	AxisAlignedBox3 parentBox = leaf->GetBoundingBox();
382
383	// update stats: we have two new leaves
384	mBvhStats.nodes += 2;
385
386	if (tData.mDepth > mBvhStats.maxDepth)
387	{
388	mBvhStats.maxDepth = tData.mDepth;
389	}
390
391	// add the new nodes to the tree
392	BvhInterior *node = new BvhInterior(parentBox, leaf->GetParent());
393
394
395	//////////////////
396	//-- create front and back leaf
397
398	AxisAlignedBox3 fbox = EvalBoundingBox(sc.mFrontObjects, &parentBox);
399	AxisAlignedBox3 bbox = EvalBoundingBox(sc.mBackObjects, &parentBox);
400
401	BvhLeaf *back = new BvhLeaf(bbox, node, (int)sc.mBackObjects.size());
402	BvhLeaf *front = new BvhLeaf(fbox, node, (int)sc.mFrontObjects.size());
403
404	BvhInterior *parent = leaf->GetParent();
405
406	// replace a link from node's parent
407	if (parent)
408	{
409	parent->ReplaceChildLink(leaf, node);
410	node->SetParent(parent);
411	}
412	else // no parent => this node is the root
413	{
414	mRoot = node;
415	}
416
417	// and setup child links
418	node->SetupChildLinks(front, back);
419
420	++ mBvhStats.splits;
421
422
423	////////////////////////////////////////
424	//-- fill front and back traversal data with the new values
425
426	frontData.mDepth = backData.mDepth = tData.mDepth + 1;
427
428	frontData.mNode = front;
429	backData.mNode = back;
430
431	back->mObjects = sc.mBackObjects;
432	front->mObjects = sc.mFrontObjects;
433
434	// if the number of rays is too low, no assumptions can be made
435	// (=> switch to surface area heuristics?)
436	frontData.mNumRays = CountRays(sc.mFrontObjects);
437	backData.mNumRays = CountRays(sc.mBackObjects);
438
439	AssociateObjectsWithLeaf(back);
440	AssociateObjectsWithLeaf(front);
441
442	#if PROBABILIY_IS_BV_VOLUME
443	// volume of bvh (= probability that this bvh can be seen)
444	frontData.mProbability = fbox.GetVolume();
445	backData.mProbability = bbox.GetVolume();
446	#else
447	// compute probability of this node being visible,
448	// i.e., volume of the view cells that can see this node
449	frontData.mProbability = EvalViewCellsVolume(sc.mFrontObjects);
450	backData.mProbability = EvalViewCellsVolume(sc.mBackObjects);
451	#endif
452
453	// how often was max cost ratio missed in this branch?
454	frontData.mMaxCostMisses = sc.GetMaxCostMisses();
455	backData.mMaxCostMisses = sc.GetMaxCostMisses();
456
457	// set the time stamp so the order of traversal can be reconstructed
458	node->SetTimeStamp(mHierarchyManager->mTimeStamp ++);
459
460	// assign the objects in sorted order
461	if (mUseGlobalSorting)
462	{
463	AssignSortedObjects(sc, frontData, backData);
464	}
465
466	// return the new interior node
467	return node;
468	}
469
470
471	BvhNode *BvHierarchy::Subdivide(SplitQueue &tQueue,
472	SubdivisionCandidate *splitCandidate,
473	const bool globalCriteriaMet)
474	{
475	BvhSubdivisionCandidate *sc =
476	dynamic_cast<BvhSubdivisionCandidate *>(splitCandidate);
477	BvhTraversalData &tData = sc->mParentData;
478
479	BvhNode *currentNode = tData.mNode;
480
481	if (!LocalTerminationCriteriaMet(tData) && !globalCriteriaMet)
482	{
483	//////////////
484	//-- continue subdivision
485
486	BvhTraversalData tFrontData;
487	BvhTraversalData tBackData;
488
489	// create new interior node and two leaf node
490	currentNode = SubdivideNode(*sc, tFrontData, tBackData);
491
492	// decrease the weighted average cost of the subdivisoin
493	mTotalCost -= sc->GetRenderCostDecrease();
494	mPvsEntries += sc->GetPvsEntriesIncr();
495
496	// subdivision statistics
497	if (1) PrintSubdivisionStats(*sc);
498
499
500	///////////////////////////
501	//-- push the new split candidates on the queue
502
503	BvhSubdivisionCandidate *frontCandidate =
504	new BvhSubdivisionCandidate(tFrontData);
505	BvhSubdivisionCandidate *backCandidate =
506	new BvhSubdivisionCandidate(tBackData);
507
508	EvalSubdivisionCandidate(*frontCandidate);
509	EvalSubdivisionCandidate(*backCandidate);
510
511	// cross reference
512	tFrontData.mNode->SetSubdivisionCandidate(frontCandidate);
513	tBackData.mNode->SetSubdivisionCandidate(backCandidate);
514
515	//cout << "f: " << frontCandidate->GetPriority() << " b: " << backCandidate->GetPriority() << endl;
516	tQueue.Push(frontCandidate);
517	tQueue.Push(backCandidate);
518	}
519
520	/////////////////////////////////
521	//-- node is a leaf => terminate traversal
522
523	if (currentNode->IsLeaf())
524	{
525	/////////////////////
526	//-- store additional info
527	EvaluateLeafStats(tData);
528
529	// this leaf is no candidate for splitting anymore
530	// => detach subdivision candidate
531	tData.mNode->SetSubdivisionCandidate(NULL);
532	// detach node so we don't delete it with the traversal data
533	tData.mNode = NULL;
534	}
535
536	return currentNode;
537	}
538
539
540	float BvHierarchy::EvalPriority(const BvhSubdivisionCandidate &splitCandidate,
541	const float renderCostDecr,
542	const float oldRenderCost) const
543	{
544	float priority;
545
546	if (mIsInitialSubdivision)
547	{
548	priority = (float)-splitCandidate.mParentData.mDepth;
549	return priority;
550	}
551
552	BvhLeaf *leaf = splitCandidate.mParentData.mNode;
553
554	// surface area heuristics is used when there is
555	// no view space subdivision available.
556	// In order to have some prioritized traversal,
557	// we use this formula instead
558	if (mHierarchyManager->GetViewSpaceSubdivisionType() ==
559	HierarchyManager::NO_VIEWSPACE_SUBDIV)
560	{
561	priority = EvalSahCost(leaf);
562	}
563	else
564	{
565	// take render cost of node into account
566	// otherwise danger of being stuck in a local minimum!
567	const float factor = mRenderCostDecreaseWeight;
568
569	priority = factor * renderCostDecr + (1.0f - factor) * oldRenderCost;
570
571	if (mHierarchyManager->mConsiderMemory)
572	{
573	priority /= ((float)splitCandidate.GetPvsEntriesIncr() + mMemoryConst);
574	}
575	}
576
577	// hack: don't allow empty splits to be taken
578	if (splitCandidate.mFrontObjects.empty() \|\| splitCandidate.mBackObjects.empty())
579	priority = 0;
580
581	return priority;
582	}
583
584
585	void BvHierarchy::EvalSubdivisionCandidate(BvhSubdivisionCandidate &splitCandidate,
586	bool computeSplitPlane)
587	{
588	if (computeSplitPlane)
589	{
590	const bool sufficientSamples =
591	splitCandidate.mParentData.mNumRays > mMinRaysForVisibility;
592
593	const bool useVisibiliyBasedHeuristics =
594	mUseSah &&
595	(mHierarchyManager->GetViewSpaceSubdivisionType() ==
596	HierarchyManager::KD_BASED_VIEWSPACE_SUBDIV) &&
597	sufficientSamples;
598
599	// compute best object partition
600	const float ratio = SelectObjectPartition(splitCandidate.mParentData,
601	splitCandidate.mFrontObjects,
602	splitCandidate.mBackObjects,
603	useVisibiliyBasedHeuristics);
604
605	// cost ratio violated?
606	const bool maxCostRatioViolated = mTermMaxCostRatio < ratio;
607	const int previousMisses = splitCandidate.mParentData.mMaxCostMisses;
608
609	splitCandidate.SetMaxCostMisses(
610	maxCostRatioViolated ? previousMisses + 1 : previousMisses);
611	}
612
613	BvhLeaf *leaf = splitCandidate.mParentData.mNode;
614
615	const float oldProp = EvalViewCellsVolume(leaf->mObjects);
616	const float oldRenderCost = EvalRenderCost(leaf->mObjects);
617
618	// compute global decrease in render cost
619	const float newRenderCost = EvalRenderCost(splitCandidate.mFrontObjects) +
620	EvalRenderCost(splitCandidate.mBackObjects);
621
622	const float renderCostDecr = oldRenderCost - newRenderCost;
623
624	splitCandidate.SetRenderCostDecrease(renderCostDecr);
625
626	// increase in pvs entries
627	const int pvsEntriesIncr = EvalPvsEntriesIncr(splitCandidate);
628	splitCandidate.SetPvsEntriesIncr(pvsEntriesIncr);
629
630	#ifdef GTP_DEBUG
631	Debug << "old render cost: " << oldRenderCost << endl;
632	Debug << "new render cost: " << newRenderCost << endl;
633	Debug << "render cost decrease: " << renderCostDecr << endl;
634	#endif
635
636	const float priority = EvalPriority(splitCandidate,
637	oldRenderCost,
638	renderCostDecr);
639
640	// compute global decrease in render cost
641	splitCandidate.SetPriority(priority);
642	}
643
644
645	int BvHierarchy::EvalPvsEntriesIncr(BvhSubdivisionCandidate &splitCandidate) const
646	{
647	const int oldPvsSize = CountViewCells(splitCandidate.mParentData.mNode->mObjects);
648
649	const int fPvsSize = CountViewCells(splitCandidate.mFrontObjects);
650	const int bPvsSize = CountViewCells(splitCandidate.mBackObjects);
651
652	return fPvsSize + bPvsSize - oldPvsSize;
653	}
654
655
656	inline bool BvHierarchy::LocalTerminationCriteriaMet(const BvhTraversalData &tData) const
657	{
658	const bool terminationCriteriaMet =
659	(0
660	\|\| ((int)tData.mNode->mObjects.size() <= 1)//mTermMinObjects)
661	//\|\| (data.mProbability <= mTermMinProbability)
662	//\|\| (data.mNumRays <= mTermMinRays)
663	);
664
665	#ifdef _DEBUG
666	if (terminationCriteriaMet)
667	{
668	cout << "bvh local termination criteria met:" << endl;
669	cout << "objects: " << tData.mNode->mObjects.size() << " " << mTermMinObjects << endl;
670	}
671	#endif
672	return terminationCriteriaMet;
673	}
674
675
676	inline bool BvHierarchy::GlobalTerminationCriteriaMet(const BvhTraversalData &data) const
677	{
678	// note: tracking for global cost termination
679	// does not make much sense for interleaved vsp / osp partition
680	// as it is the responsibility of the hierarchy manager
681
682	const bool terminationCriteriaMet =
683	(0
684	\|\| (mBvhStats.Leaves() >= mTermMaxLeaves)
685	//\|\| (mBvhStats.mGlobalCostMisses >= mTermGlobalCostMissTolerance)
686	//\|\| mOutOfMemory
687	);
688
689	#ifdef GTP_DEBUG
690	if (terminationCriteriaMet)
691	{
692	Debug << "bvh global termination criteria met:" << endl;
693	Debug << "cost misses: " << mBvhStats.mGlobalCostMisses << " " << mTermGlobalCostMissTolerance << endl;
694	Debug << "leaves: " << mBvhStats.Leaves() << " " << mTermMaxLeaves << endl;
695	}
696	#endif
697	return terminationCriteriaMet;
698	}
699
700
701	void BvHierarchy::EvaluateLeafStats(const BvhTraversalData &data)
702	{
703	// the node became a leaf -> evaluate stats for leafs
704	BvhLeaf *leaf = data.mNode;
705
706	++ mCreatedLeaves;
707
708
709	if (data.mProbability <= mTermMinProbability)
710	{
711	++ mBvhStats.minProbabilityNodes;
712	}
713
714	////////////////////////////////////////////
715	// depth related stuff
716
717	if (data.mDepth < mBvhStats.minDepth)
718	{
719	mBvhStats.minDepth = data.mDepth;
720	}
721
722	if (data.mDepth >= mTermMaxDepth)
723	{
724	++ mBvhStats.maxDepthNodes;
725	}
726
727	// accumulate depth to compute average depth
728	mBvhStats.accumDepth += data.mDepth;
729
730
731	////////////////////////////////////////////
732	// objects related stuff
733
734	// note: the sum should alwaysbe total number of objects for bvh
735	mBvhStats.objectRefs += (int)leaf->mObjects.size();
736
737	if ((int)leaf->mObjects.size() <= mTermMinObjects)
738	{
739	++ mBvhStats.minObjectsNodes;
740	}
741
742	if (leaf->mObjects.empty())
743	{
744	++ mBvhStats.emptyNodes;
745	}
746
747	if ((int)leaf->mObjects.size() > mBvhStats.maxObjectRefs)
748	{
749	mBvhStats.maxObjectRefs = (int)leaf->mObjects.size();
750	}
751
752	if ((int)leaf->mObjects.size() < mBvhStats.minObjectRefs)
753	{
754	mBvhStats.minObjectRefs = (int)leaf->mObjects.size();
755	}
756
757	////////////////////////////////////////////
758	// ray related stuff
759
760	// note: this number should always accumulate to the total number of rays
761	mBvhStats.rayRefs += data.mNumRays;
762
763	if (data.mNumRays <= mTermMinRays)
764	{
765	++ mBvhStats.minRaysNodes;
766	}
767
768	if (data.mNumRays > mBvhStats.maxRayRefs)
769	{
770	mBvhStats.maxRayRefs = data.mNumRays;
771	}
772
773	if (data.mNumRays < mBvhStats.minRayRefs)
774	{
775	mBvhStats.minRayRefs = data.mNumRays;
776	}
777
778	#ifdef _DEBUG
779	cout << "depth: " << data.mDepth << " objects: " << (int)leaf->mObjects.size()
780	<< " rays: " << data.mNumRays << " rays / objects "
781	<< (float)data.mNumRays / (float)leaf->mObjects.size() << endl;
782	#endif
783	}
784
785
786	#if 0
787
788	/// compute object boundaries using spatial mid split
789	float BvHierarchy::EvalLocalObjectPartition(const BvhTraversalData &tData,
790	const int axis,
791	ObjectContainer &objectsFront,
792	ObjectContainer &objectsBack)
793	{
794	const float maxBox = tData.mBoundingBox.Max(axis);
795	const float minBox = tData.mBoundingBox.Min(axis);
796
797	float midPoint = (maxBox + minBox) * 0.5f;
798
799	ObjectContainer::const_iterator oit, oit_end = tData.mNode->mObjects.end();
800
801	for (oit = tData.mNode->mObjects.begin(); oit != oit_end; ++ oit)
802	{
803	Intersectable obj = oit;
804	const AxisAlignedBox3 box = obj->GetBox();
805
806	const float objMid = (box.Max(axis) + box.Min(axis)) * 0.5f;
807
808	// object mailed => belongs to back objects
809	if (objMid < midPoint)
810	{
811	objectsBack.push_back(obj);
812	}
813	else
814	{
815	objectsFront.push_back(obj);
816	}
817	}
818
819	const float oldRenderCost = EvalRenderCost(tData.mNode->mObjects);
820	const float newRenderCost = EvalRenderCost(objectsFront) * EvalRenderCost(objectsBack);
821
822	const float ratio = newRenderCost / oldRenderCost;
823	return ratio;
824	}
825
826	#else
827
828	/// compute object partition by getting balanced objects on the left and right side
829	float BvHierarchy::EvalLocalObjectPartition(const BvhTraversalData &tData,
830	const int axis,
831	ObjectContainer &objectsFront,
832	ObjectContainer &objectsBack)
833	{
834	PrepareLocalSubdivisionCandidates(tData, axis);
835
836	SortableEntryContainer::const_iterator cit, cit_end = mSubdivisionCandidates->end();
837
838	int i = 0;
839	const int border = (int)tData.mNode->mObjects.size() / 2;
840
841	for (cit = mSubdivisionCandidates->begin(); cit != cit_end; ++ cit, ++ i)
842	{
843	Intersectable obj = (cit).mObject;
844
845	// object mailed => belongs to back objects
846	if (i < border)
847	{
848	objectsBack.push_back(obj);
849	}
850	else
851	{
852	objectsFront.push_back(obj);
853	}
854	}
855
856	#if 1
857	const float cost = (tData.mNode->GetBoundingBox().Size().DrivingAxis() == axis) ? -1.0f : 0.0f;
858	#else
859	const float oldRenderCost = EvalRenderCost(tData.mNode->mObjects);
860	const float newRenderCost = EvalRenderCost(objectsFront) * EvalRenderCost(objectsBack);
861
862	const float cost = newRenderCost / oldRenderCost;
863	#endif
864
865	return cost;
866	}
867	#endif
868
869	#if 1
870
871	float BvHierarchy::EvalSah(const BvhTraversalData &tData,
872	const int axis,
873	ObjectContainer &objectsFront,
874	ObjectContainer &objectsBack)
875	{
876	// go through the lists, count the number of objects left and right
877	// and evaluate the following cost funcion:
878	// C = ct_div_ci + (ol + or) / queries
879	PrepareLocalSubdivisionCandidates(tData, axis);
880
881	const float totalRenderCost = EvalAbsCost(tData.mNode->mObjects);
882	float objectsLeft = 0, objectsRight = totalRenderCost;
883
884	const AxisAlignedBox3 nodeBbox = tData.mNode->GetBoundingBox();
885	const float boxArea = nodeBbox.SurfaceArea();
886
887	float minSum = 1e20f;
888
889	float minBorder = nodeBbox.Max(axis);
890	float maxBorder = nodeBbox.Min(axis);
891
892	float areaLeft = 0, areaRight = 0;
893
894	SortableEntryContainer::const_iterator currentPos =
895	mSubdivisionCandidates->begin();
896
897	vector<float> bordersRight;
898
899	// we keep track of both borders of the bounding boxes =>
900	// store the events in descending order
901
902	bordersRight.resize(mSubdivisionCandidates->size());
903
904	SortableEntryContainer::reverse_iterator rcit =
905	mSubdivisionCandidates->rbegin(), rcit_end =
906	mSubdivisionCandidates->rend();
907
908	vector<float>::reverse_iterator rbit = bordersRight.rbegin();
909
910	for (; rcit != rcit_end; ++ rcit, ++ rbit)
911	{
912	Intersectable obj = (rcit).mObject;
913	const AxisAlignedBox3 obox = obj->GetBox();
914
915	if (obox.Min(axis) < minBorder)
916	{
917	minBorder = obox.Min(axis);
918	}
919
920	(*rbit) = minBorder;
921	}
922
923	// temporary surface areas
924	float al = 0;
925	float ar = boxArea;
926
927	vector<float>::const_iterator bit = bordersRight.begin();
928	SortableEntryContainer::const_iterator cit, cit_end = mSubdivisionCandidates->end();
929
930	for (cit = mSubdivisionCandidates->begin(); cit != cit_end; ++ cit, ++ bit)
931	{
932	Intersectable obj = (cit).mObject;
933
934	const float renderCost = mViewCellsManager->EvalRenderCost(obj);
935
936	objectsLeft += renderCost;
937	objectsRight -= renderCost;
938
939	const AxisAlignedBox3 obox = obj->GetBox();
940
941	// the borders of the bounding boxes have changed
942	if (obox.Max(axis) > maxBorder)
943	{
944	maxBorder = obox.Max(axis);
945	}
946
947	minBorder = (*bit);
948
949	AxisAlignedBox3 lbox = nodeBbox;
950	AxisAlignedBox3 rbox = nodeBbox;
951
952	lbox.SetMax(axis, maxBorder);
953	rbox.SetMin(axis, minBorder);
954
955	al = lbox.SurfaceArea();
956	ar = rbox.SurfaceArea();
957
958	const bool noValidSplit = ((objectsLeft <= Limits::Small) \|\| (objectsRight <= Limits::Small));
959	const float sum = noValidSplit ? 1e25 : objectsLeft * al + objectsRight * ar;
960
961	/cout << "pos=" << (cit).mPos << "\t q=(" << objectsLeft << "," << objectsRight <<")\t r=("
962	<< lbox.SurfaceArea() << "," << rbox.SurfaceArea() << ")" << endl;
963	cout << "minborder: " << minBorder << " maxborder: " << maxBorder << endl;
964	cout << "cost= " << sum << endl;*/
965
966	if (sum < minSum)
967	{
968	minSum = sum;
969	areaLeft = al;
970	areaRight = ar;
971
972	// objects belong to left side now
973	for (; currentPos != (cit + 1); ++ currentPos);
974	}
975	}
976
977	////////////
978	//-- assign object to front and back volume
979
980	// belongs to back bv
981	for (cit = mSubdivisionCandidates->begin(); cit != currentPos; ++ cit)
982	objectsBack.push_back((*cit).mObject);
983
984	// belongs to front bv
985	for (cit = currentPos; cit != cit_end; ++ cit)
986	objectsFront.push_back((*cit).mObject);
987
988	float newCost = minSum / boxArea;
989	float ratio = newCost / totalRenderCost;
990
991	#ifdef GTP_DEBUG
992	cout << "\n\nobjects=(" << (int)objectsBack.size() << "," << (int)objectsFront.size() << " of "
993	<< (int)tData.mNode->mObjects.size() << ")\t area=("
994	<< areaLeft << ", " << areaRight << ", " << boxArea << ")" << endl
995	<< "cost= " << newCost << " oldCost=" << totalRenderCost / boxArea << endl;
996	#endif
997
998	return ratio;
999	}
1000
1001	#else
1002
1003	float BvHierarchy::EvalSah(const BvhTraversalData &tData,
1004	const int axis,
1005	ObjectContainer &objectsFront,
1006	ObjectContainer &objectsBack)
1007	{
1008	// go through the lists, count the number of objects left and right
1009	// and evaluate the following cost funcion:
1010	// C = ct_div_ci + (ol + or) / queries
1011	PrepareLocalSubdivisionCandidates(tData, axis);
1012
1013	const float totalRenderCost = EvalAbsCost(tData.mNode->mObjects);
1014	float objectsLeft = 0, objectsRight = totalRenderCost;
1015
1016	const AxisAlignedBox3 nodeBbox = tData.mNode->GetBoundingBox();
1017
1018	const float minBox = nodeBbox.Min(axis);
1019	const float maxBox = nodeBbox.Max(axis);
1020	const float boxArea = nodeBbox.SurfaceArea();
1021
1022	float minSum = 1e20f;
1023
1024	Vector3 minBorder = nodeBbox.Max();
1025	Vector3 maxBorder = nodeBbox.Min();
1026
1027	float areaLeft = 0, areaRight = 0;
1028
1029	SortableEntryContainer::const_iterator currentPos =
1030	mSubdivisionCandidates->begin();
1031
1032	vector<Vector3> bordersRight;
1033
1034	// we keep track of both borders of the bounding boxes =>
1035	// store the events in descending order
1036	bordersRight.resize(mSubdivisionCandidates->size());
1037
1038	SortableEntryContainer::reverse_iterator rcit =
1039	mSubdivisionCandidates->rbegin(), rcit_end =
1040	mSubdivisionCandidates->rend();
1041
1042	vector<Vector3>::reverse_iterator rbit = bordersRight.rbegin();
1043
1044	for (; rcit != rcit_end; ++ rcit, ++ rbit)
1045	{
1046	Intersectable obj = (rcit).mObject;
1047	const AxisAlignedBox3 obox = obj->GetBox();
1048
1049	for (int i = 0; i < 3; ++ i)
1050	{
1051	if (obox.Min(i) < minBorder[i])
1052	{
1053	minBorder[i] = obox.Min(i);
1054	}
1055	}
1056
1057	(*rbit) = minBorder;
1058	}
1059
1060	// temporary surface areas
1061	float al = 0;
1062	float ar = boxArea;
1063
1064	vector<Vector3>::const_iterator bit = bordersRight.begin();
1065	SortableEntryContainer::const_iterator cit, cit_end =
1066	mSubdivisionCandidates->end();
1067
1068	for (cit = mSubdivisionCandidates->begin(); cit != cit_end; ++ cit, ++ bit)
1069	{
1070	Intersectable obj = (cit).mObject;
1071
1072	const float renderCost = mViewCellsManager->EvalRenderCost(obj);
1073
1074	objectsLeft += renderCost;
1075	objectsRight -= renderCost;
1076
1077	const AxisAlignedBox3 obox = obj->GetBox();
1078
1079	AxisAlignedBox3 lbox = nodeBbox;
1080	AxisAlignedBox3 rbox = nodeBbox;
1081
1082	// the borders of the left bounding box have changed
1083	for (int i = 0; i < 3; ++ i)
1084	{
1085	if (obox.Max(i) > maxBorder[i])
1086	{
1087	maxBorder[i] = obox.Max(i);
1088	}
1089	}
1090
1091	minBorder = (*bit);
1092
1093	lbox.SetMax(maxBorder);
1094	rbox.SetMin(minBorder);
1095
1096	al = lbox.SurfaceArea();
1097	ar = rbox.SurfaceArea();
1098
1099	const bool noValidSplit = ((objectsLeft <= Limits::Small) \|\| (objectsRight <= Limits::Small));
1100	const float sum = noValidSplit ? 1e25 : objectsLeft * al + objectsRight * ar;
1101
1102	/cout << "pos=" << (cit).mPos << "\t q=(" << objectsLeft << "," << objectsRight <<")\t r=("
1103	<< lbox.SurfaceArea() << "," << rbox.SurfaceArea() << ")" << endl;
1104	cout << "minborder: " << minBorder << " maxborder: " << maxBorder << endl;
1105	cout << "cost= " << sum << endl;*/
1106
1107	if (sum < minSum)
1108	{
1109	minSum = sum;
1110	areaLeft = al;
1111	areaRight = ar;
1112
1113	// objects belong to left side now
1114	for (; currentPos != (cit + 1); ++ currentPos);
1115	}
1116	}
1117
1118	/////////////
1119	//-- assign object to front and back volume
1120
1121	// belongs to back bv
1122	for (cit = mSubdivisionCandidates->begin(); cit != currentPos; ++ cit)
1123	objectsBack.push_back((*cit).mObject);
1124
1125	// belongs to front bv
1126	for (cit = currentPos; cit != cit_end; ++ cit)
1127	objectsFront.push_back((*cit).mObject);
1128
1129	float newCost = minSum / boxArea;
1130	float ratio = newCost / totalRenderCost;
1131
1132	#ifdef GTP_DEBUG
1133	cout << "\n\nobjects=(" << (int)objectsBack.size() << "," << (int)objectsFront.size() << " of "
1134	<< (int)tData.mNode->mObjects.size() << ")\t area=("
1135	<< areaLeft << ", " << areaRight << ", " << boxArea << ")" << endl
1136	<< "cost= " << newCost << " oldCost=" << totalRenderCost / boxArea << endl;
1137	#endif
1138
1139	return ratio;
1140	}
1141
1142	#endif
1143
1144	static bool PrepareOutput(const int axis,
1145	const int leaves,
1146	ofstream &sumStats,
1147	ofstream &vollStats,
1148	ofstream &volrStats)
1149	{
1150	if ((axis == 0) && (leaves > 0) && (leaves < 90))
1151	{
1152	char str[64];
1153	sprintf(str, "tmp/bvh_heur_sum-%04d.log", leaves);
1154	sumStats.open(str);
1155	sprintf(str, "tmp/bvh_heur_voll-%04d.log", leaves);
1156	vollStats.open(str);
1157	sprintf(str, "tmp/bvh_heur_volr-%04d.log", leaves);
1158	volrStats.open(str);
1159	}
1160
1161	return sumStats.is_open() && vollStats.is_open() && volrStats.is_open();
1162	}
1163
1164
1165	static void PrintHeuristics(const float objectsRight,
1166	const float sum,
1167	const float volLeft,
1168	const float volRight,
1169	const float viewSpaceVol,
1170	ofstream &sumStats,
1171	ofstream &vollStats,
1172	ofstream &volrStats)
1173	{
1174	sumStats
1175	<< "#Position\n" << objectsRight << endl
1176	<< "#Sum\n" << sum / viewSpaceVol << endl
1177	<< "#Vol\n" << (volLeft + volRight) / viewSpaceVol << endl;
1178
1179	vollStats
1180	<< "#Position\n" << objectsRight << endl
1181	<< "#Vol\n" << volLeft / viewSpaceVol << endl;
1182
1183	volrStats
1184	<< "#Position\n" << objectsRight << endl
1185	<< "#Vol\n" << volRight / viewSpaceVol << endl;
1186	}
1187
1188
1189	float BvHierarchy::EvalLocalCostHeuristics(const BvhTraversalData &tData,
1190	const int axis,
1191	ObjectContainer &objectsFront,
1192	ObjectContainer &objectsBack)
1193	{
1194	/////////////////////////////////////////////
1195	//-- go through the lists, count the number of objects
1196	//-- left and right and evaluate the cost funcion
1197
1198	// prepare the heuristics by setting mailboxes and counters
1199	const float totalVol = PrepareHeuristics(tData, axis);
1200
1201	// local helper variables
1202	float volLeft = 0;
1203	float volRight = totalVol;
1204
1205	const float nTotalObjects = EvalAbsCost(tData.mNode->mObjects);
1206	float nObjectsLeft = 0;
1207	float nObjectsRight = nTotalObjects;
1208
1209	const float viewSpaceVol =
1210	mViewCellsManager->GetViewSpaceBox().GetVolume();
1211
1212	SortableEntryContainer::const_iterator backObjectsStart =
1213	mSubdivisionCandidates->begin();
1214
1215	/////////////////////////////////
1216	//-- the parameters for the current optimum
1217
1218	float volBack = volLeft;
1219	float volFront = volRight;
1220	float newRenderCost = nTotalObjects * totalVol;
1221
1222	#ifdef GTP_DEBUG
1223	ofstream sumStats;
1224	ofstream vollStats;
1225	ofstream volrStats;
1226
1227	const bool printStats = PrepareOutput(axis,
1228	mBvhStats.Leaves(),
1229	sumStats,
1230	vollStats,
1231	volrStats);
1232	#endif
1233
1234	///////////////////////
1235	//-- the sweep heuristics
1236	//-- traverse through events and find best split plane
1237
1238	SortableEntryContainer::const_iterator cit,
1239	cit_end = cit_end = mSubdivisionCandidates->end();
1240
1241	for (cit = mSubdivisionCandidates->begin(); cit != cit_end; ++ cit)
1242	{
1243	Intersectable object = (cit).mObject;
1244
1245	// evaluate change in l and r volume
1246	// voll = view cells that see only left node (i.e., left pvs)
1247	// volr = view cells that see only right node (i.e., right pvs)
1248	EvalHeuristicsContribution(object, volLeft, volRight);
1249
1250	const float rc = mViewCellsManager->EvalRenderCost(object);
1251
1252	nObjectsLeft += rc;
1253	nObjectsRight -= rc;
1254
1255	// split is only valid if #objects on left and right is not zero
1256	const bool noValidSplit = ((nObjectsLeft <= Limits::Small) \|\|
1257	(nObjectsRight <= Limits::Small));
1258
1259	// the heuristics
1260	const float sum = noValidSplit ?
1261	1e25 : volLeft * (float)nObjectsLeft + volRight * (float)nObjectsRight;
1262
1263	#ifdef GTP_DEBUG
1264	if (printStats)
1265	{
1266	PrintHeuristics(nObjectsRight, sum, volLeft, volRight, viewSpaceVol,
1267	sumStats, vollStats, volrStats);
1268	}
1269	#endif
1270
1271	if (sum < newRenderCost)
1272	{
1273	newRenderCost = sum;
1274
1275	volBack = volLeft;
1276	volFront = volRight;
1277
1278	// objects belongs to left side now
1279	for (; backObjectsStart != (cit + 1); ++ backObjectsStart);
1280	}
1281	}
1282
1283	////////////////////////////////////////
1284	//-- assign object to front and back volume
1285
1286	// belongs to back bv
1287	for (cit = mSubdivisionCandidates->begin(); cit != backObjectsStart; ++ cit)
1288	{
1289	objectsBack.push_back((*cit).mObject);
1290	}
1291	// belongs to front bv
1292	for (cit = backObjectsStart; cit != cit_end; ++ cit)
1293	{
1294	objectsFront.push_back((*cit).mObject);
1295	}
1296
1297	// render cost of the old parent
1298	const float oldRenderCost = (float)nTotalObjects * totalVol + Limits::Small;
1299	// the relative cost ratio
1300	const float ratio = newRenderCost / oldRenderCost;
1301
1302	#ifdef GTP_DEBUG
1303	Debug << "\n§§§§ bvh eval const decrease §§§§" << endl
1304	<< "back pvs: " << (int)objectsBack.size() << " front pvs: "
1305	<< (int)objectsFront.size() << " total pvs: " << nTotalObjects << endl
1306	<< "back p: " << volBack / viewSpaceVol << " front p "
1307	<< volFront / viewSpaceVol << " p: " << totalVol / viewSpaceVol << endl
1308	<< "old rc: " << oldRenderCost / viewSpaceVol << " new rc: "
1309	<< newRenderCost / viewSpaceVol << endl
1310	<< "render cost decrease: "
1311	<< oldRenderCost / viewSpaceVol - newRenderCost / viewSpaceVol << endl;
1312	#endif
1313
1314	return ratio;
1315	}
1316
1317
1318	void BvHierarchy::PrepareLocalSubdivisionCandidates(const BvhTraversalData &tData,
1319	const int axis)
1320	{
1321	//-- insert object queries
1322	ObjectContainer *objects = mUseGlobalSorting ?
1323	tData.mSortedObjects[axis] : &tData.mNode->mObjects;
1324
1325	CreateLocalSubdivisionCandidates(*objects, &mSubdivisionCandidates, !mUseGlobalSorting, axis);
1326	}
1327
1328
1329	void BvHierarchy::CreateLocalSubdivisionCandidates(const ObjectContainer &objects,
1330	SortableEntryContainer **subdivisionCandidates,
1331	const bool sort,
1332	const int axis)
1333	{
1334	(*subdivisionCandidates)->clear();
1335
1336	// compute requested size and look if subdivision candidate has to be recomputed
1337	const int requestedSize = (int)objects.size() * 2;
1338
1339	// creates a sorted split candidates array
1340	if ((*subdivisionCandidates)->capacity() > 500000 &&
1341	requestedSize < (int)((*subdivisionCandidates)->capacity() / 10) )
1342	{
1343	delete (*subdivisionCandidates);
1344	(*subdivisionCandidates) = new SortableEntryContainer;
1345	}
1346
1347	(*subdivisionCandidates)->reserve(requestedSize);
1348
1349	ObjectContainer::const_iterator oit, oit_end = objects.end();
1350
1351	for (oit = objects.begin(); oit < oit_end; ++ oit)
1352	{
1353	Intersectable object = oit;
1354	const AxisAlignedBox3 &box = object->GetBox();
1355	const float midPt = (box.Min(axis) + box.Max(axis)) * 0.5f;
1356
1357	(*subdivisionCandidates)->push_back(SortableEntry(object, midPt));
1358	}
1359
1360	if (sort)
1361	{ // no presorted candidate list
1362	stable_sort((subdivisionCandidates)->begin(), (subdivisionCandidates)->end());
1363	}
1364	}
1365
1366
1367	const BvhStatistics &BvHierarchy::GetStatistics() const
1368	{
1369	return mBvhStats;
1370	}
1371
1372
1373	float BvHierarchy::PrepareHeuristics(const BvhTraversalData &tData,
1374	const int axis)
1375	{
1376	BvhLeaf *leaf = tData.mNode;
1377	float vol = 0;
1378
1379	// sort so we can use a sweep from right to left
1380	PrepareLocalSubdivisionCandidates(tData, axis);
1381
1382	// collect and mark the view cells as belonging to front pvs
1383	ViewCellContainer viewCells;
1384
1385	const int numRays = CollectViewCells(tData.mNode->mObjects, viewCells, true, true);
1386	//cout << "number of rays: " << numRays << endl;
1387
1388	ViewCellContainer::const_iterator vit, vit_end = viewCells.end();
1389	for (vit = viewCells.begin(); vit != vit_end; ++ vit)
1390	{
1391	#if USE_VOLUMES_FOR_HEURISTICS
1392	const float volIncr = (*vit)->GetVolume();
1393	#else
1394	const float volIncr = 1.0f;
1395	#endif
1396	vol += volIncr;
1397	}
1398
1399	// we will mail view cells switching to the back side
1400	ViewCell::NewMail();
1401
1402	return vol;
1403	}
1404
1405	///////////////////////////////////////////////////////////
1406
1407
1408	void BvHierarchy::EvalHeuristicsContribution(Intersectable *obj,
1409	float &volLeft,
1410	float &volRight)
1411	{
1412	// collect all view cells associated with this objects
1413	// (also multiple times, if they are pierced by several rays)
1414	ViewCellContainer viewCells;
1415	const bool useMailboxing = false;
1416
1417	CollectViewCells(obj, viewCells, useMailboxing, false, true);
1418
1419	// classify view cells and compute volume contri accordingly
1420	// possible view cell classifications:
1421	// view cell mailed => view cell can be seen from left child node
1422	// view cell counter > 0 view cell can be seen from right child node
1423	// combined: view cell volume belongs to both nodes
1424	ViewCellContainer::const_iterator vit, vit_end = viewCells.end();
1425
1426	for (vit = viewCells.begin(); vit != vit_end; ++ vit)
1427	{
1428	// view cells can also be seen from left child node
1429	ViewCell viewCell = vit;
1430	#if USE_VOLUMES_FOR_HEURISTICS
1431	const float vol = viewCell->GetVolume();
1432	#else
1433	const float vol = 1.0f;
1434	#endif
1435	if (!viewCell->Mailed())
1436	{
1437	viewCell->Mail();
1438	// we now see view cell from both nodes
1439	// => add volume to left node
1440	volLeft += vol;
1441	}
1442
1443	// last reference into the right node
1444	if (-- viewCell->mCounter == 0)
1445	{
1446	// view cell was previously seen from both nodes =>
1447	// remove volume from right node
1448	volRight -= vol;
1449	}
1450	}
1451	}
1452
1453
1454	void BvHierarchy::SetViewCellsManager(ViewCellsManager *vcm)
1455	{
1456	mViewCellsManager = vcm;
1457	}
1458
1459
1460	AxisAlignedBox3 BvHierarchy::GetBoundingBox() const
1461	{
1462	return mBoundingBox;
1463	}
1464
1465
1466	float BvHierarchy::SelectObjectPartition(const BvhTraversalData &tData,
1467	ObjectContainer &frontObjects,
1468	ObjectContainer &backObjects,
1469	bool useVisibilityBasedHeuristics)
1470	{
1471	if (mIsInitialSubdivision)
1472	{
1473	ApplyInitialSplit(tData, frontObjects, backObjects);
1474	return 0;
1475	}
1476
1477	ObjectContainer nFrontObjects[3];
1478	ObjectContainer nBackObjects[3];
1479	float nCostRatio[3];
1480
1481	int sAxis = 0;
1482	int bestAxis = -1;
1483
1484	if (mOnlyDrivingAxis)
1485	{
1486	const AxisAlignedBox3 box = tData.mNode->GetBoundingBox();
1487	sAxis = box.Size().DrivingAxis();
1488	}
1489
1490	// only use a subset of the rays for visibility based heuristics
1491	if (mUseCostHeuristics && useVisibilityBasedHeuristics)
1492	{
1493	VssRayContainer rays;
1494	// maximal 2 objects share the same ray
1495	rays.reserve(tData.mNumRays * 2);
1496	CollectRays(tData.mNode->mObjects, rays);
1497
1498	const float prop = (float)mMaxTests / (float)tData.mNumRays;
1499
1500	VssRay::NewMail();
1501
1502	VssRayContainer::const_iterator rit, rit_end = rays.end();
1503
1504	int nRays = 0;
1505
1506	for (rit = rays.begin(); rit != rit_end; ++ rit)
1507	{
1508	if ((mMaxTests >= (int)rays.size()) \|\| (Random(1.0f) < prop))
1509	{
1510	(*rit)->Mail();
1511	++ nRays;
1512	}
1513	}
1514	}
1515
1516	////////////////////////////////////
1517	//-- evaluate split cost for all three axis
1518
1519	for (int axis = 0; axis < 3; ++ axis)
1520	{
1521	if (!mOnlyDrivingAxis \|\| (axis == sAxis))
1522	{
1523	if (mUseCostHeuristics)
1524	{
1525	//////////////////////////////////
1526	//-- split objects using heuristics
1527
1528	if (useVisibilityBasedHeuristics)
1529	{
1530	///////////
1531	//-- heuristics using objects weighted by view cells volume
1532	nCostRatio[axis] =
1533	EvalLocalCostHeuristics(tData,
1534	axis,
1535	nFrontObjects[axis],
1536	nBackObjects[axis]);
1537	}
1538	else
1539	{
1540	//////////////////
1541	//-- view cells not constructed yet => use surface area heuristic
1542	nCostRatio[axis] = EvalSah(tData,
1543	axis,
1544	nFrontObjects[axis],
1545	nBackObjects[axis]);
1546	}
1547	}
1548	else
1549	{
1550	//-- split objects using some simple criteria
1551	nCostRatio[axis] =
1552	EvalLocalObjectPartition(tData, axis, nFrontObjects[axis], nBackObjects[axis]);
1553	}
1554
1555	// no good results for degenerate axis split
1556	if (1 &&
1557	(tData.mNode->GetBoundingBox().Size(axis) < 0.0001))//Limits::Small))
1558	{
1559	nCostRatio[axis] += 9999;
1560	}
1561
1562	if ((bestAxis == -1) \|\| (nCostRatio[axis] < nCostRatio[bestAxis]))
1563	{
1564	bestAxis = axis;
1565	}
1566	}
1567	}
1568
1569	////////////////
1570	//-- assign values
1571
1572	frontObjects = nFrontObjects[bestAxis];
1573	backObjects = nBackObjects[bestAxis];
1574
1575	//cout << "val: " << nCostRatio[bestAxis] << " axis: " << bestAxis << endl;
1576	return nCostRatio[bestAxis];
1577	}
1578
1579
1580	int BvHierarchy::AssociateObjectsWithRays(const VssRayContainer &rays) const
1581	{
1582	int nRays = 0;
1583	VssRayContainer::const_iterator rit, rit_end = rays.end();
1584
1585	VssRay::NewMail();
1586
1587	for (rit = rays.begin(); rit != rays.end(); ++ rit)
1588	{
1589	VssRay ray = (rit);
1590
1591	if (ray->mTerminationObject)
1592	{
1593	ray->mTerminationObject->GetOrCreateRays()->push_back(ray);
1594	if (!ray->Mailed())
1595	{
1596	ray->Mail();
1597	++ nRays;
1598	}
1599	}
1600
1601	#if COUNT_ORIGIN_OBJECTS
1602
1603	if (ray->mOriginObject)
1604	{
1605	ray->mOriginObject->GetOrCreateRays()->push_back(ray);
1606
1607	if (!ray->Mailed())
1608	{
1609	ray->Mail();
1610	++ nRays;
1611	}
1612	}
1613	#endif
1614	}
1615
1616	return nRays;
1617	}
1618
1619
1620	void BvHierarchy::PrintSubdivisionStats(const SubdivisionCandidate &sc)
1621	{
1622	const float costDecr = sc.GetRenderCostDecrease();
1623
1624	mSubdivisionStats
1625	<< "#Leaves\n" << mBvhStats.Leaves() << endl
1626	<< "#RenderCostDecrease\n" << costDecr << endl
1627	<< "#TotalRenderCost\n" << mTotalCost << endl
1628	<< "#EntriesInPvs\n" << mPvsEntries << endl;
1629	}
1630
1631
1632	void BvHierarchy::CollectRays(const ObjectContainer &objects,
1633	VssRayContainer &rays) const
1634	{
1635	VssRay::NewMail();
1636	ObjectContainer::const_iterator oit, oit_end = objects.end();
1637
1638	// evaluate reverse pvs and view cell volume on left and right cell
1639	// note: should I take all leaf objects or rather the objects hit by rays?
1640	for (oit = objects.begin(); oit != oit_end; ++ oit)
1641	{
1642	Intersectable obj = oit;
1643	VssRayContainer::const_iterator rit, rit_end = obj->GetOrCreateRays()->end();
1644
1645	for (rit = obj->GetOrCreateRays()->begin(); rit < rit_end; ++ rit)
1646	{
1647	VssRay ray = (rit);
1648
1649	if (!ray->Mailed())
1650	{
1651	ray->Mail();
1652	rays.push_back(ray);
1653	}
1654	}
1655	}
1656	}
1657
1658
1659	float BvHierarchy::EvalAbsCost(const ObjectContainer &objects)// const
1660	{
1661	#if USE_BETTER_RENDERCOST_EST
1662	ObjectContainer::const_iterator oit, oit_end = objects.end();
1663
1664	for (oit = objects.begin(); oit != oit_end; ++ oit)
1665	{
1666	objRenderCost += ViewCellsManager::GetRendercost(*oit);
1667	}
1668	#else
1669	return (float)objects.size();
1670	#endif
1671	}
1672
1673
1674	float BvHierarchy::EvalSahCost(BvhLeaf *leaf) const
1675	{
1676	////////////////
1677	//-- surface area heuristics
1678	if (leaf->mObjects.empty())
1679	return 0.0f;
1680
1681	const AxisAlignedBox3 box = GetBoundingBox(leaf);
1682	const float area = box.SurfaceArea();
1683	const float viewSpaceArea = mViewCellsManager->GetViewSpaceBox().SurfaceArea();
1684
1685	return EvalAbsCost(leaf->mObjects) * area / viewSpaceArea;
1686	}
1687
1688
1689	float BvHierarchy::EvalRenderCost(const ObjectContainer &objects) const
1690	{
1691	///////////////
1692	//-- render cost heuristics
1693
1694	const float viewSpaceVol = mViewCellsManager->GetViewSpaceBox().GetVolume();
1695
1696	// probability that view point lies in a view cell which sees this node
1697	const float p = EvalViewCellsVolume(objects) / viewSpaceVol;
1698	const float objRenderCost = EvalAbsCost(objects);
1699
1700	return objRenderCost * p;
1701	}
1702
1703
1704	AxisAlignedBox3 BvHierarchy::EvalBoundingBox(const ObjectContainer &objects,
1705	const AxisAlignedBox3 *parentBox) const
1706	{
1707	// if there are no objects in this box, box size is set to parent box size.
1708	// Question: Invalidate box instead?
1709	if (parentBox && objects.empty())
1710	return *parentBox;
1711
1712	AxisAlignedBox3 box;
1713	box.Initialize();
1714
1715	ObjectContainer::const_iterator oit, oit_end = objects.end();
1716
1717	for (oit = objects.begin(); oit != oit_end; ++ oit)
1718	{
1719	Intersectable obj = oit;
1720	// grow bounding box to include all objects
1721	box.Include(obj->GetBox());
1722	}
1723
1724	return box;
1725	}
1726
1727
1728	void BvHierarchy::CollectLeaves(BvhNode root, vector<BvhLeaf > &leaves) const
1729	{
1730	stack<BvhNode *> nodeStack;
1731	nodeStack.push(root);
1732
1733	while (!nodeStack.empty())
1734	{
1735	BvhNode *node = nodeStack.top();
1736	nodeStack.pop();
1737
1738	if (node->IsLeaf())
1739	{
1740	BvhLeaf leaf = (BvhLeaf )node;
1741	leaves.push_back(leaf);
1742	}
1743	else
1744	{
1745	BvhInterior interior = (BvhInterior )node;
1746
1747	nodeStack.push(interior->GetBack());
1748	nodeStack.push(interior->GetFront());
1749	}
1750	}
1751	}
1752
1753
1754	AxisAlignedBox3 BvHierarchy::GetBoundingBox(BvhNode *node) const
1755	{
1756	return node->GetBoundingBox();
1757	}
1758
1759
1760	int BvHierarchy::CollectViewCells(const ObjectContainer &objects,
1761	ViewCellContainer &viewCells,
1762	const bool setCounter,
1763	const bool onlyMailedRays) const
1764	{
1765	ViewCell::NewMail();
1766	ObjectContainer::const_iterator oit, oit_end = objects.end();
1767
1768	int numRays = 0;
1769	// loop through all object and collect view cell pvs of this node
1770	for (oit = objects.begin(); oit != oit_end; ++ oit)
1771	{
1772	// always use only mailed objects
1773	numRays += CollectViewCells(*oit, viewCells, true, setCounter, onlyMailedRays);
1774	}
1775
1776	return numRays;
1777	}
1778
1779
1780	int BvHierarchy::CollectViewCells(Intersectable *obj,
1781	ViewCellContainer &viewCells,
1782	const bool useMailBoxing,
1783	const bool setCounter,
1784	const bool onlyMailedRays) const
1785	{
1786	VssRayContainer::const_iterator rit, rit_end = obj->GetOrCreateRays()->end();
1787
1788	int numRays = 0;
1789
1790	for (rit = obj->GetOrCreateRays()->begin(); rit < rit_end; ++ rit)
1791	{
1792	VssRay ray = (rit);
1793
1794	if (onlyMailedRays && !ray->Mailed())
1795	{//if (onlyMailedRays)cout << "u";
1796	continue;
1797	}//else if (onlyMailedRays) cout << "z";
1798
1799	//ray->Mail();
1800	++ numRays;
1801
1802	ViewCellContainer tmpViewCells;
1803	mHierarchyManager->mVspTree->GetViewCells(*ray, tmpViewCells);
1804
1805	// matt: probably slow to allocate memory for view cells every time
1806	ViewCellContainer::const_iterator vit, vit_end = tmpViewCells.end();
1807
1808	for (vit = tmpViewCells.begin(); vit != vit_end; ++ vit)
1809	{
1810	ViewCell vc = vit;
1811
1812	// store view cells
1813	if (!useMailBoxing \|\| !vc->Mailed())
1814	{
1815	if (useMailBoxing)
1816	{
1817	vc->Mail();
1818	if (setCounter)
1819	{
1820	vc->mCounter = 0;
1821	}
1822	}
1823	viewCells.push_back(vc);
1824	}
1825
1826	if (setCounter)
1827	{
1828	++ vc->mCounter;
1829	}
1830	}
1831	}
1832
1833	return numRays;
1834	}
1835
1836
1837	int BvHierarchy::CountViewCells(Intersectable *obj) const
1838	{
1839	int result = 0;
1840
1841	VssRayContainer::const_iterator rit, rit_end = obj->GetOrCreateRays()->end();
1842
1843	for (rit = obj->GetOrCreateRays()->begin(); rit < rit_end; ++ rit)
1844	{
1845	VssRay ray = (rit);
1846	ViewCellContainer tmpViewCells;
1847
1848	mHierarchyManager->mVspTree->GetViewCells(*ray, tmpViewCells);
1849
1850	ViewCellContainer::const_iterator vit, vit_end = tmpViewCells.end();
1851	for (vit = tmpViewCells.begin(); vit != vit_end; ++ vit)
1852	{
1853	ViewCell vc = vit;
1854
1855	// store view cells
1856	if (!vc->Mailed())
1857	{
1858	vc->Mail();
1859	++ result;
1860	}
1861	}
1862	}
1863
1864	return result;
1865	}
1866
1867
1868	int BvHierarchy::CountViewCells(const ObjectContainer &objects) const
1869	{
1870	int nViewCells = 0;
1871	ViewCell::NewMail();
1872	ObjectContainer::const_iterator oit, oit_end = objects.end();
1873
1874	// loop through all object and collect view cell pvs of this node
1875	for (oit = objects.begin(); oit != oit_end; ++ oit)
1876	{
1877	nViewCells += CountViewCells(*oit);
1878	}
1879
1880	return nViewCells;
1881	}
1882
1883
1884	void BvHierarchy::CollectDirtyCandidates(BvhSubdivisionCandidate *sc,
1885	vector<SubdivisionCandidate *> &dirtyList,
1886	const bool onlyUnmailed)
1887	{
1888	BvhTraversalData &tData = sc->mParentData;
1889	BvhLeaf *node = tData.mNode;
1890
1891	ViewCellContainer viewCells;
1892	ViewCell::NewMail();
1893	int numRays = CollectViewCells(node->mObjects, viewCells, false, false);
1894
1895	if (0) cout << "collected " << (int)viewCells.size() << " dirty candidates" << endl;
1896
1897	// split candidates handling
1898	// these view cells are thrown into dirty list
1899	ViewCellContainer::const_iterator vit, vit_end = viewCells.end();
1900
1901	for (vit = viewCells.begin(); vit != vit_end; ++ vit)
1902	{
1903	VspViewCell vc = dynamic_cast<VspViewCell >(*vit);
1904	VspLeaf *leaf = vc->mLeaves[0];
1905
1906	SubdivisionCandidate *candidate = leaf->GetSubdivisionCandidate();
1907
1908	// is this leaf still a split candidate?
1909	if (candidate && (!onlyUnmailed \|\| !candidate->Mailed()))
1910	{
1911	candidate->Mail();
1912	candidate->SetDirty(true);
1913	dirtyList.push_back(candidate);
1914	}
1915	}
1916	}
1917
1918
1919	BvhNode *BvHierarchy::GetRoot() const
1920	{
1921	return mRoot;
1922	}
1923
1924
1925	bool BvHierarchy::IsObjectInLeaf(BvhLeaf leaf, Intersectable object) const
1926	{
1927	ObjectContainer::const_iterator oit =
1928	lower_bound(leaf->mObjects.begin(), leaf->mObjects.end(), object, ilt);
1929
1930	// objects sorted by id
1931	if ((oit != leaf->mObjects.end()) && ((*oit)->GetId() == object->GetId()))
1932	{
1933	return true;
1934	}
1935	else
1936	{
1937	return false;
1938	}
1939	}
1940
1941
1942	BvhLeaf BvHierarchy::GetLeaf(Intersectable object, BvhNode *node) const
1943	{
1944	// rather use the simple version
1945	if (!object)
1946	return NULL;
1947	return object->mBvhLeaf;
1948
1949	///////////////////////////////////////
1950	// start from root of tree
1951
1952	if (node == NULL)
1953	node = mRoot;
1954
1955	vector<BvhLeaf *> leaves;
1956
1957	stack<BvhNode *> nodeStack;
1958	nodeStack.push(node);
1959
1960	BvhLeaf *leaf = NULL;
1961
1962	while (!nodeStack.empty())
1963	{
1964	BvhNode *node = nodeStack.top();
1965	nodeStack.pop();
1966
1967	if (node->IsLeaf())
1968	{
1969	leaf = dynamic_cast<BvhLeaf *>(node);
1970
1971	if (IsObjectInLeaf(leaf, object))
1972	{
1973	return leaf;
1974	}
1975	}
1976	else
1977	{
1978	// find point
1979	BvhInterior interior = dynamic_cast<BvhInterior >(node);
1980
1981	if (interior->GetBack()->GetBoundingBox().Includes(object->GetBox()))
1982	{
1983	nodeStack.push(interior->GetBack());
1984	}
1985
1986	// search both sides as we are using bounding volumes
1987	if (interior->GetFront()->GetBoundingBox().Includes(object->GetBox()))
1988	{
1989	nodeStack.push(interior->GetFront());
1990	}
1991	}
1992	}
1993
1994	return leaf;
1995	}
1996
1997
1998	bool BvHierarchy::Export(OUT_STREAM &stream)
1999	{
2000	ExportNode(mRoot, stream);
2001
2002	return true;
2003	}
2004
2005
2006	void BvHierarchy::ExportObjects(BvhLeaf *leaf, OUT_STREAM &stream)
2007	{
2008	ObjectContainer::const_iterator oit, oit_end = leaf->mObjects.end();
2009	for (oit = leaf->mObjects.begin(); oit != oit_end; ++ oit)
2010	{
2011	stream << (*oit)->GetId() << " ";
2012	}
2013	}
2014
2015
2016	void BvHierarchy::ExportNode(BvhNode *node, OUT_STREAM &stream)
2017	{
2018	if (node->IsLeaf())
2019	{
2020	BvhLeaf leaf = dynamic_cast<BvhLeaf >(node);
2021	const AxisAlignedBox3 box = leaf->GetBoundingBox();
2022	stream << "<Leaf"
2023	<< " min=\"" << box.Min().x << " " << box.Min().y << " " << box.Min().z << "\""
2024	<< " max=\"" << box.Max().x << " " << box.Max().y << " " << box.Max().z << "\""
2025	<< " objects=\"";
2026
2027	//-- export objects
2028	ExportObjects(leaf, stream);
2029
2030	stream << "\" />" << endl;
2031	}
2032	else
2033	{
2034	BvhInterior interior = dynamic_cast<BvhInterior >(node);
2035	const AxisAlignedBox3 box = interior->GetBoundingBox();
2036
2037	stream << "<Interior"
2038	<< " min=\"" << box.Min().x << " " << box.Min().y << " " << box.Min().z << "\""
2039	<< " max=\"" << box.Max().x << " " << box.Max().y << " " << box.Max().z
2040	<< "\">" << endl;
2041
2042	ExportNode(interior->GetBack(), stream);
2043	ExportNode(interior->GetFront(), stream);
2044
2045	stream << "</Interior>" << endl;
2046	}
2047	}
2048
2049
2050	float BvHierarchy::EvalViewCellsVolume(const ObjectContainer &objects) const
2051	{
2052	float vol = 0;
2053
2054	ViewCellContainer viewCells;
2055
2056	// we have to account for all view cells that can
2057	// be seen from the objects
2058	int numRays = CollectViewCells(objects, viewCells, false, false);
2059
2060	ViewCellContainer::const_iterator vit, vit_end = viewCells.end();
2061
2062	for (vit = viewCells.begin(); vit != vit_end; ++ vit)
2063	{
2064	vol += (*vit)->GetVolume();
2065	}
2066
2067	return vol;
2068	}
2069
2070
2071	void BvHierarchy::Initialise(const ObjectContainer &objects)
2072	{
2073	AxisAlignedBox3 box = EvalBoundingBox(objects);
2074
2075	///////
2076	//-- create new root
2077
2078	BvhLeaf *bvhleaf = new BvhLeaf(box, NULL, (int)objects.size());
2079	bvhleaf->mObjects = objects;
2080	mRoot = bvhleaf;
2081
2082	// compute bounding box from objects
2083	mBoundingBox = mRoot->GetBoundingBox();
2084
2085	// associate root with current objects
2086	AssociateObjectsWithLeaf(bvhleaf);
2087	}
2088
2089
2090	/*
2091	Mesh *BvHierarchy::MergeLeafToMesh()
2092	{
2093	vector<BvhLeaf *> leaves;
2094	CollectLeaves(leaves);
2095
2096	vector<BvhLeaf *>::const_iterator lit, lit_end = leaves.end();
2097
2098	for (lit = leaves.begin(); lit != lit_end; ++ lit)
2099	{
2100	Mesh mesh = MergeLeafToMesh(lit);
2101	}
2102	}*/
2103
2104
2105	void BvHierarchy::PrepareConstruction(SplitQueue &tQueue,
2106	const VssRayContainer &sampleRays,
2107	const ObjectContainer &objects)
2108	{
2109	///////////////////////////////////////
2110	//-- we assume that we have objects sorted by their id =>
2111	//-- we don't have to sort them here and an binary search
2112	//-- for identifying if a object is in a leaf.
2113
2114	mBvhStats.Reset();
2115	mBvhStats.Start();
2116	mBvhStats.nodes = 1;
2117
2118	// store pointer to this tree
2119	BvhSubdivisionCandidate::sBvHierarchy = this;
2120
2121	// root and bounding box was already constructed
2122	BvhLeaf bvhLeaf = dynamic_cast<BvhLeaf >(mRoot);
2123
2124	// multiply termination criterium for comparison,
2125	// so it can be set between zero and one and
2126	// no division is necessary during traversal
2127
2128	#if PROBABILIY_IS_BV_VOLUME
2129	mTermMinProbability *= mBoundingBox.GetVolume();
2130	// probability that bounding volume is seen
2131	const float prop = GetBoundingBox().GetVolume();
2132	#else
2133	mTermMinProbability *= mVspTree->GetBoundingBox().GetVolume();
2134	// probability that volume is "seen" from the view cells
2135	const float prop = EvalViewCellsVolume(objects);
2136	#endif
2137
2138	// only rays intersecting objects in node are interesting
2139	const int nRays = AssociateObjectsWithRays(sampleRays);
2140	//cout << "using " << nRays << " of " << (int)sampleRays.size() << " rays" << endl;
2141
2142	// create bvh traversal data
2143	BvhTraversalData oData(bvhLeaf, 0, prop, nRays);
2144
2145	// create sorted object lists for the first data
2146	if (mUseGlobalSorting)
2147	{
2148	AssignInitialSortedObjectList(oData, objects);
2149	}
2150
2151
2152	///////////////////
2153	//-- add first candidate for object space partition
2154
2155	BvhSubdivisionCandidate *oSubdivisionCandidate =
2156	new BvhSubdivisionCandidate(oData);
2157
2158	bvhLeaf->SetSubdivisionCandidate(oSubdivisionCandidate);
2159
2160	mTotalCost = EvalRenderCost(objects);
2161	mPvsEntries = CountViewCells(objects);
2162
2163	if (mApplyInitialPartition)
2164	{
2165	vector<SubdivisionCandidate *> candidateContainer;
2166
2167	mIsInitialSubdivision = true;
2168
2169	// evaluate priority
2170	EvalSubdivisionCandidate(*oSubdivisionCandidate);
2171	PrintSubdivisionStats(*oSubdivisionCandidate);
2172
2173	ApplyInitialSubdivision(oSubdivisionCandidate, candidateContainer);
2174
2175	mIsInitialSubdivision = false;
2176
2177	vector<SubdivisionCandidate *>::const_iterator cit, cit_end = candidateContainer.end();
2178
2179	for (cit = candidateContainer.begin(); cit != cit_end; ++ cit)
2180	{
2181	// reevaluate priority
2182	EvalSubdivisionCandidate(*oSubdivisionCandidate);
2183	tQueue.Push(*cit);
2184	}
2185	}
2186	else
2187	{
2188	// evaluate priority
2189	EvalSubdivisionCandidate(*oSubdivisionCandidate);
2190	PrintSubdivisionStats(*oSubdivisionCandidate);
2191
2192	tQueue.Push(oSubdivisionCandidate);
2193	}
2194
2195	cout << "!!size: " << GetStatistics().Leaves() << endl;
2196	}
2197
2198
2199	void BvHierarchy::AssignInitialSortedObjectList(BvhTraversalData &tData,
2200	const ObjectContainer &objects)
2201	{
2202	// we sort the objects as a preprocess so they don't have
2203	// to be sorted for each split
2204	for (int i = 0; i < 3; ++ i)
2205	{
2206	SortableEntryContainer *sortedObjects = new SortableEntryContainer();
2207
2208	CreateLocalSubdivisionCandidates(objects,
2209	&sortedObjects,
2210	true,
2211	i);
2212
2213	// copy list into traversal data list
2214	tData.mSortedObjects[i] = new ObjectContainer();
2215	tData.mSortedObjects[i]->reserve((int)objects.size());
2216
2217	SortableEntryContainer::const_iterator oit, oit_end = sortedObjects->end();
2218
2219	for (oit = sortedObjects->begin(); oit != oit_end; ++ oit)
2220	{
2221	tData.mSortedObjects[i]->push_back((*oit).mObject);
2222	}
2223
2224	delete sortedObjects;
2225	}
2226
2227	// last sorted list: by size
2228	tData.mSortedObjects[3] = new ObjectContainer();
2229	tData.mSortedObjects[3]->reserve((int)objects.size());
2230
2231	*(tData.mSortedObjects[3]) = objects;
2232	stable_sort(tData.mSortedObjects[3]->begin(), tData.mSortedObjects[3]->end(), smallerSize);
2233	}
2234
2235
2236	void BvHierarchy::AssignSortedObjects(const BvhSubdivisionCandidate &sc,
2237	BvhTraversalData &frontData,
2238	BvhTraversalData &backData)
2239	{
2240	Intersectable::NewMail();
2241
2242	// we sorted the objects as a preprocess so they don't have
2243	// to be sorted for each split
2244	ObjectContainer::const_iterator fit, fit_end = sc.mFrontObjects.end();
2245
2246	for (fit = sc.mFrontObjects.begin(); fit != fit_end; ++ fit)
2247	{
2248	(*fit)->Mail();
2249	}
2250
2251	for (int i = 0; i < 4; ++ i)
2252	{
2253	frontData.mSortedObjects[i] = new ObjectContainer();
2254	backData.mSortedObjects[i] = new ObjectContainer();
2255
2256	frontData.mSortedObjects[i]->reserve((int)sc.mFrontObjects.size());
2257	backData.mSortedObjects[i]->reserve((int)sc.mFrontObjects.size());
2258
2259	ObjectContainer::const_iterator oit, oit_end = sc.mParentData.mSortedObjects[i]->end();
2260
2261	for (oit = sc.mParentData.mSortedObjects[i]->begin(); oit != oit_end; ++ oit)
2262	{
2263	if ((*oit)->Mailed())
2264	{
2265	frontData.mSortedObjects[i]->push_back(*oit);
2266	}
2267	else
2268	{
2269	backData.mSortedObjects[i]->push_back(*oit);
2270	}
2271	}
2272	}
2273	}
2274
2275
2276	void BvHierarchy::Reset(SplitQueue &tQueue,
2277	const VssRayContainer &sampleRays,
2278	const ObjectContainer &objects)
2279	{
2280	// reset stats
2281	mBvhStats.Reset();
2282	mBvhStats.Start();
2283	mBvhStats.nodes = 1;
2284
2285	// reset root
2286	DEL_PTR(mRoot);
2287
2288	BvhLeaf *bvhleaf = new BvhLeaf(mBoundingBox, NULL, (int)objects.size());
2289	bvhleaf->mObjects = objects;
2290	mRoot = bvhleaf;
2291
2292	#if PROBABILIY_IS_BV_VOLUME
2293	mTermMinProbability *= mBoundingBox.GetVolume();
2294	// probability that bounding volume is seen
2295	const float prop = GetBoundingBox().GetVolume();
2296	#else
2297	mTermMinProbability *= mVspTree->GetBoundingBox().GetVolume();
2298	// probability that volume is "seen" from the view cells
2299	const float prop = EvalViewCellsVolume(objects);
2300	#endif
2301
2302	const int nRays = CountRays(objects);
2303	BvhLeaf bvhLeaf = dynamic_cast<BvhLeaf >(mRoot);
2304
2305	// create bvh traversal data
2306	BvhTraversalData oData(bvhLeaf, 0, prop, nRays);
2307
2308	AssignInitialSortedObjectList(oData, objects);
2309
2310
2311	///////////////////
2312	//-- add first candidate for object space partition
2313
2314	BvhSubdivisionCandidate *oSubdivisionCandidate =
2315	new BvhSubdivisionCandidate(oData);
2316
2317	EvalSubdivisionCandidate(*oSubdivisionCandidate);
2318	bvhLeaf->SetSubdivisionCandidate(oSubdivisionCandidate);
2319
2320	const float viewSpaceVol = mViewCellsManager->GetViewSpaceBox().GetVolume();
2321	mTotalCost = (float)objects.size() * prop / viewSpaceVol;
2322
2323	PrintSubdivisionStats(*oSubdivisionCandidate);
2324
2325	tQueue.Push(oSubdivisionCandidate);
2326	}
2327
2328
2329	void BvhStatistics::Print(ostream &app) const
2330	{
2331	app << "=========== BvHierarchy statistics ===============\n";
2332
2333	app << setprecision(4);
2334
2335	app << "#N_CTIME ( Construction time [s] )\n" << Time() << " \n";
2336
2337	app << "#N_NODES ( Number of nodes )\n" << nodes << "\n";
2338
2339	app << "#N_INTERIORS ( Number of interior nodes )\n" << Interior() << "\n";
2340
2341	app << "#N_LEAVES ( Number of leaves )\n" << Leaves() << "\n";
2342
2343	app << "#AXIS_ALIGNED_SPLITS (number of axis aligned splits)\n" << splits << endl;
2344
2345	app << "#N_MAXCOSTNODES ( Percentage of leaves with terminated because of max cost ratio )\n"
2346	<< maxCostNodes * 100 / (double)Leaves() << endl;
2347
2348	app << "#N_PMINPROBABILITYLEAVES ( Percentage of leaves with mininum probability )\n"
2349	<< minProbabilityNodes * 100 / (double)Leaves() << endl;
2350
2351
2352	//////////////////////////////////////////////////
2353
2354	app << "#N_PMAXDEPTHLEAVES ( Percentage of leaves at maximum depth )\n"
2355	<< maxDepthNodes * 100 / (double)Leaves() << endl;
2356
2357	app << "#N_PMAXDEPTH ( Maximal reached depth )\n" << maxDepth << endl;
2358
2359	app << "#N_PMINDEPTH ( Minimal reached depth )\n" << minDepth << endl;
2360
2361	app << "#AVGDEPTH ( average depth )\n" << AvgDepth() << endl;
2362
2363
2364	////////////////////////////////////////////////////////
2365
2366	app << "#N_PMINOBJECTSLEAVES ( Percentage of leaves with mininum objects )\n"
2367	<< minObjectsNodes * 100 / (double)Leaves() << endl;
2368
2369	app << "#N_MAXOBJECTREFS ( Max number of object refs / leaf )\n" << maxObjectRefs << "\n";
2370
2371	app << "#N_MINOBJECTREFS ( Min number of object refs / leaf )\n" << minObjectRefs << "\n";
2372
2373	app << "#N_EMPTYLEAFS ( Empty leafs )\n" << emptyNodes << "\n";
2374
2375	app << "#N_PAVGOBJECTSLEAVES ( average object refs / leaf)\n" << AvgObjectRefs() << endl;
2376
2377
2378	////////////////////////////////////////////////////////
2379
2380	app << "#N_PMINRAYSLEAVES ( Percentage of leaves with mininum rays )\n"
2381	<< minRaysNodes * 100 / (double)Leaves() << endl;
2382
2383	app << "#N_MAXRAYREFS ( Max number of ray refs / leaf )\n" << maxRayRefs << "\n";
2384
2385	app << "#N_MINRAYREFS ( Min number of ray refs / leaf )\n" << minRayRefs << "\n";
2386
2387	app << "#N_PAVGRAYLEAVES ( average ray refs / leaf )\n" << AvgRayRefs() << endl;
2388
2389	app << "#N_PAVGRAYCONTRIBLEAVES ( Average ray contribution)\n" <<
2390	rayRefs / (double)objectRefs << endl;
2391
2392	app << "#N_PMAXRAYCONTRIBLEAVES ( Percentage of leaves with maximal ray contribution )\n"<<
2393	maxRayContriNodes * 100 / (double)Leaves() << endl;
2394
2395	app << "#N_PGLOBALCOSTMISSES ( Global cost misses )\n" << mGlobalCostMisses << endl;
2396
2397	app << "========== END OF BvHierarchy statistics ==========\n";
2398	}
2399
2400
2401	// TODO: return memory usage in MB
2402	float BvHierarchy::GetMemUsage() const
2403	{
2404	return (float)(sizeof(BvHierarchy)
2405	+ mBvhStats.Leaves() * sizeof(BvhLeaf)
2406	+ mBvhStats.Interior() * sizeof(BvhInterior)
2407	) / float(1024 * 1024);
2408	}
2409
2410
2411	void BvHierarchy::SetActive(BvhNode *node) const
2412	{
2413	vector<BvhLeaf *> leaves;
2414
2415	// sets the pointers to the currently active view cells
2416	CollectLeaves(node, leaves);
2417	vector<BvhLeaf *>::const_iterator lit, lit_end = leaves.end();
2418
2419	for (lit = leaves.begin(); lit != lit_end; ++ lit)
2420	{
2421	(*lit)->SetActiveNode(node);
2422	}
2423	}
2424
2425
2426	BvhNode *BvHierarchy::SubdivideAndCopy(SplitQueue &tQueue,
2427	SubdivisionCandidate *splitCandidate)
2428	{
2429	BvhSubdivisionCandidate *sc =
2430	dynamic_cast<BvhSubdivisionCandidate *>(splitCandidate);
2431	BvhTraversalData &tData = sc->mParentData;
2432
2433	BvhNode *currentNode = tData.mNode;
2434	BvhNode oldNode = (BvhNode )splitCandidate->mEvaluationHack;
2435
2436	if (!oldNode->IsLeaf())
2437	{
2438	//////////////
2439	//-- continue subdivision
2440
2441	BvhTraversalData tFrontData;
2442	BvhTraversalData tBackData;
2443
2444	BvhInterior oldInterior = dynamic_cast<BvhInterior >(oldNode);
2445
2446	sc->mFrontObjects.clear();
2447	sc->mBackObjects.clear();
2448
2449	oldInterior->GetFront()->CollectObjects(sc->mFrontObjects);
2450	oldInterior->GetBack()->CollectObjects(sc->mBackObjects);
2451
2452	// evaluate the changes in render cost and pvs entries
2453	EvalSubdivisionCandidate(*sc, false);
2454
2455	// create new interior node and two leaf node
2456	currentNode = SubdivideNode(*sc, tFrontData, tBackData);
2457
2458	//oldNode->mRenderCostDecr += sc->GetRenderCostDecrease();
2459	//oldNode->mPvsEntriesIncr += sc->GetPvsEntriesIncr();
2460
2461	//oldNode->mRenderCostDecr = sc->GetRenderCostDecrease();
2462	//oldNode->mPvsEntriesIncr = sc->GetPvsEntriesIncr();
2463
2464	///////////////////////////
2465	//-- push the new split candidates on the queue
2466
2467	BvhSubdivisionCandidate *frontCandidate = new BvhSubdivisionCandidate(tFrontData);
2468	BvhSubdivisionCandidate *backCandidate = new BvhSubdivisionCandidate(tBackData);
2469
2470	frontCandidate->SetPriority((float)-oldInterior->GetFront()->GetTimeStamp());
2471	backCandidate->SetPriority((float)-oldInterior->GetBack()->GetTimeStamp());
2472
2473	frontCandidate->mEvaluationHack = oldInterior->GetFront();
2474	backCandidate->mEvaluationHack = oldInterior->GetBack();
2475
2476	// cross reference
2477	tFrontData.mNode->SetSubdivisionCandidate(frontCandidate);
2478	tBackData.mNode->SetSubdivisionCandidate(backCandidate);
2479
2480	//cout << "f: " << frontCandidate->GetPriority() << " b: " << backCandidate->GetPriority() << endl;
2481	tQueue.Push(frontCandidate);
2482	tQueue.Push(backCandidate);
2483	}
2484
2485	/////////////////////////////////
2486	//-- node is a leaf => terminate traversal
2487
2488	if (currentNode->IsLeaf())
2489	{
2490	// this leaf is no candidate for splitting anymore
2491	// => detach subdivision candidate
2492	tData.mNode->SetSubdivisionCandidate(NULL);
2493	// detach node so we don't delete it with the traversal data
2494	tData.mNode = NULL;
2495	}
2496
2497	return currentNode;
2498	}
2499
2500
2501	void BvHierarchy::CollectObjects(const AxisAlignedBox3 &box, ObjectContainer &objects)
2502	{
2503	stack<BvhNode *> nodeStack;
2504
2505	nodeStack.push(mRoot);
2506
2507	while (!nodeStack.empty())
2508	{
2509	BvhNode *node = nodeStack.top();
2510
2511	nodeStack.pop();
2512
2513	if (node->IsLeaf())
2514	{
2515	BvhLeaf leaf = (BvhLeaf )node;
2516	if (Overlap(box, leaf->GetBoundingBox())) {
2517	Intersectable *object = leaf;
2518	if (!object->Mailed()) {
2519	object->Mail();
2520	objects.push_back(object);
2521	}
2522	}
2523	}
2524	else
2525	{
2526	BvhInterior interior = (BvhInterior )node;
2527
2528	if (Overlap(box, interior->GetBoundingBox()))
2529	nodeStack.push(interior->GetFront());
2530
2531	if (Overlap(box, interior->GetBoundingBox()))
2532	nodeStack.push(interior->GetBack());
2533	}
2534	}
2535	}
2536
2537
2538	void BvHierarchy::ApplyInitialSubdivision(SubdivisionCandidate *firstCandidate,
2539	vector<SubdivisionCandidate *> &candidateContainer)
2540	{
2541	SplitQueue tempQueue;
2542	tempQueue.Push(firstCandidate);
2543
2544	while (!tempQueue.Empty())
2545	{
2546	SubdivisionCandidate *candidate = tempQueue.Top();
2547	tempQueue.Pop();
2548
2549	BvhSubdivisionCandidate *bsc =
2550	dynamic_cast<BvhSubdivisionCandidate *>(candidate);
2551
2552	if (!InitialTerminationCriteriaMet(bsc->mParentData))
2553	{
2554	const bool globalCriteriaMet = GlobalTerminationCriteriaMet(bsc->mParentData);
2555
2556	BvhNode *node = Subdivide(tempQueue, bsc, globalCriteriaMet);
2557
2558	// not needed anymore
2559	delete bsc;
2560	}
2561	else // initial preprocessing finished for this candidate
2562	{
2563	// add to candidate container
2564	candidateContainer.push_back(bsc);
2565	}
2566	}
2567	}
2568
2569
2570	void BvHierarchy::ApplyInitialSplit(const BvhTraversalData &tData,
2571	ObjectContainer &frontObjects,
2572	ObjectContainer &backObjects)
2573	{
2574	ObjectContainer *objects = tData.mSortedObjects[3];
2575
2576	ObjectContainer::const_iterator oit, oit_end = objects->end();
2577
2578	float maxAreaDiff = -1.0f;
2579
2580	ObjectContainer::const_iterator backObjectsStart = objects->begin();
2581
2582	for (oit = objects->begin(); oit != (objects->end() - 1); ++ oit)
2583	{
2584	Intersectable objS = oit;
2585	Intersectable objL = (oit + 1);
2586
2587	const float areaDiff =
2588	objL->GetBox().SurfaceArea() - objS->GetBox().SurfaceArea();
2589
2590	if (areaDiff > maxAreaDiff)
2591	{
2592	maxAreaDiff = areaDiff;
2593	backObjectsStart = oit + 1;
2594	}
2595	}
2596
2597	// belongs to back bv
2598	for (oit = objects->begin(); oit != backObjectsStart; ++ oit)
2599	{
2600	frontObjects.push_back(*oit);
2601	}
2602
2603	// belongs to front bv
2604	for (oit = backObjectsStart; oit != oit_end; ++ oit)
2605	{
2606	backObjects.push_back(*oit);
2607	}
2608
2609	cout << "front: " << (int)frontObjects.size() << " back: " << (int)backObjects.size() << " "
2610	<< backObjects.front()->GetBox().SurfaceArea() - frontObjects.back()->GetBox().SurfaceArea() << endl;
2611	}
2612
2613
2614	inline static float AreaRatio(Intersectable smallObj, Intersectable largeObj)
2615	{
2616	const float areaSmall = smallObj->GetBox().SurfaceArea();
2617	const float areaLarge = largeObj->GetBox().SurfaceArea();
2618
2619	return areaLarge / (areaLarge - areaSmall + Limits::Small);
2620	}
2621
2622
2623	bool BvHierarchy::InitialTerminationCriteriaMet(const BvhTraversalData &tData) const
2624	{
2625	const bool terminationCriteriaMet =
2626	(0
2627	\|\| ((int)tData.mNode->mObjects.size() < mInitialMinObjects)
2628	\|\| (tData.mNode->mObjects.back()->GetBox().SurfaceArea() < mInitialMinArea)
2629	\|\| (AreaRatio(tData.mNode->mObjects.front(), tData.mNode->mObjects.back()) > mInitialMaxAreaRatio)
2630	);
2631
2632	cout << "criteria met: " << terminationCriteriaMet << " " << (int)tData.mNode->mObjects.size() << " " << mInitialMinObjects << endl
2633	<< AreaRatio(tData.mNode->mObjects.front(), tData.mNode->mObjects.back()) << " " << mInitialMaxAreaRatio << endl
2634	<< tData.mNode->mObjects.back()->GetBox().SurfaceArea() << " " << mInitialMinArea << endl;
2635
2636	return terminationCriteriaMet;
2637	}
2638
2639
2640	// HACK
2641	float BvHierarchy::GetRenderCostIncrementially(BvhNode *node) const
2642	{
2643	if (node->mRenderCost < 0)
2644	{
2645	//cout <<"p";
2646	if (node->IsLeaf())
2647	{
2648	BvhLeaf leaf = dynamic_cast<BvhLeaf >(node);
2649	node->mRenderCost = EvalAbsCost(leaf->mObjects);
2650	}
2651	else
2652	{
2653	BvhInterior interior = dynamic_cast<BvhInterior >(node);
2654
2655	node->mRenderCost = GetRenderCostIncrementially(interior->GetFront()) +
2656	GetRenderCostIncrementially(interior->GetBack());
2657	}
2658	}
2659
2660	return node->mRenderCost;
2661	}
2662
2663
2664	}

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

Download in other formats: