--- a/src/kdtree.cc	Thu Nov 22 21:46:09 2007 +0100
+++ b/src/kdtree.cc	Fri Nov 23 01:24:33 2007 +0100
@@ -22,18 +22,24 @@
 
 void KdNode::subdivide(BBox bbox, int depth)
 {
+	if (depth >= 10 || shapes.size() <= 2)
+	{
+		leaf = true;
+		return;
+	}
+
 	// choose split axis
 	axis = 0;
-	if (bbox.R.y - bbox.L.y > bbox.R.x - bbox.L.x)
+	if (bbox.h() > bbox.w() && bbox.h() > bbox.d())
 		axis = 1;
-	if (bbox.R.z - bbox.L.z > bbox.R.y - bbox.L.y)
+	if (bbox.d() > bbox.w() && bbox.d() > bbox.h())
 		axis = 2;
 
 	// *** find optimal split position
 	SortableShapeList sslist(shapes, axis);
 	sort(sslist.begin(), sslist.end());
 
-	SplitList splitlist = SplitList();
+	SplitList splitlist;
 
 	SortableShapeList::iterator sh;
 	for (sh = sslist.begin(); sh != sslist.end(); sh++)
@@ -42,6 +48,8 @@
 		splitlist.push_back(SplitPos(sh->bbox.R.cell[axis]));
 	}
 	sort(splitlist.begin(), splitlist.end());
+	// remove duplicate splits
+	splitlist.resize(unique(splitlist.begin(), splitlist.end()) - splitlist.begin());
 
 	// find all posible splits
 	SplitList::iterator spl;
@@ -51,7 +59,10 @@
 		for (sh = sslist.begin(), rest = sslist.size(); sh != sslist.end(); sh++, rest--)
 		{
 			if (sh->hasMark())
+			{
+				spl->lnum++;
 				continue;
+			}
 			// if shape is completely contained in split plane
 			if (spl->pos == sh->bbox.L.cell[axis] == sh->bbox.R.cell[axis])
 			{
@@ -69,12 +80,13 @@
 					else
 						spl->lnum++;
 				}
+				sh->setMark();
 			} else
 			// if shape is on left side of split plane
 			if (sh->bbox.R.cell[axis] <= spl->pos)
 			{
+				spl->lnum++;
 				sh->setMark();
-				spl->lnum++;
 			} else
 			// if shape occupies both sides of split plane
 			if (sh->bbox.L.cell[axis] < spl->pos && sh->bbox.R.cell[axis] > spl->pos)
@@ -105,8 +117,8 @@
 		lbb.R.cell[axis] = spl->pos;
 		rbb.L.cell[axis] = spl->pos;
 		float SAL = 2*(lbb.w()*lbb.h() + lbb.w()*lbb.d() + lbb.h()*lbb.d());
-        float SAR = 2*(rbb.w()*rbb.h() + rbb.w()*rbb.d() + rbb.h()*rbb.d());
-        float splitcost = K + SAL/SAV*(K+spl->lnum) + SAR/SAV*(K+spl->rnum);
+		float SAR = 2*(rbb.w()*rbb.h() + rbb.w()*rbb.d() + rbb.h()*rbb.d());
+		float splitcost = K + SAL/SAV*(K+spl->lnum) + SAR/SAV*(K+spl->rnum);
 
 		if (splitcost < cost)
 		{
@@ -119,6 +131,29 @@
 	if (leaf)
 		return;
 
+#if 1
+// export kd-tree as .obj for visualization
+// this would be hard to reconstruct later
+	static ofstream F("kdtree.obj");
+	Vector3 v;
+	static int f=0;
+	v.cell[axis] = splitpos->pos;
+	v.cell[(axis+1)%3] = bbox.L.cell[(axis+1)%3];
+	v.cell[(axis+2)%3] = bbox.L.cell[(axis+2)%3];
+	F << "v " << v.x << " " << v.y << " " << v.z << endl;
+	v.cell[(axis+1)%3] = bbox.L.cell[(axis+1)%3];
+	v.cell[(axis+2)%3] = bbox.R.cell[(axis+2)%3];
+	F << "v " << v.x << " " << v.y << " " << v.z << endl;
+	v.cell[(axis+1)%3] = bbox.R.cell[(axis+1)%3];
+	v.cell[(axis+2)%3] = bbox.R.cell[(axis+2)%3];
+	F << "v " << v.x << " " << v.y << " " << v.z << endl;
+	v.cell[(axis+1)%3] = bbox.R.cell[(axis+1)%3];
+	v.cell[(axis+2)%3] = bbox.L.cell[(axis+2)%3];
+	F << "v " << v.x << " " << v.y << " " << v.z << endl;
+	F << "f " << f+1 << " " << f+2 << " " << f+3 << " " << f+4 << endl;
+	f += 4;
+#endif
+
 	split = splitpos->pos;
 	float lnum = splitpos->lnum;
 	float rnum = splitpos->rnum;
@@ -184,9 +219,28 @@
 	children[1].subdivide(rbb, depth+1);
 }
 
-void KdTree::build()
+void KdTree::optimize()
 {
 	root = new KdNode();
 	root->shapes = shapes;
 	root->subdivide(bbox, 0);
+	built = true;
 }
+
+void KdTree::save(ostream &str, KdNode *node)
+{
+	if (!built)
+		return;
+	if (node == NULL)
+		node = root;
+	if (node->isLeaf())
+		str << "(leaf: " << node->shapes.size() << " shapes)";
+	else
+	{
+		str << "(split " << (char)('x'+node->getAxis()) << " at " << node->getSplit() << "; L=";
+		save(str, node->getLeftChild());
+		str << "; R=";
+		save(str, node->getRightChild());
+		str << ";)";
+	}
+}