Merge branch 'master' into boost_filesystem

Conflicts:
	tests/CMakeLists.txt

10 files changed, 226 insertions, 87 deletions

diff --git a/src/CSGTermEvaluator.cc b/src/CSGTermEvaluator.cc
index 1aedfec..fc76d56 100644
--- a/src/CSGTermEvaluator.cc
+++ b/src/CSGTermEvaluator.cc
@@ -48,11 +48,11 @@ void CSGTermEvaluator::applyToChildren(const AbstractNode &node, CSGTermEvaluato
 			t1 = t2;
 		} else if (t2 && t1) {
 			if (op == CSGT_UNION) {
-				t1.reset(new CSGTerm(CSGTerm::TYPE_UNION, t1, t2));
+				t1 = CSGTerm::createCSGTerm(CSGTerm::TYPE_UNION, t1, t2);
 			} else if (op == CSGT_DIFFERENCE) {
-				t1.reset(new CSGTerm(CSGTerm::TYPE_DIFFERENCE, t1, t2));
+				t1 = CSGTerm::createCSGTerm(CSGTerm::TYPE_DIFFERENCE, t1, t2);
 			} else if (op == CSGT_INTERSECTION) {
-				t1.reset(new CSGTerm(CSGTerm::TYPE_INTERSECTION, t1, t2));
+				t1 = CSGTerm::createCSGTerm(CSGTerm::TYPE_INTERSECTION, t1, t2);
 			}
 		}
 	}
diff --git a/src/PolySetCGALEvaluator.cc b/src/PolySetCGALEvaluator.cc
index 7a9566b..3285b46 100644
--- a/src/PolySetCGALEvaluator.cc
+++ b/src/PolySetCGALEvaluator.cc
@@ -396,7 +396,7 @@ PolySet *PolySetCGALEvaluator::evaluatePolySet(const RenderNode &node)
 	CGAL_Nef_polyhedron N = this->cgalevaluator.evaluateCGALMesh(node);
 	PolySet *ps = NULL;
 	if (!N.empty()) {
-		if (!N.p3->is_simple()) {
+		if (N.dim == 3 && !N.p3->is_simple()) {
 			PRINTF("WARNING: Body of render() isn't valid 2-manifold!");
 		}
 		else {
diff --git a/src/context.cc b/src/context.cc
index 176bf8d..df884de 100644
--- a/src/context.cc
+++ b/src/context.cc
@@ -187,7 +187,7 @@ void register_builtin(Context &ctx)
 	ctx.functions_p = &Builtins::instance()->functions();
 	ctx.modules_p = &Builtins::instance()->modules();
 	ctx.set_variable("$fn", Value(0.0));
-	ctx.set_variable("$fs", Value(1.0));
+	ctx.set_variable("$fs", Value(2.0));
 	ctx.set_variable("$fa", Value(12.0));
 	ctx.set_variable("$t", Value(0.0));
 	
diff --git a/src/csgterm.cc b/src/csgterm.cc
index b21a20c..56fcbb5 100644
--- a/src/csgterm.cc
+++ b/src/csgterm.cc
@@ -26,6 +26,7 @@
 
 #include "csgterm.h"
 #include "polyset.h"
+#include "linalg.h"
 #include <sstream>
 
 /*!
@@ -47,101 +48,187 @@
  */
 
 
+shared_ptr<CSGTerm> CSGTerm::createCSGTerm(type_e type, shared_ptr<CSGTerm> left, shared_ptr<CSGTerm> right)
+{
+	if (type != TYPE_PRIMITIVE) {
+		// In case we're creating a CSG terms from a pruned tree, left/right can be NULL
+		if (!right) {
+			if (type == TYPE_UNION || type == TYPE_DIFFERENCE) return left;
+			else return right;
+		}
+		if (!left) {
+			if (type == TYPE_UNION) return right;
+			else return left;
+		}
+	}
+
+  // Pruning the tree. For details, see:
+  // http://www.cc.gatech.edu/~turk/my_papers/pxpl_csg.pdf
+	const BoundingBox &leftbox = left->getBoundingBox();
+	const BoundingBox &rightbox = right->getBoundingBox();
+	if (type == TYPE_INTERSECTION) {
+		BoundingBox newbox(leftbox.min().cwise().max(rightbox.min()),
+											 leftbox.max().cwise().min(rightbox.max()));
+		if (newbox.isNull()) {
+			return shared_ptr<CSGTerm>(); // Prune entire product
+		}
+	}
+	else if (type == TYPE_DIFFERENCE) {
+		BoundingBox newbox(leftbox.min().cwise().max(rightbox.min()),
+											 leftbox.max().cwise().min(rightbox.max()));
+		if (newbox.isNull()) {
+			return left; // Prune the negative component
+		}
+	}
+
+	return shared_ptr<CSGTerm>(new CSGTerm(type, left, right));
+}
+
+shared_ptr<CSGTerm> CSGTerm::createCSGTerm(type_e type, CSGTerm *left, CSGTerm *right)
+{
+	return createCSGTerm(type, shared_ptr<CSGTerm>(left), shared_ptr<CSGTerm>(right));
+}
+
 CSGTerm::CSGTerm(const shared_ptr<PolySet> &polyset, const Transform3d &matrix, const double color[4], const std::string &label)
-	: type(TYPE_PRIMITIVE), polyset(polyset), label(label)
+	: type(TYPE_PRIMITIVE), polyset(polyset), label(label), m(matrix)
 {
-	this->m = matrix;
 	for (int i = 0; i < 4; i++) this->color[i] = color[i];
+	initBoundingBox();
 }
 
 CSGTerm::CSGTerm(type_e type, shared_ptr<CSGTerm> left, shared_ptr<CSGTerm> right)
-	: type(type), left(left), right(right)
+	: type(type), left(left), right(right), m(Transform3d::Identity())
 {
+	initBoundingBox();
 }
 
 CSGTerm::CSGTerm(type_e type, CSGTerm *left, CSGTerm *right)
-	: type(type), left(left), right(right)
+	: type(type), left(left), right(right), m(Transform3d::Identity())
 {
+	initBoundingBox();
 }
 
 CSGTerm::~CSGTerm()
 {
 }
 
+void CSGTerm::initBoundingBox()
+{
+	if (this->type == TYPE_PRIMITIVE) {
+		this->bbox = this->m * this->polyset->getBoundingBox();
+	}
+	else {
+		const BoundingBox &leftbox = this->left->getBoundingBox();
+		const BoundingBox &rightbox = this->right->getBoundingBox();
+		switch (this->type) {
+		case TYPE_UNION:
+			this->bbox = this->m * BoundingBox(leftbox.min().cwise().min(rightbox.min()), 
+																				 leftbox.max().cwise().max(rightbox.max()));
+			break;
+		case TYPE_INTERSECTION:
+			this->bbox = this->m * BoundingBox(leftbox.min().cwise().max(rightbox.min()),
+																				 leftbox.max().cwise().min(rightbox.max()));
+			break;
+		case TYPE_DIFFERENCE:
+			this->bbox = this->m * leftbox;
+			break;
+		case TYPE_PRIMITIVE:
+			break;
+		default:
+			assert(false);
+		}
+	}
+}
 
-shared_ptr<CSGTerm> CSGTerm::normalize(shared_ptr<CSGTerm> &term)
+shared_ptr<CSGTerm> CSGTerm::normalize(shared_ptr<CSGTerm> term)
 {
 	// This function implements the CSG normalization
-	// Reference: Florian Kirsch, Juergen Doeller,
-	// OpenCSG: A Library for Image-Based CSG Rendering,
-	// University of Potsdam, Hasso-Plattner-Institute, Germany
-	// http://www.opencsg.org/data/csg_freenix2005_paper.pdf
-
-	if (term->type == TYPE_PRIMITIVE) return term;
-
-	shared_ptr<CSGTerm> x = normalize(term->left);
-	shared_ptr<CSGTerm> y = normalize(term->right);
-
-	shared_ptr<CSGTerm> t1(term);
-	if (x != term->left || y != term->right) t1.reset(new CSGTerm(term->type, x, y));
+  // Reference:
+	// Goldfeather, J., Molnar, S., Turk, G., and Fuchs, H. Near
+	// Realtime CSG Rendering Using Tree Normalization and Geometric
+	// Pruning. IEEE Computer Graphics and Applications, 9(3):20-28,
+	// 1989.
+  // http://www.cc.gatech.edu/~turk/my_papers/pxpl_csg.pdf
+
+	if (term->type == TYPE_PRIMITIVE) {
+		return term;
+	}
 
-	shared_ptr<CSGTerm> t2;
-	while (1) {
-		t2 = normalize_tail(t1);
-		if (t1 == t2)	break;
-		t1 = t2;
+	do {
+		while (term && normalize_tail(term)) { }
+		if (!term || term->type == TYPE_PRIMITIVE) return term;
+		term->left = normalize(term->left);
+	} while (term->type != TYPE_UNION &&
+					 (term->right->type != TYPE_PRIMITIVE || term->left->type == TYPE_UNION));
+	term->right = normalize(term->right);
+
+	// FIXME: Do we need to take into account any transformation of item here?
+	if (!term->right) {
+		if (term->type == TYPE_UNION || term->type == TYPE_DIFFERENCE) return term->left;
+		else return term->right;
+	}
+	if (!term->left) {
+		if (term->type == TYPE_UNION) return term->right;
+		else return term->left;
 	}
 
-	return t1;
+	return term;
 }
 
-shared_ptr<CSGTerm> CSGTerm::normalize_tail(shared_ptr<CSGTerm> &term)
+bool CSGTerm::normalize_tail(shared_ptr<CSGTerm> &term)
 {
+	if (term->type == TYPE_UNION || term->type == TYPE_PRIMITIVE) return false;
+
 	// Part A: The 'x . (y . z)' expressions
 
 	shared_ptr<CSGTerm> x = term->left;
 	shared_ptr<CSGTerm> y = term->right->left;
 	shared_ptr<CSGTerm> z = term->right->right;
 
-	CSGTerm *result = NULL;
+	shared_ptr<CSGTerm> result = term;
 
 	// 1.  x - (y + z) -> (x - y) - z
 	if (term->type == TYPE_DIFFERENCE && term->right->type == TYPE_UNION) {
-		result = new CSGTerm(TYPE_DIFFERENCE, 
-												 shared_ptr<CSGTerm>(new CSGTerm(TYPE_DIFFERENCE, x, y)),
+		term = createCSGTerm(TYPE_DIFFERENCE, 
+												 createCSGTerm(TYPE_DIFFERENCE, x, y),
 												 z);
+		return true;
 	}
 	// 2.  x * (y + z) -> (x * y) + (x * z)
 	else if (term->type == TYPE_INTERSECTION && term->right->type == TYPE_UNION) {
-		result = new CSGTerm(TYPE_UNION, 
-												 new CSGTerm(TYPE_INTERSECTION, x, y), 
-												 new CSGTerm(TYPE_INTERSECTION, x, z));
+		term = createCSGTerm(TYPE_UNION, 
+												 createCSGTerm(TYPE_INTERSECTION, x, y), 
+												 createCSGTerm(TYPE_INTERSECTION, x, z));
+		return true;
 	}
 	// 3.  x - (y * z) -> (x - y) + (x - z)
 	else if (term->type == TYPE_DIFFERENCE && term->right->type == TYPE_INTERSECTION) {
-		result = new CSGTerm(TYPE_UNION, 
-												 new CSGTerm(TYPE_DIFFERENCE, x, y), 
-												 new CSGTerm(TYPE_DIFFERENCE, x, z));
+		term = createCSGTerm(TYPE_UNION, 
+												 createCSGTerm(TYPE_DIFFERENCE, x, y), 
+												 createCSGTerm(TYPE_DIFFERENCE, x, z));
+		return true;
 	}
 	// 4.  x * (y * z) -> (x * y) * z
 	else if (term->type == TYPE_INTERSECTION && term->right->type == TYPE_INTERSECTION) {
-		result = new CSGTerm(TYPE_INTERSECTION, 
-												 shared_ptr<CSGTerm>(new CSGTerm(TYPE_INTERSECTION, x, y)),
+		term = createCSGTerm(TYPE_INTERSECTION, 
+												 createCSGTerm(TYPE_INTERSECTION, x, y),
 												 z);
+		return true;
 	}
 	// 5.  x - (y - z) -> (x - y) + (x * z)
 	else if (term->type == TYPE_DIFFERENCE && term->right->type == TYPE_DIFFERENCE) {
-		result = new CSGTerm(TYPE_UNION, 
-												 new CSGTerm(TYPE_DIFFERENCE, x, y), 
-												 new CSGTerm(TYPE_INTERSECTION, x, z));
+		term = createCSGTerm(TYPE_UNION, 
+												 createCSGTerm(TYPE_DIFFERENCE, x, y), 
+												 createCSGTerm(TYPE_INTERSECTION, x, z));
+		return true;
 	}
 	// 6.  x * (y - z) -> (x * y) - z
 	else if (term->type == TYPE_INTERSECTION && term->right->type == TYPE_DIFFERENCE) {
-		result = new CSGTerm(TYPE_DIFFERENCE, 
-												 shared_ptr<CSGTerm>(new CSGTerm(TYPE_INTERSECTION, x, y)),
+		term = createCSGTerm(TYPE_DIFFERENCE, 
+												 createCSGTerm(TYPE_INTERSECTION, x, y),
 												 z);
+		return true;
 	}
-	if (result) return shared_ptr<CSGTerm>(result);
 
 	// Part B: The '(x . y) . z' expressions
 
@@ -151,26 +238,27 @@ shared_ptr<CSGTerm> CSGTerm::normalize_tail(shared_ptr<CSGTerm> &term)
 
 	// 7. (x - y) * z  -> (x * z) - y
 	if (term->left->type == TYPE_DIFFERENCE && term->type == TYPE_INTERSECTION) {
-		result = new CSGTerm(TYPE_DIFFERENCE, 
-												 shared_ptr<CSGTerm>(new CSGTerm(TYPE_INTERSECTION, x, z)), 
+		term = createCSGTerm(TYPE_DIFFERENCE, 
+												 createCSGTerm(TYPE_INTERSECTION, x, z), 
 												 y);
+		return true;
 	}
 	// 8. (x + y) - z  -> (x - z) + (y - z)
 	else if (term->left->type == TYPE_UNION && term->type == TYPE_DIFFERENCE) {
-		result = new CSGTerm(TYPE_UNION, 
-												 new CSGTerm(TYPE_DIFFERENCE, x, z), 
-												 new CSGTerm(TYPE_DIFFERENCE, y, z));
+		term = createCSGTerm(TYPE_UNION, 
+												 createCSGTerm(TYPE_DIFFERENCE, x, z), 
+												 createCSGTerm(TYPE_DIFFERENCE, y, z));
+		return true;
 	}
 	// 9. (x + y) * z  -> (x * z) + (y * z)
 	else if (term->left->type == TYPE_UNION && term->type == TYPE_INTERSECTION) {
-		result = new CSGTerm(TYPE_UNION, 
-												 new CSGTerm(TYPE_INTERSECTION, x, z), 
-												 new CSGTerm(TYPE_INTERSECTION, y, z));
+		term = createCSGTerm(TYPE_UNION, 
+												 createCSGTerm(TYPE_INTERSECTION, x, z), 
+												 createCSGTerm(TYPE_INTERSECTION, y, z));
+		return true;
 	}
 
-	if (result) return shared_ptr<CSGTerm>(result);
-	
-	return term;
+	return false;
 }
 
 std::string CSGTerm::dump()
@@ -239,11 +327,7 @@ BoundingBox CSGChain::getBoundingBox() const
 		if (types[i] != CSGTerm::TYPE_DIFFERENCE) {
 			BoundingBox psbox = polysets[i]->getBoundingBox();
 			if (!psbox.isNull()) {
-				Eigen::Transform3d t;
-				// Column-major vs. Row-major
-				t = matrices[i];
-				bbox.extend(t * psbox.min());
-				bbox.extend(t * psbox.max());
+				bbox.extend(matrices[i] * psbox);
 			}
 		}
 	}
diff --git a/src/csgterm.h b/src/csgterm.h
index 1895839..4930349 100644
--- a/src/csgterm.h
+++ b/src/csgterm.h
@@ -18,23 +18,35 @@ public:
 		TYPE_DIFFERENCE
 	};
 
+	static shared_ptr<CSGTerm> createCSGTerm(type_e type, shared_ptr<CSGTerm> left, shared_ptr<CSGTerm> right);
+	static shared_ptr<CSGTerm> createCSGTerm(type_e type, CSGTerm *left, CSGTerm *right);
+
 	type_e type;
 	shared_ptr<PolySet> polyset;
 	std::string label;
 	shared_ptr<CSGTerm> left;
 	shared_ptr<CSGTerm> right;
-	Transform3d m;
-	double color[4];
+	BoundingBox bbox;
 
 	CSGTerm(const shared_ptr<PolySet> &polyset, const Transform3d &matrix, const double color[4], const std::string &label);
-	CSGTerm(type_e type, shared_ptr<CSGTerm> left, shared_ptr<CSGTerm> right);
-	CSGTerm(type_e type, CSGTerm *left, CSGTerm *right);
 	~CSGTerm();
 
-	static shared_ptr<CSGTerm> normalize(shared_ptr<CSGTerm> &term);
-	static shared_ptr<CSGTerm> normalize_tail(shared_ptr<CSGTerm> &term);
+	const BoundingBox &getBoundingBox() const { return this->bbox; }
+
+	static shared_ptr<CSGTerm> normalize(shared_ptr<CSGTerm> term);
+	static bool normalize_tail(shared_ptr<CSGTerm> &term);
 
 	std::string dump();
+private:
+	CSGTerm(type_e type, shared_ptr<CSGTerm> left, shared_ptr<CSGTerm> right);
+	CSGTerm(type_e type, CSGTerm *left, CSGTerm *right);
+
+	void initBoundingBox();
+
+	Transform3d m;
+	double color[4];
+
+	friend class CSGChain;
 };
 
 class CSGChain
diff --git a/src/export.cc b/src/export.cc
index 6c427dd..99bce98 100644
--- a/src/export.cc
+++ b/src/export.cc
@@ -84,17 +84,27 @@ void export_stl(CGAL_Nef_polyhedron *root_N, std::ostream &output, QProgressDial
 			stream << x3 << " " << y3 << " " << z3;
 			std::string vs3 = stream.str();
 			if (vs1 != vs2 && vs1 != vs3 && vs2 != vs3) {
-				
-				double nx = (y1-y2)*(z1-z3) - (z1-z2)*(y1-y3);
-				double ny = (z1-z2)*(x1-x3) - (x1-x2)*(z1-z3);
-				double nz = (x1-x2)*(y1-y3) - (y1-y2)*(x1-x3);
+				// The above condition ensures that vs1-vs2, vs1-vs3, and their cross
+				// product are non-zero. Floating point arithmetic may however truncate
+				// small values to 0. This can be avoided by first scaling the components
+				// of vs1-vs2 and vs1-vs3. This has no effect on the resulting unit
+				// normal vector.
+				double dn[6] = { x1-x2, y1-y2, z1-z2, x1-x3, y1-y3, z1-z3 };
+				double maxdn = 0;
+				int i;
+				for (i = 0; i < 6; ++i) {
+					double dx = dn[i];
+					if (dx < 0) dx = -dx;
+					if (dx > maxdn) maxdn = dx;
+				}
+				for (i = 0; i < 6; ++i) dn[i] /= maxdn;
+				double nx = dn[1]*dn[5] - dn[2]*dn[4];
+				double ny = dn[2]*dn[3] - dn[0]*dn[5];
+				double nz = dn[0]*dn[4] - dn[1]*dn[3];
 				double nlength = sqrt(nx*nx + ny*ny + nz*nz);
-				// Avoid generating normals for polygons with zero area
-				double eps = 0.000001;
-				if (nlength < eps) nlength = 1.0;
-				output << "  facet normal " 
-							 << nx / nlength << " " 
-							 << ny / nlength << " " 
+				output << "  facet normal "
+							 << nx / nlength << " "
+							 << ny / nlength << " "
 							 << nz / nlength << "\n";
 				output << "    outer loop\n";
 				output << "      vertex " << vs1 << "\n";
@@ -164,7 +174,7 @@ void export_dxf(CGAL_Nef_polyhedron *root_N, std::ostream &output, QProgressDial
 						 << y2 << "\n";
 		}
 	}
-	
+
 	output << "  0\n"
 				 << "ENDSEC\n";
 
diff --git a/src/glview.cc b/src/glview.cc
index c96fe01..d9f6bb5 100644
--- a/src/glview.cc
+++ b/src/glview.cc
@@ -28,6 +28,7 @@
 #include "Preferences.h"
 #include "renderer.h"
 #include "rendersettings.h"
+#include "linalg.h"
 
 #include <QApplication>
 #include <QWheelEvent>
@@ -388,12 +389,12 @@ void GLView::paintGL()
 
 	gluLookAt(0.0, -viewer_distance, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 1.0);
 
-	glTranslated(object_trans_x, object_trans_y, object_trans_z);
-
 	glRotated(object_rot_x, 1.0, 0.0, 0.0);
 	glRotated(object_rot_y, 0.0, 1.0, 0.0);
 	glRotated(object_rot_z, 0.0, 0.0, 1.0);
 
+	glTranslated(object_trans_x, object_trans_y, object_trans_z);
+
   // FIXME: Crosshairs and axes are lighted, this doesn't make sense and causes them
   // to change color based on view orientation.
 	if (showcrosshairs)
@@ -500,6 +501,7 @@ void GLView::paintGL()
 		// FIXME: This was an attempt to keep contrast with background, but is suboptimal
 		// (e.g. nearly invisible against a gray background).
 		int r,g,b;
+		r=g=b=0;
 //		bgcol.getRgb(&r, &g, &b);
 		glColor3d((255.0-r)/255.0, (255.0-g)/255.0, (255.0-b)/255.0);
 		glBegin(GL_LINES);
@@ -557,7 +559,6 @@ void GLView::mousePressEvent(QMouseEvent *event)
 	setFocus();
 }
 
-
 void GLView::normalizeAngle(GLdouble& angle)
 {
 	while(angle < 0)
@@ -594,8 +595,37 @@ void GLView::mouseMoveEvent(QMouseEvent *event)
 			if ((QApplication::keyboardModifiers() & Qt::ShiftModifier) != 0) {
 				viewer_distance += (GLdouble)dy;
 			} else {
-				object_trans_x += dx;
-				object_trans_z -= dy;
+
+      double mx = +(dx) * viewer_distance/1000;
+      double my = -(dy) * viewer_distance/1000;
+
+			Matrix3d aax, aay, aaz, tm3;
+			aax = Eigen::AngleAxisd(-(object_rot_x/180) * M_PI, Vector3d::UnitX());
+			aay = Eigen::AngleAxisd(-(object_rot_y/180) * M_PI, Vector3d::UnitY());
+			aaz = Eigen::AngleAxisd(-(object_rot_z/180) * M_PI, Vector3d::UnitZ());
+			tm3 = Matrix3d::Identity();
+			tm3 = aaz * (aay * (aax * tm3));
+
+			Matrix4d tm;
+			tm = Matrix4d::Identity();
+			for (int i=0;i<3;i++) for (int j=0;j<3;j++) tm(j,i)=tm3(j,i);
+
+			Matrix4d vec;
+			vec <<
+        0,  0,  0,  mx,
+        0,  0,  0,  0,
+        0,  0,  0,  my,
+        0,  0,  0,  1
+			;
+			if ((QApplication::keyboardModifiers() & Qt::ShiftModifier) != 0) {
+        vec(0,3) = 0;
+        vec(1,3) = my;
+        vec(2,3) = 0;
+      }
+      tm = tm * vec;
+      object_trans_x += tm(0,3);
+      object_trans_y += tm(1,3);
+      object_trans_z += tm(2,3);
 			}
 		}
 		updateGL();
diff --git a/src/linalg.h b/src/linalg.h
index e20d8d8..c1a14d1 100644
--- a/src/linalg.h
+++ b/src/linalg.h
@@ -10,6 +10,9 @@ using Eigen::Vector3d;
 typedef Eigen::AlignedBox<double, 3> BoundingBox;
 using Eigen::Matrix3f;
 using Eigen::Matrix3d;
+using Eigen::Matrix4d;
 using Eigen::Transform3d;
 
+BoundingBox operator*(const Transform3d &m, const BoundingBox &box);
+
 #endif
diff --git a/src/primitives.cc b/src/primitives.cc
index b3fa45f..466a0c7 100644
--- a/src/primitives.cc
+++ b/src/primitives.cc
@@ -245,7 +245,7 @@ int get_fragments_from_r(double r, double fn, double fs, double fa)
 	if (r < GRID_FINE) return 0;
 	if (fn > 0.0)
 		return (int)fn;
-	return (int)ceil(fmax(fmin(360.0 / fa, r*M_PI / fs), 5));
+	return (int)ceil(fmax(fmin(360.0 / fa, r*2*M_PI / fs), 5));
 }
 
 struct point2d {
diff --git a/src/surface.cc b/src/surface.cc
index 2af09dc..fe1c6aa 100644
--- a/src/surface.cc
+++ b/src/surface.cc
@@ -142,8 +142,8 @@ PolySet *SurfaceNode::evaluate_polyset(class PolySetEvaluator *) const
 
 	p->convexity = convexity;
 
-	double ox = center ? -columns/2.0 : 0;
-	double oy = center ? -lines/2.0 : 0;
+	double ox = center ? -(columns-1)/2.0 : 0;
+	double oy = center ? -(lines-1)/2.0 : 0;
 
 	for (int i = 1; i < lines; i++)
 	for (int j = 1; j < columns; j++)