Tools

-----
More improvements to boolean tools.  The main change (although very little in
code) is changing fuzzy comparisons of floating point values.  For testing, a
new define is added in intern/boolop/intern/BOP_MathUtils.h called
VAR_EPSILON, which enables better comparisons.  This is turned on by default;
undefining it will revert to using the previous comparisons.  The downside of
these new comparisons is a loss in speed, but the resulting meshes are more
likely to be manifold (although still not always).

The other changes include speed improvements based on profiling results and
fixes for the improper creation of triangular faces with only two vertices.

8 files changed, 194 insertions, 136 deletions

diff --git a/intern/boolop/intern/BOP_Face.cpp b/intern/boolop/intern/BOP_Face.cpp
index 80c917f2838..ebe34237d4f 100644
--- a/intern/boolop/intern/BOP_Face.cpp
+++ b/intern/boolop/intern/BOP_Face.cpp
@@ -44,6 +44,8 @@ BOP_Face::BOP_Face(MT_Plane3 plane, BOP_Index originalFace)
 	m_plane        = plane;
 	m_tag          = UNCLASSIFIED;
 	m_originalFace = originalFace;
+	m_split        = 0;
+	m_bbox         = NULL;
 }
 
 /**
@@ -197,6 +199,14 @@ bool BOP_Face3::getNextVertex(BOP_Index v, BOP_Index &w)
  */
 void BOP_Face3::replaceVertexIndex(BOP_Index oldIndex, BOP_Index newIndex)
 {
+	/* if the old index really exists, and new index also exists already,
+	 * don't create an edge with both vertices == newIndex */
+
+	if( (m_indexs[0] == oldIndex || m_indexs[1] == oldIndex || m_indexs[2] == oldIndex) &&
+			(m_indexs[0] == newIndex || m_indexs[1] == newIndex || m_indexs[2] == newIndex) ) {
+		setTAG(BROKEN);
+	}
+
 	if (m_indexs[0] == oldIndex) m_indexs[0] = newIndex;
 	else if (m_indexs[1] == oldIndex) m_indexs[1] = newIndex;
 	else if (m_indexs[2] == oldIndex) m_indexs[2] = newIndex;
diff --git a/intern/boolop/intern/BOP_Face.h b/intern/boolop/intern/BOP_Face.h
index 7db5ab1fe5c..1d854ec00ca 100644
--- a/intern/boolop/intern/BOP_Face.h
+++ b/intern/boolop/intern/BOP_Face.h
@@ -34,6 +34,7 @@
 #include "BOP_Tag.h"
 #include "MT_Plane3.h"
 #include "BOP_Indexs.h"
+#include "BOP_BBox.h"
 #include <iostream>
 #include <vector>
 using namespace std;
@@ -53,10 +54,12 @@ private:
 protected:
 	BOP_Index    m_indexs[4];
 	unsigned int m_size;
+	unsigned int m_split;
+	BOP_BBox     *m_bbox;
 
 public:
 	BOP_Face(MT_Plane3 plane, BOP_Index originalFace);
-	virtual ~BOP_Face(){};
+	virtual ~BOP_Face(){if (m_bbox) delete m_bbox;};
 	inline MT_Plane3 getPlane() const {return m_plane;};
 	inline void setPlane(const MT_Plane3 plane) {m_plane = plane;};
 	inline BOP_TAG getTAG() const {return m_tag;};
@@ -65,7 +68,15 @@ public:
 	inline void setOriginalFace(const BOP_Index originalFace) {m_originalFace=originalFace;};
 	inline BOP_Index getVertex(unsigned int i) const {return m_indexs[i];};
 	inline void setVertex(const BOP_Index idx, const BOP_Index i) {m_indexs[idx]=i;};
+	inline unsigned int getSplit() const {return m_split;};
+	inline void setSplit(const unsigned int i) {m_split=i;};
+
 	void invert();
+	inline void setBBox(const MT_Point3& p1,const MT_Point3& p2,const MT_Point3& p3) {
+    m_bbox = new BOP_BBox(p1, p2, p3);};
+	inline BOP_BBox *getBBox() {return m_bbox;};
+	inline void freeBBox(){if (m_bbox!=NULL) {delete m_bbox; m_bbox=NULL;} };
+
 	inline unsigned int size() const {return m_size;};
 	
 	virtual bool getEdgeIndex(BOP_Index v1, BOP_Index v2, unsigned int &e) = 0;
diff --git a/intern/boolop/intern/BOP_Face2Face.cpp b/intern/boolop/intern/BOP_Face2Face.cpp
index 9c8ee9f7787..ef67e5dd24b 100644
--- a/intern/boolop/intern/BOP_Face2Face.cpp
+++ b/intern/boolop/intern/BOP_Face2Face.cpp
@@ -147,9 +147,18 @@ void BOP_mergeVertexs(BOP_Mesh *mesh, unsigned int firstFace);
  * @param mesh mesh that contains the faces, edges and vertices
  * @param facesA set of faces from object A
  * @param facesB set of faces from object B
+ *
+ * Two optimizations were added here:
+ *   1) keep the bounding box for a face once it's created; this is
+ *      especially important for B faces, since they were being created and
+ *      recreated over and over
+ *   2) associate a "split" index in the faceB vector with each A face; when
+ *      an A face is split, we will not need to recheck any B faces have
+ *      already been checked against that original A face
  */
+
 void BOP_Face2Face(BOP_Mesh *mesh, BOP_Faces *facesA, BOP_Faces *facesB)
-{		
+{
 	for(unsigned int idxFaceA=0;idxFaceA<facesA->size();idxFaceA++) {
 		BOP_Face *faceA = (*facesA)[idxFaceA];
 		MT_Plane3 planeA = faceA->getPlane();
@@ -157,23 +166,33 @@ void BOP_Face2Face(BOP_Mesh *mesh, BOP_Faces *facesA, BOP_Faces *facesB)
 		MT_Point3 p2 = mesh->getVertex(faceA->getVertex(1))->getPoint();
 		MT_Point3 p3 = mesh->getVertex(faceA->getVertex(2))->getPoint();
 
-		BOP_BBox boxA(p1,p2,p3);
-	
-		for(unsigned int idxFaceB=0;
+	/* get (or create) bounding box for face A */
+		if( faceA->getBBox() == NULL )
+        	faceA->setBBox(p1,p2,p3);
+		BOP_BBox *boxA = faceA->getBBox();
+
+	/* start checking B faces with the previously stored split index */
+
+		for(unsigned int idxFaceB=faceA->getSplit();
 			idxFaceB<facesB->size() && (faceA->getTAG() != BROKEN) && (faceA->getTAG() != PHANTOM);) {
 			BOP_Face *faceB = (*facesB)[idxFaceB];
+			faceA->setSplit(idxFaceB);
 			if ((faceB->getTAG() != BROKEN) && (faceB->getTAG() != PHANTOM)) {
-			  BOP_BBox boxB(mesh->getVertex(faceB->getVertex(0))->getPoint(),
-					mesh->getVertex(faceB->getVertex(1))->getPoint(),
-					mesh->getVertex(faceB->getVertex(2))->getPoint());
-			  if (boxA.intersect(boxB)) {
 
+	/* get (or create) bounding box for face B */
+				if( faceB->getBBox() == NULL )
+        			faceB->setBBox(mesh->getVertex(faceB->getVertex(0))->getPoint(),
+                    mesh->getVertex(faceB->getVertex(1))->getPoint(),
+                    mesh->getVertex(faceB->getVertex(2))->getPoint());
+			  BOP_BBox *boxB = faceB->getBBox();
+
+			  if (boxA->intersect(*boxB)) {
 			    MT_Plane3 planeB = faceB->getPlane();
 			    if (BOP_containsPoint(planeB,p1) && 
 				BOP_containsPoint(planeB,p2) && 
 				BOP_containsPoint(planeB,p3)) {
 			      if (BOP_orientation(planeB,planeA)>0) {
-				BOP_intersectCoplanarFaces(mesh,facesB,faceA,faceB,false);
+				    BOP_intersectCoplanarFaces(mesh,facesB,faceA,faceB,false);
 			      }
 			    }
 			    else {
@@ -181,10 +200,7 @@ void BOP_Face2Face(BOP_Mesh *mesh, BOP_Faces *facesA, BOP_Faces *facesB)
 			    }
 			  }			  
 			}
-			if (faceB->getTAG()==BROKEN){
-			  facesB->erase(facesB->begin()+idxFaceB);
-			}else
-			  idxFaceB++;
+			idxFaceB++;
 		}
 	}
 	
@@ -388,56 +404,11 @@ void BOP_mergeVertexs(BOP_Mesh *mesh, unsigned int firstFace)
 	for(unsigned int idxFace = firstFace; idxFace < numFaces; idxFace++) {
 		BOP_Face *face = mesh->getFace(idxFace);
 		if ((face->getTAG() != BROKEN) && (face->getTAG() != PHANTOM)) {
-			BOP_Index v1 = face->getVertex(0);
-			BOP_Index v2 = face->getVertex(1);
-			BOP_Index v3 = face->getVertex(2);
-			MT_Point3 vertex1 = mesh->getVertex(v1)->getPoint();
-			MT_Point3 vertex2 = mesh->getVertex(v2)->getPoint();
-			MT_Point3 vertex3 = mesh->getVertex(v3)->getPoint();
-			int dist12 = BOP_comp(vertex1,vertex2);
-			int dist13 = BOP_comp(vertex1,vertex3);
-			int dist23 = BOP_comp(vertex2,vertex3);
-
-			if (dist12 == 0) {
-				if (dist13 == 0) {
-					// v1 ~= v2 , v1 ~= v3 , v2 ~= v3
-					mesh->replaceVertexIndex(v2,v1);
-					mesh->replaceVertexIndex(v3,v1);		
-				}
-				else {
-					if (dist23 == 0) {
-						mesh->replaceVertexIndex(v1,v2);
-						mesh->replaceVertexIndex(v3,v2);
-					}
-					else {
-						mesh->replaceVertexIndex(v1,v2);
-					}
-				}
-			}
-			else {
-				if (dist13 == 0) {
-					// v1 ~= v3
-					if (dist23 == 0) {
-						mesh->replaceVertexIndex(v1,v3);
-						mesh->replaceVertexIndex(v2,v3);
-					}
-					else {
-						mesh->replaceVertexIndex(v1,v3);
-					}
-				}
-				else {
-					if (dist23 == 0) {
-						// v2 ~= v3
-						mesh->replaceVertexIndex(v2,v3);
-					} else {
-						// all differents
-						if (BOP_collinear(vertex1,vertex2,vertex3)) {
-							// collinear triangle 
-							face->setTAG(PHANTOM);
-						}
-					}
-				}
-			}
+			MT_Point3 vertex1 = mesh->getVertex(face->getVertex(0))->getPoint();
+			MT_Point3 vertex2 = mesh->getVertex(face->getVertex(1))->getPoint();
+			MT_Point3 vertex3 = mesh->getVertex(face->getVertex(2))->getPoint();
+			if (BOP_collinear(vertex1,vertex2,vertex3)) // collinear triangle 
+				face->setTAG(PHANTOM);
 		}
 	}
 }
@@ -510,7 +481,7 @@ void BOP_mergeSort(MT_Point3 *points, unsigned int *face, unsigned int &size, bo
 	if (size == 2) {
 
 		// Trivial case, only test the merge ...
-		if (BOP_comp(0,points[0].distance(points[1]))==0) {
+		if (BOP_fuzzyZero(points[0].distance(points[1]))) {
 			face[0] = 3;
 			size--;
 		}
@@ -595,8 +566,8 @@ void BOP_mergeSort(MT_Point3 *points, unsigned int *face, unsigned int &size, bo
 		// Merge data
 		MT_Scalar d1 = sortedPoints[1].distance(sortedPoints[0]);
 		MT_Scalar d2 = sortedPoints[1].distance(sortedPoints[2]);
-		if (BOP_comp(0,d1)==0 && sortedFaces[1] != sortedFaces[0]) {
-			if (BOP_comp(0,d2)==0 && sortedFaces[1] != sortedFaces[2])  {
+		if (BOP_fuzzyZero(d1) && sortedFaces[1] != sortedFaces[0]) {
+			if (BOP_fuzzyZero(d2) && sortedFaces[1] != sortedFaces[2])  {
 				if (d1 < d2) {
 					// merge 0 and 1
 					sortedFaces[0] = 3;
@@ -608,7 +579,7 @@ void BOP_mergeSort(MT_Point3 *points, unsigned int *face, unsigned int &size, bo
 					if (size == 3) {
 						// merge 1 and 2 ???
 						d1 = sortedPoints[1].distance(sortedPoints[2]);
-						if (BOP_comp(0,d1)==0 && sortedFaces[1] != sortedFaces[2])  {
+						if (BOP_fuzzyZero(d1) && sortedFaces[1] != sortedFaces[2])  {
 							// merge!
 							sortedFaces[1] = 3;
 							size--;
@@ -636,7 +607,7 @@ void BOP_mergeSort(MT_Point3 *points, unsigned int *face, unsigned int &size, bo
 				if (size == 3) {
 					// merge 1 i 2 ???
 					d1 = sortedPoints[1].distance(sortedPoints[2]);
-					if (BOP_comp(0,d1)==0 && sortedFaces[1] != sortedFaces[2])  {
+					if (BOP_fuzzyZero(d1) && sortedFaces[1] != sortedFaces[2])  {
 						// merge!
 						sortedFaces[1] = 3;
 						size--;
@@ -645,7 +616,7 @@ void BOP_mergeSort(MT_Point3 *points, unsigned int *face, unsigned int &size, bo
 			}     
 		}
 		else {
-			if (BOP_comp(0,d2)==0 && sortedFaces[1] != sortedFaces[2])  {
+			if (BOP_fuzzyZero(d2) && sortedFaces[1] != sortedFaces[2])  {
 				// merge 1 and 2
 				sortedFaces[1] = 3;
 				for(i = 2; i<size-1;i++) {
@@ -656,7 +627,7 @@ void BOP_mergeSort(MT_Point3 *points, unsigned int *face, unsigned int &size, bo
 			}
 			else if (size == 4) {
 				d1 = sortedPoints[2].distance(sortedPoints[3]);
-				if (BOP_comp(0,d1)==0 && sortedFaces[2] != sortedFaces[3])  {
+				if (BOP_fuzzyZero(d1) && sortedFaces[2] != sortedFaces[3])  {
 					// merge 2 and 3
 					sortedFaces[2] = 3;
 					size--;
diff --git a/intern/boolop/intern/BOP_MathUtils.cpp b/intern/boolop/intern/BOP_MathUtils.cpp
index 251bbb9e138..e0d96b465ee 100644
--- a/intern/boolop/intern/BOP_MathUtils.cpp
+++ b/intern/boolop/intern/BOP_MathUtils.cpp
@@ -43,9 +43,38 @@ using namespace std;
  */
 int BOP_comp(const MT_Scalar A, const MT_Scalar B)
 {
+#ifndef VAR_EPSILON
 	if (A >= B + BOP_EPSILON) return 1;
 	else if (B >= A + BOP_EPSILON) return -1;
 	else return 0; 
+#else
+	int expA, expB;
+	float mant;
+	frexp(A, &expA);	/* get exponents of each number */
+	frexp(B, &expB);
+
+	if(expA < expB)		/* find the larger exponent */
+		expA = expB;
+	mant = frexp((A-B), &expB);	/* get exponent of the difference */
+	/* mantissa will only be zero is (A-B) is really zero; otherwise, also
+	 * also allow a "reasonably" small exponent or "reasonably large"
+	 * difference in exponents to be considers "close to zero" */
+	if( mant == 0 || expB < -30 || expA - expB > 31) return 0;
+	else if( mant > 0) return 1;
+	else return -1;
+#endif
+}
+
+/**
+ * Compares a scalar with EPSILON accuracy.
+ * @param A scalar
+ * @return 1 if A > 0, -1 if A < 0, 0 otherwise
+ */
+int BOP_comp0(const MT_Scalar A)
+{
+	if (A >= BOP_EPSILON) return 1;
+	else if (0 >= A + BOP_EPSILON) return -1;
+	else return 0; 
 }
 
 /**
@@ -56,6 +85,7 @@ int BOP_comp(const MT_Scalar A, const MT_Scalar B)
  */
 int BOP_comp(const MT_Tuple3& A, const MT_Tuple3& B)
 {
+#ifndef VAR_EPSILON
 	if (A.x() >= (B.x() + BOP_EPSILON)) return 1;
 	else if (B.x() >= (A.x() + BOP_EPSILON)) return -1;
 	else if (A.y() >= (B.y() + BOP_EPSILON)) return 1;
@@ -63,6 +93,13 @@ int BOP_comp(const MT_Tuple3& A, const MT_Tuple3& B)
 	else if (A.z() >= (B.z() + BOP_EPSILON)) return 1;
 	else if (B.z() >= (A.z() + BOP_EPSILON)) return -1;
 	else return 0;
+#else
+	int result = BOP_comp(A.x(), B.x());
+	if (result != 0) return result;
+	result = BOP_comp(A.y(), B.y());
+	if (result != 0) return result;
+	return BOP_comp(A.z(), B.z());
+#endif
 }
 
 /**
@@ -116,17 +153,20 @@ bool BOP_between(const MT_Point3& p1, const MT_Point3& p2, const MT_Point3& p3)
  */
 bool BOP_collinear(const MT_Point3& p1, const MT_Point3& p2, const MT_Point3& p3)
 {
+	if( BOP_comp(p1,p2) == 0 || BOP_comp(p2,p3) == 0 ) return true;
+
 	MT_Vector3 v1 = p2 - p1;
 	MT_Vector3 v2 = p3 - p2;
 
 	/* normalize vectors before taking their cross product, so its length 
      * has some actual meaning */
+	// if(MT_fuzzyZero(v1.length()) || MT_fuzzyZero(v2.length())) return true;
 	v1.normalize();	
 	v2.normalize();
 
 	MT_Vector3 w = v1.cross(v2);
 	
-	return (BOP_comp(w.x(),0.0) == 0) && (BOP_comp(w.y(),0.0) == 0) && (BOP_comp(w.z(),0.0) == 0);
+	return (BOP_fuzzyZero(w.x()) && BOP_fuzzyZero(w.y()) && BOP_fuzzyZero(w.z()));
 }
 
 /**
@@ -192,17 +232,17 @@ bool BOP_intersect(const MT_Vector3& vL1, const MT_Point3& pL1, const MT_Vector3
 	MT_Scalar t = -1;
 	MT_Scalar den = (vL1.y()*vL2.x() - vL1.x() * vL2.y());
 	
-	if (BOP_comp(den,0)) {
+	if (!BOP_fuzzyZero(den)) {
 		t =  (pL2.y()*vL1.x() - vL1.y()*pL2.x() + pL1.x()*vL1.y() - pL1.y()*vL1.x()) / den ;
 	}
 	else {
 		den = (vL1.y()*vL2.z() - vL1.z() * vL2.y());
-		if (BOP_comp(den,0)) {
+		if (!BOP_fuzzyZero(den)) {
 			t =  (pL2.y()*vL1.z() - vL1.y()*pL2.z() + pL1.z()*vL1.y() - pL1.y()*vL1.z()) / den ;
 		}
 		else {
 			den = (vL1.x()*vL2.z() - vL1.z() * vL2.x());
-			if (BOP_comp(den,0)) {
+			if (!BOP_fuzzyZero(den)) {
 				t =  (pL2.x()*vL1.z() - vL1.x()*pL2.z() + pL1.z()*vL1.x() - pL1.x()*vL1.z()) / den ;
 			}
 			else {
@@ -326,7 +366,7 @@ MT_Scalar BOP_orientation(const MT_Plane3& p1, const MT_Plane3& p2)
 int BOP_classify(const MT_Point3& p, const MT_Plane3& plane)
 {
 	// Compare plane - point distance with zero
-	return BOP_comp(plane.signedDistance(p),0);
+	return BOP_comp0(plane.signedDistance(p));
 }
 
 /**
@@ -368,7 +408,7 @@ MT_Point3 BOP_intersectPlane(const MT_Plane3& plane, const MT_Point3& p1, const
  */
 bool BOP_containsPoint(const MT_Plane3& plane, const MT_Point3& point)
 {
-	return BOP_comp(plane.signedDistance(point),0) == 0;
+	return BOP_fuzzyZero(plane.signedDistance(point));
 }
 
 /**
@@ -424,8 +464,8 @@ MT_Scalar BOP_EpsilonDistance(const MT_Point3& p0, const MT_Point3& p1, const MT
 	MT_Scalar d1 = p0.distance(q);
 	MT_Scalar d;
 	
-	if (BOP_comp(d0,0)==0) d = 1.0;
-	else if (BOP_comp(d1,0)==0) d = 0.0;
+	if (BOP_fuzzyZero(d0)) d = 1.0;
+	else if (BOP_fuzzyZero(d1)) d = 0.0;
 	else d = d1 / d0;
 	return d;
 }
diff --git a/intern/boolop/intern/BOP_MathUtils.h b/intern/boolop/intern/BOP_MathUtils.h
index 61458bd8a78..3cff2d6a23b 100644
--- a/intern/boolop/intern/BOP_MathUtils.h
+++ b/intern/boolop/intern/BOP_MathUtils.h
@@ -1,4 +1,7 @@
 /**
+ *
+ * $Id$
+ *
  * ***** BEGIN GPL/BL DUAL LICENSE BLOCK *****
  *
  * This program is free software; you can redistribute it and/or
@@ -23,7 +26,7 @@
  *
  * The Original Code is: all of this file.
  *
- * Contributor(s): none yet.
+ * Contributor(s): Marc Freixas, Ken Hughes
  *
  * ***** END GPL/BL DUAL LICENSE BLOCK *****
  */
@@ -36,16 +39,23 @@
 #include "MT_Point3.h"
 #include "MT_Plane3.h"
 
+/* define this to give better precision comparisons */
+#define VAR_EPSILON
+
+#ifndef VAR_EPSILON
 const MT_Scalar BOP_EPSILON(1.0e-5);
+#else
+const MT_Scalar BOP_EPSILON(9.3132257461547852e-10);	/* ~= 2**-30 */
+#endif
 
 inline int BOP_sign(MT_Scalar x) {
     return x < 0.0 ? -1 : x > 0.0 ? 1 : 0;
 }
 inline MT_Scalar BOP_abs(MT_Scalar x) { return fabs(x); }
-inline bool      BOP_fuzzyZero(MT_Scalar x) { return BOP_abs(x) < BOP_EPSILON; }
-
 int BOP_comp(const MT_Scalar A, const MT_Scalar B);
 int BOP_comp(const MT_Tuple3& A, const MT_Tuple3& B);
+int BOP_comp0(const MT_Scalar A);
+inline bool BOP_fuzzyZero(MT_Scalar x) { return BOP_comp0(x) == 0; }
 int BOP_exactComp(const MT_Scalar A, const MT_Scalar B);
 int BOP_exactComp(const MT_Tuple3& A, const MT_Tuple3& B);
 bool BOP_between(const MT_Point3& p1, const MT_Point3& p2, const MT_Point3& p3);
diff --git a/intern/boolop/intern/BOP_Merge.cpp b/intern/boolop/intern/BOP_Merge.cpp
index 5839d38181a..c2b1a2db2b7 100644
--- a/intern/boolop/intern/BOP_Merge.cpp
+++ b/intern/boolop/intern/BOP_Merge.cpp
@@ -60,32 +60,6 @@ void BOP_Merge::mergeFaces(BOP_Mesh *m, BOP_Index v)
 	// Merge faces
 	mergeFaces();
 
-	/*
-	 * HACK: somehow triangular faces are being created with two vertices the
-	 * same.  If it's happening in BOP_Mesh::replaceVertexIndex() we should
-	 * be catching it, so either it's not happening there or we aren't
-	 * catching it (duh).  Until we figure this out, this hack cleans things.
-	 * 
-	 * Test for any invalid faces: if any two vertices are the same of a
-	 * triangle, the face is broken.  Further, I don't believe it's possible
-	 * to have any quads at this point, so if we find one send a message
-	 * to stdout.
-	 */
-
-	BOP_Faces faces = m_mesh->getFaces();
-	const BOP_IT_Faces ifacesIEnd = (faces.end());
-	for(BOP_IT_Faces faceI=faces.begin();faceI!=ifacesIEnd;faceI++) {
-		if ((*faceI)->getTAG() != BROKEN ) {
-			BOP_Index i1 = (*faceI)->getVertex(0);
-			BOP_Index i2 = (*faceI)->getVertex(1);
-			BOP_Index i3 = (*faceI)->getVertex(2);
-			if ( (*faceI)->size() == 4)
-				cout << "BOP_Merge::mergeFaces found a quad: this is an error" << endl;
-			if (i1 == i2 || i2 == i3 || i3 == i1 )
-				(*faceI)->setTAG(BROKEN);
-		}
-	}
-
 	do {
 		// Add quads ...
 		cont = createQuads();
diff --git a/intern/boolop/intern/BOP_Mesh.cpp b/intern/boolop/intern/BOP_Mesh.cpp
index 3b194ef72d4..0b70cc61533 100644
--- a/intern/boolop/intern/BOP_Mesh.cpp
+++ b/intern/boolop/intern/BOP_Mesh.cpp
@@ -681,6 +681,23 @@ unsigned int BOP_Mesh::getNumFaces(BOP_TAG tag)
 }
 
 /**
+ * Marks faces which bad edges as BROKEN (invalid face, no further processing).
+ * @param edge edge which is being replaced
+ * @param mesh mesh containing faces
+ */
+
+static void removeBrokenFaces( BOP_Edge *edge, BOP_Mesh *mesh )
+{
+	BOP_Faces m_faces = mesh->getFaces();
+
+	BOP_Indexs edgeFaces = edge->getFaces();
+	const BOP_IT_Indexs edgeFacesEnd = edgeFaces.end();
+	for(BOP_IT_Indexs idxFace=edgeFaces.begin();idxFace!=edgeFacesEnd;
+			   idxFace++)
+		m_faces[*idxFace]->setTAG(BROKEN);
+}
+
+/**
  * Replaces a vertex index.
  * @param oldIndex old vertex index
  * @param newIndex new vertex index
@@ -695,9 +712,6 @@ BOP_Index BOP_Mesh::replaceVertexIndex(BOP_Index oldIndex, BOP_Index newIndex)
 	BOP_Vertex *newVertex = m_vertexs[newIndex];
 	BOP_Indexs oldEdges = oldVertex->getEdges();
 
-	BOP_Index edgeIdx=0;
-	bool found = false;
-
 	// Update faces to the newIndex
 	BOP_IT_Indexs oldEdgesEnd = oldEdges.end();
 	for(oldEdgeIndex=oldEdges.begin();oldEdgeIndex!=oldEdgesEnd;
@@ -706,14 +720,9 @@ BOP_Index BOP_Mesh::replaceVertexIndex(BOP_Index oldIndex, BOP_Index newIndex)
 		if ((edge->getVertex1()==oldIndex && edge->getVertex2()==newIndex) ||
 			(edge->getVertex2()==oldIndex && edge->getVertex1()==newIndex)) {
 			// Remove old edge  ==> set edge faces to BROKEN      
-			BOP_Indexs edgeFaces = edge->getFaces();
-			const BOP_IT_Indexs edgeFacesEnd = edgeFaces.end();
-			for(BOP_IT_Indexs idxFace=edgeFaces.begin();idxFace!=edgeFacesEnd;
-				   idxFace++) {
-				m_faces[*idxFace]->setTAG(BROKEN);
-			}
-			edgeIdx = *oldEdgeIndex;
-			found = true;
+			removeBrokenFaces( edge, this );
+			oldVertex->removeEdge(*oldEdgeIndex);
+			newVertex->removeEdge(*oldEdgeIndex);
 		}
 		else {
 			BOP_Indexs faces = edge->getFaces();
@@ -724,10 +733,6 @@ BOP_Index BOP_Mesh::replaceVertexIndex(BOP_Index oldIndex, BOP_Index newIndex)
 			}
 		}
 	} 
-	if (found) {
-		oldVertex->removeEdge(edgeIdx);
-		newVertex->removeEdge(edgeIdx);
-	}
 
 	oldEdgesEnd = oldEdges.end();
 	for(oldEdgeIndex=oldEdges.begin();oldEdgeIndex!=oldEdgesEnd;
@@ -739,6 +744,10 @@ BOP_Index BOP_Mesh::replaceVertexIndex(BOP_Index oldIndex, BOP_Index newIndex)
 		v1 = (v1==oldIndex?edge->getVertex2():v1);      
 		if ((edge2 = getEdge(newIndex,v1)) == NULL) {
 			edge->replaceVertexIndex(oldIndex,newIndex);
+			if ( edge->getVertex1() == edge->getVertex2() ) {
+				removeBrokenFaces( edge, this );
+				oldVertex->removeEdge(*oldEdgeIndex);
+			}
 #ifdef HASH
 			rehashVertex(oldIndex,newIndex,v1);
 #endif
@@ -754,6 +763,11 @@ BOP_Index BOP_Mesh::replaceVertexIndex(BOP_Index oldIndex, BOP_Index newIndex)
 			BOP_Vertex *oppositeVertex = m_vertexs[v1];
 			oppositeVertex->removeEdge(*oldEdgeIndex);
 			edge->replaceVertexIndex(oldIndex,newIndex);
+			if ( edge->getVertex1() == edge->getVertex2() ) {
+				removeBrokenFaces( edge, this );
+				oldVertex->removeEdge(*oldEdgeIndex);
+				newVertex->removeEdge(*oldEdgeIndex);
+			}
 #ifdef HASH
 			rehashVertex(oldIndex,newIndex,v1);
 #endif
@@ -1063,3 +1077,4 @@ void BOP_Mesh::updatePlanes()
 	  face->setPlane(plane);
 	}
 }
+
diff --git a/intern/boolop/intern/BOP_Triangulator.cpp b/intern/boolop/intern/BOP_Triangulator.cpp
index ad23e69d8c0..fd7b3154195 100644
--- a/intern/boolop/intern/BOP_Triangulator.cpp
+++ b/intern/boolop/intern/BOP_Triangulator.cpp
@@ -90,8 +90,11 @@ void BOP_triangulateA(BOP_Mesh *mesh, BOP_Faces *faces, BOP_Face * face, BOP_Ind
 	
 	BOP_addFace(mesh, faces, face1, face->getTAG());
 	BOP_addFace(mesh, faces, face2, face->getTAG());
+	face1->setSplit(face->getSplit());
+	face2->setSplit(face->getSplit());
 	
 	face->setTAG(BROKEN);
+	face->freeBBox();
 }
 
 /**
@@ -122,7 +125,11 @@ void BOP_triangulateB(BOP_Mesh* mesh, BOP_Faces* faces, BOP_Face* face, BOP_Inde
 	BOP_addFace(mesh,faces,face1,face->getTAG());
 	BOP_addFace(mesh,faces,face2,face->getTAG());
 	BOP_addFace(mesh,faces,face3,face->getTAG());
+	face1->setSplit(face->getSplit());
+	face2->setSplit(face->getSplit());
+	face3->setSplit(face->getSplit());
 	face->setTAG(BROKEN);
+	face->freeBBox();
 }
 
 
@@ -180,26 +187,33 @@ void BOP_triangulateC_split(BOP_Mesh* mesh, BOP_Faces* faces, BOP_Face* face,
 {
 	BOP_Index v = BOP_getTriangleVertex(mesh, v1, v2, v4, v5);
 	BOP_Index w = (v == v4 ? v5 : v4);
+	BOP_Face *face1 = new BOP_Face3(v1, v, w, face->getPlane(), face->getOriginalFace());
+	BOP_Face *face2 = new BOP_Face3(v1, v2, v, face->getPlane(), face->getOriginalFace());
+	BOP_Face *face3 = new BOP_Face3(v1, w, v3, face->getPlane(), face->getOriginalFace());
 	
 	// v1 v w defines the nice triangle in the correct order
 	// v1 v2 v defines one lateral triangle in the correct order
 	// v1 w v3 defines the other lateral triangle in the correct order
 	// w v v2 v3 defines the quad in the correct order
 	
-	BOP_addFace(mesh, faces, new BOP_Face3(v1, v, w, face->getPlane(), 
-				face->getOriginalFace()), face->getTAG());
-	BOP_addFace(mesh, faces, new BOP_Face3(v1, v2, v, face->getPlane(), 
-				face->getOriginalFace()), face->getTAG());
-	BOP_addFace(mesh, faces, new BOP_Face3(v1, w, v3, face->getPlane(), 
-				face->getOriginalFace()), face->getTAG());
+	BOP_addFace(mesh, faces, face1, face->getTAG());
+	BOP_addFace(mesh, faces, face2, face->getTAG());
+	BOP_addFace(mesh, faces, face3, face->getTAG());
+
+	face1->setSplit(face->getSplit());
+	face2->setSplit(face->getSplit());
+	face3->setSplit(face->getSplit());
 	
 	BOP_Face *faces45[2];
 	
 	BOP_splitQuad(mesh, face->getPlane(), v2, v3, w, v, faces45, face->getOriginalFace());
 	BOP_addFace(mesh, faces, faces45[0], face->getTAG());
 	BOP_addFace(mesh, faces, faces45[1], face->getTAG());
+	faces45[0]->setSplit(face->getSplit());
+	faces45[1]->setSplit(face->getSplit());
 	
 	face->setTAG(BROKEN); 
+	face->freeBBox();
 }
 
 
@@ -254,15 +268,19 @@ void BOP_triangulateD_split(BOP_Mesh* mesh, BOP_Faces* faces, BOP_Face* face,
 {
 	BOP_Index v = BOP_getNearestVertex(mesh, v1, v4, v5);
 	BOP_Index w = (v == v4 ? v5 : v4);
+	BOP_Face *face1 = new BOP_Face3(v1, v, v3, face->getPlane(), face->getOriginalFace());
+	BOP_Face *face2 = new BOP_Face3(v, w, v3, face->getPlane(), face->getOriginalFace());
+	BOP_Face *face3 = new BOP_Face3(w, v2, v3, face->getPlane(), face->getOriginalFace());
 	
-	BOP_addFace(mesh, faces, new BOP_Face3(v1, v, v3, face->getPlane(), 
-				face->getOriginalFace()), face->getTAG());
-	BOP_addFace(mesh, faces, new BOP_Face3(v, w, v3, face->getPlane(), 
-				face->getOriginalFace()), face->getTAG());
-	BOP_addFace(mesh, faces, new BOP_Face3(w, v2, v3, face->getPlane(), 
-				face->getOriginalFace()), face->getTAG());
-	
+	BOP_addFace(mesh, faces, face1, face->getTAG());
+	BOP_addFace(mesh, faces, face2, face->getTAG());
+	BOP_addFace(mesh, faces, face3, face->getTAG());
+	face1->setSplit(face->getSplit());
+	face2->setSplit(face->getSplit());
+	face3->setSplit(face->getSplit());
+
 	face->setTAG(BROKEN);
+	face->freeBBox();
 }
 
 
@@ -328,7 +346,11 @@ void BOP_triangulateE(BOP_Mesh* mesh, BOP_Faces* faces, BOP_Face* face,
 	BOP_addFace(mesh, faces, face1, face->getTAG());
 	BOP_addFace(mesh, faces, faces23[0], face->getTAG());
 	BOP_addFace(mesh, faces, faces23[1], face->getTAG());
+	face1->setSplit(face->getSplit());
+	faces23[0]->setSplit(face->getSplit());
+	faces23[1]->setSplit(face->getSplit());
 	face->setTAG(BROKEN);
+	face->freeBBox();
 }
 
 /**
@@ -380,8 +402,13 @@ void BOP_triangulateF(BOP_Mesh* mesh, BOP_Faces* faces, BOP_Face* face,
 	BOP_addFace(mesh, faces, faces12[1], face->getTAG());
 	BOP_addFace(mesh, faces, faces34[0], face->getTAG());
 	BOP_addFace(mesh, faces, faces34[1], face->getTAG());
+	faces12[0]->setSplit(face->getSplit());
+	faces12[1]->setSplit(face->getSplit());
+	faces34[0]->setSplit(face->getSplit());
+	faces34[1]->setSplit(face->getSplit());
 	
 	face->setTAG(BROKEN);
+	face->freeBBox();
 }
 
 /**