EigenMesh3D raycaster should now be able to pick a correct intersection on the object or inside a hole

src/libslic3r/SLA/Common.cpp

@@ -183,18 +183,6 @@ void BoxIndex::foreach(std::function<void (const BoxIndexEl &)> fn)
 }
 
 
-namespace {
-// Iterates over hits and holes and returns the true hit, possibly
-// on the inside of a hole. Free function so it can return igl::Hit
-// without including igl in a header.
-igl::Hit filter_hits(const std::vector<EigenMesh3D::hit_result>& hits,
-                     const std::vector<DrainHole>& holes)
-{
-    return igl::Hit();
-}
-
-} // namespace
-
 /* ****************************************************************************
  * EigenMesh3D implementation
  * ****************************************************************************/
@@ -280,22 +268,27 @@ EigenMesh3D::query_ray_hit(const Vec3d &s,
                            const std::vector<DrainHole>* holes
                           ) const
 {
+    assert(is_approx(dir.norm(), 1.));
     igl::Hit hit;
     hit.t = std::numeric_limits<float>::infinity();
 
-    if (! holes)
+    if (! holes) {
         m_aabb->intersect_ray(m_V, m_F, s, dir, hit);
-    else
-        hit = filter_hits(query_ray_hits(s, dir), *holes);
-    
         hit_result ret(*this);
         ret.m_t = double(hit.t);
         ret.m_dir = dir;
         ret.m_source = s;
-    if(!std::isinf(hit.t) && !std::isnan(hit.t)) ret.m_face_id = hit.id;
+        if(!std::isinf(hit.t) && !std::isnan(hit.t))
+            ret.m_normal = this->normal_by_face_id(hit.id);
 
         return ret;
     }
+    else {
+        // If there are holes, the hit_results will be made by
+        // query_ray_hits (object) and filter_hits (holes):
+        return filter_hits(query_ray_hits(s, dir), *holes);
+    }
+}
 
 std::vector<EigenMesh3D::hit_result>
 EigenMesh3D::query_ray_hits(const Vec3d &s, const Vec3d &dir) const
@@ -316,12 +309,91 @@ EigenMesh3D::query_ray_hits(const Vec3d &s, const Vec3d &dir) const
         outs.back().m_dir = dir;
         outs.back().m_source = s;
         if(!std::isinf(hit.t) && !std::isnan(hit.t))
-            outs.back().m_face_id = hit.id;
+            outs.back().m_normal = this->normal_by_face_id(hit.id);
     }
 
     return outs;
 }
 
+EigenMesh3D::hit_result EigenMesh3D::filter_hits(
+                     const std::vector<EigenMesh3D::hit_result>& object_hits,
+                     const std::vector<DrainHole>& holes) const
+{
+    hit_result out(*this);
+    out.m_t = std::nan("");
+
+    if (! holes.empty() && ! object_hits.empty()) {
+        Vec3d s = object_hits.front().source();
+        Vec3d dir = object_hits.front().direction();
+
+        struct HoleHit {
+            HoleHit(float t_p, const Vec3d& normal_p, bool entry_p) :
+                t(t_p), normal(normal_p), entry(entry_p) {}
+            float t;
+            Vec3d normal;
+            bool entry;
+        };
+        std::vector<HoleHit> hole_isects;
+
+        // Collect hits on all holes, preserve information about entry/exit
+        for (const sla::DrainHole& hole : holes) {
+            std::array<std::pair<float, Vec3d>, 2> isects;
+            if (hole.get_intersections(s.cast<float>(),
+                                       dir.cast<float>(), isects)) {
+                hole_isects.emplace_back(isects[0].first, isects[0].second, true);
+                hole_isects.emplace_back(isects[1].first, isects[1].second, false);
+            }
+        }
+        // Holes can intersect each other, sort the hits by t
+        std::sort(hole_isects.begin(), hole_isects.end(),
+                  [](const HoleHit& a, const HoleHit& b) { return a.t < b.t; });
+
+        // Now inspect the intersections with object and holes, keep track how
+        // deep are we nested in mesh/holes and pick the correct intersection
+        int hole_nested = 0;
+        int object_nested = 0;
+
+        bool is_hole = false;
+        bool is_entry = false;
+        const HoleHit* next_hole_hit = &hole_isects.front();
+        const hit_result* next_mesh_hit = &object_hits.front();
+
+        while (next_hole_hit || next_mesh_hit) {
+            if (next_hole_hit && next_mesh_hit) // still have hole and obj hits
+                is_hole = (next_hole_hit->t < next_mesh_hit->m_t);
+            else
+                is_hole = next_hole_hit; // one or the other ran out
+
+            // Is this entry or exit hit?
+            is_entry = is_hole ? next_hole_hit->entry : ! next_mesh_hit->is_inside();
+
+            if (! is_hole && is_entry && hole_nested == 0) {
+                // This mesh point is the one we seek
+                return *next_mesh_hit;
+            }
+            if (is_hole && ! is_entry && object_nested != 0) {
+                // This holehit is the one we seek
+                out.m_t = next_hole_hit->t;
+                out.m_normal = next_hole_hit->normal;
+                out.m_source = s;
+                out.m_dir = dir;
+                return out;
+            }
+
+            hole_nested += (is_hole ? (is_entry ? 1 : -1) : 0);
+            object_nested += (! is_hole ? (is_entry ? 1 : -1) : 0);
+
+            // Advance the pointer
+            if (is_hole && next_hole_hit++ == &hole_isects.back())
+                next_hole_hit = nullptr;
+            if (! is_hole && next_mesh_hit++ == &object_hits.back())
+                next_mesh_hit = nullptr;
+        }
+    }
+
+    return out;
+}
+
 #ifdef SLIC3R_SLA_NEEDS_WINDTREE
 EigenMesh3D::si_result EigenMesh3D::signed_distance(const Vec3d &p) const {
     double sign = 0; double sqdst = 0; int i = 0;  Vec3d c;

src/libslic3r/SLA/EigenMesh3D.hpp

@@ -42,10 +42,10 @@ public:
     // Result of a raycast
     class hit_result {
         double m_t = std::nan("");
-        int m_face_id = -1;
         const EigenMesh3D *m_mesh = nullptr;
         Vec3d m_dir;
         Vec3d m_source;
+        Vec3d m_normal = Vec3d::Zero();
         friend class EigenMesh3D;
         
         // A valid object of this class can only be obtained from
@@ -60,26 +60,23 @@ public:
         
         inline double distance() const { return m_t; }
         inline const Vec3d& direction() const { return m_dir; }
+        inline const Vec3d& source() const { return m_source; }
         inline Vec3d position() const { return m_source + m_dir * m_t; }
-        inline int face() const { return m_face_id; }
         inline bool is_valid() const { return m_mesh != nullptr; }
+        inline bool is_hit() const { return m_normal != Vec3d::Zero(); }
         
         // Hit_result can decay into a double as the hit distance.
         inline operator double() const { return distance(); }
 
-        inline Vec3d normal() const {
+        inline const Vec3d& normal() const {
-            if(m_face_id < 0 || !is_valid()) return {};
+            if(!is_valid())
-            auto trindex    = m_mesh->m_F.row(m_face_id);
+                throw std::runtime_error("EigenMesh3D::hit_result::normal() "
-            const Vec3d& p1 = m_mesh->V().row(trindex(0));
+                                         "called on invalid object.");
-            const Vec3d& p2 = m_mesh->V().row(trindex(1));
+            return m_normal;
-            const Vec3d& p3 = m_mesh->V().row(trindex(2));
-            Eigen::Vector3d U = p2 - p1;
-            Eigen::Vector3d V = p3 - p1;
-            return U.cross(V).normalized();
         }
         
-        inline bool is_inside() {
+        inline bool is_inside() const {
-            return m_face_id >= 0 && normal().dot(m_dir) > 0;
+            return normal().dot(m_dir) > 0;
         }
     };
     
@@ -91,6 +88,11 @@ public:
     // Casts a ray on the mesh and returns all hits
     std::vector<hit_result> query_ray_hits(const Vec3d &s, const Vec3d &dir) const;
 
+    // Iterates over hits and holes and returns the true hit, possibly
+    // on the inside of a hole.
+    hit_result filter_hits(const std::vector<EigenMesh3D::hit_result>& obj_hits,
+                           const std::vector<DrainHole>& holes) const;
+
     class si_result {
         double m_value;
         int m_fidx;
@@ -123,6 +125,16 @@ public:
         Vec3d c;
         return squared_distance(p, i, c);
     }
+
+    Vec3d normal_by_face_id(int face_id) const {
+        auto trindex    = F().row(face_id);
+        const Vec3d& p1 = V().row(trindex(0));
+        const Vec3d& p2 = V().row(trindex(1));
+        const Vec3d& p3 = V().row(trindex(2));
+        Eigen::Vector3d U = p2 - p1;
+        Eigen::Vector3d V = p3 - p1;
+        return U.cross(V).normalized();
+    }
 };
 
 // Calculate the normals for the selected points (from 'points' set) on the

src/libslic3r/SLA/Hollowing.cpp

@@ -141,22 +141,26 @@ bool DrainHole::is_inside(const Vec3f& pt) const
 }
 
 
-// Given a line s+dir*t, return parameter t of intersections  with the hole.
+// Given a line s+dir*t, find parameter t of intersections with the hole
-// If there is no intersection, returns nan.
+// and the normal (points inside the hole). Outputs through out reference,
-std::pair<float, float> DrainHole::get_intersections(const Vec3f& s,
+// returns true if two intersections were found.
-                                                     const Vec3f& dir) const
+bool DrainHole::get_intersections(const Vec3f& s, const Vec3f& dir,
+                                  std::array<std::pair<float, Vec3d>, 2>& out)
+                                  const
 {
     assert(is_approx(normal.norm(), 1.f));
     const Eigen::ParametrizedLine<float, 3> ray(s, dir.normalized());
 
-    std::pair<float, float> out(std::nan(""), std::nan(""));
+    for (size_t i=0; i<2; ++i)
+        out[i] = std::make_pair(std::nan(""), Vec3d::Zero());
+
     const float sqr_radius = pow(radius, 2.f);
 
     // first check a bounding sphere of the hole:
     Vec3f center = pos+normal*height/2.f;
     float sqr_dist_limit = pow(height/2.f, 2.f) + sqr_radius ;
     if (ray.squaredDistance(center) > sqr_dist_limit)
-        return out;
+        return false;
 
     // The line intersects the bounding sphere, look for intersections with
     // bases of the cylinder.
@@ -170,7 +174,8 @@ std::pair<float, float> DrainHole::get_intersections(const Vec3f& s,
             Vec3f intersection = ray.intersectionPoint(base);
             // Only accept the point if it is inside the cylinder base.
             if ((cylinder_center-intersection).squaredNorm() < sqr_radius) {
-                (found ? out.second : out.first) = ray.intersectionParameter(base);
+                out[found].first = ray.intersectionParameter(base);
+                out[found].second = (i==0 ? 1. : -1.) * normal.cast<double>();
                 ++found;
             }
         }
@@ -200,11 +205,13 @@ std::pair<float, float> DrainHole::get_intersections(const Vec3f& s,
         for (int i=-1; i<=1 && found !=2; i+=2) {
             Vec3f isect = closest + i*dist * projected_ray.direction();
             float par = (isect-proj_origin).norm() / par_scale;
+            Vec3d hit_normal = (pos-isect).normalized().cast<double>();
             isect = ray.pointAt(par);
             // check that the intersection is between the base planes:
             float vert_dist = base.signedDistance(isect);
             if (vert_dist > 0.f && vert_dist < height) {
-                (found ? out.second : out.first) = par;
+                out[found].first = par;
+                out[found].second = hit_normal;
                 ++found;
             }
         }
@@ -212,14 +219,14 @@ std::pair<float, float> DrainHole::get_intersections(const Vec3f& s,
 
     // If only one intersection was found, it is some corner case,
     // no intersection will be returned:
-    if (found != 0)
+    if (found != 2)
-        return std::pair<float, float>(std::nan(""), std::nan(""));
+        return false;
 
     // Sort the intersections:
-    if (out.first > out.second)
+    if (out[0].first > out[1].first)
-        std::swap(out.first, out.second);
+        std::swap(out[0], out[1]);
 
-    return out;
+    return true;
 }
 
 }} // namespace Slic3r::sla

src/libslic3r/SLA/Hollowing.hpp

@@ -39,9 +39,8 @@ struct DrainHole
 
     bool is_inside(const Vec3f& pt) const;
 
-    std::pair<float, float> get_intersections(const Vec3f& s,
+    bool get_intersections(const Vec3f& s, const Vec3f& dir,
-                                              const Vec3f& dir
+                           std::array<std::pair<float, Vec3d>, 2>& out) const;
-                                             ) const;
     
     template<class Archive> inline void serialize(Archive &ar)
     {

src/libslic3r/SLA/SupportPointGenerator.cpp

@@ -83,8 +83,8 @@ void SupportPointGenerator::project_onto_mesh(std::vector<sla::SupportPoint>& po
                 sla::EigenMesh3D::hit_result hit_up   = m_emesh.query_ray_hit(p.cast<double>(), Vec3d(0., 0., 1.));
                 sla::EigenMesh3D::hit_result hit_down = m_emesh.query_ray_hit(p.cast<double>(), Vec3d(0., 0., -1.));
 
-                bool up   = hit_up.face() != -1;
+                bool up   = hit_up.is_hit();
-                bool down = hit_down.face() != -1;
+                bool down = hit_down.is_hit();
 
                 if (!up && !down)
                     continue;