Context Navigation

source: GTP/trunk/Lib/Vis/Preprocessing/src/VspKdTree.cpp @ 667

Revision 667, 59.5 KB checked in by mattausch, 19 years ago (diff)
test version build script for testing. code is set uo for testing also

Line
1
2	// ================================================================
3	// $Id: lsds_kdtree.cpp,v 1.18 2005/04/16 09:34:21 bittner Exp $
4	// ****************************************************************
5	/**
6	The KD tree based LSDS
7	*/
8	// Initial coding by
9	/**
10	@author Jiri Bittner
11	*/
12
13	// Standard headers
14	#include <stack>
15	#include <queue>
16	#include <algorithm>
17	#include <fstream>
18	#include <string>
19
20	#include "VspKdTree.h"
21
22	#include "Environment.h"
23	#include "VssRay.h"
24	#include "Intersectable.h"
25	#include "Ray.h"
26	#include "RayInfo.h"
27	#include "ViewCellsManager.h"
28	#include "ViewCellBsp.h"
29
30	// Static variables
31	int VspKdLeaf::sMailId = 0;
32	int VspKdMergeCandidate::sMaxPvsSize = 150;
33	float VspKdMergeCandidate::sOverallCost = Limits::Small;
34	#define USE_FIXEDPOINT_T 0
35
36	/// Adds object to the pvs of the front and back node
37	inline void AddObject2Pvs(Intersectable *object,
38	const int side,
39	int &pvsBack,
40	int &pvsFront)
41	{
42	if (!object)
43	return;
44
45	if (side <= 0)
46	{
47	if (!object->Mailed() && !object->Mailed(2))
48	{
49	++ pvsBack;
50
51	if (object->Mailed(1))
52	object->Mail(2);
53	else
54	object->Mail();
55	}
56	}
57
58	if (side >= 0)
59	{
60	if (!object->Mailed(1) && !object->Mailed(2))
61	{
62	++ pvsFront;
63
64	if (object->Mailed())
65	object->Mail(2);
66	else
67	object->Mail(1);
68	}
69	}
70	}
71
72	/**************************************************************/
73	/* class VspKdNode implementation */
74	/**************************************************************/
75
76	// Inline constructor
77	VspKdNode::VspKdNode(VspKdNode *p):
78	mParent(p), mAxis(-1), mDepth(p ? p->mDepth + 1 : 0)
79	{}
80
81	VspKdNode::~VspKdNode()
82	{
83	};
84
85	inline VspKdNode *VspKdNode::GetParent() const
86	{
87	return mParent;
88	}
89
90	inline void VspKdNode::SetParent(VspKdNode *p)
91	{
92	mParent = p;
93	}
94
95	bool VspKdNode::IsLeaf() const
96	{
97	return mAxis == -1;
98	}
99
100	int VspKdNode::GetAccessTime()
101	{
102	return 0x7FFFFFF;
103	}
104
105	/**************************************************************/
106	/* VspKdInterior implementation */
107	/**************************************************************/
108
109	VspKdInterior::VspKdInterior(VspKdInterior *p):
110	VspKdNode(p), mBack(NULL), mFront(NULL), mAccesses(0), mLastAccessTime(-1)
111	{
112	}
113
114	int VspKdInterior::GetAccessTime()
115	{
116	return mLastAccessTime;
117	}
118
119	void VspKdInterior::SetupChildLinks(VspKdNode b, VspKdNode f)
120	{
121	mBack = b;
122	mFront = f;
123	b->SetParent(this);
124	f->SetParent(this);
125	}
126
127	void VspKdInterior::ReplaceChildLink(VspKdNode *oldChild,
128	VspKdNode *newChild)
129	{
130	if (mBack == oldChild)
131	mBack = newChild;
132	else
133	mFront = newChild;
134	}
135
136	int VspKdInterior::Type() const
137	{
138	return EInterior;
139	}
140
141	VspKdInterior::~VspKdInterior()
142	{
143	DEL_PTR(mBack);
144	DEL_PTR(mFront);
145	}
146
147	void VspKdInterior::Print(ostream &s) const
148	{
149	switch (mAxis)
150	{
151	case 0:
152	s << "x "; break;
153	case 1:
154	s << "y "; break;
155	case 2:
156	s << "z "; break;
157	}
158
159	s << mPosition << " ";
160
161	mBack->Print(s);
162	mFront->Print(s);
163	}
164
165	int VspKdInterior::ComputeRayIntersection(const RayInfo &rayData, float &t)
166	{
167	return rayData.ComputeRayIntersection(mAxis, mPosition, t);
168	}
169
170	VspKdNode *VspKdInterior::GetBack() const
171	{
172	return mBack;
173	}
174
175	VspKdNode *VspKdInterior::GetFront() const
176	{
177	return mFront;
178	}
179
180
181	/*******************************************************************/
182	/* class VspKdLeaf implementation */
183	/*******************************************************************/
184
185
186	VspKdLeaf::VspKdLeaf(VspKdNode *p,
187	const int nRays,
188	const int maxCostMisses):
189	VspKdNode(p), mRays(), mPvsSize(0), mValidPvs(false), mViewCell(NULL),
190	mMaxCostMisses(maxCostMisses), mPvs(NULL), mVolume(0)
191	{
192	mRays.reserve(nRays);
193	}
194
195	VspKdLeaf::~VspKdLeaf()
196	{
197	DEL_PTR(mPvs);
198	}
199
200
201	int VspKdLeaf::GetMaxCostMisses()
202	{
203	return mMaxCostMisses;
204	}
205
206
207	int VspKdLeaf::Type() const
208	{
209	return ELeaf;
210	}
211
212	void VspKdLeaf::Print(ostream &s) const
213	{
214	s << endl << "L: r = " << (int)mRays.size() << endl;
215	};
216
217	void VspKdLeaf::AddRay(const RayInfo &data)
218	{
219	mValidPvs = false;
220	mRays.push_back(data);
221	data.mRay->Ref();
222	}
223
224	int VspKdLeaf::GetPvsSize() const
225	{
226	return mPvsSize;
227	}
228
229	void VspKdLeaf::SetPvsSize(const int s)
230	{
231	mPvsSize = s;
232	}
233
234	void VspKdLeaf::Mail()
235	{
236	mMailbox = sMailId;
237	}
238
239	void VspKdLeaf::NewMail()
240	{
241	++ sMailId;
242	}
243
244	bool VspKdLeaf::Mailed() const
245	{
246	return mMailbox == sMailId;
247	}
248
249	bool VspKdLeaf::Mailed(const int mail) const
250	{
251	return mMailbox == mail;
252	}
253
254	int VspKdLeaf::GetMailbox() const
255	{
256	return mMailbox;
257	}
258
259	float VspKdLeaf::GetAvgRayContribution() const
260	{
261	return GetPvsSize() / ((float)mRays.size() + Limits::Small);
262	}
263
264	float VspKdLeaf::GetSqrRayContribution() const
265	{
266	return sqr(GetAvgRayContribution());
267	}
268
269	RayInfoContainer &VspKdLeaf::GetRays()
270	{
271	return mRays;
272	}
273
274	void VspKdLeaf::ExtractPvs(ObjectContainer &objects) const
275	{
276	RayInfoContainer::const_iterator it, it_end = mRays.end();
277
278	for (it = mRays.begin(); it != it_end; ++ it)
279	{
280	if ((*it).mRay->mTerminationObject)
281	objects.push_back((*it).mRay->mTerminationObject);
282	if ((*it).mRay->mOriginObject)
283	objects.push_back((*it).mRay->mOriginObject);
284	}
285	}
286
287	void VspKdLeaf::GetRays(VssRayContainer &rays)
288	{
289	RayInfoContainer::const_iterator it, it_end = mRays.end();
290
291	for (it = mRays.begin(); it != mRays.end(); ++ it)
292	rays.push_back((*it).mRay);
293	}
294
295	VspKdViewCell *VspKdLeaf::GetViewCell()
296	{
297	return mViewCell;
298	}
299
300	void VspKdLeaf::SetViewCell(VspKdViewCell *viewCell)
301	{
302	mViewCell = viewCell;
303	}
304
305
306	void VspKdLeaf::UpdatePvsSize()
307	{
308	if (!mValidPvs)
309	{
310	Intersectable::NewMail();
311	int pvsSize = 0;
312	for(RayInfoContainer::iterator ri = mRays.begin();
313	ri != mRays.end(); ++ ri)
314	{
315	if ((*ri).mRay->IsActive())
316	{
317	Intersectable *object;
318	#if BIDIRECTIONAL_RAY
319	object = (*ri).mRay->mOriginObject;
320
321	if (object && !object->Mailed())
322	{
323	++ pvsSize;
324	object->Mail();
325	}
326	#endif
327	object = (*ri).mRay->mTerminationObject;
328	if (object && !object->Mailed())
329	{
330	++ pvsSize;
331	object->Mail();
332	}
333	}
334	}
335
336	mPvsSize = pvsSize;
337	mValidPvs = true;
338	}
339	}
340
341
342	/***************************************************************/
343	/* class VspKdIntermediate implementation */
344	/***************************************************************/
345
346
347	VspKdIntermediate::VspKdIntermediate(VspKdInterior *p):
348	VspKdNode(p), mBack(NULL), mFront(NULL), mAccesses(0), mLastAccessTime(-1)
349	{
350	}
351
352
353	VspKdIntermediate::~VspKdIntermediate()
354	{
355	DEL_PTR(mFront);
356	DEL_PTR(mBack);
357	}
358
359
360	int VspKdIntermediate::Type() const
361	{
362	return EIntermediate;
363	}
364
365
366	VspKdLeaf *VspKdIntermediate::GetBack() const
367	{
368	return mBack;
369	}
370
371
372	VspKdLeaf *VspKdIntermediate::GetFront() const
373	{
374	return mFront;
375	}
376
377
378	void VspKdIntermediate::Print(ostream &s) const
379	{
380	s << endl << mSplitPlane << endl;
381	}
382
383
384	void VspKdIntermediate::SetupChildLinks(VspKdLeaf b, VspKdLeaf f)
385	{
386	mBack = b;
387	mFront = f;
388
389	b->SetParent(this);
390	f->SetParent(this);
391	}
392
393
394	int VspKdIntermediate::ComputeRayIntersection(const RayInfo &rayData, float &t)
395	{
396	return rayData.ComputeRayIntersection(mSplitPlane, t);
397	}
398
399
400	/*************************************************************/
401	/* class VspKdTree implementation */
402	/*************************************************************/
403
404
405	VspKdTree::VspKdTree(): mOnlyDrivingAxis(false)
406	{
407	environment->GetIntValue("VspKdTree.Termination.maxDepth", mTermMaxDepth);
408	environment->GetIntValue("VspKdTree.Termination.minPvs", mTermMinPvs);
409	environment->GetIntValue("VspKdTree.Termination.minRays", mTermMinRays);
410	environment->GetFloatValue("VspKdTree.Termination.maxRayContribution", mTermMaxRayContribution);
411	environment->GetFloatValue("VspKdTree.Termination.maxCostRatio", mTermMaxCostRatio);
412	environment->GetFloatValue("VspKdTree.Termination.minSize", mTermMinSize);
413
414	environment->GetIntValue("VspKdTree.Termination.missTolerance", mTermMissTolerance);
415
416	mTermMinSize = sqr(mTermMinSize);
417
418	environment->GetFloatValue("VspKdTree.epsilon", mEpsilon);
419	environment->GetFloatValue("VspKdTree.ct_div_ci", mCtDivCi);
420
421	environment->GetFloatValue("VspKdTree.maxTotalMemory", mMaxTotalMemory);
422	environment->GetFloatValue("VspKdTree.maxStaticMemory", mMaxStaticMemory);
423
424	environment->GetIntValue("VspKdTree.accessTimeThreshold", mAccessTimeThreshold);
425	environment->GetIntValue("VspKdTree.minCollapseDepth", mMinCollapseDepth);
426
427	//environment->GetIntValue("ViewCells.maxViewCells", mMaxViewCells);
428
429	environment->GetIntValue("VspKdTree.PostProcess.minViewCells", mMergeMinViewCells);
430	environment->GetFloatValue("VspKdTree.PostProcess.maxCostRatio", mMergeMaxCostRatio);
431	environment->GetIntValue("VspKdTree.PostProcess.maxPvsSize",
432	VspKdMergeCandidate::sMaxPvsSize);
433
434	environment->GetBoolValue("VspKdTree.splitUseOnlyDrivingAxis", mOnlyDrivingAxis);
435	environment->GetIntValue("VspKdTree.Termination.maxViewCells", mMaxViewCells);
436
437	//-- output
438	Debug << "======= vsp kd tree options ========" << endl;
439	Debug << "max depth: "<< mTermMaxDepth << endl;
440	Debug << "min pvs: "<< mTermMinPvs << endl;
441	Debug << "min rays: "<< mTermMinRays << endl;
442	Debug << "max ray contribution: "<< mTermMaxRayContribution << endl;
443	Debug << "max cost ratio: " << mTermMaxCostRatio << endl;
444	Debug << "min size: " << mTermMinSize << endl;
445
446	Debug << "split type: ";
447
448	//-- split type
449	char sname[128];
450	environment->GetStringValue("VspKdTree.splitType", sname);
451	string name(sname);
452
453	if (name.compare("regular") == 0)
454	{
455	Debug << "regular split" << endl;
456	splitType = ESplitRegular;
457	}
458	else
459	{
460	if (name.compare("heuristic") == 0)
461	{
462	Debug << "heuristic split" << endl;
463	splitType = ESplitHeuristic;
464	}
465	else
466	{
467	cerr << "Invalid VspKdTree split type " << name << endl;
468	exit(1);
469	}
470	}
471
472	mRoot = NULL;
473	mSplitCandidates = new vector<SortableEntry>;
474	}
475
476
477	VspKdTree::~VspKdTree()
478	{
479	DEL_PTR(mRoot);
480	DEL_PTR(mSplitCandidates);
481	}
482
483	void VspKdStatistics::Print(ostream &app) const
484	{
485	app << "===== VspKdTree statistics ===============\n";
486
487	app << "#N_RAYS ( Number of rays )\n"
488	<< rays << endl;
489
490	app << "#N_INITPVS ( Initial PVS size )\n"
491	<< initialPvsSize << endl;
492
493	app << "#N_NODES ( Number of nodes )\n" << nodes << "\n";
494
495	app << "#N_LEAVES ( Number of leaves )\n" << Leaves() << "\n";
496
497	app << "#N_SPLITS ( Number of splits in axes x y z)\n";
498
499	for (int i = 0; i < 3; ++ i)
500	app << splits[i] << " ";
501	app << endl;
502
503	app << "#N_RAYREFS ( Number of rayRefs )\n" <<
504	rayRefs << "\n";
505
506	app << "#N_RAYRAYREFS ( Number of rayRefs / ray )\n" <<
507	rayRefs / (double)rays << "\n";
508
509	app << "#N_LEAFRAYREFS ( Number of rayRefs / leaf )\n" <<
510	rayRefs / (double)Leaves() << "\n";
511
512	app << "#N_MAXRAYREFS ( Max number of rayRefs / leaf )\n" <<
513	maxRayRefs << "\n";
514
515	// app << setprecision(4);
516
517	app << "#N_PMAXDEPTHLEAVES ( Percentage of leaves at maxdepth )\n" <<
518	maxDepthNodes * 100 / (double)Leaves() << endl;
519
520	app << "#N_PMINPVSLEAVES ( Percentage of leaves with mininimal PVS )\n" <<
521	minPvsNodes * 100 / (double)Leaves() << endl;
522
523	app << "#N_PMINRAYSLEAVES ( Percentage of leaves with minimal number of rays)\n" <<
524	minRaysNodes * 100 / (double)Leaves() << endl;
525
526	app << "#N_PMINSIZELEAVES ( Percentage of leaves with minSize )\n"<<
527	minSizeNodes * 100 / (double)Leaves() << endl;
528
529	app << "#N_PMAXRAYCONTRIBLEAVES ( Percentage of leaves with maximal ray contribution )\n" <<
530	maxRayContribNodes * 100 / (double)Leaves() << endl;
531
532	app << "#N_PMAXRCOSTLEAVES ( Percentage of leaves with maximal cost ratio )\n" <<
533	maxCostNodes * 100 / (double)Leaves() << endl;
534
535	app << "#N_ADDED_RAYREFS ( Number of dynamically added ray references )\n"<<
536	addedRayRefs << endl;
537
538	app << "#N_REMOVED_RAYREFS ( Number of dynamically removed ray references )\n"<<
539	removedRayRefs << endl;
540
541	// app << setprecision(4);
542
543	app << "#N_( Maximal PVS size / leaf)\n"
544	<< maxPvsSize << endl;
545
546	app << "#N_( Average PVS size / leaf)\n"
547	<< (double) accPvsSize / (double)Leaves() << endl;
548
549	app << "#N_CTIME ( Construction time [s] )\n"
550	<< Time() << " \n";
551
552	app << "===== END OF VspKdTree statistics ==========\n";
553	}
554
555
556	void VspKdTree::CollectViewCells(ViewCellContainer &viewCells) const
557	{
558	stack<VspKdNode *> nodeStack;
559	nodeStack.push(mRoot);
560
561	ViewCell::NewMail();
562
563	while (!nodeStack.empty())
564	{
565	VspKdNode *node = nodeStack.top();
566	nodeStack.pop();
567
568	if (node->IsLeaf())
569	{
570	ViewCell viewCell = dynamic_cast<VspKdLeaf >(node)->mViewCell;
571
572	if (!viewCell->Mailed())
573	{
574	viewCell->Mail();
575	viewCells.push_back(viewCell);
576	}
577	}
578	else
579	{
580	VspKdInterior interior = dynamic_cast<VspKdInterior >(node);
581
582	nodeStack.push(interior->GetFront());
583	nodeStack.push(interior->GetBack());
584	}
585	}
586	}
587
588	void VspKdTree::Construct(const VssRayContainer &rays,
589	AxisAlignedBox3 *forcedBoundingBox)
590	{
591	mStat.Start();
592
593	mMaxMemory = mMaxStaticMemory;
594	DEL_PTR(mRoot);
595
596	mRoot = new VspKdLeaf(NULL, (int)rays.size());
597
598	// first construct a leaf that will get subdivided
599	VspKdLeaf leaf = dynamic_cast<VspKdLeaf >(mRoot);
600
601	mStat.nodes = 1;
602	mBox.Initialize();
603
604	//-- compute bounding box
605	if (forcedBoundingBox)
606	mBox = *forcedBoundingBox;
607	else
608	for (VssRayContainer::const_iterator ri = rays.begin(); ri != rays.end(); ++ ri)
609	{
610	if ((*ri)->mOriginObject)
611	mBox.Include((*ri)->GetOrigin());
612	if ((*ri)->mTerminationObject)
613	mBox.Include((*ri)->GetTermination());
614	}
615
616	Debug << "bbox = " << mBox << endl;
617
618	for (VssRayContainer::const_iterator ri = rays.begin(); ri != rays.end(); ++ ri)
619	{
620	//leaf->AddRay(RayInfo(*ri));
621	VssRay ray = ri;
622	float minT, maxT;
623
624	// TODO: not very efficient to implictly cast between rays types!
625	if (mBox.GetRaySegment(*ray, minT, maxT))
626	{
627	float len = ray->Length();
628
629	if (!len)
630	len = Limits::Small;
631
632	leaf->AddRay(RayInfo(ray, minT / len, maxT / len));
633	}
634	}
635
636	mStat.rays = (int)leaf->mRays.size();
637	leaf->UpdatePvsSize();
638
639	mStat.initialPvsSize = leaf->GetPvsSize();
640
641	// Subdivide();
642	mRoot = Subdivide(TraversalData(leaf, mBox, 0));
643
644	if (mSplitCandidates)
645	{
646	// force release of this vector
647	delete mSplitCandidates;
648	mSplitCandidates = new vector<SortableEntry>;
649	}
650
651	mStat.Stop();
652
653	mStat.Print(cout);
654	Debug << "#Total memory=" << GetMemUsage() << endl;
655	}
656
657
658
659	VspKdNode *VspKdTree::Subdivide(const TraversalData &tdata)
660	{
661	VspKdNode *result = NULL;
662
663	priority_queue<TraversalData> tStack;
664	//stack<TraversalData> tStack;
665
666	tStack.push(tdata);
667
668	AxisAlignedBox3 backBox;
669	AxisAlignedBox3 frontBox;
670
671	int lastMem = 0;
672
673	while (!tStack.empty())
674	{
675	float mem = GetMemUsage();
676
677	if (lastMem / 10 != ((int)mem) / 10)
678	{
679	Debug << mem << " MB" << endl;
680	}
681	lastMem = (int)mem;
682
683	if (1 && mem > mMaxMemory)
684	{
685	Debug << "memory limit reached: " << mem << endl;
686	// count statistics on unprocessed leafs
687	while (!tStack.empty())
688	{
689	EvaluateLeafStats(tStack.top());
690	tStack.pop();
691	}
692	break;
693	}
694
695	TraversalData data = tStack.top();
696	tStack.pop();
697
698	VspKdNode node = SubdivideNode(dynamic_cast<VspKdLeaf >(data.mNode),
699	data.mBox, backBox, frontBox);
700	if (result == NULL)
701	result = node;
702
703	if (!node->IsLeaf())
704	{
705	VspKdInterior interior = dynamic_cast<VspKdInterior >(node);
706
707	// push the children on the stack
708	tStack.push(TraversalData(interior->GetBack(), backBox, data.mDepth + 1));
709	tStack.push(TraversalData(interior->GetFront(), frontBox, data.mDepth + 1));
710	}
711	else
712	{
713	EvaluateLeafStats(data);
714	}
715	}
716
717	return result;
718	}
719
720
721	// returns selected plane for subdivision
722	int VspKdTree::SelectPlane(VspKdLeaf *leaf,
723	const AxisAlignedBox3 &box,
724	float &position,
725	int &raysBack,
726	int &raysFront,
727	int &pvsBack,
728	int &pvsFront)
729	{
730	int axis = 0;
731	float costRatio = 0;
732
733	if (splitType == ESplitRegular)
734	{
735	costRatio = BestCostRatioRegular(leaf,
736	box,
737	axis,
738	position,
739	raysBack,
740	raysFront,
741	pvsBack,
742	pvsFront);
743	}
744	else if (splitType == ESplitHeuristic)
745	{
746	costRatio = BestCostRatioHeuristic(leaf,
747	box,
748	axis,
749	position,
750	raysBack,
751	raysFront,
752	pvsBack,
753	pvsFront);
754	}
755	else
756	{
757	cerr << "VspKdTree: Unknown split heuristics\n";
758	exit(1);
759	}
760
761	if (costRatio > mTermMaxCostRatio)
762	{
763	//Debug << "Too big cost ratio " << costRatio << endl;
764	++ leaf->mMaxCostMisses;
765	if (leaf->mMaxCostMisses > mTermMissTolerance)
766	{
767	++ mStat.maxCostNodes;
768	return -1;
769	}
770	}
771
772	if (0)
773	Debug << "pvs=" << leaf->mPvsSize
774	<< " rays=" << (int)leaf->mRays.size()
775	<< " rc=" << leaf->GetAvgRayContribution()
776	<< " axis=" << axis << endl;
777
778	return axis;
779	}
780
781
782
783	float VspKdTree::EvalCostRatio(VspKdLeaf *leaf,
784	const AxisAlignedBox3 &box,
785	const int axis,
786	const float position,
787	int &raysBack,
788	int &raysFront,
789	int &pvsBack,
790	int &pvsFront)
791	{
792	raysBack = 0;
793	raysFront = 0;
794	pvsFront = 0;
795	pvsBack = 0;
796
797	Intersectable::NewMail(3);
798
799	// eval pvs size
800	const int pvsSize = leaf->GetPvsSize();
801
802	// this is the main ray classification loop!
803	for(RayInfoContainer::iterator ri = leaf->mRays.begin();
804	ri != leaf->mRays.end(); ++ ri)
805	{
806	if (!(*ri).mRay->IsActive())
807	continue;
808
809	// determine the side of this ray with respect to the plane
810	float t;
811	int side = (*ri).ComputeRayIntersection(axis, position, t);
812
813	if (side <= 0)
814	++ raysBack;
815
816	if (side >= 0)
817	++ raysFront;
818
819	AddObject2Pvs((*ri).mRay->mTerminationObject, side, pvsBack, pvsFront);
820	}
821
822	//-- only one of these options should be one
823
824	if (0) //-- only pvs
825	{
826	const float sum = float(pvsBack + pvsFront);
827	const float oldCost = (float)pvsSize;
828
829	return sum / oldCost;
830	}
831
832	//-- pvs + probability heuristics
833	float pBack, pFront, pOverall;
834
835	if (1)
836	{
837	// box length substitute for probability
838	const float minBox = box.Min(axis);
839	const float maxBox = box.Max(axis);
840
841	pBack = position - minBox;
842	pFront = maxBox - position;
843	pOverall = maxBox - minBox;
844	}
845
846	if (0) //-- area substitute for probability
847	{
848	pOverall = box.SurfaceArea();
849
850	const bool useMidSplit = true;
851	//const bool useMidSplit = false;
852
853	if (!useMidSplit)
854	{
855	Vector3 pos = box.Max(); pos[axis] = position;
856	pBack = AxisAlignedBox3(box.Min(), pos).SurfaceArea();
857
858	pos = box.Min(); pos[axis] = position;
859	pFront = AxisAlignedBox3(pos, box.Max()).SurfaceArea();
860	}
861	else
862	{
863	//-- simplified computation for mid split
864	const int axis2 = (axis + 1) % 3;
865	const int axis3 = (axis + 2) % 3;
866
867	const float faceArea =
868	(box.Max(axis2) - box.Min(axis2)) *
869	(box.Max(axis3) - box.Min(axis3));
870
871	pBack = pFront = pOverall * 0.5f + faceArea;
872	}
873	}
874
875	//-- ray density substitute for probability
876	if (0)
877	{
878	pBack = (float)raysBack;
879	pFront = (float)raysFront;
880	pOverall = (float)leaf->mRays.size();
881	}
882
883	//Debug << axis << " " << pvsSize << " " << pvsBack << " " << pvsFront << endl;
884	//Debug << pFront << " " << pBack << " " << pOverall << endl;
885
886	// float sum = raysBack(position - minBox) + raysFront(maxBox - position);
887	const float newCost = pvsBack * pBack + pvsFront * pFront;
888	// float oldCost = leaf->mRays.size();
889	const float oldCost = (float)pvsSize * pOverall;
890
891	return (mCtDivCi + newCost) / oldCost;
892	}
893
894
895	float VspKdTree::BestCostRatioRegular(VspKdLeaf *leaf,
896	const AxisAlignedBox3 &box,
897	int &axis,
898	float &position,
899	int &raysBack,
900	int &raysFront,
901	int &pvsBack,
902	int &pvsFront)
903	{
904	int nRaysBack[3], nRaysFront[3];
905	int nPvsBack[3], nPvsFront[3];
906
907	float nPosition[3];
908	float nCostRatio[3];
909	int bestAxis = -1;
910
911	//AxisAlignedBox3 sBox = GetBBox(leaf);
912	const int sAxis = box.Size().DrivingAxis();
913
914	for (axis = 0; axis < 3; ++ axis)
915	{
916	if (!mOnlyDrivingAxis \|\| axis == sAxis)
917	{
918
919	nPosition[axis] = (box.Min()[axis] + box.Max()[axis]) * 0.5f;
920
921	nCostRatio[axis] = EvalCostRatio(leaf,
922	box,
923	axis,
924	nPosition[axis],
925	nRaysBack[axis],
926	nRaysFront[axis],
927	nPvsBack[axis],
928	nPvsFront[axis]);
929	if (bestAxis == -1)
930	{
931	bestAxis = axis;
932	}
933	else if (nCostRatio[axis] < nCostRatio[bestAxis])
934	{
935	/*Debug << "pvs front " << nPvsBack[axis]
936	<< " pvs back " << nPvsFront[axis]
937	<< " overall pvs " << leaf->GetPvsSize() << endl;*/
938	bestAxis = axis;
939	}
940
941	}
942	}
943
944	//-- assign best axis
945	axis = bestAxis;
946	position = nPosition[bestAxis];
947
948	raysBack = nRaysBack[bestAxis];
949	raysFront = nRaysFront[bestAxis];
950
951	pvsBack = nPvsBack[bestAxis];
952	pvsFront = nPvsFront[bestAxis];
953
954	return nCostRatio[bestAxis];
955	}
956
957	float VspKdTree::BestCostRatioHeuristic(VspKdLeaf *leaf,
958	const AxisAlignedBox3 &box,
959	int &axis,
960	float &position,
961	int &raysBack,
962	int &raysFront,
963	int &pvsBack,
964	int &pvsFront)
965	{
966	//AxisAlignedBox3 box = GetBBox(leaf);
967	// AxisAlignedBox3 dirBox = GetDirBBox(node);
968
969	axis = box.Size().DrivingAxis();
970
971	SortSplitCandidates(leaf, axis);
972
973	// go through the lists, count the number of objects left and right
974	// and evaluate the following cost funcion:
975	// C = ct_div_ci + (qlrl + qrrr)/queries
976
977	int rl = 0, rr = (int)leaf->mRays.size();
978	int pl = 0, pr = leaf->GetPvsSize();
979
980	float minBox = box.Min(axis);
981	float maxBox = box.Max(axis);
982	float sizeBox = maxBox - minBox;
983
984	float minBand = minBox + 0.1f*(maxBox - minBox);
985	float maxBand = minBox + 0.9f*(maxBox - minBox);
986
987	float sum = rr*sizeBox;
988	float minSum = 1e20f;
989
990	Intersectable::NewMail();
991
992	// set all object as belonging to the fron pvs
993	for(RayInfoContainer::iterator ri = leaf->mRays.begin();
994	ri != leaf->mRays.end(); ++ ri)
995	{
996	if ((*ri).mRay->IsActive())
997	{
998	Intersectable object = (ri).mRay->mTerminationObject;
999
1000	if (object)
1001	if (!object->Mailed())
1002	{
1003	object->Mail();
1004	object->mCounter = 1;
1005	}
1006	else
1007	++ object->mCounter;
1008	}
1009	}
1010
1011	Intersectable::NewMail();
1012
1013	for (vector<SortableEntry>::const_iterator ci = mSplitCandidates->begin();
1014	ci < mSplitCandidates->end(); ++ ci)
1015	{
1016	VssRay *ray;
1017
1018	switch ((*ci).type)
1019	{
1020	case SortableEntry::ERayMin:
1021	{
1022	++ rl;
1023	ray = (*ci).ray;
1024	Intersectable *object = ray->mTerminationObject;
1025	if (object && !object->Mailed())
1026	{
1027	object->Mail();
1028	++ pl;
1029	}
1030	break;
1031	}
1032	case SortableEntry::ERayMax:
1033	{
1034	-- rr;
1035
1036	ray = (*ci).ray;
1037	Intersectable *object = ray->mTerminationObject;
1038
1039	if (object)
1040	{
1041	if (-- object->mCounter == 0)
1042	-- pr;
1043	}
1044
1045	break;
1046	}
1047	}
1048
1049	if ((ci).value > minBand && (ci).value < maxBand)
1050	{
1051	sum = pl((ci).value - minBox) + pr(maxBox - (ci).value);
1052
1053	// cout<<"pos="<<(*ci).value<<"\t q=("<<ql<<","<<qr<<")\t r=("<<rl<<","<<rr<<")"<<endl;
1054	// cout<<"cost= "<<sum<<endl;
1055
1056	if (sum < minSum)
1057	{
1058	minSum = sum;
1059	position = (*ci).value;
1060
1061	raysBack = rl;
1062	raysFront = rr;
1063
1064	pvsBack = pl;
1065	pvsFront = pr;
1066
1067	}
1068	}
1069	}
1070
1071	float oldCost = (float)leaf->GetPvsSize();
1072	float newCost = mCtDivCi + minSum / sizeBox;
1073	float ratio = newCost / oldCost;
1074
1075	//Debug << "costRatio=" << ratio << " pos=" << position << " t=" << (position - minBox) / (maxBox - minBox)
1076	// <<"\t q=(" << queriesBack << "," << queriesFront << ")\t r=(" << raysBack << "," << raysFront << ")" << endl;
1077
1078	return ratio;
1079	}
1080
1081	void VspKdTree::SortSplitCandidates(VspKdLeaf *node,
1082	const int axis)
1083	{
1084	mSplitCandidates->clear();
1085
1086	int requestedSize = 2 * (int)(node->mRays.size());
1087	// creates a sorted split candidates array
1088	if (mSplitCandidates->capacity() > 500000 &&
1089	requestedSize < (int)(mSplitCandidates->capacity()/10) )
1090	{
1091	DEL_PTR(mSplitCandidates);
1092	mSplitCandidates = new vector<SortableEntry>;
1093	}
1094
1095	mSplitCandidates->reserve(requestedSize);
1096
1097	// insert all queries
1098	for(RayInfoContainer::const_iterator ri = node->mRays.begin();
1099	ri < node->mRays.end(); ++ ri)
1100	{
1101	bool positive = (*ri).mRay->HasPosDir(axis);
1102	mSplitCandidates->push_back(SortableEntry(positive ? SortableEntry::ERayMin : SortableEntry::ERayMax,
1103	(ri).ExtrapOrigin(axis), (ri).mRay));
1104
1105	mSplitCandidates->push_back(SortableEntry(positive ? SortableEntry::ERayMax : SortableEntry::ERayMin,
1106	(ri).ExtrapTermination(axis), (ri).mRay));
1107	}
1108
1109	stable_sort(mSplitCandidates->begin(), mSplitCandidates->end());
1110	}
1111
1112
1113	void VspKdTree::EvaluateLeafStats(const TraversalData &data)
1114	{
1115	// the node became a leaf -> evaluate stats for leafs
1116	VspKdLeaf leaf = dynamic_cast<VspKdLeaf >(data.mNode);
1117
1118	if (leaf->GetPvsSize() > mStat.maxPvsSize)
1119	mStat.maxPvsSize = leaf->GetPvsSize();
1120
1121	if ((int)(leaf->mRays.size()) > mStat.maxRayRefs)
1122	mStat.maxRayRefs = (int)leaf->mRays.size();
1123
1124	mStat.rays += (int)leaf->mRays.size();
1125
1126	if (data.mDepth >= mTermMaxDepth)
1127	++ mStat.maxDepthNodes;
1128
1129	if (leaf->GetPvsSize() < mTermMinPvs)
1130	++ mStat.minPvsNodes;
1131
1132	mStat.accPvsSize += leaf->GetPvsSize();
1133
1134	if ((int)leaf->GetRays().size() < mTermMinRays)
1135	++ mStat.minRaysNodes;
1136
1137	if (leaf->GetAvgRayContribution() > mTermMaxRayContribution)
1138	++ mStat.maxRayContribNodes;
1139
1140	if (SqrMagnitude(data.mBox.Size()) <= mTermMinSize)
1141	++ mStat.minSizeNodes;
1142	}
1143
1144
1145	inline bool VspKdTree::TerminationCriteriaMet(const VspKdLeaf *leaf,
1146	const AxisAlignedBox3 &box) const
1147	{
1148	return ((leaf->GetPvsSize() < mTermMinPvs) \|\|
1149	(leaf->mRays.size() < mTermMinRays) \|\|
1150	(leaf->GetAvgRayContribution() > mTermMaxRayContribution ) \|\|
1151	(leaf->mDepth >= mTermMaxDepth) \|\|
1152	(mStat.Leaves() >= mMaxViewCells) \|\|
1153	(SqrMagnitude(box.Size()) <= mTermMinSize));
1154	}
1155
1156
1157	VspKdNode VspKdTree::SubdivideNode(VspKdLeaf leaf,
1158	const AxisAlignedBox3 &box,
1159	AxisAlignedBox3 &backBBox,
1160	AxisAlignedBox3 &frontBBox)
1161	{
1162	if (TerminationCriteriaMet(leaf, box))
1163	{
1164	if (1 && leaf->mDepth >= mTermMaxDepth)
1165	{
1166	Debug << "Warning: max depth reached depth=" << (int)leaf->mDepth<<" rays=" << (int)leaf->mRays.size() << endl;
1167	Debug << "Bbox: " << GetBBox(leaf) << endl;
1168	}
1169	//Debug << "depth: " << (int)leaf->mDepth << " pvs: " << leaf->GetPvsSize() << " rays: " << leaf->mRays.size() << endl;
1170
1171	return leaf;
1172	}
1173
1174	float position;
1175	// first count ray sides
1176	int raysBack;
1177	int raysFront;
1178	int pvsBack;
1179	int pvsFront;
1180
1181	// select subdivision axis
1182	const int axis = SelectPlane(leaf, box, position, raysBack, raysFront, pvsBack, pvsFront);
1183	//Debug << "rays back=" << raysBack << " rays front=" << raysFront << " pvs back=" << pvsBack << " pvs front=" << pvsFront << endl;
1184
1185	if (axis == -1)
1186	return leaf;
1187
1188	mStat.nodes += 2;
1189	++ mStat.splits[axis];
1190
1191	// add the new nodes to the tree
1192	VspKdInterior *node =
1193	new VspKdInterior(dynamic_cast<VspKdInterior *>(leaf->mParent));
1194
1195	node->mAxis = axis;
1196	node->mPosition = position;
1197	node->mBox = box;
1198
1199	backBBox = box;
1200	frontBBox = box;
1201
1202	VspKdLeaf *back = new VspKdLeaf(node, raysBack, leaf->GetMaxCostMisses());
1203	back->SetPvsSize(pvsBack);
1204	VspKdLeaf *front = new VspKdLeaf(node, raysFront, leaf->GetMaxCostMisses());
1205	front->SetPvsSize(pvsFront);
1206
1207	// replace a link from node's parent
1208	if (leaf->mParent)
1209	{
1210	VspKdInterior parent = dynamic_cast<VspKdInterior >(leaf->mParent);
1211	parent->ReplaceChildLink(leaf, node);
1212	}
1213
1214	// and setup child links
1215	node->SetupChildLinks(back, front);
1216
1217	if (axis <= VspKdNode::SPLIT_Z)
1218	{
1219	backBBox.SetMax(axis, position);
1220	frontBBox.SetMin(axis, position);
1221
1222	for(RayInfoContainer::iterator ri = leaf->mRays.begin();
1223	ri != leaf->mRays.end(); ++ ri)
1224	{
1225	if ((*ri).mRay->IsActive())
1226	{
1227	// first unref ray from the former leaf
1228	(*ri).mRay->Unref();
1229	float t;
1230	// determine the side of this ray with respect to the plane
1231	int side = node->ComputeRayIntersection(*ri, t);
1232
1233	if (side == 0)
1234	{
1235	if ((*ri).mRay->HasPosDir(axis))
1236	{
1237	back->AddRay(RayInfo((*ri).mRay,
1238	(*ri).mMinT,
1239	t));
1240	front->AddRay(RayInfo((*ri).mRay,
1241	t,
1242	(*ri).mMaxT));
1243	}
1244	else
1245	{
1246	back->AddRay(RayInfo((*ri).mRay,
1247	t,
1248	(*ri).mMaxT));
1249	front->AddRay(RayInfo((*ri).mRay,
1250	(*ri).mMinT,
1251	t));
1252	}
1253	}
1254	else
1255	if (side == 1)
1256	front->AddRay(*ri);
1257	else
1258	back->AddRay(*ri);
1259	} else
1260	(*ri).mRay->Unref();
1261	}
1262	}
1263	else
1264	{
1265	// rays front/back
1266
1267	for (RayInfoContainer::iterator ri = leaf->mRays.begin();
1268	ri != leaf->mRays.end(); ++ ri)
1269	{
1270	if ((*ri).mRay->IsActive())
1271	{
1272	// first unref ray from the former leaf
1273	(*ri).mRay->Unref();
1274
1275	int side;
1276	if ((*ri).mRay->GetDirParametrization(axis - 3) > position)
1277	side = 1;
1278	else
1279	side = -1;
1280
1281	if (side == 1)
1282	front->AddRay(*ri);
1283	else
1284	back->AddRay(*ri);
1285	}
1286	else
1287	(*ri).mRay->Unref();
1288	}
1289	}
1290
1291	// update stats
1292	mStat.rayRefs -= (int)leaf->mRays.size();
1293	mStat.rayRefs += raysBack + raysFront;
1294
1295	DEL_PTR(leaf);
1296
1297	return node;
1298	}
1299
1300	int VspKdTree::ReleaseMemory(const int time)
1301	{
1302	stack<VspKdNode *> tstack;
1303
1304	// find a node in the tree which subtree will be collapsed
1305	int maxAccessTime = time - mAccessTimeThreshold;
1306	int released = 0;
1307
1308	tstack.push(mRoot);
1309
1310	while (!tstack.empty())
1311	{
1312	VspKdNode *node = tstack.top();
1313	tstack.pop();
1314
1315	if (!node->IsLeaf())
1316	{
1317	VspKdInterior in = dynamic_cast<VspKdInterior >(node);
1318	// cout<<"depth="<<(int)in->depth<<" time="<<in->lastAccessTime<<endl;
1319	if (in->mDepth >= mMinCollapseDepth && in->mLastAccessTime <= maxAccessTime)
1320	{
1321	released = CollapseSubtree(node, time);
1322	break;
1323	}
1324	if (in->GetBack()->GetAccessTime() < in->GetFront()->GetAccessTime())
1325	{
1326	tstack.push(in->GetFront());
1327	tstack.push(in->GetBack());
1328	}
1329	else
1330	{
1331	tstack.push(in->GetBack());
1332	tstack.push(in->GetFront());
1333	}
1334	}
1335	}
1336
1337	while (tstack.empty())
1338	{
1339	// could find node to collaps...
1340	// cout<<"Could not find a node to release "<<endl;
1341	break;
1342	}
1343
1344	return released;
1345	}
1346
1347
1348	VspKdNode VspKdTree::SubdivideLeaf(VspKdLeaf leaf,
1349	const float sizeThreshold)
1350	{
1351	VspKdNode *node = leaf;
1352
1353	AxisAlignedBox3 leafBBox = GetBBox(leaf);
1354
1355	static int pass = 0;
1356	++ pass;
1357
1358	// check if we should perform a dynamic subdivision of the leaf
1359	if (// leaf->mRays.size() > (unsigned)termMinCost &&
1360	(leaf->GetPvsSize() >= mTermMinPvs) &&
1361	(SqrMagnitude(leafBBox.Size()) > sizeThreshold) )
1362	{
1363	// memory check and release
1364	if (GetMemUsage() > mMaxTotalMemory)
1365	ReleaseMemory(pass);
1366
1367	AxisAlignedBox3 backBBox, frontBBox;
1368
1369	// subdivide the node
1370	node = SubdivideNode(leaf, leafBBox, backBBox, frontBBox);
1371	}
1372	return node;
1373	}
1374
1375
1376
1377	void VspKdTree::UpdateRays(VssRayContainer &remove,
1378	VssRayContainer &add)
1379	{
1380	VspKdLeaf::NewMail();
1381
1382	// schedule rays for removal
1383	for(VssRayContainer::const_iterator ri = remove.begin();
1384	ri != remove.end(); ++ ri)
1385	{
1386	(*ri)->ScheduleForRemoval();
1387	}
1388
1389	int inactive = 0;
1390
1391	for (VssRayContainer::const_iterator ri = remove.begin(); ri != remove.end(); ++ ri)
1392	{
1393	if ((*ri)->ScheduledForRemoval())
1394	// RemoveRay(*ri, NULL, false);
1395	// !!! BUG - with true it does not work correctly - aggreated delete
1396	RemoveRay(*ri, NULL, true);
1397	else
1398	++ inactive;
1399	}
1400
1401	// cout<<"all/inactive"<<remove.size()<<"/"<<inactive<<endl;
1402	for (VssRayContainer::const_iterator ri = add.begin(); ri != add.end(); ++ ri)
1403	{
1404	AddRay(*ri);
1405	}
1406	}
1407
1408
1409	void VspKdTree::RemoveRay(VssRay *ray,
1410	vector<VspKdLeaf > affectedLeaves,
1411	const bool removeAllScheduledRays)
1412	{
1413	stack<RayTraversalData> tstack;
1414
1415	tstack.push(RayTraversalData(mRoot, RayInfo(ray)));
1416
1417	RayTraversalData data;
1418
1419	// cout<<"Number of ray refs = "<<ray->RefCount()<<endl;
1420	while (!tstack.empty())
1421	{
1422	data = tstack.top();
1423	tstack.pop();
1424
1425	if (!data.mNode->IsLeaf())
1426	{
1427	// split the set of rays in two groups intersecting the
1428	// two subtrees
1429	TraverseInternalNode(data, tstack);
1430	}
1431	else
1432	{
1433	// remove the ray from the leaf
1434	// find the ray in the leaf and swap it with the last ray...
1435	VspKdLeaf leaf = (VspKdLeaf )data.mNode;
1436
1437	if (!leaf->Mailed())
1438	{
1439	leaf->Mail();
1440
1441	if (affectedLeaves)
1442	affectedLeaves->push_back(leaf);
1443
1444	if (removeAllScheduledRays)
1445	{
1446	int tail = (int)leaf->mRays.size() - 1;
1447
1448	for (int i=0; i < (int)(leaf->mRays.size()); ++ i)
1449	{
1450	if (leaf->mRays[i].mRay->ScheduledForRemoval())
1451	{
1452	// find a ray to replace it with
1453	while (tail >= i && leaf->mRays[tail].mRay->ScheduledForRemoval())
1454	{
1455	++ mStat.removedRayRefs;
1456	leaf->mRays[tail].mRay->Unref();
1457	leaf->mRays.pop_back();
1458
1459	-- tail;
1460	}
1461
1462	if (tail < i)
1463	break;
1464
1465	++ mStat.removedRayRefs;
1466
1467	leaf->mRays[i].mRay->Unref();
1468	leaf->mRays[i] = leaf->mRays[tail];
1469	leaf->mRays.pop_back();
1470	-- tail;
1471	}
1472	}
1473	}
1474	}
1475
1476	if (!removeAllScheduledRays)
1477	for (int i=0; i < (int)leaf->mRays.size(); i++)
1478	{
1479	if (leaf->mRays[i].mRay == ray)
1480	{
1481	++ mStat.removedRayRefs;
1482	ray->Unref();
1483	leaf->mRays[i] = leaf->mRays[leaf->mRays.size() - 1];
1484	leaf->mRays.pop_back();
1485	// check this ray again
1486	break;
1487	}
1488	}
1489	}
1490	}
1491
1492	if (ray->RefCount() != 0)
1493	{
1494	cerr << "Error: Number of remaining refs = " << ray->RefCount() << endl;
1495	exit(1);
1496	}
1497
1498	}
1499
1500	void VspKdTree::AddRay(VssRay *ray)
1501	{
1502	stack<RayTraversalData> tstack;
1503
1504	tstack.push(RayTraversalData(mRoot, RayInfo(ray)));
1505
1506	RayTraversalData data;
1507
1508	while (!tstack.empty())
1509	{
1510	data = tstack.top();
1511	tstack.pop();
1512
1513	if (!data.mNode->IsLeaf())
1514	{
1515	TraverseInternalNode(data, tstack);
1516	}
1517	else
1518	{
1519	// remove the ray from the leaf
1520	// find the ray in the leaf and swap it with the last ray
1521	VspKdLeaf leaf = dynamic_cast<VspKdLeaf >(data.mNode);
1522
1523	leaf->AddRay(data.mRayData);
1524	++ mStat.addedRayRefs;
1525	}
1526	}
1527	}
1528
1529	void VspKdTree::TraverseInternalNode(RayTraversalData &data,
1530	stack<RayTraversalData> &tstack)
1531	{
1532	VspKdInterior in = dynamic_cast<VspKdInterior >(data.mNode);
1533
1534	if (in->mAxis <= VspKdNode::SPLIT_Z)
1535	{
1536	// determine the side of this ray with respect to the plane
1537	float t;
1538	const int side = in->ComputeRayIntersection(data.mRayData, t);
1539
1540	if (side == 0)
1541	{
1542	if (data.mRayData.mRay->HasPosDir(in->mAxis))
1543	{
1544	tstack.push(RayTraversalData(in->GetBack(),
1545	RayInfo(data.mRayData.mRay, data.mRayData.mMinT, t)));
1546
1547	tstack.push(RayTraversalData(in->GetFront(),
1548	RayInfo(data.mRayData.mRay, t, data.mRayData.mMaxT)));
1549
1550	}
1551	else
1552	{
1553	tstack.push(RayTraversalData(in->GetBack(),
1554	RayInfo(data.mRayData.mRay,
1555	t,
1556	data.mRayData.mMaxT)));
1557	tstack.push(RayTraversalData(in->GetFront(),
1558	RayInfo(data.mRayData.mRay,
1559	data.mRayData.mMinT,
1560	t)));
1561	}
1562	}
1563	else
1564	if (side == 1)
1565	tstack.push(RayTraversalData(in->GetFront(), data.mRayData));
1566	else
1567	tstack.push(RayTraversalData(in->GetBack(), data.mRayData));
1568	}
1569	else
1570	{
1571	// directional split
1572	if (data.mRayData.mRay->GetDirParametrization(in->mAxis - 3) > in->mPosition)
1573	tstack.push(RayTraversalData(in->GetFront(), data.mRayData));
1574	else
1575	tstack.push(RayTraversalData(in->GetBack(), data.mRayData));
1576	}
1577	}
1578
1579
1580	int VspKdTree::CollapseSubtree(VspKdNode *sroot, const int time)
1581	{
1582	// first count all rays in the subtree
1583	// use mail 1 for this purpose
1584	stack<VspKdNode *> tstack;
1585
1586	int rayCount = 0;
1587	int totalRayCount = 0;
1588	int collapsedNodes = 0;
1589
1590	#if DEBUG_COLLAPSE
1591	cout << "Collapsing subtree" << endl;
1592	cout << "accessTime=" << sroot->GetAccessTime() << endl;
1593	cout << "depth=" << (int)sroot->depth << endl;
1594	#endif
1595
1596	// tstat.collapsedSubtrees++;
1597	// tstat.collapseDepths += (int)sroot->depth;
1598	// tstat.collapseAccessTimes += time - sroot->GetAccessTime();
1599	tstack.push(sroot);
1600	VssRay::NewMail();
1601
1602	while (!tstack.empty())
1603	{
1604	++ collapsedNodes;
1605
1606	VspKdNode *node = tstack.top();
1607	tstack.pop();
1608	if (node->IsLeaf())
1609	{
1610	VspKdLeaf leaf = (VspKdLeaf ) node;
1611
1612	for(RayInfoContainer::iterator ri = leaf->mRays.begin();
1613	ri != leaf->mRays.end(); ++ ri)
1614	{
1615	++ totalRayCount;
1616
1617	if ((ri).mRay->IsActive() && !(ri).mRay->Mailed())
1618	{
1619	(*ri).mRay->Mail();
1620	++ rayCount;
1621	}
1622	}
1623	}
1624	else
1625	{
1626	tstack.push(dynamic_cast<VspKdInterior *>(node)->GetFront());
1627	tstack.push(dynamic_cast<VspKdInterior *>(node)->GetBack());
1628	}
1629	}
1630
1631	VssRay::NewMail();
1632
1633	// create a new node that will hold the rays
1634	VspKdLeaf *newLeaf = new VspKdLeaf(sroot->mParent, rayCount);
1635
1636	VspKdInterior *parent =
1637	dynamic_cast<VspKdInterior *>(newLeaf->mParent);
1638
1639	if (parent)
1640	{
1641	parent->ReplaceChildLink(sroot, newLeaf);
1642	}
1643	tstack.push(sroot);
1644
1645	while (!tstack.empty())
1646	{
1647	VspKdNode *node = tstack.top();
1648	tstack.pop();
1649
1650	if (node->IsLeaf())
1651	{
1652	VspKdLeaf leaf = dynamic_cast<VspKdLeaf >(node);
1653
1654	for(RayInfoContainer::iterator ri = leaf->mRays.begin();
1655	ri != leaf->mRays.end(); ++ ri)
1656	{
1657	// unref this ray from the old node
1658	if ((*ri).mRay->IsActive())
1659	{
1660	(*ri).mRay->Unref();
1661	if (!(*ri).mRay->Mailed())
1662	{
1663	(*ri).mRay->Mail();
1664	newLeaf->AddRay(*ri);
1665	}
1666	}
1667	else
1668	(*ri).mRay->Unref();
1669	}
1670	}
1671	else
1672	{
1673	VspKdInterior *interior =
1674	dynamic_cast<VspKdInterior *>(node);
1675	tstack.push(interior->GetBack());
1676	tstack.push(interior->GetFront());
1677	}
1678	}
1679
1680	// delete the node and all its children
1681	DEL_PTR(sroot);
1682
1683	// for(VspKdNode::SRayContainer::iterator ri = newleaf->mRays.begin();
1684	// ri != newleaf->mRays.end(); ++ ri)
1685	// (*ri).ray->UnMail(2);
1686
1687	#if DEBUG_COLLAPSE
1688	cout << "Total memory before=" << GetMemUsage() << endl;
1689	#endif
1690
1691	mStat.nodes -= collapsedNodes - 1;
1692	mStat.rayRefs -= totalRayCount - rayCount;
1693
1694	#if DEBUG_COLLAPSE
1695	cout << "collapsed nodes" << collapsedNodes << endl;
1696	cout << "collapsed rays" << totalRayCount - rayCount << endl;
1697	cout << "Total memory after=" << GetMemUsage() << endl;
1698	cout << "================================" << endl;
1699	#endif
1700
1701	// tstat.collapsedNodes += collapsedNodes;
1702	// tstat.collapsedRays += totalRayCount - rayCount;
1703	return totalRayCount - rayCount;
1704	}
1705
1706
1707	int VspKdTree::GetPvsSize(VspKdNode *node, const AxisAlignedBox3 &box) const
1708	{
1709	stack<VspKdNode *> tstack;
1710	tstack.push(mRoot);
1711
1712	Intersectable::NewMail();
1713	int pvsSize = 0;
1714
1715	while (!tstack.empty())
1716	{
1717	VspKdNode *node = tstack.top();
1718	tstack.pop();
1719
1720	if (node->IsLeaf())
1721	{
1722	VspKdLeaf leaf = (VspKdLeaf )node;
1723	for (RayInfoContainer::iterator ri = leaf->mRays.begin();
1724	ri != leaf->mRays.end(); ++ ri)
1725	{
1726	if ((*ri).mRay->IsActive())
1727	{
1728	Intersectable *object;
1729	#if BIDIRECTIONAL_RAY
1730	object = (*ri).mRay->mOriginObject;
1731
1732	if (object && !object->Mailed())
1733	{
1734	++ pvsSize;
1735	object->Mail();
1736	}
1737	#endif
1738	object = (*ri).mRay->mTerminationObject;
1739	if (object && !object->Mailed())
1740	{
1741	++ pvsSize;
1742	object->Mail();
1743	}
1744	}
1745	}
1746	}
1747	else
1748	{
1749	VspKdInterior in = dynamic_cast<VspKdInterior >(node);
1750
1751	if (in->mAxis < 3)
1752	{
1753	if (box.Max(in->mAxis) >= in->mPosition)
1754	tstack.push(in->GetFront());
1755
1756	if (box.Min(in->mAxis) <= in->mPosition)
1757	tstack.push(in->GetBack());
1758	}
1759	else
1760	{
1761	// both nodes for directional splits
1762	tstack.push(in->GetFront());
1763	tstack.push(in->GetBack());
1764	}
1765	}
1766	}
1767
1768	return pvsSize;
1769	}
1770
1771	void VspKdTree::GetRayContributionStatistics(float &minRayContribution,
1772	float &maxRayContribution,
1773	float &avgRayContribution)
1774	{
1775	stack<VspKdNode *> tstack;
1776	tstack.push(mRoot);
1777
1778	minRayContribution = 1.0f;
1779	maxRayContribution = 0.0f;
1780
1781	float sumRayContribution = 0.0f;
1782	int leaves = 0;
1783
1784	while (!tstack.empty())
1785	{
1786	VspKdNode *node = tstack.top();
1787	tstack.pop();
1788	if (node->IsLeaf())
1789	{
1790	leaves++;
1791	VspKdLeaf leaf = (VspKdLeaf )node;
1792	float c = leaf->GetAvgRayContribution();
1793
1794	if (c > maxRayContribution)
1795	maxRayContribution = c;
1796	if (c < minRayContribution)
1797	minRayContribution = c;
1798	sumRayContribution += c;
1799
1800	}
1801	else
1802	{
1803	VspKdInterior in = (VspKdInterior )node;
1804
1805	tstack.push(in->GetFront());
1806	tstack.push(in->GetBack());
1807	}
1808	}
1809
1810	Debug << "sum=" << sumRayContribution << endl;
1811	Debug << "leaves=" << leaves << endl;
1812	avgRayContribution = sumRayContribution / (float)leaves;
1813	}
1814
1815
1816	int VspKdTree::GenerateRays(const float ratioPerLeaf,
1817	SimpleRayContainer &rays)
1818	{
1819	stack<VspKdNode *> tstack;
1820	tstack.push(mRoot);
1821
1822	while (!tstack.empty())
1823	{
1824	VspKdNode *node = tstack.top();
1825	tstack.pop();
1826
1827	if (node->IsLeaf())
1828	{
1829	VspKdLeaf leaf = dynamic_cast<VspKdLeaf >(node);
1830
1831	float c = leaf->GetAvgRayContribution();
1832	int num = (int)(c*ratioPerLeaf + 0.5);
1833	// cout<<num<<" ";
1834
1835	for (int i=0; i < num; i++)
1836	{
1837	Vector3 origin = GetBBox(leaf).GetRandomPoint();
1838	/*Vector3 dirVector = GetDirBBox(leaf).GetRandomPoint();
1839	Vector3 direction = Vector3(sin(dirVector.x), sin(dirVector.y), cos(dirVector.x));
1840	//cout<<"dir vector.x="<<dirVector.x<<"direction'.x="<<atan2(direction.x, direction.y)<<endl;
1841	rays.push_back(SimpleRay(origin, direction));*/
1842	}
1843
1844	}
1845	else
1846	{
1847	VspKdInterior *in =
1848	dynamic_cast<VspKdInterior *>(node);
1849
1850	tstack.push(in->GetFront());
1851	tstack.push(in->GetBack());
1852	}
1853	}
1854
1855	return (int)rays.size();
1856	}
1857
1858
1859	float VspKdTree::GetAvgPvsSize()
1860	{
1861	stack<VspKdNode *> tstack;
1862	tstack.push(mRoot);
1863
1864	int sumPvs = 0;
1865	int leaves = 0;
1866
1867	while (!tstack.empty())
1868	{
1869	VspKdNode *node = tstack.top();
1870	tstack.pop();
1871
1872	if (node->IsLeaf())
1873	{
1874	VspKdLeaf leaf = (VspKdLeaf )node;
1875	// update pvs size
1876	leaf->UpdatePvsSize();
1877	sumPvs += leaf->GetPvsSize();
1878	leaves++;
1879	}
1880	else
1881	{
1882	VspKdInterior in = (VspKdInterior )node;
1883
1884	tstack.push(in->GetFront());
1885	tstack.push(in->GetBack());
1886	}
1887	}
1888
1889	return sumPvs / (float)leaves;
1890	}
1891
1892	VspKdNode *VspKdTree::GetRoot() const
1893	{
1894	return mRoot;
1895	}
1896
1897	AxisAlignedBox3 VspKdTree::GetBBox(VspKdNode *node) const
1898	{
1899	if (node->mParent == NULL)
1900	return mBox;
1901
1902	if (!node->IsLeaf())
1903	{
1904	if (node->Type() == VspKdNode::EIntermediate)
1905	return (dynamic_cast<VspKdInterior *>(node))->mBox;
1906	else
1907	return (dynamic_cast<VspKdIntermediate *>(node))->mBox;
1908	}
1909
1910	VspKdInterior parent = dynamic_cast<VspKdInterior >(node->mParent);
1911
1912	if (parent->mAxis >= 3)
1913	return parent->mBox;
1914
1915	AxisAlignedBox3 box(parent->mBox);
1916	if (parent->GetFront() == node)
1917	box.SetMin(parent->mAxis, parent->mPosition);
1918	else
1919	box.SetMax(parent->mAxis, parent->mPosition);
1920
1921	return box;
1922	}
1923
1924	int VspKdTree::GetRootPvsSize() const
1925	{
1926	return GetPvsSize(mRoot, mBox);
1927	}
1928
1929	const VspKdStatistics &VspKdTree::GetStatistics() const
1930	{
1931	return mStat;
1932	}
1933
1934	void VspKdTree::AddRays(VssRayContainer &add)
1935	{
1936	VssRayContainer remove;
1937	UpdateRays(remove, add);
1938	}
1939
1940	// return memory usage in MB
1941	float VspKdTree::GetMemUsage() const
1942	{
1943	return
1944	(sizeof(VspKdTree) +
1945	mStat.Leaves() * sizeof(VspKdLeaf) +
1946	mStat.Interior() * sizeof(VspKdInterior) +
1947	mStat.rayRefs * sizeof(RayInfo)) / (1024.0f * 1024.0f);
1948	}
1949
1950	float VspKdTree::GetRayMemUsage() const
1951	{
1952	return mStat.rays * (sizeof(VssRay))/(1024.0f * 1024.0f);
1953	}
1954
1955
1956	int VspKdTree::ComputePvsSize(VspKdNode *node,
1957	const RayInfoContainer &globalRays) const
1958	{
1959	int pvsSize = 0;
1960
1961	const AxisAlignedBox3 box = GetBBox(node);
1962
1963	RayInfoContainer::const_iterator it, it_end = globalRays.end();
1964
1965	Intersectable::NewMail();
1966
1967	// warning: implicit conversion from VssRay to Ray
1968	for (it = globalRays.begin(); it != globalRays.end(); ++ it)
1969	{
1970	VssRay ray = (it).mRay;
1971
1972	// ray intersects node bounding box
1973	if (box.GetMinMaxT(*ray, NULL, NULL))
1974	{
1975	if (ray->mTerminationObject && !ray->mTerminationObject->Mailed())
1976	{
1977	ray->mTerminationObject->Mail();
1978	++ pvsSize;
1979	}
1980	}
1981	}
1982	return pvsSize;
1983	}
1984
1985
1986	void VspKdTree::CollectLeaves(vector<VspKdLeaf *> &leaves) const
1987	{
1988	stack<VspKdNode *> nodeStack;
1989	nodeStack.push(mRoot);
1990
1991	while (!nodeStack.empty())
1992	{
1993	VspKdNode *node = nodeStack.top();
1994
1995	nodeStack.pop();
1996
1997	if (node->IsLeaf())
1998	{
1999	VspKdLeaf leaf = dynamic_cast<VspKdLeaf >(node);
2000	leaves.push_back(leaf);
2001	}
2002	else
2003	{
2004	VspKdInterior interior = dynamic_cast<VspKdInterior >(node);
2005
2006	nodeStack.push(interior->GetBack());
2007	nodeStack.push(interior->GetFront());
2008	}
2009	}
2010	}
2011
2012
2013	int VspKdTree::FindNeighbors(VspKdLeaf *n,
2014	vector<VspKdLeaf *> &neighbors,
2015	bool onlyUnmailed)
2016	{
2017	stack<VspKdNode *> nodeStack;
2018
2019	nodeStack.push(mRoot);
2020
2021	AxisAlignedBox3 box = GetBBox(n);
2022
2023	while (!nodeStack.empty())
2024	{
2025	VspKdNode *node = nodeStack.top();
2026	nodeStack.pop();
2027
2028	if (node->IsLeaf())
2029	{
2030	VspKdLeaf leaf = dynamic_cast<VspKdLeaf >(node);
2031
2032	if (leaf != n && (!onlyUnmailed \|\| !leaf->Mailed()))
2033	neighbors.push_back(leaf);
2034	}
2035	else
2036	{
2037	VspKdInterior interior = dynamic_cast<VspKdInterior >(node);
2038
2039	if (interior->mPosition > box.Max(interior->mAxis))
2040	nodeStack.push(interior->mBack);
2041	else
2042	{
2043	if (interior->mPosition < box.Min(interior->mAxis))
2044	nodeStack.push(interior->mFront);
2045	else
2046	{
2047	// random decision
2048	nodeStack.push(interior->mBack);
2049	nodeStack.push(interior->mFront);
2050	}
2051	}
2052	}
2053	}
2054
2055	return (int)neighbors.size();
2056	}
2057
2058
2059	void VspKdTree::SetViewCellsManager(ViewCellsManager *vcm)
2060	{
2061	mViewCellsManager = vcm;
2062	}
2063
2064
2065	int VspKdTree::CastLineSegment(const Vector3 &origin,
2066	const Vector3 &termination,
2067	vector<ViewCell *> &viewcells)
2068	{
2069	int hits = 0;
2070
2071	float mint = 0.0f, maxt = 1.0f;
2072	const Vector3 dir = termination - origin;
2073
2074	stack<LineTraversalData> tStack;
2075
2076	Intersectable::NewMail();
2077
2078	Vector3 entp = origin;
2079	Vector3 extp = termination;
2080
2081	VspKdNode *node = mRoot;
2082	VspKdNode *farChild;
2083
2084	float position;
2085	int axis;
2086
2087	while (1)
2088	{
2089	if (!node->IsLeaf())
2090	{
2091	VspKdInterior in = dynamic_cast<VspKdInterior >(node);
2092	position = in->mPosition;
2093	axis = in->mAxis;
2094
2095	if (entp[axis] <= position)
2096	{
2097	if (extp[axis] <= position)
2098	{
2099	node = in->mBack;
2100	// cases N1,N2,N3,P5,Z2,Z3
2101	continue;
2102	} else
2103	{
2104	// case N4
2105	node = in->mBack;
2106	farChild = in->mFront;
2107	}
2108	}
2109	else
2110	{
2111	if (position <= extp[axis])
2112	{
2113	node = in->mFront;
2114	// cases P1,P2,P3,N5,Z1
2115	continue;
2116	}
2117	else
2118	{
2119	node = in->mFront;
2120	farChild = in->mBack;
2121	// case P4
2122	}
2123	}
2124
2125	// $$ modification 3.5.2004 - hints from Kamil Ghais
2126	// case N4 or P4
2127	float tdist = (position - origin[axis]) / dir[axis];
2128	tStack.push(LineTraversalData(farChild, extp, maxt)); //TODO
2129	extp = origin + dir * tdist;
2130	maxt = tdist;
2131	}
2132	else
2133	{
2134	// compute intersection with all objects in this leaf
2135	VspKdLeaf leaf = dynamic_cast<VspKdLeaf >(node);
2136	ViewCell *vc = leaf->GetViewCell();
2137
2138	if (!vc->Mailed())
2139	{
2140	vc->Mail();
2141	viewcells.push_back(vc);
2142	++ hits;
2143	}
2144
2145	if (0) //matt TODO: REMOVE LATER
2146	leaf->mRays.push_back(RayInfo(new VssRay(origin, termination, NULL, NULL, 0)));
2147
2148	// get the next node from the stack
2149	if (tStack.empty())
2150	break;
2151
2152	entp = extp;
2153	mint = maxt;
2154
2155	LineTraversalData &s = tStack.top();
2156	node = s.mNode;
2157	extp = s.mExitPoint;
2158	maxt = s.mMaxT;
2159	tStack.pop();
2160	}
2161	}
2162
2163	return hits;
2164	}
2165
2166
2167	void VspKdTree::GenerateViewCell(VspKdLeaf *leaf)
2168	{
2169	RayInfoContainer::const_iterator it, it_end = leaf->GetRays().end();
2170
2171	VspKdViewCell vc = dynamic_cast<VspKdViewCell >(mViewCellsManager->GenerateViewCell());
2172	leaf->SetViewCell(vc);
2173
2174	vc->SetVolume(GetBBox(leaf).GetVolume());
2175	vc->SetArea(GetBBox(leaf).SurfaceArea());
2176
2177	vc->mLeaf = leaf;
2178
2179	for (it = leaf->GetRays().begin(); it != it_end; ++ it)
2180	{
2181	VssRay ray = (it).mRay;
2182
2183	if (ray->mTerminationObject)
2184	vc->GetPvs().AddSample(ray->mTerminationObject, ray->mPdf);
2185
2186	if (ray->mOriginObject)
2187	vc->GetPvs().AddSample(ray->mOriginObject, ray->mPdf);
2188	}
2189	}
2190
2191
2192	bool VspKdTree::MergeViewCells(VspKdLeaf l1, VspKdLeaf l2)
2193	{
2194	//-- change pointer to view cells of all leaves associated
2195	//-- with the previous view cells
2196	VspKdViewCell *fVc = l1->GetViewCell();
2197	VspKdViewCell *bVc = l2->GetViewCell();
2198
2199	VspKdViewCell vc = dynamic_cast<VspKdViewCell >(
2200	mViewCellsManager->MergeViewCells(fVc, bVc));
2201
2202	// if merge was unsuccessful
2203	if (!vc) return false;
2204	// TODO
2205	#if VC_HISTORY
2206	// set new size of view cell
2207	vc->SetVolume(fVc->GetVolume() + bVc->GetVolume());
2208	// new area
2209	vc->SetArea(fVc->GetArea() + bVc->GetArea());
2210
2211	vector<VspKdLeaf *> fLeaves = fVc->mLeaves;
2212	vector<VspKdLeaf *> bLeaves = bVc->mLeaves;
2213
2214	vector<VspKdLeaf *>::const_iterator it;
2215
2216	//-- change view cells of all the other leaves the view cell belongs to
2217	for (it = fLeaves.begin(); it != fLeaves.end(); ++ it)
2218	{
2219	(*it)->SetViewCell(vc);
2220	vc->mLeaves.push_back(*it);
2221	}
2222
2223	for (it = bLeaves.begin(); it != bLeaves.end(); ++ it)
2224	{
2225	(*it)->SetViewCell(vc);
2226	vc->mLeaves.push_back(*it);
2227	}
2228	#endif
2229	vc->Mail();
2230
2231	// clean up old view cells
2232	DEL_PTR(fVc);
2233	DEL_PTR(bVc);
2234
2235	return true;
2236	}
2237
2238	void VspKdTree::CollectMergeCandidates()
2239	{
2240	VspKdMergeCandidate::sOverallCost = 0;
2241	vector<VspKdLeaf *> leaves;
2242
2243	// collect the leaves, e.g., the "voxels" that will build the view cells
2244	CollectLeaves(leaves);
2245
2246	VspKdLeaf::NewMail();
2247
2248	vector<VspKdLeaf *>::const_iterator it, it_end = leaves.end();
2249
2250	// generate view cells
2251	for (it = leaves.begin(); it != it_end; ++ it)
2252	GenerateViewCell(*it);
2253
2254	// find merge candidates and push them into queue
2255	for (it = leaves.begin(); it != it_end; ++ it)
2256	{
2257	VspKdLeaf leaf = it;
2258	ViewCell *vc = leaf->GetViewCell();
2259	// no leaf is part of two merge candidates
2260	leaf->Mail();
2261	VspKdMergeCandidate::sOverallCost +=
2262	vc->GetVolume() * vc->GetPvs().GetSize();
2263	vector<VspKdLeaf *> neighbors;
2264	FindNeighbors(leaf, neighbors, true);
2265
2266	vector<VspKdLeaf *>::const_iterator nit,
2267	nit_end = neighbors.end();
2268
2269	for (nit = neighbors.begin(); nit != nit_end; ++ nit)
2270	{
2271	mMergeQueue.push(VspKdMergeCandidate(leaf, *nit));
2272	}
2273	}
2274	}
2275
2276
2277	int VspKdTree::MergeViewCells(const VssRayContainer &rays)
2278	{
2279	CollectMergeCandidates();
2280
2281	int merged = 0;
2282
2283	int vcSize = mStat.Leaves();
2284	// use priority queue to merge leaves
2285	while (!mMergeQueue.empty() && (vcSize > mMergeMinViewCells) &&
2286	(mMergeQueue.top().GetMergeCost() <
2287	mMergeMaxCostRatio * VspKdMergeCandidate::sOverallCost))
2288	{
2289	//Debug << "mergecost: " << mergeQueue.top().GetMergeCost() /
2290	//VspKdMergeCandidate::sOverallCost << " " << mMergeMaxCostRatio << endl;
2291	VspKdMergeCandidate mc = mMergeQueue.top();
2292	mMergeQueue.pop();
2293
2294	// both view cells equal!
2295	if (mc.GetLeaf1()->GetViewCell() == mc.GetLeaf2()->GetViewCell())
2296	continue;
2297
2298	if (mc.Valid())
2299	{
2300	ViewCell::NewMail();
2301	MergeViewCells(mc.GetLeaf1(), mc.GetLeaf2());
2302
2303	++ merged;
2304	-- vcSize;
2305	// increase absolute merge cost
2306	VspKdMergeCandidate::sOverallCost += mc.GetMergeCost();
2307	}
2308	// merge candidate not valid, because one of the leaves was already
2309	// merged with another one
2310	else
2311	{
2312	// validate and reinsert into queue
2313	mc.SetValid();
2314	mMergeQueue.push(mc);
2315	//Debug << "validate " << mc.GetMergeCost() << endl;
2316	}
2317	}
2318
2319	Debug << "merged " << merged << " of " << mStat.Leaves() << " leaves" << endl;
2320
2321	//TODO: should return sample contributions
2322	return merged;
2323	}
2324
2325
2326	void VspKdTree::RepairViewCellsLeafLists()
2327	{
2328	// list not valid anymore => clear
2329	stack<VspKdNode *> nodeStack;
2330	nodeStack.push(mRoot);
2331
2332	ViewCell::NewMail();
2333
2334	while (!nodeStack.empty())
2335	{
2336	VspKdNode *node = nodeStack.top();
2337	nodeStack.pop();
2338
2339	if (node->IsLeaf())
2340	{
2341	VspKdLeaf leaf = dynamic_cast<VspKdLeaf >(node);
2342
2343	VspKdViewCell *viewCell = leaf->GetViewCell();
2344	#if VC_HISTORY
2345	if (!viewCell->Mailed())
2346	{
2347	viewCell->mLeaves.clear();
2348	viewCell->Mail();
2349	}
2350
2351	viewCell->mLeaves.push_back(leaf);
2352	#endif
2353	}
2354	else
2355	{
2356	VspKdInterior interior = dynamic_cast<VspKdInterior >(node);
2357
2358	nodeStack.push(interior->GetFront());
2359	nodeStack.push(interior->GetBack());
2360	}
2361	}
2362	}
2363
2364
2365	VspKdNode VspKdTree::CollapseTree(VspKdNode node, int &collapsed)
2366	{
2367	if (node->IsLeaf())
2368	return node;
2369
2370	VspKdInterior interior = dynamic_cast<VspKdInterior >(node);
2371
2372	VspKdNode *front = CollapseTree(interior->GetFront(), collapsed);
2373	VspKdNode *back = CollapseTree(interior->GetBack(), collapsed);
2374
2375	if (front->IsLeaf() && back->IsLeaf())
2376	{
2377	VspKdLeaf frontLeaf = dynamic_cast<VspKdLeaf >(front);
2378	VspKdLeaf backLeaf = dynamic_cast<VspKdLeaf >(back);
2379
2380	//-- collapse tree
2381	if (frontLeaf->GetViewCell() == backLeaf->GetViewCell())
2382	{
2383	VspKdViewCell *vc = frontLeaf->GetViewCell();
2384
2385	VspKdLeaf *leaf = new VspKdLeaf(interior->GetParent(), 0);
2386	leaf->SetViewCell(vc);
2387	//leaf->SetTreeValid(frontLeaf->TreeValid());
2388
2389	// replace a link from node's parent
2390	if (leaf->mParent)
2391	{
2392	// replace a link from node's parent
2393	VspKdInterior *parent =
2394	dynamic_cast<VspKdInterior *>(leaf->mParent);
2395	parent->ReplaceChildLink(node, leaf);
2396	}
2397	++ collapsed;
2398	delete interior;
2399
2400	return leaf;
2401	}
2402	}
2403
2404	return node;
2405	}
2406
2407
2408	int VspKdTree::CollapseTree()
2409	{
2410	int collapsed = 0;
2411	CollapseTree(mRoot, collapsed);
2412	// revalidate leaves
2413	RepairViewCellsLeafLists();
2414
2415	return collapsed;
2416	}
2417
2418
2419	int VspKdTree::RefineViewCells(const VssRayContainer &rays)
2420	{
2421	int shuffled = 0;
2422
2423	Debug << "refining " << (int)mMergeQueue.size() << " candidates " << endl;
2424	BspLeaf::NewMail();
2425
2426	// Use priority queue of remaining leaf pairs
2427	// These candidates either share the same view cells or
2428	// are border leaves which share a boundary.
2429	// We test if they can be shuffled, i.e.,
2430	// either one leaf is made part of one view cell or the other
2431	// leaf is made part of the other view cell. It is tested if the
2432	// remaining view cells are "better" than the old ones.
2433	while (!mMergeQueue.empty())
2434	{
2435	VspKdMergeCandidate mc = mMergeQueue.top();
2436	mMergeQueue.pop();
2437
2438	// both view cells equal or already shuffled
2439	if ((mc.GetLeaf1()->GetViewCell() == mc.GetLeaf2()->GetViewCell()) \|\|
2440	(mc.GetLeaf1()->Mailed()) \|\| (mc.GetLeaf2()->Mailed()))
2441	continue;
2442
2443	// candidate for shuffling
2444	const bool wasShuffled = false;
2445	//ShuffleLeaves(mc.GetLeaf1(), mc.GetLeaf2());
2446	// matt: restore
2447	//-- stats
2448	if (wasShuffled)
2449	++ shuffled;
2450	}
2451
2452	return shuffled;
2453	}
2454
2455
2456	inline int AddedPvsSize(ObjectPvs pvs1, const ObjectPvs &pvs2)
2457	{
2458	return pvs1.AddPvs(pvs2);
2459	}
2460
2461
2462	inline int SubtractedPvsSize(ObjectPvs pvs1, const ObjectPvs &pvs2)
2463	{
2464	return pvs1.SubtractPvs(pvs2);
2465	}
2466
2467
2468	float EvalShuffleCost(VspKdLeaf leaf, VspKdViewCell vc1, VspKdViewCell *vc2)
2469	{
2470	const int pvs1 = SubtractedPvsSize(vc1->GetPvs(), *leaf->mPvs);
2471	const int pvs2 = AddedPvsSize(vc2->GetPvs(), *leaf->mPvs);
2472
2473	const float vol1 = vc1->GetVolume() - leaf->mVolume;
2474	const float vol2 = vc2->GetVolume() + leaf->mVolume;
2475
2476	const float cost1 = pvs1 * vol1;
2477	const float cost2 = pvs2 * vol2;
2478
2479	return cost1 + cost2;
2480	}
2481
2482
2483	void VspKdTree::ShuffleLeaf(VspKdLeaf *leaf,
2484	VspKdViewCell *vc1,
2485	VspKdViewCell *vc2) const
2486	{
2487	// compute new pvs and area
2488	// TODO
2489	#if VC_HISTORY
2490	vc1->GetPvs().SubtractPvs(*leaf->mPvs);
2491	vc2->GetPvs().AddPvs(*leaf->mPvs);
2492
2493	vc1->SetArea(vc1->GetVolume() - leaf->mVolume);
2494	vc2->SetArea(vc2->GetVolume() + leaf->mVolume);
2495
2496	/// add to second view cell
2497	vc2->mLeaves.push_back(leaf);
2498
2499	// erase leaf from old view cell
2500	vector<VspKdLeaf *>::iterator it = vc1->mLeaves.begin();
2501
2502	for (; *it != leaf; ++ it);
2503	vc1->mLeaves.erase(it);
2504
2505	leaf->SetViewCell(vc2); // finally change view cell
2506	#endif
2507	}
2508
2509
2510	bool VspKdTree::ShuffleLeaves(VspKdLeaf leaf1, VspKdLeaf leaf2) const
2511	{
2512	VspKdViewCell *vc1 = leaf1->GetViewCell();
2513	VspKdViewCell *vc2 = leaf2->GetViewCell();
2514
2515	const float cost1 = vc1->GetPvs().GetSize() * vc1->GetArea();
2516	const float cost2 = vc2->GetPvs().GetSize() * vc2->GetArea();
2517
2518	const float oldCost = cost1 + cost2;
2519
2520	float shuffledCost1 = Limits::Infinity;
2521	float shuffledCost2 = Limits::Infinity;
2522
2523	// the view cell should not be empty after the shuffle
2524	//todo
2525	/* if (vc1->mLeaves.size() > 1)
2526	shuffledCost1 = EvalShuffleCost(leaf1, vc1, vc2);
2527	if (vc2->mLeaves.size() > 1)
2528	shuffledCost2 = EvalShuffleCost(leaf2, vc2, vc1);
2529
2530	*/ // shuffling unsuccessful
2531	if ((oldCost <= shuffledCost1) && (oldCost <= shuffledCost2))
2532	return false;
2533
2534	if (shuffledCost1 < shuffledCost2)
2535	{
2536	ShuffleLeaf(leaf1, vc1, vc2);
2537	leaf1->Mail();
2538	}
2539	else
2540	{
2541	ShuffleLeaf(leaf2, vc2, vc1);
2542	leaf2->Mail();
2543	}
2544
2545	return true;
2546	}
2547
2548
2549
2550
2551	/*********************************************************************/
2552	/* VspKdMergeCandidate implementation */
2553	/*********************************************************************/
2554
2555
2556	VspKdMergeCandidate::VspKdMergeCandidate(VspKdLeaf l1, VspKdLeaf l2):
2557	mMergeCost(0),
2558	mLeaf1(l1),
2559	mLeaf2(l2),
2560	mLeaf1Id(l1->GetViewCell()->mMailbox),
2561	mLeaf2Id(l2->GetViewCell()->mMailbox)
2562	{
2563	EvalMergeCost();
2564	}
2565
2566	float VspKdMergeCandidate::GetLeaf1Cost() const
2567	{
2568	VspKdViewCell *vc = mLeaf1->GetViewCell();
2569	return vc->GetPvs().GetSize() * vc->GetVolume();
2570	}
2571
2572	float VspKdMergeCandidate::GetLeaf2Cost() const
2573	{
2574	VspKdViewCell *vc = mLeaf2->GetViewCell();
2575	return vc->GetPvs().GetSize() * vc->GetVolume();
2576	}
2577
2578	void VspKdMergeCandidate::EvalMergeCost()
2579	{
2580	//-- compute pvs difference
2581	VspKdViewCell *vc1 = mLeaf1->GetViewCell();
2582	VspKdViewCell *vc2 = mLeaf2->GetViewCell();
2583
2584	const int diff1 = vc1->GetPvs().Diff(vc2->GetPvs());
2585	const int vcPvs = diff1 + vc1->GetPvs().GetSize();
2586
2587	//-- compute ratio of old cost
2588	//-- (added size of left and right view cell times pvs size)
2589	//-- to new rendering cost (size of merged view cell times pvs size)
2590	const float oldCost = GetLeaf1Cost() + GetLeaf2Cost();
2591
2592	const float newCost =
2593	(float)vcPvs * (vc1->GetVolume() + vc2->GetVolume());
2594
2595	mMergeCost = newCost - oldCost;
2596
2597	// if (vcPvs > sMaxPvsSize) // strong penalty if pvs size too large
2598	// mMergeCost += 1.0;
2599	}
2600
2601	void VspKdMergeCandidate::SetLeaf1(VspKdLeaf *l)
2602	{
2603	mLeaf1 = l;
2604	}
2605
2606	void VspKdMergeCandidate::SetLeaf2(VspKdLeaf *l)
2607	{
2608	mLeaf2 = l;
2609	}
2610
2611
2612	VspKdLeaf *VspKdMergeCandidate::GetLeaf1()
2613	{
2614	return mLeaf1;
2615	}
2616
2617
2618	VspKdLeaf *VspKdMergeCandidate::GetLeaf2()
2619	{
2620	return mLeaf2;
2621	}
2622
2623
2624	bool VspKdMergeCandidate::Valid() const
2625	{
2626	return
2627	(mLeaf1->GetViewCell()->mMailbox == mLeaf1Id) &&
2628	(mLeaf2->GetViewCell()->mMailbox == mLeaf2Id);
2629	}
2630
2631
2632	float VspKdMergeCandidate::GetMergeCost() const
2633	{
2634	return mMergeCost;
2635	}
2636
2637
2638	void VspKdMergeCandidate::SetValid()
2639	{
2640	mLeaf1Id = mLeaf1->GetViewCell()->mMailbox;
2641	mLeaf2Id = mLeaf2->GetViewCell()->mMailbox;
2642
2643	EvalMergeCost();
2644	}

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

Download in other formats: