faster sutherland-hodgman iteration, added missing partial specialization for integer division involving vectors

include/mitsuba/core/vector.h

@@ -98,10 +98,8 @@ public:
 	/// Divide the vector by the given scalar
 	TVector2 operator/(T f) const {
 #ifdef MTS_DEBUG
-		if (f == 0) {
+		if (f == 0)
 			SLog(EWarn, "Vector2: Division by zero!");
-			//exit(-1);
-		}
 #endif
 		T recip = (T) 1 / f;
 		return TVector2(x * recip, y * recip);
@@ -109,10 +107,8 @@ public:
 
 	TVector2 &operator/=(T f) {
 #ifdef MTS_DEBUG
-		if (f == 0) {
+		if (f == 0)
 			SLog(EWarn, "Vector2: Division by zero!");
-			//exit(-1);
-		}
 #endif
 		T recip = (T) 1 / f;
 		x *= recip; y *= recip; 
@@ -184,6 +180,25 @@ template <typename T> inline TVector2<T> normalize(const TVector2<T> &v) {
 	return v / v.length();
 }
 
+template <> inline TVector2<int> TVector2<int>::operator/(int s) const {
+#ifdef MTS_DEBUG
+	if (s == 0) 
+		SLog(EWarn, "Vector2i: Division by zero!");
+#endif
+	return TVector2(x/s, y/s);
+}
+
+template <> inline TVector2<int> &TVector2<int>::operator/=(int s) {
+#ifdef MTS_DEBUG
+	if (s == 0) 
+		SLog(EWarn, "Vector2i: Division by zero!");
+#endif
+
+	x /= s;
+	y /= s;
+	return *this;
+}
+
 /**
  * \headerfile mitsuba/core/vector.h mitsuba/mitsuba.h
  * \brief Parameterizable three-dimensional vector data structure
@@ -260,10 +275,8 @@ public:
 	/// Divide the vector by the given scalar
 	TVector3 operator/(T f) const {
 #ifdef MTS_DEBUG
-		if (f == 0) {
+		if (f == 0)
 			SLog(EWarn, "Vector3: Division by zero!");
-			//exit(-1);
-		}
 #endif
 		T recip = (T) 1 / f;
 		return TVector3(x * recip, y * recip, z * recip);
@@ -271,10 +284,8 @@ public:
 
 	TVector3 &operator/=(T f) {
 #ifdef MTS_DEBUG
-		if (f == 0) {
+		if (f == 0)
 			SLog(EWarn, "Vector3: Division by zero!");
-			//exit(-1);
-		}
 #endif
 		T recip = (T) 1 / f;
 		x *= recip; y *= recip; z *= recip;
@@ -356,6 +367,27 @@ template <typename T> inline TVector3<T> normalize(const TVector3<T> &v) {
 	return v / v.length();
 }
 
+template <> inline TVector3<int> TVector3<int>::operator/(int s) const {
+#ifdef MTS_DEBUG
+	if (s == 0) 
+		SLog(EWarn, "Vector3i: Division by zero!");
+#endif
+	return TVector3(x/s, y/s, z/s);
+}
+
+template <> inline TVector3<int> &TVector3<int>::operator/=(int s) {
+#ifdef MTS_DEBUG
+	if (s == 0) 
+		SLog(EWarn, "Vector3i: Division by zero!");
+#endif
+
+	x /= s;
+	y /= s;
+	z /= s;
+	return *this;
+}
+
+
 /**
  * \headerfile mitsuba/core/vector.h mitsuba/mitsuba.h
  * \brief Parameterizable four-dimensional vector data structure
@@ -433,10 +465,8 @@ public:
 	/// Divide the vector by the given scalar
 	TVector4 operator/(T f) const {
 #ifdef MTS_DEBUG
-		if (f == 0) {
+		if (f == 0)
 			SLog(EWarn, "Vector4: Division by zero!");
-			//exit(-1);
-		}
 #endif
 		T recip = (T) 1 / f;
 		return TVector4(x * recip, y * recip, z * recip, w * recip);
@@ -444,10 +474,8 @@ public:
 
 	TVector4 &operator/=(T f) {
 #ifdef MTS_DEBUG
-		if (f == 0) {
+		if (f == 0)
 			SLog(EWarn, "Vector4: Division by zero!");
-			//exit(-1);
-		}
 #endif
 		T recip = (T) 1 / f;
 		x *= recip; y *= recip; z *= recip; w *= recip;
@@ -521,6 +549,27 @@ template <typename T> inline TVector4<T> normalize(const TVector4<T> &v) {
 	return v / v.length();
 }
 
+template <> inline TVector4<int> TVector4<int>::operator/(int s) const {
+#ifdef MTS_DEBUG
+	if (s == 0) 
+		SLog(EWarn, "Vector4i: Division by zero!");
+#endif
+	return TVector4(x/s, y/s, z/s, w/s);
+}
+
+template <> inline TVector4<int> &TVector4<int>::operator/=(int s) {
+#ifdef MTS_DEBUG
+	if (s == 0) 
+		SLog(EWarn, "Vector4i: Division by zero!");
+#endif
+
+	x /= s;
+	y /= s;
+	z /= s;
+	w /= s;
+	return *this;
+}
+
 MTS_NAMESPACE_END
 
 #endif /* __VECTOR_H */

src/libcore/triangle.cpp

@@ -96,65 +96,70 @@ AABB Triangle::getAABB(const Vertex *buffer) const {
 	return aabb;
 }
 
-inline void sutherlandHodgman(std::vector<Point> &vertices, int axis, 
+#define MAX_VERTS 10
-	Float splitPos, bool isMinimum) {
-	int vertexCount = (int) vertices.size();
-	if (vertexCount < 2)
-		return;
 
-	Point cur = vertices[0];
+static int sutherlandHodgman(Point *input, int inCount, Point *output, int axis, 
-	Float sign = isMinimum ? 1.0f : -1.0f, 
+		Float splitPos, bool isMinimum) {
-	      distance = sign * (cur[axis] - splitPos);
+	if (inCount < 3)
-	bool curIsInside = (distance >= 0);
+		return 0;
 
-	for (int i=0; i<vertexCount; i++) {
+	Point cur         = input[0];
-		Point next = vertices[(i+1)%vertexCount];
+	Float sign        = isMinimum ? 1.0f : -1.0f;
+	Float distance    = sign * (cur[axis] - splitPos);
+	bool  curIsInside = (distance >= 0);
+	int   outCount    = 0;
+
+	for (int i=0; i<inCount; ++i) {
+		Point next = input[(i+1)%inCount];
 		distance = sign * (next[axis] - splitPos);
 		bool nextIsInside = (distance >= 0);
 
 		if (curIsInside && nextIsInside) {
 			/* Both this and the next vertex are inside, add to the list */
-			vertices.push_back(next);
+			SAssertEx(outCount + 1 < MAX_VERTS, "Overflow in sutherlandHodgman()!");
+			output[outCount++] = next;
 		} else if (curIsInside && !nextIsInside) {
 			/* Going outside -- add the intersection */
 			Float t = (splitPos - cur[axis]) / (next[axis] - cur[axis]);
-			vertices.push_back(cur + (next - cur) * t);
+			SAssertEx(outCount + 1 < MAX_VERTS, "Overflow in sutherlandHodgman()!");
+			output[outCount++] = cur + (next - cur) * t;
 		} else if (!curIsInside && nextIsInside) {
 			/* Coming back inside -- add the intersection + next vertex */
 			Float t = (splitPos - cur[axis]) / (next[axis] - cur[axis]);
-			vertices.push_back(cur + (next - cur) * t);
+			SAssertEx(outCount + 2 < MAX_VERTS, "Overflow in sutherlandHodgman()!");
-			vertices.push_back(next);
+			output[outCount++] = cur + (next - cur) * t;
+			output[outCount++] = next;
 		} else {
 			/* Entirely outside - do not add anything */
 		}
 		cur = next;
 		curIsInside = nextIsInside;
 	}
-	vertices.erase(vertices.begin(), vertices.begin() + vertexCount);
+	return outCount;
 }
 
 AABB Triangle::getClippedAABB(const Vertex *buffer, const AABB &aabb) const {
-	std::vector<Point> vertices;
 	/* Reserve room for some additional vertices */
-	vertices.reserve(8);
+	Point vertices1[MAX_VERTS], vertices2[MAX_VERTS];
+	int nVertices = 3;
+
 	for (int i=0; i<3; ++i)
-		vertices.push_back(buffer[idx[i]].v);
+		vertices1[i] = buffer[idx[i]].v;
 
 	for (int axis=0; axis<3; ++axis) {
-		sutherlandHodgman(vertices, axis, aabb.min[axis], true);
+		nVertices = sutherlandHodgman(vertices1, nVertices, vertices2, axis, aabb.min[axis], true);
-		sutherlandHodgman(vertices, axis, aabb.max[axis], false);
+		nVertices = sutherlandHodgman(vertices2, nVertices, vertices1, axis, aabb.max[axis], false);
 	}
 
 	AABB result;
-	for (unsigned int i=0; i<vertices.size(); ++i)
+	for (int i=0; i<nVertices; ++i)
-		result.expandBy(vertices[i]);
+		result.expandBy(vertices1[i]);
 
 	/* Cover up some numerical imprecisions */
-
+	for (int i=0; i<3; ++i) {
-	for (int i=0; i<3; ++i)
 		result.min[i] -= Epsilon * std::abs(result.min[i]);
-	for (int i=0; i<3; ++i)
 		result.max[i] += Epsilon * std::abs(result.max[i]);
+	}
 
 	result.clip(aabb);