mirror of https://github.com/CGAL/cgal
Merge eec0ae3bf5 into 26a5fc70e4
This commit is contained in:
commit
a08fec3149
|
|
@ -712,6 +712,11 @@ public:
|
|||
// Check if we have a single intersection point.
|
||||
const Point_2* ip = std::get_if<Point_2>(&*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));
|
||||
|
|
|
|||
|
|
@ -406,6 +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)
|
||||
|
|
|
|||
|
|
@ -0,0 +1,228 @@
|
|||
|
||||
/*!
|
||||
\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 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.
|
||||
|
||||
\cgalRefines{AosTraits_2}
|
||||
|
||||
\cgalHasModelsBegin
|
||||
\cgalHasModels{CGAL::Float_snap_rounding_traits_2<Kernel>}
|
||||
\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();
|
||||
|
||||
/*!
|
||||
|
||||
*/
|
||||
Converter_to_exact converter_to_exact_object();
|
||||
|
||||
/*!
|
||||
|
||||
*/
|
||||
Converter_from_exact converter_from_exact_object();
|
||||
|
||||
|
||||
/// @}
|
||||
|
||||
}; /* end FloatSnapRoundingTraits_2 */
|
||||
|
||||
|
||||
namespace FSRTraits_2{
|
||||
|
||||
/*!
|
||||
\ingroup PkgSnapRounding2Concepts
|
||||
\cgalConcept
|
||||
\cgalHasModelsBegin
|
||||
\cgalHasModelsBare{\link FloatSnapRoundingTraits_2::Construct_round_point_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 PkgSnapRounding2Concepts
|
||||
\cgalConcept
|
||||
\cgalHasModelsBegin
|
||||
\cgalHasModelsBare{\link FloatSnapRoundingTraits_2::Squared_round_bound_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 PkgSnapRounding2Concepts
|
||||
\cgalConcept
|
||||
\cgalHasModelsBegin
|
||||
\cgalHasModelsBare{\link FloatSnapRoundingTraits_2::Converter_to_exact `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 PkgSnapRounding2Concepts
|
||||
\cgalConcept
|
||||
\cgalHasModelsBegin
|
||||
\cgalHasModelsBare{\link FloatSnapRoundingTraits_2::Converter_from_exact `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 */
|
||||
|
|
@ -1,3 +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
|
||||
|
|
|
|||
|
|
@ -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()
|
||||
|
|
|
|||
|
|
@ -0,0 +1,68 @@
|
|||
#include <CGAL/Exact_predicates_exact_constructions_kernel.h>
|
||||
#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
|
||||
#include <CGAL/Float_snap_rounding_2.h>
|
||||
#include <CGAL/Arr_segment_traits_2.h>
|
||||
#include <CGAL/Random.h>
|
||||
|
||||
typedef CGAL::Exact_predicates_exact_constructions_kernel Kernel;
|
||||
typedef CGAL::Exact_predicates_inexact_constructions_kernel Epick;
|
||||
typedef Kernel::Segment_2 Segment_2;
|
||||
typedef Kernel::Point_2 Point_2;
|
||||
typedef Kernel::FT FT;
|
||||
typedef std::vector<Point_2 > Polyline_2;
|
||||
|
||||
int main(int argc, char *argv[])
|
||||
{
|
||||
std::vector< Segment_2 > segs;
|
||||
|
||||
if(argc>1){
|
||||
std::cout << "Read segments in " << argv[1] << std::endl;
|
||||
std::ifstream in(argv[1]);
|
||||
|
||||
int n;
|
||||
in >> n;
|
||||
segs.reserve(n);
|
||||
for (int i=0; i<n; ++i)
|
||||
{
|
||||
double x1,x2,y1,y2;
|
||||
in >> 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));
|
||||
|
||||
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 << "Computes the intersections and snaps the segments" << std::endl;
|
||||
std::vector< Segment_2> 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;
|
||||
|
||||
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;
|
||||
}
|
||||
|
|
@ -0,0 +1,830 @@
|
|||
// 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
|
||||
|
||||
#ifdef DOUBLE_2D_SNAP_VERBOSE
|
||||
#include <iostream>
|
||||
#endif
|
||||
|
||||
#include <CGAL/Surface_sweep_2_algorithms.h>
|
||||
|
||||
#include <CGAL/Float_snap_rounding_traits_2.h>
|
||||
|
||||
#include <CGAL/intersection_2.h>
|
||||
#include <set>
|
||||
#include <vector>
|
||||
|
||||
#include <CGAL/Named_function_parameters.h>
|
||||
|
||||
#include <CGAL/mutex.h>
|
||||
|
||||
namespace CGAL {
|
||||
|
||||
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 <class Concurrency_tag=Sequential_tag, class Traits, class PointsRange , class PolylineRange>
|
||||
void snap_rounding_scan(PointsRange &pts, PolylineRange &polylines, const Traits &traits){
|
||||
using Point_2 = typename Traits::Point_2;
|
||||
using Segment_2 = typename Traits::Segment_2;
|
||||
|
||||
using Polyline = std::remove_cv_t<typename std::iterator_traits<typename PolylineRange::iterator>::value_type>;
|
||||
|
||||
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;
|
||||
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();
|
||||
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();
|
||||
|
||||
// Precompute round bounds for points
|
||||
std::vector<double> round_bound_pts;
|
||||
round_bound_pts.reserve(pts.size());
|
||||
for(std::size_t i=0; i<pts.size(); ++i)
|
||||
round_bound_pts.emplace_back(round_bound(pts[i]));
|
||||
|
||||
enum class EVENT_TYPE{INSERT, MIDDLE, REMOVE};
|
||||
struct Event{
|
||||
Event(size_t pi_, size_t li_, EVENT_TYPE t_):pi(pi_),li(li_),type(t_){}
|
||||
size_t pi;
|
||||
size_t li;
|
||||
EVENT_TYPE type;
|
||||
};
|
||||
|
||||
auto event_comparator=[&](Event a, Event b){
|
||||
const Point_2 &pa=pts[a.pi];
|
||||
const Point_2 &pb=pts[b.pi];
|
||||
|
||||
// We sort event lexicographically
|
||||
if(a.pi!=b.pi)
|
||||
return less_xy_2(pa,pb);
|
||||
|
||||
// If two events use the same point, we place remove event first for a smaller y_order range
|
||||
if(a.type != b.type)
|
||||
return a.type > b.type;
|
||||
// Only to have a complete order
|
||||
return a.li < b.li;
|
||||
};
|
||||
|
||||
auto pi_below_li=[&](size_t li, std::pair<size_t, size_t> 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=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];
|
||||
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
|
||||
|
||||
// Input polylines are supposed going from left to right
|
||||
std::vector<Event> event_queue;
|
||||
event_queue.reserve(2*polylines.size());
|
||||
for(size_t li=0; li<polylines.size(); ++li){
|
||||
Polyline &pl=polylines[li];
|
||||
event_queue.emplace_back(pl.front(), li, EVENT_TYPE::INSERT);
|
||||
for(size_t pi=1; pi<pl.size()-1; ++pi){
|
||||
event_queue.emplace_back(pl[pi], li, EVENT_TYPE::MIDDLE);
|
||||
assert(less_xy_2(pts[pl[pi-1]], pts[pl[pi]]));
|
||||
}
|
||||
event_queue.emplace_back(pl.back(), li, EVENT_TYPE::REMOVE);
|
||||
assert(less_xy_2(pts[pl.front()], pts[pl.back()]));
|
||||
}
|
||||
#ifdef DOUBLE_2D_SNAP_VERBOSE
|
||||
std::cout << "Sort the event queue along x direction" << std::endl;
|
||||
#endif
|
||||
std::sort(event_queue.begin(), event_queue.end(), event_comparator);
|
||||
#ifdef DOUBLE_2D_SNAP_VERBOSE
|
||||
std::cout << "Goes through Event" << std::endl;
|
||||
#endif
|
||||
// TODO Vector is suboptimal for arbitrary insertion (Actually negligeable in the running time)
|
||||
// Track order of the lines along y along the events
|
||||
std::vector< size_t > y_order;
|
||||
y_order.push_back(event_queue[0].li);
|
||||
for(size_t i=1; i<event_queue.size();++i)
|
||||
{
|
||||
Event event = event_queue[i];
|
||||
#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
|
||||
// Find the position of the event along the y direction
|
||||
std::vector<size_t>::iterator pos_it = y_order.begin();
|
||||
if(!y_order.empty())
|
||||
pos_it = std::lower_bound(y_order.begin(), y_order.end(), std::pair<size_t, size_t>(event.pi, event.li), pi_below_li);
|
||||
if(event.type==EVENT_TYPE::REMOVE){
|
||||
assert(*pos_it==event.li);
|
||||
pos_it = y_order.erase(pos_it);
|
||||
}
|
||||
|
||||
if(!y_order.empty()){
|
||||
auto close_event=[&](std::size_t pi, std::size_t li){
|
||||
if((pi==polylines[li].front()) || (pi==polylines[li].back()))
|
||||
return true;
|
||||
|
||||
const Point_2 &p = pts[pi];
|
||||
Polyline &pl=polylines[li];
|
||||
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];
|
||||
bound = (std::max)(bound, round_bound_pts[pl.front()]);
|
||||
bound = (std::max)(bound, round_bound_pts[pl.back()]);
|
||||
bound*=4;
|
||||
|
||||
if(possibly(csq_dist_2(p, seg, bound)!=CGAL::LARGER))
|
||||
{
|
||||
if constexpr (std::is_same_v<Exact_predicates_exact_constructions_kernel, typename Traits::Exact_type>)
|
||||
{
|
||||
internal::Evaluate<typename Traits::FT> evaluate;
|
||||
// We refine the pts to reduce the rounding shift and check again
|
||||
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
|
||||
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]]);
|
||||
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)
|
||||
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;
|
||||
};
|
||||
|
||||
if(event.pi != event_queue[i-1].pi)
|
||||
{
|
||||
// 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)){
|
||||
if((pi!=polylines[*above].front()) && (pi!=polylines[*above].back()))
|
||||
pi=pts.size()-1;
|
||||
++above;
|
||||
}
|
||||
|
||||
// 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;
|
||||
}
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
if(event.type==EVENT_TYPE::INSERT){
|
||||
y_order.insert(pos_it, event.li);
|
||||
}
|
||||
}
|
||||
|
||||
// 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 less_y_2(pts[i],pts[j]);
|
||||
};
|
||||
std::vector< std::size_t > 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 <class Concurrency_tag=Sequential_tag, class Traits, class PointsRange , class PolylineRange>
|
||||
void merge_duplicate_points_in_polylines(PointsRange &pts, PolylineRange &polylines, const Traits &traits)
|
||||
{
|
||||
using Point_2 = typename Traits::Point_2;
|
||||
using Polyline = std::remove_cv_t<typename std::iterator_traits<typename PolylineRange::iterator>::value_type>;
|
||||
|
||||
using Less_xy_2 = typename Traits::Less_xy_2;
|
||||
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]);
|
||||
};
|
||||
|
||||
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<Point_2> new_pts;
|
||||
std::vector<std::size_t> 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 (Polyline& polyline : polylines) {
|
||||
std::vector<std::size_t> updated_polyline;
|
||||
for (std::size_t i=0; i<polyline.size(); ++i) {
|
||||
std::size_t new_pi=old_to_new_index[polyline[i]];
|
||||
if(i==0 || (new_pi!=updated_polyline[updated_polyline.size()-1]))
|
||||
updated_polyline.push_back(new_pi);
|
||||
assert(new_pi<pts.size());
|
||||
assert(pts[new_pi]==new_pts[polyline[i]]);
|
||||
}
|
||||
std::swap(polyline, updated_polyline);
|
||||
}
|
||||
}
|
||||
|
||||
// Some points may have collapsed on a vertical segment, we subdivide these vertical segments accordingly
|
||||
template <class Concurrency_tag=Sequential_tag, class Traits, class PointsRange , class PolylineRange>
|
||||
void snap_post_process(PointsRange &pts, PolylineRange &polylines, const Traits &traits)
|
||||
{
|
||||
using Polyline = std::remove_cv_t<typename std::iterator_traits<typename PolylineRange::iterator>::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<std::size_t, decltype(Less_indexes_xy_2)> 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<std::size_t, decltype(Less_indexes_xy_2)>::iterator;
|
||||
|
||||
for(Polyline &poly: polylines){
|
||||
std::vector<std::size_t> 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()){
|
||||
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);
|
||||
}
|
||||
}
|
||||
|
||||
template <class Concurrency_tag=Sequential_tag, class Traits, class PointsRange , class PolylineRange>
|
||||
void double_snap_rounding_2_disjoint(PointsRange &pts, PolylineRange &polylines, const Traits &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 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.
|
||||
*
|
||||
* 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 <class InputIterator , class OutputContainer, class NamedParameters = parameters::Default_named_parameters>
|
||||
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<internal_np::concurrency_tag_t,
|
||||
NamedParameters,
|
||||
Sequential_tag>::type;
|
||||
|
||||
using InputKernel = typename Kernel_traits<std::remove_cv_t<typename std::iterator_traits<InputIterator>::value_type>>::Kernel;
|
||||
using DefaultTraits = Float_snap_rounding_traits_2<InputKernel>;
|
||||
using Traits = typename internal_np::Lookup_named_param_def<internal_np::geom_traits_t,
|
||||
NamedParameters,
|
||||
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<Segment_2>::iterator;
|
||||
|
||||
using Polyline = std::remove_cv_t<typename std::iterator_traits<typename OutputContainer::iterator>::value_type>;
|
||||
|
||||
using Less_xy_2 = typename Traits::Less_xy_2;
|
||||
|
||||
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<Segment_2> convert_input;
|
||||
for(InputIterator it=begin; it!=end; ++it)
|
||||
convert_input.push_back(to_exact(*it));
|
||||
std::vector<Segment_2> 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));
|
||||
|
||||
#ifdef DOUBLE_2D_SNAP_VERBOSE
|
||||
std::cout << "Change format to range of points and indexes" << std::endl;
|
||||
#endif
|
||||
std::set<Point_2> unique_point_set;
|
||||
std::map<Point_2, int> point_to_index;
|
||||
std::vector<Point_2> 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()];
|
||||
if(Less_xy_2()(it->source(), it->target()))
|
||||
polylines.push_back({index1, index2});
|
||||
else
|
||||
polylines.push_back({index2, index1});
|
||||
}
|
||||
|
||||
// Main algorithm
|
||||
internal::double_snap_rounding_2_disjoint<Concurrency_tag, Traits>(pts, polylines, traits);
|
||||
|
||||
#ifdef DOUBLE_2D_SNAP_VERBOSE
|
||||
std::cout << "Build output" << std::endl;
|
||||
#endif
|
||||
|
||||
// Output polylines
|
||||
for(auto &poly: polylines){
|
||||
Polyline new_line;
|
||||
for(std::size_t pi: poly)
|
||||
new_line.push_back(from_exact(pts[pi]));
|
||||
out.push_back(new_line);
|
||||
}
|
||||
|
||||
return out.end();
|
||||
}
|
||||
|
||||
/**
|
||||
* ingroup PkgSnapRounding2Ref
|
||||
*
|
||||
* 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
|
||||
* @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 <class InputIterator , class OutputIterator, class NamedParameters = parameters::Default_named_parameters>
|
||||
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<internal_np::concurrency_tag_t,
|
||||
NamedParameters,
|
||||
Sequential_tag>::type;
|
||||
|
||||
using InputKernel = typename Kernel_traits<std::remove_cv_t<typename std::iterator_traits<InputIterator>::value_type>>::Kernel;
|
||||
using DefaultTraits = Float_snap_rounding_traits_2<InputKernel>;
|
||||
using Traits = typename internal_np::Lookup_named_param_def<internal_np::geom_traits_t,
|
||||
NamedParameters,
|
||||
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 SegmentRange = std::vector<Segment_2>;
|
||||
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;
|
||||
|
||||
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();
|
||||
|
||||
Less_xy_2 less_xy_2 = traits.less_xy_2_object();
|
||||
Construct_segment_2 segment_2 = traits.construct_segment_2_object();
|
||||
|
||||
std::vector<Segment_2> convert_input;
|
||||
for(InputIterator it=begin; it!=end; ++it)
|
||||
if(it->source()!=it->target())
|
||||
convert_input.push_back(to_exact(*it));
|
||||
std::vector<Segment_2> 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));
|
||||
|
||||
#ifdef DOUBLE_2D_SNAP_VERBOSE
|
||||
std::cout << "Change format to range of points and indexes" << std::endl;
|
||||
#endif
|
||||
std::set<Point_2, Less_xy_2> unique_point_set(less_xy_2);
|
||||
std::map<Point_2, std::size_t> point_to_index;
|
||||
std::vector<Point_2> pts;
|
||||
std::vector< std::vector< std::size_t> > polylines;
|
||||
|
||||
// Transform range of the segments in the range of points and polyline of indexes
|
||||
for(SegmentRangeIterator 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(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()];
|
||||
if(less_xy_2(it->source(), it->target()))
|
||||
polylines.push_back({index1, index2});
|
||||
else
|
||||
polylines.push_back({index2, index1});
|
||||
}
|
||||
|
||||
|
||||
// Main algorithm
|
||||
internal::double_snap_rounding_2_disjoint<Concurrency_tag>(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<std::size_t,std::size_t> > set_out_segs;
|
||||
for(auto &poly: polylines)
|
||||
for(std::size_t i=1; i<poly.size(); ++i)
|
||||
set_out_segs.emplace((std::min)(poly[i-1],poly[i]),(std::max)(poly[i-1],poly[i]));
|
||||
for(auto &pair: set_out_segs){
|
||||
*out++=from_exact(segment_2(pts[pair.first], pts[pair.second]));
|
||||
assert(pts[pair.first]!=pts[pair.second]);
|
||||
}
|
||||
|
||||
return out;
|
||||
}
|
||||
|
||||
/**
|
||||
* 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`.
|
||||
* Any polygon is guarantee to remain a Polygon in the output but may present pinched section or/and common vertices or segments with
|
||||
* other polygons.
|
||||
*
|
||||
* @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
|
||||
*
|
||||
* \cgalParamNBegin{compute_intersection}
|
||||
* \cgalParamDescription{Enable intersection computation between the polygons before performing snapping.}
|
||||
* \cgalParamType{boolean}
|
||||
* \cgalParamDefault{false}
|
||||
* \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
|
||||
* @warning The convex property of the polygons is not necessarly preserved
|
||||
*/
|
||||
template <class InputIterator, class OutputIterator, class NamedParameters = parameters::Default_named_parameters>
|
||||
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<internal_np::concurrency_tag_t,
|
||||
NamedParameters,
|
||||
Sequential_tag>::type;
|
||||
|
||||
using Polygon_2 = typename std::iterator_traits<InputIterator>::value_type;
|
||||
using InputKernel = typename Kernel_traits<typename Polygon_2::Point_2>::Kernel;
|
||||
using DefaultTraits = Float_snap_rounding_traits_2<InputKernel>;
|
||||
using Traits = typename internal_np::Lookup_named_param_def<internal_np::geom_traits_t,
|
||||
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;
|
||||
|
||||
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::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();
|
||||
|
||||
#ifdef DOUBLE_2D_SNAP_VERBOSE
|
||||
std::cout << "Change format to range of points and indexes" << std::endl;
|
||||
#endif
|
||||
std::vector<Point_2> pts;
|
||||
std::vector< std::vector< std::size_t> > polylines;
|
||||
std::vector<size_t> polygon_index;
|
||||
|
||||
if(compute_intersections){
|
||||
// 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));
|
||||
for(InputIterator it=begin; it!=end; ++it){
|
||||
size_t index_start=polylines.size();
|
||||
polygon_index.push_back(polylines.size());
|
||||
for(const typename Polygon_2::Point_2 &p: it->vertices())
|
||||
pts.push_back(to_exact(p));
|
||||
for(size_t i=0; i<P.size()-1; ++i)
|
||||
if(less_xy_2(pts[i], pts[i+1]))
|
||||
polylines.push_back({i, i+1});
|
||||
else
|
||||
polylines.push_back({i+1, i});
|
||||
polylines.push_back({pts.size()-1,index_start});
|
||||
assert(pts.size()==polylines.size());
|
||||
}
|
||||
polygon_index.push_back(polylines.size());
|
||||
}
|
||||
|
||||
// Main algorithm
|
||||
internal::double_snap_rounding_2_disjoint<Concurrency_tag>(pts, polylines, traits);
|
||||
|
||||
#ifdef DOUBLE_2D_SNAP_VERBOSE
|
||||
std::cout << "Build output" << std::endl;
|
||||
#endif
|
||||
|
||||
// Output a range of segments while removing duplicate ones
|
||||
for(std::size_t input_ind=0; input_ind<std::size_t(std::distance(begin,end)); ++input_ind){
|
||||
for(std::size_t pl_ind=polygon_index[input_ind]; pl_ind<polygon_index[input_ind+1]; ++pl_ind){
|
||||
std::vector<std::size_t> &poly = polylines[pl_ind];
|
||||
Polygon_2 P;
|
||||
for(std::size_t i=1; i<poly.size(); ++i)
|
||||
P.push_back(from_exact(pts[poly[i]]));
|
||||
}
|
||||
*out++=P;
|
||||
}
|
||||
}
|
||||
|
||||
|
||||
/**
|
||||
* 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.
|
||||
*
|
||||
* \tparam Polygon_2 model of CGAL::Polygon_2
|
||||
* \tparam NamedParameters a sequence of \ref bgl_namedparameters "Named Parameters"
|
||||
* \param P the input polygon
|
||||
* \param out the output polygon
|
||||
* \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{a multiplier of the error value for boundary edges to preserve the boundaries}
|
||||
* \cgalParamType{double}
|
||||
* \cgalParamDefault{100}
|
||||
* \cgalParamNEnd
|
||||
* @warning The convex property is not necessarly preserved
|
||||
*/
|
||||
template <class Polygon_2, class NamedParameters = parameters::Default_named_parameters>
|
||||
void compute_snapped_polygon_2(const Polygon_2 &P, Polygon_2 &out, const NamedParameters &np = parameters::default_values())
|
||||
{
|
||||
std::array<Polygon_2, 1> vec({P});
|
||||
std::vector<Polygon_2> out_vec;
|
||||
compute_snapped_polygons_2(vec.begin(), vec.end(), std::back_inserter(out_vec), np);
|
||||
out = out_vec[0];
|
||||
}
|
||||
|
||||
|
||||
|
||||
} //namespace CGAL
|
||||
|
||||
#endif
|
||||
|
|
@ -0,0 +1,116 @@
|
|||
// 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 <CGAL/Arr_segment_traits_2.h>
|
||||
|
||||
#include <CGAL/Exact_predicates_exact_constructions_kernel.h>
|
||||
#include <CGAL/Cartesian_converter.h>
|
||||
|
||||
#include <CGAL/Named_function_parameters.h>
|
||||
|
||||
#include <CGAL/mutex.h>
|
||||
|
||||
namespace CGAL {
|
||||
|
||||
/*!
|
||||
\ingroup PkgSnapRounding2Ref
|
||||
|
||||
The class `Float_snap_rounding_traits_2<Kernel>` 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<Gmpq>`.
|
||||
|
||||
\cgalModels{SnapRoundingTraits_2}
|
||||
|
||||
*/
|
||||
template<typename Input_Kernel, typename Exact_Kernel = Exact_predicates_exact_constructions_kernel>
|
||||
struct Float_snap_rounding_traits_2: Arr_segment_traits_2<Exact_Kernel>{
|
||||
typedef Arr_segment_traits_2<Exact_Kernel> 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;
|
||||
|
||||
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<Input_Kernel, Exact_Kernel> Converter_to_exact;
|
||||
typedef Cartesian_converter<Exact_Kernel, Input_Kernel> 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<double>::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
|
||||
|
|
@ -59,6 +59,10 @@ public:
|
|||
public:
|
||||
void operator()(const Point_2& p, NT pixel_size, NT &x, NT &y)
|
||||
{
|
||||
// TODO a proper and efficient double_ceil was done for snap_rounding_polygon_soup
|
||||
#if 0
|
||||
x = double_ceil((x / pixel_size) - 0.5) * pixel_size;
|
||||
#else
|
||||
NT x_tmp = p.x() / pixel_size;
|
||||
NT y_tmp = p.y() / pixel_size;
|
||||
|
||||
|
|
@ -81,6 +85,7 @@ public:
|
|||
|
||||
x = x_floor * pixel_size + pixel_size / NT(2.0);
|
||||
y = y_floor * pixel_size + pixel_size / NT(2.0);
|
||||
#endif
|
||||
}
|
||||
};
|
||||
|
||||
|
|
|
|||
|
|
@ -7,6 +7,8 @@ project(Snap_rounding_2_Tests)
|
|||
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")
|
||||
|
|
|
|||
|
|
@ -0,0 +1,331 @@
|
|||
#include <CGAL/Exact_predicates_exact_constructions_kernel.h>
|
||||
#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
|
||||
|
||||
#include <CGAL/Float_snap_rounding_traits_2.h>
|
||||
#include <CGAL/Float_snap_rounding_2.h>
|
||||
#include <CGAL/Arr_segment_traits_2.h>
|
||||
#include <CGAL/Surface_sweep_2_algorithms.h>
|
||||
|
||||
#include <CGAL/Boolean_set_operations_2.h>
|
||||
|
||||
#include <CGAL/Random.h>
|
||||
|
||||
typedef CGAL::Exact_predicates_exact_constructions_kernel Kernel;
|
||||
typedef CGAL::Float_snap_rounding_traits_2<Kernel> Traits_2;
|
||||
typedef Kernel::Segment_2 Segment_2;
|
||||
typedef Kernel::Point_2 Point_2;
|
||||
typedef Kernel::Vector_2 Vector_2;
|
||||
typedef std::vector<Point_2 > Polyline_2;
|
||||
typedef std::vector<Polyline_2> Polyline_range_2;
|
||||
typedef Kernel::FT FT;
|
||||
|
||||
typedef CGAL::Arr_segment_traits_2<Kernel> Arr_segment_traits_2;
|
||||
typedef Arr_segment_traits_2::Curve_2 Curve_2;
|
||||
|
||||
typedef CGAL::Polygon_2<Kernel> Polygon_2;
|
||||
typedef CGAL::Polygon_with_holes_2<Kernel> Polygon_with_holes_2;
|
||||
typedef std::vector<Polygon_with_holes_2> Pwh_vec_2;
|
||||
|
||||
#ifdef BENCH_AND_VERBOSE_FLOAT_SNAP_ROUNDING_2
|
||||
#include <CGAL/Real_timer.h>
|
||||
#include <CGAL/Cartesian_converter.h>
|
||||
#include <CGAL/Simple_cartesian.h>
|
||||
typedef CGAL::Simple_cartesian<double> Naive;
|
||||
typedef CGAL::Cartesian_converter<Kernel,Naive> EK_to_IK;
|
||||
typedef CGAL::Cartesian_converter<Naive, Kernel> IK_to_EK;
|
||||
#endif
|
||||
|
||||
#ifdef COMPARE_WITH_INTEGER_SNAP_ROUNDING_2
|
||||
#include <CGAL/Snap_rounding_traits_2.h>
|
||||
#include <CGAL/Snap_rounding_2.h>
|
||||
typedef CGAL::Snap_rounding_traits_2<Kernel> SnapTraits;
|
||||
#endif
|
||||
|
||||
//Biggest double with ulp smaller than an integer
|
||||
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){
|
||||
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)));
|
||||
}
|
||||
|
||||
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<class Iterator, class OutRange>
|
||||
void compute_subcurves_and_naive_snap(Iterator begin, Iterator end, OutRange &out){
|
||||
EK_to_IK to_inexact;
|
||||
IK_to_EK to_exact;
|
||||
std::vector<Curve_2> segs;
|
||||
for(Iterator it=begin; it!=end; ++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;
|
||||
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);
|
||||
}
|
||||
#endif
|
||||
|
||||
void test(const std::vector<Segment_2> &segs){
|
||||
std::vector<Segment_2> 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()) <<std::endl;
|
||||
out.clear();
|
||||
t.reset();
|
||||
t.start();
|
||||
#endif
|
||||
CGAL::compute_snapped_subcurves_2(segs.begin(), segs.end(), std::back_inserter(out));
|
||||
std::vector<Segment_2> out2(out.begin(), out.end());
|
||||
#ifdef BENCH_AND_VERBOSE_FLOAT_SNAP_ROUNDING_2
|
||||
t.stop();
|
||||
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<SnapTraits>(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
|
||||
#ifdef BENCH_AND_VERBOSE_FLOAT_SNAP_ROUNDING_2
|
||||
std::cout << "\n";
|
||||
#endif
|
||||
}
|
||||
|
||||
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<Segment_2> segs;
|
||||
for(size_t i=0; i!=nb_segments; ++i)
|
||||
segs.emplace_back(random_point(r), random_point(r));
|
||||
|
||||
test(segs);
|
||||
}
|
||||
|
||||
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<Polygon_2> polygons;
|
||||
for(size_t i=0; i!=nb_polygons; ++i){
|
||||
Polygon_2 poly;
|
||||
for(size_t j=0; j!=nb_pts; ++j)
|
||||
poly.push_back(Point_2(random_point(r), random_point(r)));
|
||||
polygons.emplace_back(poly);
|
||||
}
|
||||
|
||||
std::vector<Polygon_2> out;
|
||||
CGAL::compute_snapped_polygons_2(polygons.begin(), polygons.end(), std::back_inserter(out));
|
||||
|
||||
std::vector<Segment_2> segs;
|
||||
for(const Polygon_2 &poly: out)
|
||||
for(size_t i=1; i<poly.size(); ++i)
|
||||
segs.emplace_back(poly[i-1],poly[i]);
|
||||
assert(!CGAL::do_curves_intersect(segs.begin(), segs.end()));
|
||||
}
|
||||
|
||||
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<Segment_2> segs;
|
||||
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));
|
||||
}
|
||||
|
||||
test(segs);
|
||||
}
|
||||
|
||||
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);
|
||||
#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;
|
||||
};
|
||||
|
||||
Polygon_2 scene=random_rotated_square();
|
||||
Polygon_2 snap_scene;
|
||||
Pwh_vec_2 out_intersection;
|
||||
|
||||
for(size_t i=0; i<nb_iterations; ++i){
|
||||
out_intersection.clear();
|
||||
CGAL::intersection(random_rotated_square(), scene, std::back_inserter(out_intersection));
|
||||
assert(out_intersection.size()==1 && out_intersection[0].number_of_holes()==0);
|
||||
#ifdef BENCH_AND_VERBOSE_FLOAT_SNAP_ROUNDING_2
|
||||
CGAL::Real_timer t;
|
||||
t.start();
|
||||
#endif
|
||||
compute_snapped_polygon_2(out_intersection[0].outer_boundary(), snap_scene);
|
||||
#ifdef BENCH_AND_VERBOSE_FLOAT_SNAP_ROUNDING_2
|
||||
t.stop();
|
||||
std::cout << "Iteration " << i << std::endl;
|
||||
std::cout << "Polygon size: " << out_intersection[0].outer_boundary().size()
|
||||
<< " , Snapped polygon size: " << snap_scene.size() << std::endl;
|
||||
std::cout << "is convex: " << snap_scene.is_convex() << std::endl;
|
||||
std::cout << "Running time: " << t.time() << std::endl;
|
||||
#endif
|
||||
scene=snap_scene;
|
||||
}
|
||||
}
|
||||
|
||||
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 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();
|
||||
|
||||
// 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));
|
||||
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);
|
||||
|
||||
// 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));
|
||||
test(segs);
|
||||
segs.clear();
|
||||
|
||||
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);
|
||||
segs.clear();
|
||||
}
|
||||
|
||||
void test_float_snap_rounding(){
|
||||
#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;
|
||||
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);
|
||||
}
|
||||
|
||||
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);
|
||||
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);
|
||||
}
|
||||
|
|
@ -0,0 +1,42 @@
|
|||
// #define DOUBLE_2D_SNAP_VERBOSE
|
||||
// #define DOUBLE_2D_SNAP_FULL_VERBOSE
|
||||
#define BENCH_AND_VERBOSE_FLOAT_SNAP_ROUNDING_2
|
||||
|
||||
#include <CGAL/Exact_rational.h>
|
||||
#include <CGAL/Cartesian.h>
|
||||
#include <CGAL/Cartesian_converter.h>
|
||||
|
||||
#include <CGAL/Float_snap_rounding_traits_2.h>
|
||||
#include <CGAL/Float_snap_rounding_2.h>
|
||||
|
||||
typedef CGAL::Exact_predicates_exact_constructions_kernel Epeck;
|
||||
typedef CGAL::Cartesian<CGAL::Exact_rational> Rational_Kernel;
|
||||
typedef CGAL::Float_snap_rounding_traits_2<Rational_Kernel, Rational_Kernel> Rational_Traits;
|
||||
|
||||
template<class Traits>
|
||||
struct Test{
|
||||
typedef typename Traits::Segment_2 Segment_2;
|
||||
typedef typename Traits::Point_2 Point_2;
|
||||
|
||||
void fix_test(){
|
||||
std::vector<Segment_2> 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<Segment_2> out;
|
||||
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()));
|
||||
}
|
||||
};
|
||||
|
||||
int main(/*int argc,char *argv[]*/)
|
||||
{
|
||||
Test<Rational_Traits>().fix_test();
|
||||
return(0);
|
||||
}
|
||||
|
|
@ -6,25 +6,56 @@
|
|||
|
||||
#include <CGAL/Exact_predicates_exact_constructions_kernel.h>
|
||||
#include <CGAL/Arr_segment_traits_2.h>
|
||||
#include <CGAL/Arr_polyline_traits_2.h>
|
||||
#include <CGAL/Surface_sweep_2_algorithms.h>
|
||||
|
||||
typedef CGAL::Exact_predicates_exact_constructions_kernel Kernel;
|
||||
typedef Kernel::Point_2 Point_2;
|
||||
typedef CGAL::Arr_segment_traits_2<Kernel> 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<Kernel>;
|
||||
using Traits = CGAL::Arr_polyline_traits_2<Segment_traits>;
|
||||
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<Point> 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<Segment> 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<Point_2> 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<Point_2>(std::cout, "\n"));
|
||||
|
||||
// Compute the non-intersecting sub-segments induced by the input segments.
|
||||
std::list<Segment_2> sub_segs;
|
||||
std::list<My_polyline> 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;
|
||||
}
|
||||
|
|
|
|||
|
|
@ -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.
|
||||
|
|
|
|||
Loading…
Reference in New Issue