Separate coplanarity from faceID (#991)

* separated faceID

* clean up tests

* fix AsOriginal

* reenable CreateFaces and keep mesh relation intact

* fix compilation

* updates docs

bindings/python/manifold3d.cpp

@@ -351,9 +351,7 @@ NB_MODULE(manifold3d, m) {
           "Get the surface area of the manifold\n This is clamped to zero for "
           "a given face if they are within the Precision().")
       .def("original_id", &Manifold::OriginalID, manifold__original_id)
-      .def("as_original", &Manifold::AsOriginal,
-           nb::arg("property_tolerance") = nb::list(),
-           manifold__as_original__property_tolerance)
+      .def("as_original", &Manifold::AsOriginal, manifold__as_original)
       .def("is_empty", &Manifold::IsEmpty, manifold__is_empty)
       .def("decompose", &Manifold::Decompose, manifold__decompose)
       .def("split", &Manifold::Split, nb::arg("cutter"),

bindings/wasm/bindings.cpp

@@ -168,7 +168,7 @@ EMSCRIPTEN_BINDINGS(whatever) {
       .function("calculateCurvature", &Manifold::CalculateCurvature)
       .function("_CalculateNormals", &Manifold::CalculateNormals)
       .function("originalID", &Manifold::OriginalID)
-      .function("_AsOriginal", &Manifold::AsOriginal);
+      .function("asOriginal", &Manifold::AsOriginal);
 
   // Manifold Static Methods
   function("_Cube", &Manifold::Cube);

bindings/wasm/bindings.js

@@ -228,14 +228,6 @@ Module.setup = function() {
     return this._CalculateNormals(normalIdx, minSharpAngle);
   };
 
-  Module.Manifold.prototype.asOriginal = function(propertyTolerance = []) {
-    const tol = new Module.Vector_f64();
-    toVec(tol, propertyTolerance);
-    const result = this._AsOriginal(tol);
-    tol.delete();
-    return result
-  };
-
   Module.Manifold.prototype.setProperties = function(numProp, func) {
     const oldNumProp = this.numProp();
     const wasmFuncPtr = addFunction(function(newPtr, vec3Ptr, oldPtr) {

include/manifold/manifold.h

@@ -78,13 +78,10 @@ struct MeshGLP {
   /// This matrix is stored in column-major order and the length of the overall
   /// vector is 12 * runOriginalID.size().
   std::vector<Precision> runTransform;
-  /// Optional: Length NumTri, contains an ID of the source face this triangle
-  /// comes from. When auto-generated, this ID will be a triangle index into the
-  /// original mesh. All neighboring coplanar triangles from that input mesh
-  /// will refer to a single triangle of that group as the faceID. When
-  /// supplying faceIDs, ensure that triangles with the same ID are in fact
-  /// coplanar and have consistent properties (within some tolerance) or the
-  /// output will be surprising.
+  /// Optional: Length NumTri, contains the source face ID this
+  /// triangle comes from. When auto-generated, this ID will be a triangle index
+  /// into the original mesh. This index/ID is purely for external use (e.g.
+  /// recreating polygonal faces) and will not affect Manifold's algorithms.
   std::vector<I> faceID;
   /// Optional: The X-Y-Z-W weighted tangent vectors for smooth Refine(). If
   /// non-empty, must be exactly four times as long as Mesh.triVerts. Indexed
@@ -245,7 +242,7 @@ class Manifold {
    */
   ///@{
   int OriginalID() const;
-  Manifold AsOriginal(const std::vector<double>& propertyTolerance = {}) const;
+  Manifold AsOriginal() const;
   static uint32_t ReserveIDs(uint32_t);
   ///@}
 

src/impl.cpp

@@ -93,14 +93,11 @@ struct UpdateMeshID {
 
 struct CoplanarEdge {
   VecView<std::pair<int, int>> face2face;
-  VecView<std::pair<int, int>> vert2vert;
   VecView<double> triArea;
   VecView<const Halfedge> halfedge;
   VecView<const vec3> vertPos;
   VecView<const TriRef> triRef;
   VecView<const ivec3> triProp;
-  VecView<const double> prop;
-  VecView<const double> propTol;
   const int numProp;
   const double precision;
 
@@ -118,19 +115,9 @@ struct CoplanarEdge {
     const int jointNum = edge.pairedHalfedge - 3 * pairFace;
 
     if (numProp > 0) {
-      const int prop0 = triProp[edgeFace][baseNum];
-      const int prop1 = triProp[pairFace][jointNum == 2 ? 0 : jointNum + 1];
-      bool propEqual = true;
-      for (int p = 0; p < numProp; ++p) {
-        if (std::abs(prop[numProp * prop0 + p] - prop[numProp * prop1 + p]) >
-            propTol[p]) {
-          propEqual = false;
-          break;
-        }
-      }
-      if (propEqual) {
-        vert2vert[edgeIdx] = std::make_pair(prop0, prop1);
-      }
+      if (triProp[edgeFace][baseNum] != triProp[pairFace][Next3(jointNum)] ||
+          triProp[edgeFace][Next3(baseNum)] != triProp[pairFace][jointNum])
+        return;
     }
 
     if (!edge.IsForward()) return;
@@ -158,32 +145,6 @@ struct CoplanarEdge {
     // Only operate on coplanar triangles
     if (volume > std::max(area, areaPair) * precision) return;
 
-    // Check property linearity
-    if (area > 0) {
-      normal /= area;
-      for (int i = 0; i < numProp; ++i) {
-        const double scale = precision / propTol[i];
-
-        const double baseProp = prop[numProp * triProp[edgeFace][baseNum] + i];
-        const double jointProp =
-            prop[numProp * triProp[pairFace][jointNum] + i];
-        const double edgeProp = prop[numProp * triProp[edgeFace][edgeNum] + i];
-        const double pairProp = prop[numProp * triProp[pairFace][pairNum] + i];
-
-        const vec3 iJointVec =
-            jointVec + normal * scale * (jointProp - baseProp);
-        const vec3 iEdgeVec = edgeVec + normal * scale * (edgeProp - baseProp);
-        const vec3 iPairVec = pairVec + normal * scale * (pairProp - baseProp);
-
-        vec3 cross = la::cross(iJointVec, iEdgeVec);
-        const double areaP = std::max(
-            la::length(cross), la::length(la::cross(iPairVec, iJointVec)));
-        const double volumeP = std::abs(la::dot(cross, iPairVec));
-        // Only operate on consistent triangles
-        if (volumeP > areaP * precision) return;
-      }
-    }
-
     face2face[edgeIdx] = std::make_pair(edgeFace, pairFace);
   }
 };
@@ -339,40 +300,67 @@ void Manifold::Impl::RemoveUnreferencedVerts() {
            });
 }
 
-void Manifold::Impl::InitializeOriginal() {
+void Manifold::Impl::InitializeOriginal(bool keepFaceID) {
   const int meshID = ReserveIDs(1);
   meshRelation_.originalID = meshID;
   auto& triRef = meshRelation_.triRef;
   triRef.resize(NumTri());
   for_each_n(autoPolicy(NumTri(), 1e5), countAt(0), NumTri(),
-             [meshID, &triRef](const int tri) {
-               triRef[tri] = {meshID, meshID, tri};
+             [meshID, keepFaceID, &triRef](const int tri) {
+               triRef[tri] = {meshID, meshID, tri,
+                              keepFaceID ? triRef[tri].faceID : tri};
              });
   meshRelation_.meshIDtransform.clear();
   meshRelation_.meshIDtransform[meshID] = {meshID};
 }
 
-void Manifold::Impl::CreateFaces(const std::vector<double>& propertyTolerance) {
+void Manifold::Impl::CreateFaces() {
   ZoneScoped;
-  constexpr double kDefaultPropTolerance = 1e-5;
-  Vec<double> propertyToleranceD =
-      propertyTolerance.empty()
-          ? Vec<double>(meshRelation_.numProp, kDefaultPropTolerance)
-          : propertyTolerance;
-
   Vec<std::pair<int, int>> face2face(halfedge_.size(), {-1, -1});
   Vec<std::pair<int, int>> vert2vert(halfedge_.size(), {-1, -1});
   Vec<double> triArea(NumTri());
-  for_each_n(autoPolicy(halfedge_.size(), 1e4), countAt(0), halfedge_.size(),
-             CoplanarEdge({face2face, vert2vert, triArea, halfedge_, vertPos_,
-                           meshRelation_.triRef, meshRelation_.triProperties,
-                           meshRelation_.properties, propertyToleranceD,
-                           meshRelation_.numProp, precision_}));
 
-  if (meshRelation_.triProperties.size() > 0) {
+  const size_t numProp = NumProp();
+  if (numProp > 0) {
+    for_each_n(
+        autoPolicy(halfedge_.size(), 1e4), countAt(0), halfedge_.size(),
+        [&vert2vert, numProp, this](const int edgeIdx) {
+          const Halfedge edge = halfedge_[edgeIdx];
+          const Halfedge pair = halfedge_[edge.pairedHalfedge];
+          const int edgeFace = edgeIdx / 3;
+          const int pairFace = edge.pairedHalfedge / 3;
+
+          if (meshRelation_.triRef[edgeFace].meshID !=
+              meshRelation_.triRef[pairFace].meshID)
+            return;
+
+          const int baseNum = edgeIdx - 3 * edgeFace;
+          const int jointNum = edge.pairedHalfedge - 3 * pairFace;
+
+          const int prop0 = meshRelation_.triProperties[edgeFace][baseNum];
+          const int prop1 =
+              meshRelation_
+                  .triProperties[pairFace][jointNum == 2 ? 0 : jointNum + 1];
+          bool propEqual = true;
+          for (size_t p = 0; p < numProp; ++p) {
+            if (meshRelation_.properties[numProp * prop0 + p] !=
+                meshRelation_.properties[numProp * prop1 + p]) {
+              propEqual = false;
+              break;
+            }
+          }
+          if (propEqual) {
+            vert2vert[edgeIdx] = std::make_pair(prop0, prop1);
+          }
+        });
     DedupePropVerts(meshRelation_.triProperties, vert2vert);
   }
 
+  for_each_n(autoPolicy(halfedge_.size(), 1e4), countAt(0), halfedge_.size(),
+             CoplanarEdge({face2face, triArea, halfedge_, vertPos_,
+                           meshRelation_.triRef, meshRelation_.triProperties,
+                           meshRelation_.numProp, precision_}));
+
   std::vector<int> components;
   const int numComponent = GetLabels(components, face2face, NumTri());
 
@@ -394,7 +382,7 @@ void Manifold::Impl::CreateFaces(const std::vector<double>& propertyTolerance) {
   for (size_t tri = 0; tri < NumTri(); ++tri) {
     const int referenceTri = comp2tri[components[tri]];
     if (referenceTri >= 0) {
-      triRef[tri].tri = referenceTri;
+      triRef[tri].faceID = referenceTri;
     }
   }
 }
@@ -511,8 +499,9 @@ void Manifold::Impl::WarpBatch(std::function<void(VecView<vec3>)> warpFunc) {
   faceNormal_.resize(0);  // force recalculation of triNormal
   CalculateNormals();
   SetPrecision();
-  InitializeOriginal();
   Finish();
+  CreateFaces();
+  meshRelation_.originalID = -1;
 }
 
 Manifold::Impl Manifold::Impl::Transform(const mat3x4& transform_) const {

src/impl.h

@@ -147,6 +147,7 @@ struct Manifold::Impl {
           ref.meshID = meshID;
           ref.originalID = originalID;
           ref.tri = meshGL.faceID.empty() ? tri : meshGL.faceID[tri];
+          ref.faceID = tri;
         }
 
         if (meshGL.runTransform.empty()) {
@@ -209,6 +210,8 @@ struct Manifold::Impl {
       InitializeOriginal();
     }
 
+    CreateFaces();
+
     SimplifyTopology();
     Finish();
 
@@ -240,9 +243,9 @@ struct Manifold::Impl {
     } while (current != halfedge);
   }
 
-  void CreateFaces(const std::vector<double>& propertyTolerance = {});
+  void CreateFaces();
   void RemoveUnreferencedVerts();
-  void InitializeOriginal();
+  void InitializeOriginal(bool keepFaceID = false);
   void CreateHalfedges(const Vec<ivec3>& triVerts);
   void CalculateNormals();
   void IncrementMeshIDs();

src/manifold.cpp

@@ -282,13 +282,6 @@ Manifold::Manifold(const MeshGL& meshGL)
  * runs.
  *
  * @param meshGL The input MeshGL.
- * @param propertyTolerance A vector of precision values for each property
- * beyond position. If specified, the propertyTolerance vector must have size =
- * numProp - 3. This is the amount of interpolation error allowed before two
- * neighboring triangles are considered to be on a property boundary edge.
- * Property boundary edges will be retained across operations even if the
- * triangles are coplanar. Defaults to 1e-5, which works well for most
- * properties in the [-1, 1] range.
  */
 Manifold::Manifold(const MeshGL64& meshGL64)
     : pNode_(std::make_shared<CsgLeafNode>(std::make_shared<Impl>(meshGL64))) {}
@@ -418,21 +411,12 @@ int Manifold::OriginalID() const {
 }
 
 /**
- * This function condenses all coplanar faces in the relation, and
- * collapses those edges. In the process the relation to ancestor meshes is lost
- * and this new Manifold is marked an original. Properties are preserved, so if
- * they do not match across an edge, that edge will be kept.
- *
- * @param propertyTolerance A vector of precision values for each property
- * beyond position. If specified, the propertyTolerance vector must have size =
- * numProp - 3. This is the amount of interpolation error allowed before two
- * neighboring triangles are considered to be on a property boundary edge.
- * Property boundary edges will be retained across operations even if the
- * triangles are coplanar. Defaults to 1e-5, which works well for most
- * single-precision properties in the [-1, 1] range.
+ * This removes all relations (originalID, faceID, transform) to ancestor meshes
+ * and this new Manifold is marked an original. It also collapses colinear edges
+ * - these don't get collapsed at boundaries where originalID changes, so the
+ * reset may allow flat faces to be further simplified.
  */
-Manifold Manifold::AsOriginal(
-    const std::vector<double>& propertyTolerance) const {
+Manifold Manifold::AsOriginal() const {
   auto oldImpl = GetCsgLeafNode().GetImpl();
   if (oldImpl->status_ != Error::NoError) {
     auto newImpl = std::make_shared<Impl>();
@@ -441,9 +425,10 @@ Manifold Manifold::AsOriginal(
   }
   auto newImpl = std::make_shared<Impl>(*oldImpl);
   newImpl->InitializeOriginal();
-  newImpl->CreateFaces(propertyTolerance);
+  newImpl->CreateFaces();
   newImpl->SimplifyTopology();
   newImpl->Finish();
+  newImpl->InitializeOriginal(true);
   return Manifold(std::make_shared<CsgLeafNode>(newImpl));
 }
 

src/quickhull.cpp

@@ -851,6 +851,7 @@ void Manifold::Impl::Hull(VecView<vec3> vertPos) {
   CalculateNormals();
   InitializeOriginal();
   Finish();
+  CreateFaces();
 }
 
 }  // namespace manifold

src/shared.h

@@ -143,12 +143,14 @@ struct TriRef {
   /// The OriginalID of the mesh this triangle came from. This ID is ideal for
   /// reapplying properties like UV coordinates to the output mesh.
   int originalID;
-  /// The triangle index of the original triangle this was part of:
-  /// Mesh.triVerts[tri].
+  /// Probably the triangle index of the original triangle this was part of:
+  /// Mesh.triVerts[tri], but it's an input, so just pass it along unchanged.
   int tri;
+  /// Triangles with the same face ID are coplanar.
+  int faceID;
 
   bool SameFace(const TriRef& other) const {
-    return meshID == other.meshID && tri == other.tri;
+    return meshID == other.meshID && faceID == other.faceID;
   }
 };
 
@@ -214,7 +216,8 @@ inline std::ostream& operator<<(std::ostream& stream, const Barycentric& bary) {
 
 inline std::ostream& operator<<(std::ostream& stream, const TriRef& ref) {
   return stream << "meshID: " << ref.meshID
-                << ", originalID: " << ref.originalID << ", tri: " << ref.tri;
+                << ", originalID: " << ref.originalID << ", tri: " << ref.tri
+                << ", faceID: " << ref.faceID;
 }
 #endif
 }  // namespace manifold

src/smoothing.cpp

@@ -986,13 +986,12 @@ void Manifold::Impl::Refine(std::function<int(vec3, vec4, vec4)> edgeDivisions,
   if (old.halfedgeTangent_.size() == old.halfedge_.size()) {
     for_each_n(autoPolicy(NumTri(), 1e4), countAt(0), NumVert(),
                InterpTri({vertPos_, vertBary, &old}));
-    // Make original since the subdivided faces have been warped into
-    // being non-coplanar, and hence not being related to the original faces.
-    InitializeOriginal();
   }
 
   halfedgeTangent_.resize(0);
   Finish();
+  CreateFaces();
+  meshRelation_.originalID = -1;
 }
 
 }  // namespace manifold

test/boolean_complex_test.cpp

@@ -28,9 +28,8 @@ using namespace manifold;
  */
 
 TEST(BooleanComplex, Sphere) {
-  Manifold sphere = Manifold::Sphere(1.0, 12);
-  MeshGL sphereGL = WithPositionColors(sphere);
-  sphere = Manifold(sphereGL);
+  Manifold sphere = WithPositionColors(Manifold::Sphere(1.0, 12));
+  MeshGL sphereGL = sphere.GetMeshGL();
 
   Manifold sphere2 = sphere.Translate(vec3(0.5));
   Manifold result = sphere - sphere2;
@@ -52,8 +51,8 @@ TEST(BooleanComplex, Sphere) {
 }
 
 TEST(BooleanComplex, MeshRelation) {
-  MeshGL gyroidMeshGL = WithPositionColors(Gyroid());
-  Manifold gyroid(gyroidMeshGL);
+  Manifold gyroid = WithPositionColors(Gyroid()).AsOriginal();
+  MeshGL gyroidMeshGL = gyroid.GetMeshGL();
 
   Manifold gyroid2 = gyroid.Translate(vec3(2.0));