Context Navigation

source: trunk/VUT/GtpVisibilityPreprocessor/src/VspKdTree.cpp @ 582

Revision 582, 60.3 KB checked in by mattausch, 18 years ago (diff)
fixed bug in mergueue to find root of merge and sort out doube view cells

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

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

Download in other formats: