From 173ae7b0914b7a1ace81e262553e6b0f40caa13c Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?L=C3=A9o=20Valque?= Date: Fri, 7 Mar 2025 16:27:32 +0100 Subject: [PATCH 01/41] Creation of float_snap_rounding_2.h and its example .cpp --- .../Snap_rounding_2/float_snap_rounding.cpp | 71 +++++++ .../include/CGAL/Float_snap_rounding_2.h | 182 ++++++++++++++++++ .../include/CGAL/Snap_rounding_traits_2.h | 5 + 3 files changed, 258 insertions(+) create mode 100644 Snap_rounding_2/examples/Snap_rounding_2/float_snap_rounding.cpp create mode 100644 Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h diff --git a/Snap_rounding_2/examples/Snap_rounding_2/float_snap_rounding.cpp b/Snap_rounding_2/examples/Snap_rounding_2/float_snap_rounding.cpp new file mode 100644 index 00000000000..2cf4203bdb8 --- /dev/null +++ b/Snap_rounding_2/examples/Snap_rounding_2/float_snap_rounding.cpp @@ -0,0 +1,71 @@ +#include +#include +#include + +typedef CGAL::Exact_predicates_exact_constructions_kernel Kernel; +typedef Kernel::Segment_2 Segment_2; +typedef Kernel::Point_2 Point_2; +typedef Kernel::FT FT; + +int main() +{ + std::vector pts; + std::vector< std::vector > segs; + + Kernel::FT eps(std::pow(2, -60)); + Kernel::FT one(1.); + Kernel::FT two(2.); + + pts.emplace_back(one-eps, one); + pts.emplace_back(-one+eps, -one+2*eps); + pts.emplace_back(eps/2, eps/2); + pts.emplace_back(one, -one); + pts.emplace_back(0, two-eps/2); + pts.emplace_back(two, 0); + pts.emplace_back(-two+eps, -FT(4.)); + pts.emplace_back(two, two); + pts.emplace_back(-two, two); + + segs.push_back(std::vector()); + segs[0].push_back(0); + segs[0].push_back(1); + segs.push_back(std::vector()); + segs[1].push_back(2); + segs[1].push_back(3); + segs.push_back(std::vector()); + segs[2].push_back(4); + segs[2].push_back(5); + segs.push_back(std::vector()); + segs[3].push_back(4); + segs[3].push_back(6); + segs.push_back(std::vector()); + segs[4].push_back(7); + segs[4].push_back(8); + + std::cout << "Input" << std::endl; + std::cout << "Points" << std::endl; + for(auto &p: pts) + std::cout << p << std::endl; + std::cout << "Segments:" << std::endl; + for(auto &seg: segs){ + for(auto pi : seg) + std::cout << pi << " "; + std::cout << std::endl; + } + std::cout << "\n\n" << std::endl; + + CGAL::float_snap_rounding_2(pts, segs); + + std::cout << "Output" << std::endl; + std::cout << "Points" << std::endl; + for(auto &p: pts) + std::cout << p << std::endl; + std::cout << "Segments:" << std::endl; + for(auto &seg: segs){ + for(auto pi : seg) + std::cout << pi << " "; + std::cout << std::endl; + } + + return(0); +} diff --git a/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h b/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h new file mode 100644 index 00000000000..efb25b2efca --- /dev/null +++ b/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h @@ -0,0 +1,182 @@ +// Copyright (c) 2025 Geometry Factory. +// All rights reserved. +// +// This file is part of CGAL (www.cgal.org). +// +// $URL$ +// $Id$ +// SPDX-License-Identifier: GPL-3.0-or-later OR LicenseRef-Commercial +// +// +// author(s) : Léo Valque + +#ifndef CGAL_FLOAT_SNAP_ROUNDING_2_H +#define CGAL_FLOAT_SNAP_ROUNDING_2_H + +#include + + +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include + +namespace CGAL { +namespace Box_intersection_d { + +template +class Box_with_index_d: public Box_d< NT_, N, ID_NONE> { +protected: + size_t m_index; +public: + typedef Box_d< NT_, N, ID_NONE> Base; + typedef NT_ NT; + typedef size_t ID; + + Box_with_index_d() {} + Box_with_index_d( ID i) : m_index(i) {} + Box_with_index_d( bool complete, ID i): Base(complete), m_index(i) {} + Box_with_index_d(NT l[N], NT h[N], ID i) : Base( l, h), m_index(i) {} + Box_with_index_d( const Bbox_2& b, ID i) : Base( b), m_index(i) {} + Box_with_index_d( const Bbox_3& b, ID i) : Base( b), m_index(i) {} + ID index() const { return m_index; } + ID id() const { return m_index; } +}; + +// Generic template signature of boxes, specialized for ID_FROM_HANDLE policy +#if 0 +template +class Box_with_index_d : public Box_d< NT_, N, IdPolicy> { +protected: + size_t m_index; +public: + typedef Box_d< NT_, N, IdPolicy> Base; + typedef NT_ NT; + typedef size_t ID; + + Box_with_index_d() {} + Box_with_index_d( ID i) : m_index(i) {} + Box_with_index_d( bool complete, ID i): Base(complete), m_index(i) {} + Box_with_index_d(NT l[N], NT h[N], ID i) : Base( l, h), m_index(i) {} + Box_with_index_d( const Bbox_2& b, ID i) : Base( b), m_index(i) {} + Box_with_index_d( const Bbox_3& b, ID i) : Base( b), m_index(i) {} + ID index() const { return m_index; } +}; + +// Specialization for ID_FROM_INDEX policy +template +class Box_with_index_d + : public Box_d< NT_, N, ID_NONE> { +protected: + size_t m_index; +public: + typedef Box_d< NT_, N, ID_NONE> Base; + typedef NT_ NT; + typedef size_t ID; + + Box_with_index_d() {} + Box_with_index_d( ID i) : m_index(i) {} + Box_with_index_d( bool complete, ID i): Base(complete), m_index(i) {} + Box_with_index_d(NT l[N], NT h[N], ID i) : Base( l, h), m_index(i) {} + Box_with_index_d( const Bbox_2& b, ID i) : Base( b), m_index(i) {} + Box_with_index_d( const Bbox_3& b, ID i) : Base( b), m_index(i) {} + ID index() const { return m_index; } + ID id() const { return m_index; } +}; +#endif +} + +template +void float_snap_rounding_2(PointRange& points, + PolylineRange& segments) +{ + using Point_2 = std::remove_cv_t::value_type>; + using Kernel = typename Kernel_traits::Kernel; + using Segment_2 = typename Kernel::Segment_2; + using Vector_2 = typename Kernel::Vector_2; + using PBox = CGAL::Box_intersection_d::Box_with_index_d; + using SBox = CGAL::Box_intersection_d::Box_with_index_d; + + auto round_bound=[](Point_2& p){ + return std::pow(p.x().approx().sup()-p.x().approx().inf(),2)+std::pow(p.y().approx().sup()-p.y().approx().inf(), 2); + }; + + + //Kd-Tree to exhibits pairs + std::vector points_boxes; + std::vector segs_boxes; + for(size_t i=0; i Date: Mon, 17 Mar 2025 20:00:11 +0100 Subject: [PATCH 02/41] transform to same API as existing snap_rounding, remove duplicate, efficient subdivision of segments by vertices --- .../Snap_rounding_2/float_snap_rounding.cpp | 69 ++--- .../include/CGAL/Float_snap_rounding_2.h | 239 +++++++++++++++--- .../include/CGAL/Surface_sweep_2.h | 2 +- .../CGAL/Surface_sweep_2/Subcurves_visitor.h | 1 + 4 files changed, 223 insertions(+), 88 deletions(-) diff --git a/Snap_rounding_2/examples/Snap_rounding_2/float_snap_rounding.cpp b/Snap_rounding_2/examples/Snap_rounding_2/float_snap_rounding.cpp index 2cf4203bdb8..4c4aecd6b3a 100644 --- a/Snap_rounding_2/examples/Snap_rounding_2/float_snap_rounding.cpp +++ b/Snap_rounding_2/examples/Snap_rounding_2/float_snap_rounding.cpp @@ -1,71 +1,46 @@ #include #include #include +#include +#include -typedef CGAL::Exact_predicates_exact_constructions_kernel Kernel; -typedef Kernel::Segment_2 Segment_2; +typedef CGAL::Exact_predicates_exact_constructions_kernel Kernel; +typedef CGAL::Arr_segment_traits_2 Traits_2; +typedef Traits_2::Curve_2 Segment_2; typedef Kernel::Point_2 Point_2; typedef Kernel::FT FT; int main() { std::vector pts; - std::vector< std::vector > segs; + std::vector< Segment_2 > segs; - Kernel::FT eps(std::pow(2, -60)); - Kernel::FT one(1.); - Kernel::FT two(2.); + Kernel::FT e(std::pow(2, -60)); - pts.emplace_back(one-eps, one); - pts.emplace_back(-one+eps, -one+2*eps); - pts.emplace_back(eps/2, eps/2); - pts.emplace_back(one, -one); - pts.emplace_back(0, two-eps/2); - pts.emplace_back(two, 0); - pts.emplace_back(-two+eps, -FT(4.)); - pts.emplace_back(two, two); - pts.emplace_back(-two, two); - - segs.push_back(std::vector()); - segs[0].push_back(0); - segs[0].push_back(1); - segs.push_back(std::vector()); - segs[1].push_back(2); - segs[1].push_back(3); - segs.push_back(std::vector()); - segs[2].push_back(4); - segs[2].push_back(5); - segs.push_back(std::vector()); - segs[3].push_back(4); - segs[3].push_back(6); - segs.push_back(std::vector()); - segs[4].push_back(7); - segs[4].push_back(8); + segs.emplace_back(Point_2(1-e, 1), Point_2(-1-e, -1+2*e)); + segs.emplace_back(Point_2(e/2, e/2), Point_2(1, -1)); + segs.emplace_back(Point_2(0, 2-e/2), Point_2(2, 0)); + segs.emplace_back(Point_2(0, 2-e/2), Point_2(-2+e, -4)); + segs.emplace_back(Point_2(-2, 2), Point_2(2, 2)); + segs.emplace_back(Point_2(7, 7), Point_2(7+e, 7+e)); + segs.emplace_back(Point_2(5, 7-e), Point_2(9, 7-e)); std::cout << "Input" << std::endl; - std::cout << "Points" << std::endl; - for(auto &p: pts) - std::cout << p << std::endl; - std::cout << "Segments:" << std::endl; for(auto &seg: segs){ - for(auto pi : seg) - std::cout << pi << " "; - std::cout << std::endl; + std::cout << seg.source() << " -> " << seg.target() << std::endl; } std::cout << "\n\n" << std::endl; - CGAL::float_snap_rounding_2(pts, segs); + std::vector< std::vector > out; + CGAL::double_snap_rounding_2(segs.begin(), segs.end(), out); std::cout << "Output" << std::endl; - std::cout << "Points" << std::endl; - for(auto &p: pts) - std::cout << p << std::endl; - std::cout << "Segments:" << std::endl; - for(auto &seg: segs){ - for(auto pi : seg) - std::cout << pi << " "; + for(auto &poly: out){ + std::cout << poly[0]; + for(size_t i=1; i!=poly.size(); ++i) + std::cout << " -> " << poly[i]; std::cout << std::endl; } - return(0); + return 0; } diff --git a/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h b/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h index efb25b2efca..4f60d5a9b13 100644 --- a/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h +++ b/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h @@ -13,13 +13,11 @@ #ifndef CGAL_FLOAT_SNAP_ROUNDING_2_H #define CGAL_FLOAT_SNAP_ROUNDING_2_H -#include - #include #include #include -#include +// #include #include #include #include @@ -96,85 +94,246 @@ public: #endif } -template -void float_snap_rounding_2(PointRange& points, - PolylineRange& segments) +template +void double_snap_rounding_2_disjoint(PointsRange &pts, + PolylineRange &polylines) { - using Point_2 = std::remove_cv_t::value_type>; + using Point_2 = std::remove_cv_t::value_type>; using Kernel = typename Kernel_traits::Kernel; - using Segment_2 = typename Kernel::Segment_2; - using Vector_2 = typename Kernel::Vector_2; + using Comparison_result = typename Kernel::Comparison_result; + using Segment_2 = typename Kernel::Segment_2; + using Vector_2 = typename Kernel::Vector_2; using PBox = CGAL::Box_intersection_d::Box_with_index_d; using SBox = CGAL::Box_intersection_d::Box_with_index_d; - auto round_bound=[](Point_2& p){ - return std::pow(p.x().approx().sup()-p.x().approx().inf(),2)+std::pow(p.y().approx().sup()-p.y().approx().inf(), 2); + auto comp_by_x=[&](size_t ai, size_t bi){ + return compare(pts[ai].x(),pts[bi].x())==SMALLER; + }; + auto comp_by_y=[&](size_t ai, size_t bi){ + return compare(pts[ai].y(),pts[bi].y())==SMALLER; + }; + auto comp_by_x_first=[&](size_t ai, size_t bi){ + Comparison_result res=compare(pts[ai].x(),pts[bi].x()); + if(res==EQUAL) + return compare(pts[ai].y(),pts[bi].y())==SMALLER; + return res==SMALLER; + }; + auto comp_by_y_first=[&](size_t ai, size_t bi){ + Comparison_result res=compare(pts[ai].y(),pts[bi].y()); + if(res==EQUAL) + return compare(pts[ai].x(),pts[bi].x())==SMALLER; + return res==SMALLER; }; + // Compute the order of the points along the 2 axis + // Sorted the points may perform exact computations and thus refine the intervals of the coordinates values + // This refine ensures that the order of the points will be preserved by the rounding + using Iterator_set_x = typename std::set::iterator; + using Iterator_set_y = typename std::set::iterator; + std::set p_sort_by_x(comp_by_x_first); + std::set p_sort_by_y(comp_by_y_first); + for(size_t i=0; i!=pts.size(); ++i) + { + p_sort_by_x.insert(i); + p_sort_by_y.insert(i); + } + // std::vector p_sort_by_y(pts.size(),0); + // std::iota(p_sort_by_x.begin(),p_sort_by_x.end(),0); + // std::iota(p_sort_by_y.begin(),p_sort_by_y.end(),0); + // std::sort(p_sort_by_x.begin(),p_sort_by_x.end(),comp_by_x); + // std::sort(p_sort_by_y.begin(),p_sort_by_y.end(),comp_by_y); + //Kd-Tree to exhibits pairs std::vector points_boxes; std::vector segs_boxes; - for(size_t i=0; i unique_points(p_sort_by_x.begin(),p_sort_by_x.end()); + std::sort(unique_points.begin(),unique_points.end(),comp_by_x_first); + std::vector new_pts; + std::vector old_to_new_index(pts.size()); + for(size_t i=0; i!=pts.size(); ++i){ + if(i==0 || (pts[unique_points[i]]!=pts[unique_points[i-1]])) + new_pts.push_back(pts[unique_points[i]]); + old_to_new_index[unique_points[i]]=new_pts.size()-1; + } + + std::swap(pts, new_pts); + // Update the polylines by remapping the old indices to new indices + for (auto& polyline : polylines) { + std::vector updated_polyline; + for (size_t i=0; i vec_sort_by_x_first(p_sort_by_x.begin(),p_sort_by_x.end()); + // std::vector< size_t > vec_sort_by_y_first(p_sort_by_y.begin(),p_sort_by_y.end()); + // std::sort(vec_sort_by_x_first.begin(),vec_sort_by_x_first.end(),comp_by_x_first); + // std::sort(vec_sort_by_y_first.begin(),vec_sort_by_y_first.end(),comp_by_y_first); + for(auto &poly: polylines){ + std::vector updated_polyline; + updated_polyline.push_back(poly[0]); + for(size_t i=1; i!=poly.size(); ++i){ + if(pts[poly[i-1]].x()==pts[poly[i]].x()){ + Iterator_set_x start, end; + if(comp_by_x_first(poly[i-1],poly[i])){ + start=p_sort_by_x.upper_bound(poly[i-1]); + end=p_sort_by_x.lower_bound(poly[i]); + } else { + start=p_sort_by_x.upper_bound(poly[i]); + end=p_sort_by_x.lower_bound(poly[i-1]); + } + for(auto it=start; it!=end; ++it){ + updated_polyline.push_back(*it); + } + } + if(pts[poly[i-1]].y()==pts[poly[i]].y()){ + Iterator_set_y start, end; + if(comp_by_y_first(poly[i-1],poly[i])){ + start=p_sort_by_y.upper_bound(poly[i-1]); + end=p_sort_by_y.lower_bound(poly[i]); + } else { + start=p_sort_by_y.upper_bound(poly[i]); + end=p_sort_by_y.lower_bound(poly[i-1]); + } + for(auto it=start; it!=end; ++it){ + updated_polyline.push_back(*it); + } + } + updated_polyline.push_back(poly[i]); + } + std::swap(poly, updated_polyline); + } + // compute_subcurves(input_begin, input_end, polylines); +} + +template +typename OutputContainer::iterator double_snap_rounding_2(InputIterator input_begin, + InputIterator input_end, + OutputContainer& output) +{ + using Segment_2 = std::remove_cv_t::value_type>; + using Polyline = std::remove_cv_t::value_type>; + using Point_2 = typename Default_arr_traits::Traits::Point_2; + + std::vector segs; + compute_subcurves(input_begin, input_end, std::back_inserter(segs)); + + std::set unique_point_set; + std::map point_to_index; + std::vector pts; + std::vector< std::vector< size_t> > polylines; + + // Transform range of the segments in the range of points and polyline of indexes + for(typename std::vector::iterator it=segs.begin(); it!=segs.end(); ++it) + { + const Point_2& p1 = it->source(); + const Point_2& p2 = it->target(); + + // Check and insert the first endpoint if it's not already added + if (unique_point_set.find(p1) == unique_point_set.end()) { + unique_point_set.insert(p1); + pts.push_back(p1); + point_to_index[p1] = pts.size() - 1; + } + if (unique_point_set.find(p2) == unique_point_set.end()) { + unique_point_set.insert(p2); + pts.push_back(p2); + point_to_index[p2] = pts.size() - 1; + } + } + + for(InputIterator it=input_begin; it!=input_end; ++it) + { + size_t index1 = point_to_index[it->source()]; + size_t index2 = point_to_index[it->target()]; + polylines.push_back({index1, index2}); + } + + double_snap_rounding_2_disjoint(pts, polylines); + for(auto &poly: polylines){ + Polyline new_line; + for(size_t pi: poly){ + new_line.push_back(pts[pi]); + } + output.push_back(new_line); + } + + return output.begin(); } } //namespace CGAL diff --git a/Surface_sweep_2/include/CGAL/Surface_sweep_2.h b/Surface_sweep_2/include/CGAL/Surface_sweep_2.h index 2c007c5955a..93d46ecc8a9 100644 --- a/Surface_sweep_2/include/CGAL/Surface_sweep_2.h +++ b/Surface_sweep_2/include/CGAL/Surface_sweep_2.h @@ -198,7 +198,7 @@ protected: /*! Compute intersections between the two given curves. * If the two curves intersect, create a new event (or use the event that - * already exits in the intersection point) and insert the curves to the + * already exists in the intersection point) and insert the curves to the * event. * \param curve1 The first curve. * \param curve2 The second curve. diff --git a/Surface_sweep_2/include/CGAL/Surface_sweep_2/Subcurves_visitor.h b/Surface_sweep_2/include/CGAL/Surface_sweep_2/Subcurves_visitor.h index 9a660c211fd..e47eb16e3b5 100644 --- a/Surface_sweep_2/include/CGAL/Surface_sweep_2/Subcurves_visitor.h +++ b/Surface_sweep_2/include/CGAL/Surface_sweep_2/Subcurves_visitor.h @@ -76,6 +76,7 @@ public: m_overlapping(overlapping) {} + template void sweep(CurveIterator begin, CurveIterator end) { From d229e9fabd97f9e4b66d30f11df280ceb9b9ddaf Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?L=C3=A9o=20Valque?= Date: Fri, 21 Mar 2025 17:02:34 +0100 Subject: [PATCH 03/41] Test of float snap rounding 2 and intense debug --- .../Snap_rounding_2/float_snap_rounding.cpp | 5 +- .../include/CGAL/Float_snap_rounding_2.h | 174 +++++++++++++++--- .../test/Snap_rounding_2/CMakeLists.txt | 17 +- 3 files changed, 164 insertions(+), 32 deletions(-) diff --git a/Snap_rounding_2/examples/Snap_rounding_2/float_snap_rounding.cpp b/Snap_rounding_2/examples/Snap_rounding_2/float_snap_rounding.cpp index 4c4aecd6b3a..4866de1edf0 100644 --- a/Snap_rounding_2/examples/Snap_rounding_2/float_snap_rounding.cpp +++ b/Snap_rounding_2/examples/Snap_rounding_2/float_snap_rounding.cpp @@ -9,13 +9,14 @@ typedef CGAL::Arr_segment_traits_2 Traits_2; typedef Traits_2::Curve_2 Segment_2; typedef Kernel::Point_2 Point_2; typedef Kernel::FT FT; +typedef std::vector Polyline_2; int main() { std::vector pts; std::vector< Segment_2 > segs; - Kernel::FT e(std::pow(2, -60)); + FT e(std::pow(2, -60)); segs.emplace_back(Point_2(1-e, 1), Point_2(-1-e, -1+2*e)); segs.emplace_back(Point_2(e/2, e/2), Point_2(1, -1)); @@ -31,7 +32,7 @@ int main() } std::cout << "\n\n" << std::endl; - std::vector< std::vector > out; + std::vector< Polyline_2 > out; CGAL::double_snap_rounding_2(segs.begin(), segs.end(), out); std::cout << "Output" << std::endl; diff --git a/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h b/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h index 4f60d5a9b13..7720e709dd9 100644 --- a/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h +++ b/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h @@ -125,10 +125,12 @@ void double_snap_rounding_2_disjoint(PointsRange &pts, return res==SMALLER; }; - +#ifdef DOUBLE_2D_SNAP_VERBOSE + std::cout << "Sort the input points" << std::endl; +#endif // Compute the order of the points along the 2 axis // Sorted the points may perform exact computations and thus refine the intervals of the coordinates values - // This refine ensures that the order of the points will be preserved by the rounding + // This refine ensures that the order of the points will be preserved when rounded using Iterator_set_x = typename std::set::iterator; using Iterator_set_y = typename std::set::iterator; std::set p_sort_by_x(comp_by_x_first); @@ -138,11 +140,6 @@ void double_snap_rounding_2_disjoint(PointsRange &pts, p_sort_by_x.insert(i); p_sort_by_y.insert(i); } - // std::vector p_sort_by_y(pts.size(),0); - // std::iota(p_sort_by_x.begin(),p_sort_by_x.end(),0); - // std::iota(p_sort_by_y.begin(),p_sort_by_y.end(),0); - // std::sort(p_sort_by_x.begin(),p_sort_by_x.end(),comp_by_x); - // std::sort(p_sort_by_y.begin(),p_sort_by_y.end(),comp_by_y); //Kd-Tree to exhibits pairs std::vector points_boxes; @@ -152,8 +149,10 @@ void double_snap_rounding_2_disjoint(PointsRange &pts, for(size_t i=0; i::infinity())-p.x().approx().inf(),2)+ + std::pow(std::nextafter(p.y().approx().sup(),std::numeric_limits::infinity())-p.y().approx().inf(),2); }; auto callback=[&](PBox &bp, SBox &bseg){ size_t pi=bp.index(); @@ -161,35 +160,54 @@ void double_snap_rounding_2_disjoint(PointsRange &pts, size_t si1=polylines[bseg.index()][0]; size_t si2=polylines[bseg.index()][1]; - // the point is a vertex of the segment if((pi==si1) || (pi==si2)) return; Point_2& p= pts[bp.index()]; Segment_2 seg(pts[polylines[bseg.index()][0]], pts[polylines[bseg.index()][1]]); - double round_bound_s=(std::max)(round_bound(pts[si1]), round_bound(pts[si2])); + // (A+B)^2 <= 4*max(A^2,B^2) and we take some margin + double bound=16*(std::max)(round_bound(pts[pi]), (std::max)(round_bound(pts[si1]), round_bound(pts[si2]))); - if(possibly(Kernel().compare_squared_distance_2_object()(p, seg, round_bound_s)!=CGAL::LARGER)){ + // If the segment is closed to the vertex and we subdivide it at same x coordinate of that vertex + if(possibly(Kernel().compare_squared_distance_2_object()(p, seg, bound)!=CGAL::LARGER) && + compare(seg.source().x(),p.x())!=compare(seg.target().x(),p.x())) + { pts.emplace_back(p.x(), seg.supporting_line().y_at_x(p.x())); auto pair=p_sort_by_x.insert(pts.size()-1); if(pair.second) p_sort_by_y.insert(pts.size()-1); else - pts.pop_back(); + pts.pop_back(); // Remove the new point if it is already exist polylines[si].emplace_back(*pair.first); } }; - CGAL::box_intersection_d(points_boxes.begin(), points_boxes.end(), segs_boxes.begin(), segs_boxes.end(), callback); +#ifdef DOUBLE_2D_SNAP_VERBOSE + std::cout << "Exhibit pairs of possible intersections" << std::endl; +#endif - //sort new vertices + do{ + size_t size_before=pts.size(); + CGAL::box_intersection_d(points_boxes.begin(), points_boxes.end(), segs_boxes.begin(), segs_boxes.end(), callback); + points_boxes.clear(); + // The new vertices may intersect a segment when rounded, we repeat until they are not new vertices +#ifdef DOUBLE_2D_SNAP_VERBOSE + std::cout << points_boxes.size()-size_before << " subdivisions performed" << std::endl; +#endif + for(size_t i=size_before; i unique_points(p_sort_by_x.begin(),p_sort_by_x.end()); std::sort(unique_points.begin(),unique_points.end(),comp_by_x_first); std::vector new_pts; @@ -220,7 +245,6 @@ void double_snap_rounding_2_disjoint(PointsRange &pts, } std::swap(pts, new_pts); - // Update the polylines by remapping the old indices to new indices for (auto& polyline : polylines) { std::vector updated_polyline; for (size_t i=0; i vec_sort_by_x_first(p_sort_by_x.begin(),p_sort_by_x.end()); - // std::vector< size_t > vec_sort_by_y_first(p_sort_by_y.begin(),p_sort_by_y.end()); - // std::sort(vec_sort_by_x_first.begin(),vec_sort_by_x_first.end(),comp_by_x_first); - // std::sort(vec_sort_by_y_first.begin(),vec_sort_by_y_first.end(),comp_by_y_first); + for(auto &poly: polylines){ std::vector updated_polyline; updated_polyline.push_back(poly[0]); for(size_t i=1; i!=poly.size(); ++i){ if(pts[poly[i-1]].x()==pts[poly[i]].x()){ Iterator_set_x start, end; + // Get all vertices between the two endpoints along x order if(comp_by_x_first(poly[i-1],poly[i])){ start=p_sort_by_x.upper_bound(poly[i-1]); end=p_sort_by_x.lower_bound(poly[i]); @@ -257,12 +283,14 @@ void double_snap_rounding_2_disjoint(PointsRange &pts, start=p_sort_by_x.upper_bound(poly[i]); end=p_sort_by_x.lower_bound(poly[i-1]); } + // Add all endpoints between them to the polyline for(auto it=start; it!=end; ++it){ updated_polyline.push_back(*it); } } if(pts[poly[i-1]].y()==pts[poly[i]].y()){ Iterator_set_y start, end; + // Get all vertices between the two endpoints along y order if(comp_by_y_first(poly[i-1],poly[i])){ start=p_sort_by_y.upper_bound(poly[i-1]); end=p_sort_by_y.lower_bound(poly[i]); @@ -270,6 +298,7 @@ void double_snap_rounding_2_disjoint(PointsRange &pts, start=p_sort_by_y.upper_bound(poly[i]); end=p_sort_by_y.lower_bound(poly[i-1]); } + // Add all endpoints between them to the polyline for(auto it=start; it!=end; ++it){ updated_polyline.push_back(*it); } @@ -278,9 +307,12 @@ void double_snap_rounding_2_disjoint(PointsRange &pts, } std::swap(poly, updated_polyline); } - // compute_subcurves(input_begin, input_end, polylines); } + +/* +TODO doc +*/ template typename OutputContainer::iterator double_snap_rounding_2(InputIterator input_begin, InputIterator input_end, @@ -290,9 +322,15 @@ typename OutputContainer::iterator double_snap_rounding_2(InputIterator input_ using Polyline = std::remove_cv_t::value_type>; using Point_2 = typename Default_arr_traits::Traits::Point_2; - std::vector segs; + std::vector segs(input_begin, input_end); +#ifdef DOUBLE_2D_SNAP_VERBOSE + std::cout << "Solved intersections" << std::endl; +#endif compute_subcurves(input_begin, input_end, std::back_inserter(segs)); +#ifdef DOUBLE_2D_SNAP_VERBOSE + std::cout << "Change format to range of points and indexes" << std::endl; +#endif std::set unique_point_set; std::map point_to_index; std::vector pts; @@ -317,25 +355,103 @@ typename OutputContainer::iterator double_snap_rounding_2(InputIterator input_ } } - for(InputIterator it=input_begin; it!=input_end; ++it) + for(InputIterator it=segs.begin(); it!=segs.end(); ++it) { size_t index1 = point_to_index[it->source()]; size_t index2 = point_to_index[it->target()]; polylines.push_back({index1, index2}); } + // Main algorithm double_snap_rounding_2_disjoint(pts, polylines); + +#ifdef DOUBLE_2D_SNAP_VERBOSE + std::cout << "Build output" << std::endl; +#endif + + // Output polylines + output.clear(); for(auto &poly: polylines){ Polyline new_line; - for(size_t pi: poly){ - new_line.push_back(pts[pi]); - } + for(size_t pi: poly) + new_line.push_back(pts[pi]); output.push_back(new_line); } return output.begin(); } +/* +TODO doc +*/ +template +typename OutputContainer::iterator compute_snap_subcurves_2(InputIterator input_begin, + InputIterator input_end, + OutputContainer& output) +{ + using Segment_2 = std::remove_cv_t::value_type>; + using Point_2 = typename Default_arr_traits::Traits::Point_2; + + std::vector segs; +#ifdef DOUBLE_2D_SNAP_VERBOSE + std::cout << "Solved intersections" << std::endl; +#endif + + compute_subcurves(input_begin, input_end, std::back_inserter(segs)); + +#ifdef DOUBLE_2D_SNAP_VERBOSE + std::cout << "Change format to range of points and indexes" << std::endl; +#endif + std::set unique_point_set; + std::map point_to_index; + std::vector pts; + std::vector< std::vector< size_t> > polylines; + + for(typename std::vector::iterator it=segs.begin(); it!=segs.end(); ++it) + { + const Point_2& p1 = it->source(); + const Point_2& p2 = it->target(); + + // Check and insert the endpoints if they are not already added + if (unique_point_set.find(p1) == unique_point_set.end()) { + unique_point_set.insert(p1); + pts.push_back(p1); + point_to_index[p1] = pts.size() - 1; + } + if (unique_point_set.find(p2) == unique_point_set.end()) { + unique_point_set.insert(p2); + pts.push_back(p2); + point_to_index[p2] = pts.size() - 1; + } + } + + for(InputIterator it=segs.begin(); it!=segs.end(); ++it) + { + size_t index1 = point_to_index[it->source()]; + size_t index2 = point_to_index[it->target()]; + polylines.push_back({index1, index2}); + } + + // Main algorithm + double_snap_rounding_2_disjoint(pts, polylines); + +#ifdef DOUBLE_2D_SNAP_VERBOSE + std::cout << "Build output" << std::endl; +#endif + + // Output a range of segments while removing duplicate ones + std::set< std::pair > set_out_segs; + output.clear(); + for(auto &poly: polylines){ + for(size_t i=1; i Date: Mon, 24 Mar 2025 08:50:10 +0100 Subject: [PATCH 04/41] Box_with_index has known explicit ID --- Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h b/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h index 7720e709dd9..ca2b949f657 100644 --- a/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h +++ b/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h @@ -33,11 +33,11 @@ namespace CGAL { namespace Box_intersection_d { template -class Box_with_index_d: public Box_d< NT_, N, ID_NONE> { +class Box_with_index_d: public Box_d< NT_, N, ID_EXPLICIT> { protected: size_t m_index; public: - typedef Box_d< NT_, N, ID_NONE> Base; + typedef Box_d< NT_, N, ID_EXPLICIT> Base; typedef NT_ NT; typedef size_t ID; @@ -48,7 +48,7 @@ public: Box_with_index_d( const Bbox_2& b, ID i) : Base( b), m_index(i) {} Box_with_index_d( const Bbox_3& b, ID i) : Base( b), m_index(i) {} ID index() const { return m_index; } - ID id() const { return m_index; } + // ID id() const { return m_index; } }; // Generic template signature of boxes, specialized for ID_FROM_HANDLE policy From a35acd9ba6c66455d2b9d9ad041d82d26914732f Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?L=C3=A9o=20Valque?= Date: Mon, 24 Mar 2025 08:55:04 +0100 Subject: [PATCH 05/41] remove blank line unvolontary added to Subcruves_visitor --- Surface_sweep_2/include/CGAL/Surface_sweep_2/Subcurves_visitor.h | 1 - 1 file changed, 1 deletion(-) diff --git a/Surface_sweep_2/include/CGAL/Surface_sweep_2/Subcurves_visitor.h b/Surface_sweep_2/include/CGAL/Surface_sweep_2/Subcurves_visitor.h index e47eb16e3b5..9a660c211fd 100644 --- a/Surface_sweep_2/include/CGAL/Surface_sweep_2/Subcurves_visitor.h +++ b/Surface_sweep_2/include/CGAL/Surface_sweep_2/Subcurves_visitor.h @@ -76,7 +76,6 @@ public: m_overlapping(overlapping) {} - template void sweep(CurveIterator begin, CurveIterator end) { From ad3bde0bf497cce2e542439828fdfaa9071a96ae Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?L=C3=A9o=20Valque?= Date: Mon, 24 Mar 2025 09:36:35 +0100 Subject: [PATCH 06/41] some cleaning --- .../Snap_rounding_2/float_snap_rounding.cpp | 4 +- .../include/CGAL/Float_snap_rounding_2.h | 27 ++--- .../test_float_snap_rounding_2.cpp | 104 ++++++++++++++++++ 3 files changed, 117 insertions(+), 18 deletions(-) create mode 100644 Snap_rounding_2/test/Snap_rounding_2/test_float_snap_rounding_2.cpp diff --git a/Snap_rounding_2/examples/Snap_rounding_2/float_snap_rounding.cpp b/Snap_rounding_2/examples/Snap_rounding_2/float_snap_rounding.cpp index 4866de1edf0..21ff49bc960 100644 --- a/Snap_rounding_2/examples/Snap_rounding_2/float_snap_rounding.cpp +++ b/Snap_rounding_2/examples/Snap_rounding_2/float_snap_rounding.cpp @@ -1,8 +1,8 @@ -#include +// #define DOUBLE_2D_SNAP_VERBOSE + #include #include #include -#include typedef CGAL::Exact_predicates_exact_constructions_kernel Kernel; typedef CGAL::Arr_segment_traits_2 Traits_2; diff --git a/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h b/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h index ca2b949f657..db594ce8034 100644 --- a/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h +++ b/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h @@ -15,19 +15,12 @@ #include -#include -#include -// #include #include #include #include #include -#include #include -#include -#include -#include -#include +#include namespace CGAL { namespace Box_intersection_d { @@ -112,6 +105,7 @@ void double_snap_rounding_2_disjoint(PointsRange &pts, auto comp_by_y=[&](size_t ai, size_t bi){ return compare(pts[ai].y(),pts[bi].y())==SMALLER; }; + // Total order is required to using set auto comp_by_x_first=[&](size_t ai, size_t bi){ Comparison_result res=compare(pts[ai].x(),pts[bi].x()); if(res==EQUAL) @@ -131,6 +125,7 @@ void double_snap_rounding_2_disjoint(PointsRange &pts, // Compute the order of the points along the 2 axis // Sorted the points may perform exact computations and thus refine the intervals of the coordinates values // This refine ensures that the order of the points will be preserved when rounded + // However, except for this reason these sets are unused using Iterator_set_x = typename std::set::iterator; using Iterator_set_y = typename std::set::iterator; std::set p_sort_by_x(comp_by_x_first); @@ -141,7 +136,7 @@ void double_snap_rounding_2_disjoint(PointsRange &pts, p_sort_by_y.insert(i); } - //Kd-Tree to exhibits pairs + //Prepare boxes for box_intersection_d std::vector points_boxes; std::vector segs_boxes; for(size_t i=0; i::infinity())-p.x().approx().inf(),2)+ std::pow(std::nextafter(p.y().approx().sup(),std::numeric_limits::infinity())-p.y().approx().inf(),2); }; + + // Callback used for box_intersection_d auto callback=[&](PBox &bp, SBox &bseg){ size_t pi=bp.index(); size_t si=bseg.index(); @@ -166,10 +163,10 @@ void double_snap_rounding_2_disjoint(PointsRange &pts, Point_2& p= pts[bp.index()]; Segment_2 seg(pts[polylines[bseg.index()][0]], pts[polylines[bseg.index()][1]]); - // (A+B)^2 <= 4*max(A^2,B^2) and we take some margin + // (A+B)^2 <= 4*max(A^2,B^2), we take some margin double bound=16*(std::max)(round_bound(pts[pi]), (std::max)(round_bound(pts[si1]), round_bound(pts[si2]))); - // If the segment is closed to the vertex and we subdivide it at same x coordinate of that vertex + // If the segment is closed to the vertex, we subdivide it at same x coordinate that this vertex if(possibly(Kernel().compare_squared_distance_2_object()(p, seg, bound)!=CGAL::LARGER) && compare(seg.source().x(),p.x())!=compare(seg.target().x(),p.x())) { @@ -179,7 +176,7 @@ void double_snap_rounding_2_disjoint(PointsRange &pts, p_sort_by_y.insert(pts.size()-1); else pts.pop_back(); // Remove the new point if it is already exist - polylines[si].emplace_back(*pair.first); + polylines[si].emplace_back(*pair.first); // Some duplicates maybe introduced, it will be removed later } }; @@ -191,7 +188,7 @@ void double_snap_rounding_2_disjoint(PointsRange &pts, size_t size_before=pts.size(); CGAL::box_intersection_d(points_boxes.begin(), points_boxes.end(), segs_boxes.begin(), segs_boxes.end(), callback); points_boxes.clear(); - // The new vertices may intersect a segment when rounded, we repeat until they are not new vertices + // The new vertices may intersect another segment when rounded, we repeat until they are not new vertices #ifdef DOUBLE_2D_SNAP_VERBOSE std::cout << points_boxes.size()-size_before << " subdivisions performed" << std::endl; #endif @@ -260,7 +257,7 @@ void double_snap_rounding_2_disjoint(PointsRange &pts, // The algorithm prevents the a vertex that goes through a segment but a vertex may lie on an horizontal/vertical segments after rounding std::cout << "Subdivide horizontal and vertical segments with vertices on them" << std::endl; #endif - //The order may have changed (Example: (1,1)<(1,2)<(1+e,1)) + //The order may have changed (Example: (1,1)<(1,2)<(1+e,1)), we recompute it p_sort_by_x.clear(); p_sort_by_y.clear(); for(size_t i=0; i!=pts.size(); ++i) @@ -342,7 +339,6 @@ typename OutputContainer::iterator double_snap_rounding_2(InputIterator input_ const Point_2& p1 = it->source(); const Point_2& p2 = it->target(); - // Check and insert the first endpoint if it's not already added if (unique_point_set.find(p1) == unique_point_set.end()) { unique_point_set.insert(p1); pts.push_back(p1); @@ -412,7 +408,6 @@ typename OutputContainer::iterator compute_snap_subcurves_2(InputIterator inpu const Point_2& p1 = it->source(); const Point_2& p2 = it->target(); - // Check and insert the endpoints if they are not already added if (unique_point_set.find(p1) == unique_point_set.end()) { unique_point_set.insert(p1); pts.push_back(p1); diff --git a/Snap_rounding_2/test/Snap_rounding_2/test_float_snap_rounding_2.cpp b/Snap_rounding_2/test/Snap_rounding_2/test_float_snap_rounding_2.cpp new file mode 100644 index 00000000000..190047bb59d --- /dev/null +++ b/Snap_rounding_2/test/Snap_rounding_2/test_float_snap_rounding_2.cpp @@ -0,0 +1,104 @@ +// #define DOUBLE_2D_SNAP_VERBOSE + +#include +#include +#include +#include + +#include +#include + +typedef CGAL::Exact_predicates_exact_constructions_kernel Kernel; +typedef CGAL::Arr_segment_traits_2 Traits_2; +typedef Traits_2::Curve_2 Segment_2; +typedef Kernel::Point_2 Point_2; +typedef Kernel::Vector_2 Vector_2; +typedef std::vector Polyline_2; +typedef Kernel::FT FT; + +Point_2 random_point(CGAL::Random& r, double m = -1, double M = 1) +{ + return Point_2(FT(r.get_double(m, M)), FT(r.get_double(m, M))); +} + +Point_2 random_noise_point(CGAL::Random& r, double m, double M, Vector_2 v) +{ + return Point_2(FT(r.get_double(m, M) + CGAL::to_double(v.x())), FT(r.get_double(m, M) + CGAL::to_double(v.y()))); +} + +void test_fully_random(CGAL::Random &r, size_t nb_segments){ + std::cout << "Test fully random" << std::endl; + std::vector segs; + std::vector out; + + for(size_t i=0; i!=nb_segments; ++i) + segs.emplace_back(random_point(r), random_point(r)); + + CGAL::Real_timer t; + t.start(); + CGAL::compute_snap_subcurves_2(segs.begin(), segs.end(), out); + t.stop(); + + std::cout << t.time() << "sec" << std::endl; + assert(!CGAL::do_curves_intersect(out.begin(), out.end())); +} + +void test_almost_indentical_segments(CGAL::Random &r, size_t nb_segments, Vector_2 source, Vector_2 target){ + // Difficult test + std::cout << "Test with almost identical segments from " << source << " to " << target << std::endl; + std::vector segs; + std::vector out; + + double eps=std::pow(2,-45); + for(size_t i=0; i!=nb_segments; ++i){ + segs.emplace_back(random_noise_point(r, -eps, eps, source), random_noise_point(r, -eps, eps, target)); + } + + CGAL::Real_timer t; + t.start(); + CGAL::compute_snap_subcurves_2(segs.begin(), segs.end(), out); + t.stop(); + + std::cout << t.time() << "sec" << std::endl; + assert(!CGAL::do_curves_intersect(out.begin(), out.end())); +} + +void test_multi_almost_indentical_segments(CGAL::Random &r, size_t nb_segments){ + for(double x1=-1; x1<=1; ++x1) + for(double y1=-1; y1<=1; ++y1) + for(double x2=x1; x2<=1; ++x2) + for(double y2=y1; y2<=1; ++y2) + if(Vector_2(x1,y1)!=Vector_2(x2,y2)) + test_almost_indentical_segments(r, nb_segments, Vector_2(x1,y1), Vector_2(x2,y2)); +} + +void test_data(){ + std::vector segs; + std::vector out; + + segs.emplace_back(Point_2( 0.353288537902858912, 0.0229645738916604246 ), + Point_2( 0.353619283885096169, 0.0230228644565691754 )); + segs.emplace_back(Point_2( 0.353619283885096169, 0.0230228644565691754 ), + Point_2( 0.354002104897765901, 0.0214915804089364851 )); + segs.emplace_back(Point_2( 0.354002104897765901, 0.0214915804089364851 ), + Point_2( 0.354002107499122531, 0.0214915700035099161 )); + segs.emplace_back(Point_2( 0.354002107499122531, 0.0214915700035099161 ), + Point_2( 0.353288537902858912, 0.0229645738916604246 )); + segs.emplace_back(Point_2( 0.354002107499122531, 0.0214915700035098953 ), + Point_2( 0.354002104897765901, 0.0214915804089364851 )); + + CGAL::compute_snap_subcurves_2(segs.begin(), segs.end(), out); + assert(!CGAL::do_curves_intersect(out.begin(), out.end())); +} + +int main(int argc,char *argv[]) +{ + CGAL::Random rp; + CGAL::Random r(argc==1?rp.get_seed():std::stoi(argv[1])); + std::cout << "random seed = " << r.get_seed() << std::endl; + std::cout << std::setprecision(17); + test_data(); + test_fully_random(r,2000); + test_multi_almost_indentical_segments(r,200); + return(0); +} From 717cd99c5c2f8e490009a4fc4355069951414a9d Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?L=C3=A9o=20Valque?= Date: Wed, 2 Apr 2025 09:33:39 +0200 Subject: [PATCH 07/41] delete specific test from data --- .../test_float_snap_rounding_2.cpp | 19 ------------------- 1 file changed, 19 deletions(-) diff --git a/Snap_rounding_2/test/Snap_rounding_2/test_float_snap_rounding_2.cpp b/Snap_rounding_2/test/Snap_rounding_2/test_float_snap_rounding_2.cpp index 190047bb59d..a575eb07a9b 100644 --- a/Snap_rounding_2/test/Snap_rounding_2/test_float_snap_rounding_2.cpp +++ b/Snap_rounding_2/test/Snap_rounding_2/test_float_snap_rounding_2.cpp @@ -72,25 +72,6 @@ void test_multi_almost_indentical_segments(CGAL::Random &r, size_t nb_segments){ test_almost_indentical_segments(r, nb_segments, Vector_2(x1,y1), Vector_2(x2,y2)); } -void test_data(){ - std::vector segs; - std::vector out; - - segs.emplace_back(Point_2( 0.353288537902858912, 0.0229645738916604246 ), - Point_2( 0.353619283885096169, 0.0230228644565691754 )); - segs.emplace_back(Point_2( 0.353619283885096169, 0.0230228644565691754 ), - Point_2( 0.354002104897765901, 0.0214915804089364851 )); - segs.emplace_back(Point_2( 0.354002104897765901, 0.0214915804089364851 ), - Point_2( 0.354002107499122531, 0.0214915700035099161 )); - segs.emplace_back(Point_2( 0.354002107499122531, 0.0214915700035099161 ), - Point_2( 0.353288537902858912, 0.0229645738916604246 )); - segs.emplace_back(Point_2( 0.354002107499122531, 0.0214915700035098953 ), - Point_2( 0.354002104897765901, 0.0214915804089364851 )); - - CGAL::compute_snap_subcurves_2(segs.begin(), segs.end(), out); - assert(!CGAL::do_curves_intersect(out.begin(), out.end())); -} - int main(int argc,char *argv[]) { CGAL::Random rp; From 1672576a87fb7e54d6685716d22c130414db0c68 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?L=C3=A9o=20Valque?= Date: Wed, 2 Apr 2025 09:35:41 +0200 Subject: [PATCH 08/41] hide include iostream behind debug macro --- Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h b/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h index db594ce8034..5ab413d20eb 100644 --- a/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h +++ b/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h @@ -13,8 +13,9 @@ #ifndef CGAL_FLOAT_SNAP_ROUNDING_2_H #define CGAL_FLOAT_SNAP_ROUNDING_2_H - +#ifdef DOUBLE_2D_SNAP_VERBOSE #include +#endif #include #include #include From 829aba2842ee35a9c527e0a298472a4ee2c0db00 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?L=C3=A9o=20Valque?= Date: Wed, 2 Apr 2025 09:36:52 +0200 Subject: [PATCH 09/41] remove unused specialization of BOx_with_index --- .../include/CGAL/Float_snap_rounding_2.h | 42 ------------------- 1 file changed, 42 deletions(-) diff --git a/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h b/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h index 5ab413d20eb..7be03b379f8 100644 --- a/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h +++ b/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h @@ -42,50 +42,8 @@ public: Box_with_index_d( const Bbox_2& b, ID i) : Base( b), m_index(i) {} Box_with_index_d( const Bbox_3& b, ID i) : Base( b), m_index(i) {} ID index() const { return m_index; } - // ID id() const { return m_index; } }; -// Generic template signature of boxes, specialized for ID_FROM_HANDLE policy -#if 0 -template -class Box_with_index_d : public Box_d< NT_, N, IdPolicy> { -protected: - size_t m_index; -public: - typedef Box_d< NT_, N, IdPolicy> Base; - typedef NT_ NT; - typedef size_t ID; - - Box_with_index_d() {} - Box_with_index_d( ID i) : m_index(i) {} - Box_with_index_d( bool complete, ID i): Base(complete), m_index(i) {} - Box_with_index_d(NT l[N], NT h[N], ID i) : Base( l, h), m_index(i) {} - Box_with_index_d( const Bbox_2& b, ID i) : Base( b), m_index(i) {} - Box_with_index_d( const Bbox_3& b, ID i) : Base( b), m_index(i) {} - ID index() const { return m_index; } -}; - -// Specialization for ID_FROM_INDEX policy -template -class Box_with_index_d - : public Box_d< NT_, N, ID_NONE> { -protected: - size_t m_index; -public: - typedef Box_d< NT_, N, ID_NONE> Base; - typedef NT_ NT; - typedef size_t ID; - - Box_with_index_d() {} - Box_with_index_d( ID i) : m_index(i) {} - Box_with_index_d( bool complete, ID i): Base(complete), m_index(i) {} - Box_with_index_d(NT l[N], NT h[N], ID i) : Base( l, h), m_index(i) {} - Box_with_index_d( const Bbox_2& b, ID i) : Base( b), m_index(i) {} - Box_with_index_d( const Bbox_3& b, ID i) : Base( b), m_index(i) {} - ID index() const { return m_index; } - ID id() const { return m_index; } -}; -#endif } template From c9df11e375151e418ff8700a4d24e856aec82710 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?L=C3=A9o=20Valque?= Date: Wed, 2 Apr 2025 09:38:06 +0200 Subject: [PATCH 10/41] replace size_t per std::size_t --- .../include/CGAL/Float_snap_rounding_2.h | 83 +++++++++---------- 1 file changed, 41 insertions(+), 42 deletions(-) diff --git a/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h b/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h index 7be03b379f8..a7cdbb9e420 100644 --- a/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h +++ b/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h @@ -29,11 +29,11 @@ namespace Box_intersection_d { template class Box_with_index_d: public Box_d< NT_, N, ID_EXPLICIT> { protected: - size_t m_index; + std::size_t m_index; public: typedef Box_d< NT_, N, ID_EXPLICIT> Base; typedef NT_ NT; - typedef size_t ID; + typedef std::size_t ID; Box_with_index_d() {} Box_with_index_d( ID i) : m_index(i) {} @@ -47,8 +47,7 @@ public: } template -void double_snap_rounding_2_disjoint(PointsRange &pts, - PolylineRange &polylines) +void double_snap_rounding_2_disjoint(PointsRange &pts, PolylineRange &polylines) { using Point_2 = std::remove_cv_t::value_type>; using Kernel = typename Kernel_traits::Kernel; @@ -58,20 +57,20 @@ void double_snap_rounding_2_disjoint(PointsRange &pts, using PBox = CGAL::Box_intersection_d::Box_with_index_d; using SBox = CGAL::Box_intersection_d::Box_with_index_d; - auto comp_by_x=[&](size_t ai, size_t bi){ + auto comp_by_x=[&](std::size_t ai, std::size_t bi){ return compare(pts[ai].x(),pts[bi].x())==SMALLER; }; - auto comp_by_y=[&](size_t ai, size_t bi){ + auto comp_by_y=[&](std::size_t ai, std::size_t bi){ return compare(pts[ai].y(),pts[bi].y())==SMALLER; }; // Total order is required to using set - auto comp_by_x_first=[&](size_t ai, size_t bi){ + auto comp_by_x_first=[&](std::size_t ai, std::size_t bi){ Comparison_result res=compare(pts[ai].x(),pts[bi].x()); if(res==EQUAL) return compare(pts[ai].y(),pts[bi].y())==SMALLER; return res==SMALLER; }; - auto comp_by_y_first=[&](size_t ai, size_t bi){ + auto comp_by_y_first=[&](std::size_t ai, std::size_t bi){ Comparison_result res=compare(pts[ai].y(),pts[bi].y()); if(res==EQUAL) return compare(pts[ai].x(),pts[bi].x())==SMALLER; @@ -85,11 +84,11 @@ void double_snap_rounding_2_disjoint(PointsRange &pts, // Sorted the points may perform exact computations and thus refine the intervals of the coordinates values // This refine ensures that the order of the points will be preserved when rounded // However, except for this reason these sets are unused - using Iterator_set_x = typename std::set::iterator; - using Iterator_set_y = typename std::set::iterator; - std::set p_sort_by_x(comp_by_x_first); - std::set p_sort_by_y(comp_by_y_first); - for(size_t i=0; i!=pts.size(); ++i) + using Iterator_set_x = typename std::set::iterator; + using Iterator_set_y = typename std::set::iterator; + std::set p_sort_by_x(comp_by_x_first); + std::set p_sort_by_y(comp_by_y_first); + for(std::size_t i=0; i!=pts.size(); ++i) { p_sort_by_x.insert(i); p_sort_by_y.insert(i); @@ -98,9 +97,9 @@ void double_snap_rounding_2_disjoint(PointsRange &pts, //Prepare boxes for box_intersection_d std::vector points_boxes; std::vector segs_boxes; - for(size_t i=0; i unique_points(p_sort_by_x.begin(),p_sort_by_x.end()); + std::vector< std::size_t > unique_points(p_sort_by_x.begin(),p_sort_by_x.end()); std::sort(unique_points.begin(),unique_points.end(),comp_by_x_first); std::vector new_pts; - std::vector old_to_new_index(pts.size()); - for(size_t i=0; i!=pts.size(); ++i){ + std::vector old_to_new_index(pts.size()); + for(std::size_t i=0; i!=pts.size(); ++i){ if(i==0 || (pts[unique_points[i]]!=pts[unique_points[i-1]])) new_pts.push_back(pts[unique_points[i]]); old_to_new_index[unique_points[i]]=new_pts.size()-1; @@ -202,9 +201,9 @@ void double_snap_rounding_2_disjoint(PointsRange &pts, std::swap(pts, new_pts); for (auto& polyline : polylines) { - std::vector updated_polyline; - for (size_t i=0; i updated_polyline; + for (std::size_t i=0; i updated_polyline; + std::vector updated_polyline; updated_polyline.push_back(poly[0]); - for(size_t i=1; i!=poly.size(); ++i){ + for(std::size_t i=1; i!=poly.size(); ++i){ if(pts[poly[i-1]].x()==pts[poly[i]].x()){ Iterator_set_x start, end; // Get all vertices between the two endpoints along x order @@ -290,7 +289,7 @@ typename OutputContainer::iterator double_snap_rounding_2(InputIterator input_ std::set unique_point_set; std::map point_to_index; std::vector pts; - std::vector< std::vector< size_t> > polylines; + std::vector< std::vector< std::size_t> > polylines; // Transform range of the segments in the range of points and polyline of indexes for(typename std::vector::iterator it=segs.begin(); it!=segs.end(); ++it) @@ -312,8 +311,8 @@ typename OutputContainer::iterator double_snap_rounding_2(InputIterator input_ for(InputIterator it=segs.begin(); it!=segs.end(); ++it) { - size_t index1 = point_to_index[it->source()]; - size_t index2 = point_to_index[it->target()]; + std::size_t index1 = point_to_index[it->source()]; + std::size_t index2 = point_to_index[it->target()]; polylines.push_back({index1, index2}); } @@ -328,7 +327,7 @@ typename OutputContainer::iterator double_snap_rounding_2(InputIterator input_ output.clear(); for(auto &poly: polylines){ Polyline new_line; - for(size_t pi: poly) + for(std::size_t pi: poly) new_line.push_back(pts[pi]); output.push_back(new_line); } @@ -360,7 +359,7 @@ typename OutputContainer::iterator compute_snap_subcurves_2(InputIterator inpu std::set unique_point_set; std::map point_to_index; std::vector pts; - std::vector< std::vector< size_t> > polylines; + std::vector< std::vector< std::size_t> > polylines; for(typename std::vector::iterator it=segs.begin(); it!=segs.end(); ++it) { @@ -381,8 +380,8 @@ typename OutputContainer::iterator compute_snap_subcurves_2(InputIterator inpu for(InputIterator it=segs.begin(); it!=segs.end(); ++it) { - size_t index1 = point_to_index[it->source()]; - size_t index2 = point_to_index[it->target()]; + std::size_t index1 = point_to_index[it->source()]; + std::size_t index2 = point_to_index[it->target()]; polylines.push_back({index1, index2}); } @@ -394,10 +393,10 @@ typename OutputContainer::iterator compute_snap_subcurves_2(InputIterator inpu #endif // Output a range of segments while removing duplicate ones - std::set< std::pair > set_out_segs; + std::set< std::pair > set_out_segs; output.clear(); for(auto &poly: polylines){ - for(size_t i=1; i Date: Wed, 2 Apr 2025 09:44:18 +0200 Subject: [PATCH 11/41] add a reserve in the iteration of the snap --- Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h b/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h index a7cdbb9e420..501693affbf 100644 --- a/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h +++ b/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h @@ -146,9 +146,10 @@ void double_snap_rounding_2_disjoint(PointsRange &pts, PolylineRange &polylines) std::size_t size_before=pts.size(); CGAL::box_intersection_d(points_boxes.begin(), points_boxes.end(), segs_boxes.begin(), segs_boxes.end(), callback); points_boxes.clear(); + points_boxes.reserve(pts.size()-size_before); // The new vertices may intersect another segment when rounded, we repeat until they are not new vertices #ifdef DOUBLE_2D_SNAP_VERBOSE - std::cout << points_boxes.size()-size_before << " subdivisions performed" << std::endl; + std::cout << pts.size()-size_before << " subdivisions performed" << std::endl; #endif for(std::size_t i=size_before; i Date: Tue, 8 Apr 2025 17:45:03 +0200 Subject: [PATCH 12/41] Write float_snap_rounding_traits --- .../include/CGAL/Float_snap_rounding_2.h | 155 +++++++++++++----- .../test_float_snap_rounding_2.cpp | 1 - 2 files changed, 110 insertions(+), 46 deletions(-) diff --git a/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h b/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h index 501693affbf..2107d3bbab9 100644 --- a/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h +++ b/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h @@ -24,6 +24,59 @@ #include namespace CGAL { + +template +struct Float_snap_rounding_traits_2{ + typedef Kernel_ Kernel; + + typedef typename Kernel::FT FT; + typedef typename Kernel::Point_2 Point_2; + typedef typename Kernel::Segment_2 Segment_2; + typedef typename Kernel::Vector_2 Vector_2; + + typedef typename Kernel::Less_x_2 Less_x_2; + typedef typename Kernel::Less_y_2 Less_y_2; + typedef typename Kernel::Less_xy_2 Less_xy_2; + typedef typename Kernel::Less_yx_2 Less_yx_2; + typedef typename Kernel::Equal_2 Equal_2; + + typedef typename Kernel::Construct_source_2 Construct_source_2; + typedef typename Kernel::Construct_target_2 Construct_target_2; + + struct Squared_round_bound_2{ + double operator()(const FT &x){ + double b=std::nextafter(to_interval(x).second - to_interval(x).first, std::numeric_limits::infinity()); + return b*b; + } + double operator()(const Point_2 &p){ + return (*this)(p.x())+(*this)(p.y()); + } + }; + + Squared_round_bound_2 squared_round_bound_2_object(){ + return Squared_round_bound_2(); + } + + typedef typename Kernel::Compare_squared_distance_2 Compare_squared_distance_2; + Compare_squared_distance_2 compare_squared_distance_2_object(){ + return Kernel::squared_distance_2_object(); + } + + // typedef typename Base_kernel::Segment_2 Segment_2; + // typedef typename Base_kernel::Iso_rectangle_2 Iso_rectangle_2; + // typedef typename Base_kernel::Vector_2 Vector_2; + // typedef typename Base_kernel::Line_2 Line_2; + // typedef typename Base_kernel::Aff_transformation_2 Aff_transformation_2; + // typedef typename Base_kernel::Direction_2 Direction_2; + // typedef typename Base_kernel::Construct_vertex_2 Construct_vertex_2 ; + // typedef typename Base_kernel::Construct_segment_2 Construct_segment_2 ; + // typedef typename Base_kernel::Construct_iso_rectangle_2 Construct_iso_rectangle_2; + // typedef typename Base_kernel::Compare_y_2 Compare_y_2; + + + +}; + namespace Box_intersection_d { template @@ -46,35 +99,48 @@ public: } -template +template void double_snap_rounding_2_disjoint(PointsRange &pts, PolylineRange &polylines) { - using Point_2 = std::remove_cv_t::value_type>; - using Kernel = typename Kernel_traits::Kernel; + using Kernel = typename Traits::Kernel; + using Point_2 = typename Traits::Point_2; + using Segment_2 = typename Traits::Segment_2; + using Vector_2 = typename Traits::Vector_2; + using Comparison_result = typename Kernel::Comparison_result; - using Segment_2 = typename Kernel::Segment_2; - using Vector_2 = typename Kernel::Vector_2; + using PBox = CGAL::Box_intersection_d::Box_with_index_d; using SBox = CGAL::Box_intersection_d::Box_with_index_d; - auto comp_by_x=[&](std::size_t ai, std::size_t bi){ - return compare(pts[ai].x(),pts[bi].x())==SMALLER; + using Less_x_2 = typename Traits::Less_x_2; + using Less_y_2 = typename Traits::Less_y_2; + using Less_xy_2 = typename Traits::Less_xy_2; + using Less_yx_2 = typename Traits::Less_yx_2; + using Equal_2 = typename Traits::Equal_2; + + using Construct_source_2 = typename Traits::Construct_source_2; + using Construct_target_2 = typename Traits::Construct_target_2; + + using Compare_squared_distance_2 = typename Traits::Compare_squared_distance_2; + using Squared_round_bound_2 = typename Traits::Squared_round_bound_2; + + Compare_squared_distance_2 csq_dist_2; + Squared_round_bound_2 round_bound; + Construct_source_2 source; + Construct_target_2 target; + Equal_2 equal; + + auto Less_indexes_x_2=[&](size_t i, size_t j){ + return Less_x_2()(pts[i], pts[j]); }; - auto comp_by_y=[&](std::size_t ai, std::size_t bi){ - return compare(pts[ai].y(),pts[bi].y())==SMALLER; + auto Less_indexes_y_2=[&](size_t i, size_t j){ + return Less_y_2()(pts[i], pts[j]); }; - // Total order is required to using set - auto comp_by_x_first=[&](std::size_t ai, std::size_t bi){ - Comparison_result res=compare(pts[ai].x(),pts[bi].x()); - if(res==EQUAL) - return compare(pts[ai].y(),pts[bi].y())==SMALLER; - return res==SMALLER; + auto Less_indexes_xy_2=[&](size_t i, size_t j){ + return Less_xy_2()(pts[i], pts[j]); }; - auto comp_by_y_first=[&](std::size_t ai, std::size_t bi){ - Comparison_result res=compare(pts[ai].y(),pts[bi].y()); - if(res==EQUAL) - return compare(pts[ai].x(),pts[bi].x())==SMALLER; - return res==SMALLER; + auto Less_indexes_yx_2=[&](size_t i, size_t j){ + return Less_yx_2()(pts[i], pts[j]); }; #ifdef DOUBLE_2D_SNAP_VERBOSE @@ -84,10 +150,10 @@ void double_snap_rounding_2_disjoint(PointsRange &pts, PolylineRange &polylines) // Sorted the points may perform exact computations and thus refine the intervals of the coordinates values // This refine ensures that the order of the points will be preserved when rounded // However, except for this reason these sets are unused - using Iterator_set_x = typename std::set::iterator; - using Iterator_set_y = typename std::set::iterator; - std::set p_sort_by_x(comp_by_x_first); - std::set p_sort_by_y(comp_by_y_first); + using Iterator_set_x = typename std::set::iterator; + using Iterator_set_y = typename std::set::iterator; + std::set p_sort_by_x(Less_indexes_xy_2); + std::set p_sort_by_y(Less_indexes_yx_2); for(std::size_t i=0; i!=pts.size(); ++i) { p_sort_by_x.insert(i); @@ -100,20 +166,14 @@ void double_snap_rounding_2_disjoint(PointsRange &pts, PolylineRange &polylines) for(std::size_t i=0; i::infinity())-p.x().approx().inf(),2)+ - std::pow(std::nextafter(p.y().approx().sup(),std::numeric_limits::infinity())-p.y().approx().inf(),2); - }; + segs_boxes.emplace_back(pts[polylines[i][0]].bbox()+pts[polylines[i][polylines[i].size()-1]].bbox(),i); // Callback used for box_intersection_d auto callback=[&](PBox &bp, SBox &bseg){ std::size_t pi=bp.index(); std::size_t si=bseg.index(); std::size_t si1=polylines[bseg.index()][0]; - std::size_t si2=polylines[bseg.index()][1]; + std::size_t si2=polylines[bseg.index()][polylines[bseg.index()].size()-1]; if((pi==si1) || (pi==si2)) return; @@ -122,11 +182,14 @@ void double_snap_rounding_2_disjoint(PointsRange &pts, PolylineRange &polylines) Segment_2 seg(pts[polylines[bseg.index()][0]], pts[polylines[bseg.index()][1]]); // (A+B)^2 <= 4*max(A^2,B^2), we take some margin - double bound=16*(std::max)(round_bound(pts[pi]), (std::max)(round_bound(pts[si1]), round_bound(pts[si2]))); + double bound = round_bound(pts[pi]); + for(size_t i=0; i unique_points(p_sort_by_x.begin(),p_sort_by_x.end()); - std::sort(unique_points.begin(),unique_points.end(),comp_by_x_first); + std::sort(unique_points.begin(),unique_points.end(),Less_indexes_xy_2); std::vector new_pts; std::vector old_to_new_index(pts.size()); for(std::size_t i=0; i!=pts.size(); ++i){ - if(i==0 || (pts[unique_points[i]]!=pts[unique_points[i-1]])) + if(i==0 || !equal(pts[unique_points[i]],pts[unique_points[i-1]])) new_pts.push_back(pts[unique_points[i]]); old_to_new_index[unique_points[i]]=new_pts.size()-1; } @@ -232,7 +295,7 @@ void double_snap_rounding_2_disjoint(PointsRange &pts, PolylineRange &polylines) if(pts[poly[i-1]].x()==pts[poly[i]].x()){ Iterator_set_x start, end; // Get all vertices between the two endpoints along x order - if(comp_by_x_first(poly[i-1],poly[i])){ + if(Less_indexes_xy_2(poly[i-1],poly[i])){ start=p_sort_by_x.upper_bound(poly[i-1]); end=p_sort_by_x.lower_bound(poly[i]); } else { @@ -247,7 +310,7 @@ void double_snap_rounding_2_disjoint(PointsRange &pts, PolylineRange &polylines) if(pts[poly[i-1]].y()==pts[poly[i]].y()){ Iterator_set_y start, end; // Get all vertices between the two endpoints along y order - if(comp_by_y_first(poly[i-1],poly[i])){ + if(Less_indexes_yx_2(poly[i-1],poly[i])){ start=p_sort_by_y.upper_bound(poly[i-1]); end=p_sort_by_y.lower_bound(poly[i]); } else { @@ -277,6 +340,7 @@ typename OutputContainer::iterator double_snap_rounding_2(InputIterator input_ using Segment_2 = std::remove_cv_t::value_type>; using Polyline = std::remove_cv_t::value_type>; using Point_2 = typename Default_arr_traits::Traits::Point_2; + using Kernel = typename Kernel_traits::Kernel; std::vector segs(input_begin, input_end); #ifdef DOUBLE_2D_SNAP_VERBOSE @@ -318,7 +382,7 @@ typename OutputContainer::iterator double_snap_rounding_2(InputIterator input_ } // Main algorithm - double_snap_rounding_2_disjoint(pts, polylines); + double_snap_rounding_2_disjoint >(pts, polylines); #ifdef DOUBLE_2D_SNAP_VERBOSE std::cout << "Build output" << std::endl; @@ -341,11 +405,12 @@ TODO doc */ template typename OutputContainer::iterator compute_snap_subcurves_2(InputIterator input_begin, - InputIterator input_end, - OutputContainer& output) + InputIterator input_end, + OutputContainer& output) { using Segment_2 = std::remove_cv_t::value_type>; using Point_2 = typename Default_arr_traits::Traits::Point_2; + using Kernel = typename Kernel_traits::Kernel; std::vector segs; #ifdef DOUBLE_2D_SNAP_VERBOSE @@ -387,7 +452,7 @@ typename OutputContainer::iterator compute_snap_subcurves_2(InputIterator inpu } // Main algorithm - double_snap_rounding_2_disjoint(pts, polylines); + double_snap_rounding_2_disjoint >(pts, polylines); #ifdef DOUBLE_2D_SNAP_VERBOSE std::cout << "Build output" << std::endl; diff --git a/Snap_rounding_2/test/Snap_rounding_2/test_float_snap_rounding_2.cpp b/Snap_rounding_2/test/Snap_rounding_2/test_float_snap_rounding_2.cpp index a575eb07a9b..316205ab8df 100644 --- a/Snap_rounding_2/test/Snap_rounding_2/test_float_snap_rounding_2.cpp +++ b/Snap_rounding_2/test/Snap_rounding_2/test_float_snap_rounding_2.cpp @@ -78,7 +78,6 @@ int main(int argc,char *argv[]) CGAL::Random r(argc==1?rp.get_seed():std::stoi(argv[1])); std::cout << "random seed = " << r.get_seed() << std::endl; std::cout << std::setprecision(17); - test_data(); test_fully_random(r,2000); test_multi_almost_indentical_segments(r,200); return(0); From 4f89066f485144b2e63704fe415eb902728c1978 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?L=C3=A9o=20Valque?= Date: Tue, 8 Apr 2025 18:25:33 +0200 Subject: [PATCH 13/41] add round in the traits --- .../include/CGAL/Float_snap_rounding_2.h | 47 ++++++++++++------- .../test_float_snap_rounding_2.cpp | 4 +- 2 files changed, 33 insertions(+), 18 deletions(-) diff --git a/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h b/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h index 2107d3bbab9..3a54fe5cef8 100644 --- a/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h +++ b/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h @@ -53,6 +53,15 @@ struct Float_snap_rounding_traits_2{ } }; + struct Round_2{ + double operator()(const FT &x){ + return to_double(x); + } + Point_2 operator()(const Point_2 &p){ + return Point_2((*this)(p.x()),(*this)(p.y())); + } + }; + Squared_round_bound_2 squared_round_bound_2_object(){ return Squared_round_bound_2(); } @@ -109,6 +118,8 @@ void double_snap_rounding_2_disjoint(PointsRange &pts, PolylineRange &polylines) using Comparison_result = typename Kernel::Comparison_result; + using Polyline = std::remove_cv_t::value_type>; + using PBox = CGAL::Box_intersection_d::Box_with_index_d; using SBox = CGAL::Box_intersection_d::Box_with_index_d; @@ -123,12 +134,14 @@ void double_snap_rounding_2_disjoint(PointsRange &pts, PolylineRange &polylines) using Compare_squared_distance_2 = typename Traits::Compare_squared_distance_2; using Squared_round_bound_2 = typename Traits::Squared_round_bound_2; + using Round_2 = typename Traits::Round_2; Compare_squared_distance_2 csq_dist_2; Squared_round_bound_2 round_bound; Construct_source_2 source; Construct_target_2 target; Equal_2 equal; + Round_2 round; auto Less_indexes_x_2=[&](size_t i, size_t j){ return Less_x_2()(pts[i], pts[j]); @@ -172,19 +185,21 @@ void double_snap_rounding_2_disjoint(PointsRange &pts, PolylineRange &polylines) auto callback=[&](PBox &bp, SBox &bseg){ std::size_t pi=bp.index(); std::size_t si=bseg.index(); - std::size_t si1=polylines[bseg.index()][0]; - std::size_t si2=polylines[bseg.index()][polylines[bseg.index()].size()-1]; + + Polyline& pl=polylines[bseg.index()]; + std::size_t si1=pl[0]; + std::size_t si2=pl[pl.size()-1]; if((pi==si1) || (pi==si2)) return; Point_2& p= pts[bp.index()]; - Segment_2 seg(pts[polylines[bseg.index()][0]], pts[polylines[bseg.index()][1]]); + Segment_2 seg(pts[pl[0]], pts[pl[1]]); // (A+B)^2 <= 4*max(A^2,B^2), we take some margin double bound = round_bound(pts[pi]); - for(size_t i=0; i -typename OutputContainer::iterator compute_snap_subcurves_2(InputIterator input_begin, - InputIterator input_end, - OutputContainer& output) +typename OutputContainer::iterator compute_snapped_subcurves_2(InputIterator input_begin, + InputIterator input_end, + OutputContainer& output) { using Segment_2 = std::remove_cv_t::value_type>; using Point_2 = typename Default_arr_traits::Traits::Point_2; diff --git a/Snap_rounding_2/test/Snap_rounding_2/test_float_snap_rounding_2.cpp b/Snap_rounding_2/test/Snap_rounding_2/test_float_snap_rounding_2.cpp index 316205ab8df..e670767c1b0 100644 --- a/Snap_rounding_2/test/Snap_rounding_2/test_float_snap_rounding_2.cpp +++ b/Snap_rounding_2/test/Snap_rounding_2/test_float_snap_rounding_2.cpp @@ -36,7 +36,7 @@ void test_fully_random(CGAL::Random &r, size_t nb_segments){ CGAL::Real_timer t; t.start(); - CGAL::compute_snap_subcurves_2(segs.begin(), segs.end(), out); + CGAL::compute_snapped_subcurves_2(segs.begin(), segs.end(), out); t.stop(); std::cout << t.time() << "sec" << std::endl; @@ -56,7 +56,7 @@ void test_almost_indentical_segments(CGAL::Random &r, size_t nb_segments, Vector CGAL::Real_timer t; t.start(); - CGAL::compute_snap_subcurves_2(segs.begin(), segs.end(), out); + CGAL::compute_snapped_subcurves_2(segs.begin(), segs.end(), out); t.stop(); std::cout << t.time() << "sec" << std::endl; From 09684719e58b0499ebb8c53b73d97b1183a6d0c3 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?L=C3=A9o=20Valque?= Date: Fri, 29 Aug 2025 16:08:47 +0200 Subject: [PATCH 14/41] Improve running time of compute_snap_rounding_2 --- .../include/CGAL/Float_snap_rounding_2.h | 78 +++++++++++++++---- 1 file changed, 61 insertions(+), 17 deletions(-) diff --git a/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h b/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h index 3a54fe5cef8..20b31594443 100644 --- a/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h +++ b/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h @@ -23,6 +23,10 @@ #include #include +typedef CGAL::Exact_predicates_inexact_constructions_kernel Epeck; + +#include + namespace CGAL { template @@ -33,6 +37,7 @@ struct Float_snap_rounding_traits_2{ typedef typename Kernel::Point_2 Point_2; typedef typename Kernel::Segment_2 Segment_2; typedef typename Kernel::Vector_2 Vector_2; + typedef typename Kernel::Line_2 Line_2; typedef typename Kernel::Less_x_2 Less_x_2; typedef typename Kernel::Less_y_2 Less_y_2; @@ -108,13 +113,14 @@ public: } -template +template void double_snap_rounding_2_disjoint(PointsRange &pts, PolylineRange &polylines) { using Kernel = typename Traits::Kernel; using Point_2 = typename Traits::Point_2; using Segment_2 = typename Traits::Segment_2; using Vector_2 = typename Traits::Vector_2; + using Line_2 = typename Traits::Line_2; using Comparison_result = typename Kernel::Comparison_result; @@ -157,6 +163,7 @@ void double_snap_rounding_2_disjoint(PointsRange &pts, PolylineRange &polylines) }; #ifdef DOUBLE_2D_SNAP_VERBOSE + std::cout << "Nb of points: " << pts.size() << " , nb of polylines: " << polylines.size() << std::endl; std::cout << "Sort the input points" << std::endl; #endif // Compute the order of the points along the 2 axis @@ -173,14 +180,22 @@ void double_snap_rounding_2_disjoint(PointsRange &pts, PolylineRange &polylines) p_sort_by_y.insert(i); } + const double max_coordinate=std::max(std::max(to_double(pts[*p_sort_by_x.begin()].x()), to_double(pts[*(--p_sort_by_x.end())].x())), + std::max(to_double(pts[*p_sort_by_y.begin()].y()), to_double(pts[*(--p_sort_by_y.end())].y()))); + const double global_bound=max_coordinate*std::pow(2, -20); + //Prepare boxes for box_intersection_d std::vector points_boxes; std::vector segs_boxes; - for(std::size_t i=0; i round_bound_pts; + for(std::size_t i=0; i(points_boxes.begin(), points_boxes.end(), segs_boxes.begin(), segs_boxes.end(), callback); points_boxes.clear(); points_boxes.reserve(pts.size()-size_before); // The new vertices may intersect another segment when rounded, we repeat until they are not new vertices #ifdef DOUBLE_2D_SNAP_VERBOSE + std::cout << nb_calls << std::endl; + std::cout << nb_tests << std::endl; + std::cout << nb_execute << std::endl; + std::cout << nb_already_exists << std::endl; std::cout << pts.size()-size_before << " subdivisions performed" << std::endl; #endif for(std::size_t i=size_before; i +template typename OutputContainer::iterator double_snap_rounding_2(InputIterator input_begin, InputIterator input_end, OutputContainer& output) @@ -397,7 +433,7 @@ typename OutputContainer::iterator double_snap_rounding_2(InputIterator input_ } // Main algorithm - double_snap_rounding_2_disjoint >(pts, polylines); + double_snap_rounding_2_disjoint >(pts, polylines); #ifdef DOUBLE_2D_SNAP_VERBOSE std::cout << "Build output" << std::endl; @@ -415,10 +451,17 @@ typename OutputContainer::iterator double_snap_rounding_2(InputIterator input_ return output.begin(); } -/* -TODO doc +/** +* ingroup +* +* Given a range of segments, compute rounded subsegments that are pairwise disjoint in their interior, as induced by the input curves. +* +* @tparam Concurrency_tag That template parameter enables to choose whether the algorithm is to be run in +* parallel, if CGAL::Parallel_tag is specified and CGAL has been linked with the Intel TBB library, or sequentially, if CGAL::Sequential_tag - the default value - is specified. +* @tparam InputIterator iterator of a segment range +* @tparam OutputContainer inserter of a segment range */ -template +template typename OutputContainer::iterator compute_snapped_subcurves_2(InputIterator input_begin, InputIterator input_end, OutputContainer& output) @@ -427,11 +470,10 @@ typename OutputContainer::iterator compute_snapped_subcurves_2(InputIterator i using Point_2 = typename Default_arr_traits::Traits::Point_2; using Kernel = typename Kernel_traits::Kernel; - std::vector segs; + std::vector segs(input_begin, input_end); #ifdef DOUBLE_2D_SNAP_VERBOSE std::cout << "Solved intersections" << std::endl; #endif - compute_subcurves(input_begin, input_end, std::back_inserter(segs)); #ifdef DOUBLE_2D_SNAP_VERBOSE @@ -442,6 +484,7 @@ typename OutputContainer::iterator compute_snapped_subcurves_2(InputIterator i std::vector pts; std::vector< std::vector< std::size_t> > polylines; + // Transform range of the segments in the range of points and polyline of indexes for(typename std::vector::iterator it=segs.begin(); it!=segs.end(); ++it) { const Point_2& p1 = it->source(); @@ -466,8 +509,9 @@ typename OutputContainer::iterator compute_snapped_subcurves_2(InputIterator i polylines.push_back({index1, index2}); } + // Main algorithm - double_snap_rounding_2_disjoint >(pts, polylines); + double_snap_rounding_2_disjoint >(pts, polylines); #ifdef DOUBLE_2D_SNAP_VERBOSE std::cout << "Build output" << std::endl; From 9f7f0002b01b75d14f4c3994abb492536c2c7ed1 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?L=C3=A9o=20Valque?= Date: Fri, 29 Aug 2025 16:10:01 +0200 Subject: [PATCH 15/41] float_snap_rounding example --- .../examples/Snap_rounding_2/CMakeLists.txt | 8 +++ .../Snap_rounding_2/float_snap_rounding.cpp | 70 ++++++++++++------- 2 files changed, 54 insertions(+), 24 deletions(-) diff --git a/Snap_rounding_2/examples/Snap_rounding_2/CMakeLists.txt b/Snap_rounding_2/examples/Snap_rounding_2/CMakeLists.txt index e09fa21f5d5..315389de4a2 100644 --- a/Snap_rounding_2/examples/Snap_rounding_2/CMakeLists.txt +++ b/Snap_rounding_2/examples/Snap_rounding_2/CMakeLists.txt @@ -14,3 +14,11 @@ file( foreach(cppfile ${cppfiles}) create_single_source_cgal_program("${cppfile}") endforeach() + +find_package(TBB QUIET) +include(CGAL_TBB_support) +if(TARGET CGAL::TBB_support) + target_link_libraries(float_snap_rounding PRIVATE CGAL::TBB_support) +else() + message(STATUS "NOTICE: Intel TBB was not found. Sequential code will be used.") +endif() diff --git a/Snap_rounding_2/examples/Snap_rounding_2/float_snap_rounding.cpp b/Snap_rounding_2/examples/Snap_rounding_2/float_snap_rounding.cpp index 21ff49bc960..1783d0b570f 100644 --- a/Snap_rounding_2/examples/Snap_rounding_2/float_snap_rounding.cpp +++ b/Snap_rounding_2/examples/Snap_rounding_2/float_snap_rounding.cpp @@ -1,46 +1,68 @@ -// #define DOUBLE_2D_SNAP_VERBOSE - #include +#include #include #include +#include typedef CGAL::Exact_predicates_exact_constructions_kernel Kernel; +typedef CGAL::Exact_predicates_inexact_constructions_kernel Epick; typedef CGAL::Arr_segment_traits_2 Traits_2; typedef Traits_2::Curve_2 Segment_2; typedef Kernel::Point_2 Point_2; typedef Kernel::FT FT; typedef std::vector Polyline_2; -int main() +int main(int argc, char *argv[]) { - std::vector pts; std::vector< Segment_2 > segs; - FT e(std::pow(2, -60)); + if(argc>1){ + std::cout << "Read segments in " << argv[1] << std::endl; + std::ifstream in(argv[1]); - segs.emplace_back(Point_2(1-e, 1), Point_2(-1-e, -1+2*e)); - segs.emplace_back(Point_2(e/2, e/2), Point_2(1, -1)); - segs.emplace_back(Point_2(0, 2-e/2), Point_2(2, 0)); - segs.emplace_back(Point_2(0, 2-e/2), Point_2(-2+e, -4)); - segs.emplace_back(Point_2(-2, 2), Point_2(2, 2)); - segs.emplace_back(Point_2(7, 7), Point_2(7+e, 7+e)); - segs.emplace_back(Point_2(5, 7-e), Point_2(9, 7-e)); + int n; + in >> n; + segs.reserve(n); + for (int i=0; i> x1 >> y1 >> x2 >> y2; + if(Point_2(x1,y1)!=Point_2(x2,y2)) + segs.emplace_back(Point_2(x1, y1), Point_2(x2, y2)); + } + } else { + //Example with a non-trivial rounding + FT e(std::pow(2, -60)); - std::cout << "Input" << std::endl; - for(auto &seg: segs){ - std::cout << seg.source() << " -> " << seg.target() << std::endl; + segs.emplace_back(Point_2(1-e, 1), Point_2(-1-e, -1+2*e)); + segs.emplace_back(Point_2(e/2, e/2), Point_2(1, -1)); + segs.emplace_back(Point_2(0, 2-e/2), Point_2(2, 0)); + segs.emplace_back(Point_2(0, 2-e/2), Point_2(-2+e, -4)); + segs.emplace_back(Point_2(-2, 2), Point_2(2, 2)); + segs.emplace_back(Point_2(7, 7), Point_2(7+e, 7+e)); + segs.emplace_back(Point_2(5, 7-e), Point_2(9, 7-e)); } - std::cout << "\n\n" << std::endl; - std::vector< Polyline_2 > out; - CGAL::double_snap_rounding_2(segs.begin(), segs.end(), out); + std::cout << "Computes the intersections and snaps the segments" << std::endl; + std::vector< Segment_2> out; + CGAL::compute_snapped_subcurves_2(segs.begin(), segs.end(), out); + std::cout << "Does the output intersect: " << CGAL::do_curves_intersect(out.begin(), out.end()) << std::endl; + std::cout << "Size of the output: " << out.size() << std::endl; - std::cout << "Output" << std::endl; - for(auto &poly: out){ - std::cout << poly[0]; - for(size_t i=1; i!=poly.size(); ++i) - std::cout << " -> " << poly[i]; - std::cout << std::endl; + std::string out_path=(argc>2)?argv[2]:"out.segs"; + std::cout << "Write the outputs in " << out_path << std::endl; + std::ofstream outf(out_path); + CGAL::Random r; + + outf << std::setprecision(17); + outf << out.size() << std::endl; + for (auto &seg: out) + { + double x1 = CGAL::to_double(seg.source().x()); + double y1 = CGAL::to_double(seg.source().y()); + double x2 = CGAL::to_double(seg.target().x()); + double y2 = CGAL::to_double(seg.target().y()); + outf << x1 << " " << y1 << " " << x2 << " " << y2 << std::endl; } return 0; From 39d3ade9cd29abb4c5ee7ffd7251ce3f4f781179 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?L=C3=A9o=20Valque?= Date: Fri, 10 Oct 2025 17:59:53 +0200 Subject: [PATCH 16/41] rewrite test of Float_snap_rounding_2 --- .../include/CGAL/Arr_segment_traits_2.h | 5 + .../Snap_rounding_2/float_snap_rounding.cpp | 13 +- .../include/CGAL/Float_snap_rounding_2.h | 202 ++++++++++++++---- .../test_float_snap_rounding_2.cpp | 194 +++++++++++++++-- 4 files changed, 349 insertions(+), 65 deletions(-) diff --git a/Arrangement_on_surface_2/include/CGAL/Arr_segment_traits_2.h b/Arrangement_on_surface_2/include/CGAL/Arr_segment_traits_2.h index 0e2d4819b1b..8de90e1f4ef 100644 --- a/Arrangement_on_surface_2/include/CGAL/Arr_segment_traits_2.h +++ b/Arrangement_on_surface_2/include/CGAL/Arr_segment_traits_2.h @@ -728,6 +728,11 @@ public: // Check if we have a single intersection point. const Point_2* ip = std::get_if(&*res); if (ip != nullptr) { + if(!(cv1.is_vertical() ? + m_traits.is_in_y_range_2_object()(cv1, *ip) : + m_traits.is_in_x_range_2_object()(cv1, *ip))){ + std::cout << cv1 << std::endl; + } CGAL_assertion(cv1.is_vertical() ? m_traits.is_in_y_range_2_object()(cv1, *ip) : m_traits.is_in_x_range_2_object()(cv1, *ip)); diff --git a/Snap_rounding_2/examples/Snap_rounding_2/float_snap_rounding.cpp b/Snap_rounding_2/examples/Snap_rounding_2/float_snap_rounding.cpp index 1783d0b570f..36eba6c5371 100644 --- a/Snap_rounding_2/examples/Snap_rounding_2/float_snap_rounding.cpp +++ b/Snap_rounding_2/examples/Snap_rounding_2/float_snap_rounding.cpp @@ -1,6 +1,6 @@ #include #include -#include +// #include #include #include @@ -8,9 +8,9 @@ typedef CGAL::Exact_predicates_exact_constructions_kernel Kernel; typedef CGAL::Exact_predicates_inexact_constructions_kernel Epick; typedef CGAL::Arr_segment_traits_2 Traits_2; typedef Traits_2::Curve_2 Segment_2; -typedef Kernel::Point_2 Point_2; -typedef Kernel::FT FT; -typedef std::vector Polyline_2; +typedef Kernel::Point_2 Point_2; +typedef Kernel::FT FT; +typedef std::vector Polyline_2; int main(int argc, char *argv[]) { @@ -45,8 +45,9 @@ int main(int argc, char *argv[]) std::cout << "Computes the intersections and snaps the segments" << std::endl; std::vector< Segment_2> out; - CGAL::compute_snapped_subcurves_2(segs.begin(), segs.end(), out); - std::cout << "Does the output intersect: " << CGAL::do_curves_intersect(out.begin(), out.end()) << std::endl; + CGAL::Kernel_traits::value_type>>::Kernel::toto toto; + // CGAL::compute_snapped_subcurves_2(segs.begin(), segs.end(), out); + // std::cout << "Does the output intersect: " << CGAL::do_curves_intersect(out.begin(), out.end()) << std::endl; std::cout << "Size of the output: " << out.size() << std::endl; std::string out_path=(argc>2)?argv[2]:"out.segs"; diff --git a/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h b/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h index 20b31594443..db38ad48ae3 100644 --- a/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h +++ b/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h @@ -18,26 +18,34 @@ #endif #include #include +#include +#include #include +#include +#include +#include #include #include #include -typedef CGAL::Exact_predicates_inexact_constructions_kernel Epeck; +using Epeck=CGAL::Exact_predicates_exact_constructions_kernel; #include namespace CGAL { -template +template struct Float_snap_rounding_traits_2{ - typedef Kernel_ Kernel; + typedef Exact_Kernel Kernel; + + typedef Arr_segment_traits_2 Traits_2; typedef typename Kernel::FT FT; typedef typename Kernel::Point_2 Point_2; - typedef typename Kernel::Segment_2 Segment_2; + // typedef typename Kernel::Segment_2 Segment_2; + typedef typename Traits_2::Segment_2 Segment_2; typedef typename Kernel::Vector_2 Vector_2; - typedef typename Kernel::Line_2 Line_2; + typedef typename Kernel::Line_2 Line_2; typedef typename Kernel::Less_x_2 Less_x_2; typedef typename Kernel::Less_y_2 Less_y_2; @@ -48,32 +56,51 @@ struct Float_snap_rounding_traits_2{ typedef typename Kernel::Construct_source_2 Construct_source_2; typedef typename Kernel::Construct_target_2 Construct_target_2; + typedef Cartesian_converter Converter_in; + typedef Cartesian_converter Converter_out; + struct Squared_round_bound_2{ - double operator()(const FT &x){ + double operator()(const FT &x) const{ double b=std::nextafter(to_interval(x).second - to_interval(x).first, std::numeric_limits::infinity()); return b*b; } - double operator()(const Point_2 &p){ + double operator()(const Point_2 &p) const{ return (*this)(p.x())+(*this)(p.y()); } }; struct Round_2{ - double operator()(const FT &x){ + double operator()(const FT &x) const{ return to_double(x); } - Point_2 operator()(const Point_2 &p){ + Point_2 operator()(const Point_2 &p) const{ return Point_2((*this)(p.x()),(*this)(p.y())); } }; - Squared_round_bound_2 squared_round_bound_2_object(){ + Converter_in converter_to_exact_object() const{ + return Converter_in(); + } + + Converter_out converter_from_exact_object() const{ + return Converter_out(); + } + + Construct_source_2 construct_source_2_object() const{ + return Construct_source_2(); + } + + Construct_target_2 construct_target_2_object() const{ + return Construct_target_2(); + } + + Squared_round_bound_2 squared_round_bound_2_object() const{ return Squared_round_bound_2(); } typedef typename Kernel::Compare_squared_distance_2 Compare_squared_distance_2; - Compare_squared_distance_2 compare_squared_distance_2_object(){ - return Kernel::squared_distance_2_object(); + Compare_squared_distance_2 compare_squared_distance_2_object() const{ + return Kernel().compare_squared_distance_2_object(); } // typedef typename Base_kernel::Segment_2 Segment_2; @@ -114,7 +141,7 @@ public: } template -void double_snap_rounding_2_disjoint(PointsRange &pts, PolylineRange &polylines) +void double_snap_rounding_2_disjoint(PointsRange &pts, PolylineRange &polylines, const Traits &traits) { using Kernel = typename Traits::Kernel; using Point_2 = typename Traits::Point_2; @@ -142,10 +169,10 @@ void double_snap_rounding_2_disjoint(PointsRange &pts, PolylineRange &polylines) using Squared_round_bound_2 = typename Traits::Squared_round_bound_2; using Round_2 = typename Traits::Round_2; - Compare_squared_distance_2 csq_dist_2; - Squared_round_bound_2 round_bound; - Construct_source_2 source; - Construct_target_2 target; + Compare_squared_distance_2 csq_dist_2 = traits.compare_squared_distance_2_object(); + Squared_round_bound_2 round_bound = traits.squared_round_bound_2_object(); + Construct_source_2 source = traits.construct_source_2_object(); + Construct_target_2 target = traits.construct_target_2_object(); Equal_2 equal; Round_2 round; @@ -376,12 +403,14 @@ void double_snap_rounding_2_disjoint(PointsRange &pts, PolylineRange &polylines) /** * ingroup * -* Given a range of segments, compute rounded polylines that are pairwise disjoint in their interior, as induced by the input curves. +* Given a range of segments, compute rounded subsegments that are pairwise disjoint in their interiors, based on the input curves. +* The output is a range of polyline with each polyline corresponding to an input segment. * -* @tparam Concurrency_tag enables to choose whether the algorithm is to be run in -* parallel, if CGAL::Parallel_tag is specified and CGAL has been linked with the Intel TBB library, or sequentially, if CGAL::Sequential_tag - the default value - is specified. -* @tparam InputIterator iterator of a segment range -* @tparam OutputContainer inserter of a segment range +* @tparam Concurrency_tag allows choosing whether the algorithm runs in parallel or sequentially. +* If `CGAL::Parallel_tag` is specified and CGAL is linked with the Intel TBB library, the algorithm will run in parallel. +* Otherwise, if `CGAL::Sequential_tag` is specified (the default), the algorithm will run sequentially. +* @tparam InputIterator iterator over a range of `Segment_2` +* @tparam OutputContainer inserter over a range of `Polyline`. `Polyline` must be a type that provides a `push_back(Point_2)` function. */ template typename OutputContainer::iterator double_snap_rounding_2(InputIterator input_begin, @@ -460,21 +489,30 @@ typename OutputContainer::iterator double_snap_rounding_2(InputIterator input_ * parallel, if CGAL::Parallel_tag is specified and CGAL has been linked with the Intel TBB library, or sequentially, if CGAL::Sequential_tag - the default value - is specified. * @tparam InputIterator iterator of a segment range * @tparam OutputContainer inserter of a segment range +* @tparam The exact kernel needed for computation (Epeck by default) */ -template -typename OutputContainer::iterator compute_snapped_subcurves_2(InputIterator input_begin, - InputIterator input_end, - OutputContainer& output) +template ::value_type>>::Kernel> > +typename OutputContainer::iterator compute_snapped_subcurves_2(InputIterator input_begin, + InputIterator input_end, + OutputContainer& output, + const Traits& traits=Traits()) { - using Segment_2 = std::remove_cv_t::value_type>; - using Point_2 = typename Default_arr_traits::Traits::Point_2; - using Kernel = typename Kernel_traits::Kernel; + using Point_2 = typename Traits::Point_2; + using Segment_2 = typename Traits::Segment_2; + using I2E = typename Traits::Converter_in; + using E2O = typename Traits::Converter_out; + using VectorIterator = typename std::vector::iterator; + I2E converter_to_exact=traits.converter_to_exact_object(); + E2O converter_from_exact=traits.converter_from_exact_object(); - std::vector segs(input_begin, input_end); + std::vector convert_input; + for(InputIterator it=input_begin; it!=input_end; ++it) + convert_input.push_back(Segment_2(converter_to_exact(*it))); + std::vector segs; #ifdef DOUBLE_2D_SNAP_VERBOSE std::cout << "Solved intersections" << std::endl; #endif - compute_subcurves(input_begin, input_end, std::back_inserter(segs)); + compute_subcurves(convert_input.begin(), convert_input.end(), std::back_inserter(segs)); #ifdef DOUBLE_2D_SNAP_VERBOSE std::cout << "Change format to range of points and indexes" << std::endl; @@ -485,7 +523,7 @@ typename OutputContainer::iterator compute_snapped_subcurves_2(InputIterator i std::vector< std::vector< std::size_t> > polylines; // Transform range of the segments in the range of points and polyline of indexes - for(typename std::vector::iterator it=segs.begin(); it!=segs.end(); ++it) + for(VectorIterator it=segs.begin(); it!=segs.end(); ++it) { const Point_2& p1 = it->source(); const Point_2& p2 = it->target(); @@ -502,7 +540,7 @@ typename OutputContainer::iterator compute_snapped_subcurves_2(InputIterator i } } - for(InputIterator it=segs.begin(); it!=segs.end(); ++it) + for(VectorIterator it=segs.begin(); it!=segs.end(); ++it) { std::size_t index1 = point_to_index[it->source()]; std::size_t index2 = point_to_index[it->target()]; @@ -511,7 +549,7 @@ typename OutputContainer::iterator compute_snapped_subcurves_2(InputIterator i // Main algorithm - double_snap_rounding_2_disjoint >(pts, polylines); + double_snap_rounding_2_disjoint(pts, polylines, traits); #ifdef DOUBLE_2D_SNAP_VERBOSE std::cout << "Build output" << std::endl; @@ -524,8 +562,100 @@ typename OutputContainer::iterator compute_snapped_subcurves_2(InputIterator i for(std::size_t i=1; i::value_type>>::Kernel> > +typename OutputContainer::iterator snap_polygons_2(InputIterator input_begin, + InputIterator input_end, + OutputContainer& output, + const Traits& traits=Traits()) +{ + using Point_2 = typename Traits::Point_2; + using Segment_2 = typename Traits::Segment_2; + using I2E = typename Traits::Converter_in; + using E2O = typename Traits::Converter_out; + using VectorIterator = typename std::vector::iterator; + I2E converter_to_exact=traits.converter_to_exact_object(); + E2O converter_from_exact=traits.converter_from_exact_object(); + + std::vector convert_input; + for(InputIterator it=input_begin; it!=input_end; ++it) + convert_input.push_back(Segment_2(converter_to_exact(*it))); + std::vector segs; +#ifdef DOUBLE_2D_SNAP_VERBOSE + std::cout << "Solved intersections" << std::endl; +#endif + compute_subcurves(convert_input.begin(), convert_input.end(), std::back_inserter(segs)); + +#ifdef DOUBLE_2D_SNAP_VERBOSE + std::cout << "Change format to range of points and indexes" << std::endl; +#endif + std::set unique_point_set; + std::map point_to_index; + std::vector pts; + std::vector< std::vector< std::size_t> > polylines; + + // Transform range of the segments in the range of points and polyline of indexes + for(VectorIterator it=segs.begin(); it!=segs.end(); ++it) + { + const Point_2& p1 = it->source(); + const Point_2& p2 = it->target(); + + if (unique_point_set.find(p1) == unique_point_set.end()) { + unique_point_set.insert(p1); + pts.push_back(p1); + point_to_index[p1] = pts.size() - 1; + } + if (unique_point_set.find(p2) == unique_point_set.end()) { + unique_point_set.insert(p2); + pts.push_back(p2); + point_to_index[p2] = pts.size() - 1; + } + } + + for(VectorIterator it=segs.begin(); it!=segs.end(); ++it) + { + std::size_t index1 = point_to_index[it->source()]; + std::size_t index2 = point_to_index[it->target()]; + polylines.push_back({index1, index2}); + } + + + // Main algorithm + double_snap_rounding_2_disjoint(pts, polylines, traits); + +#ifdef DOUBLE_2D_SNAP_VERBOSE + std::cout << "Build output" << std::endl; +#endif + + // Output a range of segments while removing duplicate ones + std::set< std::pair > set_out_segs; + output.clear(); + for(auto &poly: polylines){ + for(std::size_t i=1; i +#include #include #include #include +#include + +#include + #include #include typedef CGAL::Exact_predicates_exact_constructions_kernel Kernel; -typedef CGAL::Arr_segment_traits_2 Traits_2; -typedef Traits_2::Curve_2 Segment_2; +typedef CGAL::Float_snap_rounding_traits_2 Traits_2; +typedef Kernel::Segment_2 Segment_2; typedef Kernel::Point_2 Point_2; typedef Kernel::Vector_2 Vector_2; typedef std::vector Polyline_2; +typedef std::vector Polyline_range_2; typedef Kernel::FT FT; +typedef CGAL::Arr_segment_traits_2 Arr_segment_traits_2; +typedef Arr_segment_traits_2::Curve_2 Curve_2; + +typedef CGAL::Snap_rounding_traits_2 SnapTraits; + +typedef CGAL::Polygon_2 Polygon_2; +typedef CGAL::Polygon_with_holes_2 Polygon_with_holes_2; +typedef std::vector Pwh_vec_2; + +typedef CGAL::Cartesian Naive; +typedef CGAL::Cartesian_converter EK_to_IK; +typedef CGAL::Cartesian_converter IK_to_EK; + +//Biggest double with ulp smaller than an integer +constexpr double maxFloat = std::pow(2,23); +constexpr double maxDouble = std::pow(2,52); + +Segment_2 scale_segment(const Segment_2 &s){ + auto scale_value=[](FT x){ + return maxDouble * CGAL::to_double(x); + }; + return Segment_2(Point_2(scale_value(s.start().x()), scale_value(s.start().y())), + Point_2(scale_value(s.end().x()), scale_value(s.end().y()))); +} + Point_2 random_point(CGAL::Random& r, double m = -1, double M = 1) { return Point_2(FT(r.get_double(m, M)), FT(r.get_double(m, M))); @@ -26,41 +57,135 @@ Point_2 random_noise_point(CGAL::Random& r, double m, double M, Vector_2 v) return Point_2(FT(r.get_double(m, M) + CGAL::to_double(v.x())), FT(r.get_double(m, M) + CGAL::to_double(v.y()))); } +template +void compute_subcurves_and_naive_snap(Iterator begin, Iterator end, OutRange &out){ + EK_to_IK to_inexact; + IK_to_EK to_exact; + std::vector segs; + for(Iterator it=begin; it!=end; ++it) + segs.emplace_back(*it); + CGAL::compute_subcurves(segs.begin(), segs.end(), std::back_inserter(out)); + //Naive snap + OutRange out2; + for(Segment_2 &seg: out){ + Segment_2 s = to_exact(to_inexact(seg)); + if(!s.is_degenerate()) + out2.emplace_back(s); + } + std::swap(out,out2); +} + +void test(const std::vector &segs){ + std::vector out; +#ifdef BENCH_AND_VERBOSE_FLOAT_SNAP_ROUNDING_2 + CGAL::Real_timer t; + std::cout << "Input size: " << segs.size() << std::endl; + t.start(); + compute_subcurves_and_naive_snap(segs.begin(), segs.end(), out); + t.stop(); + std::cout << "Naive snap size: " << out.size() << " ,running time: " << t.time() << " and output do intersect: " << CGAL::do_curves_intersect(out.begin(), out.end()) < segs; - std::vector out; - for(size_t i=0; i!=nb_segments; ++i) segs.emplace_back(random_point(r), random_point(r)); - CGAL::Real_timer t; - t.start(); - CGAL::compute_snapped_subcurves_2(segs.begin(), segs.end(), out); - t.stop(); - - std::cout << t.time() << "sec" << std::endl; - assert(!CGAL::do_curves_intersect(out.begin(), out.end())); + test(segs); + // CGAL::Real_timer t; + // Polyline_range_2 output_list; + // t.start(); + // CGAL::snap_rounding_2(segs.begin(), segs.end(), output_list, 1./maxFloat); + // t.stop(); + // std::cout << "Running time with integer 2D snap: " << t.time() << "sec" << std::endl; + // t.reset(); } void test_almost_indentical_segments(CGAL::Random &r, size_t nb_segments, Vector_2 source, Vector_2 target){ // Difficult test std::cout << "Test with almost identical segments from " << source << " to " << target << std::endl; std::vector segs; - std::vector out; - - double eps=std::pow(2,-45); + double eps=std::pow(2,-40); for(size_t i=0; i!=nb_segments; ++i){ segs.emplace_back(random_noise_point(r, -eps, eps, source), random_noise_point(r, -eps, eps, target)); } - CGAL::Real_timer t; - t.start(); - CGAL::compute_snapped_subcurves_2(segs.begin(), segs.end(), out); - t.stop(); + test(segs); +} + +void test_iterative_square_intersection(CGAL::Random &r, size_t nb_iterations){ + auto add_random_rotated_square=[&](std::vector &segs){ + double theta=r.get_double(0, CGAL_PI/2); + double cos_t = std::cos(theta); + double sin_t = std::sin(theta); + Point_2 a( cos_t, sin_t); + Point_2 b( sin_t,-cos_t); + Point_2 c(-cos_t,-sin_t); + Point_2 d(-sin_t, cos_t); + segs.emplace_back(a,b); + segs.emplace_back(b,c); + segs.emplace_back(c,d); + segs.emplace_back(d,a); + }; + + std::vector segs; + std::vector out; + std::vector arr_segs; + + CGAL::Real_timer t; + for(size_t i=0; i segs; + FT e(std::pow(2, -60)); + segs.emplace_back(Point_2(0, 1+e), Point_2(2, 1)); + segs.emplace_back(Point_2(1, 1), Point_2(1, 0)); + + test(segs); + segs.clear(); + + segs.emplace_back(Point_2(1-e, 1), Point_2(-1-e, -1+2*e)); + segs.emplace_back(Point_2(e/2, e/2), Point_2(1, -1)); + segs.emplace_back(Point_2(0, 2-e/2), Point_2(2, 0)); + segs.emplace_back(Point_2(0, 2-e/2), Point_2(-2+e, -4)); + segs.emplace_back(Point_2(-2, 2), Point_2(2, 2)); + segs.emplace_back(Point_2(7, 7), Point_2(7+e, 7+e)); + segs.emplace_back(Point_2(5, 7-e), Point_2(9, 7-e)); + test(segs); +} + int main(int argc,char *argv[]) { CGAL::Random rp; CGAL::Random r(argc==1?rp.get_seed():std::stoi(argv[1])); std::cout << "random seed = " << r.get_seed() << std::endl; std::cout << std::setprecision(17); - test_fully_random(r,2000); - test_multi_almost_indentical_segments(r,200); + fix_test(); + test_fully_random(r,1000); + test_multi_almost_indentical_segments(r,100); + test_iterative_square_intersection(r,500); return(0); } From 64922a63bbb24a3f25906d3838fe41b6ffe5dbd1 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?L=C3=A9o=20Valque?= Date: Wed, 15 Oct 2025 17:10:29 +0200 Subject: [PATCH 17/41] Clean Float_2D_snap --- .../include/CGAL/Float_snap_rounding_2.h | 140 +++++++----------- .../test_float_snap_rounding_2.cpp | 73 +++++---- 2 files changed, 101 insertions(+), 112 deletions(-) diff --git a/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h b/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h index db38ad48ae3..dc2ff684d8a 100644 --- a/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h +++ b/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h @@ -16,6 +16,7 @@ #ifdef DOUBLE_2D_SNAP_VERBOSE #include #endif + #include #include #include @@ -28,37 +29,33 @@ #include #include -using Epeck=CGAL::Exact_predicates_exact_constructions_kernel; - #include namespace CGAL { -template -struct Float_snap_rounding_traits_2{ - typedef Exact_Kernel Kernel; +template +struct Float_snap_rounding_traits_2: Arr_segment_traits_2{ + typedef Arr_segment_traits_2 Base; - typedef Arr_segment_traits_2 Traits_2; + typedef typename Base::FT FT; + typedef typename Base::Point_2 Point_2; + typedef typename Base::Segment_2 Segment_2; + typedef typename Base::Vector_2 Vector_2; + typedef typename Base::Line_2 Line_2; - typedef typename Kernel::FT FT; - typedef typename Kernel::Point_2 Point_2; - // typedef typename Kernel::Segment_2 Segment_2; - typedef typename Traits_2::Segment_2 Segment_2; - typedef typename Kernel::Vector_2 Vector_2; - typedef typename Kernel::Line_2 Line_2; + typedef typename Base::Less_x_2 Less_x_2; + typedef typename Base::Less_y_2 Less_y_2; + typedef typename Base::Less_xy_2 Less_xy_2; + typedef typename Base::Less_yx_2 Less_yx_2; + typedef typename Base::Equal_2 Equal_2; - typedef typename Kernel::Less_x_2 Less_x_2; - typedef typename Kernel::Less_y_2 Less_y_2; - typedef typename Kernel::Less_xy_2 Less_xy_2; - typedef typename Kernel::Less_yx_2 Less_yx_2; - typedef typename Kernel::Equal_2 Equal_2; - - typedef typename Kernel::Construct_source_2 Construct_source_2; - typedef typename Kernel::Construct_target_2 Construct_target_2; + typedef typename Base::Construct_source_2 Construct_source_2; + typedef typename Base::Construct_target_2 Construct_target_2; typedef Cartesian_converter Converter_in; typedef Cartesian_converter Converter_out; + // Return an upperbound of the squared distance between a point and its rounded value struct Squared_round_bound_2{ double operator()(const FT &x) const{ double b=std::nextafter(to_interval(x).second - to_interval(x).first, std::numeric_limits::infinity()); @@ -86,40 +83,18 @@ struct Float_snap_rounding_traits_2{ return Converter_out(); } - Construct_source_2 construct_source_2_object() const{ - return Construct_source_2(); - } - - Construct_target_2 construct_target_2_object() const{ - return Construct_target_2(); - } - Squared_round_bound_2 squared_round_bound_2_object() const{ return Squared_round_bound_2(); } - typedef typename Kernel::Compare_squared_distance_2 Compare_squared_distance_2; - Compare_squared_distance_2 compare_squared_distance_2_object() const{ - return Kernel().compare_squared_distance_2_object(); + Round_2 round_2_object() const{ + return Round_2(); } - - // typedef typename Base_kernel::Segment_2 Segment_2; - // typedef typename Base_kernel::Iso_rectangle_2 Iso_rectangle_2; - // typedef typename Base_kernel::Vector_2 Vector_2; - // typedef typename Base_kernel::Line_2 Line_2; - // typedef typename Base_kernel::Aff_transformation_2 Aff_transformation_2; - // typedef typename Base_kernel::Direction_2 Direction_2; - // typedef typename Base_kernel::Construct_vertex_2 Construct_vertex_2 ; - // typedef typename Base_kernel::Construct_segment_2 Construct_segment_2 ; - // typedef typename Base_kernel::Construct_iso_rectangle_2 Construct_iso_rectangle_2; - // typedef typename Base_kernel::Compare_y_2 Compare_y_2; - - - }; namespace Box_intersection_d { +// Since we used only box_self_intersection_d, we may not have intersection of two distinct boxes with same index template class Box_with_index_d: public Box_d< NT_, N, ID_EXPLICIT> { protected: @@ -143,49 +118,48 @@ public: template void double_snap_rounding_2_disjoint(PointsRange &pts, PolylineRange &polylines, const Traits &traits) { - using Kernel = typename Traits::Kernel; - using Point_2 = typename Traits::Point_2; + using Point_2 = typename Traits::Point_2; using Segment_2 = typename Traits::Segment_2; using Vector_2 = typename Traits::Vector_2; - using Line_2 = typename Traits::Line_2; - - using Comparison_result = typename Kernel::Comparison_result; + using Line_2 = typename Traits::Line_2; using Polyline = std::remove_cv_t::value_type>; - using PBox = CGAL::Box_intersection_d::Box_with_index_d; - using SBox = CGAL::Box_intersection_d::Box_with_index_d; + using PointsRangeIterator = typename PointsRange::iterator; + using PolylineRangeIterator = typename PolylineRange::iterator; - using Less_x_2 = typename Traits::Less_x_2; - using Less_y_2 = typename Traits::Less_y_2; + using PBox = CGAL::Box_intersection_d::Box_with_info_d; + using SBox = CGAL::Box_intersection_d::Box_with_handle_d; + + using Less_x_2 = typename Traits::Less_x_2; + using Less_y_2 = typename Traits::Less_y_2; using Less_xy_2 = typename Traits::Less_xy_2; using Less_yx_2 = typename Traits::Less_yx_2; - using Equal_2 = typename Traits::Equal_2; - - using Construct_source_2 = typename Traits::Construct_source_2; - using Construct_target_2 = typename Traits::Construct_target_2; + using Equal_2 = typename Traits::Equal_2; + using Round_2 = typename Traits::Round_2; + using Construct_source_2 = typename Traits::Construct_source_2; + using Construct_target_2 = typename Traits::Construct_target_2; using Compare_squared_distance_2 = typename Traits::Compare_squared_distance_2; - using Squared_round_bound_2 = typename Traits::Squared_round_bound_2; - using Round_2 = typename Traits::Round_2; + using Squared_round_bound_2 = typename Traits::Squared_round_bound_2; Compare_squared_distance_2 csq_dist_2 = traits.compare_squared_distance_2_object(); Squared_round_bound_2 round_bound = traits.squared_round_bound_2_object(); Construct_source_2 source = traits.construct_source_2_object(); Construct_target_2 target = traits.construct_target_2_object(); - Equal_2 equal; - Round_2 round; + Equal_2 equal = traits.equal_2_object(); + Round_2 round = traits.round_2_object(); - auto Less_indexes_x_2=[&](size_t i, size_t j){ + auto Less_indexes_x_2=[&](std::size_t i, std::size_t j){ return Less_x_2()(pts[i], pts[j]); }; - auto Less_indexes_y_2=[&](size_t i, size_t j){ + auto Less_indexes_y_2=[&](std::size_t i, std::size_t j){ return Less_y_2()(pts[i], pts[j]); }; - auto Less_indexes_xy_2=[&](size_t i, size_t j){ + auto Less_indexes_xy_2=[&](std::size_t i, std::size_t j){ return Less_xy_2()(pts[i], pts[j]); }; - auto Less_indexes_yx_2=[&](size_t i, size_t j){ + auto Less_indexes_yx_2=[&](std::size_t i, std::size_t j){ return Less_yx_2()(pts[i], pts[j]); }; @@ -212,30 +186,34 @@ void double_snap_rounding_2_disjoint(PointsRange &pts, PolylineRange &polylines, const double global_bound=max_coordinate*std::pow(2, -20); //Prepare boxes for box_intersection_d + std::vector new_points; std::vector points_boxes; std::vector segs_boxes; + + //We create that vector to avoid multiple computations std::vector round_bound_pts; + for(std::size_t i=0; isize()-1]].bbox(),it); CGAL_MUTEX mutex_callback; // Callback used for box_intersection_d auto callback=[&](PBox &bp, SBox &bseg){ - std::size_t pi=bp.index(); - std::size_t si=bseg.index(); + std::size_t pi=bp.info(); - Polyline& pl=polylines[bseg.index()]; + Polyline& pl=*(bseg.handle()); std::size_t si1=pl[0]; std::size_t si2=pl[pl.size()-1]; + //Check if p is one endpoint of the segment if((pi==si1) || (pi==si2)) return; - Point_2& p= pts[bp.index()]; + Point_2& p= pts[pi]; // Early exit for better running time if(certainly(csq_dist_2(p, Line_2(pts[pl[0]], pts[pl[1]]), global_bound)==LARGER)) @@ -244,10 +222,8 @@ void double_snap_rounding_2_disjoint(PointsRange &pts, PolylineRange &polylines, Segment_2 seg(pts[pl[0]], pts[pl[1]]); // (A+B)^2 <= 4*max(A^2,B^2), we take some margin - // double bound = round_bound(pts[pi]); double bound=round_bound_pts[pi]; - for(size_t i=0; i(points_boxes.begin(), points_boxes.end(), segs_boxes.begin(), segs_boxes.end(), callback); points_boxes.clear(); points_boxes.reserve(pts.size()-size_before); - // The new vertices may intersect another segment when rounded, we repeat until they are not new vertices #ifdef DOUBLE_2D_SNAP_VERBOSE - std::cout << nb_calls << std::endl; - std::cout << nb_tests << std::endl; - std::cout << nb_execute << std::endl; - std::cout << nb_already_exists << std::endl; - std::cout << pts.size()-size_before << " subdivisions performed" << std::endl; + ++turn_nb; + std::cout << "Turn " << turn_nb << ": " << pts.size()-size_before << " subdivisions performed" << std::endl; #endif + // The new vertices may intersect another segment when rounded, we repeat until they are not new vertices for(std::size_t i=size_before; i #include @@ -91,9 +93,18 @@ void test(const std::vector &segs){ CGAL::compute_snapped_subcurves_2(segs.begin(), segs.end(), out); #ifdef BENCH_AND_VERBOSE_FLOAT_SNAP_ROUNDING_2 t.stop(); - std::cout << "Formal snap size: " << out.size() << " ,running time: " << t.time() << "\n" << std::endl; + std::cout << "Formal snap size: " << out.size() << " ,running time: " << t.time() << std::endl; #endif assert(!CGAL::do_curves_intersect(out.begin(), out.end())); +#ifdef COMPARE_WITH_INTEGER_SNAP_ROUNDING_2 + Polyline_range_2 output_list; + t.reset(); + t.start(); + CGAL::snap_rounding_2(segs.begin(), segs.end(), output_list, 1./maxDouble); + t.stop(); + std::cout << "Running time with integer 2D snap (scaled at 10^15): " << t.time() << std::endl; +#endif + std::cout << "\n"; } void test_fully_random(CGAL::Random &r, size_t nb_segments){ @@ -103,13 +114,6 @@ void test_fully_random(CGAL::Random &r, size_t nb_segments){ segs.emplace_back(random_point(r), random_point(r)); test(segs); - // CGAL::Real_timer t; - // Polyline_range_2 output_list; - // t.start(); - // CGAL::snap_rounding_2(segs.begin(), segs.end(), output_list, 1./maxFloat); - // t.stop(); - // std::cout << "Running time with integer 2D snap: " << t.time() << "sec" << std::endl; - // t.reset(); } void test_almost_indentical_segments(CGAL::Random &r, size_t nb_segments, Vector_2 source, Vector_2 target){ @@ -163,30 +167,30 @@ void test_iterative_square_intersection(CGAL::Random &r, size_t nb_iterations){ } -void test_iterative_square_intersection(CGAL::Random &r, size_t nb_iterations){ - auto random_rotated_square=[&](){ - double theta=r.get_double(0, CGAL_PI/2); - double cos_t = std::cos(theta); - double sin_t = std::sin(theta); - Polygon_2 P; - P.push_back( cos_t, sin_t); - P.push_back( sin_t,-cos_t); - P.push_back(-cos_t,-sin_t); - P.push_back(-sin_t, cos_t); - return P; - }; +// void test_iterative_square_intersection(CGAL::Random &r, size_t nb_iterations){ +// auto random_rotated_square=[&](){ +// double theta=r.get_double(0, CGAL_PI/2); +// double cos_t = std::cos(theta); +// double sin_t = std::sin(theta); +// Polygon_2 P; +// P.push_back( cos_t, sin_t); +// P.push_back( sin_t,-cos_t); +// P.push_back(-cos_t,-sin_t); +// P.push_back(-sin_t, cos_t); +// return P; +// }; - Pwh_vec_2 scene, out_intersection; +// Pwh_vec_2 scene, out_intersection; - CGAL::intersection(random_rotated_square, random_rotated_square, scene); +// CGAL::intersection(random_rotated_square, random_rotated_square, scene); - for(size_t i=0; i segs; + FT e(std::pow(2, -60)); + segs.emplace_back(Point_2(0, 0), Point_2(1, 1)); + segs.emplace_back(Point_2(0.5+e, 0.5), Point_2(1, -1)); + segs.emplace_back(Point_2(0.5-e, 0.5), Point_2(-1, 3)); + + test(segs); +} + int main(int argc,char *argv[]) { CGAL::Random rp; CGAL::Random r(argc==1?rp.get_seed():std::stoi(argv[1])); std::cout << "random seed = " << r.get_seed() << std::endl; std::cout << std::setprecision(17); + test_box_intersection(); fix_test(); - test_fully_random(r,1000); - test_multi_almost_indentical_segments(r,100); - test_iterative_square_intersection(r,500); + // test_fully_random(r,1000); + // test_multi_almost_indentical_segments(r,100); + // test_iterative_square_intersection(r,500); return(0); } From c380ad2fb663f287e09d64f34612a3bafc1da55c Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?L=C3=A9o=20Valque?= Date: Thu, 16 Oct 2025 15:14:42 +0200 Subject: [PATCH 18/41] Snap of polygons --- .../include/CGAL/Float_snap_rounding_2.h | 100 ++++++++--------- .../test_float_snap_rounding_2.cpp | 103 ++++++------------ 2 files changed, 79 insertions(+), 124 deletions(-) diff --git a/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h b/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h index dc2ff684d8a..e535fb0c6c7 100644 --- a/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h +++ b/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h @@ -476,12 +476,12 @@ typename OutputContainer::iterator compute_snapped_subcurves_2(InputIterator using I2E = typename Traits::Converter_in; using E2O = typename Traits::Converter_out; using VectorIterator = typename std::vector::iterator; - I2E converter_to_exact=traits.converter_to_exact_object(); - E2O converter_from_exact=traits.converter_from_exact_object(); + I2E to_exact=traits.converter_to_exact_object(); + E2O from_exact=traits.converter_from_exact_object(); std::vector convert_input; for(InputIterator it=input_begin; it!=input_end; ++it) - convert_input.push_back(Segment_2(converter_to_exact(*it))); + convert_input.push_back(Segment_2(to_exact(*it))); std::vector segs; #ifdef DOUBLE_2D_SNAP_VERBOSE std::cout << "Solved intersections" << std::endl; @@ -537,7 +537,7 @@ typename OutputContainer::iterator compute_snapped_subcurves_2(InputIterator set_out_segs.emplace((std::min)(poly[i-1],poly[i]),(std::max)(poly[i-1],poly[i])); } for(auto &pair: set_out_segs){ - output.emplace_back(converter_from_exact(pts[pair.first]), converter_from_exact(pts[pair.second])); + output.emplace_back(from_exact(pts[pair.first]), from_exact(pts[pair.second])); assert(pts[pair.first]!=pts[pair.second]); } @@ -547,71 +547,67 @@ typename OutputContainer::iterator compute_snapped_subcurves_2(InputIterator /** * ingroup * -* Given a range of segments, compute rounded subsegments that are pairwise disjoint in their interior, as induced by the input curves. +* Given a Polygon_2, compute rounded segments that are pairwise disjoint in their interior, as induced by the input polygon. +* The output is guarantee to be a Polygon but may present pinched section. * * @tparam Concurrency_tag That template parameter enables to choose whether the algorithm is to be run in * parallel, if CGAL::Parallel_tag is specified and CGAL has been linked with the Intel TBB library, or sequentially, if CGAL::Sequential_tag - the default value - is specified. * @tparam InputIterator iterator of a segment range * @tparam OutputContainer inserter of a segment range * @tparam The exact kernel needed for computation (Epeck by default) +* +* @warning The convex property is not necessarly preserved */ -template ::value_type>>::Kernel> > -typename OutputContainer::iterator snap_polygons_2(InputIterator input_begin, - InputIterator input_end, - OutputContainer& output, - const Traits& traits=Traits()) +template ::Kernel> > +void snap_polygons_2(const Polygon_2 &P, + Polygon_2 &out, + const Traits& traits=Traits(), + bool check_duplicates = false) { using Point_2 = typename Traits::Point_2; using Segment_2 = typename Traits::Segment_2; using I2E = typename Traits::Converter_in; using E2O = typename Traits::Converter_out; using VectorIterator = typename std::vector::iterator; - I2E converter_to_exact=traits.converter_to_exact_object(); - E2O converter_from_exact=traits.converter_from_exact_object(); - - std::vector convert_input; - for(InputIterator it=input_begin; it!=input_end; ++it) - convert_input.push_back(Segment_2(converter_to_exact(*it))); - std::vector segs; -#ifdef DOUBLE_2D_SNAP_VERBOSE - std::cout << "Solved intersections" << std::endl; -#endif - compute_subcurves(convert_input.begin(), convert_input.end(), std::back_inserter(segs)); + I2E to_exact=traits.converter_to_exact_object(); + E2O from_exact=traits.converter_from_exact_object(); #ifdef DOUBLE_2D_SNAP_VERBOSE std::cout << "Change format to range of points and indexes" << std::endl; #endif - std::set unique_point_set; - std::map point_to_index; std::vector pts; std::vector< std::vector< std::size_t> > polylines; - // Transform range of the segments in the range of points and polyline of indexes - for(VectorIterator it=segs.begin(); it!=segs.end(); ++it) - { - const Point_2& p1 = it->source(); - const Point_2& p2 = it->target(); + if(check_duplicates){ + std::set unique_point_set; + std::map point_to_index; - if (unique_point_set.find(p1) == unique_point_set.end()) { - unique_point_set.insert(p1); - pts.push_back(p1); - point_to_index[p1] = pts.size() - 1; + // Transform the polygon in a range of points and polylines of indexes + for(const typename Polygon_2::Point_2 &p_: P.vertices()) + { + Point_2 p=to_exact(p_); + if (unique_point_set.find(p) == unique_point_set.end()) { + unique_point_set.insert(p); + pts.push_back(p); + point_to_index[p] = pts.size() - 1; + } } - if (unique_point_set.find(p2) == unique_point_set.end()) { - unique_point_set.insert(p2); - pts.push_back(p2); - point_to_index[p2] = pts.size() - 1; + + for(const typename Polygon_2::Segment_2 &s: P.edges()){ + Point_2 p1=to_exact(s.source()); + Point_2 p2=to_exact(s.target()); + std::size_t index1 = point_to_index[p1]; + std::size_t index2 = point_to_index[p2]; + polylines.push_back({index1, index2}); } + } else { + for(const typename Polygon_2::Point_2 &p: P.vertices()) + pts.push_back(to_exact(p)); + for(size_t i=0; isource()]; - std::size_t index2 = point_to_index[it->target()]; - polylines.push_back({index1, index2}); - } - - // Main algorithm double_snap_rounding_2_disjoint(pts, polylines, traits); @@ -620,20 +616,14 @@ typename OutputContainer::iterator snap_polygons_2(InputIterator input_begin, #endif // Output a range of segments while removing duplicate ones - std::set< std::pair > set_out_segs; - output.clear(); - for(auto &poly: polylines){ + out.clear(); + for(auto &poly: polylines) for(std::size_t i=1; i #include @@ -129,69 +127,47 @@ void test_almost_indentical_segments(CGAL::Random &r, size_t nb_segments, Vector } void test_iterative_square_intersection(CGAL::Random &r, size_t nb_iterations){ - auto add_random_rotated_square=[&](std::vector &segs){ + auto random_rotated_square=[&](){ double theta=r.get_double(0, CGAL_PI/2); - double cos_t = std::cos(theta); - double sin_t = std::sin(theta); - Point_2 a( cos_t, sin_t); - Point_2 b( sin_t,-cos_t); - Point_2 c(-cos_t,-sin_t); - Point_2 d(-sin_t, cos_t); - segs.emplace_back(a,b); - segs.emplace_back(b,c); - segs.emplace_back(c,d); - segs.emplace_back(d,a); +#ifdef BENCH_AND_VERBOSE_FLOAT_SNAP_ROUNDING_2 + std::cout << "Angle: " << theta << std::endl; +#endif + FT cos_t(std::cos(theta)); + FT sin_t(std::sin(theta)); + Polygon_2 P; + P.push_back(Point_2( cos_t, sin_t)); + P.push_back(Point_2(-sin_t, cos_t)); + P.push_back(Point_2(-cos_t,-sin_t)); + P.push_back(Point_2( sin_t,-cos_t)); + return P; }; - std::vector segs; - std::vector out; - std::vector arr_segs; + Polygon_2 scene=random_rotated_square(); + Polygon_2 snap_scene; + Pwh_vec_2 out_intersection; - CGAL::Real_timer t; for(size_t i=0; i segs; - FT e(std::pow(2, -60)); - segs.emplace_back(Point_2(0, 0), Point_2(1, 1)); - segs.emplace_back(Point_2(0.5+e, 0.5), Point_2(1, -1)); - segs.emplace_back(Point_2(0.5-e, 0.5), Point_2(-1, 3)); - - test(segs); -} - int main(int argc,char *argv[]) { CGAL::Random rp; CGAL::Random r(argc==1?rp.get_seed():std::stoi(argv[1])); std::cout << "random seed = " << r.get_seed() << std::endl; std::cout << std::setprecision(17); - test_box_intersection(); fix_test(); // test_fully_random(r,1000); // test_multi_almost_indentical_segments(r,100); - // test_iterative_square_intersection(r,500); + test_iterative_square_intersection(r,2000); return(0); } From 443ed4801cd846fd45ace5051c5c5ce266fea347 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?L=C3=A9o=20Valque?= Date: Thu, 16 Oct 2025 18:43:54 +0200 Subject: [PATCH 19/41] NamedParameter --- .../internal/parameters_interface.h | 2 + .../include/CGAL/Float_snap_rounding_2.h | 186 +++++++++++++++--- .../test_float_snap_rounding_2.cpp | 22 ++- 3 files changed, 178 insertions(+), 32 deletions(-) diff --git a/STL_Extension/include/CGAL/STL_Extension/internal/parameters_interface.h b/STL_Extension/include/CGAL/STL_Extension/internal/parameters_interface.h index 73c54858bf2..f70fb981631 100644 --- a/STL_Extension/include/CGAL/STL_Extension/internal/parameters_interface.h +++ b/STL_Extension/include/CGAL/STL_Extension/internal/parameters_interface.h @@ -409,3 +409,5 @@ CGAL_add_named_parameter(maximum_height_t, maximum_height, maximum_height) // List of named parameters used in the package 'Constrained_triangulation_3' CGAL_add_named_parameter(plc_face_id_t, plc_face_id, plc_face_id) CGAL_add_named_parameter(with_plc_face_id_t, with_plc_face_id, with_plc_face_id) +//List of named parameters used in Snap_rounding_2 +CGAL_add_named_parameter(do_intersection_computation_t, do_intersection_computation, do_intersection_computation) diff --git a/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h b/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h index e535fb0c6c7..7d45f9436c0 100644 --- a/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h +++ b/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h @@ -29,6 +29,8 @@ #include #include +#include + #include namespace CGAL { @@ -386,15 +388,37 @@ void double_snap_rounding_2_disjoint(PointsRange &pts, PolylineRange &polylines, * @tparam InputIterator iterator over a range of `Segment_2` * @tparam OutputContainer inserter over a range of `Polyline`. `Polyline` must be a type that provides a `push_back(Point_2)` function. */ -template +template typename OutputContainer::iterator double_snap_rounding_2(InputIterator input_begin, InputIterator input_end, - OutputContainer& output) + OutputContainer& output, + const NamedParameters &np = parameters::default_values()) { - using Segment_2 = std::remove_cv_t::value_type>; + using Concurrency_tag = typename internal_np::Lookup_named_param_def::type; + + using InputKernel = typename Kernel_traits::value_type>>::Kernel; + using DefaultTraits = Float_snap_rounding_traits_2; + using Traits = typename internal_np::Lookup_named_param_def::type; + + using Point_2 = typename Traits::Point_2; + using Segment_2 = typename Traits::Segment_2; + using I2E = typename Traits::Converter_in; + using E2O = typename Traits::Converter_out; + using VectorIterator = typename std::vector::iterator; + using Polyline = std::remove_cv_t::value_type>; - using Point_2 = typename Default_arr_traits::Traits::Point_2; - using Kernel = typename Kernel_traits::Kernel; + + using parameters::choose_parameter; + using parameters::get_parameter; + + const Traits &traits = choose_parameter(get_parameter(np, internal_np::geom_traits), Traits()); + + I2E to_exact=traits.converter_to_exact_object(); + E2O from_exact=traits.converter_from_exact_object(); std::vector segs(input_begin, input_end); #ifdef DOUBLE_2D_SNAP_VERBOSE @@ -436,7 +460,7 @@ typename OutputContainer::iterator double_snap_rounding_2(InputIterator input_ } // Main algorithm - double_snap_rounding_2_disjoint >(pts, polylines); + double_snap_rounding_2_disjoint(pts, polylines, traits); #ifdef DOUBLE_2D_SNAP_VERBOSE std::cout << "Build output" << std::endl; @@ -465,17 +489,33 @@ typename OutputContainer::iterator double_snap_rounding_2(InputIterator input_ * @tparam OutputContainer inserter of a segment range * @tparam The exact kernel needed for computation (Epeck by default) */ -template ::value_type>>::Kernel> > +template typename OutputContainer::iterator compute_snapped_subcurves_2(InputIterator input_begin, InputIterator input_end, OutputContainer& output, - const Traits& traits=Traits()) + const NamedParameters &np = parameters::default_values()) { + using Concurrency_tag = typename internal_np::Lookup_named_param_def::type; + + using InputKernel = typename Kernel_traits::value_type>>::Kernel; + using DefaultTraits = Float_snap_rounding_traits_2; + using Traits = typename internal_np::Lookup_named_param_def::type; + using Point_2 = typename Traits::Point_2; using Segment_2 = typename Traits::Segment_2; using I2E = typename Traits::Converter_in; using E2O = typename Traits::Converter_out; using VectorIterator = typename std::vector::iterator; + + using parameters::choose_parameter; + using parameters::get_parameter; + + const Traits &traits = choose_parameter(get_parameter(np, internal_np::geom_traits), Traits()); + I2E to_exact=traits.converter_to_exact_object(); E2O from_exact=traits.converter_from_exact_object(); @@ -547,28 +587,62 @@ typename OutputContainer::iterator compute_snapped_subcurves_2(InputIterator /** * ingroup * -* Given a Polygon_2, compute rounded segments that are pairwise disjoint in their interior, as induced by the input polygon. -* The output is guarantee to be a Polygon but may present pinched section. +* Given a range of Polygon_2, compute rounded segments that are pairwise disjoint in their interior, as induced by the input polygons. +* Any polygon is guarantee to remain a Polygon in the output but may present pinched section and common vertices or segments with +* other polygons. * -* @tparam Concurrency_tag That template parameter enables to choose whether the algorithm is to be run in -* parallel, if CGAL::Parallel_tag is specified and CGAL has been linked with the Intel TBB library, or sequentially, if CGAL::Sequential_tag - the default value - is specified. -* @tparam InputIterator iterator of a segment range -* @tparam OutputContainer inserter of a segment range -* @tparam The exact kernel needed for computation (Epeck by default) +* @tparam InputIterator iterator of a CGAL::Polygon_2 range +* @tparam OutputContainer inserter of a CGAL::Polygon_2 range +* @tparam NamedParameters a sequence of \ref bgl_namedparameters "Named Parameters" +* \param begin,end the input polygon range +* \param out the output inserter +* \param np an optional sequence of \ref bgl_namedparameters "Named Parameters" among the ones listed below +* \cgalNamedParamsBegin +* \cgalParamNBegin{concurrency_tag} +* \cgalParamDescription{That template parameter enables to choose whether the algorithm is to be run in parallel, if CGAL::Parallel_tag +* is specified and CGAL has been linked with the Intel TBB library, or sequentially, otherwise.} +* \cgalParamType{CGAL::Concurrency_tag} +* \cgalParamDefault{CGAL::Sequential_tag} +* \cgalParamNEnd * -* @warning The convex property is not necessarly preserved +* \cgalParamNBegin{geom_traits} +* \cgalParamDescription{TODO} +* \cgalParamType{double} +* \cgalParamDefault{100} +* \cgalParamNEnd +* @warning The convex property of the polygons is not necessarly preserved */ -template ::Kernel> > -void snap_polygons_2(const Polygon_2 &P, - Polygon_2 &out, - const Traits& traits=Traits(), - bool check_duplicates = false) +template +void snap_polygons_2(InputIterator begin, + InputIterator end, + OutputContainer out, + const NamedParameters &np = parameters::default_values()) { + using Concurrency_tag = typename internal_np::Lookup_named_param_def::type; + + using Polygon_2 = typename std::iterator_traits::value_type; + using InputKernel = typename Kernel_traits::Kernel; + using DefaultTraits = Float_snap_rounding_traits_2; + using Traits = typename internal_np::Lookup_named_param_def::type; + using Point_2 = typename Traits::Point_2; using Segment_2 = typename Traits::Segment_2; using I2E = typename Traits::Converter_in; using E2O = typename Traits::Converter_out; using VectorIterator = typename std::vector::iterator; + + using parameters::choose_parameter; + using parameters::get_parameter; + + Polygon_2 P=*begin; + + const Traits &traits = choose_parameter(get_parameter(np, internal_np::geom_traits), Traits()); + const bool compute_intersections = choose_parameter(get_parameter(np, internal_np::do_intersection_computation), false); + I2E to_exact=traits.converter_to_exact_object(); E2O from_exact=traits.converter_from_exact_object(); @@ -577,8 +651,9 @@ void snap_polygons_2(const Polygon_2 &P, #endif std::vector pts; std::vector< std::vector< std::size_t> > polylines; + std::vector polygon_index; - if(check_duplicates){ + if(compute_intersections){ std::set unique_point_set; std::map point_to_index; @@ -601,11 +676,19 @@ void snap_polygons_2(const Polygon_2 &P, polylines.push_back({index1, index2}); } } else { - for(const typename Polygon_2::Point_2 &p: P.vertices()) - pts.push_back(to_exact(p)); - for(size_t i=0; ivertices()) + pts.push_back(to_exact(p)); + for(size_t i=0; i &poly = polylines[pl_ind]; + Polygon_2 P; + for(std::size_t i=1; i +void snap_polygon_2(const Polygon_2 &P, Polygon_2 &out, const NamedParameters &np = parameters::default_values()) +{ + std::array vec({P}); + std::vector out_vec; + snap_polygons_2(vec.begin(), vec.end(), std::back_inserter(out_vec), np); + out = out_vec[0]; } diff --git a/Snap_rounding_2/test/Snap_rounding_2/test_float_snap_rounding_2.cpp b/Snap_rounding_2/test/Snap_rounding_2/test_float_snap_rounding_2.cpp index b0f824a6ea8..c664e2f1029 100644 --- a/Snap_rounding_2/test/Snap_rounding_2/test_float_snap_rounding_2.cpp +++ b/Snap_rounding_2/test/Snap_rounding_2/test_float_snap_rounding_2.cpp @@ -1,3 +1,4 @@ +#define DOUBLE_2D_SNAP_VERBOSE #define BENCH_AND_VERBOSE_FLOAT_SNAP_ROUNDING_2 #include @@ -154,7 +155,7 @@ void test_iterative_square_intersection(CGAL::Random &r, size_t nb_iterations){ CGAL::Real_timer t; t.start(); #endif - snap_polygons_2(out_intersection[0].outer_boundary(), snap_scene); + snap_polygon_2(out_intersection[0].outer_boundary(), snap_scene); // snap_scene=out_intersection[0].outer_boundary(); #ifdef BENCH_AND_VERBOSE_FLOAT_SNAP_ROUNDING_2 t.stop(); @@ -198,6 +199,24 @@ void fix_test(){ test(segs); } +void test_float_snap_rounding(){ + std::cout << "Fix tests" << std::endl; + + std::vector< Segment_2 > segs; + std::vector< Polyline_2 > out; + FT e(std::pow(2, -60)); + segs.emplace_back(Point_2(1-e, 1), Point_2(-1-e, -1+2*e)); + segs.emplace_back(Point_2(e/2, e/2), Point_2(1, -1)); + segs.emplace_back(Point_2(0, 2-e/2), Point_2(2, 0)); + segs.emplace_back(Point_2(0, 2-e/2), Point_2(-2+e, -4)); + segs.emplace_back(Point_2(-2, 2), Point_2(2, 2)); + segs.emplace_back(Point_2(7, 7), Point_2(7+e, 7+e)); + segs.emplace_back(Point_2(5, 7-e), Point_2(9, 7-e)); + + double_snap_rounding_2(segs.begin(), segs.end(), out); + assert(out.size()==segs.size()); +} + int main(int argc,char *argv[]) { CGAL::Random rp; @@ -205,6 +224,7 @@ int main(int argc,char *argv[]) std::cout << "random seed = " << r.get_seed() << std::endl; std::cout << std::setprecision(17); fix_test(); + test_float_snap_rounding(); // test_fully_random(r,1000); // test_multi_almost_indentical_segments(r,100); test_iterative_square_intersection(r,2000); From 3df3151cbb5f916ca0eb53cc2ef576af4fb8214f Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?L=C3=A9o=20Valque?= Date: Thu, 16 Oct 2025 18:52:47 +0200 Subject: [PATCH 20/41] Rewrite float_snap_rounding example --- .../examples/Snap_rounding_2/float_snap_rounding.cpp | 10 ++++------ 1 file changed, 4 insertions(+), 6 deletions(-) diff --git a/Snap_rounding_2/examples/Snap_rounding_2/float_snap_rounding.cpp b/Snap_rounding_2/examples/Snap_rounding_2/float_snap_rounding.cpp index 36eba6c5371..0a7a58d5063 100644 --- a/Snap_rounding_2/examples/Snap_rounding_2/float_snap_rounding.cpp +++ b/Snap_rounding_2/examples/Snap_rounding_2/float_snap_rounding.cpp @@ -1,13 +1,12 @@ #include #include -// #include +#include #include #include typedef CGAL::Exact_predicates_exact_constructions_kernel Kernel; typedef CGAL::Exact_predicates_inexact_constructions_kernel Epick; -typedef CGAL::Arr_segment_traits_2 Traits_2; -typedef Traits_2::Curve_2 Segment_2; +typedef Kernel::Segment_2 Segment_2; typedef Kernel::Point_2 Point_2; typedef Kernel::FT FT; typedef std::vector Polyline_2; @@ -45,9 +44,8 @@ int main(int argc, char *argv[]) std::cout << "Computes the intersections and snaps the segments" << std::endl; std::vector< Segment_2> out; - CGAL::Kernel_traits::value_type>>::Kernel::toto toto; - // CGAL::compute_snapped_subcurves_2(segs.begin(), segs.end(), out); - // std::cout << "Does the output intersect: " << CGAL::do_curves_intersect(out.begin(), out.end()) << std::endl; + CGAL::compute_snapped_subcurves_2(segs.begin(), segs.end(), out); + std::cout << "Does the output intersect: " << CGAL::do_curves_intersect(out.begin(), out.end()) << std::endl; std::cout << "Size of the output: " << out.size() << std::endl; std::string out_path=(argc>2)?argv[2]:"out.segs"; From cfa69d01cf1e2110dd7cddf451a4602214d0efdc Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?L=C3=A9o=20Valque?= Date: Fri, 17 Oct 2025 14:44:39 +0200 Subject: [PATCH 21/41] Reintroduce flag NO_TESTING to test_snap_rounding_2 --- .../test/Snap_rounding_2/CMakeLists.txt | 18 ++---------------- 1 file changed, 2 insertions(+), 16 deletions(-) diff --git a/Snap_rounding_2/test/Snap_rounding_2/CMakeLists.txt b/Snap_rounding_2/test/Snap_rounding_2/CMakeLists.txt index 4a86b77fee0..7645a98dbd8 100644 --- a/Snap_rounding_2/test/Snap_rounding_2/CMakeLists.txt +++ b/Snap_rounding_2/test/Snap_rounding_2/CMakeLists.txt @@ -6,22 +6,8 @@ project(Snap_rounding_2_Tests) find_package(CGAL REQUIRED) -if(MSVC) - # Turn off a VC++ warning on a potential division by zero - # in Cartesian_kernel/include/CGAL/constructions/kernel_ftC3.h - # where CGAL_assume() does not help - set(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} /wd4723") -endif() - -# create a target per cppfile -file( - GLOB cppfiles - RELATIVE ${CMAKE_CURRENT_SOURCE_DIR} - ${CMAKE_CURRENT_SOURCE_DIR}/*.cpp) -foreach(cppfile ${cppfiles}) - create_single_source_cgal_program("${cppfile}") -endforeach() - +create_single_source_cgal_program(test_snap_rounding_2.cpp NO_TESTING) +create_single_source_cgal_program(test_float_snap_rounding_2.cpp) function(add_Snap_rounding_tests name) set(data_dir "data") From e8ea2811bae3334643d0e6b27bde81bd5879ece9 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?L=C3=A9o=20Valque?= Date: Thu, 23 Oct 2025 11:26:07 +0200 Subject: [PATCH 22/41] Use a one way scan instead of box_intersection_d --- .../include/CGAL/Float_snap_rounding_2.h | 398 +++++++++++++++++- .../test_float_snap_rounding_2.cpp | 70 ++- 2 files changed, 454 insertions(+), 14 deletions(-) diff --git a/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h b/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h index 7d45f9436c0..c6b5a9414ff 100644 --- a/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h +++ b/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h @@ -117,9 +117,374 @@ public: } +template +void scan(PointsRange &pts, PolylineRange &polylines, const Traits &traits){ + using Point_2 = typename Traits::Point_2; + using Segment_2 = typename Traits::Segment_2; + using Line_2 = typename Traits::Line_2; + + using Polyline = std::remove_cv_t::value_type>; + + using PointsRangeIterator = typename PointsRange::iterator; + using PolylineRangeIterator = typename PolylineRange::iterator; + + using PBox = CGAL::Box_intersection_d::Box_with_info_d; + using SBox = CGAL::Box_intersection_d::Box_with_handle_d; + + using Less_x_2 = typename Traits::Less_x_2; + using Less_y_2 = typename Traits::Less_y_2; + using Less_xy_2 = typename Traits::Less_xy_2; + using Less_yx_2 = typename Traits::Less_yx_2; + using Equal_2 = typename Traits::Equal_2; + using Round_2 = typename Traits::Round_2; + + using Construct_source_2 = typename Traits::Construct_source_2; + using Construct_target_2 = typename Traits::Construct_target_2; + using Compare_squared_distance_2 = typename Traits::Compare_squared_distance_2; + using Squared_round_bound_2 = typename Traits::Squared_round_bound_2; + + Compare_squared_distance_2 csq_dist_2 = traits.compare_squared_distance_2_object(); + Squared_round_bound_2 round_bound = traits.squared_round_bound_2_object(); + Construct_source_2 source = traits.construct_source_2_object(); + Construct_target_2 target = traits.construct_target_2_object(); + Equal_2 equal = traits.equal_2_object(); + Round_2 round = traits.round_2_object(); + + // Precompute round bounds for points + std::vector round_bound_pts; + round_bound_pts.reserve(pts.size()); + for(std::size_t i=0; i b.type; + // We compare with y only to have a complete order + return Less_y_2()(pa,pb); + }; + + auto pi_below_li=[&](size_t li, std::pair pair_pi_li){ + std::size_t ind_pi_support = pair_pi_li.second; + std::size_t pi = pair_pi_li.first; + if(ind_pi_support==li) return false; + + const Polyline &pl=polylines[li]; + Orientation ori=orientation(pts[pl.front()], pts[pl.back()], pts[pi]); + // When collinear, pick an other vertex of the line of pi to decide + if(ori==COLLINEAR){ + const Polyline &pi_support = polylines[ind_pi_support]; + std::size_t other_point = (pi_support.front()!=pi) ? pi_support.front() : pi_support.back(); + ori=orientation(pts[pl.front()], pts[pl.back()], pts[other_point]); + assert(ori!=COLLINEAR); + } + return ori==RIGHT_TURN; + }; + +#ifdef DOUBLE_2D_SNAP_FULL_VERBOSE + std::cout << "Input points" << std::endl; + std::size_t i=0; + for(const Point_2 &p: pts) + std::cout << i++ << ": " << p << std::endl; + std::cout << std::endl; + std::cout << "Input polylines" << std::endl; + i=0; + for(const Polyline &pl: polylines){ + std::cout << i++ << ":"; + for(std::size_t pi: pl) + std::cout << " " << pi; + std::cout << std::endl; + } + std::cout << std::endl; +#endif + +#ifdef DOUBLE_2D_SNAP_VERBOSE + std::cout << "Create the Event Queue" << std::endl; +#endif + + std::vector event_queue; + event_queue.reserve(2*polylines.size()); + for(size_t li=0; li y_order; + for(Event event: event_queue){ + // Find the position of the event along the y direction + std::vector::iterator pos_it = y_order.begin(); + if(!y_order.empty()) + pos_it = std::lower_bound(y_order.begin(), y_order.end(), std::pair(event.pi, event.li), pi_below_li); +#ifdef DOUBLE_2D_SNAP_FULL_VERBOSE + std::cout << ((event.type==EVENT_TYPE::INSERT)?"Insert": "Remove") << " Event at point " << event.pi << " on line "<< event.li << std::endl; +#endif + if(event.type==EVENT_TYPE::REMOVE){ + assert(*pos_it==event.li); + y_order.erase(pos_it); + } + + if(!y_order.empty()){ + auto close_event=[&](std::size_t pi, std::size_t li){ + // std::cout << pi << " " << li << std::endl; + if((pi==polylines[li].front()) || (pi==polylines[li].back())) + return false; + + const Point_2 &p = pts[pi]; + Polyline &pl=polylines[li]; + Segment_2 seg(pts[pl.front()], pts[pl.back()]); + + // (A+B)^2 <= 4*max(A^2,B^2), we take some margin + double bound=round_bound_pts[pi]; + for(std::size_t i=0; i indices(pts.size(),0); + std::iota(indices.begin(),indices.end(),0); + std::sort(indices.begin(),indices.end(),sort_pi); + +#ifdef DOUBLE_2D_SNAP_FULL_VERBOSE + std::cout << "Output points" << std::endl; + i=0; + for(const Point_2 &p: pts) + std::cout << i++ << ": " << p << std::endl; + std::cout << std::endl; + std::cout << "Output polylines" << std::endl; + i=0; + for(const Polyline &pl: polylines){ + std::cout << i++ << ":"; + for(std::size_t pi: pl) + std::cout << " " << pi; + std::cout << std::endl; + } + std::cout << std::endl; +#endif +} + +template +void round_and_post_process(PointsRange &pts, PolylineRange &polylines, const Traits &traits) +{ + using Point_2 = typename Traits::Point_2; + using Segment_2 = typename Traits::Segment_2; + using Vector_2 = typename Traits::Vector_2; + using Line_2 = typename Traits::Line_2; + + using Polyline = std::remove_cv_t::value_type>; + + using PointsRangeIterator = typename PointsRange::iterator; + using PolylineRangeIterator = typename PolylineRange::iterator; + + using PBox = CGAL::Box_intersection_d::Box_with_info_d; + using SBox = CGAL::Box_intersection_d::Box_with_handle_d; + + using Less_x_2 = typename Traits::Less_x_2; + using Less_y_2 = typename Traits::Less_y_2; + using Less_xy_2 = typename Traits::Less_xy_2; + using Less_yx_2 = typename Traits::Less_yx_2; + using Equal_2 = typename Traits::Equal_2; + using Round_2 = typename Traits::Round_2; + + using Construct_source_2 = typename Traits::Construct_source_2; + using Construct_target_2 = typename Traits::Construct_target_2; + using Compare_squared_distance_2 = typename Traits::Compare_squared_distance_2; + using Squared_round_bound_2 = typename Traits::Squared_round_bound_2; + + Compare_squared_distance_2 csq_dist_2 = traits.compare_squared_distance_2_object(); + Squared_round_bound_2 round_bound = traits.squared_round_bound_2_object(); + Construct_source_2 source = traits.construct_source_2_object(); + Construct_target_2 target = traits.construct_target_2_object(); + Equal_2 equal = traits.equal_2_object(); + Round_2 round = traits.round_2_object(); + + auto Less_indexes_x_2=[&](std::size_t i, std::size_t j){ + return Less_x_2()(pts[i], pts[j]); + }; + auto Less_indexes_y_2=[&](std::size_t i, std::size_t j){ + return Less_y_2()(pts[i], pts[j]); + }; + auto Less_indexes_xy_2=[&](std::size_t i, std::size_t j){ + return Less_xy_2()(pts[i], pts[j]); + }; + auto Less_indexes_yx_2=[&](std::size_t i, std::size_t j){ + return Less_yx_2()(pts[i], pts[j]); + }; + +#ifdef DOUBLE_2D_SNAP_VERBOSE + std::cout << "Round" << std::endl; +#endif + for(auto &p: pts) + p=round(p); + +#ifdef DOUBLE_2D_SNAP_VERBOSE + std::cout << "Remove duplicate points from collapsing" << std::endl; +#endif + + std::vector< std::size_t > unique_points(pts.size()); + std::iota(unique_points.begin(), unique_points.end(), 0); + std::sort(unique_points.begin(), unique_points.end(), Less_indexes_xy_2); + std::vector new_pts; + std::vector old_to_new_index(pts.size()); + for(std::size_t i=0; i!=pts.size(); ++i){ + if(i==0 || !equal(pts[unique_points[i]],pts[unique_points[i-1]])) + new_pts.push_back(pts[unique_points[i]]); + old_to_new_index[unique_points[i]]=new_pts.size()-1; + } + + std::swap(pts, new_pts); + size_t i=0; + for (Polyline& polyline : polylines) { + std::vector updated_polyline; + for (std::size_t i=0; i p_sort_by_x(Less_indexes_xy_2); + for(std::size_t i=0; i!=pts.size(); ++i) + p_sort_by_x.insert(i); + + using Iterator_set_x = typename std::set::iterator; + + i=0; + for(Polyline &poly: polylines){ + std::vector updated_polyline; + updated_polyline.push_back(poly.front()); + for(std::size_t i=1; i!=poly.size(); ++i){ + if(pts[poly[i-1]].x()==pts[poly[i]].x()){ + // updated_polyline.clear(); + // std::swap(poly, updated_polyline); + // break; + Iterator_set_x start, end; + // Get all vertices between the two endpoints along x order + if(Less_indexes_xy_2(poly[i-1],poly[i])){ + start=p_sort_by_x.upper_bound(poly[i-1]); + end=p_sort_by_x.lower_bound(poly[i]); + for(auto it=start; it!=end; ++it){ + updated_polyline.push_back(*it); + } + } else { + start=p_sort_by_x.upper_bound(poly[i]); + end=p_sort_by_x.lower_bound(poly[i-1]); + for(auto it=end; it!=start;){ + --it; + updated_polyline.push_back(*it); + } + } + } + updated_polyline.push_back(poly[i]); + } + std::swap(poly, updated_polyline); + ++i; + } +} + template void double_snap_rounding_2_disjoint(PointsRange &pts, PolylineRange &polylines, const Traits &traits) { + scan(pts, polylines, traits); + round_and_post_process(pts, polylines, traits); + return; + using Point_2 = typename Traits::Point_2; using Segment_2 = typename Traits::Segment_2; using Vector_2 = typename Traits::Vector_2; @@ -412,6 +777,8 @@ typename OutputContainer::iterator double_snap_rounding_2(InputIterator input_ using Polyline = std::remove_cv_t::value_type>; + using Less_xy_2 = typename Traits::Less_xy_2; + using parameters::choose_parameter; using parameters::get_parameter; @@ -420,11 +787,15 @@ typename OutputContainer::iterator double_snap_rounding_2(InputIterator input_ I2E to_exact=traits.converter_to_exact_object(); E2O from_exact=traits.converter_from_exact_object(); - std::vector segs(input_begin, input_end); + std::vector convert_input; + for(InputIterator it=input_begin; it!=input_end; ++it) + convert_input.push_back(Segment_2(to_exact(*it))); + std::vector segs; #ifdef DOUBLE_2D_SNAP_VERBOSE std::cout << "Solved intersections" << std::endl; + std::cout << "do intersect? " << do_curves_intersect(convert_input.begin(), convert_input.end()) << std::endl; #endif - compute_subcurves(input_begin, input_end, std::back_inserter(segs)); + compute_subcurves(convert_input.begin(), convert_input.end(), std::back_inserter(segs)); #ifdef DOUBLE_2D_SNAP_VERBOSE std::cout << "Change format to range of points and indexes" << std::endl; @@ -435,7 +806,7 @@ typename OutputContainer::iterator double_snap_rounding_2(InputIterator input_ std::vector< std::vector< std::size_t> > polylines; // Transform range of the segments in the range of points and polyline of indexes - for(typename std::vector::iterator it=segs.begin(); it!=segs.end(); ++it) + for(VectorIterator it=segs.begin(); it!=segs.end(); ++it) { const Point_2& p1 = it->source(); const Point_2& p2 = it->target(); @@ -452,11 +823,14 @@ typename OutputContainer::iterator double_snap_rounding_2(InputIterator input_ } } - for(InputIterator it=segs.begin(); it!=segs.end(); ++it) + for(VectorIterator it=segs.begin(); it!=segs.end(); ++it) { std::size_t index1 = point_to_index[it->source()]; std::size_t index2 = point_to_index[it->target()]; - polylines.push_back({index1, index2}); + if(Less_xy_2()(it->source(), it->target())) + polylines.push_back({index1, index2}); + else + polylines.push_back({index2, index1}); } // Main algorithm @@ -511,6 +885,8 @@ typename OutputContainer::iterator compute_snapped_subcurves_2(InputIterator using E2O = typename Traits::Converter_out; using VectorIterator = typename std::vector::iterator; + using Less_xy_2 = typename Traits::Less_xy_2; + using parameters::choose_parameter; using parameters::get_parameter; @@ -525,6 +901,7 @@ typename OutputContainer::iterator compute_snapped_subcurves_2(InputIterator std::vector segs; #ifdef DOUBLE_2D_SNAP_VERBOSE std::cout << "Solved intersections" << std::endl; + std::cout << "do intersect? " << do_curves_intersect(convert_input.begin(), convert_input.end()) << std::endl; #endif compute_subcurves(convert_input.begin(), convert_input.end(), std::back_inserter(segs)); @@ -558,7 +935,10 @@ typename OutputContainer::iterator compute_snapped_subcurves_2(InputIterator { std::size_t index1 = point_to_index[it->source()]; std::size_t index2 = point_to_index[it->target()]; - polylines.push_back({index1, index2}); + if(Less_xy_2()(it->source(), it->target())) + polylines.push_back({index1, index2}); + else + polylines.push_back({index2, index1}); } @@ -572,9 +952,15 @@ typename OutputContainer::iterator compute_snapped_subcurves_2(InputIterator // Output a range of segments while removing duplicate ones std::set< std::pair > set_out_segs; output.clear(); + // for(auto &poly: polylines){ + // for(std::size_t i=1; i @@ -90,6 +91,9 @@ void test(const std::vector &segs){ t.start(); #endif CGAL::compute_snapped_subcurves_2(segs.begin(), segs.end(), out); + std::vector out2(out.begin(), out.end()); + // out.clear(); + // CGAL::compute_snapped_subcurves_2(out2.begin(), out2.end(), out); #ifdef BENCH_AND_VERBOSE_FLOAT_SNAP_ROUNDING_2 t.stop(); std::cout << "Formal snap size: " << out.size() << " ,running time: " << t.time() << std::endl; @@ -180,15 +184,35 @@ void test_multi_almost_indentical_segments(CGAL::Random &r, size_t nb_segments){ void fix_test(){ std::cout << "Fix tests" << std::endl; - std::vector< Segment_2 > segs; FT e(std::pow(2, -60)); - segs.emplace_back(Point_2(0, 1+e), Point_2(2, 1)); - segs.emplace_back(Point_2(1, 1), Point_2(1, 0)); + segs.emplace_back(Point_2(0, 1), Point_2(1, 1)); + segs.emplace_back(Point_2(1, 2), Point_2(2, 2)); test(segs); segs.clear(); + // segs.emplace_back(Point_2(0, 1), Point_2(2, 1)); + // segs.emplace_back(Point_2(1, 0), Point_2(1, 0)); + // test(segs); + // segs.clear(); + + segs.emplace_back(Point_2(0, 0), Point_2(2, 0)); + segs.emplace_back(Point_2(2, 0), Point_2(1, 1)); + segs.emplace_back(Point_2(1, 1), Point_2(3, 1)); + test(segs); + segs.clear(); + + segs.emplace_back(Point_2(0, 1+e), Point_2(2, 1)); + segs.emplace_back(Point_2(1, 1), Point_2(1, 0)); + test(segs); + segs.clear(); + + // segs.emplace_back(Point_2(0, 0), Point_2(1-e, 0)); + // segs.emplace_back(Point_2(1, 1), Point_2(1, -1)); + // test(segs); + // segs.clear(); + segs.emplace_back(Point_2(1-e, 1), Point_2(-1-e, -1+2*e)); segs.emplace_back(Point_2(e/2, e/2), Point_2(1, -1)); segs.emplace_back(Point_2(0, 2-e/2), Point_2(2, 0)); @@ -199,6 +223,35 @@ void fix_test(){ test(segs); } +void big_test(){ + std::cout << "Fix tests" << std::endl; + std::vector< Segment_2 > segs; + segs.emplace_back(Point_2(0.99999999999992173,-0.9999999999991368),Point_2(0.99999999999987266,-6.016984285966173e-13)); + segs.emplace_back(Point_2(1.000000000000494,-1.0000000000002354),Point_2(0.99999999999950062,8.0391743541535603e-13)); + segs.emplace_back(Point_2(1.0000000000008489,-0.99999999999951794),Point_2(0.99999999999926925,-2.3642280455149055e-13)); + + test(segs); +} + +void big_test_2(){ + std::cout << "Fix tests" << std::endl; + std::vector< Segment_2 > segs; + segs.emplace_back(Point_2(1.0000000000008309,9.0061990813884068e-13), Point_2(0.99999999999981259,0.99999999999938394)); + segs.emplace_back(Point_2(1.0000000000005094,-3.1951552372595695e-13), Point_2(1.0000000000004887,1.0000000000006302)); + segs.emplace_back(Point_2(1.0000000000006171,-7.3034227828590228e-13), Point_2(1.0000000000004181,0.99999999999913269)); + segs.emplace_back(Point_2(1.0000000000008491,6.6531536071947998e-14), Point_2(0.9999999999993564,1.0000000000006497)); + segs.emplace_back(Point_2(1.0000000000006859,3.4155608842536406e-13), Point_2(0.99999999999990252,0.99999999999978451)); + segs.emplace_back(Point_2(1.0000000000007121,-1.6363168568197086e-13), Point_2(1.0000000000004778,0.99999999999927969)); + segs.emplace_back(Point_2(1.000000000000753,8.3769047200168284e-13), Point_2(1.0000000000003773,1.0000000000001725)); + segs.emplace_back(Point_2(1.000000000000574,6.6944274898354285e-13), Point_2(1.000000000000427,0.99999999999998845)); + segs.emplace_back(Point_2(1.0000000000006153,-2.798953111315494e-13), Point_2(1.0000000000003044,1.0000000000008289)); + segs.emplace_back(Point_2(1.0000000000008387,5.8214218612051864e-13), Point_2(1.0000000000002918,0.99999999999943212)); + segs.emplace_back(Point_2(1.0000000000007745,-3.9231414307222458e-13), Point_2(1.000000000000294,1.0000000000003408)); + + test(segs); + +} + void test_float_snap_rounding(){ std::cout << "Fix tests" << std::endl; @@ -214,7 +267,7 @@ void test_float_snap_rounding(){ segs.emplace_back(Point_2(5, 7-e), Point_2(9, 7-e)); double_snap_rounding_2(segs.begin(), segs.end(), out); - assert(out.size()==segs.size()); + // assert(out.size()==segs.size()); } int main(int argc,char *argv[]) @@ -223,10 +276,11 @@ int main(int argc,char *argv[]) CGAL::Random r(argc==1?rp.get_seed():std::stoi(argv[1])); std::cout << "random seed = " << r.get_seed() << std::endl; std::cout << std::setprecision(17); + big_test_2(); fix_test(); test_float_snap_rounding(); - // test_fully_random(r,1000); - // test_multi_almost_indentical_segments(r,100); - test_iterative_square_intersection(r,2000); + test_fully_random(r,1000); + test_multi_almost_indentical_segments(r,100); + // test_iterative_square_intersection(r,2000); return(0); } From 134330ebe85bda429ce5aada2de1adce4c02579f Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?L=C3=A9o=20Valque?= Date: Thu, 23 Oct 2025 18:18:50 +0200 Subject: [PATCH 23/41] bounding box filter and combinatorial to test to improve running time --- .../include/CGAL/Float_snap_rounding_2.h | 167 +++++++----------- .../test_float_snap_rounding_2.cpp | 11 +- 2 files changed, 71 insertions(+), 107 deletions(-) diff --git a/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h b/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h index c6b5a9414ff..3e278e8b0ac 100644 --- a/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h +++ b/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h @@ -120,6 +120,7 @@ public: template void scan(PointsRange &pts, PolylineRange &polylines, const Traits &traits){ using Point_2 = typename Traits::Point_2; + using Vector_2 = typename Traits::Vector_2; using Segment_2 = typename Traits::Segment_2; using Line_2 = typename Traits::Line_2; @@ -167,7 +168,7 @@ void scan(PointsRange &pts, PolylineRange &polylines, const Traits &traits){ auto event_comparator=[&](Event a, Event b){ const Point_2 &pa=pts[a.pi]; const Point_2 &pb=pts[b.pi]; - if(pa.x()!=pb.x()) + if(a.pi!=b.pi && pa.x()!=pb.x()) return Less_x_2()(pa,pb); // We want remove previous line before to insert new ones if(a.type!=b.type) @@ -176,8 +177,8 @@ void scan(PointsRange &pts, PolylineRange &polylines, const Traits &traits){ return a.type < b.type; else return a.type > b.type; - // We compare with y only to have a complete order - return Less_y_2()(pa,pb); + // We compare pi only to have a complete order + return a.pi < b.pi; }; auto pi_below_li=[&](size_t li, std::pair pair_pi_li){ @@ -186,7 +187,10 @@ void scan(PointsRange &pts, PolylineRange &polylines, const Traits &traits){ if(ind_pi_support==li) return false; const Polyline &pl=polylines[li]; - Orientation ori=orientation(pts[pl.front()], pts[pl.back()], pts[pi]); + + Orientation ori=COLLINEAR; + if(pl.front()!=pi && pl.back()!=pi) + ori=orientation(pts[pl.front()], pts[pl.back()], pts[pi]); // When collinear, pick an other vertex of the line of pi to decide if(ori==COLLINEAR){ const Polyline &pi_support = polylines[ind_pi_support]; @@ -260,19 +264,26 @@ void scan(PointsRange &pts, PolylineRange &polylines, const Traits &traits){ Polyline &pl=polylines[li]; Segment_2 seg(pts[pl.front()], pts[pl.back()]); - // (A+B)^2 <= 4*max(A^2,B^2), we take some margin + if(certainly(csq_dist_2(p, Line_2(pts[pl[0]], pts[pl[1]]), 1.)==LARGER)) + return false; + + // (A+B)^2 <= 4*max(A^2,B^2) double bound=round_bound_pts[pi]; for(std::size_t i=0; i= bb_seg.ymin() && bb_seg.ymax() >= bb_p.ymin() && possibly(csq_dist_2(p, seg, bound)!=CGAL::LARGER)) + // if(bb_p.ymax()+bound >= bb_seg.ymin() && bb_seg.ymax()+bound >= bb_p.ymin() && possibly(csq_dist_2(p, seg, bound)!=CGAL::LARGER)) { // TODO check duplicates for(std::size_t i: pl) if(pts[i].x()==pts[pi].x()) return false; - pts.emplace_back(p.x(), seg.supporting_line().y_at_x(p.x())); + pts.emplace_back(pts[event.pi].x(), seg.supporting_line().y_at_x(pts[event.pi].x())); + // pts.emplace_back(p.x(), seg.supporting_line().y_at_x(p.x())); std::size_t new_pi=pts.size()-1; round_bound_pts.emplace_back(round_bound(pts[new_pi])); // We insert it on pl before the last vertex @@ -355,7 +366,6 @@ void round_and_post_process(PointsRange &pts, PolylineRange &polylines, const Tr { using Point_2 = typename Traits::Point_2; using Segment_2 = typename Traits::Segment_2; - using Vector_2 = typename Traits::Vector_2; using Line_2 = typename Traits::Line_2; using Polyline = std::remove_cv_t::value_type>; @@ -385,12 +395,6 @@ void round_and_post_process(PointsRange &pts, PolylineRange &polylines, const Tr Equal_2 equal = traits.equal_2_object(); Round_2 round = traits.round_2_object(); - auto Less_indexes_x_2=[&](std::size_t i, std::size_t j){ - return Less_x_2()(pts[i], pts[j]); - }; - auto Less_indexes_y_2=[&](std::size_t i, std::size_t j){ - return Less_y_2()(pts[i], pts[j]); - }; auto Less_indexes_xy_2=[&](std::size_t i, std::size_t j){ return Less_xy_2()(pts[i], pts[j]); }; @@ -476,14 +480,54 @@ void round_and_post_process(PointsRange &pts, PolylineRange &polylines, const Tr std::swap(poly, updated_polyline); ++i; } + + //y + std::set p_sort_by_y(Less_indexes_yx_2); + for(std::size_t i=0; i!=pts.size(); ++i) + p_sort_by_y.insert(i); + + using Iterator_set_y = typename std::set::iterator; + + i=0; + for(Polyline &poly: polylines){ + std::vector updated_polyline; + updated_polyline.push_back(poly.front()); + for(std::size_t i=1; i!=poly.size(); ++i){ + if(pts[poly[i-1]].y()==pts[poly[i]].y()){ + // updated_polyline.clear(); + // std::swap(poly, updated_polyline); + // break; + Iterator_set_x start, end; + // Get all vertices between the two endpoints along x order + if(Less_indexes_yx_2(poly[i-1],poly[i])){ + start=p_sort_by_y.upper_bound(poly[i-1]); + end=p_sort_by_y.lower_bound(poly[i]); + for(auto it=start; it!=end; ++it){ + updated_polyline.push_back(*it); + } + } else { + start=p_sort_by_y.upper_bound(poly[i]); + end=p_sort_by_y.lower_bound(poly[i-1]); + for(auto it=end; it!=start;){ + --it; + updated_polyline.push_back(*it); + } + } + } + updated_polyline.push_back(poly[i]); + } + std::swap(poly, updated_polyline); + ++i; + } } template void double_snap_rounding_2_disjoint(PointsRange &pts, PolylineRange &polylines, const Traits &traits) { - scan(pts, polylines, traits); - round_and_post_process(pts, polylines, traits); - return; + // scan(pts, polylines, traits); + // round_and_post_process(pts, polylines, traits); + // return; + std::cout << "_________OLD___________" << std::endl; using Point_2 = typename Traits::Point_2; using Segment_2 = typename Traits::Segment_2; @@ -652,92 +696,7 @@ void double_snap_rounding_2_disjoint(PointsRange &pts, PolylineRange &polylines, CGAL_assertion((pl[0]==ps) && (pl[pl.size()-1]==pt)); } -#ifdef DOUBLE_2D_SNAP_VERBOSE - std::cout << "Round" << std::endl; -#endif - for(auto &p: pts) - p=round(p); - - //The order of the points on an axis must be preserved for the correctness of the algorithm - CGAL_assertion(std::is_sorted(p_sort_by_x.begin(),p_sort_by_x.end(),Less_indexes_x_2)); - CGAL_assertion(std::is_sorted(p_sort_by_y.begin(),p_sort_by_y.end(),Less_indexes_y_2)); - -#ifdef DOUBLE_2D_SNAP_VERBOSE - std::cout << "Remove duplicate points" << std::endl; -#endif - std::vector< std::size_t > unique_points(p_sort_by_x.begin(),p_sort_by_x.end()); - std::sort(unique_points.begin(),unique_points.end(),Less_indexes_xy_2); - std::vector new_pts; - std::vector old_to_new_index(pts.size()); - for(std::size_t i=0; i!=pts.size(); ++i){ - if(i==0 || !equal(pts[unique_points[i]],pts[unique_points[i-1]])) - new_pts.push_back(pts[unique_points[i]]); - old_to_new_index[unique_points[i]]=new_pts.size()-1; - } - - std::swap(pts, new_pts); - for (auto& polyline : polylines) { - std::vector updated_polyline; - for (std::size_t i=0; i updated_polyline; - updated_polyline.push_back(poly[0]); - for(std::size_t i=1; i!=poly.size(); ++i){ - if(pts[poly[i-1]].x()==pts[poly[i]].x()){ - Iterator_set_x start, end; - // Get all vertices between the two endpoints along x order - if(Less_indexes_xy_2(poly[i-1],poly[i])){ - start=p_sort_by_x.upper_bound(poly[i-1]); - end=p_sort_by_x.lower_bound(poly[i]); - } else { - start=p_sort_by_x.upper_bound(poly[i]); - end=p_sort_by_x.lower_bound(poly[i-1]); - } - // Add all endpoints between them to the polyline - for(auto it=start; it!=end; ++it){ - updated_polyline.push_back(*it); - } - } - if(pts[poly[i-1]].y()==pts[poly[i]].y()){ - Iterator_set_y start, end; - // Get all vertices between the two endpoints along y order - if(Less_indexes_yx_2(poly[i-1],poly[i])){ - start=p_sort_by_y.upper_bound(poly[i-1]); - end=p_sort_by_y.lower_bound(poly[i]); - } else { - start=p_sort_by_y.upper_bound(poly[i]); - end=p_sort_by_y.lower_bound(poly[i-1]); - } - // Add all endpoints between them to the polyline - for(auto it=start; it!=end; ++it){ - updated_polyline.push_back(*it); - } - } - updated_polyline.push_back(poly[i]); - } - std::swap(poly, updated_polyline); - } + round_and_post_process(pts, polylines, traits); } diff --git a/Snap_rounding_2/test/Snap_rounding_2/test_float_snap_rounding_2.cpp b/Snap_rounding_2/test/Snap_rounding_2/test_float_snap_rounding_2.cpp index d90e15f932e..3e838d43e7e 100644 --- a/Snap_rounding_2/test/Snap_rounding_2/test_float_snap_rounding_2.cpp +++ b/Snap_rounding_2/test/Snap_rounding_2/test_float_snap_rounding_2.cpp @@ -99,6 +99,10 @@ void test(const std::vector &segs){ std::cout << "Formal snap size: " << out.size() << " ,running time: " << t.time() << std::endl; #endif assert(!CGAL::do_curves_intersect(out.begin(), out.end())); + if(CGAL::do_curves_intersect(out.begin(), out.end())){ + std::cout << "INTERSECTING OUTPUT!!" << std::endl; + exit(1); + } #ifdef COMPARE_WITH_INTEGER_SNAP_ROUNDING_2 Polyline_range_2 output_list; t.reset(); @@ -276,10 +280,11 @@ int main(int argc,char *argv[]) CGAL::Random r(argc==1?rp.get_seed():std::stoi(argv[1])); std::cout << "random seed = " << r.get_seed() << std::endl; std::cout << std::setprecision(17); - big_test_2(); - fix_test(); + // big_test_2(); + test_almost_indentical_segments(r, 50, Vector_2(1,1), Vector_2(-1,-1)); + // fix_test(); test_float_snap_rounding(); - test_fully_random(r,1000); + // test_fully_random(r,1000); test_multi_almost_indentical_segments(r,100); // test_iterative_square_intersection(r,2000); return(0); From f82659bba797be6d8fce4d44acdf0c61ca761baf Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?L=C3=A9o=20Valque?= Date: Mon, 27 Oct 2025 18:34:42 +0100 Subject: [PATCH 24/41] cleaning and optimize Float_2D_snap --- .../include/CGAL/Float_snap_rounding_2.h | 328 ++++-------------- .../test_float_snap_rounding_2.cpp | 22 +- 2 files changed, 84 insertions(+), 266 deletions(-) diff --git a/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h b/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h index 3e278e8b0ac..83211669730 100644 --- a/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h +++ b/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h @@ -94,33 +94,10 @@ struct Float_snap_rounding_traits_2: Arr_segment_traits_2{ } }; -namespace Box_intersection_d { - -// Since we used only box_self_intersection_d, we may not have intersection of two distinct boxes with same index -template -class Box_with_index_d: public Box_d< NT_, N, ID_EXPLICIT> { -protected: - std::size_t m_index; -public: - typedef Box_d< NT_, N, ID_EXPLICIT> Base; - typedef NT_ NT; - typedef std::size_t ID; - - Box_with_index_d() {} - Box_with_index_d( ID i) : m_index(i) {} - Box_with_index_d( bool complete, ID i): Base(complete), m_index(i) {} - Box_with_index_d(NT l[N], NT h[N], ID i) : Base( l, h), m_index(i) {} - Box_with_index_d( const Bbox_2& b, ID i) : Base( b), m_index(i) {} - Box_with_index_d( const Bbox_3& b, ID i) : Base( b), m_index(i) {} - ID index() const { return m_index; } -}; - -} - template -void scan(PointsRange &pts, PolylineRange &polylines, const Traits &traits){ +void one_way_scan(PointsRange &pts, PolylineRange &polylines, const Traits &traits){ using Point_2 = typename Traits::Point_2; - using Vector_2 = typename Traits::Vector_2; + using FT = typename Traits::FT; using Segment_2 = typename Traits::Segment_2; using Line_2 = typename Traits::Line_2; @@ -168,17 +145,15 @@ void scan(PointsRange &pts, PolylineRange &polylines, const Traits &traits){ auto event_comparator=[&](Event a, Event b){ const Point_2 &pa=pts[a.pi]; const Point_2 &pb=pts[b.pi]; - if(a.pi!=b.pi && pa.x()!=pb.x()) - return Less_x_2()(pa,pb); - // We want remove previous line before to insert new ones - if(a.type!=b.type) - // Except of course if it is the same line (implying that line is vertical), we have to insert it before removing it - if(a.li==b.li) - return a.type < b.type; + if(a.pi!=b.pi) + if(pa!=pb) + return Less_xy_2()(pa,pb); else - return a.type > b.type; - // We compare pi only to have a complete order - return a.pi < b.pi; + return a.pi < b.pi; + + if(a.li!=b.li) + return a.type > b.type; + return a.type < b.type; }; auto pi_below_li=[&](size_t li, std::pair pair_pi_li){ @@ -241,14 +216,17 @@ void scan(PointsRange &pts, PolylineRange &polylines, const Traits &traits){ // TODO Vector is suboptimal arbitrary insertion, looking or redblack tree of skip list // Track order of the lines along y along the events std::vector< size_t > y_order; - for(Event event: event_queue){ + y_order.push_back(event_queue[0].li); + for(size_t i=1; i::iterator pos_it = y_order.begin(); if(!y_order.empty()) pos_it = std::lower_bound(y_order.begin(), y_order.end(), std::pair(event.pi, event.li), pi_below_li); -#ifdef DOUBLE_2D_SNAP_FULL_VERBOSE - std::cout << ((event.type==EVENT_TYPE::INSERT)?"Insert": "Remove") << " Event at point " << event.pi << " on line "<< event.li << std::endl; -#endif if(event.type==EVENT_TYPE::REMOVE){ assert(*pos_it==event.li); y_order.erase(pos_it); @@ -256,34 +234,44 @@ void scan(PointsRange &pts, PolylineRange &polylines, const Traits &traits){ if(!y_order.empty()){ auto close_event=[&](std::size_t pi, std::size_t li){ - // std::cout << pi << " " << li << std::endl; if((pi==polylines[li].front()) || (pi==polylines[li].back())) - return false; + return true; const Point_2 &p = pts[pi]; Polyline &pl=polylines[li]; Segment_2 seg(pts[pl.front()], pts[pl.back()]); - if(certainly(csq_dist_2(p, Line_2(pts[pl[0]], pts[pl[1]]), 1.)==LARGER)) - return false; - // (A+B)^2 <= 4*max(A^2,B^2) double bound=round_bound_pts[pi]; - for(std::size_t i=0; i= bb_seg.ymin() && bb_seg.ymax() >= bb_p.ymin() && possibly(csq_dist_2(p, seg, bound)!=CGAL::LARGER)) - // if(bb_p.ymax()+bound >= bb_seg.ymin() && bb_seg.ymax()+bound >= bb_p.ymin() && possibly(csq_dist_2(p, seg, bound)!=CGAL::LARGER)) + if(possibly(csq_dist_2(p, seg, bound)!=CGAL::LARGER)) { - // TODO check duplicates + pts[pi].exact(); + pts[pl[pl.size()-2]].exact(); + pts[pl.back()].exact(); + pts[pl.front()].exact(); + round_bound_pts[pi]=round_bound(pts[pi]); + round_bound_pts[pl[pl.size()-2]]=round_bound(pts[pl[pl.size()-2]]); + round_bound_pts[pl.back()]=round_bound(pts[pl.back()]); + round_bound_pts[pl.front()]=round_bound(pts[pl.front()]); + bound=round_bound_pts[pi]; + bound = (std::max)(bound, round_bound_pts[pl[pl.size()-2]]); + bound = (std::max)(bound, round_bound_pts[pl.back()]); + bound*=4; + if(csq_dist_2(p, seg, bound)==CGAL::LARGER){ + return false; + } + + // Check duplicates for(std::size_t i: pl) if(pts[i].x()==pts[pi].x()) return false; - pts.emplace_back(pts[event.pi].x(), seg.supporting_line().y_at_x(pts[event.pi].x())); - // pts.emplace_back(p.x(), seg.supporting_line().y_at_x(p.x())); + + FT y= (seg.supporting_line().y_at_x(pts[event.pi].x())); + pts.emplace_back(pts[event.pi].x(), y); std::size_t new_pi=pts.size()-1; round_bound_pts.emplace_back(round_bound(pts[new_pi])); // We insert it on pl before the last vertex @@ -301,32 +289,30 @@ void scan(PointsRange &pts, PolylineRange &polylines, const Traits &traits){ return false; }; - // Look above segments and creates a point if there too close to safe snap - auto above = pos_it; - size_t pi=event.pi; - while(above!=y_order.end() && close_event(pi,*above)){ - ++above; - pi=pts.size()-1; - } - - // same with below segments - auto below = pos_it; - pi=event.pi; - if(below!=y_order.begin()){ - --below; - while(close_event(pi,*(below))){ - if(below==y_order.begin()) - break; - pi=pts.size()-1; - --below; + if(event.pi != event_queue[i-1].pi) + { + // Look above segments and creates a point if there too close to safe snap + auto above = pos_it; + size_t pi=event.pi; + while(above!=y_order.end() && close_event(pi,*above)){ + if((pi!=polylines[*above].front()) && (pi!=polylines[*above].back())) + pi=pts.size()-1; + ++above; } - } - while(close_event(pi,*(below))){ - if(below==y_order.begin()) - break; - pi=pts.size()-1; - --below; + // same with below segments + auto below = pos_it; + pi=event.pi; + if(below!=y_order.begin()){ + --below; + while(close_event(pi,*(below))){ + if(below==y_order.begin()) + break; + if((pi!=polylines[*below].front()) && (pi!=polylines[*below].back())) + pi=pts.size()-1; + --below; + } + } } } @@ -455,9 +441,6 @@ void round_and_post_process(PointsRange &pts, PolylineRange &polylines, const Tr updated_polyline.push_back(poly.front()); for(std::size_t i=1; i!=poly.size(); ++i){ if(pts[poly[i-1]].x()==pts[poly[i]].x()){ - // updated_polyline.clear(); - // std::swap(poly, updated_polyline); - // break; Iterator_set_x start, end; // Get all vertices between the two endpoints along x order if(Less_indexes_xy_2(poly[i-1],poly[i])){ @@ -494,9 +477,6 @@ void round_and_post_process(PointsRange &pts, PolylineRange &polylines, const Tr updated_polyline.push_back(poly.front()); for(std::size_t i=1; i!=poly.size(); ++i){ if(pts[poly[i-1]].y()==pts[poly[i]].y()){ - // updated_polyline.clear(); - // std::swap(poly, updated_polyline); - // break; Iterator_set_x start, end; // Get all vertices between the two endpoints along x order if(Less_indexes_yx_2(poly[i-1],poly[i])){ @@ -524,178 +504,7 @@ void round_and_post_process(PointsRange &pts, PolylineRange &polylines, const Tr template void double_snap_rounding_2_disjoint(PointsRange &pts, PolylineRange &polylines, const Traits &traits) { - // scan(pts, polylines, traits); - // round_and_post_process(pts, polylines, traits); - // return; - std::cout << "_________OLD___________" << std::endl; - - using Point_2 = typename Traits::Point_2; - using Segment_2 = typename Traits::Segment_2; - using Vector_2 = typename Traits::Vector_2; - using Line_2 = typename Traits::Line_2; - - using Polyline = std::remove_cv_t::value_type>; - - using PointsRangeIterator = typename PointsRange::iterator; - using PolylineRangeIterator = typename PolylineRange::iterator; - - using PBox = CGAL::Box_intersection_d::Box_with_info_d; - using SBox = CGAL::Box_intersection_d::Box_with_handle_d; - - using Less_x_2 = typename Traits::Less_x_2; - using Less_y_2 = typename Traits::Less_y_2; - using Less_xy_2 = typename Traits::Less_xy_2; - using Less_yx_2 = typename Traits::Less_yx_2; - using Equal_2 = typename Traits::Equal_2; - using Round_2 = typename Traits::Round_2; - - using Construct_source_2 = typename Traits::Construct_source_2; - using Construct_target_2 = typename Traits::Construct_target_2; - using Compare_squared_distance_2 = typename Traits::Compare_squared_distance_2; - using Squared_round_bound_2 = typename Traits::Squared_round_bound_2; - - Compare_squared_distance_2 csq_dist_2 = traits.compare_squared_distance_2_object(); - Squared_round_bound_2 round_bound = traits.squared_round_bound_2_object(); - Construct_source_2 source = traits.construct_source_2_object(); - Construct_target_2 target = traits.construct_target_2_object(); - Equal_2 equal = traits.equal_2_object(); - Round_2 round = traits.round_2_object(); - - auto Less_indexes_x_2=[&](std::size_t i, std::size_t j){ - return Less_x_2()(pts[i], pts[j]); - }; - auto Less_indexes_y_2=[&](std::size_t i, std::size_t j){ - return Less_y_2()(pts[i], pts[j]); - }; - auto Less_indexes_xy_2=[&](std::size_t i, std::size_t j){ - return Less_xy_2()(pts[i], pts[j]); - }; - auto Less_indexes_yx_2=[&](std::size_t i, std::size_t j){ - return Less_yx_2()(pts[i], pts[j]); - }; - -#ifdef DOUBLE_2D_SNAP_VERBOSE - std::cout << "Nb of points: " << pts.size() << " , nb of polylines: " << polylines.size() << std::endl; - std::cout << "Sort the input points" << std::endl; -#endif - // Compute the order of the points along the 2 axis - // Sorted the points may perform exact computations and thus refine the intervals of the coordinates values - // This refine ensures that the order of the points will be preserved when rounded - // However, except for this reason these sets are unused - using Iterator_set_x = typename std::set::iterator; - using Iterator_set_y = typename std::set::iterator; - std::set p_sort_by_x(Less_indexes_xy_2); - std::set p_sort_by_y(Less_indexes_yx_2); - for(std::size_t i=0; i!=pts.size(); ++i) - { - p_sort_by_x.insert(i); - p_sort_by_y.insert(i); - } - - const double max_coordinate=std::max(std::max(to_double(pts[*p_sort_by_x.begin()].x()), to_double(pts[*(--p_sort_by_x.end())].x())), - std::max(to_double(pts[*p_sort_by_y.begin()].y()), to_double(pts[*(--p_sort_by_y.end())].y()))); - const double global_bound=max_coordinate*std::pow(2, -20); - - //Prepare boxes for box_intersection_d - std::vector new_points; - std::vector points_boxes; - std::vector segs_boxes; - - //We create that vector to avoid multiple computations - std::vector round_bound_pts; - - for(std::size_t i=0; isize()-1]].bbox(),it); - - CGAL_MUTEX mutex_callback; - // Callback used for box_intersection_d - auto callback=[&](PBox &bp, SBox &bseg){ - std::size_t pi=bp.info(); - - Polyline& pl=*(bseg.handle()); - std::size_t si1=pl[0]; - std::size_t si2=pl[pl.size()-1]; - - //Check if p is one endpoint of the segment - if((pi==si1) || (pi==si2)) - return; - - Point_2& p= pts[pi]; - - // Early exit for better running time - if(certainly(csq_dist_2(p, Line_2(pts[pl[0]], pts[pl[1]]), global_bound)==LARGER)) - return; - - Segment_2 seg(pts[pl[0]], pts[pl[1]]); - - // (A+B)^2 <= 4*max(A^2,B^2), we take some margin - double bound=round_bound_pts[pi]; - for(std::size_t i=0; i(points_boxes.begin(), points_boxes.end(), segs_boxes.begin(), segs_boxes.end(), callback); - points_boxes.clear(); - points_boxes.reserve(pts.size()-size_before); -#ifdef DOUBLE_2D_SNAP_VERBOSE - ++turn_nb; - std::cout << "Turn " << turn_nb << ": " << pts.size()-size_before << " subdivisions performed" << std::endl; -#endif - // The new vertices may intersect another segment when rounded, we repeat until they are not new vertices - for(std::size_t i=size_before; i > set_out_segs; output.clear(); - // for(auto &poly: polylines){ - // for(std::size_t i=1; i &segs){ assert(!CGAL::do_curves_intersect(out.begin(), out.end())); if(CGAL::do_curves_intersect(out.begin(), out.end())){ std::cout << "INTERSECTING OUTPUT!!" << std::endl; + std::vector bad_pts; + std::vector out2; + for(auto seg: out) + out2.emplace_back(seg.source(), seg.target()); + CGAL::compute_intersection_points(out2.begin(), out2.end(), std::back_inserter(bad_pts)); + std::cout << bad_pts[0] << std::endl; + for(size_t i=0; i event_queue; event_queue.reserve(2*polylines.size()); for(size_t li=0; li y_order; y_order.push_back(event_queue[0].li); for(size_t i=1; i::iterator pos_it = y_order.begin(); if(!y_order.empty()) @@ -249,10 +243,13 @@ void one_way_scan(PointsRange &pts, PolylineRange &polylines, const Traits &trai if(possibly(csq_dist_2(p, seg, bound)!=CGAL::LARGER)) { + // We refine the pts to reduce the rounding shift and check again pts[pi].exact(); pts[pl[pl.size()-2]].exact(); + // The two following call of exact act act on seg variable since they appear in its DAG pts[pl.back()].exact(); pts[pl.front()].exact(); + // Update the bounds round_bound_pts[pi]=round_bound(pts[pi]); round_bound_pts[pl[pl.size()-2]]=round_bound(pts[pl[pl.size()-2]]); round_bound_pts[pl.back()]=round_bound(pts[pl.back()]); @@ -261,37 +258,38 @@ void one_way_scan(PointsRange &pts, PolylineRange &polylines, const Traits &trai bound = (std::max)(bound, round_bound_pts[pl[pl.size()-2]]); bound = (std::max)(bound, round_bound_pts[pl.back()]); bound*=4; - if(csq_dist_2(p, seg, bound)==CGAL::LARGER){ - return false; - } - // Check duplicates - for(std::size_t i: pl) - if(pts[i].x()==pts[pi].x()) - return false; + // Check if the point and the segment are still too closed for a safe rounding + if(csq_dist_2(p, seg, bound)!=CGAL::LARGER){ + // Check if segment was not already subdivided by another point + for(std::size_t i: pl) + if(pts[i].x()==pts[pi].x()) + return false; - FT y= (seg.supporting_line().y_at_x(pts[event.pi].x())); - pts.emplace_back(pts[event.pi].x(), y); - std::size_t new_pi=pts.size()-1; - round_bound_pts.emplace_back(round_bound(pts[new_pi])); - // We insert it on pl before the last vertex - pl.insert(pl.end()-1, new_pi); + // Create a point on seg at the same x coordinate than p + FT y= (seg.supporting_line().y_at_x(pts[event.pi].x())); + pts.emplace_back(pts[event.pi].x(), y); + std::size_t new_pi=pts.size()-1; + round_bound_pts.emplace_back(round_bound(pts[new_pi])); + // We insert it on pl before the last vertex + pl.insert(pl.end()-1, new_pi); #ifdef DOUBLE_2D_SNAP_FULL_VERBOSE - std::cout << "Create point " << new_pi << " on " << li << " due to proximity with " << pi << "_____________________________" << std::endl; - std::cout << new_pi <<": " << pts[new_pi] << std::endl; - std::cout << li << ":"; - for(std::size_t i: pl) - std::cout << " " << i; - std::cout << std::endl; + std::cout << "Create point " << new_pi << " on " << li << " due to proximity with " << pi << "_____________________________" << std::endl; + std::cout << new_pi <<": " << pts[new_pi] << std::endl; + std::cout << li << ":"; + for(std::size_t i: pl) + std::cout << " " << i; + std::cout << std::endl; #endif - return true; + return true; + } } return false; }; if(event.pi != event_queue[i-1].pi) { - // Look above segments and creates a point if there too close to safe snap + // Look above segments and creates a point if there too close for a safe rounding auto above = pos_it; size_t pi=event.pi; while(above!=y_order.end() && close_event(pi,*above)){ @@ -322,6 +320,7 @@ void one_way_scan(PointsRange &pts, PolylineRange &polylines, const Traits &trai } // We sort the point by y to ensure the order of the point is preserved + // Round value of y coordinates are indirectly modified in the case of a filter failure auto sort_pi=[&](std::size_t i, std::size_t j){ return pts[i].y() < pts[j].y(); }; @@ -348,55 +347,21 @@ void one_way_scan(PointsRange &pts, PolylineRange &polylines, const Traits &trai } template -void round_and_post_process(PointsRange &pts, PolylineRange &polylines, const Traits &traits) +std::size_t merge_duplicate_points_in_polylines(PointsRange &pts, PolylineRange &polylines, const Traits &traits) { - using Point_2 = typename Traits::Point_2; - using Segment_2 = typename Traits::Segment_2; - using Line_2 = typename Traits::Line_2; - + using P_ID = std::size_t; + using Point_2 = typename Traits::Point_2; using Polyline = std::remove_cv_t::value_type>; - using PointsRangeIterator = typename PointsRange::iterator; - using PolylineRangeIterator = typename PolylineRange::iterator; - - using PBox = CGAL::Box_intersection_d::Box_with_info_d; - using SBox = CGAL::Box_intersection_d::Box_with_handle_d; - - using Less_x_2 = typename Traits::Less_x_2; - using Less_y_2 = typename Traits::Less_y_2; using Less_xy_2 = typename Traits::Less_xy_2; - using Less_yx_2 = typename Traits::Less_yx_2; - using Equal_2 = typename Traits::Equal_2; - using Round_2 = typename Traits::Round_2; + using Unique_point_container = std::map; - using Construct_source_2 = typename Traits::Construct_source_2; - using Construct_target_2 = typename Traits::Construct_target_2; - using Compare_squared_distance_2 = typename Traits::Compare_squared_distance_2; - using Squared_round_bound_2 = typename Traits::Squared_round_bound_2; - - Compare_squared_distance_2 csq_dist_2 = traits.compare_squared_distance_2_object(); - Squared_round_bound_2 round_bound = traits.squared_round_bound_2_object(); - Construct_source_2 source = traits.construct_source_2_object(); - Construct_target_2 target = traits.construct_target_2_object(); + using Equal_2 = typename Traits::Equal_2; Equal_2 equal = traits.equal_2_object(); - Round_2 round = traits.round_2_object(); - auto Less_indexes_xy_2=[&](std::size_t i, std::size_t j){ - return Less_xy_2()(pts[i], pts[j]); - }; - auto Less_indexes_yx_2=[&](std::size_t i, std::size_t j){ - return Less_yx_2()(pts[i], pts[j]); - }; - -#ifdef DOUBLE_2D_SNAP_VERBOSE - std::cout << "Round" << std::endl; -#endif - for(auto &p: pts) - p=round(p); - -#ifdef DOUBLE_2D_SNAP_VERBOSE - std::cout << "Remove duplicate points from collapsing" << std::endl; -#endif + auto Less_indexes_xy_2=[&](std::size_t i, std::size_t j){ + return Less_xy_2()(pts[i], pts[j]); + }; std::vector< std::size_t > unique_points(pts.size()); std::iota(unique_points.begin(), unique_points.end(), 0); @@ -423,11 +388,72 @@ void round_and_post_process(PointsRange &pts, PolylineRange &polylines, const Tr ++i; std::swap(polyline, updated_polyline); } + return 0; -#ifdef DOUBLE_2D_SNAP_VERBOSE - // The algorithm prevents the a vertex that goes through a segment but a vertex may lie on an horizontal/vertical segments after rounding - std::cout << "Subdivide vertical segments with vertices on them" << std::endl; -#endif +// const std::size_t ini_points_n = pts.size(); +// std::vector point_index(ini_points_n, 0); + +// Unique_point_container point_to_id(traits.less_xy_2_object()); + +// std::vector unique_points; +// unique_points.reserve(ini_points_n); + +// for(std::size_t i=0; i is_insert_successful = +// point_to_id.insert(std::make_pair(pts[i], unique_points.size())); + +// // #ifdef DOUBLE_2D_SNAP_FULL_VERBOSE +// if(!is_insert_successful.second) +// std::cout << "points[" << i << "] = " << pts[i] << " was already encountered" << std::endl; +// // #endif +// std::size_t id = is_insert_successful.first->second; + +// if(id == unique_points.size()) +// unique_points.push_back(pts[i]); +// point_index[i] = id; +// } + +// if(unique_points.size() != ini_points_n) +// { +// for(P_ID pl_index=0, end=polylines.size(); pl_index!=end; ++pl_index) +// { +// Polyline& pl = polylines[pl_index]; +// for(std::size_t i=0, pl_size = pl.size(); i 0) +// std::cout << "Removed (merged) " << removed_points_n << " duplicate points" << std::endl; +// // #endif +// return removed_points_n; +} + +template +void snap_post_process(PointsRange &pts, PolylineRange &polylines, const Traits &traits) +{ + using Polyline = std::remove_cv_t::value_type>; + using Less_xy_2 = typename Traits::Less_xy_2; + + Less_xy_2 less_xy_2 = traits.less_xy_2_object(); + + auto Less_indexes_xy_2=[&](std::size_t i, std::size_t j){ + return less_xy_2(pts[i], pts[j]); + }; std::set p_sort_by_x(Less_indexes_xy_2); for(std::size_t i=0; i!=pts.size(); ++i) @@ -435,7 +461,7 @@ void round_and_post_process(PointsRange &pts, PolylineRange &polylines, const Tr using Iterator_set_x = typename std::set::iterator; - i=0; + std::size_t i=0; for(Polyline &poly: polylines){ std::vector updated_polyline; updated_polyline.push_back(poly.front()); @@ -463,51 +489,36 @@ void round_and_post_process(PointsRange &pts, PolylineRange &polylines, const Tr std::swap(poly, updated_polyline); ++i; } - - //y - std::set p_sort_by_y(Less_indexes_yx_2); - for(std::size_t i=0; i!=pts.size(); ++i) - p_sort_by_y.insert(i); - - using Iterator_set_y = typename std::set::iterator; - - i=0; - for(Polyline &poly: polylines){ - std::vector updated_polyline; - updated_polyline.push_back(poly.front()); - for(std::size_t i=1; i!=poly.size(); ++i){ - if(pts[poly[i-1]].y()==pts[poly[i]].y()){ - Iterator_set_x start, end; - // Get all vertices between the two endpoints along x order - if(Less_indexes_yx_2(poly[i-1],poly[i])){ - start=p_sort_by_y.upper_bound(poly[i-1]); - end=p_sort_by_y.lower_bound(poly[i]); - for(auto it=start; it!=end; ++it){ - updated_polyline.push_back(*it); - } - } else { - start=p_sort_by_y.upper_bound(poly[i]); - end=p_sort_by_y.lower_bound(poly[i-1]); - for(auto it=end; it!=start;){ - --it; - updated_polyline.push_back(*it); - } - } - } - updated_polyline.push_back(poly[i]); - } - std::swap(poly, updated_polyline); - ++i; - } } template void double_snap_rounding_2_disjoint(PointsRange &pts, PolylineRange &polylines, const Traits &traits) { - one_way_scan(pts, polylines, traits); - round_and_post_process(pts, polylines, traits); + using Point_2 = typename Traits::Point_2; + using Construct_round_point_2 = typename Traits::Construct_round_point_2; + + Construct_round_point_2 round = traits.construct_round_point_2_object(); + + snap_rounding_scan(pts, polylines, traits); + +#ifdef DOUBLE_2D_SNAP_VERBOSE + std::cout << "Round" << std::endl; +#endif + for(Point_2 &p: pts) + p=round(p); +#ifdef DOUBLE_2D_SNAP_VERBOSE + std::cout << "Merging points that are collapsed together" << std::endl; +#endif + merge_duplicate_points_in_polylines(pts, polylines, traits); +#ifdef DOUBLE_2D_SNAP_VERBOSE + // The algorithm prevents the a vertex that goes through a segment but a vertex may lie on an horizontal/vertical segments after rounding + std::cout << "Subdivide vertical segments with vertices on them" << std::endl; +#endif + + snap_post_process(pts, polylines, traits); } +} // end of namespace internal /** * ingroup @@ -602,7 +613,7 @@ typename OutputContainer::iterator double_snap_rounding_2(InputIterator input_ } // Main algorithm - double_snap_rounding_2_disjoint(pts, polylines, traits); + internal::double_snap_rounding_2_disjoint(pts, polylines, traits); #ifdef DOUBLE_2D_SNAP_VERBOSE std::cout << "Build output" << std::endl; @@ -651,7 +662,9 @@ typename OutputContainer::iterator compute_snapped_subcurves_2(InputIterator using Segment_2 = typename Traits::Segment_2; using I2E = typename Traits::Converter_in; using E2O = typename Traits::Converter_out; - using VectorIterator = typename std::vector::iterator; + + using SegmentRange = std::vector; + using SegmentRangeIterator = typename SegmentRange::iterator; using Less_xy_2 = typename Traits::Less_xy_2; @@ -682,7 +695,7 @@ typename OutputContainer::iterator compute_snapped_subcurves_2(InputIterator std::vector< std::vector< std::size_t> > polylines; // Transform range of the segments in the range of points and polyline of indexes - for(VectorIterator it=segs.begin(); it!=segs.end(); ++it) + for(SegmentRangeIterator it=segs.begin(); it!=segs.end(); ++it) { const Point_2& p1 = it->source(); const Point_2& p2 = it->target(); @@ -699,7 +712,7 @@ typename OutputContainer::iterator compute_snapped_subcurves_2(InputIterator } } - for(VectorIterator it=segs.begin(); it!=segs.end(); ++it) + for(SegmentRangeIterator it=segs.begin(); it!=segs.end(); ++it) { std::size_t index1 = point_to_index[it->source()]; std::size_t index2 = point_to_index[it->target()]; @@ -711,7 +724,7 @@ typename OutputContainer::iterator compute_snapped_subcurves_2(InputIterator // Main algorithm - double_snap_rounding_2_disjoint(pts, polylines, traits); + internal::double_snap_rounding_2_disjoint(pts, polylines, traits); #ifdef DOUBLE_2D_SNAP_VERBOSE std::cout << "Build output" << std::endl; @@ -839,7 +852,7 @@ void snap_polygons_2(InputIterator begin, } // Main algorithm - double_snap_rounding_2_disjoint(pts, polylines, traits); + internal::double_snap_rounding_2_disjoint(pts, polylines, traits); #ifdef DOUBLE_2D_SNAP_VERBOSE std::cout << "Build output" << std::endl; diff --git a/Snap_rounding_2/test/Snap_rounding_2/test_float_snap_rounding_2.cpp b/Snap_rounding_2/test/Snap_rounding_2/test_float_snap_rounding_2.cpp index 361e3a39e2c..f221633add3 100644 --- a/Snap_rounding_2/test/Snap_rounding_2/test_float_snap_rounding_2.cpp +++ b/Snap_rounding_2/test/Snap_rounding_2/test_float_snap_rounding_2.cpp @@ -296,9 +296,10 @@ int main(int argc,char *argv[]) CGAL::Random r(argc==1?rp.get_seed():std::stoi(argv[1])); std::cout << "random seed = " << r.get_seed() << std::endl; std::cout << std::setprecision(17); - // big_test_2(); + big_test(); + big_test_2(); // test_almost_indentical_segments(r, 50, Vector_2(1,1), Vector_2(-1,-1)); - // fix_test(); + fix_test(); test_float_snap_rounding(); test_fully_random(r,1000); test_multi_almost_indentical_segments(r,200); From 58caea293adc75d453203d88e4f899a9c4c4d966 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?L=C3=A9o=20Valque?= Date: Tue, 28 Oct 2025 18:00:06 +0100 Subject: [PATCH 26/41] More cleaning of Float_snap_rounding_2 --- .../include/CGAL/Float_snap_rounding_2.h | 75 +++---------------- 1 file changed, 11 insertions(+), 64 deletions(-) diff --git a/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h b/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h index 934e7d09b44..4dd5b77c63c 100644 --- a/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h +++ b/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h @@ -347,21 +347,18 @@ void snap_rounding_scan(PointsRange &pts, PolylineRange &polylines, const Traits } template -std::size_t merge_duplicate_points_in_polylines(PointsRange &pts, PolylineRange &polylines, const Traits &traits) +void merge_duplicate_points_in_polylines(PointsRange &pts, PolylineRange &polylines, const Traits &traits) { - using P_ID = std::size_t; - using Point_2 = typename Traits::Point_2; + using Point_2 = typename Traits::Point_2; using Polyline = std::remove_cv_t::value_type>; using Less_xy_2 = typename Traits::Less_xy_2; - using Unique_point_container = std::map; - using Equal_2 = typename Traits::Equal_2; Equal_2 equal = traits.equal_2_object(); - auto Less_indexes_xy_2=[&](std::size_t i, std::size_t j){ - return Less_xy_2()(pts[i], pts[j]); - }; + auto Less_indexes_xy_2=[&](std::size_t i, std::size_t j){ + return Less_xy_2()(pts[i], pts[j]); + }; std::vector< std::size_t > unique_points(pts.size()); std::iota(unique_points.begin(), unique_points.end(), 0); @@ -388,61 +385,9 @@ std::size_t merge_duplicate_points_in_polylines(PointsRange &pts, PolylineRange ++i; std::swap(polyline, updated_polyline); } - return 0; - -// const std::size_t ini_points_n = pts.size(); -// std::vector point_index(ini_points_n, 0); - -// Unique_point_container point_to_id(traits.less_xy_2_object()); - -// std::vector unique_points; -// unique_points.reserve(ini_points_n); - -// for(std::size_t i=0; i is_insert_successful = -// point_to_id.insert(std::make_pair(pts[i], unique_points.size())); - -// // #ifdef DOUBLE_2D_SNAP_FULL_VERBOSE -// if(!is_insert_successful.second) -// std::cout << "points[" << i << "] = " << pts[i] << " was already encountered" << std::endl; -// // #endif -// std::size_t id = is_insert_successful.first->second; - -// if(id == unique_points.size()) -// unique_points.push_back(pts[i]); -// point_index[i] = id; -// } - -// if(unique_points.size() != ini_points_n) -// { -// for(P_ID pl_index=0, end=polylines.size(); pl_index!=end; ++pl_index) -// { -// Polyline& pl = polylines[pl_index]; -// for(std::size_t i=0, pl_size = pl.size(); i 0) -// std::cout << "Removed (merged) " << removed_points_n << " duplicate points" << std::endl; -// // #endif -// return removed_points_n; } +// Some points may have collapsed on a vertical segment, we subdivide these vertical segments accordingly template void snap_post_process(PointsRange &pts, PolylineRange &polylines, const Traits &traits) { @@ -676,6 +621,8 @@ typename OutputContainer::iterator compute_snapped_subcurves_2(InputIterator I2E to_exact=traits.converter_to_exact_object(); E2O from_exact=traits.converter_from_exact_object(); + Less_xy_2 less_xy_2 = traits.less_xy_2_object(); + std::vector convert_input; for(InputIterator it=input_begin; it!=input_end; ++it) convert_input.push_back(Segment_2(to_exact(*it))); @@ -689,8 +636,8 @@ typename OutputContainer::iterator compute_snapped_subcurves_2(InputIterator #ifdef DOUBLE_2D_SNAP_VERBOSE std::cout << "Change format to range of points and indexes" << std::endl; #endif - std::set unique_point_set; - std::map point_to_index; + std::set unique_point_set(less_xy_2); + std::map point_to_index; std::vector pts; std::vector< std::vector< std::size_t> > polylines; @@ -716,7 +663,7 @@ typename OutputContainer::iterator compute_snapped_subcurves_2(InputIterator { std::size_t index1 = point_to_index[it->source()]; std::size_t index2 = point_to_index[it->target()]; - if(Less_xy_2()(it->source(), it->target())) + if(less_xy_2(it->source(), it->target())) polylines.push_back({index1, index2}); else polylines.push_back({index2, index1}); From ec7b77b40adc6cc40ed68406decd337a626b28ea Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?L=C3=A9o=20Valque?= Date: Mon, 3 Nov 2025 17:52:27 +0100 Subject: [PATCH 27/41] Enhance documentation of FLoat snap rounding 2 --- .../internal/parameters_interface.h | 5 +- .../Concepts/FloatSnapRoundingTraits_2.h | 218 +++++++++++++++++ .../Snap_rounding_2/PackageDescription.txt | 2 + .../include/CGAL/Float_snap_rounding_2.h | 231 ++++++++---------- .../CGAL/Float_snap_rounding_traits_2.h | 114 +++++++++ .../test_float_snap_rounding_2.cpp | 24 +- 6 files changed, 449 insertions(+), 145 deletions(-) create mode 100644 Snap_rounding_2/doc/Snap_rounding_2/Concepts/FloatSnapRoundingTraits_2.h create mode 100644 Snap_rounding_2/include/CGAL/Float_snap_rounding_traits_2.h diff --git a/STL_Extension/include/CGAL/STL_Extension/internal/parameters_interface.h b/STL_Extension/include/CGAL/STL_Extension/internal/parameters_interface.h index f70fb981631..d99ac998b0e 100644 --- a/STL_Extension/include/CGAL/STL_Extension/internal/parameters_interface.h +++ b/STL_Extension/include/CGAL/STL_Extension/internal/parameters_interface.h @@ -406,8 +406,9 @@ CGAL_add_named_parameter(do_not_modify_geometry_t, do_not_modify_geometry, do_no CGAL_add_named_parameter_with_compatibility_ref_only(angles_param_t, angles_param, angles) CGAL_add_named_parameter(maximum_height_t, maximum_height, maximum_height) +//List of named parameters used in Snap_rounding_2 +CGAL_add_named_parameter(compute_intersection_t, compute_intersection, compute_intersection) + // List of named parameters used in the package 'Constrained_triangulation_3' CGAL_add_named_parameter(plc_face_id_t, plc_face_id, plc_face_id) CGAL_add_named_parameter(with_plc_face_id_t, with_plc_face_id, with_plc_face_id) -//List of named parameters used in Snap_rounding_2 -CGAL_add_named_parameter(do_intersection_computation_t, do_intersection_computation, do_intersection_computation) diff --git a/Snap_rounding_2/doc/Snap_rounding_2/Concepts/FloatSnapRoundingTraits_2.h b/Snap_rounding_2/doc/Snap_rounding_2/Concepts/FloatSnapRoundingTraits_2.h new file mode 100644 index 00000000000..ef4e9c54be6 --- /dev/null +++ b/Snap_rounding_2/doc/Snap_rounding_2/Concepts/FloatSnapRoundingTraits_2.h @@ -0,0 +1,218 @@ + +/*! +\ingroup PkgFloatSnapRounding2Concepts +\cgalConcept + +The concept `FloatSnapRoundingTraits_2` lists the set of requirements that must be fulfilled by +an instance of the `Traits` template-parameter of +the free functions \ref CGAL::snap_rounding_ 2() `CGAL::snap_rounding_2()`. +The list includes the nested types of the geometric primitives used in this class and +some function object types for the required predicates on those primitives. + +\cgalRefines{AosTraits_2} + +\cgalHasModelsBegin +\cgalHasModels{CGAL::Float_snap_rounding_traits_2} +\cgalHasModelsEnd +*/ + +class FloatSnapRoundingTraits_2 { +public: + +/// \name Types +/// @{ + +/*! +The number type. This type must fulfill the requirements on +`FieldNumberType` +*/ +typedef unspecified_type FT; + +/*! +Models the concept `ArrTraits::Point_2`. +*/ +typedef unspecified_type Point_2; + +/*! +Models the concept `ArrTraits::XMonotoneCurve_2`. +*/ +typedef unspecified_type Segment_2; + +/// @} + +/// \name Functor Types +/// @{ + +/*! +Models the concept `Kernel::ConstructSource_2`. +*/ +typedef unspecified_type Construct_source_2; + +/*! +Models the concept `Kernel::ConstructTarget_2`. +*/ +typedef unspecified_type Construct_target_2; + +/*! +Models the concept `Kernel::ConstructSegment_2`. +*/ +typedef unspecified_type Construct_segment_2; + +/*! +Models the concept `Kernel::LessXY_2`. +*/ +typedef unspecified_type Less_xy_2; + +/*! +Models the concept `Kernel::LessY_2`. +*/ +typedef unspecified_type Less_y_2; + +/*! +Models the concept `Kernel::Equal_2`. +*/ +typedef unspecified_type Equal_2; + +/*! +Models the concept `FSRTraits_2::ConstructRoundPoint_2`. +*/ +typedef unspecified_type Construct_round_point_2; + +/*! +Models the concept `FSRTraits_2::SquaredRoundBound_2`. +*/ +typedef unspecified_type Squared_round_bound_2; + +/*! +Models the concept `FSRTraits_2::ConverterToExact`. +*/ +typedef unspecified_type Converter_to_exact; + +/*! +Models the concept `FSRTraits_2::ConverterFromExact`. +*/ +typedef unspecified_type Converter_from_exact; + +/// @} + +/// \name Accessing Functor Objects +/// @{ + +/*! + +*/ +Construct_source_2 construct_source_2_object(); + +/*! + +*/ +Construct_target_2 construct_target_2_object(); + +/*! + +*/ +Construct_segment_2 construct_segment_2_object(); + +/*! + +*/ +Less_xy_2 less_xy_2_object(); + +/*! + +*/ +Less_y_2 less_y_2_object(); + +/*! + +*/ +Construct_round_point_2 construct_round_point_2_object(); + +/*! + +*/ +Squared_round_bound_2 squared_round_bound_2_object(); + + +/// @} + +}; /* end FloatSnapRoundingTraits_2 */ + + +namespace FSRTraits_2{ + +/*! + \ingroup PkgFloatSnapRounding2Concepts + \cgalConcept + \cgalHasModelsBegin + \cgalHasModelsBare{\link FloatSnapRoundingTraits_2::ConstructRoundPoint_2 `Float_snap_rounding_traits_2::Construct_round_point_2` \endlink} + \cgalHasModelsEnd +*/ +class ConstructRoundPoint_2 +{ + public: + + /*! + Given a point, construct its rounded version + */ + Point_2 operator()(Point_2 p); +}; + +/*! + \ingroup PkgFloatSnapRounding2Concepts + \cgalConcept + \cgalHasModelsBegin + \cgalHasModelsBare{\link FloatSnapRoundingTraits_2::SquaredRoundBound_2 `Float_snap_rounding_traits_2::Squared_round_bound_2` \endlink} + \cgalHasModelsEnd +*/ +class SquaredRoundBound_2 +{ + public: + + /*! + Given a point, compute an upper bound of the squared distance between its exact coordinates and its rounded coordinates + */ + double operator()(Point_2 p); +}; + +/*! + \ingroup PkgFloatSnapRounding2Concepts + \cgalConcept + \cgalHasModelsBegin + \cgalHasModelsBare{\link FloatSnapRoundingTraits_2::ConverterToExact `Float_snap_rounding_traits_2::Convertex_to_exact` \endlink} + \cgalHasModelsEnd +*/ +class ConverterToExact +{ + public: + typedef unspecified_type InputPoint_2; + typedef unspecified_type InputSegment_2; + /*! + Convert input points (i.e. segments) into points (i.e. segments) of the type of the traits + */ + Point_2 operator()(InputPoint_2 p); + Segment_2 operator()(InputSegment_2 p); +}; + +/*! + \ingroup PkgFloatSnapRounding2Concepts + \cgalConcept + \cgalHasModelsBegin + \cgalHasModelsBare{\link FloatSnapRoundingTraits_2::ConverterToExact `Float_snap_rounding_traits_2::Convertex_to_exact` \endlink} + \cgalHasModelsEnd +*/ +class ConverterFromExact +{ + public: + typedef unspecified_type OutputPoint_2; + typedef unspecified_type OutputSegment_2; + /*! + Convert points (i.e. segments) of the type of the traits into points (i.e. segments) of the output + */ + OutputPoint_2 operator()(Point_2 p); + OutputSegment_2 operator()(Segment_2 p); +}; + + + +} /* end of namespace SRTraits_2 */ diff --git a/Snap_rounding_2/doc/Snap_rounding_2/PackageDescription.txt b/Snap_rounding_2/doc/Snap_rounding_2/PackageDescription.txt index 9426f5b4808..b67d6079f80 100644 --- a/Snap_rounding_2/doc/Snap_rounding_2/PackageDescription.txt +++ b/Snap_rounding_2/doc/Snap_rounding_2/PackageDescription.txt @@ -1,5 +1,7 @@ /// \defgroup PkgSnapRounding2Ref Reference Manual /// \defgroup PkgSnapRounding2Concepts Concepts +/// \defgroup PkgFloatSnapRounding2Ref Reference Manual +/// \defgroup PkgFloatSnapRounding2Concepts Concepts /// \ingroup PkgSnapRounding2Ref /*! \addtogroup PkgSnapRounding2Ref diff --git a/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h b/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h index 4dd5b77c63c..6f58e0e4e16 100644 --- a/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h +++ b/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h @@ -17,13 +17,11 @@ #include #endif -#include #include -#include -#include -#include -#include -#include + +#include + +#include #include #include @@ -33,70 +31,15 @@ namespace CGAL { -template -struct Float_snap_rounding_traits_2: Arr_segment_traits_2{ - typedef Arr_segment_traits_2 Base; - - typedef typename Base::FT FT; - typedef typename Base::Point_2 Point_2; - typedef typename Base::Segment_2 Segment_2; - typedef typename Base::Vector_2 Vector_2; - typedef typename Base::Line_2 Line_2; - - typedef typename Base::Less_x_2 Less_x_2; - typedef typename Base::Less_y_2 Less_y_2; - typedef typename Base::Less_xy_2 Less_xy_2; - typedef typename Base::Less_yx_2 Less_yx_2; - typedef typename Base::Equal_2 Equal_2; - - typedef typename Base::Construct_source_2 Construct_source_2; - typedef typename Base::Construct_target_2 Construct_target_2; - - typedef Cartesian_converter Converter_in; - typedef Cartesian_converter Converter_out; - - // Return an upperbound of the squared distance between a point and its rounded value - struct Squared_round_bound_2{ - double operator()(const FT &x) const{ - double b=std::nextafter(to_interval(x).second - to_interval(x).first, std::numeric_limits::infinity()); - return b*b; - } - double operator()(const Point_2 &p) const{ - return (*this)(p.x())+(*this)(p.y()); - } - }; - - struct Construct_round_point_2{ - double operator()(const FT &x) const{ - return to_double(x); - } - Point_2 operator()(const Point_2 &p) const{ - return Point_2((*this)(p.x()),(*this)(p.y())); - } - }; - - Converter_in converter_to_exact_object() const{ - return Converter_in(); - } - - Converter_out converter_from_exact_object() const{ - return Converter_out(); - } - - Squared_round_bound_2 squared_round_bound_2_object() const{ - return Squared_round_bound_2(); - } - - Construct_round_point_2 construct_round_point_2_object() const{ - return Construct_round_point_2(); - } -}; - namespace internal{ + /* + Scan the vertices from left to right while maintening the y order of the segments. + Subdivide the segments if there are too close to a vertex + */ template void snap_rounding_scan(PointsRange &pts, PolylineRange &polylines, const Traits &traits){ - using Point_2 = typename Traits::Point_2; using FT = typename Traits::FT; + using Point_2 = typename Traits::Point_2; using Segment_2 = typename Traits::Segment_2; using Polyline = std::remove_cv_t::value_type>; @@ -105,14 +48,17 @@ void snap_rounding_scan(PointsRange &pts, PolylineRange &polylines, const Traits using PolylineRangeIterator = typename PolylineRange::iterator; using Less_xy_2 = typename Traits::Less_xy_2; - using Less_yx_2 = typename Traits::Less_yx_2; - using Equal_2 = typename Traits::Equal_2; + using Less_y_2 = typename Traits::Less_y_2; + using Construct_segment_2 = typename Traits::Construct_segment_2; + using Construct_point_at_x_on_segment_2 = typename Traits::Construct_point_at_x_on_segment_2; using Compare_squared_distance_2 = typename Traits::Compare_squared_distance_2; using Squared_round_bound_2 = typename Traits::Squared_round_bound_2; Less_xy_2 less_xy_2 = traits.less_xy_2_object(); - Equal_2 equal = traits.equal_2_object(); + Less_y_2 less_y_2 = traits.less_y_2_object(); + Construct_segment_2 segment_2 = traits.construct_segment_2_object(); + Construct_point_at_x_on_segment_2 point_at_x = traits.construct_point_at_x_on_segment_2_object(); Compare_squared_distance_2 csq_dist_2 = traits.compare_squared_distance_2_object(); Squared_round_bound_2 round_bound = traits.squared_round_bound_2_object(); @@ -187,7 +133,7 @@ void snap_rounding_scan(PointsRange &pts, PolylineRange &polylines, const Traits std::cout << "Create the Event Queue" << std::endl; #endif - // Inout polylines are supposed going from left to right + // Input polylines are supposed going from left to right std::vector event_queue; event_queue.reserve(2*polylines.size()); for(size_t li=0; li y_order; y_order.push_back(event_queue[0].li); @@ -233,7 +179,7 @@ void snap_rounding_scan(PointsRange &pts, PolylineRange &polylines, const Traits const Point_2 &p = pts[pi]; Polyline &pl=polylines[li]; - Segment_2 seg(pts[pl.front()], pts[pl.back()]); + Segment_2 seg=segment_2(pts[pl.front()], pts[pl.back()]); // (A+B)^2 <= 4*max(A^2,B^2) double bound=round_bound_pts[pi]; @@ -246,7 +192,7 @@ void snap_rounding_scan(PointsRange &pts, PolylineRange &polylines, const Traits // We refine the pts to reduce the rounding shift and check again pts[pi].exact(); pts[pl[pl.size()-2]].exact(); - // The two following call of exact act act on seg variable since they appear in its DAG + // The two following call of exact act on seg variables since they appear in its DAG pts[pl.back()].exact(); pts[pl.front()].exact(); // Update the bounds @@ -267,8 +213,7 @@ void snap_rounding_scan(PointsRange &pts, PolylineRange &polylines, const Traits return false; // Create a point on seg at the same x coordinate than p - FT y= (seg.supporting_line().y_at_x(pts[event.pi].x())); - pts.emplace_back(pts[event.pi].x(), y); + pts.push_back(point_at_x(seg, pts[event.pi].x())); std::size_t new_pi=pts.size()-1; round_bound_pts.emplace_back(round_bound(pts[new_pi])); // We insert it on pl before the last vertex @@ -322,7 +267,7 @@ void snap_rounding_scan(PointsRange &pts, PolylineRange &polylines, const Traits // We sort the point by y to ensure the order of the point is preserved // Round value of y coordinates are indirectly modified in the case of a filter failure auto sort_pi=[&](std::size_t i, std::size_t j){ - return pts[i].y() < pts[j].y(); + return less_y_2(pts[i],pts[j]); }; std::vector< std::size_t > indices(pts.size(),0); std::iota(indices.begin(),indices.end(),0); @@ -445,7 +390,6 @@ void double_snap_rounding_2_disjoint(PointsRange &pts, PolylineRange &polylines, Construct_round_point_2 round = traits.construct_round_point_2_object(); snap_rounding_scan(pts, polylines, traits); - #ifdef DOUBLE_2D_SNAP_VERBOSE std::cout << "Round" << std::endl; #endif @@ -459,7 +403,6 @@ void double_snap_rounding_2_disjoint(PointsRange &pts, PolylineRange &polylines, // The algorithm prevents the a vertex that goes through a segment but a vertex may lie on an horizontal/vertical segments after rounding std::cout << "Subdivide vertical segments with vertices on them" << std::endl; #endif - snap_post_process(pts, polylines, traits); } @@ -471,16 +414,36 @@ void double_snap_rounding_2_disjoint(PointsRange &pts, PolylineRange &polylines, * Given a range of segments, compute rounded subsegments that are pairwise disjoint in their interiors, based on the input curves. * The output is a range of polyline with each polyline corresponding to an input segment. * +* TODO Currently compute_subcurves have no visitor to obtain a polyline per input segment, thus +* currently the output contain one polyline per ubsegment of the arrangement +* * @tparam Concurrency_tag allows choosing whether the algorithm runs in parallel or sequentially. * If `CGAL::Parallel_tag` is specified and CGAL is linked with the Intel TBB library, the algorithm will run in parallel. * Otherwise, if `CGAL::Sequential_tag` is specified (the default), the algorithm will run sequentially. * @tparam InputIterator iterator over a range of `Segment_2` * @tparam OutputContainer inserter over a range of `Polyline`. `Polyline` must be a type that provides a `push_back(Point_2)` function. +* @tparam NamedParameters a sequence of \ref bgl_namedparameters "Named Parameters" +* \param begin,end the input segment range +* \param out the output container +* \param np an optional sequence of \ref bgl_namedparameters "Named Parameters" among the ones listed below +* \cgalNamedParamsBegin +* \cgalParamNBegin{concurrency_tag} +* \cgalParamDescription{That template parameter enables to choose whether the algorithm is to be run in parallel, if CGAL::Parallel_tag +* is specified and CGAL has been linked with the Intel TBB library, or sequentially, otherwise.} +* \cgalParamType{CGAL::Concurrency_tag} +* \cgalParamDefault{CGAL::Sequential_tag} +* \cgalParamNEnd +* +* \cgalParamNBegin{geom_traits} +* \cgalParamDescription{an instance of a geometric traits class} +* \cgalParamType{The traits class must respect the concept of `FloatSnapRoundingTraits_2`} +* \cgalParamDefault{an instance of `Float_snap_rounding_traits_2`} +* \cgalParamNEnd */ template -typename OutputContainer::iterator double_snap_rounding_2(InputIterator input_begin, - InputIterator input_end, - OutputContainer& output, +typename OutputContainer::iterator double_snap_rounding_2(InputIterator begin, + InputIterator end, + OutputContainer &out, const NamedParameters &np = parameters::default_values()) { using Concurrency_tag = typename internal_np::Lookup_named_param_def::iterator; using Polyline = std::remove_cv_t::value_type>; @@ -512,8 +475,8 @@ typename OutputContainer::iterator double_snap_rounding_2(InputIterator input_ E2O from_exact=traits.converter_from_exact_object(); std::vector convert_input; - for(InputIterator it=input_begin; it!=input_end; ++it) - convert_input.push_back(Segment_2(to_exact(*it))); + for(InputIterator it=begin; it!=end; ++it) + convert_input.push_back(to_exact(*it)); std::vector segs; #ifdef DOUBLE_2D_SNAP_VERBOSE std::cout << "Solved intersections" << std::endl; @@ -565,15 +528,14 @@ typename OutputContainer::iterator double_snap_rounding_2(InputIterator input_ #endif // Output polylines - output.clear(); for(auto &poly: polylines){ Polyline new_line; for(std::size_t pi: poly) new_line.push_back(pts[pi]); - output.push_back(new_line); + out.push_back(new_line); } - return output.begin(); + return out.end(); } /** @@ -585,13 +547,29 @@ typename OutputContainer::iterator double_snap_rounding_2(InputIterator input_ * parallel, if CGAL::Parallel_tag is specified and CGAL has been linked with the Intel TBB library, or sequentially, if CGAL::Sequential_tag - the default value - is specified. * @tparam InputIterator iterator of a segment range * @tparam OutputContainer inserter of a segment range -* @tparam The exact kernel needed for computation (Epeck by default) +* @tparam NamedParameters a sequence of \ref bgl_namedparameters "Named Parameters" +* \param begin,end the input segment range +* \param out the output inserter +* \param np an optional sequence of \ref bgl_namedparameters "Named Parameters" among the ones listed below +* \cgalNamedParamsBegin +* \cgalParamNBegin{concurrency_tag} +* \cgalParamDescription{That template parameter enables to choose whether the algorithm is to be run in parallel, if CGAL::Parallel_tag +* is specified and CGAL has been linked with the Intel TBB library, or sequentially, otherwise.} +* \cgalParamType{CGAL::Concurrency_tag} +* \cgalParamDefault{CGAL::Sequential_tag} +* \cgalParamNEnd +* +* \cgalParamNBegin{geom_traits} +* \cgalParamDescription{an instance of a geometric traits class} +* \cgalParamType{The traits class must respect the concept of `FloatSnapRoundingTraits_2`} +* \cgalParamDefault{an instance of `Float_snap_rounding_traits_2`} +* \cgalParamNEnd */ -template -typename OutputContainer::iterator compute_snapped_subcurves_2(InputIterator input_begin, - InputIterator input_end, - OutputContainer& output, - const NamedParameters &np = parameters::default_values()) +template +OutputIterator compute_snapped_subcurves_2(InputIterator begin, + InputIterator end, + OutputIterator out, + const NamedParameters &np = parameters::default_values()) { using Concurrency_tag = typename internal_np::Lookup_named_param_def; using SegmentRangeIterator = typename SegmentRange::iterator; using Less_xy_2 = typename Traits::Less_xy_2; + using Construct_segment_2 = typename Traits::Construct_segment_2; using parameters::choose_parameter; using parameters::get_parameter; @@ -622,10 +601,11 @@ typename OutputContainer::iterator compute_snapped_subcurves_2(InputIterator E2O from_exact=traits.converter_from_exact_object(); Less_xy_2 less_xy_2 = traits.less_xy_2_object(); + Construct_segment_2 segment_2 = traits.construct_segment_2_object(); std::vector convert_input; - for(InputIterator it=input_begin; it!=input_end; ++it) - convert_input.push_back(Segment_2(to_exact(*it))); + for(InputIterator it=begin; it!=end; ++it) + convert_input.push_back(to_exact(*it)); std::vector segs; #ifdef DOUBLE_2D_SNAP_VERBOSE std::cout << "Solved intersections" << std::endl; @@ -679,23 +659,23 @@ typename OutputContainer::iterator compute_snapped_subcurves_2(InputIterator // Output a range of segments while removing duplicate ones std::set< std::pair > set_out_segs; - output.clear(); for(auto &poly: polylines) for(std::size_t i=1; i +template void snap_polygons_2(InputIterator begin, InputIterator end, - OutputContainer out, + OutputIterator out, const NamedParameters &np = parameters::default_values()) { using Concurrency_tag = typename internal_np::Lookup_named_param_def::iterator; using parameters::choose_parameter; @@ -748,7 +734,7 @@ void snap_polygons_2(InputIterator begin, Polygon_2 P=*begin; const Traits &traits = choose_parameter(get_parameter(np, internal_np::geom_traits), Traits()); - const bool compute_intersections = choose_parameter(get_parameter(np, internal_np::do_intersection_computation), false); + const bool compute_intersections = choose_parameter(get_parameter(np, internal_np::compute_intersection), false); I2E to_exact=traits.converter_to_exact_object(); E2O from_exact=traits.converter_from_exact_object(); @@ -761,27 +747,8 @@ void snap_polygons_2(InputIterator begin, std::vector polygon_index; if(compute_intersections){ - std::set unique_point_set; - std::map point_to_index; - - // Transform the polygon in a range of points and polylines of indexes - for(const typename Polygon_2::Point_2 &p_: P.vertices()) - { - Point_2 p=to_exact(p_); - if (unique_point_set.find(p) == unique_point_set.end()) { - unique_point_set.insert(p); - pts.push_back(p); - point_to_index[p] = pts.size() - 1; - } - } - - for(const typename Polygon_2::Segment_2 &s: P.edges()){ - Point_2 p1=to_exact(s.source()); - Point_2 p2=to_exact(s.target()); - std::size_t index1 = point_to_index[p1]; - std::size_t index2 = point_to_index[p2]; - polylines.push_back({index1, index2}); - } + // TODO Need to compute_subcurves that output polylines to track the polygons + assert(0); } else { // Index of the segments that introduced a new polygon polygon_index.reserve(std::distance(begin, end)); diff --git a/Snap_rounding_2/include/CGAL/Float_snap_rounding_traits_2.h b/Snap_rounding_2/include/CGAL/Float_snap_rounding_traits_2.h new file mode 100644 index 00000000000..2a466b1f3e8 --- /dev/null +++ b/Snap_rounding_2/include/CGAL/Float_snap_rounding_traits_2.h @@ -0,0 +1,114 @@ +// Copyright (c) 2025 Geometry Factory. +// All rights reserved. +// +// This file is part of CGAL (www.cgal.org). +// +// $URL$ +// $Id$ +// SPDX-License-Identifier: GPL-3.0-or-later OR LicenseRef-Commercial +// +// +// author(s) : Léo Valque + +#ifndef CGAL_FLOAT_SNAP_ROUNDING_TRAITS_2_H +#define CGAL_FLOAT_SNAP_ROUNDING_TRAITS_2_H + +#include + +#include +#include + +#include + +#include + +namespace CGAL { + +/*! +\ingroup PkgFloatSnapRounding2Ref + +The class `Float_snap_rounding_traits_2` is a model of the +`FloatSnapRoundingTraits_2` concept, and is the only traits class supplied +with the package. +This class should be instantiated with a kernel with points and segments that are convertible to +`Exact_predicates_exact_constructions_kernel`, the user must specified to used another exact geometric kernel that conforms +to the \cgal kernel-concept, such as `Cartesian`. + +\cgalModels{SnapRoundingTraits_2} + +*/ +template +struct Float_snap_rounding_traits_2: Arr_segment_traits_2{ + typedef Arr_segment_traits_2 Base; + + typedef typename Base::FT FT; + typedef typename Base::Point_2 Point_2; + typedef typename Base::Segment_2 Segment_2; + + typedef typename Base::Less_xy_2 Less_xy_2; + typedef typename Base::Less_y_2 Less_y_2; + typedef typename Base::Equal_2 Equal_2; + + typedef typename Base::Construct_point_2 Construct_point_2; + typedef typename Base::Construct_source_2 Construct_source_2; + typedef typename Base::Construct_target_2 Construct_target_2; + typedef typename Base::Construct_segment_2 Construct_segment_2; + + typedef Cartesian_converter Converter_to_exact; + typedef Cartesian_converter Converter_from_exact; + + // Return an upperbound of the squared distance between a point and its rounded value + struct Squared_round_bound_2{ + double operator()(const FT &x) const{ + double b=std::nextafter(to_interval(x).second - to_interval(x).first, std::numeric_limits::infinity()); + return b*b; + } + double operator()(const Point_2 &p) const{ + return (*this)(p.x())+(*this)(p.y()); + } + double operator()(const Segment_2 &seg) const{ + return (std::max)( (*this)(Base().construct_source_2_object()(seg)), + (*this)(Base().construct_target_2_object()(seg))); + } + }; + + struct Construct_round_point_2{ + double operator()(const FT &x) const{ + return to_double(x); + } + Point_2 operator()(const Point_2 &p) const{ + return Point_2((*this)(p.x()),(*this)(p.y())); + } + }; + + struct Construct_point_at_x_on_segment_2{ + Point_2 operator()(const Segment_2 &seg, const FT &x) const{ + FT y= (seg.supporting_line().y_at_x(x)); + return Base().construct_point_2_object()(x, y); + } + }; + + Converter_to_exact converter_to_exact_object() const{ + return Converter_to_exact(); + } + + Converter_from_exact converter_from_exact_object() const{ + return Converter_from_exact(); + } + + Squared_round_bound_2 squared_round_bound_2_object() const{ + return Squared_round_bound_2(); + } + + Construct_round_point_2 construct_round_point_2_object() const{ + return Construct_round_point_2(); + } + + Construct_point_at_x_on_segment_2 construct_point_at_x_on_segment_2_object() const{ + return Construct_point_at_x_on_segment_2(); + } +}; + +} //namespace CGAL + +#endif diff --git a/Snap_rounding_2/test/Snap_rounding_2/test_float_snap_rounding_2.cpp b/Snap_rounding_2/test/Snap_rounding_2/test_float_snap_rounding_2.cpp index f221633add3..249fbf2de46 100644 --- a/Snap_rounding_2/test/Snap_rounding_2/test_float_snap_rounding_2.cpp +++ b/Snap_rounding_2/test/Snap_rounding_2/test_float_snap_rounding_2.cpp @@ -27,15 +27,19 @@ typedef Kernel::FT FT; typedef CGAL::Arr_segment_traits_2 Arr_segment_traits_2; typedef Arr_segment_traits_2::Curve_2 Curve_2; -typedef CGAL::Snap_rounding_traits_2 SnapTraits; +typedef CGAL::Polygon_2 Polygon_2; +typedef CGAL::Polygon_with_holes_2 Polygon_with_holes_2; +typedef std::vector Pwh_vec_2; -typedef CGAL::Polygon_2 Polygon_2; -typedef CGAL::Polygon_with_holes_2 Polygon_with_holes_2; -typedef std::vector Pwh_vec_2; +typedef CGAL::Cartesian Naive; +typedef CGAL::Cartesian_converter EK_to_IK; +typedef CGAL::Cartesian_converter IK_to_EK; -typedef CGAL::Cartesian Naive; -typedef CGAL::Cartesian_converter EK_to_IK; -typedef CGAL::Cartesian_converter IK_to_EK; +#ifdef COMPARE_WITH_INTEGER_SNAP_ROUNDING_2 +#include +#include +typedef CGAL::Snap_rounding_traits_2 SnapTraits; +#endif //Biggest double with ulp smaller than an integer constexpr double maxFloat = std::pow(2,23); @@ -90,10 +94,8 @@ void test(const std::vector &segs){ t.reset(); t.start(); #endif - CGAL::compute_snapped_subcurves_2(segs.begin(), segs.end(), out); + CGAL::compute_snapped_subcurves_2(segs.begin(), segs.end(), std::back_inserter(out)); std::vector out2(out.begin(), out.end()); - // out.clear(); - // CGAL::compute_snapped_subcurves_2(out2.begin(), out2.end(), out); #ifdef BENCH_AND_VERBOSE_FLOAT_SNAP_ROUNDING_2 t.stop(); std::cout << "Formal snap size: " << out.size() << " ,running time: " << t.time() << std::endl; @@ -123,7 +125,7 @@ void test(const std::vector &segs){ Polyline_range_2 output_list; t.reset(); t.start(); - CGAL::snap_rounding_2(segs.begin(), segs.end(), output_list, 1./maxDouble); + CGAL::snap_rounding_2(segs.begin(), segs.end(), std::back_inserter(output_list), 1./maxDouble); t.stop(); std::cout << "Running time with integer 2D snap (scaled at 10^15): " << t.time() << std::endl; #endif From fc83e5477b9b64bdc1c586cfc7f68ed2bd9eb7f1 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?L=C3=A9o=20Valque?= Date: Mon, 3 Nov 2025 18:21:23 +0100 Subject: [PATCH 28/41] change reference --- Snap_rounding_2/doc/Snap_rounding_2/PackageDescription.txt | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/Snap_rounding_2/doc/Snap_rounding_2/PackageDescription.txt b/Snap_rounding_2/doc/Snap_rounding_2/PackageDescription.txt index b67d6079f80..94641a384c2 100644 --- a/Snap_rounding_2/doc/Snap_rounding_2/PackageDescription.txt +++ b/Snap_rounding_2/doc/Snap_rounding_2/PackageDescription.txt @@ -1,7 +1,7 @@ /// \defgroup PkgSnapRounding2Ref Reference Manual /// \defgroup PkgSnapRounding2Concepts Concepts -/// \defgroup PkgFloatSnapRounding2Ref Reference Manual -/// \defgroup PkgFloatSnapRounding2Concepts Concepts +/// \defgroup PkgFloatSnapRounding2Ref Reference Manual* +/// \defgroup PkgFloatSnapRounding2Concepts Concepts* /// \ingroup PkgSnapRounding2Ref /*! \addtogroup PkgSnapRounding2Ref From 6918848e70eab8e4d552eea1a9b3cd5e0bc9d4e7 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?L=C3=A9o=20Valque?= Date: Thu, 6 Nov 2025 14:03:07 +0100 Subject: [PATCH 29/41] Fix error of the testsuite --- .../Concepts/FloatSnapRoundingTraits_2.h | 2 +- .../include/CGAL/Float_snap_rounding_2.h | 99 +++++++++---------- .../CGAL/Float_snap_rounding_traits_2.h | 2 + .../test/Snap_rounding_2/CMakeLists.txt | 1 + .../test_float_snap_rounding_2.cpp | 8 +- 5 files changed, 55 insertions(+), 57 deletions(-) diff --git a/Snap_rounding_2/doc/Snap_rounding_2/Concepts/FloatSnapRoundingTraits_2.h b/Snap_rounding_2/doc/Snap_rounding_2/Concepts/FloatSnapRoundingTraits_2.h index ef4e9c54be6..2f971628cc6 100644 --- a/Snap_rounding_2/doc/Snap_rounding_2/Concepts/FloatSnapRoundingTraits_2.h +++ b/Snap_rounding_2/doc/Snap_rounding_2/Concepts/FloatSnapRoundingTraits_2.h @@ -5,7 +5,7 @@ The concept `FloatSnapRoundingTraits_2` lists the set of requirements that must be fulfilled by an instance of the `Traits` template-parameter of -the free functions \ref CGAL::snap_rounding_ 2() `CGAL::snap_rounding_2()`. +the free functions \ref `CGAL::snap_polygons_2()`, `CGAL::compute_snapped_subcurves_2()` and `CGAL::double_snap_rounding_2()`. The list includes the nested types of the geometric primitives used in this class and some function object types for the required predicates on those primitives. diff --git a/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h b/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h index 6f58e0e4e16..e0468e25858 100644 --- a/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h +++ b/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h @@ -38,15 +38,11 @@ namespace internal{ */ template void snap_rounding_scan(PointsRange &pts, PolylineRange &polylines, const Traits &traits){ - using FT = typename Traits::FT; using Point_2 = typename Traits::Point_2; using Segment_2 = typename Traits::Segment_2; using Polyline = std::remove_cv_t::value_type>; - using PointsRangeIterator = typename PointsRange::iterator; - using PolylineRangeIterator = typename PolylineRange::iterator; - using Less_xy_2 = typename Traits::Less_xy_2; using Less_y_2 = typename Traits::Less_y_2; using Construct_segment_2 = typename Traits::Construct_segment_2; @@ -141,10 +137,10 @@ void snap_rounding_scan(PointsRange &pts, PolylineRange &polylines, const Traits event_queue.emplace_back(pl.front(), li, EVENT_TYPE::INSERT); for(size_t pi=1; pi){ + // We refine the pts to reduce the rounding shift and check again + pts[pi].exact(); + pts[pl[pl.size()-2]].exact(); + // The two following call of exact act on seg variables since they appear in its DAG + pts[pl.back()].exact(); + pts[pl.front()].exact(); + // Update the bounds + round_bound_pts[pi]=round_bound(pts[pi]); + round_bound_pts[pl[pl.size()-2]]=round_bound(pts[pl[pl.size()-2]]); + round_bound_pts[pl.back()]=round_bound(pts[pl.back()]); + round_bound_pts[pl.front()]=round_bound(pts[pl.front()]); + bound=round_bound_pts[pi]; + bound = (std::max)(bound, round_bound_pts[pl[pl.size()-2]]); + bound = (std::max)(bound, round_bound_pts[pl.back()]); + bound*=4; - // Check if the point and the segment are still too closed for a safe rounding - if(csq_dist_2(p, seg, bound)!=CGAL::LARGER){ - // Check if segment was not already subdivided by another point - for(std::size_t i: pl) - if(pts[i].x()==pts[pi].x()) - return false; - - // Create a point on seg at the same x coordinate than p - pts.push_back(point_at_x(seg, pts[event.pi].x())); - std::size_t new_pi=pts.size()-1; - round_bound_pts.emplace_back(round_bound(pts[new_pi])); - // We insert it on pl before the last vertex - pl.insert(pl.end()-1, new_pi); -#ifdef DOUBLE_2D_SNAP_FULL_VERBOSE - std::cout << "Create point " << new_pi << " on " << li << " due to proximity with " << pi << "_____________________________" << std::endl; - std::cout << new_pi <<": " << pts[new_pi] << std::endl; - std::cout << li << ":"; - for(std::size_t i: pl) - std::cout << " " << i; - std::cout << std::endl; -#endif - return true; + // Check if the point and the segment are still too closed for a safe rounding + if(csq_dist_2(p, seg, bound)==CGAL::LARGER) + return false; } + + // Check if segment was not already subdivided by another point + for(std::size_t i: pl) + if(pts[i].x()==pts[pi].x()) + return false; + + // Create a point on seg at the same x coordinate than p + pts.push_back(point_at_x(seg, pts[event.pi].x())); + std::size_t new_pi=pts.size()-1; + round_bound_pts.emplace_back(round_bound(pts[new_pi])); + // We insert it on pl before the last vertex + pl.insert(pl.end()-1, new_pi); +#ifdef DOUBLE_2D_SNAP_FULL_VERBOSE + std::cout << "Create point " << new_pi << " on " << li << " due to proximity with " << pi << "_____________________________" << std::endl; + std::cout << new_pi <<": " << pts[new_pi] << std::endl; + std::cout << li << ":"; + for(std::size_t i: pl) + std::cout << " " << i; + std::cout << std::endl; +#endif + return true; } return false; }; @@ -317,7 +316,6 @@ void merge_duplicate_points_in_polylines(PointsRange &pts, PolylineRange &polyli } std::swap(pts, new_pts); - size_t i=0; for (Polyline& polyline : polylines) { std::vector updated_polyline; for (std::size_t i=0; i::iterator; - std::size_t i=0; for(Polyline &poly: polylines){ std::vector updated_polyline; updated_polyline.push_back(poly.front()); @@ -377,7 +373,6 @@ void snap_post_process(PointsRange &pts, PolylineRange &polylines, const Traits updated_polyline.push_back(poly[i]); } std::swap(poly, updated_polyline); - ++i; } } @@ -531,7 +526,7 @@ typename OutputContainer::iterator double_snap_rounding_2(InputIterator begin, for(auto &poly: polylines){ Polyline new_line; for(std::size_t pi: poly) - new_line.push_back(pts[pi]); + new_line.push_back(from_exact(pts[pi])); out.push_back(new_line); } @@ -723,10 +718,8 @@ void snap_polygons_2(InputIterator begin, DefaultTraits>::type; using Point_2 = typename Traits::Point_2; - using Segment_2 = typename Traits::Segment_2; using I2E = typename Traits::Converter_to_exact; using E2O = typename Traits::Converter_from_exact; - using VectorIterator = typename std::vector::iterator; using parameters::choose_parameter; using parameters::get_parameter; @@ -773,9 +766,9 @@ void snap_polygons_2(InputIterator begin, #endif // Output a range of segments while removing duplicate ones - for(size_t input_ind=0; input_ind &poly = polylines[pl_ind]; + for(std::size_t input_ind=0; input_ind &poly = polylines[pl_ind]; Polygon_2 P; for(std::size_t i=1; i{ typedef Arr_segment_traits_2 Base; + typedef Exact_Kernel Exact_type; + typedef typename Base::FT FT; typedef typename Base::Point_2 Point_2; typedef typename Base::Segment_2 Segment_2; diff --git a/Snap_rounding_2/test/Snap_rounding_2/CMakeLists.txt b/Snap_rounding_2/test/Snap_rounding_2/CMakeLists.txt index 7645a98dbd8..94bc16ced5c 100644 --- a/Snap_rounding_2/test/Snap_rounding_2/CMakeLists.txt +++ b/Snap_rounding_2/test/Snap_rounding_2/CMakeLists.txt @@ -8,6 +8,7 @@ find_package(CGAL REQUIRED) create_single_source_cgal_program(test_snap_rounding_2.cpp NO_TESTING) create_single_source_cgal_program(test_float_snap_rounding_2.cpp) +create_single_source_cgal_program(test_float_snap_rounding_2_concept.cpp) function(add_Snap_rounding_tests name) set(data_dir "data") diff --git a/Snap_rounding_2/test/Snap_rounding_2/test_float_snap_rounding_2.cpp b/Snap_rounding_2/test/Snap_rounding_2/test_float_snap_rounding_2.cpp index 249fbf2de46..e8791f43f56 100644 --- a/Snap_rounding_2/test/Snap_rounding_2/test_float_snap_rounding_2.cpp +++ b/Snap_rounding_2/test/Snap_rounding_2/test_float_snap_rounding_2.cpp @@ -4,6 +4,8 @@ #include #include + +#include #include #include #include @@ -38,12 +40,12 @@ typedef CGAL::Cartesian_converter IK_to_EK; #ifdef COMPARE_WITH_INTEGER_SNAP_ROUNDING_2 #include #include -typedef CGAL::Snap_rounding_traits_2 SnapTraits; +typedef CGAL::Snap_rounding_traits_2 SnapTraits; #endif //Biggest double with ulp smaller than an integer -constexpr double maxFloat = std::pow(2,23); -constexpr double maxDouble = std::pow(2,52); +const double maxFloat = std::pow(2,23); +const double maxDouble = std::pow(2,52); Segment_2 scale_segment(const Segment_2 &s){ auto scale_value=[](FT x){ From ec71e0189f17d6c2e5f61c5788614eb692931172 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?L=C3=A9o=20Valque?= Date: Thu, 6 Nov 2025 16:31:52 +0100 Subject: [PATCH 30/41] Fix documentation --- .../Concepts/FloatSnapRoundingTraits_2.h | 12 ++++++------ Snap_rounding_2/doc/Snap_rounding_2/Doxyfile.in | 2 ++ .../doc/Snap_rounding_2/PackageDescription.txt | 2 -- Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h | 6 +++--- .../include/CGAL/Float_snap_rounding_traits_2.h | 2 +- 5 files changed, 12 insertions(+), 12 deletions(-) diff --git a/Snap_rounding_2/doc/Snap_rounding_2/Concepts/FloatSnapRoundingTraits_2.h b/Snap_rounding_2/doc/Snap_rounding_2/Concepts/FloatSnapRoundingTraits_2.h index 2f971628cc6..ec76d7e5e2b 100644 --- a/Snap_rounding_2/doc/Snap_rounding_2/Concepts/FloatSnapRoundingTraits_2.h +++ b/Snap_rounding_2/doc/Snap_rounding_2/Concepts/FloatSnapRoundingTraits_2.h @@ -1,11 +1,11 @@ /*! -\ingroup PkgFloatSnapRounding2Concepts +\ingroup PkgSnapRounding2Concepts \cgalConcept The concept `FloatSnapRoundingTraits_2` lists the set of requirements that must be fulfilled by an instance of the `Traits` template-parameter of -the free functions \ref `CGAL::snap_polygons_2()`, `CGAL::compute_snapped_subcurves_2()` and `CGAL::double_snap_rounding_2()`. +the free functions \ref CGAL::snap_polygons_2(), \ref CGAL::compute_snapped_subcurves_2() and \ref CGAL::double_snap_rounding_2(). The list includes the nested types of the geometric primitives used in this class and some function object types for the required predicates on those primitives. @@ -142,7 +142,7 @@ Squared_round_bound_2 squared_round_bound_2_object(); namespace FSRTraits_2{ /*! - \ingroup PkgFloatSnapRounding2Concepts + \ingroup PkgSnapRounding2Concepts \cgalConcept \cgalHasModelsBegin \cgalHasModelsBare{\link FloatSnapRoundingTraits_2::ConstructRoundPoint_2 `Float_snap_rounding_traits_2::Construct_round_point_2` \endlink} @@ -159,7 +159,7 @@ class ConstructRoundPoint_2 }; /*! - \ingroup PkgFloatSnapRounding2Concepts + \ingroup PkgSnapRounding2Concepts \cgalConcept \cgalHasModelsBegin \cgalHasModelsBare{\link FloatSnapRoundingTraits_2::SquaredRoundBound_2 `Float_snap_rounding_traits_2::Squared_round_bound_2` \endlink} @@ -176,7 +176,7 @@ class SquaredRoundBound_2 }; /*! - \ingroup PkgFloatSnapRounding2Concepts + \ingroup PkgSnapRounding2Concepts \cgalConcept \cgalHasModelsBegin \cgalHasModelsBare{\link FloatSnapRoundingTraits_2::ConverterToExact `Float_snap_rounding_traits_2::Convertex_to_exact` \endlink} @@ -195,7 +195,7 @@ class ConverterToExact }; /*! - \ingroup PkgFloatSnapRounding2Concepts + \ingroup PkgSnapRounding2Concepts \cgalConcept \cgalHasModelsBegin \cgalHasModelsBare{\link FloatSnapRoundingTraits_2::ConverterToExact `Float_snap_rounding_traits_2::Convertex_to_exact` \endlink} diff --git a/Snap_rounding_2/doc/Snap_rounding_2/Doxyfile.in b/Snap_rounding_2/doc/Snap_rounding_2/Doxyfile.in index 7c99d29c3b1..2d03578eb31 100644 --- a/Snap_rounding_2/doc/Snap_rounding_2/Doxyfile.in +++ b/Snap_rounding_2/doc/Snap_rounding_2/Doxyfile.in @@ -1,3 +1,5 @@ @INCLUDE = ${CGAL_DOC_PACKAGE_DEFAULTS} PROJECT_NAME = "CGAL ${CGAL_DOC_VERSION} - 2D Snap Rounding" +INPUT += ${CGAL_PACKAGE_INCLUDE_DIR}/CGAL/Float_snap_rounding_2.h \ + ${CGAL_PACKAGE_INCLUDE_DIR}/CGAL/Float_snap_rounding_traits_2.h diff --git a/Snap_rounding_2/doc/Snap_rounding_2/PackageDescription.txt b/Snap_rounding_2/doc/Snap_rounding_2/PackageDescription.txt index 94641a384c2..9426f5b4808 100644 --- a/Snap_rounding_2/doc/Snap_rounding_2/PackageDescription.txt +++ b/Snap_rounding_2/doc/Snap_rounding_2/PackageDescription.txt @@ -1,7 +1,5 @@ /// \defgroup PkgSnapRounding2Ref Reference Manual /// \defgroup PkgSnapRounding2Concepts Concepts -/// \defgroup PkgFloatSnapRounding2Ref Reference Manual* -/// \defgroup PkgFloatSnapRounding2Concepts Concepts* /// \ingroup PkgSnapRounding2Ref /*! \addtogroup PkgSnapRounding2Ref diff --git a/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h b/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h index e0468e25858..72dd129ec49 100644 --- a/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h +++ b/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h @@ -404,7 +404,7 @@ void double_snap_rounding_2_disjoint(PointsRange &pts, PolylineRange &polylines, } // end of namespace internal /** -* ingroup +* ingroup PkgSnapRounding2Ref * * Given a range of segments, compute rounded subsegments that are pairwise disjoint in their interiors, based on the input curves. * The output is a range of polyline with each polyline corresponding to an input segment. @@ -534,7 +534,7 @@ typename OutputContainer::iterator double_snap_rounding_2(InputIterator begin, } /** -* ingroup +* ingroup PkgSnapRounding2Ref * * Given a range of segments, compute rounded subsegments that are pairwise disjoint in their interior, as induced by the input curves. * @@ -666,7 +666,7 @@ OutputIterator compute_snapped_subcurves_2(InputIterator begin, } /** -* ingroup +* ingroup PkgSnapRounding2Ref * * Given a range of `Polygon_2`, compute rounded polygons such that their segments are either equal either disjoint in their interior, as induced by the input polygons. * The polygons are intended to be non-intersecting, unless the named parameter `compute_intersection` is set to `true`. diff --git a/Snap_rounding_2/include/CGAL/Float_snap_rounding_traits_2.h b/Snap_rounding_2/include/CGAL/Float_snap_rounding_traits_2.h index ab48747ae32..1d61398de05 100644 --- a/Snap_rounding_2/include/CGAL/Float_snap_rounding_traits_2.h +++ b/Snap_rounding_2/include/CGAL/Float_snap_rounding_traits_2.h @@ -25,7 +25,7 @@ namespace CGAL { /*! -\ingroup PkgFloatSnapRounding2Ref +\ingroup PkgSnapRounding2Ref The class `Float_snap_rounding_traits_2` is a model of the `FloatSnapRoundingTraits_2` concept, and is the only traits class supplied From cb5893ae92554273a6724a0ceb713f720a6a93cc Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?L=C3=A9o=20Valque?= Date: Thu, 6 Nov 2025 16:59:25 +0100 Subject: [PATCH 31/41] Fix documentation --- .../Concepts/FloatSnapRoundingTraits_2.h | 20 ++++++++++++++----- .../doc/Snap_rounding_2/Doxyfile.in | 6 ++++++ 2 files changed, 21 insertions(+), 5 deletions(-) diff --git a/Snap_rounding_2/doc/Snap_rounding_2/Concepts/FloatSnapRoundingTraits_2.h b/Snap_rounding_2/doc/Snap_rounding_2/Concepts/FloatSnapRoundingTraits_2.h index ec76d7e5e2b..e067d0e6d18 100644 --- a/Snap_rounding_2/doc/Snap_rounding_2/Concepts/FloatSnapRoundingTraits_2.h +++ b/Snap_rounding_2/doc/Snap_rounding_2/Concepts/FloatSnapRoundingTraits_2.h @@ -5,7 +5,7 @@ The concept `FloatSnapRoundingTraits_2` lists the set of requirements that must be fulfilled by an instance of the `Traits` template-parameter of -the free functions \ref CGAL::snap_polygons_2(), \ref CGAL::compute_snapped_subcurves_2() and \ref CGAL::double_snap_rounding_2(). +the free functions of TODO. The list includes the nested types of the geometric primitives used in this class and some function object types for the required predicates on those primitives. @@ -133,6 +133,16 @@ Construct_round_point_2 construct_round_point_2_object(); */ Squared_round_bound_2 squared_round_bound_2_object(); +/*! + +*/ +Converter_to_exact converter_to_exact_object(); + +/*! + +*/ +Converter_from_exact converter_from_exact_object(); + /// @} @@ -145,7 +155,7 @@ namespace FSRTraits_2{ \ingroup PkgSnapRounding2Concepts \cgalConcept \cgalHasModelsBegin - \cgalHasModelsBare{\link FloatSnapRoundingTraits_2::ConstructRoundPoint_2 `Float_snap_rounding_traits_2::Construct_round_point_2` \endlink} + \cgalHasModelsBare{\link FloatSnapRoundingTraits_2::Construct_round_point_2 `Float_snap_rounding_traits_2::Construct_round_point_2` \endlink} \cgalHasModelsEnd */ class ConstructRoundPoint_2 @@ -162,7 +172,7 @@ class ConstructRoundPoint_2 \ingroup PkgSnapRounding2Concepts \cgalConcept \cgalHasModelsBegin - \cgalHasModelsBare{\link FloatSnapRoundingTraits_2::SquaredRoundBound_2 `Float_snap_rounding_traits_2::Squared_round_bound_2` \endlink} + \cgalHasModelsBare{\link FloatSnapRoundingTraits_2::Squared_round_bound_2 `Float_snap_rounding_traits_2::Squared_round_bound_2` \endlink} \cgalHasModelsEnd */ class SquaredRoundBound_2 @@ -179,7 +189,7 @@ class SquaredRoundBound_2 \ingroup PkgSnapRounding2Concepts \cgalConcept \cgalHasModelsBegin - \cgalHasModelsBare{\link FloatSnapRoundingTraits_2::ConverterToExact `Float_snap_rounding_traits_2::Convertex_to_exact` \endlink} + \cgalHasModelsBare{\link FloatSnapRoundingTraits_2::Converter_to_exact `Float_snap_rounding_traits_2::Convertex_to_exact` \endlink} \cgalHasModelsEnd */ class ConverterToExact @@ -198,7 +208,7 @@ class ConverterToExact \ingroup PkgSnapRounding2Concepts \cgalConcept \cgalHasModelsBegin - \cgalHasModelsBare{\link FloatSnapRoundingTraits_2::ConverterToExact `Float_snap_rounding_traits_2::Convertex_to_exact` \endlink} + \cgalHasModelsBare{\link FloatSnapRoundingTraits_2::Converter_from_exact `Float_snap_rounding_traits_2::Convertex_to_exact` \endlink} \cgalHasModelsEnd */ class ConverterFromExact diff --git a/Snap_rounding_2/doc/Snap_rounding_2/Doxyfile.in b/Snap_rounding_2/doc/Snap_rounding_2/Doxyfile.in index 2d03578eb31..0e85668587a 100644 --- a/Snap_rounding_2/doc/Snap_rounding_2/Doxyfile.in +++ b/Snap_rounding_2/doc/Snap_rounding_2/Doxyfile.in @@ -1,5 +1,11 @@ @INCLUDE = ${CGAL_DOC_PACKAGE_DEFAULTS} PROJECT_NAME = "CGAL ${CGAL_DOC_VERSION} - 2D Snap Rounding" + +# custom options for this package +HIDE_UNDOC_CLASSES = true +WARN_IF_UNDOCUMENTED = false + + INPUT += ${CGAL_PACKAGE_INCLUDE_DIR}/CGAL/Float_snap_rounding_2.h \ ${CGAL_PACKAGE_INCLUDE_DIR}/CGAL/Float_snap_rounding_traits_2.h From 92807e94f77eb3c8417dce05059cdfc3a2916f31 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?L=C3=A9o=20Valque?= Date: Thu, 6 Nov 2025 17:50:52 +0100 Subject: [PATCH 32/41] Add a test with differents Kernels --- .../test_float_snap_rounding_2_concept.cpp | 43 +++++++++++++++++++ 1 file changed, 43 insertions(+) create mode 100644 Snap_rounding_2/test/Snap_rounding_2/test_float_snap_rounding_2_concept.cpp diff --git a/Snap_rounding_2/test/Snap_rounding_2/test_float_snap_rounding_2_concept.cpp b/Snap_rounding_2/test/Snap_rounding_2/test_float_snap_rounding_2_concept.cpp new file mode 100644 index 00000000000..903633f53c4 --- /dev/null +++ b/Snap_rounding_2/test/Snap_rounding_2/test_float_snap_rounding_2_concept.cpp @@ -0,0 +1,43 @@ +// #define DOUBLE_2D_SNAP_VERBOSE +// #define DOUBLE_2D_SNAP_FULL_VERBOSE +#define BENCH_AND_VERBOSE_FLOAT_SNAP_ROUNDING_2 + +#include +#include +#include + +#include + +typedef CGAL::Cartesian Rational_Kernel; +typedef CGAL::Float_snap_rounding_traits_2 Rational_Traits; + +struct ConceptTraits{ + +}; + +template +struct Test{ + typedef typename Traits::Segment_2 Segment_2; + typedef typename Traits::Point_2 Point_2; + + void fix_test(){ + std::vector segs; + segs.emplace_back(Point_2(1, 1), Point_2(-1, -1)); + segs.emplace_back(Point_2(0, 0), Point_2(1, -1)); + segs.emplace_back(Point_2(0, 2), Point_2(2, 0)); + segs.emplace_back(Point_2(0, 2), Point_2(-2, -4)); + segs.emplace_back(Point_2(-2, 2), Point_2(2, 2)); + segs.emplace_back(Point_2(5, 7), Point_2(9, 7)); + segs.emplace_back(Point_2(1,1), Point_2(3,1)); + segs.emplace_back(Point_2(1,2), Point_2(3,0)); + std::vector out; + CGAL::compute_snapped_subcurves_2(segs.begin(), segs.end(), std::back_inserter(out), CGAL::parameters::geom_traits(Traits())); + assert(!CGAL::do_curves_intersect(out.begin(), out.end())); + } +}; + +int main(/*int argc,char *argv[]*/) +{ + Test().fix_test(); + return(0); +} From 803ef3593b9919b34bc4ea9c5d1368176bac8762 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?L=C3=A9o=20Valque?= Date: Fri, 7 Nov 2025 09:54:44 +0100 Subject: [PATCH 33/41] Fix the example --- .../examples/Snap_rounding_2/float_snap_rounding.cpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/Snap_rounding_2/examples/Snap_rounding_2/float_snap_rounding.cpp b/Snap_rounding_2/examples/Snap_rounding_2/float_snap_rounding.cpp index 0a7a58d5063..9ac611ca936 100644 --- a/Snap_rounding_2/examples/Snap_rounding_2/float_snap_rounding.cpp +++ b/Snap_rounding_2/examples/Snap_rounding_2/float_snap_rounding.cpp @@ -44,7 +44,7 @@ int main(int argc, char *argv[]) std::cout << "Computes the intersections and snaps the segments" << std::endl; std::vector< Segment_2> out; - CGAL::compute_snapped_subcurves_2(segs.begin(), segs.end(), out); + CGAL::compute_snapped_subcurves_2(segs.begin(), segs.end(), std::back_inserter(out)); std::cout << "Does the output intersect: " << CGAL::do_curves_intersect(out.begin(), out.end()) << std::endl; std::cout << "Size of the output: " << out.size() << std::endl; From 19c832b62628ca80acb4b81ff2ca6cbef87b04d6 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?L=C3=A9o=20Valque?= Date: Mon, 10 Nov 2025 15:12:13 +0100 Subject: [PATCH 34/41] Some additional test --- .../include/CGAL/Float_snap_rounding_2.h | 12 ++-- .../test_float_snap_rounding_2.cpp | 69 ++++++++----------- .../test_float_snap_rounding_2_concept.cpp | 6 +- .../examples/Surface_sweep_2/plane_sweep.cpp | 51 +++++++++++--- 4 files changed, 78 insertions(+), 60 deletions(-) diff --git a/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h b/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h index 72dd129ec49..99c1d93ce57 100644 --- a/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h +++ b/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h @@ -701,10 +701,10 @@ OutputIterator compute_snapped_subcurves_2(InputIterator begin, * @warning The convex property of the polygons is not necessarly preserved */ template -void snap_polygons_2(InputIterator begin, - InputIterator end, - OutputIterator out, - const NamedParameters &np = parameters::default_values()) +void compute_snapped_polygons_2(InputIterator begin, + InputIterator end, + OutputIterator out, + const NamedParameters &np = parameters::default_values()) { using Concurrency_tag = typename internal_np::Lookup_named_param_def -void snap_polygon_2(const Polygon_2 &P, Polygon_2 &out, const NamedParameters &np = parameters::default_values()) +void compute_snapped_polygon_2(const Polygon_2 &P, Polygon_2 &out, const NamedParameters &np = parameters::default_values()) { std::array vec({P}); std::vector out_vec; - snap_polygons_2(vec.begin(), vec.end(), std::back_inserter(out_vec), np); + compute_snapped_polygons_2(vec.begin(), vec.end(), std::back_inserter(out_vec), np); out = out_vec[0]; } diff --git a/Snap_rounding_2/test/Snap_rounding_2/test_float_snap_rounding_2.cpp b/Snap_rounding_2/test/Snap_rounding_2/test_float_snap_rounding_2.cpp index e8791f43f56..f281e16f7af 100644 --- a/Snap_rounding_2/test/Snap_rounding_2/test_float_snap_rounding_2.cpp +++ b/Snap_rounding_2/test/Snap_rounding_2/test_float_snap_rounding_2.cpp @@ -17,7 +17,7 @@ #include #include -typedef CGAL::Exact_predicates_exact_constructions_kernel Kernel; +typedef CGAL::Exact_predicates_inexact_constructions_kernel Kernel; typedef CGAL::Float_snap_rounding_traits_2 Traits_2; typedef Kernel::Segment_2 Segment_2; typedef Kernel::Point_2 Point_2; @@ -103,26 +103,6 @@ void test(const std::vector &segs){ std::cout << "Formal snap size: " << out.size() << " ,running time: " << t.time() << std::endl; #endif assert(!CGAL::do_curves_intersect(out.begin(), out.end())); - if(CGAL::do_curves_intersect(out.begin(), out.end())){ - std::cout << "INTERSECTING OUTPUT!!" << std::endl; - std::vector bad_pts; - std::vector out2; - for(auto seg: out) - out2.emplace_back(seg.source(), seg.target()); - CGAL::compute_intersection_points(out2.begin(), out2.end(), std::back_inserter(bad_pts)); - std::cout << bad_pts[0] << std::endl; - for(size_t i=0; i polygons; + for(size_t i=0; i!=nb_polygons; ++i){ + Polygon_2 poly; + for(size_t j=0; j!=nb_polygons; ++j) + poly.push_back(Point_2(random_point(r), random_point(r))); + polygons.emplace_back(poly); + } + + std::vector out; + CGAL::compute_snapped_polygons_2(polygons.begin(), polygons.end(), std::back_inserter(out)); + + std::vector segs; + for(const Polygon_2 &poly: out) + for(size_t i=1; i segs; segs.emplace_back(Point_2(0.99999999999992173,-0.9999999999991368),Point_2(0.99999999999987266,-6.016984285966173e-13)); segs.emplace_back(Point_2(1.000000000000494,-1.0000000000002354),Point_2(0.99999999999950062,8.0391743541535603e-13)); segs.emplace_back(Point_2(1.0000000000008489,-0.99999999999951794),Point_2(0.99999999999926925,-2.3642280455149055e-13)); - test(segs); -} + segs.clear(); -void big_test_2(){ - std::cout << "Fix tests" << std::endl; - std::vector< Segment_2 > segs; segs.emplace_back(Point_2(1.0000000000008309,9.0061990813884068e-13), Point_2(0.99999999999981259,0.99999999999938394)); segs.emplace_back(Point_2(1.0000000000005094,-3.1951552372595695e-13), Point_2(1.0000000000004887,1.0000000000006302)); segs.emplace_back(Point_2(1.0000000000006171,-7.3034227828590228e-13), Point_2(1.0000000000004181,0.99999999999913269)); @@ -271,9 +262,8 @@ void big_test_2(){ segs.emplace_back(Point_2(1.0000000000006153,-2.798953111315494e-13), Point_2(1.0000000000003044,1.0000000000008289)); segs.emplace_back(Point_2(1.0000000000008387,5.8214218612051864e-13), Point_2(1.0000000000002918,0.99999999999943212)); segs.emplace_back(Point_2(1.0000000000007745,-3.9231414307222458e-13), Point_2(1.000000000000294,1.0000000000003408)); - test(segs); - + segs.clear(); } void test_float_snap_rounding(){ @@ -300,13 +290,12 @@ int main(int argc,char *argv[]) CGAL::Random r(argc==1?rp.get_seed():std::stoi(argv[1])); std::cout << "random seed = " << r.get_seed() << std::endl; std::cout << std::setprecision(17); - big_test(); - big_test_2(); // test_almost_indentical_segments(r, 50, Vector_2(1,1), Vector_2(-1,-1)); fix_test(); test_float_snap_rounding(); test_fully_random(r,1000); test_multi_almost_indentical_segments(r,200); + test_random_polygons(r,200,10); // test_iterative_square_intersection(r,2000); return(0); } diff --git a/Snap_rounding_2/test/Snap_rounding_2/test_float_snap_rounding_2_concept.cpp b/Snap_rounding_2/test/Snap_rounding_2/test_float_snap_rounding_2_concept.cpp index 903633f53c4..a0873feb600 100644 --- a/Snap_rounding_2/test/Snap_rounding_2/test_float_snap_rounding_2_concept.cpp +++ b/Snap_rounding_2/test/Snap_rounding_2/test_float_snap_rounding_2_concept.cpp @@ -6,15 +6,13 @@ #include #include +#include #include +typedef CGAL::Exact_predicates_exact_constructions_kernel Epeck; typedef CGAL::Cartesian Rational_Kernel; typedef CGAL::Float_snap_rounding_traits_2 Rational_Traits; -struct ConceptTraits{ - -}; - template struct Test{ typedef typename Traits::Segment_2 Segment_2; diff --git a/Surface_sweep_2/examples/Surface_sweep_2/plane_sweep.cpp b/Surface_sweep_2/examples/Surface_sweep_2/plane_sweep.cpp index 28f894cd8bf..1493b1c5f63 100644 --- a/Surface_sweep_2/examples/Surface_sweep_2/plane_sweep.cpp +++ b/Surface_sweep_2/examples/Surface_sweep_2/plane_sweep.cpp @@ -6,25 +6,56 @@ #include #include +#include #include typedef CGAL::Exact_predicates_exact_constructions_kernel Kernel; typedef Kernel::Point_2 Point_2; -typedef CGAL::Arr_segment_traits_2 Traits_2; -typedef Traits_2::Curve_2 Segment_2; +typedef Kernel::Segment_2 Segment; +typedef Kernel::Point_2 Point; +using Segment_traits = CGAL::Arr_segment_traits_2; +using Traits = CGAL::Arr_polyline_traits_2; +using My_polyline = Traits::Curve_2; int main() { + Traits traits; // Construct the input segments. - Segment_2 segments[] = {Segment_2 (Point_2 (1, 5), Point_2 (8, 5)), - Segment_2 (Point_2 (1, 1), Point_2 (8, 8)), - Segment_2 (Point_2 (3, 1), Point_2 (3, 8)), - Segment_2 (Point_2 (8, 5), Point_2 (8, 8))}; + auto polyline_construct = traits.construct_curve_2_object(); + + Point points1[5]; + points1[0] = Point(0, 0); + points1[1] = Point(2, 4); + points1[2] = Point(3, 0); + points1[3] = Point(4, 4); + points1[4] = Point(6, 0); + auto pi1 = polyline_construct(&points1[0], &points1[5]); + + std::list points2; + points2.push_back(Point(1, 3)); + points2.push_back(Point(0, 2)); + points2.push_back(Point(1, 0)); + points2.push_back(Point(2, 1)); + points2.push_back(Point(3, 0)); + points2.push_back(Point(4, 1)); + points2.push_back(Point(5, 0)); + points2.push_back(Point(6, 2)); + points2.push_back(Point(5, 3)); + points2.push_back(Point(4, 2)); + auto pi2 = polyline_construct(points2.begin(), points2.end()); + + std::vector segs; + segs.push_back(Segment(Point(0, 2), Point(1, 2))); + segs.push_back(Segment(Point(1, 2), Point(3, 6))); + segs.push_back(Segment(Point(3, 6), Point(5, 2))); + auto pi3 = polyline_construct(segs.begin(), segs.end()); + + My_polyline vec[3]={pi1, pi2, pi3}; // Compute all intersection points. std::list pts; - CGAL::compute_intersection_points(segments, segments + 4, + CGAL::compute_intersection_points(vec, vec + 3, std::back_inserter(pts)); // Print the result. @@ -33,14 +64,14 @@ int main() std::ostream_iterator(std::cout, "\n")); // Compute the non-intersecting sub-segments induced by the input segments. - std::list sub_segs; + std::list sub_segs; - CGAL::compute_subcurves(segments, segments + 4, std::back_inserter(sub_segs)); + CGAL::compute_subcurves(vec, vec + 3, std::back_inserter(sub_segs)); std::cout << "Found " << sub_segs.size() << " interior-disjoint sub-segments." << std::endl; - assert(CGAL::do_curves_intersect (segments, segments + 4)); +// assert(CGAL::do_curves_intersect (segments, segments + 4)); return 0; } From 422338c9899f0bcdd188cef34bcf31279a4335cc Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?L=C3=A9o=20Valque?= Date: Mon, 10 Nov 2025 16:58:26 +0100 Subject: [PATCH 35/41] clean test --- .../include/CGAL/Float_snap_rounding_2.h | 12 ++++- .../test_float_snap_rounding_2.cpp | 54 +++++++++++++------ 2 files changed, 48 insertions(+), 18 deletions(-) diff --git a/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h b/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h index 99c1d93ce57..f2d2b45825b 100644 --- a/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h +++ b/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h @@ -600,7 +600,8 @@ OutputIterator compute_snapped_subcurves_2(InputIterator begin, std::vector convert_input; for(InputIterator it=begin; it!=end; ++it) - convert_input.push_back(to_exact(*it)); + if(it->source()!=it->target()) + convert_input.push_back(to_exact(*it)); std::vector segs; #ifdef DOUBLE_2D_SNAP_VERBOSE std::cout << "Solved intersections" << std::endl; @@ -717,6 +718,8 @@ void compute_snapped_polygons_2(InputIterator begin, NamedParameters, DefaultTraits>::type; + using Less_xy_2 = typename Traits::Less_xy_2; + using Point_2 = typename Traits::Point_2; using I2E = typename Traits::Converter_to_exact; using E2O = typename Traits::Converter_from_exact; @@ -729,6 +732,8 @@ void compute_snapped_polygons_2(InputIterator begin, const Traits &traits = choose_parameter(get_parameter(np, internal_np::geom_traits), Traits()); const bool compute_intersections = choose_parameter(get_parameter(np, internal_np::compute_intersection), false); + Less_xy_2 less_xy_2 = traits.less_xy_2_object(); + I2E to_exact=traits.converter_to_exact_object(); E2O from_exact=traits.converter_from_exact_object(); @@ -751,7 +756,10 @@ void compute_snapped_polygons_2(InputIterator begin, for(const typename Polygon_2::Point_2 &p: it->vertices()) pts.push_back(to_exact(p)); for(size_t i=0; i #include @@ -17,7 +13,7 @@ #include #include -typedef CGAL::Exact_predicates_inexact_constructions_kernel Kernel; +typedef CGAL::Exact_predicates_exact_constructions_kernel Kernel; typedef CGAL::Float_snap_rounding_traits_2 Traits_2; typedef Kernel::Segment_2 Segment_2; typedef Kernel::Point_2 Point_2; @@ -71,7 +67,8 @@ void compute_subcurves_and_naive_snap(Iterator begin, Iterator end, OutRange &ou IK_to_EK to_exact; std::vector segs; for(Iterator it=begin; it!=end; ++it) - segs.emplace_back(*it); + if(it->source()!=it->target()) + segs.emplace_back(*it); CGAL::compute_subcurves(segs.begin(), segs.end(), std::back_inserter(out)); //Naive snap OutRange out2; @@ -111,7 +108,9 @@ void test(const std::vector &segs){ t.stop(); std::cout << "Running time with integer 2D snap (scaled at 10^15): " << t.time() << std::endl; #endif +#ifdef BENCH_AND_VERBOSE_FLOAT_SNAP_ROUNDING_2 std::cout << "\n"; +#endif } void test_fully_random(CGAL::Random &r, size_t nb_segments){ @@ -124,11 +123,11 @@ void test_fully_random(CGAL::Random &r, size_t nb_segments){ } void test_random_polygons(CGAL::Random &r, size_t nb_polygons, size_t nb_pts){ - std::cout << "Test fully random" << std::endl; + std::cout << "Test random polygons" << std::endl; std::vector polygons; for(size_t i=0; i!=nb_polygons; ++i){ Polygon_2 poly; - for(size_t j=0; j!=nb_polygons; ++j) + for(size_t j=0; j!=nb_pts; ++j) poly.push_back(Point_2(random_point(r), random_point(r))); polygons.emplace_back(poly); } @@ -145,7 +144,9 @@ void test_random_polygons(CGAL::Random &r, size_t nb_polygons, size_t nb_pts){ void test_almost_indentical_segments(CGAL::Random &r, size_t nb_segments, Vector_2 source, Vector_2 target){ // Difficult test +#ifdef BENCH_AND_VERBOSE_FLOAT_SNAP_ROUNDING_2 std::cout << "Test with almost identical segments from " << source << " to " << target << std::endl; +#endif std::vector segs; double eps=std::pow(2,-40); for(size_t i=0; i!=nb_segments; ++i){ @@ -206,36 +207,56 @@ void test_multi_almost_indentical_segments(CGAL::Random &r, size_t nb_segments){ } void fix_test(){ +#ifdef BENCH_AND_VERBOSE_FLOAT_SNAP_ROUNDING_2 std::cout << "Fix tests" << std::endl; +#endif std::vector< Segment_2 > segs; FT e(std::pow(2, -60)); + // Disjoint segments segs.emplace_back(Point_2(0, 1), Point_2(1, 1)); segs.emplace_back(Point_2(1, 2), Point_2(2, 2)); test(segs); segs.clear(); - // segs.emplace_back(Point_2(0, 1), Point_2(2, 1)); - // segs.emplace_back(Point_2(1, 0), Point_2(1, 0)); - // test(segs); - // segs.clear(); + // Degenerate segment + segs.emplace_back(Point_2(0, 1), Point_2(2, 1)); + segs.emplace_back(Point_2(1, 0), Point_2(1, 0)); + test(segs); + segs.clear(); + // Collinear segments + segs.emplace_back(Point_2(0, 0), Point_2(4, 0)); + segs.emplace_back(Point_2(2, 0), Point_2(6, 0)); + test(segs); + segs.clear(); + + // Duplicate segment + segs.emplace_back(Point_2(0, 0), Point_2(4, 0)); + segs.emplace_back(Point_2(0, 0), Point_2(4, 0)); + test(segs); + segs.clear(); + + // Polyline segs.emplace_back(Point_2(0, 0), Point_2(2, 0)); segs.emplace_back(Point_2(2, 0), Point_2(1, 1)); segs.emplace_back(Point_2(1, 1), Point_2(3, 1)); test(segs); segs.clear(); + // Horizontal intersection when round segs.emplace_back(Point_2(0, 1+e), Point_2(2, 1)); segs.emplace_back(Point_2(1, 1), Point_2(1, 0)); test(segs); segs.clear(); + // Vertical intersection when round segs.emplace_back(Point_2(0, 0), Point_2(1-e, 0)); segs.emplace_back(Point_2(1, 1), Point_2(1, -1)); test(segs); segs.clear(); + // A more complex example segs.emplace_back(Point_2(1-e, 1), Point_2(-1-e, -1+2*e)); segs.emplace_back(Point_2(e/2, e/2), Point_2(1, -1)); segs.emplace_back(Point_2(0, 2-e/2), Point_2(2, 0)); @@ -245,6 +266,7 @@ void fix_test(){ segs.emplace_back(Point_2(5, 7-e), Point_2(9, 7-e)); test(segs); + // Minimal intersection subsets of a random examples segs.emplace_back(Point_2(0.99999999999992173,-0.9999999999991368),Point_2(0.99999999999987266,-6.016984285966173e-13)); segs.emplace_back(Point_2(1.000000000000494,-1.0000000000002354),Point_2(0.99999999999950062,8.0391743541535603e-13)); segs.emplace_back(Point_2(1.0000000000008489,-0.99999999999951794),Point_2(0.99999999999926925,-2.3642280455149055e-13)); @@ -267,7 +289,9 @@ void fix_test(){ } void test_float_snap_rounding(){ - std::cout << "Fix tests" << std::endl; +#ifdef BENCH_AND_VERBOSE_FLOAT_SNAP_ROUNDING_2 + std::cout << "Test the other API" << std::endl; +#endif std::vector< Segment_2 > segs; std::vector< Polyline_2 > out; @@ -281,7 +305,6 @@ void test_float_snap_rounding(){ segs.emplace_back(Point_2(5, 7-e), Point_2(9, 7-e)); double_snap_rounding_2(segs.begin(), segs.end(), out); - // assert(out.size()==segs.size()); } int main(int argc,char *argv[]) @@ -290,12 +313,11 @@ int main(int argc,char *argv[]) CGAL::Random r(argc==1?rp.get_seed():std::stoi(argv[1])); std::cout << "random seed = " << r.get_seed() << std::endl; std::cout << std::setprecision(17); - // test_almost_indentical_segments(r, 50, Vector_2(1,1), Vector_2(-1,-1)); fix_test(); test_float_snap_rounding(); test_fully_random(r,1000); test_multi_almost_indentical_segments(r,200); - test_random_polygons(r,200,10); + // test_random_polygons(r,200,10); // test_iterative_square_intersection(r,2000); return(0); } From e63869e2b32ea4a7b484d5aa92a6e323ec3dbe83 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?L=C3=A9o=20Valque?= Date: Wed, 12 Nov 2025 10:29:11 +0100 Subject: [PATCH 36/41] Solved invalid iterator after erased with MSVC --- Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h | 2 +- .../test/Snap_rounding_2/test_float_snap_rounding_2.cpp | 4 ++++ 2 files changed, 5 insertions(+), 1 deletion(-) diff --git a/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h b/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h index f2d2b45825b..7129a3f109c 100644 --- a/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h +++ b/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h @@ -165,7 +165,7 @@ void snap_rounding_scan(PointsRange &pts, PolylineRange &polylines, const Traits pos_it = std::lower_bound(y_order.begin(), y_order.end(), std::pair(event.pi, event.li), pi_below_li); if(event.type==EVENT_TYPE::REMOVE){ assert(*pos_it==event.li); - y_order.erase(pos_it); + pos_it = y_order.erase(pos_it); } if(!y_order.empty()){ diff --git a/Snap_rounding_2/test/Snap_rounding_2/test_float_snap_rounding_2.cpp b/Snap_rounding_2/test/Snap_rounding_2/test_float_snap_rounding_2.cpp index ec673be283b..13893cc23c8 100644 --- a/Snap_rounding_2/test/Snap_rounding_2/test_float_snap_rounding_2.cpp +++ b/Snap_rounding_2/test/Snap_rounding_2/test_float_snap_rounding_2.cpp @@ -114,7 +114,9 @@ void test(const std::vector &segs){ } void test_fully_random(CGAL::Random &r, size_t nb_segments){ +#ifdef BENCH_AND_VERBOSE_FLOAT_SNAP_ROUNDING_2 std::cout << "Test fully random" << std::endl; +#endif std::vector segs; for(size_t i=0; i!=nb_segments; ++i) segs.emplace_back(random_point(r), random_point(r)); @@ -123,7 +125,9 @@ void test_fully_random(CGAL::Random &r, size_t nb_segments){ } void test_random_polygons(CGAL::Random &r, size_t nb_polygons, size_t nb_pts){ +#ifdef BENCH_AND_VERBOSE_FLOAT_SNAP_ROUNDING_2 std::cout << "Test random polygons" << std::endl; +#endif std::vector polygons; for(size_t i=0; i!=nb_polygons; ++i){ Polygon_2 poly; From b0388e4b2c8e89e73e3b489335f31c8799c26e69 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?L=C3=A9o=20Valque?= Date: Tue, 18 Nov 2025 09:41:14 +0100 Subject: [PATCH 37/41] Fix dangling reference warning of the testsuite --- .../Snap_rounding_2/test_float_snap_rounding_2_concept.cpp | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/Snap_rounding_2/test/Snap_rounding_2/test_float_snap_rounding_2_concept.cpp b/Snap_rounding_2/test/Snap_rounding_2/test_float_snap_rounding_2_concept.cpp index a0873feb600..518f8d405e4 100644 --- a/Snap_rounding_2/test/Snap_rounding_2/test_float_snap_rounding_2_concept.cpp +++ b/Snap_rounding_2/test/Snap_rounding_2/test_float_snap_rounding_2_concept.cpp @@ -29,7 +29,8 @@ struct Test{ segs.emplace_back(Point_2(1,1), Point_2(3,1)); segs.emplace_back(Point_2(1,2), Point_2(3,0)); std::vector out; - CGAL::compute_snapped_subcurves_2(segs.begin(), segs.end(), std::back_inserter(out), CGAL::parameters::geom_traits(Traits())); + Traits traits; + CGAL::compute_snapped_subcurves_2(segs.begin(), segs.end(), std::back_inserter(out), CGAL::parameters::geom_traits(traits)); assert(!CGAL::do_curves_intersect(out.begin(), out.end())); } }; From 4a2ebc50f8ae9113398dc93f955002aa68877abf Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?S=C3=A9bastien=20Loriot?= Date: Mon, 17 Nov 2025 16:58:38 +0100 Subject: [PATCH 38/41] workaround issue with MSVC2017 --- Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h | 12 +++++++----- 1 file changed, 7 insertions(+), 5 deletions(-) diff --git a/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h b/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h index 7129a3f109c..6122970d19a 100644 --- a/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h +++ b/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h @@ -185,13 +185,15 @@ void snap_rounding_scan(PointsRange &pts, PolylineRange &polylines, const Traits if(possibly(csq_dist_2(p, seg, bound)!=CGAL::LARGER)) { - if constexpr(std::is_same_v){ + if (std::is_same_v) + { + internal::Evaluate evaluate; // We refine the pts to reduce the rounding shift and check again - pts[pi].exact(); - pts[pl[pl.size()-2]].exact(); + evaluate(pts[pi]); + evaluate(pts[pl[pl.size()-2]]); // The two following call of exact act on seg variables since they appear in its DAG - pts[pl.back()].exact(); - pts[pl.front()].exact(); + evaluate(pts[pl.back()]); + evaluate(pts[pl.front()]); // Update the bounds round_bound_pts[pi]=round_bound(pts[pi]); round_bound_pts[pl[pl.size()-2]]=round_bound(pts[pl[pl.size()-2]]); From e7fe2e73a4d76fa2beae418795c84eb1d8e0e001 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?S=C3=A9bastien=20Loriot?= Date: Mon, 17 Nov 2025 18:59:57 +0100 Subject: [PATCH 39/41] push the constexpr back --- Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h b/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h index 6122970d19a..66d1630fdcc 100644 --- a/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h +++ b/Snap_rounding_2/include/CGAL/Float_snap_rounding_2.h @@ -185,7 +185,7 @@ void snap_rounding_scan(PointsRange &pts, PolylineRange &polylines, const Traits if(possibly(csq_dist_2(p, seg, bound)!=CGAL::LARGER)) { - if (std::is_same_v) + if constexpr (std::is_same_v) { internal::Evaluate evaluate; // We refine the pts to reduce the rounding shift and check again From 4c5b9174e5408b14dbabd6406ce0f535c1b8f71b Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?L=C3=A9o=20Valque?= Date: Thu, 20 Nov 2025 11:46:55 +0100 Subject: [PATCH 40/41] Use Simple_cartesian instead of Cartesian --- .../Snap_rounding_2/test_float_snap_rounding_2.cpp | 14 +++++++++----- 1 file changed, 9 insertions(+), 5 deletions(-) diff --git a/Snap_rounding_2/test/Snap_rounding_2/test_float_snap_rounding_2.cpp b/Snap_rounding_2/test/Snap_rounding_2/test_float_snap_rounding_2.cpp index 13893cc23c8..0596c647917 100644 --- a/Snap_rounding_2/test/Snap_rounding_2/test_float_snap_rounding_2.cpp +++ b/Snap_rounding_2/test/Snap_rounding_2/test_float_snap_rounding_2.cpp @@ -8,12 +8,9 @@ #include -#include - #include -#include -typedef CGAL::Exact_predicates_exact_constructions_kernel Kernel; +typedef CGAL::Exact_predicates_inexact_constructions_kernel Kernel; typedef CGAL::Float_snap_rounding_traits_2 Traits_2; typedef Kernel::Segment_2 Segment_2; typedef Kernel::Point_2 Point_2; @@ -29,9 +26,14 @@ typedef CGAL::Polygon_2 Polygon_2; typedef CGAL::Polygon_with_holes_2 Polygon_with_holes_2; typedef std::vector Pwh_vec_2; -typedef CGAL::Cartesian Naive; +#ifdef BENCH_AND_VERBOSE_FLOAT_SNAP_ROUNDING_2 +#include +#include +#include +typedef CGAL::Simple_cartesian Naive; typedef CGAL::Cartesian_converter EK_to_IK; typedef CGAL::Cartesian_converter IK_to_EK; +#endif #ifdef COMPARE_WITH_INTEGER_SNAP_ROUNDING_2 #include @@ -61,6 +63,7 @@ Point_2 random_noise_point(CGAL::Random& r, double m, double M, Vector_2 v) return Point_2(FT(r.get_double(m, M) + CGAL::to_double(v.x())), FT(r.get_double(m, M) + CGAL::to_double(v.y()))); } +#ifdef BENCH_AND_VERBOSE_FLOAT_SNAP_ROUNDING_2 template void compute_subcurves_and_naive_snap(Iterator begin, Iterator end, OutRange &out){ EK_to_IK to_inexact; @@ -79,6 +82,7 @@ void compute_subcurves_and_naive_snap(Iterator begin, Iterator end, OutRange &ou } std::swap(out,out2); } +#endif void test(const std::vector &segs){ std::vector out; From eec0ae3bf5c32bcb6208dbe7de78f543d03ebb8a Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?L=C3=A9o=20Valque?= Date: Thu, 20 Nov 2025 15:14:53 +0100 Subject: [PATCH 41/41] Forget that do_curves_intersect does not work with EPICK --- .../test/Snap_rounding_2/test_float_snap_rounding_2.cpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/Snap_rounding_2/test/Snap_rounding_2/test_float_snap_rounding_2.cpp b/Snap_rounding_2/test/Snap_rounding_2/test_float_snap_rounding_2.cpp index 0596c647917..421204fbde0 100644 --- a/Snap_rounding_2/test/Snap_rounding_2/test_float_snap_rounding_2.cpp +++ b/Snap_rounding_2/test/Snap_rounding_2/test_float_snap_rounding_2.cpp @@ -10,7 +10,7 @@ #include -typedef CGAL::Exact_predicates_inexact_constructions_kernel Kernel; +typedef CGAL::Exact_predicates_exact_constructions_kernel Kernel; typedef CGAL::Float_snap_rounding_traits_2 Traits_2; typedef Kernel::Segment_2 Segment_2; typedef Kernel::Point_2 Point_2;