Context Navigation

source: GTP/trunk/Lib/Vis/Preprocessing/src/KdTree.cpp @ 1286

Revision 1286, 30.8 KB checked in by mattausch, 18 years ago (diff)

Line
1	#include <stack>
2	#include <algorithm>
3	#include <queue>
4	#include "Environment.h"
5	#include "Mesh.h"
6	#include "KdTree.h"
7	#include "ViewCell.h"
8	#include "Beam.h"
9
10
11
12	namespace GtpVisibilityPreprocessor {
13
14	int KdNode::sMailId = 1;
15	int KdNode::sReservedMailboxes = 1;
16
17	inline static bool ilt(Intersectable obj1, Intersectable obj2)
18	{
19	return obj1->mId < obj2->mId;
20	}
21
22
23	KdNode::KdNode(KdInterior *parent):mParent(parent), mMailbox(0)
24	{
25	if (parent)
26	mDepth = parent->mDepth+1;
27	else
28	mDepth = 0;
29	}
30
31
32	KdInterior::~KdInterior()
33	{
34	// recursivly destroy children
35	DEL_PTR(mFront);
36	DEL_PTR(mBack);
37	}
38
39
40	KdTree::KdTree()
41	{
42
43
44	mRoot = new KdLeaf(NULL, 0);
45	Environment::GetSingleton()->GetIntValue("KdTree.Termination.maxNodes", mTermMaxNodes);
46	Environment::GetSingleton()->GetIntValue("KdTree.Termination.maxDepth", mTermMaxDepth);
47	Environment::GetSingleton()->GetIntValue("KdTree.Termination.minCost", mTermMinCost);
48	Environment::GetSingleton()->GetFloatValue("KdTree.Termination.maxCostRatio", mMaxCostRatio);
49	Environment::GetSingleton()->GetFloatValue("KdTree.Termination.ct_div_ci", mCt_div_ci);
50	Environment::GetSingleton()->GetFloatValue("KdTree.splitBorder", mSplitBorder);
51
52	Environment::GetSingleton()->GetBoolValue("KdTree.sahUseFaces", mSahUseFaces);
53
54	char splitType[64];
55	Environment::GetSingleton()->GetStringValue("KdTree.splitMethod", splitType);
56
57	mSplitMethod = SPLIT_SPATIAL_MEDIAN;
58	if (strcmp(splitType, "spatialMedian") == 0)
59	mSplitMethod = SPLIT_SPATIAL_MEDIAN;
60	else
61	if (strcmp(splitType, "objectMedian") == 0)
62	mSplitMethod = SPLIT_OBJECT_MEDIAN;
63	else
64	if (strcmp(splitType, "SAH") == 0)
65	mSplitMethod = SPLIT_SAH;
66	else {
67	cerr<<"Wrong kd split type "<<splitType<<endl;
68	exit(1);
69	}
70	splitCandidates = NULL;
71	}
72
73
74	KdTree::~KdTree()
75	{
76	DEL_PTR(mRoot);
77	}
78
79
80	bool
81	KdTree::Construct()
82	{
83
84	if (!splitCandidates)
85	splitCandidates = new vector<SortableEntry>;
86
87	// first construct a leaf that will get subdivide
88	KdLeaf leaf = (KdLeaf ) mRoot;
89
90	mStat.nodes = 1;
91
92	mBox.Initialize();
93	ObjectContainer::const_iterator mi;
94	for ( mi = leaf->mObjects.begin();
95	mi != leaf->mObjects.end();
96	mi++) {
97	// cout<<(*mi)->GetBox()<<endl;
98	mBox.Include((*mi)->GetBox());
99	}
100
101	cout <<"KdTree Root Box:"<<mBox<<endl;
102	mRoot = Subdivide(TraversalData(leaf, mBox, 0));
103
104	// remove the allocated array
105	delete splitCandidates;
106
107	return true;
108	}
109
110	KdNode *
111	KdTree::Subdivide(const TraversalData &tdata)
112	{
113
114	KdNode *result = NULL;
115
116	priority_queue<TraversalData> tStack;
117	// stack<STraversalData> tStack;
118
119	tStack.push(tdata);
120	AxisAlignedBox3 backBox, frontBox;
121
122	while (!tStack.empty()) {
123	// cout<<mStat.Nodes()<<" "<<mTermMaxNodes<<endl;
124	if (mStat.Nodes() > mTermMaxNodes) {
125	// if ( GetMemUsage() > maxMemory ) {
126	// count statistics on unprocessed leafs
127	while (!tStack.empty()) {
128	EvaluateLeafStats(tStack.top());
129	tStack.pop();
130	}
131	break;
132	}
133
134
135	TraversalData data = tStack.top();
136	tStack.pop();
137
138	KdNode node = SubdivideNode((KdLeaf ) data.mNode,
139	data.mBox,
140	backBox,
141	frontBox
142	);
143
144	if (result == NULL)
145	result = node;
146
147	if (!node->IsLeaf()) {
148
149	KdInterior interior = (KdInterior ) node;
150	// push the children on the stack
151	tStack.push(TraversalData(interior->mBack, backBox, data.mDepth+1));
152	tStack.push(TraversalData(interior->mFront, frontBox, data.mDepth+1));
153
154	} else {
155	EvaluateLeafStats(data);
156	}
157	}
158
159	return result;
160
161	}
162
163
164
165	bool
166	KdTree::TerminationCriteriaMet(const KdLeaf *leaf)
167	{
168	// cerr<<"\n OBJECTS="<<leaf->mObjects.size()<<endl;
169	return
170	((int)leaf->mObjects.size() <= mTermMinCost) \|\|
171	(leaf->mDepth >= mTermMaxDepth);
172
173	}
174
175
176	int
177	KdTree::SelectPlane(KdLeaf *leaf,
178	const AxisAlignedBox3 &box,
179	float &position
180	)
181	{
182	int axis = -1;
183
184	switch (mSplitMethod)
185	{
186	case SPLIT_SPATIAL_MEDIAN: {
187	axis = box.Size().DrivingAxis();
188	position = (box.Min()[axis] + box.Max()[axis])*0.5f;
189	break;
190	}
191	case SPLIT_SAH: {
192	int objectsBack, objectsFront;
193	float costRatio;
194	bool mOnlyDrivingAxis = false;
195	if (mOnlyDrivingAxis) {
196	axis = box.Size().DrivingAxis();
197	costRatio = BestCostRatio(leaf,
198	box,
199	axis,
200	position,
201	objectsBack,
202	objectsFront);
203	} else {
204	costRatio = MAX_FLOAT;
205	for (int i=0; i < 3; i++) {
206	float p;
207	float r = BestCostRatio(leaf,
208	box,
209	i,
210	p,
211	objectsBack,
212	objectsFront);
213	if (r < costRatio) {
214	costRatio = r;
215	axis = i;
216	position = p;
217	}
218	}
219	}
220
221	if (costRatio > mMaxCostRatio) {
222	//cout<<"Too big cost ratio "<<costRatio<<endl;
223	axis = -1;
224	}
225	break;
226	}
227
228	}
229	return axis;
230	}
231
232	KdNode *
233	KdTree::SubdivideNode(
234	KdLeaf *leaf,
235	const AxisAlignedBox3 &box,
236	AxisAlignedBox3 &backBBox,
237	AxisAlignedBox3 &frontBBox
238	)
239	{
240
241	if (TerminationCriteriaMet(leaf))
242	return leaf;
243
244	float position;
245
246	// select subdivision axis
247	int axis = SelectPlane( leaf, box, position );
248
249	if (axis == -1) {
250	return leaf;
251	}
252
253	mStat.nodes+=2;
254	mStat.splits[axis]++;
255
256	// add the new nodes to the tree
257	KdInterior *node = new KdInterior(leaf->mParent);
258
259	node->mAxis = axis;
260	node->mPosition = position;
261	node->mBox = box;
262
263	backBBox = box;
264	frontBBox = box;
265
266	// first count ray sides
267	int objectsBack = 0;
268	int objectsFront = 0;
269
270	backBBox.SetMax(axis, position);
271	frontBBox.SetMin(axis, position);
272
273	ObjectContainer::const_iterator mi;
274
275	for ( mi = leaf->mObjects.begin();
276	mi != leaf->mObjects.end();
277	mi++) {
278	// determine the side of this ray with respect to the plane
279	AxisAlignedBox3 box = (*mi)->GetBox();
280	if (box.Max(axis) > position )
281	objectsFront++;
282
283	if (box.Min(axis) < position )
284	objectsBack++;
285	}
286
287
288	KdLeaf *back = new KdLeaf(node, objectsBack);
289	KdLeaf *front = new KdLeaf(node, objectsFront);
290
291
292	// replace a link from node's parent
293	if ( leaf->mParent )
294	leaf->mParent->ReplaceChildLink(leaf, node);
295
296	// and setup child links
297	node->SetupChildLinks(back, front);
298
299	for (mi = leaf->mObjects.begin();
300	mi != leaf->mObjects.end();
301	mi++) {
302	// determine the side of this ray with respect to the plane
303	AxisAlignedBox3 box = (*mi)->GetBox();
304
305	// for handling multiple objects: keep track of references
306	if (leaf->IsRoot())
307	(*mi)->mReferences = 1; // initialise references at root
308
309	-- (*mi)->mReferences; // remove parent ref
310
311
312	if (box.Max(axis) >= position )
313	{
314	front->mObjects.push_back(*mi);
315	++ (*mi)->mReferences;
316	}
317
318	if (box.Min(axis) < position )
319	{
320	back->mObjects.push_back(*mi);
321	++ (*mi)->mReferences;
322	}
323
324	mStat.objectRefs -= (int)leaf->mObjects.size();
325	mStat.objectRefs += objectsBack + objectsFront;
326	}
327
328	// store objects referenced in more than one leaf
329	// for easy access
330	ProcessMultipleRefs(back);
331	ProcessMultipleRefs(front);
332
333	delete leaf;
334	return node;
335	}
336
337
338	void KdTree::ProcessMultipleRefs(KdLeaf *leaf) const
339	{
340	// find objects from multiple kd-leaves
341	ObjectContainer::const_iterator oit, oit_end = leaf->mObjects.end();
342
343	for (oit = leaf->mObjects.begin(); oit != oit_end; ++ oit)
344	{
345	Intersectable object = oit;
346
347	if (object->mReferences > 1)
348	{
349	leaf->mMultipleObjects.push_back(object);
350	}
351	}
352	}
353
354
355	void
356	KdTreeStatistics::Print(ostream &app) const
357	{
358	app << "===== KdTree statistics ===============\n";
359
360	app << "#N_NODES ( Number of nodes )\n" << nodes << "\n";
361
362	app << "#N_LEAVES ( Number of leaves )\n" << Leaves() << "\n";
363
364	app << "#N_SPLITS ( Number of splits in axes x y z dx dy dz)\n";
365	for (int i=0; i<7; i++)
366	app << splits[i] <<" ";
367	app <<endl;
368
369	app << "#N_RAYREFS ( Number of rayRefs )\n" <<
370	rayRefs << "\n";
371
372	app << "#N_RAYRAYREFS ( Number of rayRefs / ray )\n" <<
373	rayRefs/(double)rays << "\n";
374
375	app << "#N_LEAFRAYREFS ( Number of rayRefs / leaf )\n" <<
376	rayRefs/(double)Leaves() << "\n";
377
378	app << "#N_MAXOBJECTREFS ( Max number of object refs / leaf )\n" <<
379	maxObjectRefs << "\n";
380
381	app << "#N_NONEMPTYRAYREFS ( Number of rayRefs in nonEmpty leaves / non empty leaf )\n" <<
382	rayRefsNonZeroQuery/(double)(Leaves() - zeroQueryNodes) << "\n";
383
384	app << "#N_LEAFDOMAINREFS ( Number of query domain Refs / leaf )\n" <<
385	objectRefs/(double)Leaves() << "\n";
386
387	// app << setprecision(4);
388
389	app << "#N_PEMPTYLEAVES ( Percentage of leaves with zero query domains )\n"<<
390	zeroQueryNodes*100/(double)Leaves()<<endl;
391
392	app << "#N_PMAXDEPTHLEAVES ( Percentage of leaves at maxdepth )\n"<<
393	maxDepthNodes*100/(double)Leaves()<<endl;
394
395	app << "#N_PMINCOSTLEAVES ( Percentage of leaves with minCost )\n"<<
396	minCostNodes*100/(double)Leaves()<<endl;
397
398	app << "#N_ADDED_RAYREFS (Number of dynamically added ray references )\n"<<
399	addedRayRefs<<endl;
400
401	app << "#N_REMOVED_RAYREFS (Number of dynamically removed ray references )\n"<<
402	removedRayRefs<<endl;
403
404	// app << setprecision(4);
405
406	// app << "#N_CTIME ( Construction time [s] )\n"
407	// << Time() << " \n";
408
409	app << "===== END OF KdTree statistics ==========\n";
410
411	}
412
413
414
415	void
416	KdTree::EvaluateLeafStats(const TraversalData &data)
417	{
418
419	// the node became a leaf -> evaluate stats for leafs
420	KdLeaf leaf = (KdLeaf )data.mNode;
421
422	if (data.mDepth > mTermMaxDepth)
423	mStat.maxDepthNodes++;
424
425	if ( (int)(leaf->mObjects.size()) < mTermMinCost)
426	mStat.minCostNodes++;
427
428
429	if ( (int)(leaf->mObjects.size()) > mStat.maxObjectRefs)
430	mStat.maxObjectRefs = (int)leaf->mObjects.size();
431
432	}
433
434
435
436	void
437	KdTree::SortSubdivisionCandidates(
438	KdLeaf *node,
439	const int axis
440	)
441	{
442	splitCandidates->clear();
443
444	int requestedSize = 2*(int)node->mObjects.size();
445	// creates a sorted split candidates array
446	if (splitCandidates->capacity() > 500000 &&
447	requestedSize < (int)(splitCandidates->capacity()/10) ) {
448	delete splitCandidates;
449	splitCandidates = new vector<SortableEntry>;
450	}
451
452	splitCandidates->reserve(requestedSize);
453
454	// insert all queries
455	for(ObjectContainer::const_iterator mi = node->mObjects.begin();
456	mi != node->mObjects.end();
457	mi++) {
458	AxisAlignedBox3 box = (*mi)->GetBox();
459
460	splitCandidates->push_back(SortableEntry(SortableEntry::BOX_MIN,
461	box.Min(axis),
462	*mi)
463	);
464
465
466	splitCandidates->push_back(SortableEntry(SortableEntry::BOX_MAX,
467	box.Max(axis),
468	*mi)
469	);
470	}
471
472	stable_sort(splitCandidates->begin(), splitCandidates->end());
473	}
474
475
476	float
477	KdTree::BestCostRatio(
478	KdLeaf *node,
479	const AxisAlignedBox3 &box,
480	const int axis,
481	float &position,
482	int &objectsBack,
483	int &objectsFront
484	)
485	{
486
487	SortSubdivisionCandidates(node, axis);
488
489	// go through the lists, count the number of objects left and right
490	// and evaluate the following cost funcion:
491	// C = ct_div_ci + (ol + or)/queries
492
493	float totalIntersections = 0.0f;
494	vector<SortableEntry>::const_iterator ci;
495
496	for(ci = splitCandidates->begin();
497	ci < splitCandidates->end();
498	ci++)
499	if ((*ci).type == SortableEntry::BOX_MIN) {
500	totalIntersections += (*ci).intersectable->IntersectionComplexity();
501	}
502
503	float intersectionsLeft = 0;
504	float intersectionsRight = totalIntersections;
505
506	int objectsLeft = 0, objectsRight = (int)node->mObjects.size();
507
508	float minBox = box.Min(axis);
509	float maxBox = box.Max(axis);
510	float boxArea = box.SurfaceArea();
511
512	float minBand = minBox + mSplitBorder*(maxBox - minBox);
513	float maxBand = minBox + (1.0f - mSplitBorder)*(maxBox - minBox);
514
515	float minSum = 1e20f;
516
517	for(ci = splitCandidates->begin();
518	ci < splitCandidates->end();
519	ci++) {
520	switch ((*ci).type) {
521	case SortableEntry::BOX_MIN:
522	objectsLeft++;
523	intersectionsLeft += (*ci).intersectable->IntersectionComplexity();
524	break;
525	case SortableEntry::BOX_MAX:
526	objectsRight--;
527	intersectionsRight -= (*ci).intersectable->IntersectionComplexity();
528	break;
529	}
530
531	if ((ci).value > minBand && (ci).value < maxBand) {
532	AxisAlignedBox3 lbox = box;
533	AxisAlignedBox3 rbox = box;
534	lbox.SetMax(axis, (*ci).value);
535	rbox.SetMin(axis, (*ci).value);
536
537	float sum;
538	if (mSahUseFaces)
539	sum = intersectionsLeftlbox.SurfaceArea() + intersectionsRightrbox.SurfaceArea();
540	else
541	sum = objectsLeftlbox.SurfaceArea() + objectsRightrbox.SurfaceArea();
542
543	// cout<<"pos="<<(*ci).value<<"\t q=("<<ql<<","<<qr<<")\t r=("<<rl<<","<<rr<<")"<<endl;
544	// cout<<"cost= "<<sum<<endl;
545
546	if (sum < minSum) {
547	minSum = sum;
548	position = (*ci).value;
549
550	objectsBack = objectsLeft;
551	objectsFront = objectsRight;
552	}
553	}
554	}
555
556	float oldCost = mSahUseFaces ? totalIntersections : node->mObjects.size();
557	float newCost = mCt_div_ci + minSum/boxArea;
558	float ratio = newCost/oldCost;
559
560	#if 0
561	cout<<"===================="<<endl;
562	cout<<"costRatio="<<ratio<<" pos="<<position<<" t="<<(position - minBox)/(maxBox - minBox)
563	<<"\t o=("<<objectsBack<<","<<objectsFront<<")"<<endl;
564	#endif
565	return ratio;
566	}
567
568	int
569	KdTree::CastRay(
570	Ray &ray
571	)
572	{
573	int hits = 0;
574
575	stack<RayTraversalData> tStack;
576
577	float maxt = 1e6;
578	float mint = 0;
579
580	Intersectable::NewMail();
581
582	if (!mBox.GetMinMaxT(ray, &mint, &maxt))
583	return 0;
584
585	if (mint < 0)
586	mint = 0;
587
588	maxt += Limits::Threshold;
589
590	Vector3 entp = ray.Extrap(mint);
591	Vector3 extp = ray.Extrap(maxt);
592
593	KdNode *node = mRoot;
594	KdNode *farChild;
595	float position;
596	int axis;
597
598
599	while (1) {
600	if (!node->IsLeaf()) {
601	KdInterior in = (KdInterior ) node;
602	position = in->mPosition;
603	axis = in->mAxis;
604
605	if (entp[axis] <= position) {
606	if (extp[axis] <= position) {
607	node = in->mBack;
608	// cases N1,N2,N3,P5,Z2,Z3
609	continue;
610	} else {
611	// case N4
612	node = in->mBack;
613	farChild = in->mFront;
614	}
615	}
616	else {
617	if (position <= extp[axis]) {
618	node = in->mFront;
619	// cases P1,P2,P3,N5,Z1
620	continue;
621	} else {
622	node = in->mFront;
623	farChild = in->mBack;
624	// case P4
625	}
626	}
627	// $$ modification 3.5.2004 - hints from Kamil Ghais
628	// case N4 or P4
629	float tdist = (position - ray.GetLoc(axis)) / ray.GetDir(axis);
630	tStack.push(RayTraversalData(farChild, extp, maxt));
631	extp = ray.GetLoc() + ray.GetDir()*tdist;
632	maxt = tdist;
633	} else {
634	// compute intersection with all objects in this leaf
635	KdLeaf leaf = (KdLeaf ) node;
636	if (ray.mFlags & Ray::STORE_KDLEAVES)
637	ray.kdLeaves.push_back(leaf);
638
639	ObjectContainer::const_iterator mi;
640	for ( mi = leaf->mObjects.begin();
641	mi != leaf->mObjects.end();
642	mi++) {
643	Intersectable object = mi;
644	if (!object->Mailed() ) {
645	object->Mail();
646	if (ray.mFlags & Ray::STORE_TESTED_OBJECTS)
647	ray.testedObjects.push_back(object);
648
649	static int oi=1;
650	if (MeshDebug)
651	cout<<"Object "<<oi++;
652
653	hits += object->CastRay(ray);
654
655	if (MeshDebug) {
656	if (ray.intersections.size())
657	cout<<"nearest t="<<ray.intersections[0].mT<<endl;
658	else
659	cout<<"nearest t=-INF"<<endl;
660	}
661
662	}
663	}
664
665	if (hits && ray.GetType() == Ray::LOCAL_RAY)
666	if (ray.intersections[0].mT <= maxt)
667	break;
668
669	// get the next node from the stack
670	if (tStack.empty())
671	break;
672
673	entp = extp;
674	mint = maxt;
675	if (ray.GetType() == Ray::LINE_SEGMENT && mint > 1.0f)
676	break;
677
678	RayTraversalData &s = tStack.top();
679	node = s.mNode;
680	extp = s.mExitPoint;
681	maxt = s.mMaxT;
682	tStack.pop();
683	}
684	}
685	return hits;
686	}
687
688	int KdTree::CastLineSegment(const Vector3 &origin,
689	const Vector3 &termination,
690	ViewCellContainer &viewcells)
691	{
692	int hits = 0;
693
694	float mint = 0.0f, maxt = 1.0f;
695	const Vector3 dir = termination - origin;
696
697	stack<RayTraversalData> tStack;
698
699	Intersectable::NewMail();
700
701	//maxt += Limits::Threshold;
702
703	Vector3 entp = origin;
704	Vector3 extp = termination;
705
706	KdNode *node = mRoot;
707	KdNode *farChild;
708
709	float position;
710	int axis;
711
712	while (1)
713	{
714	if (!node->IsLeaf())
715	{
716	KdInterior in = dynamic_cast<KdInterior >(node);
717	position = in->mPosition;
718	axis = in->mAxis;
719
720	if (entp[axis] <= position)
721	{
722	if (extp[axis] <= position)
723	{
724	node = in->mBack;
725	// cases N1,N2,N3,P5,Z2,Z3
726	continue;
727	}
728	else
729	{
730	// case N4
731	node = in->mBack;
732	farChild = in->mFront;
733	}
734	}
735	else
736	{
737	if (position <= extp[axis])
738	{
739	node = in->mFront;
740	// cases P1,P2,P3,N5,Z1
741	continue;
742	}
743	else
744	{
745	node = in->mFront;
746	farChild = in->mBack;
747	// case P4
748	}
749	}
750
751	// $$ modification 3.5.2004 - hints from Kamil Ghais
752	// case N4 or P4
753	float tdist = (position - origin[axis]) / dir[axis];
754	//tStack.push(RayTraversalData(farChild, extp, maxt)); //TODO
755	extp = origin + dir * tdist;
756	maxt = tdist;
757	}
758	else
759	{
760	// compute intersection with all objects in this leaf
761	KdLeaf leaf = dynamic_cast<KdLeaf >(node);
762
763	// add view cell to intersections
764	ViewCell *vc = leaf->mViewCell;
765
766	if (!vc->Mailed())
767	{
768	vc->Mail();
769	viewcells.push_back(vc);
770	++ hits;
771	}
772
773	// get the next node from the stack
774	if (tStack.empty())
775	break;
776
777	entp = extp;
778	mint = maxt;
779
780	RayTraversalData &s = tStack.top();
781	node = s.mNode;
782	extp = s.mExitPoint;
783	maxt = s.mMaxT;
784	tStack.pop();
785	}
786	}
787
788	return hits;
789	}
790
791
792	void
793	KdTree::CollectObjects(const AxisAlignedBox3 &box,
794	ObjectContainer &objects)
795	{
796	stack<KdNode *> nodeStack;
797
798	nodeStack.push(mRoot);
799
800	while (!nodeStack.empty()) {
801	KdNode *node = nodeStack.top();
802	nodeStack.pop();
803	if (node->IsLeaf()) {
804	KdLeaf leaf = (KdLeaf )node;
805	for (int j=0; j < leaf->mObjects.size(); j++) {
806	Intersectable *object = leaf->mObjects[j];
807	if (!object->Mailed() && Overlap(box, object->GetBox())) {
808	object->Mail();
809	objects.push_back(object);
810	}
811	}
812	} else {
813	KdInterior interior = (KdInterior )node;
814
815	if ( box.Max()[interior->mAxis] > interior->mPosition )
816	nodeStack.push(interior->mFront);
817
818	if (box.Min()[interior->mAxis] < interior->mPosition)
819	nodeStack.push(interior->mBack);
820	}
821	}
822	}
823
824	void
825	KdTree::CollectObjects(KdNode *n, ObjectContainer &objects)
826	{
827	stack<KdNode *> nodeStack;
828
829	nodeStack.push(n);
830
831	while (!nodeStack.empty()) {
832	KdNode *node = nodeStack.top();
833	nodeStack.pop();
834	if (node->IsLeaf()) {
835	KdLeaf leaf = (KdLeaf )node;
836	for (int j=0; j < leaf->mObjects.size(); j++) {
837	Intersectable *object = leaf->mObjects[j];
838	if (!object->Mailed()) {
839	object->Mail();
840	objects.push_back(object);
841	}
842	}
843	} else {
844	KdInterior interior = (KdInterior )node;
845	nodeStack.push(interior->mFront);
846	nodeStack.push(interior->mBack);
847	}
848	}
849	}
850
851	// Find random neighbor which was not mailed
852	KdNode *
853	KdTree::FindRandomNeighbor(KdNode *n,
854	bool onlyUnmailed
855	)
856	{
857	stack<KdNode *> nodeStack;
858
859	nodeStack.push(mRoot);
860
861	AxisAlignedBox3 box = GetBox(n);
862	int mask = rand();
863
864	while (!nodeStack.empty()) {
865	KdNode *node = nodeStack.top();
866	nodeStack.pop();
867	if (node->IsLeaf()) {
868	if ( node != n && (!onlyUnmailed \|\| !node->Mailed()) )
869	return node;
870	} else {
871	KdInterior interior = (KdInterior )node;
872	if (interior->mPosition > box.Max(interior->mAxis))
873	nodeStack.push(interior->mBack);
874	else
875	if (interior->mPosition < box.Min(interior->mAxis))
876	nodeStack.push(interior->mFront);
877	else {
878	// random decision
879	if (mask&1)
880	nodeStack.push(interior->mBack);
881	else
882	nodeStack.push(interior->mFront);
883	mask = mask>>1;
884	}
885	}
886	}
887
888	return NULL;
889	}
890
891	int
892	KdTree::FindNeighbors(KdNode *n,
893	vector<KdNode *> &neighbors,
894	bool onlyUnmailed
895	)
896	{
897	stack<KdNode *> nodeStack;
898
899	nodeStack.push(mRoot);
900
901	AxisAlignedBox3 box = GetBox(n);
902
903	while (!nodeStack.empty()) {
904	KdNode *node = nodeStack.top();
905	nodeStack.pop();
906	if (node->IsLeaf()) {
907	if ( node != n && (!onlyUnmailed \|\| !node->Mailed()) )
908	neighbors.push_back(node);
909	} else {
910	KdInterior interior = (KdInterior )node;
911	if (interior->mPosition > box.Max(interior->mAxis))
912	nodeStack.push(interior->mBack);
913	else
914	if (interior->mPosition < box.Min(interior->mAxis))
915	nodeStack.push(interior->mFront);
916	else {
917	// random decision
918	nodeStack.push(interior->mBack);
919	nodeStack.push(interior->mFront);
920	}
921	}
922	}
923
924	return (int)neighbors.size();
925	}
926
927	// Find random neighbor which was not mailed
928	KdNode *
929	KdTree::GetRandomLeaf(const Plane3 &plane)
930	{
931	stack<KdNode *> nodeStack;
932
933	nodeStack.push(mRoot);
934
935	int mask = rand();
936
937	while (!nodeStack.empty()) {
938	KdNode *node = nodeStack.top();
939	nodeStack.pop();
940	if (node->IsLeaf()) {
941	return node;
942	} else {
943	KdInterior interior = (KdInterior )node;
944	KdNode *next;
945	if (GetBox(interior->mBack).Side(plane) < 0)
946	next = interior->mFront;
947	else
948	if (GetBox(interior->mFront).Side(plane) < 0)
949	next = interior->mBack;
950	else {
951	// random decision
952	if (mask&1)
953	next = interior->mBack;
954	else
955	next = interior->mFront;
956	mask = mask>>1;
957	}
958	nodeStack.push(next);
959	}
960	}
961
962
963	return NULL;
964	}
965
966	void
967	KdTree::CollectLeaves(vector<KdLeaf *> &leaves)
968	{
969	stack<KdNode *> nodeStack;
970	nodeStack.push(mRoot);
971
972	while (!nodeStack.empty()) {
973	KdNode *node = nodeStack.top();
974	nodeStack.pop();
975	if (node->IsLeaf()) {
976	KdLeaf leaf = (KdLeaf )node;
977	leaves.push_back(leaf);
978	} else {
979	KdInterior interior = (KdInterior )node;
980	nodeStack.push(interior->mBack);
981	nodeStack.push(interior->mFront);
982	}
983	}
984	}
985
986	void
987	KdTree::CreateAndCollectViewCells(ViewCellContainer &vc) const
988	{
989	stack<KdNode *> nodeStack;
990	nodeStack.push(mRoot);
991
992	while (!nodeStack.empty()) {
993	KdNode *node = nodeStack.top();
994	nodeStack.pop();
995	if (node->IsLeaf()) {
996	KdLeaf leaf = (KdLeaf )node;
997	// kdtree used as view cell container => create view cell
998	KdViewCell *viewCell = new KdViewCell();
999	leaf->mViewCell = viewCell;
1000	// push back pointer to this leaf
1001	viewCell->mLeaf = leaf;
1002	vc.push_back(viewCell);
1003	} else {
1004	KdInterior interior = (KdInterior )node;
1005	nodeStack.push(interior->mBack);
1006	nodeStack.push(interior->mFront);
1007	}
1008	}
1009	}
1010
1011	int
1012	KdTree::CollectLeafPvs()
1013	{
1014	int totalPvsSize = 0;
1015	stack<KdNode *> nodeStack;
1016
1017	nodeStack.push(mRoot);
1018
1019	while (!nodeStack.empty()) {
1020	KdNode *node = nodeStack.top();
1021	nodeStack.pop();
1022	if (node->IsLeaf()) {
1023	KdLeaf leaf = (KdLeaf )node;
1024	for (int j=0; j < leaf->mObjects.size(); j++) {
1025	Intersectable *object = leaf->mObjects[j];
1026	if (!object->Mailed()) {
1027	object->Mail();
1028	// add this node to pvs of all nodes it can see
1029	KdPvsMap::iterator ni = object->mKdPvs.mEntries.begin();
1030	for (; ni != object->mKdPvs.mEntries.end(); ni++) {
1031	KdNode node = (ni).first;
1032	// $$ JB TEMPORARY solution -> should add object PVS or explictly computed
1033	// kd tree PVS
1034	float contribution;
1035	if (leaf->mKdPvs.AddSample(node, 1.0f, contribution))
1036	totalPvsSize++;
1037	}
1038	}
1039	}
1040	} else {
1041	KdInterior interior = (KdInterior )node;
1042	nodeStack.push(interior->mFront);
1043	nodeStack.push(interior->mBack);
1044	}
1045	}
1046
1047	return totalPvsSize;
1048	}
1049
1050
1051	KdNode *
1052	KdTree::GetRandomLeaf(const bool onlyUnmailed)
1053	{
1054	stack<KdNode *> nodeStack;
1055	nodeStack.push(mRoot);
1056
1057	int mask = rand();
1058
1059	while (!nodeStack.empty()) {
1060	KdNode *node = nodeStack.top();
1061	nodeStack.pop();
1062	if (node->IsLeaf()) {
1063	if ( (!onlyUnmailed \|\| !node->Mailed()) )
1064	return node;
1065	} else {
1066	KdInterior interior = (KdInterior )node;
1067	// random decision
1068	if (mask&1)
1069	nodeStack.push(interior->mBack);
1070	else
1071	nodeStack.push(interior->mFront);
1072	mask = mask>>1;
1073	}
1074	}
1075	return NULL;
1076	}
1077
1078
1079	int
1080	KdTree::CastBeam(
1081	Beam &beam
1082	)
1083	{
1084	stack<KdNode *> nodeStack;
1085	nodeStack.push(mRoot);
1086
1087	while (!nodeStack.empty()) {
1088	KdNode *node = nodeStack.top();
1089	nodeStack.pop();
1090
1091	int side = beam.ComputeIntersection(GetBox(node));
1092	switch (side) {
1093	case -1:
1094	beam.mKdNodes.push_back(node);
1095	break;
1096	case 0:
1097	if (node->IsLeaf())
1098	beam.mKdNodes.push_back(node);
1099	else {
1100	KdInterior interior = (KdInterior )node;
1101	KdNode *first = interior->mBack;
1102	KdNode *second = interior->mFront;
1103
1104	if (interior->mAxis < 3) {
1105	// spatial split -> decide on the order of the nodes
1106	if (beam.mPlanes[0].mNormal[interior->mAxis] > 0)
1107	swap(first, second);
1108	}
1109
1110	nodeStack.push(first);
1111	nodeStack.push(second);
1112	}
1113	break;
1114	// default: cull
1115	}
1116	}
1117
1118	if (beam.mFlags & Beam::STORE_OBJECTS)
1119	{
1120	vector<KdNode *>::const_iterator it, it_end = beam.mKdNodes.end();
1121
1122	Intersectable::NewMail();
1123	for (it = beam.mKdNodes.begin(); it != it_end; ++ it)
1124	{
1125	CollectObjects(*it, beam.mObjects);
1126	}
1127	}
1128
1129	return (int)beam.mKdNodes.size();
1130	}
1131
1132
1133	#define TYPE_INTERIOR -2
1134	#define TYPE_LEAF -3
1135
1136
1137	void KdTree::ExportBinLeaf(OUT_STREAM &stream, KdLeaf *leaf)
1138	{
1139	ObjectContainer::const_iterator it, it_end = leaf->mObjects.end();
1140
1141	int type = TYPE_LEAF;
1142	int size = (int)leaf->mObjects.size();
1143
1144	stream.write(reinterpret_cast<char *>(&type), sizeof(int));
1145	stream.write(reinterpret_cast<char *>(&size), sizeof(int));
1146
1147	for (it = leaf->mObjects.begin(); it != it_end; ++ it)
1148	{
1149	Intersectable obj = it;
1150	int id = obj->mId;
1151
1152	//stream.write(reinterpret_cast<char *>(&origin), sizeof(Vector3));
1153	stream.write(reinterpret_cast<char *>(&id), sizeof(int));
1154	}
1155	}
1156
1157
1158	KdLeaf *KdTree::ImportBinLeaf(IN_STREAM &stream,
1159	KdInterior *parent,
1160	const ObjectContainer &objects)
1161	{
1162	int leafId = TYPE_LEAF;
1163	int objId = leafId;
1164	int size;
1165
1166	stream.read(reinterpret_cast<char *>(&size), sizeof(int));
1167	KdLeaf *leaf = new KdLeaf(parent, size);
1168
1169	MeshInstance dummyInst(NULL);
1170	//Debug << "objects size: " << size << endl;
1171
1172	// read object ids
1173	for (int i = 0; i < size; ++ i)
1174	{
1175	stream.read(reinterpret_cast<char *>(&objId), sizeof(int));
1176	dummyInst.SetId(objId);
1177
1178	ObjectContainer::const_iterator oit =
1179	lower_bound(objects.begin(), objects.end(), (Intersectable *)&dummyInst, ilt);
1180
1181	if ((oit != objects.end()) && ((*oit)->GetId() == objId))
1182	{
1183	leaf->mObjects.push_back(*oit);
1184	}
1185	else
1186	{
1187	Debug << "error: object with id " << objId << " does not exist" << endl;
1188	}
1189	}
1190
1191	return leaf;
1192	}
1193
1194
1195	void KdTree::ExportBinInterior(OUT_STREAM &stream, KdInterior *interior)
1196	{
1197	int interiorid = TYPE_INTERIOR;
1198	stream.write(reinterpret_cast<char *>(&interiorid), sizeof(int));
1199
1200	int axis = interior->mAxis;
1201	float pos = interior->mPosition;
1202
1203	stream.write(reinterpret_cast<char *>(&axis), sizeof(int));
1204	stream.write(reinterpret_cast<char *>(&pos), sizeof(float));
1205	}
1206
1207
1208	KdInterior KdTree::ImportBinInterior(IN_STREAM &stream, KdInterior parent)
1209	{
1210	KdInterior *interior = new KdInterior(parent);
1211
1212	int axis = interior->mAxis;
1213	float pos = interior->mPosition;
1214
1215	stream.read(reinterpret_cast<char *>(&axis), sizeof(int));
1216	stream.read(reinterpret_cast<char *>(&pos), sizeof(float));
1217
1218	interior->mAxis = axis;
1219	interior->mPosition = pos;
1220
1221	return interior;
1222	}
1223
1224
1225	bool KdTree::ExportBinTree(const string &filename)
1226	{
1227	OUT_STREAM stream(filename.c_str(), OUT_BIN_MODE);
1228
1229	//if (!stream.is_open()) return false;
1230
1231	// export binary version of mesh
1232	queue<KdNode *> tStack;
1233	tStack.push(mRoot);
1234
1235	while(!tStack.empty())
1236	{
1237	KdNode *node = tStack.front();
1238	tStack.pop();
1239
1240	if (node->IsLeaf())
1241	{
1242	//Debug << "l";
1243	ExportBinLeaf(stream, dynamic_cast<KdLeaf *>(node));
1244	}
1245	else
1246	{
1247	//Debug << "i";
1248	KdInterior interior = dynamic_cast<KdInterior >(node);
1249
1250	ExportBinInterior(stream, interior);
1251
1252	tStack.push(interior->mFront);
1253	tStack.push(interior->mBack);
1254	}
1255	}
1256
1257	return true;
1258	}
1259
1260
1261	KdNode *KdTree::LoadNextNode(IN_STREAM &stream,
1262	KdInterior *parent,
1263	const ObjectContainer &objects)
1264	{
1265	int nodeType;
1266	stream.read(reinterpret_cast<char *>(&nodeType), sizeof(int));
1267
1268	if (nodeType == TYPE_LEAF)
1269	{
1270	return ImportBinLeaf(stream, dynamic_cast<KdInterior *>(parent), objects);
1271	}
1272
1273	if (nodeType == TYPE_INTERIOR)
1274	{
1275	return ImportBinInterior(stream, dynamic_cast<KdInterior *>(parent));
1276	}
1277
1278	Debug << "error! loading failed!" << endl;
1279
1280	return NULL;
1281	}
1282
1283
1284	bool KdTree::LoadBinTree(const string &filename, ObjectContainer &objects)
1285	{
1286	// export binary version of mesh
1287	queue<TraversalData> tStack;
1288	IN_STREAM stream(filename.c_str(), IN_BIN_MODE);
1289
1290	if (!stream.is_open()) return false;
1291
1292	std::stable_sort(objects.begin(), objects.end(), ilt);
1293
1294	mBox.Initialize();
1295	ObjectContainer::const_iterator oit, oit_end = objects.end();
1296
1297	//-- compute bounding box
1298	for (oit = objects.begin(); oit != oit_end; ++ oit)
1299	{
1300	Intersectable obj = oit;
1301	// compute bounding box of view space
1302	mBox.Include(obj->GetBox());
1303	}
1304
1305	DEL_PTR(mRoot); // we make a new root
1306
1307	KdNode *node = LoadNextNode(stream, NULL, objects);
1308	mRoot = node;
1309
1310	tStack.push(TraversalData(node, mBox, 0));
1311	mStat.Reset();
1312	mStat.nodes = 1;
1313
1314	while(!tStack.empty())
1315	{
1316	TraversalData tData = tStack.front();
1317	tStack.pop();
1318
1319	KdNode *node = tData.mNode;
1320
1321	if (!node->IsLeaf())
1322	{
1323	mStat.nodes += 2;
1324
1325	//Debug << "i" ;
1326	KdInterior interior = dynamic_cast<KdInterior >(node);
1327	interior->mBox = tData.mBox;
1328
1329	KdNode *front = LoadNextNode(stream, interior, objects);
1330	KdNode *back = LoadNextNode(stream, interior, objects);
1331
1332	interior->SetupChildLinks(back, front);
1333
1334	++ mStat.splits[interior->mAxis];
1335
1336	// compute new bounding box
1337	AxisAlignedBox3 frontBox, backBox;
1338
1339	tData.mBox.Split(interior->mAxis,
1340	interior->mPosition,
1341	frontBox,
1342	backBox);
1343
1344	tStack.push(TraversalData(front, frontBox, tData.mDepth + 1));
1345	tStack.push(TraversalData(back, backBox, tData.mDepth + 1));
1346	}
1347	else
1348	{
1349	EvaluateLeafStats(tData);
1350	//Debug << "l" << endl;
1351	}
1352	}
1353	Debug << mStat << endl;
1354
1355	return true;
1356	}
1357
1358
1359	}

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

Download in other formats: