From 0516aacfe231dcdebefeee785decddbc457101c8 Mon Sep 17 00:00:00 2001
From: Wenzel Jakob <wenzel@cs.cornell.edu>
Date: Wed, 9 Feb 2011 02:13:52 +0100
Subject: [PATCH] added left-balanced kd-tree construction as another testcase

---
 include/mitsuba/core/kdtree.h | 32 ++++++++++++++++++++++++++++++++
 src/tests/test_kd.cpp         | 12 +++++++-----
 2 files changed, 39 insertions(+), 5 deletions(-)

diff --git a/include/mitsuba/core/kdtree.h b/include/mitsuba/core/kdtree.h
index 00a847bc..2914b486 100644
--- a/include/mitsuba/core/kdtree.h
+++ b/include/mitsuba/core/kdtree.h
@@ -97,6 +97,9 @@ public:
 	enum EHeuristic {
 		/// Create a balanced tree by splitting along the median
 		EBalanced = 0,
+		
+		/// Create a left-balanced tree
+		ELeftBalanced,
 
 		/**
 		 * \brief Use the sliding midpoint tree construction rule. This 
@@ -365,10 +368,39 @@ protected:
 			case EBalanced: {
 					/* Split along the median */
 					split = rangeStart + (rangeEnd-rangeStart)/2;
+					axis = m_aabb.getLargestAxis();
 					std::nth_element(rangeStart, split, rangeEnd, CoordinateOrdering(axis));
 				};
 				break;
 
+			case ELeftBalanced: {
+					size_t treeSize = rangeEnd-rangeStart;
+					/* Layer 0 contains one node */
+					size_t p = 1;
+
+					/* Traverse downwards until the first incompletely
+					   filled tree level is encountered */
+					while (2*p <= treeSize)
+						p *= 2;
+
+					/* Calculate the number of filled slots in the last level */
+					size_t remaining = treeSize - p + 1;
+
+					if (2*remaining < p) {
+						/* Case 2: The last level contains too few nodes. Remove
+						   overestimate from the left subtree node count and add
+						   the remaining nodes */
+						p = (p >> 1) + remaining;
+					}
+
+					axis = m_aabb.getLargestAxis();
+					
+					split = rangeStart + (p - 1);
+					std::nth_element(rangeStart, split, rangeEnd,
+						CoordinateOrdering(axis));
+				};
+				break;
+
 			case ESlidingMidpoint: {
 					/* Sliding midpoint rule: find a split that is close to the spatial median */
 					axis = m_aabb.getLargestAxis();
diff --git a/src/tests/test_kd.cpp b/src/tests/test_kd.cpp
index e0993d1e..2e6d34e7 100644
--- a/src/tests/test_kd.cpp
+++ b/src/tests/test_kd.cpp
@@ -130,10 +130,10 @@ public:
 	
 	void test03_pointKDTree() {
 		typedef TKDTree< BasicKDNode<Point2, Float> > KDTree2;
-		size_t nPoints = 20000;
+		size_t nPoints = 50000, nTries = 20;
 		ref<Random> random = new Random();
 
-		for (int heuristic=0; heuristic<3; ++heuristic) {
+		for (int heuristic=0; heuristic<4; ++heuristic) {
 			KDTree2 kdtree(nPoints, (KDTree2::EHeuristic) heuristic);
 
 			for (size_t i=0; i<nPoints; ++i) {
@@ -146,8 +146,10 @@ public:
 			if (heuristic == 0) {
 				Log(EInfo, "Testing the balanced kd-tree construction heuristic");
 			} else if (heuristic == 1) {
-				Log(EInfo, "Testing the sliding midpoint kd-tree construction heuristic");
+				Log(EInfo, "Testing the left-balanced kd-tree construction heuristic");
 			} else if (heuristic == 2) {
+				Log(EInfo, "Testing the sliding midpoint kd-tree construction heuristic");
+			} else if (heuristic == 3) {
 				Log(EInfo, "Testing the voxel volume kd-tree construction heuristic");
 			}
 
@@ -156,7 +158,7 @@ public:
 			Log(EInfo, "Construction time = %i ms, depth = %i", timer->getMilliseconds(), kdtree.getDepth());
 
 			for (int k=1; k<=10; ++k) {
-				size_t nTraversals = 0, nTries = 10;
+				size_t nTraversals = 0;
 				for (size_t it = 0; it < nTries; ++it) {
 					Point2 p(random->nextFloat(), random->nextFloat());
 					nTraversals += kdtree.nnSearch(p, k, results);
@@ -170,7 +172,7 @@ public:
 				}
 				Log(EInfo, "Average number of traversals for a %i-nn query = " SIZE_T_FMT, k, nTraversals / nTries);
 			}
-			size_t nTraversals = 0, nTries = 10;
+			size_t nTraversals = 0;
 			for (size_t it = 0; it < nTries; ++it) {
 				Point2 p(random->nextFloat(), random->nextFloat());
 				nTraversals += kdtree.search(p, 0.05, results);