diff --git a/src/libslic3r/Geometry/Circle.cpp b/src/libslic3r/Geometry/Circle.cpp index 4d7c38ccc..679667195 100644 --- a/src/libslic3r/Geometry/Circle.cpp +++ b/src/libslic3r/Geometry/Circle.cpp @@ -108,7 +108,7 @@ Circled circle_taubin_newton(const Vec2ds& input, size_t cycles) return out; } -Circled circle_ransac(const Vec2ds& input, size_t iterations) +Circled circle_ransac(const Vec2ds& input, size_t iterations, double* min_error) { if (input.size() < 3) return Circled::make_invalid(); @@ -132,6 +132,8 @@ Circled circle_ransac(const Vec2ds& input, size_t iterations) circle_best = c; } } + if (min_error) + *min_error = err_min; return circle_best; } diff --git a/src/libslic3r/Geometry/Circle.hpp b/src/libslic3r/Geometry/Circle.hpp index 39973d916..653102e2a 100644 --- a/src/libslic3r/Geometry/Circle.hpp +++ b/src/libslic3r/Geometry/Circle.hpp @@ -102,7 +102,7 @@ inline Vec2d circle_center_taubin_newton(const Vec2ds& input, size_t cycles = 20 Circled circle_taubin_newton(const Vec2ds& input, size_t cycles = 20); // Find circle using RANSAC randomized algorithm. -Circled circle_ransac(const Vec2ds& input, size_t iterations = 20); +Circled circle_ransac(const Vec2ds& input, size_t iterations = 20, double* min_error = nullptr); // Randomized algorithm by Emo Welzl, working with squared radii for efficiency. The returned circle radius is inflated by epsilon. template diff --git a/src/libslic3r/Measure.cpp b/src/libslic3r/Measure.cpp index 6368f30f5..08c17df47 100644 --- a/src/libslic3r/Measure.cpp +++ b/src/libslic3r/Measure.cpp @@ -8,13 +8,15 @@ #include +#define DEBUG_EXTRACT_ALL_FEATURES_AT_ONCE 0 + namespace Slic3r { namespace Measure { constexpr double feature_hover_limit = 0.5; // how close to a feature the mouse must be to highlight it -static std::pair get_center_and_radius(const std::vector& points, const Transform3d& trafo) +static std::tuple get_center_and_radius(const std::vector& points, const Transform3d& trafo, const Transform3d& trafo_inv) { Vec2ds out; double z = 0.; @@ -24,18 +26,17 @@ static std::pair get_center_and_radius(const std::vector& out.emplace_back(pt_transformed.x(), pt_transformed.y()); } - auto circle = Geometry::circle_ransac(out, 20); // FIXME: iterations? + const int iter = points.size() < 10 ? 2 : + points.size() < 100 ? 4 : + 6; - return std::make_pair(trafo.inverse() * Vec3d(circle.center.x(), circle.center.y(), z), circle.radius); + double error = std::numeric_limits::max(); + auto circle = Geometry::circle_ransac(out, iter, &error); + + return std::make_tuple(trafo.inverse() * Vec3d(circle.center.x(), circle.center.y(), z), circle.radius, error); } -static bool circle_fit_is_ok(const std::vector& pts, const Vec3d& center, double radius) -{ - for (const Vec3d& pt : pts) - if (std::abs((pt - center).norm() - radius) > 0.05) - return false; - return true; -} + static std::array orthonormal_basis(const Vec3d& v) { @@ -64,17 +65,18 @@ public: std::vector surface_features; Vec3d normal; float area; + bool features_extracted = false; }; - std::vector get_all_features() const; - std::optional get_feature(size_t face_idx, const Vec3d& point) const; - std::vector> get_planes_triangle_indices() const; - const std::vector& get_plane_features(unsigned int plane_id) const; + std::optional get_feature(size_t face_idx, const Vec3d& point); + int get_num_of_planes() const; + const std::vector& get_plane_triangle_indices(int idx) const; + const std::vector& get_plane_features(unsigned int plane_id); const TriangleMesh& get_mesh() const; private: void update_planes(); - void extract_features(); + void extract_features(int plane_idx); std::vector m_planes; std::vector m_face_to_plane; @@ -90,7 +92,13 @@ MeasuringImpl::MeasuringImpl(const indexed_triangle_set& its) : m_mesh(its) { update_planes(); - extract_features(); + + // Extracting features will be done as needed. + // To extract all planes at once, run the following: +#if DEBUG_EXTRACT_ALL_FEATURES_AT_ONCE + for (int i=0; i angles; // placed in outer scope to prevent reallocations std::vector lengths; + for (const std::vector& border : plane.borders) { + if (border.size() <= 1) + continue; - for (int i=0; i<(int)m_planes.size(); ++i) { - PlaneData& plane = m_planes[i]; - plane.surface_features.clear(); - const Vec3d& normal = plane.normal; + bool done = false; - Eigen::Quaterniond q; - q.setFromTwoVectors(plane.normal, Vec3d::UnitZ()); - Transform3d trafo = Transform3d::Identity(); - trafo.rotate(q); + if (border.size() > 4) { + const auto& [center, radius, err] = get_center_and_radius(border, trafo, trafo_inv); - for (const std::vector& border : plane.borders) { - if (border.size() <= 1) - continue; - - bool done = false; - - if (const auto& [center, radius] = get_center_and_radius(border, trafo); - (border.size()>4) && circle_fit_is_ok(border, center, radius)) { + if (err < 0.05) { // The whole border is one circle. Just add it into the list of features // and we are done. @@ -298,200 +309,190 @@ void MeasuringImpl::extract_features() done = true; } } + } - if (! done) { - // In this case, the border is not a circle and may contain circular - // segments. Try to find them and then add all remaining edges as edges. + if (! done) { + // In this case, the border is not a circle and may contain circular + // segments. Try to find them and then add all remaining edges as edges. - auto are_angles_same = [](double a, double b) { return Slic3r::is_approx(a,b,0.01); }; - auto are_lengths_same = [](double a, double b) { return Slic3r::is_approx(a,b,0.01); }; + auto are_angles_same = [](double a, double b) { return Slic3r::is_approx(a,b,0.01); }; + auto are_lengths_same = [](double a, double b) { return Slic3r::is_approx(a,b,0.01); }; - // Given an idx into border, return the index that is idx+offset position, - // while taking into account the need for wrap-around and the fact that - // the first and last point are the same. - auto offset_to_index = [border_size = int(border.size())](int idx, int offset) -> int { - assert(std::abs(offset) < border_size); - int out = idx+offset; - if (out >= border_size) - out = out - border_size; - else if (out < 0) - out = border_size + out; + // Given an idx into border, return the index that is idx+offset position, + // while taking into account the need for wrap-around and the fact that + // the first and last point are the same. + auto offset_to_index = [border_size = int(border.size())](int idx, int offset) -> int { + assert(std::abs(offset) < border_size); + int out = idx+offset; + if (out >= border_size) + out = out - border_size; + else if (out < 0) + out = border_size + out; - return out; - }; + return out; + }; - // First calculate angles at all the vertices. - angles.clear(); - lengths.clear(); - int first_different_angle_idx = 0; - for (int i=0; i M_PI) - angle = 2*M_PI - angle; + // First calculate angles at all the vertices. + angles.clear(); + lengths.clear(); + int first_different_angle_idx = 0; + for (int i=0; i M_PI) + angle = 2*M_PI - angle; - angles.push_back(angle); - lengths.push_back(v2.norm()); - if (first_different_angle_idx == 0 && angles.size() > 1) { - if (! are_angles_same(angles.back(), angles[angles.size()-2])) - first_different_angle_idx = angles.size()-1; - } + angles.push_back(angle); + lengths.push_back(v2.norm()); + if (first_different_angle_idx == 0 && angles.size() > 1) { + if (! are_angles_same(angles.back(), angles[angles.size()-2])) + first_different_angle_idx = angles.size()-1; } - assert(border.size() == angles.size()); - assert(border.size() == lengths.size()); + } + assert(border.size() == angles.size()); + assert(border.size() == lengths.size()); - // First go around the border and pick what might be circular segments. - // Save pair of indices to where such potential segments start and end. - // Also remember the length of these segments. - int start_idx = -1; - bool circle = false; - bool first_iter = true; - std::vector circles; - std::vector edges; - std::vector> circles_idxs; - //std::vector circles_lengths; - std::vector single_circle; // could be in loop-scope, but reallocations - double single_circle_length = 0.; - int first_pt_idx = offset_to_index(first_different_angle_idx, 1); - int i = first_pt_idx; - while (i != first_pt_idx || first_iter) { - if (are_angles_same(angles[i], angles[offset_to_index(i,-1)]) - && i != offset_to_index(first_pt_idx, -1) // not the last point - && i != start_idx ) { - // circle - if (! circle) { - circle = true; - single_circle.clear(); - single_circle_length = 0.; - start_idx = offset_to_index(i, -2); - single_circle = { border[start_idx], border[offset_to_index(start_idx,1)] }; - single_circle_length += lengths[offset_to_index(i, -1)]; - } + // First go around the border and pick what might be circular segments. + // Save pair of indices to where such potential segments start and end. + // Also remember the length of these segments. + int start_idx = -1; + bool circle = false; + bool first_iter = true; + std::vector circles; + std::vector edges; + std::vector> circles_idxs; + //std::vector circles_lengths; + std::vector single_circle; // could be in loop-scope, but reallocations + double single_circle_length = 0.; + int first_pt_idx = offset_to_index(first_different_angle_idx, 1); + int i = first_pt_idx; + while (i != first_pt_idx || first_iter) { + if (are_angles_same(angles[i], angles[offset_to_index(i,-1)]) + && i != offset_to_index(first_pt_idx, -1) // not the last point + && i != start_idx ) { + // circle + if (! circle) { + circle = true; + single_circle.clear(); + single_circle_length = 0.; + start_idx = offset_to_index(i, -2); + single_circle = { border[start_idx], border[offset_to_index(start_idx,1)] }; + single_circle_length += lengths[offset_to_index(i, -1)]; + } + single_circle.emplace_back(border[i]); + single_circle_length += lengths[i]; + } else { + if (circle && single_circle.size() >= 5) { // Less than 5 vertices? Not a circle. single_circle.emplace_back(border[i]); single_circle_length += lengths[i]; - } else { - if (circle && single_circle.size() >= 5) { // Less than 5 vertices? Not a circle. - single_circle.emplace_back(border[i]); - single_circle_length += lengths[i]; - bool accept_circle = true; - { - // Check that lengths of internal (!!!) edges match. - int j = offset_to_index(start_idx, 3); - while (j != i) { - if (! are_lengths_same(lengths[offset_to_index(j,-1)], lengths[j])) { - accept_circle = false; - break; - } - j = offset_to_index(j, 1); + bool accept_circle = true; + { + // Check that lengths of internal (!!!) edges match. + int j = offset_to_index(start_idx, 3); + while (j != i) { + if (! are_lengths_same(lengths[offset_to_index(j,-1)], lengths[j])) { + accept_circle = false; + break; } + j = offset_to_index(j, 1); } + } + + if (accept_circle) { + const auto& [center, radius, err] = get_center_and_radius(single_circle, trafo, trafo_inv); + + // Check that the fit went well. The tolerance is high, only to + // reject complete failures. + accept_circle &= err < 0.05; + + // If the segment subtends less than 90 degrees, throw it away. + accept_circle &= single_circle_length / radius > 0.9*M_PI/2.; if (accept_circle) { - const auto& [center, radius] = get_center_and_radius(single_circle, trafo); - - // Check that the fit went well. The tolerance is high, only to - // reject complete failures. - accept_circle &= circle_fit_is_ok(single_circle, center, radius); - - // If the segment subtends less than 90 degrees, throw it away. - accept_circle &= single_circle_length / radius > 0.9*M_PI/2.; - - if (accept_circle) { - // Add the circle and remember indices into borders. - circles_idxs.emplace_back(start_idx, i); - circles.emplace_back(SurfaceFeature(SurfaceFeatureType::Circle, center, plane.normal, std::nullopt, radius)); - } + // Add the circle and remember indices into borders. + circles_idxs.emplace_back(start_idx, i); + circles.emplace_back(SurfaceFeature(SurfaceFeatureType::Circle, center, plane.normal, std::nullopt, radius)); } } - circle = false; } - // Take care of the wrap around. - first_iter = false; + circle = false; + } + // Take care of the wrap around. + first_iter = false; + i = offset_to_index(i, 1); + } + + // We have the circles. Now go around again and pick edges, while jumping over circles. + if (circles_idxs.empty()) { + // Just add all edges. + for (int i=1; i 1 || circles_idxs.front().first != circles_idxs.front().second) { + // There is at least one circular segment. Start at its end and add edges until the start of the next one. + int i = circles_idxs.front().second; + int circle_idx = 1; + while (true) { i = offset_to_index(i, 1); - } - - // We have the circles. Now go around again and pick edges, while jumping over circles. - if (circles_idxs.empty()) { - // Just add all edges. - for (int i=1; i 1 || circles_idxs.front().first != circles_idxs.front().second) { - // There is at least one circular segment. Start at its end and add edges until the start of the next one. - int i = circles_idxs.front().second; - int circle_idx = 1; - while (true) { - i = offset_to_index(i, 1); - edges.emplace_back(SurfaceFeature(SurfaceFeatureType::Edge, border[offset_to_index(i,-1)], border[i])); - if (circle_idx < int(circles_idxs.size()) && i == circles_idxs[circle_idx].first) { - i = circles_idxs[circle_idx].second; - ++circle_idx; - } - if (i == circles_idxs.front().first) - break; + edges.emplace_back(SurfaceFeature(SurfaceFeatureType::Edge, border[offset_to_index(i,-1)], border[i])); + if (circle_idx < int(circles_idxs.size()) && i == circles_idxs[circle_idx].first) { + i = circles_idxs[circle_idx].second; + ++circle_idx; } + if (i == circles_idxs.front().first) + break; } + } - // Merge adjacent edges where needed. - assert(std::all_of(edges.begin(), edges.end(), - [](const SurfaceFeature& f) { return f.get_type() == SurfaceFeatureType::Edge; })); - for (int i=edges.size()-1; i>=0; --i) { - const auto& [first_start, first_end] = edges[i==0 ? edges.size()-1 : i-1].get_edge(); - const auto& [second_start, second_end] = edges[i].get_edge(); + // Merge adjacent edges where needed. + assert(std::all_of(edges.begin(), edges.end(), + [](const SurfaceFeature& f) { return f.get_type() == SurfaceFeatureType::Edge; })); + for (int i=edges.size()-1; i>=0; --i) { + const auto& [first_start, first_end] = edges[i==0 ? edges.size()-1 : i-1].get_edge(); + const auto& [second_start, second_end] = edges[i].get_edge(); - if (Slic3r::is_approx(first_end, second_start) - && Slic3r::is_approx((first_end-first_start).normalized().dot((second_end-second_start).normalized()), 1.)) { - // The edges have the same direction and share a point. Merge them. - edges[i==0 ? edges.size()-1 : i-1] = SurfaceFeature(SurfaceFeatureType::Edge, first_start, second_end); - edges.erase(edges.begin() + i); - } + if (Slic3r::is_approx(first_end, second_start) + && Slic3r::is_approx((first_end-first_start).normalized().dot((second_end-second_start).normalized()), 1.)) { + // The edges have the same direction and share a point. Merge them. + edges[i==0 ? edges.size()-1 : i-1] = SurfaceFeature(SurfaceFeatureType::Edge, first_start, second_end); + edges.erase(edges.begin() + i); } - - // Now move the circles and edges into the feature list for the plane. - assert(std::all_of(circles.begin(), circles.end(), [](const SurfaceFeature& f) { - return f.get_type() == SurfaceFeatureType::Circle; - })); - assert(std::all_of(edges.begin(), edges.end(), [](const SurfaceFeature& f) { - return f.get_type() == SurfaceFeatureType::Edge; - })); - plane.surface_features.insert(plane.surface_features.end(), std::make_move_iterator(circles.begin()), - std::make_move_iterator(circles.end())); - plane.surface_features.insert(plane.surface_features.end(), std::make_move_iterator(edges.begin()), - std::make_move_iterator(edges.end())); } - } - // The last surface feature is the plane itself. - Vec3d cog = Vec3d::Zero(); - size_t counter = 0; - for (const std::vector& b : plane.borders) { - for (size_t i = 1; i < b.size(); ++i) { - cog += b[i]; - ++counter; - } + // Now move the circles and edges into the feature list for the plane. + assert(std::all_of(circles.begin(), circles.end(), [](const SurfaceFeature& f) { + return f.get_type() == SurfaceFeatureType::Circle; + })); + assert(std::all_of(edges.begin(), edges.end(), [](const SurfaceFeature& f) { + return f.get_type() == SurfaceFeatureType::Edge; + })); + plane.surface_features.insert(plane.surface_features.end(), std::make_move_iterator(circles.begin()), + std::make_move_iterator(circles.end())); + plane.surface_features.insert(plane.surface_features.end(), std::make_move_iterator(edges.begin()), + std::make_move_iterator(edges.end())); } - cog /= double(counter); - plane.surface_features.emplace_back(SurfaceFeature(SurfaceFeatureType::Plane, - plane.normal, cog, std::optional(), i + 0.0001)); - - plane.borders.clear(); - plane.borders.shrink_to_fit(); } -} + // The last surface feature is the plane itself. + Vec3d cog = Vec3d::Zero(); + size_t counter = 0; + for (const std::vector& b : plane.borders) { + for (size_t i = 1; i < b.size(); ++i) { + cog += b[i]; + ++counter; + } + } + cog /= double(counter); + plane.surface_features.emplace_back(SurfaceFeature(SurfaceFeatureType::Plane, + plane.normal, cog, std::optional(), plane_idx + 0.0001)); + plane.borders.clear(); + plane.borders.shrink_to_fit(); -std::vector MeasuringImpl::get_all_features() const -{ - std::vector features; - //PlaneData& plane = m_planes[0]; - for (const PlaneData& plane : m_planes) - for (const SurfaceFeature& feature : plane.surface_features) - features.emplace_back(feature); - return features; + plane.features_extracted = true; } @@ -499,12 +500,17 @@ std::vector MeasuringImpl::get_all_features() const -std::optional MeasuringImpl::get_feature(size_t face_idx, const Vec3d& point) const + + +std::optional MeasuringImpl::get_feature(size_t face_idx, const Vec3d& point) { if (face_idx >= m_face_to_plane.size()) return std::optional(); const PlaneData& plane = m_planes[m_face_to_plane[face_idx]]; + + if (! plane.features_extracted) + extract_features(m_face_to_plane[face_idx]); size_t closest_feature_idx = size_t(-1); double min_dist = std::numeric_limits::max(); @@ -554,17 +560,24 @@ std::optional MeasuringImpl::get_feature(size_t face_idx, const -std::vector> MeasuringImpl::get_planes_triangle_indices() const +int MeasuringImpl::get_num_of_planes() const { - std::vector> out; - for (const PlaneData& plane : m_planes) - out.emplace_back(plane.facets); - return out; + return (m_planes.size()); } -const std::vector& MeasuringImpl::get_plane_features(unsigned int plane_id) const + + +const std::vector& MeasuringImpl::get_plane_triangle_indices(int idx) const +{ + assert(idx >= 0 && idx < int(m_planes.size())); + return m_planes[idx].facets; +} + +const std::vector& MeasuringImpl::get_plane_features(unsigned int plane_id) { assert(plane_id < m_planes.size()); + if (! m_planes[plane_id].features_extracted) + extract_features(plane_id); return m_planes[plane_id].surface_features; } @@ -590,11 +603,6 @@ Measuring::Measuring(const indexed_triangle_set& its) Measuring::~Measuring() {} -std::vector Measuring::get_all_features() const -{ - return priv->get_all_features(); -} - std::optional Measuring::get_feature(size_t face_idx, const Vec3d& point) const { @@ -602,10 +610,15 @@ std::optional Measuring::get_feature(size_t face_idx, const Vec3 } - -std::vector> Measuring::get_planes_triangle_indices() const +int Measuring::get_num_of_planes() const { - return priv->get_planes_triangle_indices(); + return priv->get_num_of_planes(); +} + + +const std::vector& Measuring::get_plane_triangle_indices(int idx) const +{ + return priv->get_plane_triangle_indices(idx); } const std::vector& Measuring::get_plane_features(unsigned int plane_id) const diff --git a/src/libslic3r/Measure.hpp b/src/libslic3r/Measure.hpp index 3a1a4b89d..a273135ca 100644 --- a/src/libslic3r/Measure.hpp +++ b/src/libslic3r/Measure.hpp @@ -93,20 +93,18 @@ public: // Construct the measurement object on a given its. explicit Measuring(const indexed_triangle_set& its); ~Measuring(); - - // Return a reference to a list of all features identified on the its. - // Use only for debugging. Expensive, do not call often. - std::vector get_all_features() const; + // Given a face_idx where the mouse cursor points, return a feature that // should be highlighted (if any). std::optional get_feature(size_t face_idx, const Vec3d& point) const; - // Returns a list of triangle indices for each identified plane. Each - // Plane object contains an index into this vector. Expensive, do not - // call too often. - std::vector> get_planes_triangle_indices() const; - + // Return total number of planes. + int get_num_of_planes() const; + + // Returns a list of triangle indices for given plane. + const std::vector& get_plane_triangle_indices(int idx) const; + // Returns the surface features of the plane with the given index const std::vector& get_plane_features(unsigned int plane_id) const; diff --git a/src/slic3r/GUI/Gizmos/GLGizmoMeasure.cpp b/src/slic3r/GUI/Gizmos/GLGizmoMeasure.cpp index 8c4b46826..6994bfa32 100644 --- a/src/slic3r/GUI/Gizmos/GLGizmoMeasure.cpp +++ b/src/slic3r/GUI/Gizmos/GLGizmoMeasure.cpp @@ -93,11 +93,8 @@ static std::string center_on_feature_type_as_string(Measure::SurfaceFeatureType return ret; } -static GLModel::Geometry init_plane_data(const indexed_triangle_set& its, const std::vector>& planes_triangles, int idx) +static GLModel::Geometry init_plane_data(const indexed_triangle_set& its, const std::vector& triangle_indices) { - assert(0 <= idx && idx < (int)planes_triangles.size()); - const std::vector& triangle_indices = planes_triangles[idx]; - GLModel::Geometry init_data; init_data.format = { GUI::GLModel::Geometry::EPrimitiveType::Triangles, GLModel::Geometry::EVertexLayout::P3N3 }; unsigned int i = 0; @@ -653,11 +650,10 @@ void GLGizmoMeasure::on_render() if (m_last_plane_idx != idx) { m_last_plane_idx = idx; const indexed_triangle_set& its = m_measuring->get_mesh().its; - const std::vector> planes_triangles = m_measuring->get_planes_triangle_indices(); - GLModel::Geometry init_data = init_plane_data(its, planes_triangles, idx); + const std::vector& plane_triangles = m_measuring->get_plane_triangle_indices(idx); + GLModel::Geometry init_data = init_plane_data(its, plane_triangles); m_plane.reset(); m_plane.mesh_raycaster = std::make_unique(std::make_shared(init_data.get_as_indexed_triangle_set())); - m_plane.model.init_from(std::move(init_data)); } m_raycasters.insert({ PLANE_ID, m_parent.add_raycaster_for_picking(SceneRaycaster::EType::Gizmo, PLANE_ID, *m_plane.mesh_raycaster) }); @@ -1041,10 +1037,9 @@ void GLGizmoMeasure::update_if_needed() { auto update_plane_models_cache = [this](const indexed_triangle_set& its) { m_plane_models_cache.clear(); - const std::vector> planes_triangles = m_measuring->get_planes_triangle_indices(); - for (int idx = 0; idx < (int)planes_triangles.size(); ++idx) { + for (int idx = 0; idx < m_measuring->get_num_of_planes(); ++idx) { m_plane_models_cache.emplace_back(GLModel()); - GLModel::Geometry init_data = init_plane_data(its, planes_triangles, idx); + GLModel::Geometry init_data = init_plane_data(its, m_measuring->get_plane_triangle_indices(idx)); m_plane_models_cache.back().init_from(std::move(init_data)); } };