some progress

src/boolean_result.cpp

@@ -856,6 +856,7 @@ Manifold::Impl Boolean3::Result(OpType op) const {
 
   UpdateReference(outR, inP_, inQ_, invertQ);
 
+  outR.CreateFaces();
   outR.SimplifyTopology();
   outR.RemoveUnreferencedVerts();
 

src/edge_op.cpp

@@ -33,27 +33,6 @@ bool Is01Longest(vec2 v0, vec2 v1, vec2 v2) {
   return l[0] > l[1] && l[0] > l[2];
 }
 
-struct FlagEdge {
-  VecView<const Halfedge> halfedge;
-  VecView<const TriRef> triRef;
-
-  bool operator()(int edge) const {
-    if (halfedge[edge].pairedHalfedge < 0) return false;
-    // Flag redundant edges - those where the startVert is surrounded by only
-    // two original triangles.
-    const TriRef ref0 = triRef[edge / 3];
-    int current = NextHalfedge(halfedge[edge].pairedHalfedge);
-    const TriRef ref1 = triRef[current / 3];
-    while (current != edge) {
-      current = NextHalfedge(halfedge[current].pairedHalfedge);
-      int tri = current / 3;
-      const TriRef ref = triRef[tri];
-      if (!ref.SameFace(ref0) && !ref.SameFace(ref1)) return false;
-    }
-    return true;
-  }
-};
-
 struct SwappableEdge {
   VecView<const Halfedge> halfedge;
   VecView<const vec3> vertPos;
@@ -180,10 +159,29 @@ void Manifold::Impl::SimplifyTopology() {
 
   while (1) {
     ZoneScopedN("CollapseFlaggedEdge");
+    for_each_n(policy, countAt(0_uz), nbEdges, [&bFlags, this](size_t edge) {
+      if (halfedge_[edge].pairedHalfedge < 0) {
+        bFlags[edge] = false;
+        return;
+      }
+      // Flag redundant edges - those where the startVert is surrounded by only
+      // two original triangles.
+      const TriRef ref0 = meshRelation_.triRef[edge / 3];
+      int current = NextHalfedge(halfedge_[edge].pairedHalfedge);
+      const TriRef ref1 = meshRelation_.triRef[current / 3];
+      while (current != edge) {
+        current = NextHalfedge(halfedge_[current].pairedHalfedge);
+        int tri = current / 3;
+        const TriRef ref = meshRelation_.triRef[tri];
+        if (!ref.SameFace(ref0) && !ref.SameFace(ref1)) {
+          bFlags[edge] = false;
+          return;
+        }
+      }
+      bFlags[edge] = true;
+    });
+
     numFlagged = 0;
-    FlagEdge se{halfedge_, meshRelation_.triRef};
-    for_each_n(policy, countAt(0_uz), nbEdges,
-               [&](size_t i) { bFlags[i] = se(i); });
     for (size_t i = 0; i < nbEdges; ++i) {
       if (bFlags[i]) {
         if (CollapseEdge(i, scratchBuffer)) numFlagged++;
@@ -436,27 +434,36 @@ bool Manifold::Impl::CollapseEdge(const int edge, std::vector<int>& edges) {
   // Orbit startVert
   int start = halfedge_[tri1edge[1]].pairedHalfedge;
   current = start;
-  TriRef refCheck = triRef[toRemove.pairedHalfedge / 3];
+  int tri = toRemove.pairedHalfedge / 3;
+  TriRef refCheck = triRef[tri];
+  mat2x3 projection = GetAxisAlignedProjection(faceNormal_[tri]);
   vec3 pLast = vertPos_[halfedge_[tri1edge[1]].endVert];
   while (current != tri0edge[2]) {
     current = NextHalfedge(current);
     vec3 pNext = vertPos_[halfedge_[current].endVert];
-    const int tri = current / 3;
-    const TriRef ref = triRef[tri];
-    const mat2x3 projection = GetAxisAlignedProjection(faceNormal_[tri]);
+    const TriRef ref = triRef[current / 3];
     // Don't collapse if the edge is not redundant (this may have changed due
     // to the collapse of neighbors).
     if (!ref.SameFace(refCheck)) {
-      refCheck = triRef[edge / 3];
+      tri = edge / 3;
+      refCheck = triRef[tri];
+      projection = GetAxisAlignedProjection(faceNormal_[tri]);
       if (!ref.SameFace(refCheck)) {
         return false;
+      } else {
+        // Don't collapse if the edges separating the faces are not colinear
+        // (can happen when the two faces are coplanar).
+        if (CCW(projection * pOld, projection * pLast, projection * pNew,
+                tolerance_) != 0)
+          return false;
       }
     }
 
     // Don't collapse edge if it would cause a triangle to invert.
     if (CCW(projection * pNext, projection * pLast, projection * pNew,
-            epsilon_) < 0)
+            epsilon_) < 0) {
       return false;
+    }
 
     pLast = pNext;
     current = halfedge_[current].pairedHalfedge;

src/impl.cpp

@@ -227,6 +227,7 @@ void Manifold::Impl::CreateFaces() {
     if (meshRelation_.triRef[tp.tri].faceID >= 0) continue;
 
     meshRelation_.triRef[tp.tri].faceID = tp.tri;
+    const int meshID = meshRelation_.triRef[tp.tri].meshID;
     const vec3 base = vertPos_[halfedge_[3 * tp.tri].startVert];
     const vec3 normal = faceNormal_[tp.tri];
     std::vector<int> interiorHalfedges = {3 * tp.tri, 3 * tp.tri + 1,
@@ -235,11 +236,17 @@ void Manifold::Impl::CreateFaces() {
       const int h =
           NextHalfedge(halfedge_[interiorHalfedges.back()].pairedHalfedge);
       interiorHalfedges.pop_back();
-      if (meshRelation_.triRef[h / 3].faceID >= 0) continue;
+      const int tri = h / 3;
+      TriRef& ref = meshRelation_.triRef[tri];
+      if (ref.faceID >= 0) continue;
 
       const vec3 v = vertPos_[halfedge_[h].endVert];
-      if (abs(dot(v - base, normal)) < tolerance_) {
-        meshRelation_.triRef[h / 3].faceID = tp.tri;
+
+      if (abs(dot(v - base, normal)) < tolerance_ &&
+          (ref.meshID == meshID || TriCCW(tri, tolerance_) == 0)) {
+        faceNormal_[tri] = faceNormal_[tp.tri];
+        ref.meshID = meshID;
+        ref.faceID = tp.tri;
 
         if (interiorHalfedges.empty() ||
             h != halfedge_[interiorHalfedges.back()].pairedHalfedge) {

src/impl.h

@@ -280,6 +280,7 @@ struct Manifold::Impl {
   double GetProperty(Property prop) const;
   void CalculateCurvature(int gaussianIdx, int meanIdx);
   void CalculateBBox();
+  int TriCCW(size_t tri, double epsilon) const;
   bool IsFinite() const;
   bool IsIndexInBounds(VecView<const ivec3> triVerts) const;
   void SetEpsilon(double minEpsilon = -1, bool useSingle = false);

src/properties.cpp

@@ -123,54 +123,19 @@ struct CheckHalfedges {
     return good;
   }
 };
-
-struct CheckCCW {
-  VecView<const Halfedge> halfedges;
-  VecView<const vec3> vertPos;
-  VecView<const vec3> triNormal;
-  const double tol;
-
-  bool operator()(size_t face) const {
-    if (halfedges[3 * face].pairedHalfedge < 0) return true;
-
-    const mat2x3 projection = GetAxisAlignedProjection(triNormal[face]);
-    vec2 v[3];
-    for (int i : {0, 1, 2})
-      v[i] = projection * vertPos[halfedges[3 * face + i].startVert];
-
-    int ccw = CCW(v[0], v[1], v[2], std::abs(tol));
-    bool check = tol > 0 ? ccw >= 0 : ccw == 0;
-
-#ifdef MANIFOLD_DEBUG
-    if (tol > 0 && !check) {
-      vec2 v1 = v[1] - v[0];
-      vec2 v2 = v[2] - v[0];
-      double area = v1.x * v2.y - v1.y * v2.x;
-      double base2 = std::max(la::dot(v1, v1), la::dot(v2, v2));
-      double base = std::sqrt(base2);
-      vec3 V0 = vertPos[halfedges[3 * face].startVert];
-      vec3 V1 = vertPos[halfedges[3 * face + 1].startVert];
-      vec3 V2 = vertPos[halfedges[3 * face + 2].startVert];
-      vec3 norm = la::cross(V1 - V0, V2 - V0);
-      printf(
-          "Tri %ld does not match normal, approx height = %g, base = %g\n"
-          "tol = %g, area2 = %g, base2*tol2 = %g\n"
-          "normal = %g, %g, %g\n"
-          "norm = %g, %g, %g\nverts: %d, %d, %d\n",
-          static_cast<long>(face), area / base, base, tol, area * area,
-          base2 * tol * tol, triNormal[face].x, triNormal[face].y,
-          triNormal[face].z, norm.x, norm.y, norm.z,
-          halfedges[3 * face].startVert, halfedges[3 * face + 1].startVert,
-          halfedges[3 * face + 2].startVert);
-    }
-#endif
-    return check;
-  }
-};
 }  // namespace
 
 namespace manifold {
 
+int Manifold::Impl::TriCCW(size_t tri, double tol) const {
+  if (halfedge_[3 * tri].pairedHalfedge < 0) return true;
+  const mat2x3 projection = GetAxisAlignedProjection(faceNormal_[tri]);
+  vec2 v[3];
+  for (int i : {0, 1, 2})
+    v[i] = projection * vertPos_[halfedge_[3 * tri + i].startVert];
+  return CCW(v[0], v[1], v[2], std::abs(tol));
+}
+
 /**
  * Returns true if this manifold is in fact an oriented even manifold and all of
  * the data structures are consistent.
@@ -267,18 +232,25 @@ bool Manifold::Impl::IsSelfIntersecting() const {
  */
 bool Manifold::Impl::MatchesTriNormals() const {
   if (halfedge_.size() == 0 || faceNormal_.size() != NumTri()) return true;
-  return all_of(countAt(0_uz), countAt(NumTri()),
-                CheckCCW({halfedge_, vertPos_, faceNormal_, 2 * epsilon_}));
+  return all_of(countAt(0_uz), countAt(NumTri()), [this](size_t tri) {
+    if (halfedge_[3 * tri].startVert < 0) return true;
+    const vec3 normal = faceNormal_[tri];
+    const vec3 p(dot(normal, vertPos_[halfedge_[3 * tri].startVert]),
+                 dot(normal, vertPos_[halfedge_[3 * tri + 1].startVert]),
+                 dot(normal, vertPos_[halfedge_[3 * tri + 2].startVert]));
+    return maxelem(p) - minelem(p) <= 2 * tolerance_ &&
+           TriCCW(tri, 2 * epsilon_) >= 0;
+  });
 }
 
 /**
  * Returns the number of triangles that are colinear within epsilon_.
  */
 int Manifold::Impl::NumDegenerateTris() const {
-  if (halfedge_.size() == 0 || faceNormal_.size() != NumTri()) return true;
-  return count_if(
-      countAt(0_uz), countAt(NumTri()),
-      CheckCCW({halfedge_, vertPos_, faceNormal_, -1 * epsilon_ / 2}));
+  if (halfedge_.size() == 0 || faceNormal_.size() != NumTri()) return 0;
+  return count_if(countAt(0_uz), countAt(NumTri()), [this](size_t tri) {
+    return TriCCW(tri, epsilon_ / 2) == 0;
+  });
 }
 
 double Manifold::Impl::GetProperty(Property prop) const {

test/boolean_test.cpp

@@ -42,6 +42,7 @@ TEST(Boolean, MeshGLRoundTrip) {
   Manifold result = cube + cube.Translate({1, 1, 0});
 
   ASSERT_LT(result.OriginalID(), 0);
+  EXPECT_LE(result.NumDegenerateTris(), 0);
   ExpectMeshes(result, {{18, 32}});
   RelatedGL(result, {original});
 
@@ -50,6 +51,7 @@ TEST(Boolean, MeshGLRoundTrip) {
   const Manifold result2(inGL);
 
   ASSERT_LT(result2.OriginalID(), 0);
+  EXPECT_LE(result2.NumDegenerateTris(), 0);
   ExpectMeshes(result2, {{16, 28}});
   RelatedGL(result2, {original});