diff --git a/bindings/module.cpp b/bindings/module.cpp index d31e70b..1f1f64c 100644 --- a/bindings/module.cpp +++ b/bindings/module.cpp @@ -247,6 +247,24 @@ nb::list csgTreeToPy(const std::vector>& tre // Exposed so a front end that writes its own files -- BelfrySCAD's GUI has // its own exporters -- can run the same check the CLI does rather than a // second, drifting implementation. +// Exposed for the export path: strip zero-area faces and repair the +// T-joints their removal exposes. Returns (verts, tris, report). +nb::tuple stripSliversPy(const std::vector& verts, const std::vector& tris) { + manifold::MeshGL m; + m.numProp = 3; + m.vertProperties = verts; + m.triVerts = tris; + oscadeval::SliverStripReport r; + manifold::MeshGL out = oscadeval::stripSlivers(m, r); + nb::dict rep; + rep["removed"] = r.removed; + rep["restitched"] = r.restitched; + rep["needles"] = r.needles; + rep["left_behind"] = r.leftBehind; + rep["passes"] = r.passes; + return nb::make_tuple(out.vertProperties, out.triVerts, rep); +} + nb::dict checkMeshPy(const std::vector& verts, const std::vector& tris) { manifold::MeshGL m; m.numProp = 3; @@ -708,6 +726,9 @@ NB_MODULE(_openscad_cpp_evaluator, m) { "Evaluate a .scad file; return (bodies, echoes, id_to_node, csg_tree, profile_result, dyn, dyn_explicit)."); m.def("parse_decls", &parseDecls, nb::arg("path"), "Parse a .scad file; return top-level declaration (namespace, name, start, end, line, column, origin) tuples."); + m.def("strip_slivers", &stripSliversPy, nb::arg("verts"), nb::arg("tris"), + "Remove zero-area faces and repair the T-joints their removal " + "exposes. Returns (verts, tris, report)."); m.def("check_mesh", &checkMeshPy, nb::arg("verts"), nb::arg("tris"), "Diagnose a triangle mesh against the manifoldness conditions. " "verts is a flat [x,y,z,...] list, tris a flat index list. Returns " diff --git a/include/openscad_cpp_evaluator/mesh_check.hpp b/include/openscad_cpp_evaluator/mesh_check.hpp index 0c5dfaf..52b1315 100644 --- a/include/openscad_cpp_evaluator/mesh_check.hpp +++ b/include/openscad_cpp_evaluator/mesh_check.hpp @@ -63,10 +63,11 @@ struct MeshRepairReport { size_t filledTriangles = 0; size_t splitVertices = 0; size_t unfilledHoles = 0; // boundary loops it could not close + size_t strippedSlivers = 0; // zero-area faces removed and restitched bool didAnything() const { return weldedVertices || droppedDegenerate || droppedDuplicate || reversedFaces - || filledHoles || splitVertices; + || filledHoles || splitVertices || strippedSlivers; } // "welded 12 vertices, filled 1 hole (4 triangles)" -- empty if nothing. std::string summary() const; @@ -82,4 +83,26 @@ struct MeshRepairReport { // improved as far as it got, with the remainder reported in `report`. manifold::MeshGL repairMesh(const manifold::MeshGL& mesh, MeshRepairReport& report); +struct SliverStripReport { + size_t removed = 0; // zero-area faces taken out + size_t restitched = 0; // neighbours split to close the gap they left + // Faces with two corners at one position. These need no restitching: + // their two long edges run between the same two points, so the faces on + // either side already meet once the needle is gone. + size_t needles = 0; + size_t leftBehind = 0; // slivers whose neighbour could not be found + size_t passes = 0; // removing one can expose another +}; + +// Remove zero-area faces and repair the T-joints their removal exposes. +// +// A sliver's three vertices are collinear, so one lies between the other +// two. Dropping the face leaves that middle vertex sitting on the interior +// of the neighbour's edge -- a T-joint -- and the two sides no longer share +// an edge, which reads as a hole. Splitting the neighbour at the middle +// vertex restores the match without moving any geometry. +// +// The mesh is returned unchanged if it has no slivers. +manifold::MeshGL stripSlivers(const manifold::MeshGL& mesh, SliverStripReport& report); + } // namespace oscadeval diff --git a/pyproject.toml b/pyproject.toml index 465dc77..4c63aca 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -4,7 +4,7 @@ build-backend = "scikit_build_core.build" [project] name = "openscad_cpp_evaluator" -version = "0.24.0" +version = "0.25.0" description = "C++ OpenSCAD evaluator with Python bindings" readme = "README.md" requires-python = ">=3.12" diff --git a/python/openscad_cpp_evaluator/__init__.py b/python/openscad_cpp_evaluator/__init__.py index 58808b8..f51f6b6 100644 --- a/python/openscad_cpp_evaluator/__init__.py +++ b/python/openscad_cpp_evaluator/__init__.py @@ -18,7 +18,7 @@ __all__ = [ "Evaluator", "ColoredBody", "EvalError", "ParseError", "OscObject", "parse", "to_renderable_bodies", "ManifoldCache", "CallSiteProfile", "ProfileResult", "format_csg_tree", "bodies_from_dicts", - "FastContinueSignal", "parse_ast", "parse_ast_string", "check_mesh", + "FastContinueSignal", "parse_ast", "parse_ast_string", "check_mesh", "strip_slivers", ] @@ -274,6 +274,20 @@ def parse(path: str) -> RootScope: return RootScope(decls) +def strip_slivers(verts, tris): + """Remove zero-area faces and repair the T-joints removal exposes. + + A sliver's three vertices are collinear, so one lies between the other + two. Dropping the face leaves that middle vertex on the interior of the + neighbour's edge, and the two sides no longer share an edge -- which + reads as a hole. Splitting the neighbour there restores the match + without moving any geometry. + + Returns (verts, tris, report). + """ + return _ext.strip_slivers([float(v) for v in verts], [int(t) for t in tris]) + + def check_mesh(verts, tris) -> dict: """Diagnose a triangle mesh against the manifoldness conditions. diff --git a/src/mesh_check.cpp b/src/mesh_check.cpp index 44cf77f..6577c65 100644 --- a/src/mesh_check.cpp +++ b/src/mesh_check.cpp @@ -204,6 +204,8 @@ std::string MeshRepairReport::summary() const { + (filledHoles == 1 ? " hole filled (" : " holes filled (") + std::to_string(filledTriangles) + " triangles)"); } + add(strippedSlivers, strippedSlivers == 1 ? "zero-area face stripped" + : "zero-area faces stripped"); add(unfilledHoles, unfilledHoles == 1 ? "hole left open" : "holes left open"); std::string out; for (size_t i = 0; i < parts.size(); ++i) { @@ -426,6 +428,206 @@ manifold::MeshGL repairMesh(const manifold::MeshGL& mesh, MeshRepairReport& repo out.runTransform.clear(); out.mergeFromVert.clear(); out.mergeToVert.clear(); + + // Zero-area faces last. Welding already disposes of needles -- it + // collapses the coincident pair, leaving a face that names a vertex + // twice, which step 2 drops -- but a sliver whose three corners are + // distinct and collinear survives all of the above untouched, and an + // imported STL is exactly where those turn up. + SliverStripReport strip; + out = stripSlivers(out, strip); + report.strippedSlivers = strip.removed; + return out; +} + + +namespace { + +// Which two corners of a zero-area face sit at the same position, if any. +// Such a face is a needle rather than a T-joint: its two long edges run +// between the same two points, so once it is gone the faces on either side +// already meet along one edge and nothing needs splitting. The coincident +// pair is merged so that is true by index as well as by position. +bool coincidentPair(const manifold::MeshGL& m, const std::array& f, + Vert& keep, Vert& drop) { + for (int i = 0; i < 3; ++i) { + const Vert a = f[i], b = f[(i + 1) % 3]; + if (a == b) { keep = a; drop = b; return true; } + double p[3], q[3]; + pos(m, a, p); pos(m, b, q); + if (llround(p[0] * 1e6) == llround(q[0] * 1e6) + && llround(p[1] * 1e6) == llround(q[1] * 1e6) + && llround(p[2] * 1e6) == llround(q[2] * 1e6)) { + keep = std::min(a, b); + drop = std::max(a, b); + return keep != drop; + } + } + return false; +} + +// Index of the vertex opposite the longest edge -- for three collinear +// points that is the one in the middle. +int middleOfCollinear(const manifold::MeshGL& m, const std::array& f) { + double p[3][3]; + for (int i = 0; i < 3; ++i) pos(m, f[i], p[i]); + auto d2 = [&](int a, int b) { + double s = 0; + for (int k = 0; k < 3; ++k) { const double d = p[a][k] - p[b][k]; s += d * d; } + return s; + }; + const double opp[3] = {d2(1, 2), d2(0, 2), d2(0, 1)}; + int best = 0; + for (int i = 1; i < 3; ++i) if (opp[i] > opp[best]) best = i; + return best; +} + +} // namespace + +manifold::MeshGL stripSlivers(const manifold::MeshGL& mesh, SliverStripReport& report) { + report = SliverStripReport{}; + std::vector> tris; + tris.reserve(triCount(mesh)); + for (size_t t = 0; t < triCount(mesh); ++t) { + Vert v[3]; + triVerts(mesh, t, v); + tris.push_back({v[0], v[1], v[2]}); + } + + // Removing a sliver can leave its neighbour split into pieces that are + // themselves slivers, so this repeats. Bounded because each pass must + // remove at least one face to continue. + const int kMaxPasses = 12; + for (int pass = 0; pass < kMaxPasses; ++pass) { + std::vector slivers; + for (size_t i = 0; i < tris.size(); ++i) { + const auto& f = tris[i]; + if (f[0] == f[1] || f[1] == f[2] || f[0] == f[2]) continue; // handled elsewhere + double a[3], b[3], c[3]; + pos(mesh, f[0], a); pos(mesh, f[1], b); pos(mesh, f[2], c); + const double ux = b[0] - a[0], uy = b[1] - a[1], uz = b[2] - a[2]; + const double wx = c[0] - a[0], wy = c[1] - a[1], wz = c[2] - a[2]; + const double cx = uy * wz - uz * wy, cy = uz * wx - ux * wz, cz = ux * wy - uy * wx; + if ((cx * cx + cy * cy + cz * cz) < 1e-24) slivers.push_back(i); + } + if (slivers.empty()) break; + ++report.passes; + + // Which face owns each edge, so a sliver's long-edge neighbour can + // be found. Rebuilt per pass: splitting changes it. + std::map> owners; + for (size_t i = 0; i < tris.size(); ++i) { + for (int e = 0; e < 3; ++e) { + owners[undirected(tris[i][e], tris[i][(e + 1) % 3])].push_back(i); + } + } + + std::vector isSliver(tris.size(), 0); + for (size_t i : slivers) isSliver[i] = 1; + + std::vector dead(tris.size(), 0); + std::vector> added; + std::map merge; // needle corners to fold together + for (size_t si : slivers) { + if (dead[si]) continue; + const auto f = tris[si]; + + // A needle: two corners at one point. Nothing to restitch -- + // the faces on either side already share an edge positionally, + // and merging the pair makes them share it by index too. + Vert keep = 0, drop = 0; + if (coincidentPair(mesh, f, keep, drop)) { + dead[si] = 1; + if (keep != drop) merge[drop] = keep; + ++report.removed; + ++report.needles; + continue; + } + + const int mid = middleOfCollinear(mesh, f); + const Vert m = f[mid], a = f[(mid + 1) % 3], b = f[(mid + 2) % 3]; + + // The neighbour across the long edge a-b, which is the one the + // middle vertex now sits inside. + // Prefer a neighbour that is not itself a sliver. Two slivers + // sharing their long edge are each other's only candidate, and + // splitting one into the other just moves the problem around -- + // a level-4 Menger sponge has exactly one such pair, and it was + // what stopped the last two from ever clearing. + size_t nb = SIZE_MAX, fallback = SIZE_MAX; + for (size_t cand : owners[undirected(a, b)]) { + if (cand == si || dead[cand]) continue; + if (isSliver[cand]) { + if (fallback == SIZE_MAX) fallback = cand; + continue; + } + nb = cand; + break; + } + // A non-sliver neighbour is preferred but not required. Two + // slivers sharing their long edge are each other's only + // candidate -- a level-4 Menger sponge has one such pair, and + // refusing to split into a sliver leaves them forever. Splitting + // into one still makes progress, because the halves are smaller + // and the next pass reconsiders them. + if (nb == SIZE_MAX) nb = fallback; + if (nb == SIZE_MAX) { ++report.leftBehind; continue; } + + // Split the neighbour at m, keeping its winding: the edge a-b + // appears in it in some direction, and the two pieces must + // traverse it the same way round. + const auto& n = tris[nb]; + int e = -1; + for (int i = 0; i < 3; ++i) { + const Vert x = n[i], y = n[(i + 1) % 3]; + if ((x == a && y == b) || (x == b && y == a)) { e = i; break; } + } + if (e < 0) { ++report.leftBehind; continue; } + const Vert x = n[e], y = n[(e + 1) % 3], apex = n[(e + 2) % 3]; + + dead[si] = 1; + dead[nb] = 1; + added.push_back({x, m, apex}); + added.push_back({m, y, apex}); + ++report.removed; + ++report.restitched; + } + + std::vector> next; + next.reserve(tris.size() + added.size()); + auto resolve = [&](Vert v) { + for (int hop = 0; hop < 8; ++hop) { // chains are short + auto it = merge.find(v); + if (it == merge.end()) break; + v = it->second; + } + return v; + }; + for (size_t i = 0; i < tris.size(); ++i) { + if (dead[i]) continue; + const auto& f = tris[i]; + next.push_back({resolve(f[0]), resolve(f[1]), resolve(f[2])}); + } + for (const auto& f : added) { + next.push_back({resolve(f[0]), resolve(f[1]), resolve(f[2])}); + } + if (next.size() == tris.size() && added.empty()) break; + tris.swap(next); + } + + manifold::MeshGL out = mesh; + out.triVerts.clear(); + out.triVerts.reserve(tris.size() * 3); + for (const auto& f : tris) { + out.triVerts.push_back(f[0]); + out.triVerts.push_back(f[1]); + out.triVerts.push_back(f[2]); + } + // These describe the old triangle list; a stale one is worse than none. + out.runIndex.clear(); + out.runOriginalID.clear(); + out.faceID.clear(); + out.runTransform.clear(); return out; } diff --git a/tests/test_mesh_check.cpp b/tests/test_mesh_check.cpp index 6d41e3c..2280d42 100644 --- a/tests/test_mesh_check.cpp +++ b/tests/test_mesh_check.cpp @@ -249,3 +249,130 @@ TEST(MeshRepair, ReportsWhatItDid) { repairMesh(m, r); EXPECT_NE(r.summary().find("hole filled"), std::string::npos) << r.summary(); } + +// -- sliver stripping ----------------------------------------------------- + +// The cube whose top face is fanned through a point on its diagonal, from +// the check tests above: manifold, with one zero-area face. +namespace { +manifold::MeshGL cubeWithSliver() { + return build( + {0, 0, 0, 1, 0, 0, 1, 1, 0, 0, 1, 0, + 0, 0, 1, 1, 0, 1, 1, 1, 1, 0, 1, 1, + 0.5f, 0.5f, 1}, + { + 0, 2, 1, 0, 3, 2, + 4, 5, 8, 5, 6, 8, 4, 8, 6, + 4, 6, 7, + 0, 1, 5, 0, 5, 4, + 1, 2, 6, 1, 6, 5, + 2, 3, 7, 2, 7, 6, + 3, 0, 4, 3, 4, 7, + }); +} + +double signedVolume(const manifold::MeshGL& m) { + double vol = 0; + for (size_t t = 0; t + 2 < m.triVerts.size(); t += 3) { + const float* a = &m.vertProperties[m.triVerts[t] * 3]; + const float* b = &m.vertProperties[m.triVerts[t + 1] * 3]; + const float* c = &m.vertProperties[m.triVerts[t + 2] * 3]; + vol += (a[0] * (b[1] * c[2] - b[2] * c[1]) + - a[1] * (b[0] * c[2] - b[2] * c[0]) + + a[2] * (b[0] * c[1] - b[1] * c[0])) / 6.0; + } + return vol; +} +} // namespace + +TEST(SliverStrip, RemovesTheSliverAndRestitchesTheTJoint) { + const manifold::MeshGL m = cubeWithSliver(); + ASSERT_EQ(checkMesh(m).degenerateFaces, 1u); + ASSERT_TRUE(checkMesh(m).manifold()); + + SliverStripReport r; + const manifold::MeshGL out = stripSlivers(m, r); + EXPECT_EQ(r.removed, 1u) << "sliver not removed"; + EXPECT_EQ(r.restitched, 1u) << "neighbour not split"; + + const MeshDiagnosis d = checkMesh(out); + EXPECT_EQ(d.degenerateFaces, 0u) << d.summary(); + // The point of restitching: simply deleting the face would leave the + // middle vertex sitting on the neighbour's edge, and three holes. + EXPECT_EQ(d.boundaryEdges, 0u) << "deleting alone left holes: " << d.summary(); + EXPECT_TRUE(d.manifold()) << d.summary(); +} + +TEST(SliverStrip, MovesNoGeometry) { + // A retriangulation, not a repair: the solid must be identical. + const manifold::MeshGL m = cubeWithSliver(); + SliverStripReport r; + const manifold::MeshGL out = stripSlivers(m, r); + EXPECT_NEAR(signedVolume(out), signedVolume(m), 1e-9); + // Two faces out -- the sliver and the neighbour it was stuck to -- and + // the neighbour's two halves in. The count is unchanged. + EXPECT_EQ(out.triVerts.size() / 3, m.triVerts.size() / 3); +} + +TEST(SliverStrip, LeavesACleanMeshUntouched) { + SliverStripReport r; + const manifold::MeshGL out = stripSlivers(tetra(), r); + EXPECT_EQ(r.removed, 0u); + EXPECT_EQ(r.passes, 0u) << "walked a mesh with nothing to do"; + EXPECT_EQ(out.triVerts, tetra().triVerts); +} + +// A zero-area face with two corners at one point is a needle, not a +// T-joint: its two long edges run between the same pair of points, so the +// faces on either side already meet once it is gone. Splitting a neighbour +// for it would be wrong -- there is no middle vertex to split at. +TEST(SliverStrip, ANeedleIsRemovedWithoutSplittingAnything) { + manifold::MeshGL m = tetra(); + // Vertex 4 sits exactly on vertex 1, and the needle {0,1,4} has two + // corners at that one point. + m.vertProperties.insert(m.vertProperties.end(), {1, 0, 0}); + m.triVerts.insert(m.triVerts.end(), {0u, 1u, 4u}); + ASSERT_GT(checkMesh(m).degenerateFaces, 0u); + + SliverStripReport r; + const manifold::MeshGL out = stripSlivers(m, r); + EXPECT_EQ(r.needles, 1u) << "not recognised as a needle: " << r.removed; + EXPECT_EQ(r.restitched, 0u) << "split a neighbour it did not need to"; + + const MeshDiagnosis d = checkMesh(out); + EXPECT_EQ(d.degenerateFaces, 0u) << d.summary(); + EXPECT_EQ(d.boundaryEdges, 0u) << "removal left holes: " << d.summary(); + EXPECT_TRUE(d.manifold()) << d.summary(); +} + +// import(repair=true) goes through repairMesh, so it has to deal with +// slivers too -- an imported STL is exactly where they turn up. +TEST(MeshRepair, StripsSliversAndRestitchesThem) { + const manifold::MeshGL m = cubeWithSliver(); + ASSERT_EQ(checkMesh(m).degenerateFaces, 1u); + + MeshRepairReport r; + const manifold::MeshGL out = repairMesh(m, r); + const MeshDiagnosis d = checkMesh(out); + EXPECT_EQ(d.degenerateFaces, 0u) << "sliver survived repair: " << d.summary(); + EXPECT_EQ(d.boundaryEdges, 0u) << "removal left holes: " << d.summary(); + EXPECT_TRUE(d.manifold()) << d.summary(); +} + +// A needle needs no restitching, and welding alone already handles it: +// collapsing the coincident pair turns the face into one naming a vertex +// twice, which repair drops. Worth pinning, because it is the reason +// repair coped with needles before it knew about slivers at all. +TEST(MeshRepair, WeldingAloneDisposesOfANeedle) { + manifold::MeshGL m = tetra(); + m.vertProperties.insert(m.vertProperties.end(), {1, 0, 0}); // == vertex 1 + m.triVerts.insert(m.triVerts.end(), {0u, 1u, 4u}); + + MeshRepairReport r; + const manifold::MeshGL out = repairMesh(m, r); + EXPECT_GT(r.weldedVertices, 0u) << r.summary(); + EXPECT_GT(r.droppedDegenerate, 0u) << r.summary(); + const MeshDiagnosis d = checkMesh(out); + EXPECT_EQ(d.degenerateFaces, 0u) << d.summary(); + EXPECT_TRUE(d.manifold()) << d.summary(); +}