Merge branch 'tm_openvdb_integration' into lm_hollow_gizmo

CMakeLists.txt

@@ -403,8 +403,10 @@ if(SLIC3R_STATIC)
     set(USE_BLOSC TRUE)
 endif()
 
-#find_package(OpenVDB 5.0 COMPONENTS openvdb)
+find_package(OpenVDB 5.0 COMPONENTS openvdb)
-#slic3r_remap_configs(IlmBase::Half RelWithDebInfo Release)
+if(OpenVDB_FOUND)
+    slic3r_remap_configs(IlmBase::Half RelWithDebInfo Release)
+endif()
 
 # libslic3r, PrusaSlicer GUI and the PrusaSlicer executable.
 add_subdirectory(src)

cmake/modules/FindOpenVDB.cmake

@@ -108,6 +108,18 @@ if(POLICY CMP0074)
   cmake_policy(SET CMP0074 NEW)
 endif()
 
+if(OpenVDB_FIND_QUIETLY)
+  set (_quiet "QUIET")
+else()
+  set (_quiet "")
+endif()
+
+if(OpenVDB_FIND_REQUIRED)
+  set (_required "REQUIRED")
+else()
+  set (_required "")
+endif()
+
 # Include utility functions for version information
 include(${CMAKE_CURRENT_LIST_DIR}/OpenVDBUtils.cmake)
 
@@ -146,7 +158,7 @@ set(_OPENVDB_ROOT_SEARCH_DIR "")
 
 # Additionally try and use pkconfig to find OpenVDB
 
-find_package(PkgConfig)
+find_package(PkgConfig ${_quiet} ${_required})
 pkg_check_modules(PC_OpenVDB QUIET OpenVDB)
 
 # ------------------------------------------------------------------------
@@ -250,7 +262,7 @@ OPENVDB_ABI_VERSION_FROM_PRINT(
   ABI OpenVDB_ABI
 )
 
-if(NOT OpenVDB_FIND_QUIET)
+if(NOT OpenVDB_FIND_QUIETLY)
   if(NOT OpenVDB_ABI)
     message(WARNING "Unable to determine OpenVDB ABI version from OpenVDB "
       "installation. The library major version \"${OpenVDB_MAJOR_VERSION}\" "
@@ -268,7 +280,17 @@ endif()
 
 # Add standard dependencies
 
-find_package(IlmBase COMPONENTS Half)
+macro(just_fail msg)
+  set(OpenVDB_FOUND FALSE)
+  if(OpenVDB_FIND_REQUIRED)
+    message(FATAL_ERROR msg)
+  elseif(NOT OpenVDB_FIND_QUIETLY)
+    message(ERROR msg)
+  endif()
+  return()
+endmacro()
+
+find_package(IlmBase QUIET COMPONENTS Half)
 if(NOT IlmBase_FOUND)
   pkg_check_modules(IlmBase QUIET IlmBase)
 endif()
@@ -276,20 +298,20 @@ if (IlmBase_FOUND AND NOT TARGET IlmBase::Half)
   message(STATUS "Falling back to IlmBase found by pkg-config...")
 
   find_library(IlmHalf_LIBRARY NAMES Half)
-  if(IlmHalf_LIBRARY-NOTFOUND)
+  if(IlmHalf_LIBRARY-NOTFOUND OR NOT IlmBase_INCLUDE_DIRS)
-    message(FATAL_ERROR "IlmBase::Half can not be found!")
+    just_fail("IlmBase::Half can not be found!")
   endif()
   
   add_library(IlmBase::Half UNKNOWN IMPORTED)
   set_target_properties(IlmBase::Half PROPERTIES
     IMPORTED_LOCATION "${IlmHalf_LIBRARY}"
-    INTERFACE_INCLUDE_DIRECTORIES ${IlmBase_INCLUDE_DIRS})
+    INTERFACE_INCLUDE_DIRECTORIES "${IlmBase_INCLUDE_DIRS}")
 elseif(NOT IlmBase_FOUND)
-  message(FATAL_ERROR "IlmBase::Half can not be found!")
+  just_fail("IlmBase::Half can not be found!")
 endif()
-find_package(TBB REQUIRED COMPONENTS tbb)
+find_package(TBB ${_quiet} ${_required} COMPONENTS tbb)
-find_package(ZLIB REQUIRED)
+find_package(ZLIB ${_quiet} ${_required})
-find_package(Boost REQUIRED COMPONENTS iostreams system)
+find_package(Boost ${_quiet} ${_required} COMPONENTS iostreams system )
 
 # Use GetPrerequisites to see which libraries this OpenVDB lib has linked to
 # which we can query for optional deps. This basically runs ldd/otoll/objdump
@@ -350,7 +372,7 @@ unset(_OPENVDB_PREREQUISITE_LIST)
 unset(_HAS_DEP)
 
 if(OpenVDB_USES_BLOSC)
-  find_package(Blosc )
+  find_package(Blosc QUIET)
   if(NOT Blosc_FOUND OR NOT TARGET Blosc::blosc) 
     message(STATUS "find_package could not find Blosc. Using fallback blosc search...")
     find_path(Blosc_INCLUDE_DIR blosc.h)
@@ -362,25 +384,25 @@ if(OpenVDB_USES_BLOSC)
         IMPORTED_LOCATION "${Blosc_LIBRARY}"
         INTERFACE_INCLUDE_DIRECTORIES ${Blosc_INCLUDE_DIR})
     elseif()
-      message(FATAL_ERROR "Blosc library can not be found!")
+      just_fail("Blosc library can not be found!")
     endif()
   endif()
 endif()
 
 if(OpenVDB_USES_LOG4CPLUS)
-  find_package(Log4cplus REQUIRED)
+  find_package(Log4cplus ${_quiet} ${_required})
 endif()
 
 if(OpenVDB_USES_ILM)
-  find_package(IlmBase REQUIRED)
+  find_package(IlmBase ${_quiet} ${_required})
 endif()
 
 if(OpenVDB_USES_EXR)
-  find_package(OpenEXR REQUIRED)
+  find_package(OpenEXR ${_quiet} ${_required})
 endif()
 
 if(UNIX)
-  find_package(Threads REQUIRED)
+  find_package(Threads ${_quiet} ${_required})
 endif()
 
 # Set deps. Note that the order here is important. If we're building against
@@ -481,7 +503,7 @@ foreach(COMPONENT ${OpenVDB_FIND_COMPONENTS})
   endif()
 endforeach()
 
-if(OpenVDB_FOUND AND NOT ${CMAKE_FIND_PACKAGE_NAME}_FIND_QUIETLY)
+if(OpenVDB_FOUND AND NOT OpenVDB_FIND_QUIETLY)
   message(STATUS "OpenVDB libraries: ${OpenVDB_LIBRARIES}")
 endif()
 

cmake/modules/OpenVDBUtils.cmake

@@ -125,7 +125,9 @@ function(OPENVDB_ABI_VERSION_FROM_PRINT OPENVDB_PRINT)
   cmake_parse_arguments(_VDB "QUIET" "ABI" "" ${ARGN})
 
   if(NOT EXISTS ${OPENVDB_PRINT})
+    if(NOT OpenVDB_FIND_QUIETLY)
       message(WARNING "vdb_print not found! ${OPENVDB_PRINT}")
+    endif()
     return()
   endif()
 
@@ -148,7 +150,9 @@ function(OPENVDB_ABI_VERSION_FROM_PRINT OPENVDB_PRINT)
   endif()
 
   if(${_VDB_PRINT_RETURN_STATUS})
+    if(NOT OpenVDB_FIND_QUIETLY)
       message(WARNING "vdb_print returned with status ${_VDB_PRINT_RETURN_STATUS}")
+    endif()
     return()
   endif()
 

sandboxes/openvdb/CMakeLists.txt

@@ -1,2 +1,7 @@
+if(TARGET OpenVDB::openvdb)
     add_executable(openvdb_example openvdb_example.cpp)
+
     target_link_libraries(openvdb_example libslic3r)
+    target_link_libraries(openvdb_example OpenVDB::openvdb)
+endif()
+

src/libslic3r/CMakeLists.txt

@@ -9,6 +9,11 @@ if (MINGW)
     add_compile_options(-Wa,-mbig-obj)
 endif ()
 
+set(OpenVDBUtils_SOURCES "")
+if (TARGET OpenVDB::openvdb)
+    set(OpenVDBUtils_SOURCES OpenVDBUtils.cpp OpenVDBUtils.hpp)
+endif()
+
 add_library(libslic3r STATIC
     pchheader.cpp
     pchheader.hpp
@@ -176,7 +181,9 @@ add_library(libslic3r STATIC
     MinAreaBoundingBox.cpp
     miniz_extension.hpp
     miniz_extension.cpp
+    ${OpenVDBUtils_SOURCES}
     SLA/SLACommon.hpp
+    SLA/SLACommon.cpp
     SLA/SLABoilerPlate.hpp
     SLA/SLAPad.hpp
     SLA/SLAPad.cpp
@@ -224,10 +231,13 @@ target_link_libraries(libslic3r
     qhull
     semver
     TBB::tbb
-    # OpenVDB::openvdb
     ${CMAKE_DL_LIBS}
     )
 
+if (TARGET OpenVDB::openvdb)
+    target_link_libraries(libslic3r OpenVDB::openvdb)
+endif()
+
 if(WIN32)
     target_link_libraries(libslic3r Psapi.lib)
 endif()

src/libslic3r/Format/objparser.cpp

@@ -355,6 +355,35 @@ bool objparse(const char *path, ObjData &data)
 	return true;
 }
 
+bool objparse(std::istream &stream, ObjData &data)
+{
+    try {
+        char buf[65536 * 2];
+        size_t len = 0;
+        size_t lenPrev = 0;
+        while ((len = size_t(stream.read(buf + lenPrev, 65536).gcount())) != 0) {
+            len += lenPrev;
+            size_t lastLine = 0;
+            for (size_t i = 0; i < len; ++ i)
+                if (buf[i] == '\r' || buf[i] == '\n') {
+                    buf[i] = 0;
+                    char *c = buf + lastLine;
+                    while (*c == ' ' || *c == '\t')
+                        ++ c;
+                    obj_parseline(c, data);
+                    lastLine = i + 1;
+                }
+            lenPrev = len - lastLine;
+            memmove(buf, buf + lastLine, lenPrev);
+        }
+    }
+    catch (std::bad_alloc&) {
+        printf("Out of memory\r\n");
+    }
+    
+    return true;
+}
+
 template<typename T> 
 bool savevector(FILE *pFile, const std::vector<T> &v)
 {

src/libslic3r/Format/objparser.hpp

@@ -3,6 +3,7 @@
 
 #include <string>
 #include <vector>
+#include <istream>
 
 namespace ObjParser {
 
@@ -97,6 +98,7 @@ struct ObjData {
 };
 
 extern bool objparse(const char *path, ObjData &data);
+extern bool objparse(std::istream &stream, ObjData &data);
 
 extern bool objbinsave(const char *path, const ObjData &data);
 

src/libslic3r/OpenVDBUtils.cpp

@@ -0,0 +1,121 @@
+#define NOMINMAX
+#include "OpenVDBUtils.hpp"
+#include <openvdb/tools/MeshToVolume.h>
+#include <openvdb/tools/VolumeToMesh.h>
+
+namespace Slic3r {
+
+class TriangleMeshDataAdapter {
+public:
+    const TriangleMesh &mesh;
+    
+    size_t polygonCount() const { return mesh.its.indices.size(); }
+    size_t pointCount() const   { return mesh.its.vertices.size(); }
+    size_t vertexCount(size_t) const { return 3; }
+    
+    // Return position pos in local grid index space for polygon n and vertex v
+    void getIndexSpacePoint(size_t n, size_t v, openvdb::Vec3d& pos) const;
+};
+
+class Contour3DDataAdapter {
+public:
+    const sla::Contour3D &mesh;
+    
+    size_t polygonCount() const { return mesh.faces3.size() + mesh.faces4.size(); }
+    size_t pointCount() const   { return mesh.points.size(); }
+    size_t vertexCount(size_t n) const { return n < mesh.faces3.size() ? 3 : 4; }
+    
+    // Return position pos in local grid index space for polygon n and vertex v
+    void getIndexSpacePoint(size_t n, size_t v, openvdb::Vec3d& pos) const;
+};
+
+void TriangleMeshDataAdapter::getIndexSpacePoint(size_t          n,
+                                                 size_t          v,
+                                                 openvdb::Vec3d &pos) const
+{
+    auto vidx = size_t(mesh.its.indices[n](Eigen::Index(v)));
+    Slic3r::Vec3d p = mesh.its.vertices[vidx].cast<double>();
+    pos = {p.x(), p.y(), p.z()};
+}
+
+void Contour3DDataAdapter::getIndexSpacePoint(size_t          n,
+                                              size_t          v,
+                                              openvdb::Vec3d &pos) const
+{
+    size_t vidx = 0;
+    if (n < mesh.faces3.size()) vidx = size_t(mesh.faces3[n](Eigen::Index(v)));
+    else vidx = size_t(mesh.faces4[n - mesh.faces3.size()](Eigen::Index(v)));
+    
+    Slic3r::Vec3d p = mesh.points[vidx];
+    pos = {p.x(), p.y(), p.z()};
+}
+
+openvdb::FloatGrid::Ptr meshToVolume(const TriangleMesh &mesh,
+                                     float               exteriorBandWidth,
+                                     float               interiorBandWidth,
+                                     int                 flags,
+                                     const openvdb::math::Transform &tr)
+{
+    openvdb::initialize();
+    return openvdb::tools::meshToVolume<openvdb::FloatGrid>(
+        TriangleMeshDataAdapter{mesh}, tr, exteriorBandWidth,
+        interiorBandWidth, flags);
+}
+
+// TODO: Do I need to call initialize? Seems to work without it as well but the
+// docs say it should be called ones. It does a mutex lock-unlock sequence all
+// even if was called previously.
+
+openvdb::FloatGrid::Ptr meshToVolume(const sla::Contour3D &          mesh,
+                                     float exteriorBandWidth,
+                                     float interiorBandWidth,
+                                     int flags,
+                                     const openvdb::math::Transform &tr)
+{
+    openvdb::initialize();
+    return openvdb::tools::meshToVolume<openvdb::FloatGrid>(
+        Contour3DDataAdapter{mesh}, tr, exteriorBandWidth, interiorBandWidth,
+        flags);
+}
+
+inline Vec3f to_vec3f(const openvdb::Vec3s &v) { return Vec3f{v.x(), v.y(), v.z()}; }
+inline Vec3d to_vec3d(const openvdb::Vec3s &v) { return to_vec3f(v).cast<double>(); }
+inline Vec3i to_vec3i(const openvdb::Vec3I &v) { return Vec3i{int(v[0]), int(v[1]), int(v[2])}; }
+inline Vec4i to_vec4i(const openvdb::Vec4I &v) { return Vec4i{int(v[0]), int(v[1]), int(v[2]), int(v[3])}; }
+
+template<class Grid>
+sla::Contour3D _volumeToMesh(const Grid &grid,
+                             double      isovalue,
+                             double      adaptivity,
+                             bool        relaxDisorientedTriangles)
+{
+    openvdb::initialize();
+    
+    std::vector<openvdb::Vec3s> points;
+    std::vector<openvdb::Vec3I> triangles;
+    std::vector<openvdb::Vec4I> quads;
+    
+    openvdb::tools::volumeToMesh(grid, points, triangles, quads, isovalue,
+                                 adaptivity, relaxDisorientedTriangles);
+    
+    sla::Contour3D ret;
+    ret.points.reserve(points.size());
+    ret.faces3.reserve(triangles.size());
+    ret.faces4.reserve(quads.size());
+    
+    for (auto &v : points) ret.points.emplace_back(to_vec3d(v));
+    for (auto &v : triangles) ret.faces3.emplace_back(to_vec3i(v));
+    for (auto &v : quads) ret.faces4.emplace_back(to_vec4i(v));
+    
+    return ret;
+}
+
+sla::Contour3D volumeToMesh(const openvdb::FloatGrid &grid,
+                            double                    isovalue,
+                            double                    adaptivity,
+                            bool relaxDisorientedTriangles)
+{
+    return _volumeToMesh(grid, isovalue, adaptivity, relaxDisorientedTriangles);
+}
+
+} // namespace Slic3r

src/libslic3r/OpenVDBUtils.hpp

@@ -0,0 +1,29 @@
+#ifndef OPENVDBUTILS_HPP
+#define OPENVDBUTILS_HPP
+
+#include <libslic3r/TriangleMesh.hpp>
+#include <libslic3r/SLA/SLABoilerPlate.hpp>
+#include <openvdb/openvdb.h>
+
+namespace Slic3r {
+
+openvdb::FloatGrid::Ptr meshToVolume(const TriangleMesh &mesh,
+                                     float exteriorBandWidth = 3.0f,
+                                     float interiorBandWidth = 3.0f,
+                                     int   flags             = 0,
+                                     const openvdb::math::Transform &tr = {});
+
+openvdb::FloatGrid::Ptr meshToVolume(const sla::Contour3D &mesh,
+                                     float exteriorBandWidth = 3.0f,
+                                     float interiorBandWidth = 3.0f,
+                                     int   flags             = 0,
+                                     const openvdb::math::Transform &tr = {});
+
+sla::Contour3D volumeToMesh(const openvdb::FloatGrid &grid,
+                            double                    isovalue   = 0.0,
+                            double                    adaptivity = 0.0,
+                            bool relaxDisorientedTriangles       = true);
+
+} // namespace Slic3r
+
+#endif // OPENVDBUTILS_HPP

src/libslic3r/SLA/SLABoilerPlate.hpp

@@ -12,91 +12,11 @@
 #include "SLASpatIndex.hpp"
 
 namespace Slic3r {
+
+typedef Eigen::Matrix<int, 4, 1, Eigen::DontAlign> Vec4i;
+
 namespace sla {
 
-/// Intermediate struct for a 3D mesh
-struct Contour3D {
-    Pointf3s points;
-    std::vector<Vec3i> indices;
-
-    Contour3D& merge(const Contour3D& ctr)
-    {
-        auto s3 = coord_t(points.size());
-        auto s = indices.size();
-
-        points.insert(points.end(), ctr.points.begin(), ctr.points.end());
-        indices.insert(indices.end(), ctr.indices.begin(), ctr.indices.end());
-
-        for(size_t n = s; n < indices.size(); n++) {
-            auto& idx = indices[n]; idx.x() += s3; idx.y() += s3; idx.z() += s3;
-        }
-        
-        return *this;
-    }
-
-    Contour3D& merge(const Pointf3s& triangles)
-    {
-        const size_t offs = points.size();
-        points.insert(points.end(), triangles.begin(), triangles.end());
-        indices.reserve(indices.size() + points.size() / 3);
-        
-        for(int i = int(offs); i < int(points.size()); i += 3)
-            indices.emplace_back(i, i + 1, i + 2);
-        
-        return *this;
-    }
-
-    // Write the index triangle structure to OBJ file for debugging purposes.
-    void to_obj(std::ostream& stream)
-    {
-        for(auto& p : points) {
-            stream << "v " << p.transpose() << "\n";
-        }
-
-        for(auto& f : indices) {
-            stream << "f " << (f + Vec3i(1, 1, 1)).transpose() << "\n";
-        }
-    }
-};
-
-using ClusterEl = std::vector<unsigned>;
-using ClusteredPoints = std::vector<ClusterEl>;
-
-// Clustering a set of points by the given distance.
-ClusteredPoints cluster(const std::vector<unsigned>& indices,
-                        std::function<Vec3d(unsigned)> pointfn,
-                        double dist,
-                        unsigned max_points);
-
-ClusteredPoints cluster(const PointSet& points,
-                        double dist,
-                        unsigned max_points);
-
-ClusteredPoints cluster(
-    const std::vector<unsigned>& indices,
-    std::function<Vec3d(unsigned)> pointfn,
-    std::function<bool(const PointIndexEl&, const PointIndexEl&)> predicate,
-    unsigned max_points);
-
-
-// Calculate the normals for the selected points (from 'points' set) on the
-// mesh. This will call squared distance for each point.
-PointSet normals(const PointSet& points,
-                 const EigenMesh3D& mesh,
-                 double eps = 0.05,  // min distance from edges
-                 std::function<void()> throw_on_cancel = [](){},
-                 const std::vector<unsigned>& selected_points = {});
-
-/// Mesh from an existing contour.
-inline TriangleMesh mesh(const Contour3D& ctour) {
-    return {ctour.points, ctour.indices};
-}
-
-/// Mesh from an evaporating 3D contour
-inline TriangleMesh mesh(Contour3D&& ctour) {
-    return {std::move(ctour.points), std::move(ctour.indices)};
-}
-
 }
 }
 

src/libslic3r/SLA/SLACommon.cpp

@@ -0,0 +1,147 @@
+#include "SLACommon.hpp"
+#include <libslic3r/Format/objparser.hpp>
+
+namespace Slic3r { namespace sla {
+
+Contour3D::Contour3D(const TriangleMesh &trmesh)
+{
+    points.reserve(trmesh.its.vertices.size());
+    faces3.reserve(trmesh.its.indices.size());
+    
+    for (auto &v : trmesh.its.vertices)
+        points.emplace_back(v.cast<double>());
+    
+    std::copy(trmesh.its.indices.begin(), trmesh.its.indices.end(),
+              std::back_inserter(faces3));
+}
+
+Contour3D::Contour3D(TriangleMesh &&trmesh)
+{
+    points.reserve(trmesh.its.vertices.size());
+    
+    for (auto &v : trmesh.its.vertices)
+        points.emplace_back(v.cast<double>());
+    
+    faces3.swap(trmesh.its.indices);
+}
+
+Contour3D::Contour3D(const EigenMesh3D &emesh) {
+    points.reserve(size_t(emesh.V().rows()));
+    faces3.reserve(size_t(emesh.F().rows()));
+    
+    for (int r = 0; r < emesh.V().rows(); r++)
+        points.emplace_back(emesh.V().row(r).cast<double>());
+    
+    for (int i = 0; i < emesh.F().rows(); i++)
+        faces3.emplace_back(emesh.F().row(i));
+}
+
+Contour3D &Contour3D::merge(const Contour3D &ctr)
+{
+    auto N = coord_t(points.size());
+    auto N_f3 = faces3.size();
+    auto N_f4 = faces4.size();
+    
+    points.insert(points.end(), ctr.points.begin(), ctr.points.end());
+    faces3.insert(faces3.end(), ctr.faces3.begin(), ctr.faces3.end());
+    faces4.insert(faces4.end(), ctr.faces4.begin(), ctr.faces4.end());
+    
+    for(size_t n = N_f3; n < faces3.size(); n++) {
+        auto& idx = faces3[n]; idx.x() += N; idx.y() += N; idx.z() += N;
+    }
+    
+    for(size_t n = N_f4; n < faces4.size(); n++) {
+        auto& idx = faces4[n]; for (int k = 0; k < 4; k++) idx(k) += N;
+    }        
+    
+    return *this;
+}
+
+Contour3D &Contour3D::merge(const Pointf3s &triangles)
+{
+    const size_t offs = points.size();
+    points.insert(points.end(), triangles.begin(), triangles.end());
+    faces3.reserve(faces3.size() + points.size() / 3);
+    
+    for(int i = int(offs); i < int(points.size()); i += 3)
+        faces3.emplace_back(i, i + 1, i + 2);
+    
+    return *this;
+}
+
+void Contour3D::to_obj(std::ostream &stream)
+{
+    for(auto& p : points)
+        stream << "v " << p.transpose() << "\n";
+    
+    for(auto& f : faces3) 
+        stream << "f " << (f + Vec3i(1, 1, 1)).transpose() << "\n";
+    
+    for(auto& f : faces4)
+        stream << "f " << (f + Vec4i(1, 1, 1, 1)).transpose() << "\n";
+}
+
+void Contour3D::from_obj(std::istream &stream)
+{
+    ObjParser::ObjData data;
+    ObjParser::objparse(stream, data);
+    
+    points.reserve(data.coordinates.size() / 4 + 1);
+    auto &coords = data.coordinates;
+    for (size_t i = 0; i < coords.size(); i += 4)
+        points.emplace_back(coords[i], coords[i + 1], coords[i + 2]);
+    
+    Vec3i triangle;
+    Vec4i quad;
+    size_t v = 0;
+    while(v < data.vertices.size()) {
+        size_t N = 0;
+        size_t i = v;
+        while (data.vertices[v++].coordIdx != -1) ++N;
+        
+        std::function<void(int, int)> setfn;
+        if (N < 3 || N > 4) continue;
+        else if (N == 3) setfn = [&triangle](int k, int f) { triangle(k) = f; };
+        else setfn = [&quad](int k, int f) { quad(k) = f; };
+        
+        for (size_t j = 0; j < N; ++j)
+            setfn(int(j), data.vertices[i + j].coordIdx);
+    }
+}
+
+TriangleMesh to_triangle_mesh(const Contour3D &ctour) {
+    if (ctour.faces4.empty()) return {ctour.points, ctour.faces3};
+    
+    std::vector<Vec3i> triangles;
+    
+    triangles.reserve(ctour.faces3.size() + 2 * ctour.faces4.size());
+    std::copy(ctour.faces3.begin(), ctour.faces3.end(),
+              std::back_inserter(triangles));
+    
+    for (auto &quad : ctour.faces4) {
+        triangles.emplace_back(quad(0), quad(1), quad(2));
+        triangles.emplace_back(quad(2), quad(3), quad(0));
+    }
+    
+    return {ctour.points, std::move(triangles)};
+}
+
+TriangleMesh to_triangle_mesh(Contour3D &&ctour) {
+    if (ctour.faces4.empty())
+        return {std::move(ctour.points), std::move(ctour.faces3)};
+    
+    std::vector<Vec3i> triangles;
+    
+    triangles.reserve(ctour.faces3.size() + 2 * ctour.faces4.size());
+    std::copy(ctour.faces3.begin(), ctour.faces3.end(),
+              std::back_inserter(triangles));
+    
+    for (auto &quad : ctour.faces4) {
+        triangles.emplace_back(quad(0), quad(1), quad(2));
+        triangles.emplace_back(quad(2), quad(3), quad(0));
+    }
+    
+    return {std::move(ctour.points), std::move(triangles)};
+}
+
+}} // namespace Slic3r::sla

src/libslic3r/SLA/SLACommon.hpp

@@ -5,6 +5,11 @@
 #include <vector>
 #include <Eigen/Geometry>
 
+#include "SLASpatIndex.hpp"
+
+#include <libslic3r/ExPolygon.hpp>
+#include <libslic3r/TriangleMesh.hpp>
+
 // #define SLIC3R_SLA_NEEDS_WINDTREE
 
 namespace Slic3r {
@@ -12,8 +17,7 @@ namespace Slic3r {
 // Typedefs from Point.hpp
 typedef Eigen::Matrix<float, 3, 1, Eigen::DontAlign> Vec3f;
 typedef Eigen::Matrix<double, 3, 1, Eigen::DontAlign> Vec3d;
-
+typedef Eigen::Matrix<int, 4, 1, Eigen::DontAlign> Vec4i;
-class TriangleMesh;
 
 namespace sla {
     
@@ -59,9 +63,11 @@ struct SupportPoint
 
     bool operator==(const SupportPoint &sp) const
     {
-        return (pos == sp.pos) && head_front_radius == sp.head_front_radius &&
+        float rdiff = std::abs(head_front_radius - sp.head_front_radius);
+        return (pos == sp.pos) && rdiff < float(EPSILON) &&
                is_new_island == sp.is_new_island;
     }
+    
     bool operator!=(const SupportPoint &sp) const { return !(sp == (*this)); }
 
     template<class Archive> void serialize(Archive &ar)
@@ -72,8 +78,11 @@ struct SupportPoint
 
 using SupportPoints = std::vector<SupportPoint>;
 
+struct Contour3D;
+
 /// An index-triangle structure for libIGL functions. Also serves as an
-/// alternative (raw) input format for the SLASupportTree
+/// alternative (raw) input format for the SLASupportTree.
+//  Implemented in SLASupportTreeIGL.cpp
 class EigenMesh3D {
     class AABBImpl;
 
@@ -86,6 +95,7 @@ public:
 
     EigenMesh3D(const TriangleMesh&);
     EigenMesh3D(const EigenMesh3D& other);
+    EigenMesh3D(const Contour3D &other);
     EigenMesh3D& operator=(const EigenMesh3D&);
 
     ~EigenMesh3D();
@@ -180,6 +190,63 @@ public:
 
 using PointSet = Eigen::MatrixXd;
 
+
+/// Dumb vertex mesh consisting of triangles (or) quads. Capable of merging with
+/// other meshes of this type and converting to and from other mesh formats.
+struct Contour3D {
+    Pointf3s points;
+    std::vector<Vec3i> faces3;
+    std::vector<Vec4i> faces4;
+    
+    Contour3D() = default;
+    Contour3D(const TriangleMesh &trmesh);
+    Contour3D(TriangleMesh &&trmesh);
+    Contour3D(const EigenMesh3D  &emesh);
+    
+    Contour3D& merge(const Contour3D& ctr);
+    Contour3D& merge(const Pointf3s& triangles);
+    
+    // Write the index triangle structure to OBJ file for debugging purposes.
+    void to_obj(std::ostream& stream);
+    void from_obj(std::istream &stream);
+    
+    inline bool empty() const { return points.empty() || (faces4.empty() && faces3.empty()); }
+};
+
+using ClusterEl = std::vector<unsigned>;
+using ClusteredPoints = std::vector<ClusterEl>;
+
+// Clustering a set of points by the given distance.
+ClusteredPoints cluster(const std::vector<unsigned>& indices,
+                        std::function<Vec3d(unsigned)> pointfn,
+                        double dist,
+                        unsigned max_points);
+
+ClusteredPoints cluster(const PointSet& points,
+                        double dist,
+                        unsigned max_points);
+
+ClusteredPoints cluster(
+    const std::vector<unsigned>& indices,
+    std::function<Vec3d(unsigned)> pointfn,
+    std::function<bool(const PointIndexEl&, const PointIndexEl&)> predicate,
+    unsigned max_points);
+
+
+// Calculate the normals for the selected points (from 'points' set) on the
+// mesh. This will call squared distance for each point.
+PointSet normals(const PointSet& points,
+    const EigenMesh3D& convert_mesh,
+    double eps = 0.05,  // min distance from edges
+    std::function<void()> throw_on_cancel = [](){},
+    const std::vector<unsigned>& selected_points = {});
+
+/// Mesh from an existing contour.
+TriangleMesh to_triangle_mesh(const Contour3D& ctour);
+
+/// Mesh from an evaporating 3D contour
+TriangleMesh to_triangle_mesh(Contour3D&& ctour);
+
 } // namespace sla
 } // namespace Slic3r
 

src/libslic3r/SLA/SLAPad.cpp

@@ -69,7 +69,7 @@ Contour3D walls(
 
     // Shorthand for the vertex arrays
     auto& upts = upper.points, &lpts = lower.points;
-    auto& rpts = ret.points; auto& ind = ret.indices;
+    auto& rpts = ret.points; auto& ind = ret.faces3;
 
     // If the Z levels are flipped, or the offset difference is negative, we
     // will interpret that as the triangles normals should be inverted.
@@ -677,7 +677,7 @@ void create_pad(const ExPolygons &sup_blueprint,
                 ThrowOnCancel thr)
 {
     Contour3D t = create_pad_geometry(sup_blueprint, model_blueprint, cfg, thr);
-    out.merge(mesh(std::move(t)));
+    out.merge(to_triangle_mesh(std::move(t)));
 }
 
 std::string PadConfig::validate() const

src/libslic3r/SLA/SLASupportTreeBuilder.cpp

@@ -12,7 +12,7 @@ Contour3D sphere(double rho, Portion portion, double fa) {
     if(rho <= 1e-6 && rho >= -1e-6) return ret;
     
     auto& vertices = ret.points;
-    auto& facets = ret.indices;
+    auto& facets = ret.faces3;
     
     // Algorithm:
     // Add points one-by-one to the sphere grid and form facets using relative
@@ -102,7 +102,7 @@ Contour3D cylinder(double r, double h, size_t ssteps, const Vec3d &sp)
     
     auto steps = int(ssteps);
     auto& points = ret.points;
-    auto& indices = ret.indices;
+    auto& indices = ret.faces3;
     points.reserve(2*ssteps);
     double a = 2*PI/steps;
     
@@ -211,8 +211,8 @@ Head::Head(double       r_big_mm,
         coord_t i1s1 = coord_t(idx1), i1s2 = coord_t(idx2);
         coord_t i2s1 = i1s1 + 1, i2s2 = i1s2 + 1;
         
-        mesh.indices.emplace_back(i1s1, i2s1, i2s2);
+        mesh.faces3.emplace_back(i1s1, i2s1, i2s2);
-        mesh.indices.emplace_back(i1s1, i2s2, i1s2);
+        mesh.faces3.emplace_back(i1s1, i2s2, i1s2);
     }
     
     auto i1s1 = coord_t(s1.points.size()) - coord_t(steps);
@@ -220,8 +220,8 @@ Head::Head(double       r_big_mm,
     auto i1s2 = coord_t(s1.points.size());
     auto i2s2 = coord_t(s1.points.size()) + coord_t(steps) - 1;
     
-    mesh.indices.emplace_back(i2s2, i2s1, i1s1);
+    mesh.faces3.emplace_back(i2s2, i2s1, i1s1);
-    mesh.indices.emplace_back(i1s2, i2s2, i1s1);
+    mesh.faces3.emplace_back(i1s2, i2s2, i1s1);
     
     // To simplify further processing, we translate the mesh so that the
     // last vertex of the pointing sphere (the pinpoint) will be at (0,0,0)
@@ -240,7 +240,7 @@ Pillar::Pillar(const Vec3d &jp, const Vec3d &endp, double radius, size_t st):
         // move the data.
         Contour3D body = cylinder(radius, height, st, endp);
         mesh.points.swap(body.points);
-        mesh.indices.swap(body.indices);
+        mesh.faces3.swap(body.faces3);
     }
 }
 
@@ -275,7 +275,7 @@ Pillar &Pillar::add_base(double baseheight, double radius)
     base.points.emplace_back(endpt);
     base.points.emplace_back(ep);
     
-    auto& indices = base.indices;
+    auto& indices = base.faces3;
     auto hcenter = int(base.points.size() - 1);
     auto lcenter = int(base.points.size() - 2);
     auto offs = int(steps);
@@ -466,7 +466,7 @@ const TriangleMesh &SupportTreeBuilder::merged_mesh() const
         return m_meshcache;
     }
     
-    m_meshcache = mesh(merged);
+    m_meshcache = to_triangle_mesh(merged);
     
     // The mesh will be passed by const-pointer to TriangleMeshSlicer,
     // which will need this.

src/libslic3r/SLA/SLASupportTreeIGL.cpp

@@ -228,6 +228,26 @@ EigenMesh3D::EigenMesh3D(const EigenMesh3D &other):
     m_V(other.m_V), m_F(other.m_F), m_ground_level(other.m_ground_level),
     m_aabb( new AABBImpl(*other.m_aabb) ) {}
 
+EigenMesh3D::EigenMesh3D(const Contour3D &other)
+{
+    m_V.resize(Eigen::Index(other.points.size()), 3);
+    m_F.resize(Eigen::Index(other.faces3.size() + 2 * other.faces4.size()), 3);
+    
+    for (Eigen::Index i = 0; i < Eigen::Index(other.points.size()); ++i)
+        m_V.row(i) = other.points[size_t(i)];
+    
+    for (Eigen::Index i = 0; i < Eigen::Index(other.faces3.size()); ++i)
+        m_F.row(i) = other.faces3[size_t(i)];
+    
+    size_t N = other.faces3.size() + 2 * other.faces4.size();
+    for (size_t i = other.faces3.size(); i < N; i += 2) {
+        size_t quad_idx = (i - other.faces3.size()) / 2;
+        auto & quad     = other.faces4[quad_idx];
+        m_F.row(Eigen::Index(i)) = Vec3i{quad(0), quad(1), quad(2)};
+        m_F.row(Eigen::Index(i + 1)) = Vec3i{quad(2), quad(3), quad(0)};
+    }
+}
+
 EigenMesh3D &EigenMesh3D::operator=(const EigenMesh3D &other)
 {
     m_V = other.m_V;

tests/example/example_tests_main.cpp

@@ -1,4 +1,3 @@
-#define CATCH_CONFIG_MAIN
 #include <catch_main.hpp>
 
 TEST_CASE("Is example succesful", "[example]") {

tests/libslic3r/CMakeLists.txt

@@ -1,4 +1,5 @@
 get_filename_component(_TEST_NAME ${CMAKE_CURRENT_LIST_DIR} NAME)
+
 add_executable(${_TEST_NAME}_tests 
 	${_TEST_NAME}_tests.cpp
 	test_3mf.cpp
@@ -10,6 +11,11 @@ add_executable(${_TEST_NAME}_tests
 	test_polygon.cpp
 	test_stl.cpp
 	)
+
+if (TARGET OpenVDB::openvdb)
+    target_sources(${_TEST_NAME}_tests PRIVATE test_hollowing.cpp)
+endif()
+    
 target_link_libraries(${_TEST_NAME}_tests test_common libslic3r)
 set_property(TARGET ${_TEST_NAME}_tests PROPERTY FOLDER "tests")
 

tests/libslic3r/test_hollowing.cpp

@@ -0,0 +1,78 @@
+#include <iostream>
+#include <fstream>
+#include <catch2/catch.hpp>
+
+#include "libslic3r/OpenVDBUtils.hpp"
+#include "libslic3r/Format/OBJ.hpp"
+
+#if defined(WIN32) || defined(_WIN32)
+#define PATH_SEPARATOR R"(\)"
+#else
+#define PATH_SEPARATOR R"(/)"
+#endif
+
+static Slic3r::TriangleMesh load_model(const std::string &obj_filename)
+{
+    Slic3r::TriangleMesh mesh;
+    auto fpath = TEST_DATA_DIR PATH_SEPARATOR + obj_filename;
+    Slic3r::load_obj(fpath.c_str(), &mesh);
+    return mesh;
+}
+
+static bool _check_normals(const Slic3r::sla::Contour3D &mesh)
+{
+    for (auto & face : mesh.faces3)
+    {
+        
+    }
+    
+    return false;
+}
+
+TEST_CASE("Passing OpenVDB grid conversion produce similar geometry.", "[Hollowing]")
+{
+    Slic3r::TriangleMesh in_mesh = load_model("20mm_cube.obj");
+    Slic3r::sla::Contour3D imesh = Slic3r::sla::Contour3D{in_mesh};
+    auto ptr = Slic3r::meshToVolume(imesh, {});
+    
+    REQUIRE(ptr);
+    
+    std::cout << "Grid class = " << ptr->getGridClass() << std::endl;
+    
+    Slic3r::sla::Contour3D omesh = Slic3r::volumeToMesh(*ptr, -2.9, 1.0, true);
+    
+    std::cout << "Triangle count: " << omesh.faces3.size() << std::endl;
+    std::cout << "Quad count: " << omesh.faces4.size() << std::endl;
+    
+    REQUIRE(!omesh.empty());
+    
+    SECTION("Converting to Contour3D to TriangleMesh") {
+        Slic3r::TriangleMesh msh = Slic3r::sla::to_triangle_mesh(omesh);
+        
+        msh.require_shared_vertices();
+        msh.WriteOBJFile("out_tr.obj");
+        
+        REQUIRE(msh.volume() == Approx(in_mesh.volume()));        
+    }
+    
+//    omesh.faces4.clear();
+    std::fstream outfile{"out.obj", std::ios::out};
+    omesh.to_obj(outfile);
+    
+}
+
+//TEST_CASE("Negative 3D offset should produce smaller object.", "[Hollowing]")
+//{
+//    Slic3r::sla::Contour3D imesh = Slic3r::sla::Contour3D{load_model("20mm_cube.obj")};
+//    auto ptr = Slic3r::meshToVolume(imesh, {});
+    
+//    REQUIRE(ptr);
+    
+//    Slic3r::sla::Contour3D omesh = Slic3r::volumeToMesh(*ptr, -1., 0.0, true);
+    
+//    REQUIRE(!omesh.empty());
+    
+//    imesh.merge(omesh);
+//    std::fstream merged_outfile("merged_out.obj", std::ios::out);
+//    imesh.to_obj(merged_outfile);
+//}

tests/sla_print/CMakeLists.txt

@@ -1,5 +1,5 @@
 get_filename_component(_TEST_NAME ${CMAKE_CURRENT_LIST_DIR} NAME)
-add_executable(${_TEST_NAME}_tests ${_TEST_NAME}_tests.cpp)
+add_executable(${_TEST_NAME}_tests ${_TEST_NAME}_tests_main.cpp sla_print_tests.cpp)
 target_link_libraries(${_TEST_NAME}_tests test_common libslic3r)
 set_property(TARGET ${_TEST_NAME}_tests PROPERTY FOLDER "tests")
 

tests/sla_print/sla_print_tests.cpp

@@ -1,9 +1,9 @@
-#include <catch_main.hpp>
-
 #include <unordered_set>
 #include <unordered_map>
 #include <random>
 
+#include <catch2/catch.hpp>
+
 // Debug
 #include <fstream>
 
@@ -689,3 +689,17 @@ TEST_CASE("RasterizedPolygonAreaShouldMatch", "[SLARasterOutput]") {
     
     REQUIRE(diff <= predict_error(poly, pixdim));
 }
+
+TEST_CASE("Triangle mesh conversions should be correct", "[SLAConversions]")
+{
+    sla::Contour3D cntr;
+    
+    {
+        std::fstream infile{"extruder_idler_quads.obj", std::ios::in};
+        cntr.from_obj(infile);
+    }
+    
+    
+    
+    
+}

tests/sla_print/sla_print_tests_main.cpp

@@ -0,0 +1 @@
+#include <catch_main.hpp>