diff --git a/Arrangement_on_surface_2/include/CGAL/Arr_geodesic_arc_on_sphere_traits_2.h b/Arrangement_on_surface_2/include/CGAL/Arr_geodesic_arc_on_sphere_traits_2.h index c4903e1c8164..411dda8c5e83 100644 --- a/Arrangement_on_surface_2/include/CGAL/Arr_geodesic_arc_on_sphere_traits_2.h +++ b/Arrangement_on_surface_2/include/CGAL/Arr_geodesic_arc_on_sphere_traits_2.h @@ -39,9 +39,86 @@ namespace CGAL { +/*! Represent an extended 3D direction that is used in turn to represent a + * spherical-arc endpoint. The extended data consists of two flags that + * indicate whether the point is on the x and on a y boundaries, + * respectively. + */ +template +class Arr_extended_direction_3 : public Kernel::Direction_3 { +public: + using FT = typename Kernel::FT; + using Direction_3 = typename Kernel::Direction_3; + + /*! Enumeration of discontinuity type */ + enum Location_type { + NO_BOUNDARY_LOC = 0, + MIN_BOUNDARY_LOC, + MID_BOUNDARY_LOC, + MAX_BOUNDARY_LOC + }; + +private: + using Direction_2 = typename Kernel::Direction_2; + + //! The point discontinuity type + Location_type m_location; + + inline Sign x_sign(Direction_3 d) const { return CGAL::sign(d.dx()); } + + inline Sign y_sign(Direction_3 d) const { return CGAL::sign(d.dy()); } + + inline Sign z_sign(Direction_3 d) const { return CGAL::sign(d.dz()); } + +public: + /*! Default constructor */ + Arr_extended_direction_3() : + Direction_3(0, 0, 1), + m_location(MAX_BOUNDARY_LOC) + {} + + /*! Constructor */ + Arr_extended_direction_3(const Direction_3& dir, Location_type location) : + Direction_3(dir), + m_location(location) + {} + + /*! Copy constructor */ + Arr_extended_direction_3(const Arr_extended_direction_3& other) : + Direction_3(static_cast(other)) + { m_location = other.discontinuity_type(); } + + /*! Assignment operator */ + Arr_extended_direction_3& operator=(const Arr_extended_direction_3& other) { + *(static_cast(this)) = static_cast(other); + m_location = other.discontinuity_type(); + return (*this); + } + + /*! Set the location type of the point. + */ + void set_location(Location_type location) { m_location = location; } + + /*! Obtain the location type of the point. + */ + Location_type location() const { return m_location; } + + /*! Obtain the discontinuity type of the point. + * \todo deprecate this one; use the above instead. + */ + Location_type discontinuity_type() const { return m_location; } + + bool is_no_boundary() const { return (m_location == NO_BOUNDARY_LOC); } + + bool is_min_boundary() const { return (m_location == MIN_BOUNDARY_LOC); } + + bool is_mid_boundary() const { return (m_location == MID_BOUNDARY_LOC); } + + bool is_max_boundary() const { return (m_location == MAX_BOUNDARY_LOC); } +}; + template class Arr_x_monotone_geodesic_arc_on_sphere_3; template class Arr_geodesic_arc_on_sphere_3; -template class Arr_extended_direction_3; /*! A traits class-template for constructing and maintaining arcs of great * circles embedded on spheres. It is parameterized from a (linear) geometry @@ -54,44 +131,43 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { friend class Arr_extended_direction_3; public: - typedef Kernel_ Kernel; + using Kernel = Kernel_; // Category tags: - typedef Tag_true Has_left_category; - typedef Tag_true Has_merge_category; - typedef Tag_false Has_do_intersect_category; + using Has_left_category = Tag_true; + using Has_merge_category = Tag_true; + using Has_do_intersect_category = Tag_false; - typedef Arr_identified_side_tag Left_side_category; - typedef Arr_contracted_side_tag Bottom_side_category; - typedef Arr_contracted_side_tag Top_side_category; - typedef Arr_identified_side_tag Right_side_category; + using Left_side_category = Arr_identified_side_tag; + using Bottom_side_category = Arr_contracted_side_tag; + using Top_side_category = Arr_contracted_side_tag; + using Right_side_category = Arr_identified_side_tag; - typedef std::integral_constant Zero_atan_y; + using Zero_atan_y = std::integral_constant; // Traits objects - typedef Arr_extended_direction_3 Point_2; - typedef Arr_x_monotone_geodesic_arc_on_sphere_3 X_monotone_curve_2; - typedef Arr_geodesic_arc_on_sphere_3 Curve_2; - typedef unsigned int Multiplicity; + using Point_2 = Arr_extended_direction_3; + using X_monotone_curve_2 = Arr_x_monotone_geodesic_arc_on_sphere_3; + using Curve_2 = Arr_geodesic_arc_on_sphere_3; + using Multiplicity = std::size_t; /*! Default constructor */ Arr_geodesic_arc_on_sphere_traits_2() {} protected: - typedef typename Kernel::FT FT; + using FT = typename Kernel::FT; - typedef typename Kernel::Direction_3 Direction_3; - typedef typename Kernel::Vector_3 Vector_3; - typedef typename Kernel::Direction_2 Direction_2; - typedef typename Kernel::Vector_2 Vector_2; + using Direction_3 = typename Kernel::Direction_3; + using Vector_3 = typename Kernel::Vector_3; + using Direction_2 = typename Kernel::Direction_2; + using Vector_2 = typename Kernel::Vector_2; /*! Obtain the intersection of the identification arc and the xy plane. * By default, it is the vector directed along the negative x axis * (x = -infinity). * \return the intersection of the identification arc and the xy plane. */ - inline static const Direction_2& identification_xy() - { + inline static const Direction_2& identification_xy() { static const Direction_2 d(atan_x, atan_y); return d; } @@ -101,8 +177,7 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { * (y = infinity). * \return the normal of the plane that contains the identification arc. */ - inline static const Direction_3& identification_normal() - { + inline static const Direction_3& identification_normal() { static const Direction_3 d(atan_y, -atan_x, 0); return d; } @@ -110,8 +185,7 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { /*! Obtain the 2D direction directed along the negative x axis * \return the direction directed at x = -infinity */ - inline static const Direction_2& neg_x_2() - { + inline static const Direction_2& neg_x_2() { CGAL_STATIC_THREAD_LOCAL_VARIABLE_2(Direction_2, d, -1, 0); return d; } @@ -119,8 +193,7 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { /*! Obtain the 2D direction directed along the negative y axis * \return the direction directed at y = -infinity */ - inline static const Direction_2& neg_y_2() - { + inline static const Direction_2& neg_y_2() { CGAL_STATIC_THREAD_LOCAL_VARIABLE_2(Direction_2, d, 0, -1); return d; } @@ -186,8 +259,7 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { * \param dir the direction. */ inline Oriented_side oriented_side(const Direction_3& normal, - const Direction_3 dir) const - { + const Direction_3 dir) const { FT dot = normal.vector() * dir.vector(); return CGAL::sign(dot); } @@ -197,9 +269,8 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { * \param d2 the second direction. * \return the relative orientation of d1 and d2. */ - inline Orientation orientation(const Direction_2& d1, const Direction_2& d2) - const - { + inline Orientation orientation(const Direction_2& d1, + const Direction_2& d2) const { const Kernel& kernel(*this); return kernel.orientation_2_object()(d1.vector(), d2.vector()); } @@ -209,8 +280,7 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { * \param d2 the second direction. */ inline Direction_3 construct_normal_3(const Direction_3& d1, - const Direction_3& d2) const - { + const Direction_3& d2) const { const Kernel& kernel(*this); Vector_3 v = kernel.construct_cross_product_vector_3_object()(d1.vector(), d2.vector()); @@ -224,8 +294,7 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { * \return true if dir is contained in plane; false otherwise. * \pre the plane contains the origin. */ - inline bool has_on(const Direction_3& normal, const Direction_3& dir) const - { + inline bool has_on(const Direction_3& normal, const Direction_3& dir) const { FT dot = normal.vector() * dir.vector(); return CGAL::sign(dot) == ZERO; } @@ -239,8 +308,7 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { * LARGER - v(d1) > v(d2). */ inline Comparison_result compare_y(const Direction_3& d1, - const Direction_3& d2) const - { + const Direction_3& d2) const { Vector_3 v1 = d1.vector(); Vector_3 v2 = d2.vector(); @@ -274,8 +342,7 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { * LARGER - u(d1) > u(d2). */ inline Comparison_result compare_x(const Direction_2& d1, - const Direction_2& d2) const - { + const Direction_2& d2) const { const Kernel& kernel = *this; if (kernel.equal_2_object()(d1, d2)) return EQUAL; const Direction_2& d = identification_xy(); @@ -293,8 +360,7 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { * \pre d2 does not coincide with any pole. */ inline Comparison_result compare_x(const Direction_3& d1, - const Direction_3& d2) const - { + const Direction_3& d2) const { // Compare the projections onto the xy plane: Direction_2 d1_2 = project_xy(d1); Direction_2 d2_2 = project_xy(d2); @@ -313,8 +379,7 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { * \pre d2 does not lie on the discontinuity arc. */ inline Comparison_result compare_xy(const Direction_3& d1, - const Direction_3& d2) const - { + const Direction_3& d2) const { Comparison_result res = compare_x(d1, d2); if (res == EQUAL) return compare_y(d1, d2); return res; @@ -327,9 +392,9 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { * false otherwise. * \pre point does not coincide with one of the poles */ - bool is_in_x_range(const X_monotone_curve_2& xcv, const Point_2& point) const - { - typedef Arr_geodesic_arc_on_sphere_traits_2 Traits; + bool is_in_x_range(const X_monotone_curve_2& xcv, + const Point_2& point) const { + using Traits = Arr_geodesic_arc_on_sphere_traits_2; CGAL_precondition(!point.is_min_boundary()); CGAL_precondition(!point.is_max_boundary()); @@ -358,8 +423,7 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { */ void intersection_with_identification(const X_monotone_curve_2& xcv, Direction_3& dp, - std::true_type) const - { + std::true_type) const { const Direction_3& normal = xcv.normal(); dp = (CGAL::sign(normal.dz()) == POSITIVE) ? Direction_3(-(normal.dz()), 0, normal.dx()) : @@ -371,8 +435,7 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { */ void intersection_with_identification(const X_monotone_curve_2& xcv, Direction_3& dp, - std::false_type) const - { + std::false_type) const { const Direction_3& normal = xcv.normal(); FT z((atan_x * normal.dx() + atan_y * normal.dy()) / -(normal.dz())); @@ -383,8 +446,7 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { * \param[in] cv the curve */ bool overlap_with_identification(const X_monotone_curve_2& xcv, - std::true_type) const - { + std::true_type) const { const Direction_3& normal = xcv.normal(); return ((x_sign(normal) == ZERO) && (((y_sign(normal) == NEGATIVE) && !xcv.is_directed_right()) || @@ -395,8 +457,7 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { * \param[in] cv the curve */ bool overlap_with_identification(const X_monotone_curve_2& xcv, - std::false_type) const - { + std::false_type) const { const Direction_3& normal = xcv.normal(); const Direction_3& iden_normal = identification_normal(); const Direction_2 iden_normal_xy = project_xy(iden_normal); @@ -417,9 +478,9 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { /*! A functor that constructs a point on the sphere. */ class Construct_point_2 { protected: - typedef Arr_geodesic_arc_on_sphere_traits_2 Traits; + using Traits = Arr_geodesic_arc_on_sphere_traits_2; - /*! The traits (in case it has state) */ + //! The traits (in case it has state) const Traits& m_traits; /*! Constructor @@ -436,8 +497,7 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { * \param[in] y the y coordinate * \param[in] z the z coordinate */ - Point_2 operator()(const FT& x, const FT& y, const FT& z) - { + Point_2 operator()(const FT& x, const FT& y, const FT& z) { Point_2 p; Direction_3& d(p); d = Direction_3(x, y, z); @@ -449,8 +509,7 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { * direction. * \param other the other direction */ - Point_2 operator()(const Direction_3& other) - { + Point_2 operator()(const Direction_3& other) { Point_2 p; Direction_3& d(p); d = Direction_3(other); @@ -461,8 +520,7 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { /*! Initialize a point on the sphere, * \param[in] p the point to initialize. */ - void init(Point_2& p, std::true_type) const - { + void init(Point_2& p, std::true_type) const { const Direction_3& dir = p; if (y_sign(dir) != ZERO) { p.set_location(Point_2::NO_BOUNDARY_LOC); @@ -480,8 +538,7 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { /*! Initialize a point on the sphere, * \param[in] p the point to initialize. */ - void init(Point_2& p, std::false_type) const - { + void init(Point_2& p, std::false_type) const { const Direction_3& dir = p; if ((x_sign(dir) == ZERO) && (y_sign(dir) == ZERO)) { typename Point_2::Location_type location = (z_sign(dir) == NEGATIVE) ? @@ -508,9 +565,9 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { /*! A functor that constructs an x-monotone geodesic arc on the sphere. */ class Construct_x_monotone_curve_2 { protected: - typedef Arr_geodesic_arc_on_sphere_traits_2 Traits; + using Traits = Arr_geodesic_arc_on_sphere_traits_2; - /*! The traits (in case it has state) */ + //! The traits (in case it has state) const Traits& m_traits; /*! Constructor @@ -536,9 +593,8 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { * \pre the source and target must not coincide. * \pre the source and target cannot be antipodal. */ - X_monotone_curve_2 operator()(const Point_2& source, const Point_2& target) - const - { + X_monotone_curve_2 operator()(const Point_2& source, + const Point_2& target) const { X_monotone_curve_2 xcv; xcv.set_source(source); @@ -561,8 +617,7 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { * \param plane the containing plane. * \pre the plane is not vertical */ - X_monotone_curve_2 operator()(const Direction_3& normal) const - { + X_monotone_curve_2 operator()(const Direction_3& normal) const { X_monotone_curve_2 xcv; xcv.set_normal(normal); @@ -608,12 +663,11 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { * \param target the target point. * \pre the source and target cannot be equal. */ - void init(X_monotone_curve_2& xcv) const - { + void init(X_monotone_curve_2& xcv) const { const Point_2& source = xcv.source(); const Point_2& target = xcv.target(); - typedef Arr_geodesic_arc_on_sphere_traits_2 Traits; + using Traits = Arr_geodesic_arc_on_sphere_traits_2; const Kernel& kernel(m_traits); CGAL_precondition(!kernel.equal_3_object()(Direction_3(source), @@ -701,9 +755,9 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { /*! A functor that constructs a geodesic arc on the sphere. */ class Construct_curve_2 { protected: - typedef Arr_geodesic_arc_on_sphere_traits_2 Traits; + using Traits = Arr_geodesic_arc_on_sphere_traits_2; - /*! The traits (in case it has state) */ + //! The traits (in case it has state) const Traits& m_traits; /*! Constructor @@ -728,8 +782,7 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { * \pre the source and target cannot be equal. * \pre the source and target cannot be the opoosite of each other. */ - Curve_2 operator()(const Point_2& source, const Point_2& target) - { + Curve_2 operator()(const Point_2& source, const Point_2& target) { Curve_2 cv; cv.set_source(source); cv.set_target(target); @@ -794,8 +847,7 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { Orientation orient = m_traits.orientation(s, t); const Kernel& kernel = m_traits; - typename Kernel::Counterclockwise_in_between_2 cc_in_between = - kernel.counterclockwise_in_between_2_object(); + auto cc_in_between = kernel.counterclockwise_in_between_2_object(); const Direction_2& d = Traits::identification_xy(); if (orient == LEFT_TURN) { @@ -819,8 +871,7 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { * \pre Both endpoints lie on the given plane. */ Curve_2 operator()(const Point_2& source, const Point_2& target, - const Direction_3& normal) - { + const Direction_3& normal) { Curve_2 cv; cv.set_source(source); @@ -829,7 +880,7 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { cv.set_is_degenerate(false); cv.set_is_empty(false); - typedef Arr_geodesic_arc_on_sphere_traits_2 Traits; + using Traits = Arr_geodesic_arc_on_sphere_traits_2; CGAL_precondition(m_traits.has_on(normal, source)); CGAL_precondition(m_traits.has_on(normal, target)); @@ -912,8 +963,7 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { Direction_2 s = project(source); Direction_2 t = project(target); const Direction_2& ny = Traits::neg_y_2(); - typename Kernel::Counterclockwise_in_between_2 ccib = - kernel.counterclockwise_in_between_2_object(); + auto ccib = kernel.counterclockwise_in_between_2_object(); cv.set_is_x_monotone((plane_is_positive && !ccib(ny, s, t)) || (!plane_is_positive && !ccib(ny, t, s))); @@ -929,8 +979,7 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { const Direction_2& d = Traits::identification_xy(); Direction_2 s = Traits::project_xy(source); Direction_2 t = Traits::project_xy(target); - typename Kernel::Counterclockwise_in_between_2 ccib = - kernel.counterclockwise_in_between_2_object(); + auto ccib = kernel.counterclockwise_in_between_2_object(); bool plane_is_positive = (z_sign(normal) == POSITIVE); cv.set_is_x_monotone((plane_is_positive && !ccib(d, s, t)) || (!plane_is_positive && !ccib(d, t, s))); @@ -961,9 +1010,9 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { /*! A functor that compares the x-coordinates of two directional points */ class Compare_x_2 { protected: - typedef Arr_geodesic_arc_on_sphere_traits_2 Traits; + using Traits = Arr_geodesic_arc_on_sphere_traits_2; - /*! The traits (in case it has state) */ + //! The traits (in case it has state) const Traits& m_traits; /*! Constructor @@ -983,8 +1032,7 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { * \pre p1 does not lie on the boundary. * \pre p2 does not lie on the boundary. */ - Comparison_result operator()(const Point_2& p1, const Point_2& p2) const - { + Comparison_result operator()(const Point_2& p1, const Point_2& p2) const { CGAL_precondition(p1.is_no_boundary()); CGAL_precondition(p2.is_no_boundary()); @@ -996,8 +1044,7 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { /*! Obtain the positive (north) pole * \return the positive (north) pole */ - inline static const Point_2& pos_pole() - { + inline static const Point_2& pos_pole() { static const Point_2 p(Direction_3(0, 0, 1), Point_2::MAX_BOUNDARY_LOC); return p; } @@ -1005,8 +1052,7 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { /*! Obtain the negative (south) pole * \return the negative (south) pole */ - inline static const Point_2& neg_pole() - { + inline static const Point_2& neg_pole() { static const Point_2 p(Direction_3(0, 0, -1), Point_2::MIN_BOUNDARY_LOC); return p; } @@ -1020,9 +1066,9 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { */ class Compare_xy_2 { protected: - typedef Arr_geodesic_arc_on_sphere_traits_2 Traits; + using Traits = Arr_geodesic_arc_on_sphere_traits_2; - /*! The traits (in case it has state) */ + //! The traits (in case it has state) const Traits& m_traits; /*! Constructor @@ -1044,8 +1090,7 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { * \pre p1 does not lie on the boundary. * \pre p2 does not lie on the boundary. */ - Comparison_result operator()(const Point_2& p1, const Point_2& p2) const - { + Comparison_result operator()(const Point_2& p1, const Point_2& p2) const { CGAL_precondition(p1.is_no_boundary()); CGAL_precondition(p2.is_no_boundary()); @@ -1094,8 +1139,7 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { * \return true if the curve is a vertical spherical_arc; false otherwise. * \pre the arc is not degenerate (consists of a single point) */ - bool operator()(const X_monotone_curve_2& xc) const - { + bool operator()(const X_monotone_curve_2& xc) const { CGAL_precondition(!xc.is_degenerate()); return xc.is_vertical(); } @@ -1109,9 +1153,9 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { */ class Compare_y_at_x_2 { protected: - typedef Arr_geodesic_arc_on_sphere_traits_2 Traits; + using Traits = Arr_geodesic_arc_on_sphere_traits_2; - /*! The traits (in case it has state) */ + //! The traits (in case it has state) const Traits& m_traits; /*! Constructor @@ -1132,8 +1176,7 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { * \pre p is in the x-range of xc. */ Comparison_result operator()(const Point_2& p, - const X_monotone_curve_2& xc) const - { + const X_monotone_curve_2& xc) const { CGAL_precondition(!p.is_min_boundary() && !p.is_max_boundary()); CGAL_precondition(m_traits.is_in_x_range(xc, p)); @@ -1169,9 +1212,9 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { */ class Compare_y_at_x_left_2 { protected: - typedef Arr_geodesic_arc_on_sphere_traits_2 Traits; + using Traits = Arr_geodesic_arc_on_sphere_traits_2; - /*! The traits (in case it has state) */ + //! The traits (in case it has state) const Traits& m_traits; /*! Constructor @@ -1195,8 +1238,7 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { */ Comparison_result operator()(const X_monotone_curve_2& xc1, const X_monotone_curve_2& xc2, - const Point_2& p) const - { + const Point_2& p) const { CGAL_precondition(! xc1.is_degenerate()); CGAL_precondition(! xc2.is_degenerate()); CGAL_precondition(p == xc1.right()); @@ -1227,8 +1269,7 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { auto opposite_3 = kernel.construct_opposite_direction_3_object(); Direction_3 opposite_p = opposite_3(p); if (kernel.equal_3_object()(opposite_p, Direction_3(l1)) || - kernel.equal_3_object()(opposite_p, Direction_3(l2))) - { + kernel.equal_3_object()(opposite_p, Direction_3(l2))) { Sign xsign = Traits::x_sign(p); Sign ysign = Traits::y_sign(p); Project project = (xsign == ZERO) ? @@ -1286,9 +1327,9 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { */ class Compare_y_at_x_right_2 { protected: - typedef Arr_geodesic_arc_on_sphere_traits_2 Traits; + using Traits = Arr_geodesic_arc_on_sphere_traits_2; - /*! The traits (in case it has state) */ + //! The traits (in case it has state) const Traits& m_traits; /*! Constructor @@ -1312,8 +1353,7 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { */ Comparison_result operator()(const X_monotone_curve_2& xc1, const X_monotone_curve_2& xc2, - const Point_2& p) const - { + const Point_2& p) const { CGAL_precondition(! xc1.is_degenerate()); CGAL_precondition(! xc2.is_degenerate()); @@ -1335,8 +1375,7 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { auto opposite_3 = kernel.construct_opposite_direction_3_object(); Direction_3 opposite_p = opposite_3(p); if (kernel.equal_3_object()(opposite_p, Direction_3(r1)) || - kernel.equal_3_object()(opposite_p, Direction_3(r2))) - { + kernel.equal_3_object()(opposite_p, Direction_3(r2))) { Sign xsign = Traits::x_sign(p); Sign ysign = Traits::y_sign(p); Project project = (xsign == ZERO) ? @@ -1402,9 +1441,9 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { */ class Equal_2 { protected: - typedef Arr_geodesic_arc_on_sphere_traits_2 Traits; + using Traits = Arr_geodesic_arc_on_sphere_traits_2; - /*! The traits (in case it has state) */ + //! The traits (in case it has state) const Traits& m_traits; /*! Constructor @@ -1422,10 +1461,9 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { * \return true if the two curves are the same; false otherwise. */ bool operator()(const X_monotone_curve_2& xc1, - const X_monotone_curve_2& xc2) const - { + const X_monotone_curve_2& xc2) const { const Kernel& kernel = m_traits; - typename Kernel::Equal_3 equal_3 = kernel.equal_3_object(); + auto equal_3 = kernel.equal_3_object(); if (xc1.is_full() || xc2.is_full()) { if (!xc1.is_full() || !xc2.is_full()) return false; auto opposite_3 = kernel.construct_opposite_direction_3_object(); @@ -1447,8 +1485,7 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { * \param p2 the second point. * \return true if the two point are the same; false otherwise. */ - bool operator()(const Point_2& p1, const Point_2& p2) const - { + bool operator()(const Point_2& p1, const Point_2& p2) const { const Kernel& kernel = m_traits; return kernel.equal_3_object()(Direction_3(p1), Direction_3(p2)); } @@ -1466,9 +1503,9 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { */ class Parameter_space_in_x_2 { protected: - typedef Arr_geodesic_arc_on_sphere_traits_2 Traits; + using Traits = Arr_geodesic_arc_on_sphere_traits_2; - /*! The traits (in case it has state) */ + //! The traits (in case it has state) const Traits& m_traits; /*! Constructor @@ -1498,8 +1535,7 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { * \pre xcv does not coincide with the identification */ Arr_parameter_space operator()(const X_monotone_curve_2& xcv, - Arr_curve_end ce) const - { + Arr_curve_end ce) const { CGAL_precondition(!m_traits.is_on_y_identification_2_object()(xcv)); // vertical, but not on identification! if (xcv.is_vertical()) return ARR_INTERIOR; @@ -1516,8 +1552,7 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { * \param p the point. * \return the parameter space at p. */ - Arr_parameter_space operator()(const Point_2& p) const - { + Arr_parameter_space operator()(const Point_2& p) const { CGAL_precondition(p.is_no_boundary()); CGAL_USE(p); return ARR_INTERIOR; @@ -1551,8 +1586,7 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { * maximal end. */ Arr_parameter_space operator()(const X_monotone_curve_2& xcv, - Arr_curve_end ce) const - { + Arr_curve_end ce) const { return (ce == ARR_MIN_END) ? ((xcv.left().is_min_boundary()) ? ARR_BOTTOM_BOUNDARY: ARR_INTERIOR) : ((xcv.right().is_max_boundary()) ? ARR_TOP_BOUNDARY : ARR_INTERIOR); @@ -1565,8 +1599,7 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { * \param p the point. * \return the parameter space at p. */ - Arr_parameter_space operator()(const Point_2& p) const - { + Arr_parameter_space operator()(const Point_2& p) const { return (p.is_min_boundary()) ? ARR_BOTTOM_BOUNDARY : (p.is_max_boundary()) ? ARR_TOP_BOUNDARY : ARR_INTERIOR; @@ -1583,9 +1616,9 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { */ class Compare_x_on_boundary_2 { protected: - typedef Arr_geodesic_arc_on_sphere_traits_2 Traits; + using Traits = Arr_geodesic_arc_on_sphere_traits_2; - /*! The traits (in case it has state) */ + //! The traits (in case it has state) const Traits& m_traits; /*! Constructor @@ -1613,8 +1646,8 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { */ Comparison_result operator()(const Point_2& point, const X_monotone_curve_2& xcv, - Arr_curve_end CGAL_precondition_code(ce)) const - { + Arr_curve_end CGAL_precondition_code(ce)) + const { CGAL_precondition(point.is_no_boundary()); CGAL_precondition_code (const Point_2& p2 = (ce == ARR_MIN_END) ? xcv.left() : xcv.right();); @@ -1659,8 +1692,8 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { Comparison_result operator()(const X_monotone_curve_2& xcv1, Arr_curve_end CGAL_precondition_code(ce1), const X_monotone_curve_2& xcv2, - Arr_curve_end CGAL_precondition_code(ce2)) const - { + Arr_curve_end CGAL_precondition_code(ce2)) + const { CGAL_precondition_code (const Point_2& p1 = (ce1 == ARR_MIN_END) ? xcv1.left() : xcv1.right();); CGAL_precondition(!p1.is_no_boundary()); @@ -1695,9 +1728,7 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { */ Comparison_result operator()(const Point_2& /* p1 */, const Point_2& /* p2 */) const - { - CGAL_error(); return EQUAL; - } + { CGAL_error(); return EQUAL; } }; /*! Obtain a Compare_x_on_boundary_2 function object. @@ -1710,9 +1741,9 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { */ class Compare_x_near_boundary_2 { protected: - typedef Arr_geodesic_arc_on_sphere_traits_2 Traits; + using Traits = Arr_geodesic_arc_on_sphere_traits_2; - /*! The traits (in case it has state) */ + //! The traits (in case it has state) const Traits& m_traits; /*! Constructor @@ -1746,8 +1777,8 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { CGAL_precondition_code(xcv1), const X_monotone_curve_2& CGAL_precondition_code(xcv2), - Arr_curve_end CGAL_precondition_code(ce)) const - { + Arr_curve_end CGAL_precondition_code(ce)) + const { CGAL_precondition_code (const Point_2& p1 = (ce == ARR_MIN_END) ? xcv1.left() : xcv1.right();); CGAL_precondition(!p1.is_no_boundary()); @@ -1774,9 +1805,9 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { */ class Compare_y_near_boundary_2 { protected: - typedef Arr_geodesic_arc_on_sphere_traits_2 Traits; + using Traits = Arr_geodesic_arc_on_sphere_traits_2; - /*! The traits (in case it has state) */ + //! The traits (in case it has state) const Traits& m_traits; /*! Constructor @@ -1800,8 +1831,7 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { */ Comparison_result operator()(const X_monotone_curve_2& xcv1, const X_monotone_curve_2& xcv2, - Arr_curve_end ce) const - { + Arr_curve_end ce) const { CGAL_precondition(! xcv1.is_degenerate()); CGAL_precondition(! xcv2.is_degenerate()); @@ -1918,9 +1948,9 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { */ class Is_on_y_identification_2 { protected: - typedef Arr_geodesic_arc_on_sphere_traits_2 Traits; + using Traits = Arr_geodesic_arc_on_sphere_traits_2; - /*! The traits (in case it has state) */ + //! The traits (in case it has state) const Traits& m_traits; /*! Constructor @@ -1944,8 +1974,7 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { * \return a Boolean indicating whether xcv coincides with the vertical * identification arc. */ - bool operator()(const X_monotone_curve_2& xcv) const - { + bool operator()(const X_monotone_curve_2& xcv) const { /* If the curve is not vertical and non of its endpoints lie on the * boundary, the arc itself cannot lie on the identification arc. */ @@ -1977,9 +2006,9 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { */ class Compare_y_on_boundary_2 { protected: - typedef Arr_geodesic_arc_on_sphere_traits_2 Traits; + using Traits = Arr_geodesic_arc_on_sphere_traits_2; - /*! The traits (in case it has state) */ + //! The traits (in case it has state) const Traits& m_traits; /*! Constructor @@ -2000,8 +2029,7 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { * \pre p1 lies on the vertical identification arc including the poles! * \pre p2 lies on the vertical identification arc including the poles! */ - Comparison_result operator()(const Point_2& p1, const Point_2& p2) const - { + Comparison_result operator()(const Point_2& p1, const Point_2& p2) const { // first deal with the 'degenerate' case of poles! if (p1.is_min_boundary()) { if (p2.is_min_boundary()) return EQUAL; @@ -2039,9 +2067,9 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { */ class Make_x_monotone_2 { protected: - typedef Arr_geodesic_arc_on_sphere_traits_2 Traits; + using Traits = Arr_geodesic_arc_on_sphere_traits_2; - /*! The traits (in case it has state) */ + //! The traits (in case it has state) const Traits& m_traits; /*! Constructor @@ -2063,8 +2091,8 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { */ template OutputIterator operator()(const Curve_2& c, OutputIterator oi) const { - typedef std::variant - Make_x_monotone_result; + using Make_x_monotone_result = std::variant; + // std::cout << "full: " << c.is_full() << std::endl; // std::cout << "vert: " << c.is_vertical() << std::endl; // std::cout << "xmon: " << c.is_x_monotone() << std::endl; @@ -2199,9 +2227,9 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { /*! A functor that splits an x-monotone arc at a directional point. */ class Split_2 { protected: - typedef Arr_geodesic_arc_on_sphere_traits_2 Traits; + using Traits = Arr_geodesic_arc_on_sphere_traits_2; - /*! The traits (in case it has state) */ + //! The traits (in case it has state) const Traits& m_traits; /*! Constructor @@ -2223,14 +2251,12 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { * \pre xc is not degenerate */ void operator()(const X_monotone_curve_2& xc, const Point_2& p, - X_monotone_curve_2& xc1, X_monotone_curve_2& xc2) const - { + X_monotone_curve_2& xc1, X_monotone_curve_2& xc2) const { CGAL_precondition(!xc.is_degenerate()); const Point_2& source = xc.source(); const Point_2& target = xc.target(); CGAL_precondition_code(const Kernel& kernel = m_traits); - CGAL_precondition_code - (typename Kernel::Equal_3 equal_3 = kernel.equal_3_object()); + CGAL_precondition_code(auto equal_3 = kernel.equal_3_object()); CGAL_precondition(!equal_3(Direction_3(p), Direction_3(source))); CGAL_precondition(!equal_3(Direction_3(p), Direction_3(target))); @@ -2268,9 +2294,9 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { /*! The clockwise-in-between function object */ class Clockwise_in_between_2 { protected: - typedef Arr_geodesic_arc_on_sphere_traits_2 Traits; + using Traits = Arr_geodesic_arc_on_sphere_traits_2; - /*! The traits (in case it has state) */ + //! The traits (in case it has state) const Traits& m_traits; /*! Constructor @@ -2282,8 +2308,7 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { public: bool operator()(const Direction_2& d, - const Direction_2& d1, const Direction_2& d2) const - { + const Direction_2& d1, const Direction_2& d2) const { const Kernel& kernel = m_traits; return kernel.counterclockwise_in_between_2_object()(d, d2, d1); } @@ -2318,9 +2343,8 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { bool vertical, const In_between& in_between, Project project, - OutputIterator oi) const - { - typedef std::pair Intersection_point; + OutputIterator oi) const { + using Intersection_point = std::pair; const Kernel& kernel = m_traits; typename Kernel::Equal_2 equal = kernel.equal_2_object(); @@ -2472,8 +2496,7 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { * \pre point lies in the underlying plane of xc. */ bool is_in_between(const Point_2& point, - const X_monotone_curve_2& xc) const - { + const X_monotone_curve_2& xc) const { const Kernel& kernel = m_traits; CGAL_precondition(m_traits.has_on(xc.normal(), point)); @@ -2522,9 +2545,9 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { } protected: - typedef Arr_geodesic_arc_on_sphere_traits_2 Traits; + using Traits = Arr_geodesic_arc_on_sphere_traits_2; - /*! The traits (in case it has state) */ + //! The traits (in case it has state) const Traits& m_traits; /*! Constructor @@ -2547,22 +2570,17 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { template OutputIterator operator()(const X_monotone_curve_2& xc1, const X_monotone_curve_2& xc2, - OutputIterator oi) const - { + OutputIterator oi) const { // std::cout << "xc1: " << xc1 << std::endl // << "xc2: " << xc2 << std::endl; CGAL_precondition(!xc1.is_degenerate()); CGAL_precondition(!xc2.is_degenerate()); - typedef typename Kernel::Counterclockwise_in_between_2 - Counterclockwise_in_between_2; - typedef typename Kernel::Equal_3 Equal_3; - - typedef std::pair Intersection_point; + using Intersection_point = std::pair; const Kernel& kernel = m_traits; - Equal_3 equal_3 = kernel.equal_3_object(); + auto equal_3 = kernel.equal_3_object(); const Direction_3& normal1 = xc1.normal(); const Direction_3& normal2 = xc2.normal(); @@ -2571,7 +2589,7 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { if (equal_3(normal1, normal2) || equal_3(opposite_normal1, normal2)) { // The underlying planes are the same - Counterclockwise_in_between_2 ccib = kernel.counterclockwise_in_between_2_object(); + auto ccib = kernel.counterclockwise_in_between_2_object(); auto cib = m_traits.clockwise_in_between_2_object(); if (xc1.is_vertical()) { @@ -2659,9 +2677,9 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { /*! A functor that tests whether two x-monotone arcs can be merged. */ class Are_mergeable_2 { - typedef Arr_geodesic_arc_on_sphere_traits_2 Traits; + using Traits = Arr_geodesic_arc_on_sphere_traits_2; - /*! The traits (in case it has state) */ + //! The traits (in case it has state) const Traits& m_traits; /*! Constructor @@ -2681,14 +2699,13 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { * 2. share a common endpoint that is not on the identification arc */ bool operator()(const X_monotone_curve_2& xc1, - const X_monotone_curve_2& xc2) const - { + const X_monotone_curve_2& xc2) const { if (xc1.is_empty() || xc2.is_empty()) return true; if ((xc1.is_full() || xc1.is_meridian()) && (xc2.is_full() || xc2.is_meridian())) return false; const Kernel& kernel = m_traits; - typename Kernel::Equal_3 equal = kernel.equal_3_object(); + auto equal = kernel.equal_3_object(); // Down cast to pass to kernel member functions const Direction_3& xc1_left = xc1.left(); @@ -2731,9 +2748,9 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { /*! A functor that merges two x-monotone arcs into one */ class Merge_2 { protected: - typedef Arr_geodesic_arc_on_sphere_traits_2 Traits; + using Traits = Arr_geodesic_arc_on_sphere_traits_2; - /*! The traits (in case it has state) */ + //! The traits (in case it has state) const Traits& m_traits; /*! Constructor @@ -2752,8 +2769,7 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { */ void operator()(const X_monotone_curve_2& xc1, const X_monotone_curve_2& xc2, - X_monotone_curve_2& xc) const - { + X_monotone_curve_2& xc) const { CGAL_precondition (m_traits.are_mergeable_2_object()(xc1, xc2) == true); if (xc1.is_degenerate() || xc1.is_empty()) { @@ -2767,7 +2783,7 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { } const Kernel& kernel = m_traits; - typename Kernel::Equal_3 equal = kernel.equal_3_object(); + auto equal = kernel.equal_3_object(); // Down cast to pass to kernel member functions const Direction_3& xc1_right = xc1.right(); @@ -2832,9 +2848,11 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { /// \name Functor definitions for the landmarks point-location strategy. //@{ - typedef double Approximate_number_type; - typedef CGAL::Cartesian Approximate_kernel; - typedef Approximate_kernel::Point_2 Approximate_point_2; + using Approximate_number_type = double; + using Approximate_kernel = CGAL::Cartesian; + using Approximate_point_2 = Arr_extended_direction_3; + using Approximate_kernel_vector_3 = Approximate_kernel::Vector_3; + using Approximate_kernel_direction_3 = Approximate_kernel::Direction_3; class Approximate_2 { public: @@ -2846,21 +2864,116 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { * approximation of p's y-coordinate (if i == 1). */ Approximate_number_type operator()(const Point_2& p, int i) const { - CGAL_precondition((i == 0) || (i == 1)); - return (i == 0) ? CGAL::to_double(p.x()) : CGAL::to_double(p.y()); + CGAL_precondition((i == 0) || (i == 1) || (i == 2)); + return (i == 0) ? CGAL::to_double(p.dx()) : + ((i == 1) ? CGAL::to_double(p.dy()) : CGAL::to_double(p.dz())); } /*! Obtain an approximation of a point. */ - Approximate_point_2 operator()(const Point_2& p) const - { return Approximate_point_2(operator()(p, 0), operator()(p, 1)); } + Approximate_point_2 operator()(const Point_2& p) const { + Approximate_kernel::Direction_3 dir(operator()(p, 0), operator()(p, 1), + operator()(p, 2)); + auto loc = static_cast(p.location()); + return Approximate_point_2(dir, loc); + } /*! Obtain an approximation of an \f$x\f$-monotone curve. */ template - OutputIterator operator()(const X_monotone_curve_2& /* xcv */, double /* error */, - OutputIterator /* oi */, bool /* l2r */ = true) const { - CGAL_error_msg("Not implemented yet!"); + OutputIterator operator()(const X_monotone_curve_2& xcv, + Approximate_number_type error, + OutputIterator oi, bool l2r = true) const { + const auto& s = xcv.source(); + const auto& t = xcv.target(); + const auto& n = xcv.normal(); + const auto dx = CGAL::to_double(n.dx()); + const auto dy = CGAL::to_double(n.dy()); + const auto dz = CGAL::to_double(n.dz()); + + Approximate_point_2 as, at; + Approximate_kernel_vector_3 vn; + if (xcv.is_directed_right() == l2r) { + // Get the approximate points + as = (*this)(s); + at = (*this)(t); + vn = Approximate_kernel_vector_3(dx, dy, dz); + } + else { + // Get the approximate points + as = (*this)(t); + at = (*this)(s); + vn = Approximate_kernel_vector_3(-dx, -dy, -dz); + } + + // convert the approximate points to vectors with approximate-kernel + auto vs = approximate_vector_3(as); + auto vt = approximate_vector_3(at); + + // normalize the vectors + auto normalize = [](auto& x) { x /= std::sqrt(x.squared_length()); }; + normalize(vs); + normalize(vt); + normalize(vn); + + // Define the spanning vectors of the coordinate system where we are + // going to make the approximation: + auto axis_x = vs; // x-axis will coincide with the vector from the + // origin to the normalized SOURCE-vector + auto axis_z = vn; // this will make sure that the orientation of the + // approximated curve is consistent with the curve + auto axis_y = CGAL::cross_product(axis_z, axis_x); + normalize(axis_y); + + // In this coordinate system the source has local coords (0,0), hence its + // initial angle with the X-axis is 0 degrees (radians) + // Compute the local coordinates and the angle it makes with the X-axis + Approximate_number_type theta; + if (xcv.is_full()) theta = 2.0 * CGAL_PI; + else { + auto ltx = CGAL::scalar_product(axis_x, vt); + auto lty = CGAL::scalar_product(axis_y, vt); + theta = std::atan2(lty, ltx); + if (theta < 0) + theta += 2.0 * CGAL_PI; + } + + // compute the number of divisions given the requested error + const Approximate_number_type R = 1.0; // radius is always 1 + Approximate_number_type dtheta = 2.0 * std::acos(1 - error / R); + int num_segs = std::ceil(theta / dtheta); + dtheta = theta / num_segs; + + // generate the points approximating the curve + const auto loc = Approximate_point_2::NO_BOUNDARY_LOC; + *oi++ = approximate_point_2(vs, loc); // source vector + for (int i = 1; i < num_segs; ++i) { + const Approximate_number_type angle = i * dtheta; + auto p = std::cos(angle) * axis_x + std::sin(angle) * axis_y; + *oi++ = approximate_point_2(p, loc); + } + *oi++ = approximate_point_2(vt, loc); // target vector + + return oi; + } + + private: + Approximate_kernel_vector_3 + approximate_vector_3(const Approximate_point_2& p) const + { return Approximate_kernel_vector_3(p.dx(), p.dy(), p.dz()); }; + + Approximate_kernel_vector_3 + approximate_vector_3(const Direction_3& d) const { + return Approximate_kernel_vector_3(CGAL::to_double(d.dx()), + CGAL::to_double(d.dy()), + CGAL::to_double(d.dz())); + }; + + Approximate_point_2 + approximate_point_2(const Approximate_kernel_vector_3& v, + const Approximate_point_2::Location_type loc) const { + Approximate_kernel_direction_3 d(v.x(), v.y(), v.z()); + return Approximate_point_2(d, loc); } }; @@ -2874,9 +2987,7 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { class Compare_endpoints_xy_2 { public: - - /*! - * Compare the endpoints of an $x$-monotone curve lexicographically. + /*! Compare the endpoints of an $x$-monotone curve lexicographically. * (assuming the curve has a designated source and target points). * \param xc the curve. * \return SMALLER if the curve is directed right; @@ -2908,8 +3019,7 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { #if 0 /*! Inserter for the spherical_arc class used by the traits-class */ template - friend OutputStream& operator<<(OutputStream& os, const Point_2& p) - { + friend OutputStream& operator<<(OutputStream& os, const Point_2& p) { CGAL::To_double todouble; os << static_cast(todouble(p.dx())) << ", " << static_cast(todouble(p.dy())) << ", " @@ -2920,101 +3030,20 @@ class Arr_geodesic_arc_on_sphere_traits_2 : public Kernel_ { /*! Inserter for the spherical_arc class used by the traits-class */ template friend OutputStream& operator<<(OutputStream& os, - const X_monotone_curve_2& xc) - { + const X_monotone_curve_2& xc) { os << "(" << xc.left() << "), (" << xc.right() << ")"; return os; } /*! Extractor for the spherical_arc class used by the traits-class */ template - friend InputStream& operator>>(InputStream& is, X_monotone_curve_2& arc) - { + friend InputStream& operator>>(InputStream& is, X_monotone_curve_2& arc) { CGAL_error_msg("Not implemented yet!"); return is; } #endif }; -/*! Represent an extended 3D direction that is used in turn to represent a - * spherical-arc endpoint. The extended data consists of two flags that - * indicate whether the point is on the x and on a y boundaries, - * respectively. - */ -template -class Arr_extended_direction_3 : public Kernel::Direction_3 { -public: - typedef typename Kernel::FT FT; - typedef typename Kernel::Direction_3 Direction_3; - - /*! Enumeration of discontinuity type */ - enum Location_type { - NO_BOUNDARY_LOC = 0, - MIN_BOUNDARY_LOC, - MID_BOUNDARY_LOC, - MAX_BOUNDARY_LOC - }; - -private: - typedef typename Kernel::Direction_2 Direction_2; - - /*! The point discontinuity type */ - Location_type m_location; - - inline Sign x_sign(Direction_3 d) const { return CGAL::sign(d.dx()); } - - inline Sign y_sign(Direction_3 d) const { return CGAL::sign(d.dy()); } - - inline Sign z_sign(Direction_3 d) const { return CGAL::sign(d.dz()); } - -public: - /*! Default constructor */ - Arr_extended_direction_3() : - Direction_3(0, 0, 1), - m_location(MAX_BOUNDARY_LOC) - {} - - /*! Constructor */ - Arr_extended_direction_3(const Direction_3& dir, Location_type location) : - Direction_3(dir), - m_location(location) - {} - - /*! Copy constructor */ - Arr_extended_direction_3(const Arr_extended_direction_3& other) : - Direction_3(static_cast(other)) - { m_location = other.discontinuity_type(); } - - /*! Assignment operator */ - Arr_extended_direction_3& operator=(const Arr_extended_direction_3& other) - { - *(static_cast(this)) = static_cast(other); - m_location = other.discontinuity_type(); - return (*this); - } - - /*! Set the location of the point. - */ - void set_location(Location_type location) { m_location = location; } - - /*! Obtain the location of the point. - */ - Location_type location() const { return m_location; } - - /*! Obtain the discontinuity type of the point. - * \todo deprecate this one; use the above instead. - */ - Location_type discontinuity_type() const { return m_location; } - - bool is_no_boundary() const { return (m_location == NO_BOUNDARY_LOC); } - - bool is_min_boundary() const { return (m_location == MIN_BOUNDARY_LOC); } - - bool is_mid_boundary() const { return (m_location == MID_BOUNDARY_LOC); } - - bool is_max_boundary() const { return (m_location == MAX_BOUNDARY_LOC); } -}; - /*! A Representation of an x-monotone great circular arc embedded on a sphere, * as used by the Arr_geodesic_arc_on_sphere_traits_2 traits-class * An x-monotone great circular arc cannot cross the closed hemi-circle arc of @@ -3026,38 +3055,38 @@ class Arr_extended_direction_3 : public Kernel::Direction_3 { template class Arr_x_monotone_geodesic_arc_on_sphere_3 { public: - typedef Kernel_ Kernel; - typedef typename Kernel::Direction_3 Direction_3; - typedef typename Kernel::Plane_3 Plane_3; - typedef typename Kernel::Vector_3 Vector_3; - typedef typename Kernel::Direction_2 Direction_2; + using Kernel = Kernel_; + using Direction_3 = typename Kernel::Direction_3; + using Plane_3 = typename Kernel::Plane_3; + using Vector_3 = typename Kernel::Vector_3; + using Direction_2 = typename Kernel::Direction_2; protected: // For some reason compilation under Windows fails without the qualifier - typedef CGAL::Arr_extended_direction_3 Arr_extended_direction_3; + using Arr_extended_direction_3 = CGAL::Arr_extended_direction_3; - /*! The source point of the arc. */ + //! The source point of the arc. Arr_extended_direction_3 m_source; - /*! The target point of the arc. */ + //! The target point of the arc. Arr_extended_direction_3 m_target; - /*! The direction of the plane that contains the arc. */ + //! The direction of the plane that contains the arc. Direction_3 m_normal; - /*! The arc is vertical. */ + //! The arc is vertical. bool m_is_vertical; - /*! Target (lexicographically) larger than source. */ + //! Target (lexicographically) larger than source. bool m_is_directed_right; - /*! The arc is a full circle. */ + //! The arc is a full circle. bool m_is_full; - /* The arc is degenerate - it consists of a single point. */ + //! The arc is degenerate - it consists of a single point. bool m_is_degenerate; - /*! The arc is empty. */ + //! The arc is empty. bool m_is_empty; inline Sign x_sign(Direction_3 d) const { return CGAL::sign(d.dx()); } @@ -3106,8 +3135,7 @@ class Arr_x_monotone_geodesic_arc_on_sphere_3 { * \param other the other arc */ Arr_x_monotone_geodesic_arc_on_sphere_3 - (const Arr_x_monotone_geodesic_arc_on_sphere_3& other) - { + (const Arr_x_monotone_geodesic_arc_on_sphere_3& other) { m_source = other.m_source; m_target = other.m_target; m_normal = other.m_normal; @@ -3120,8 +3148,7 @@ class Arr_x_monotone_geodesic_arc_on_sphere_3 { /*! Assignment operator */ Arr_x_monotone_geodesic_arc_on_sphere_3& operator= - (const Arr_x_monotone_geodesic_arc_on_sphere_3& other) - { + (const Arr_x_monotone_geodesic_arc_on_sphere_3& other) { m_source = other.m_source; m_target = other.m_target; m_normal = other.m_normal; @@ -3147,9 +3174,8 @@ class Arr_x_monotone_geodesic_arc_on_sphere_3 { * \param target the target point. * \pre the source and target cannot be equal. */ - void init() - { - typedef Arr_geodesic_arc_on_sphere_traits_2 Traits; + void init() { + using Traits = Arr_geodesic_arc_on_sphere_traits_2; Kernel kernel; CGAL_precondition(!kernel.equal_3_object()(Direction_3(m_source), @@ -3238,8 +3264,7 @@ class Arr_x_monotone_geodesic_arc_on_sphere_3 { m_is_directed_right(z_sign(normal) == POSITIVE), m_is_full(true), m_is_degenerate(false), - m_is_empty(false) - { + m_is_empty(false) { CGAL_precondition(z_sign(normal) != ZERO); #if (CGAL_IDENTIFICATION_XY == CGAL_X_MINUS_1_Y_0) @@ -3247,8 +3272,8 @@ class Arr_x_monotone_geodesic_arc_on_sphere_3 { Direction_3(-(normal.dz()), 0, normal.dx()) : Direction_3(normal.dz(), 0, -(normal.dx())); #else - typedef Arr_geodesic_arc_on_sphere_traits_2 Traits; - typedef typename Kernel::FT FT; + using Traits = Arr_geodesic_arc_on_sphere_traits_2; + using FT = typename Kernel::FT; const Direction_2& xy = Traits::identification_xy(); FT x = xy.dx(); @@ -3274,8 +3299,7 @@ class Arr_x_monotone_geodesic_arc_on_sphere_3 { m_is_directed_right(z_sign(normal) == POSITIVE), m_is_full(true), m_is_degenerate(false), - m_is_empty(false) - { + m_is_empty(false) { CGAL_precondition(has_on(point)); CGAL_precondition(z_sign(normal) != ZERO); #if !defined(CGAL_FULL_X_MONOTONE_GEODESIC_ARC_ON_SPHERE_IS_SUPPORTED) @@ -3300,8 +3324,7 @@ class Arr_x_monotone_geodesic_arc_on_sphere_3 { m_normal(normal), m_is_full(false), m_is_degenerate(false), - m_is_empty(false) - { + m_is_empty(false) { CGAL_precondition(has_on(source)); CGAL_precondition(has_on(target)); @@ -3411,8 +3434,7 @@ class Arr_x_monotone_geodesic_arc_on_sphere_3 { #if 0 /*! Create a bounding box for the spherical_arc */ - Bbox_2 bbox() const - { + Bbox_2 bbox() const { Kernel kernel; Segment_2 seg = kernel.construct_spherical_arc_2_object()(this->m_source, this->m_target); @@ -3421,8 +3443,7 @@ class Arr_x_monotone_geodesic_arc_on_sphere_3 { #endif /*! Flip the spherical_arc (swap it source and target) */ - Arr_x_monotone_geodesic_arc_on_sphere_3 opposite() const - { + Arr_x_monotone_geodesic_arc_on_sphere_3 opposite() const { Arr_x_monotone_geodesic_arc_on_sphere_3 opp; opp.m_source = this->m_target; opp.m_target = this->m_source; @@ -3441,8 +3462,7 @@ class Arr_x_monotone_geodesic_arc_on_sphere_3 { * \return true if dir is contained in plane; false otherwise. * \pre the plane contains the origin. */ - inline bool has_on(const Direction_3& dir) const - { + inline bool has_on(const Direction_3& dir) const { typename Kernel::FT dot = normal().vector() * dir.vector(); return CGAL::sign(dot) == ZERO; } @@ -3460,25 +3480,25 @@ template class Arr_geodesic_arc_on_sphere_3 : public Arr_x_monotone_geodesic_arc_on_sphere_3 { public: - typedef Kernel_ Kernel; + using Kernel = Kernel_; protected: - typedef Arr_x_monotone_geodesic_arc_on_sphere_3 Base; + using Base = Arr_x_monotone_geodesic_arc_on_sphere_3; public: - typedef typename Base::Plane_3 Plane_3; - typedef typename Base::Direction_3 Direction_3; - typedef typename Base::Direction_2 Direction_2; + using Plane_3 = typename Base::Plane_3; + using Direction_3 = typename Base::Direction_3; + using Direction_2 = typename Base::Direction_2; protected: // For some reason compilation under Windows fails without the qualifier - typedef CGAL::Arr_extended_direction_3 Arr_extended_direction_3; + using Arr_extended_direction_3 = CGAL::Arr_extended_direction_3; using Base::x_sign; using Base::y_sign; using Base::z_sign; - /*! Indicates whether the arc is x-monotone */ + //! Indicates whether the arc is x-monotone bool m_is_x_monotone; public: @@ -3534,8 +3554,7 @@ class Arr_geodesic_arc_on_sphere_3 : */ Arr_geodesic_arc_on_sphere_3(const Arr_extended_direction_3& source, const Arr_extended_direction_3& target, - const Direction_3& normal) - { + const Direction_3& normal) { Kernel kernel; this->set_source(source); @@ -3544,7 +3563,7 @@ class Arr_geodesic_arc_on_sphere_3 : this->set_is_degenerate(false); this->set_is_empty(false); - typedef Arr_geodesic_arc_on_sphere_traits_2 Traits; + using Traits = Arr_geodesic_arc_on_sphere_traits_2; CGAL_precondition(this->has_on(source)); CGAL_precondition(this->has_on(target)); @@ -3624,8 +3643,7 @@ class Arr_geodesic_arc_on_sphere_3 : Direction_2 s = project(source); Direction_2 t = project(target); const Direction_2& ny = Traits::neg_y_2(); - typename Kernel::Counterclockwise_in_between_2 ccib = - kernel.counterclockwise_in_between_2_object(); + auto ccib = kernel.counterclockwise_in_between_2_object(); set_is_x_monotone((plane_is_positive && !ccib(ny, s, t)) || (!plane_is_positive && !ccib(ny, t, s))); @@ -3641,8 +3659,7 @@ class Arr_geodesic_arc_on_sphere_3 : const Direction_2& d = Traits::identification_xy(); Direction_2 s = Traits::project_xy(source); Direction_2 t = Traits::project_xy(target); - typename Kernel::Counterclockwise_in_between_2 ccib = - kernel.counterclockwise_in_between_2_object(); + auto ccib = kernel.counterclockwise_in_between_2_object(); bool plane_is_positive = (z_sign(normal) == POSITIVE); set_is_x_monotone((plane_is_positive && !ccib(d, s, t)) || (!plane_is_positive && !ccib(d, t, s))); @@ -3651,8 +3668,7 @@ class Arr_geodesic_arc_on_sphere_3 : /*! Construct a full spherical_arc from a normal to a plane. * \param normal the normal to the plane containing the arc. */ - Arr_geodesic_arc_on_sphere_3(const Direction_3& normal) - { + Arr_geodesic_arc_on_sphere_3(const Direction_3& normal) { this->normal(normal); this->set_is_vertical(CGAL::sign(normal.dz()) == ZERO); this->set_is_directed_right(true); @@ -3676,8 +3692,7 @@ class Arr_geodesic_arc_on_sphere_3 : /*! Inserter for the spherical_arc class used by the traits-class */ template OutputStream& operator<<(OutputStream& os, - const Arr_extended_direction_3& ed) -{ + const Arr_extended_direction_3& ed) { #if defined(CGAL_ARR_GEODESIC_ARC_ON_SPHERE_DETAILS) os << "(" << ed.dx() << ", " << ed.dy() << ", " << ed.dz(); @@ -3700,8 +3715,7 @@ OutputStream& operator<<(OutputStream& os, template OutputStream& operator<<(OutputStream& os, - const Arr_x_monotone_geodesic_arc_on_sphere_3& arc) -{ + const Arr_x_monotone_geodesic_arc_on_sphere_3& arc) { #if defined(CGAL_ARR_GEODESIC_ARC_ON_SPHERE_DETAILS) os << "(" << "(" << arc.source() << "), (" << arc.target() << ")" @@ -3721,10 +3735,9 @@ operator<<(OutputStream& os, /*! Extractor for the spherical-arc point class used by the traits-class */ template InputStream& -operator>>(InputStream& is, Arr_extended_direction_3& point) -{ - typedef Kernel_ Kernel; - typedef Arr_extended_direction_3 Point; +operator>>(InputStream& is, Arr_extended_direction_3& point) { + using Kernel = Kernel_; + using Point = Arr_extended_direction_3; // CGAL_error_msg("Importing a geodesic point is not supported!"); typename Kernel::Direction_3 d; is >> d; @@ -3738,10 +3751,9 @@ operator>>(InputStream& is, Arr_extended_direction_3& point) template InputStream& operator>>(InputStream& is, - Arr_x_monotone_geodesic_arc_on_sphere_3& arc) -{ - typedef Kernel_ Kernel; - typedef Arr_extended_direction_3 Point; + Arr_x_monotone_geodesic_arc_on_sphere_3& arc) { + using Kernel = Kernel_; + using Point = Arr_extended_direction_3; // CGAL_error_msg("Importing a geodesic arc is not supported!\n"); diff --git a/Arrangement_on_surface_2/include/CGAL/Arrangement_2/Arrangement_zone_2_impl.h b/Arrangement_on_surface_2/include/CGAL/Arrangement_2/Arrangement_zone_2_impl.h index 5d1e19f73c73..2889f219237f 100644 --- a/Arrangement_on_surface_2/include/CGAL/Arrangement_2/Arrangement_zone_2_impl.h +++ b/Arrangement_on_surface_2/include/CGAL/Arrangement_2/Arrangement_zone_2_impl.h @@ -28,8 +28,7 @@ namespace CGAL { // template void Arrangement_zone_2:: -init_with_hint(const X_monotone_curve_2& cv, Pl_result_type obj) -{ +init_with_hint(const X_monotone_curve_2& cv, Pl_result_type obj) { #if defined(ARR_ZONE_VERBOSE) std::cout << "init_with_hint()" << std::endl; #endif @@ -80,8 +79,7 @@ init_with_hint(const X_monotone_curve_2& cv, Pl_result_type obj) // notifications for the visitor. // template -void Arrangement_zone_2::compute_zone() -{ +void Arrangement_zone_2::compute_zone() { #if defined(ARR_ZONE_VERBOSE) std::cout << "compute_zone()" << std::endl; #endif @@ -137,7 +135,7 @@ void Arrangement_zone_2::compute_zone() if (m_found_overlap) { // In this case m_cv overlaps the curve associated with m_intersect_he. // Compute the overlapping subcurve. - bool dummy; + Arr_parameter_space dummy; auto obj = _compute_next_intersection(m_intersect_he, false, dummy); m_overlap_cv = std::get(*obj); @@ -153,7 +151,7 @@ void Arrangement_zone_2::compute_zone() // overlaps the curve associated with this edge. m_intersect_he = m_arr.non_const_handle(*hh); - bool dummy; + Arr_parameter_space dummy; auto obj = _compute_next_intersection(m_intersect_he, false, dummy); m_overlap_cv = std::get(*obj); @@ -211,7 +209,7 @@ void Arrangement_zone_2::compute_zone() if (m_found_overlap) { // In this case m_cv overlaps the curve associated with m_intersect_he. // Compute the overlapping subcurve to the right of curr_v. - bool dummy; + Arr_parameter_space dummy; auto obj = _compute_next_intersection(m_intersect_he, false, dummy); m_overlap_cv = std::get(*obj); @@ -253,8 +251,7 @@ template bool Arrangement_zone_2:: do_overlap_impl(const X_monotone_curve_2& cv1, const X_monotone_curve_2& cv2, - const Point_2& p, Arr_not_all_sides_oblivious_tag) const -{ + const Point_2& p, Arr_not_all_sides_oblivious_tag) const { typename Traits_adaptor_2::Compare_y_at_x_right_2 cmp_right = m_geom_traits->compare_y_at_x_right_2_object(); @@ -307,8 +304,7 @@ do_overlap_impl(const X_monotone_curve_2& cv1, // template bool Arrangement_zone_2:: -_find_prev_around_vertex(Vertex_handle v, Halfedge_handle& he) -{ +_find_prev_around_vertex(Vertex_handle v, Halfedge_handle& he) { #if defined(ARR_ZONE_VERBOSE) std::cout << "_find_prev_around_vertex(" << v->point() << ")" << std::endl; #endif @@ -393,8 +389,7 @@ typename Arrangement_zone_2::Halfedge_handle Arrangement_zone_2:: _direct_intersecting_edge_to_right(const X_monotone_curve_2& cv_ins, const Point_2& cv_left_pt, - Halfedge_handle query_he) -{ + Halfedge_handle query_he) { #if defined(ARR_ZONE_VERBOSE) std::cout << "_direct_intersecting_edge_to_right() " << cv_left_pt << std::endl; @@ -446,8 +441,7 @@ template typename Arrangement_zone_2::Halfedge_handle Arrangement_zone_2:: _direct_intersecting_edge_to_left(const X_monotone_curve_2& cv_ins, - Halfedge_handle query_he) -{ + Halfedge_handle query_he) { #if defined(ARR_ZONE_VERBOSE) std::cout << "_direct_intersecting_edge_to_left()" << std::endl; #endif @@ -505,6 +499,74 @@ _direct_intersecting_edge_to_left(const X_monotone_curve_2& cv_ins, } } +//! Implementation for no boundary conditions. +template +bool Arrangement_zone_2:: +is_intersection_valid_impl(const Point_2& ip, + Arr_parameter_space& /* intersection_location */, + Arr_all_sides_oblivious_tag) const +{ return (m_geom_traits->compare_xy_2_object()(ip, m_left_pt) == LARGER); } + +//! Implementation for left and right identified boundaries. +template +bool Arrangement_zone_2:: +is_intersection_valid_impl(const Point_2& ip, + Arr_parameter_space& intersection_location, + Arr_has_identified_side_tag) const { + auto equal = m_geom_traits->equal_2_object(); + auto is_on_y_ident = m_geom_traits->is_on_y_identification_2_object(); + // Case 1: the curve lies on the y-identification + if (is_on_y_ident(m_cv)) { + if (equal(ip, m_left_pt)) return false; + // We set the location to be on the left as a convention + intersection_location = ARR_LEFT_BOUNDARY; + return true; + } + + // Case 2: The left-end lies on the left boundary + // (It cannot lie on the right boundary) + // If the intersection point is not equal to the left-end, it must lie to + // its right; thus valid. + if (m_left_on_boundary) return ! equal(ip, m_left_pt); + + // Case 3: The right-end lies on the right boundary + if (m_right_on_boundary && equal(ip, m_right_pt)) { + intersection_location = ARR_RIGHT_BOUNDARY; + return true; + } + + // We have a simple intersection point; + // make sure it lies to the right of m_left_pt. + return (m_geom_traits->compare_xy_2_object()(ip, m_left_pt) == LARGER); +} + +/*! Implementation for all the rest. + * It would be better to split into the various cases, which is the cartesian + * product of (contructed, closed, open) X (contructed, closed, open) + */ +template +bool Arrangement_zone_2:: +is_intersection_valid_impl(const Point_2& ip, + Arr_parameter_space& intersection_location, + Arr_boundary_cond_tag) const { + auto equal = m_geom_traits->equal_2_object(); + if (m_left_on_boundary) { + // The left-end lies on the left boundary. If the intersection point is not + // equal to the left-end, the intersection is valid, because it must lie to + // its right. + if (m_has_left_pt) return ! equal(ip, m_left_pt); + else return true; + } + if (m_has_right_pt && m_right_on_boundary && equal(ip, m_right_pt)) { + intersection_location = ARR_RIGHT_BOUNDARY; + return true; + } + + // We have a simple intersection point; + // make sure it lies to the right of m_left_pt. + return (m_geom_traits->compare_xy_2_object()(ip, m_left_pt) == LARGER); +} + //----------------------------------------------------------------------------- // Get the next intersection of cv with the given halfedge. // @@ -513,14 +575,12 @@ typename Arrangement_zone_2::Optional_intersection Arrangement_zone_2:: _compute_next_intersection(Halfedge_handle he, bool skip_first_point, - bool& intersection_on_right_boundary) -{ + Arr_parameter_space& intersection_location) { #if defined(ARR_ZONE_VERBOSE) std::cout << "_compute_next_intersection(" << he->curve() << ", " << skip_first_point << ")" << std::endl; #endif - auto equal = m_geom_traits->equal_2_object(); auto compare_xy = m_geom_traits->compare_xy_2_object(); auto ctr_min = m_geom_traits->construct_min_vertex_2_object(); auto is_closed = m_geom_traits->is_closed_2_object(); @@ -529,12 +589,12 @@ _compute_next_intersection(Halfedge_handle he, const X_monotone_curve_2* p_curve = &(he->curve()); // Try to locate the intersections with this curve in the intersections map. - Intersect_map_iterator iter = m_inter_map.find(p_curve); + auto iter = m_inter_map.find(p_curve); const Intersection_point* ip; const X_monotone_curve_2* icv; bool valid_intersection; - intersection_on_right_boundary = false; + intersection_location = ARR_INTERIOR; if (iter != m_inter_map.end()) { // The intersections with the curve have already been computed. // Retrieve the intersections list from the map. @@ -548,25 +608,9 @@ _compute_next_intersection(Halfedge_handle he, // Compare that current object with m_left_pt (if exists). ip = std::get_if(&(inter_list.front())); if (ip != nullptr) { - // We have an intersection point - if (m_left_on_boundary) { - // The left-end lies on the left boundary. If the intersection point - // is not equal to the left-end, the intersection is valid, because - // it must lie to its right. - if (m_has_left_pt) valid_intersection = ! equal(ip->first, m_left_pt); - else valid_intersection = true; - } - else if (m_has_right_pt && m_right_on_boundary && - equal(ip->first, m_right_pt)) - { - valid_intersection = true; - intersection_on_right_boundary = true; - } - else { - // We have a simple intersection point - make sure it lies to the - // right of m_left_pt. - valid_intersection = (compare_xy(ip->first, m_left_pt) == LARGER); - } + // We have an intersection point - + valid_intersection = + is_intersection_valid(ip->first, intersection_location); } else { // We have an overlapping subcurve. @@ -576,7 +620,8 @@ _compute_next_intersection(Halfedge_handle he, if (is_closed(*icv, ARR_MIN_END)) { // The curve has a valid left point - make sure it lies to the // right of m_left_pt. - valid_intersection = (compare_xy(ctr_min(*icv), m_left_pt) != SMALLER); + valid_intersection = + (compare_xy(ctr_min(*icv), m_left_pt) != SMALLER); } // In this case the overlap is not valid. else valid_intersection = false; @@ -613,23 +658,8 @@ _compute_next_intersection(Halfedge_handle he, if (ip != nullptr) { // We have an intersection point - // Check whether we need to skip the first intersection - if (is_first && skip_first_point) valid_intersection = false; - else if (m_left_on_boundary) { - // The left-end lies on the left boundary. If the intersection point - // is not equal to the left-end, the intersection is valid, because - // it must lie to its right. - if (m_has_left_pt) valid_intersection = ! equal(ip->first, m_left_pt); - else valid_intersection = true; - } - else if (m_right_on_boundary && m_has_right_pt && - equal(ip->first, m_right_pt)) - { - valid_intersection = true; - intersection_on_right_boundary = true; - } - else { - valid_intersection = (compare_xy(ip->first, m_left_pt) == LARGER); - } + valid_intersection = (is_first && skip_first_point) ? false : + is_intersection_valid(ip->first, intersection_location); } else if (m_left_on_boundary) { // The left end is on the boundary, so all overlapping curves are valid, @@ -671,13 +701,12 @@ _compute_next_intersection(Halfedge_handle he, // template void Arrangement_zone_2:: -_remove_next_intersection(Halfedge_handle he) -{ +_remove_next_intersection(Halfedge_handle he) { // Get a pointer to the curve associated with the halfedge. const X_monotone_curve_2* p_curve = &(he->curve()); // Locate the intersections with this curve in the intersections map. - Intersect_map_iterator iter = m_inter_map.find(p_curve); + auto iter = m_inter_map.find(p_curve); CGAL_assertion(iter != m_inter_map.end()); CGAL_assertion(! iter->second.empty()); @@ -692,8 +721,7 @@ _remove_next_intersection(Halfedge_handle he) template bool Arrangement_zone_2:: _is_to_left_impl(const Point_2& p, Halfedge_handle he, - Arr_not_all_sides_oblivious_tag) const -{ + Arr_not_all_sides_oblivious_tag) const { #if defined(ARR_ZONE_VERBOSE) std::cout << "_is_to_left_impl(" << p << "," << he->curve() << ")" << std::endl; @@ -708,12 +736,12 @@ _is_to_left_impl(const Point_2& p, Halfedge_handle he, if (ps_x_min == ARR_LEFT_BOUNDARY) return false; auto ps_in_y = m_geom_traits->parameter_space_in_y_2_object(); - auto ps_y = ps_in_y(he->curve(), ARR_MIN_END); - if (ps_y != ARR_INTERIOR) { + auto ps_y_min = ps_in_y(he->curve(), ARR_MIN_END); + if (ps_y_min != ARR_INTERIOR) { // Check if p is to the left of the minimal curve-end: auto cmp_x = m_geom_traits->compare_x_point_curve_end_2_object(); const auto res = cmp_x(p, he->curve(), ARR_MIN_END); - return ((res == SMALLER) || (res == EQUAL && ps_y == ARR_TOP_BOUNDARY)); + return ((res == SMALLER) || (res == EQUAL && ps_y_min == ARR_TOP_BOUNDARY)); } // In case the minimal curve-end does not have boundary conditions, simply @@ -724,13 +752,13 @@ _is_to_left_impl(const Point_2& p, Halfedge_handle he, } //----------------------------------------------------------------------------- -// Determine whether a given point lies completely to the right of a given curve. +// Determine whether a given point lies completely to the right of a given +// curve. // template bool Arrangement_zone_2:: _is_to_right_impl(const Point_2& p, Halfedge_handle he, - Arr_not_all_sides_oblivious_tag) const -{ + Arr_not_all_sides_oblivious_tag) const { #if defined(ARR_ZONE_VERBOSE) std::cout << "_is_to_right_impl(" << p << "," << he->curve() << ")" << std::endl; @@ -747,12 +775,12 @@ _is_to_right_impl(const Point_2& p, Halfedge_handle he, if (ps_x_max == ARR_LEFT_BOUNDARY) return true; auto ps_in_y = m_geom_traits->parameter_space_in_y_2_object(); - auto ps_y = ps_in_y(he->curve(), ARR_MAX_END); - if (ps_y != ARR_INTERIOR) { + auto ps_y_max = ps_in_y(he->curve(), ARR_MAX_END); + if (ps_y_max != ARR_INTERIOR) { // Check if p is to the right of the maximal curve-end: auto cmp_x = m_geom_traits->compare_x_point_curve_end_2_object(); auto res = cmp_x(p, he->curve(), ARR_MAX_END); - return ((res == LARGER) || (res == EQUAL && ps_y == ARR_BOTTOM_BOUNDARY)); + return ((res == LARGER) || (res == EQUAL && ps_y_max == ARR_BOTTOM_BOUNDARY)); } // In case the maximal curve-end does not have boundary conditions, simply @@ -768,8 +796,7 @@ _is_to_right_impl(const Point_2& p, Halfedge_handle he, template void Arrangement_zone_2:: _leftmost_intersection(Ccb_halfedge_circulator he_curr, bool on_boundary, - bool& leftmost_on_right_boundary) -{ + Arr_parameter_space& leftmost_location) { #if defined(ARR_ZONE_VERBOSE) std::cout << "_leftmost_intersection(" << he_curr->curve() << ", " << on_boundary << ")" << std::endl; @@ -793,7 +820,7 @@ _leftmost_intersection(Ccb_halfedge_circulator he_curr, bool on_boundary, // entirely to the right of m_intersect_p, its intersection with m_cv (if any) // cannot lie to the left of this point. We therefore do not need to compute // this intersection. - if (m_found_intersect && ! leftmost_on_right_boundary && + if (m_found_intersect && (leftmost_location == ARR_INTERIOR) && _is_to_left(m_intersect_p, he_curr)) return; @@ -821,13 +848,14 @@ _leftmost_intersection(Ccb_halfedge_circulator he_curr, bool on_boundary, // do not intersect. if (! left_equals_curr_endpoint && ((! m_left_on_boundary && _is_to_right(m_left_pt, he_curr)) || - ! is_in_x_range(m_cv, he_curr->curve()))) + ! is_in_x_range(m_cv, he_curr->curve()))) { return; + } // Compute the next intersection of m_cv and the current halfedge. - bool intersection_on_right_boundary; + Arr_parameter_space intersection_location; auto iobj = _compute_next_intersection(he_curr, left_equals_curr_endpoint, - intersection_on_right_boundary); + intersection_location); if (iobj) { // We have found an intersection (either a simple point or an @@ -839,8 +867,8 @@ _leftmost_intersection(Ccb_halfedge_circulator he_curr, bool on_boundary, // Found a simple intersection point. Check if it is the leftmost // intersection point so far. if (! m_found_intersect || - (! intersection_on_right_boundary && - (leftmost_on_right_boundary || + ((intersection_location != ARR_RIGHT_BOUNDARY) && + ((leftmost_location == ARR_RIGHT_BOUNDARY) || compare_xy(ip, m_intersect_p) == SMALLER))) { // Store the leftmost intersection point and the halfedge handle. @@ -848,7 +876,7 @@ _leftmost_intersection(Ccb_halfedge_circulator he_curr, bool on_boundary, m_ip_multiplicity = int_p->second; m_intersect_he = he_curr; m_found_overlap = false; - leftmost_on_right_boundary = intersection_on_right_boundary; + leftmost_location = intersection_location; } } else { @@ -880,8 +908,7 @@ _leftmost_intersection(Ccb_halfedge_circulator he_curr, bool on_boundary, // template void Arrangement_zone_2:: -_leftmost_intersection_with_face_boundary(Face_handle face, bool on_boundary) -{ +_leftmost_intersection_with_face_boundary(Face_handle face, bool on_boundary) { #if defined(ARR_ZONE_VERBOSE) std::cout << "_leftmost_intersection_with_face_boundary(" << on_boundary << ")" << std::endl; @@ -895,25 +922,23 @@ _leftmost_intersection_with_face_boundary(Face_handle face, bool on_boundary) auto compare_xy = m_geom_traits->compare_xy_2_object(); auto is_in_x_range = m_geom_traits->is_in_x_range_2_object(); - bool leftmost_on_right_boundary = false; + Arr_parameter_space leftmost_location = ARR_INTERIOR; // Traverse the face outer-boundaries; iterate through all outer CCBs. for (auto occb_it = face->outer_ccbs_begin(); - occb_it != face->outer_ccbs_end(); ++occb_it) - { + occb_it != face->outer_ccbs_end(); ++occb_it) { Ccb_halfedge_circulator he_first = *occb_it; Ccb_halfedge_circulator he_curr = he_first; - do _leftmost_intersection(he_curr, on_boundary, leftmost_on_right_boundary); + do _leftmost_intersection(he_curr, on_boundary, leftmost_location); while (++he_curr != he_first); } // Traverse the face inner-boundaries; iterate through all inner CCBs (holes). for (auto iccb_it = face->inner_ccbs_begin(); - iccb_it != face->inner_ccbs_end(); ++iccb_it) - { + iccb_it != face->inner_ccbs_end(); ++iccb_it) { Ccb_halfedge_circulator he_first = *iccb_it; Ccb_halfedge_circulator he_curr = he_first; - do _leftmost_intersection(he_curr, on_boundary, leftmost_on_right_boundary); + do _leftmost_intersection(he_curr, on_boundary, leftmost_location); while (++he_curr != he_first); } @@ -921,17 +946,17 @@ _leftmost_intersection_with_face_boundary(Face_handle face, bool on_boundary) // Traverse the isolated vertices inside the face (if there exist any), and // check whether an isolated vertex lies on the curve. - typedef typename Arrangement_2::Isolated_vertex_iterator - Isolated_vertex_iterator; - for (Isolated_vertex_iterator iv_it = face->isolated_vertices_begin(); - iv_it != face->isolated_vertices_end(); ++iv_it) - { + // MSVC17 requires an explicit (non-const) type for the iterator. + typename Arrangement_2::Isolated_vertex_iterator iv_it; + for (iv_it = face->isolated_vertices_begin(); + iv_it != face->isolated_vertices_end(); ++iv_it) { // If the isolated vertex is not in the x-range of our curve, disregard it. if (! is_in_x_range(m_cv, iv_it->point())) continue; // If we already have an intersection point, compare it to the current // isolated vertex, in order to filter unnecessary computations. - if (m_found_intersect && compare_xy(iv_it->point(), m_intersect_p) == LARGER) + if (m_found_intersect && + (compare_xy(iv_it->point(), m_intersect_p) == LARGER)) continue; // In case the isolated vertex lies on the curve, update the intersection @@ -960,8 +985,7 @@ _leftmost_intersection_with_face_boundary(Face_handle face, bool on_boundary) // template bool Arrangement_zone_2:: -_zone_in_face(Face_handle face, bool on_boundary) -{ +_zone_in_face(Face_handle face, bool on_boundary) { #if defined(ARR_ZONE_VERBOSE) std::cout << "_zone_in_face(" << on_boundary << ")" << std::endl; #endif @@ -1095,7 +1119,7 @@ _zone_in_face(Face_handle face, bool on_boundary) // Associate the intersection list of the original curve with the // right subcurve, while we can associate an empty list with the // left subcurve, as we are now done with it. - Intersect_map_iterator iter = m_inter_map.find(p_orig_curve); + auto iter = m_inter_map.find(p_orig_curve); Intersect_list empty_inter_list; m_inter_map[p_right_subcurve] = iter->second; @@ -1190,8 +1214,7 @@ _zone_in_face(Face_handle face, bool on_boundary) // curve currently associated with m_intersect_he. // template -bool Arrangement_zone_2::_zone_in_overlap() -{ +bool Arrangement_zone_2::_zone_in_overlap() { #if defined(ARR_ZONE_VERBOSE) std::cout << "_zone_in_overlap()" << std::endl; #endif diff --git a/Arrangement_on_surface_2/include/CGAL/Arrangement_zone_2.h b/Arrangement_on_surface_2/include/CGAL/Arrangement_zone_2.h index c32deb11ecac..cb0e71178337 100644 --- a/Arrangement_on_surface_2/include/CGAL/Arrangement_zone_2.h +++ b/Arrangement_on_surface_2/include/CGAL/Arrangement_zone_2.h @@ -55,65 +55,68 @@ namespace CGAL { template class Arrangement_zone_2 { public: - typedef Arrangement_ Arrangement_2; - typedef typename Arrangement_2::Geometry_traits_2 Geometry_traits_2; - typedef typename Arrangement_2::Topology_traits Topology_traits; + using Arrangement_2 = Arrangement_; + using Geometry_traits_2 = typename Arrangement_2::Geometry_traits_2; + using Topology_traits = typename Arrangement_2::Topology_traits; protected: - typedef Arr_traits_adaptor_2 Traits_adaptor_2; + using Traits_adaptor_2 = Arr_traits_adaptor_2; - typedef typename Traits_adaptor_2::Left_side_category Left_side_category; - typedef typename Traits_adaptor_2::Bottom_side_category Bottom_side_category; - typedef typename Traits_adaptor_2::Top_side_category Top_side_category; - typedef typename Traits_adaptor_2::Right_side_category Right_side_category; + using Left_side_category = typename Traits_adaptor_2::Left_side_category; + using Bottom_side_category = typename Traits_adaptor_2::Bottom_side_category; + using Top_side_category = typename Traits_adaptor_2::Top_side_category; + using Right_side_category = typename Traits_adaptor_2::Right_side_category; static_assert(Arr_sane_identified_tagging::value); + // Categories for dispatching + using Are_all_sides_oblivious_category = + typename Arr_all_sides_oblivious_category::result; + + using Left_or_right_sides_category = + typename Arr_two_sides_category::result; public: - typedef ZoneVisitor_ Visitor; + using Visitor = ZoneVisitor_; - typedef typename Arrangement_2::Vertex_handle Vertex_handle; - typedef typename Arrangement_2::Halfedge_handle Halfedge_handle; - typedef typename Arrangement_2::Face_handle Face_handle; + using Vertex_handle = typename Arrangement_2::Vertex_handle; + using Halfedge_handle = typename Arrangement_2::Halfedge_handle; + using Face_handle = typename Arrangement_2::Face_handle; - typedef std::pair Visitor_result; + using Visitor_result = std::pair; - typedef typename Geometry_traits_2::Point_2 Point_2; - typedef typename Geometry_traits_2::X_monotone_curve_2 X_monotone_curve_2; - typedef typename Geometry_traits_2::Multiplicity Multiplicity; + using Point_2 = typename Geometry_traits_2::Point_2; + using X_monotone_curve_2 = typename Geometry_traits_2::X_monotone_curve_2; + using Multiplicity = typename Geometry_traits_2::Multiplicity; protected: - typedef typename Arr_all_sides_oblivious_category::result - Are_all_sides_oblivious_category; - - typedef typename Arrangement_2::Vertex_const_handle Vertex_const_handle; - typedef typename Arrangement_2::Halfedge_const_handle Halfedge_const_handle; - typedef typename Arrangement_2::Face_const_handle Face_const_handle; + // General types + using Vertex_const_handle = typename Arrangement_2::Vertex_const_handle; + using Halfedge_const_handle = typename Arrangement_2::Halfedge_const_handle; + using Face_const_handle = typename Arrangement_2::Face_const_handle; - typedef typename Arrangement_2::Ccb_halfedge_circulator - Ccb_halfedge_circulator; + using Ccb_halfedge_circulator = + typename Arrangement_2::Ccb_halfedge_circulator; // Types used for caching intersection points: - typedef std::pair Intersection_point; - typedef std::variant - Intersection_result; - typedef std::optional Optional_intersection; - typedef std::list Intersect_list; - typedef std::map - Intersect_map; - typedef typename Intersect_map::iterator Intersect_map_iterator; + using Intersection_point = std::pair; + using Intersection_result = + std::variant; + using Optional_intersection = std::optional; + using Intersect_list = std::list; + using Intersect_map = std::map; - typedef std::set Curves_set; - typedef typename Curves_set::iterator Curves_set_iterator; + using Curves_set = std::set; + using Curves_set_iterator = typename Curves_set::iterator; - typedef Arr_point_location_result Pl_result; - typedef typename Pl_result::Type Pl_result_type; + using Pl_result = Arr_point_location_result; + using Pl_result_type = typename Pl_result::Type; // Data members: Arrangement_2& m_arr; // The associated arrangement. @@ -193,8 +196,7 @@ class Arrangement_zone_2 { * \param pl A point-location object associated with the arrangement. */ template - void init(const X_monotone_curve_2& cv, const PointLocation& pl) - { + void init(const X_monotone_curve_2& cv, const PointLocation& pl) { // Set the curve and check whether its left end has boundary conditions. m_cv = cv; @@ -267,8 +269,7 @@ class Arrangement_zone_2 { */ bool do_overlap_impl(const X_monotone_curve_2& cv1, const X_monotone_curve_2& cv2, - const Point_2& p, Arr_all_sides_oblivious_tag) const - { + const Point_2& p, Arr_all_sides_oblivious_tag) const { return m_geom_traits->compare_y_at_x_right_2_object()(cv1, cv2, p) == EQUAL; } @@ -336,7 +337,7 @@ class Arrangement_zone_2 { Optional_intersection _compute_next_intersection(Halfedge_handle he, bool skip_first_point, - bool& intersect_on_right_boundary); + Arr_parameter_space& intersection_location); /*! Remove the next intersection of m_cv with the given halfedge from the map. * \param he A handle to the halfedge. @@ -381,8 +382,7 @@ class Arrangement_zone_2 { { return (_is_to_right_impl(p, he, Are_all_sides_oblivious_category())); } bool _is_to_right_impl(const Point_2& p, Halfedge_handle he, - Arr_all_sides_oblivious_tag) const - { + Arr_all_sides_oblivious_tag) const { return (((he->direction() == ARR_LEFT_TO_RIGHT) && m_geom_traits->compare_xy_2_object()(p, he->target()->point()) == LARGER) || @@ -394,12 +394,33 @@ class Arrangement_zone_2 { bool _is_to_right_impl(const Point_2& p, Halfedge_handle he, Arr_not_all_sides_oblivious_tag) const; + /*! Check whether an intersection point is valid. A valid intersection point + * must be to the left of the left end of the curve involved. + */ + bool is_intersection_valid(const Point_2& ip, + Arr_parameter_space& intersection_location) const { + return is_intersection_valid_impl(ip, intersection_location, + Left_or_right_sides_category()); + } + + bool is_intersection_valid_impl(const Point_2& ip, + Arr_parameter_space& intersection_location, + Arr_all_sides_oblivious_tag) const; + + bool is_intersection_valid_impl(const Point_2& ip, + Arr_parameter_space& intersection_location, + Arr_has_identified_side_tag) const; + + bool is_intersection_valid_impl(const Point_2& ip, + Arr_parameter_space& intersection_location, + Arr_boundary_cond_tag) const; + /*! Compute the (lexicographically) leftmost intersection of the query * curve with a given halfedge on the boundary of a face in the arrangement. */ void _leftmost_intersection(Ccb_halfedge_circulator he_curr, bool on_boundary, - bool& leftmost_on_right_boundary); + Arr_parameter_space& leftmost_location); /*! Compute the (lexicographically) leftmost intersection of the query * curve with the boundary of a given face in the arrangement. diff --git a/Arrangement_on_surface_2/include/CGAL/draw_arrangement_2.h b/Arrangement_on_surface_2/include/CGAL/draw_arrangement_2.h index c7637f04b919..76e2eaf32832 100644 --- a/Arrangement_on_surface_2/include/CGAL/draw_arrangement_2.h +++ b/Arrangement_on_surface_2/include/CGAL/draw_arrangement_2.h @@ -29,7 +29,9 @@ #include #include +#include #include +#include namespace CGAL { @@ -46,32 +48,32 @@ struct Default_color_generator { }; // Viewer class for`< Polygon_2 -template -class Arr_2_basic_viewer_qt : public Basic_viewer_qt { - using Arr = Arrangement_2_; +class Aos_2_basic_viewer_qt : public Basic_viewer_qt { + using Aos = ArrangementOnSurface_2; using Color_generator = ColorGenerator; using Base = Basic_viewer_qt; - using Gt = typename Arr::Geometry_traits_2; - using Point = typename Arr::Point_2; - using X_monotone_curve = typename Arr::X_monotone_curve_2; - using Vertex_const_handle = typename Arr::Vertex_const_handle; - using Halfedge_const_handle = typename Arr::Halfedge_const_handle; - using Face_const_handle = typename Arr::Face_const_handle; + using Gt = typename Aos::Geometry_traits_2; + using Point = typename Aos::Point_2; + using X_monotone_curve = typename Aos::X_monotone_curve_2; + using Vertex_const_handle = typename Aos::Vertex_const_handle; + using Halfedge_const_handle = typename Aos::Halfedge_const_handle; + using Face_const_handle = typename Aos::Face_const_handle; using Ccb_halfedge_const_circulator = - typename Arr::Ccb_halfedge_const_circulator; + typename Aos::Ccb_halfedge_const_circulator; public: /// Construct the viewer. /// @param arr the arrangement to view /// @param title the title of the window - Arr_2_basic_viewer_qt(QWidget* parent, const Arr& arr, + Aos_2_basic_viewer_qt(QWidget* parent, const Aos& aos, Color_generator color_generator, const char* title = "2D Arrangement Basic Viewer", bool draw_vertices = false) : // First draw: vertices; edges, faces; multi-color; no inverse normal Base(parent, title, draw_vertices, true, true, false, false), - m_arr(arr), + m_aos(aos), m_color_generator(color_generator) { // mimic the computation of Camera::pixelGLRatio() @@ -154,32 +156,47 @@ class Arr_2_basic_viewer_qt : public Basic_viewer_qt { */ CGAL::Bbox_2 bounding_box() { CGAL::Bbox_2 bbox; - const auto* traits = this->m_arr.geometry_traits(); + const auto* traits = this->m_aos.geometry_traits(); // At this point we assume that the arrangement is not open, and thus the // bounding box is defined by the vertices. - for (auto it = m_arr.vertices_begin(); it != m_arr.vertices_end(); ++it) + for (auto it = m_aos.vertices_begin(); it != m_aos.vertices_end(); ++it) bounding_box_impl1(bbox, it->point(), *traits, 0); return bbox; } + /*! Add all faces. + */ + template + void add_faces(const Traits&) { + for (auto it = m_aos.unbounded_faces_begin(); + it != m_aos.unbounded_faces_end(); ++it) + add_face(it); + } + + /*! Add all faces. + */ + template + void + add_faces(Arr_geodesic_arc_on_sphere_traits_2 const&) + { add_face(m_aos.faces_begin()); } + /*! Add all elements to be drawn. */ void add_elements() { + // std::cout << "add_elements()\n"; // std::cout << "ratio: " << this->pixel_ratio() << std::endl; clear(); m_visited.clear(); - if (m_arr.is_empty()) return; - for (auto it = m_arr.unbounded_faces_begin(); - it != m_arr.unbounded_faces_end(); ++it) - add_face(it); + if (m_aos.is_empty()) return; + add_faces(*(this->m_aos.geometry_traits())); - // Add edges that do not separe faces. - for (auto it = m_arr.edges_begin(); it != m_arr.edges_end(); ++it) + // Add edges that do not separate faces. + for (auto it = m_aos.edges_begin(); it != m_aos.edges_end(); ++it) if (it->face() == it->twin()->face()) draw_curve(it->curve()); // Add all points - for (auto it = m_arr.vertices_begin(); it != m_arr.vertices_end(); ++it) + for (auto it = m_aos.vertices_begin(); it != m_aos.vertices_end(); ++it) draw_point(it->point()); m_visited.clear(); @@ -190,11 +207,20 @@ class Arr_2_basic_viewer_qt : public Basic_viewer_qt { double pixel_ratio() const { return m_pixel_ratio; } protected: + template + Halfedge_const_handle + find_smallest(Ccb_halfedge_const_circulator circ, + Arr_geodesic_arc_on_sphere_traits_2 const&) + { return circ; } + /*! Find the halfedge incident to the lexicographically smallest vertex * along the CCB, such that there is no other halfedge underneath. */ - Halfedge_const_handle find_smallest(Ccb_halfedge_const_circulator circ) { - const auto* traits = this->m_arr.geometry_traits(); + template + Halfedge_const_handle find_smallest(Ccb_halfedge_const_circulator circ, + const Traits&) { + // std::cout << "find_smallest()\n"; + const auto* traits = this->m_aos.geometry_traits(); auto cmp_xy = traits->compare_xy_2_object(); auto cmp_y = traits->compare_y_at_x_right_2_object(); @@ -241,6 +267,7 @@ class Arr_2_basic_viewer_qt : public Basic_viewer_qt { template void draw_approximate_region(Halfedge_const_handle curr, const Approximate& approx) { + // std::cout << "draw_approximate_region()\n"; std::vector polyline; double error(this->pixel_ratio()); bool l2r = curr->direction() == ARR_LEFT_TO_RIGHT; @@ -258,7 +285,7 @@ class Arr_2_basic_viewer_qt : public Basic_viewer_qt { */ template void draw_exact_curve(const XMonotoneCurve& curve) { - const auto* traits = this->m_arr.geometry_traits(); + const auto* traits = this->m_aos.geometry_traits(); auto ctr_min = traits->construct_min_vertex_2_object(); auto ctr_max = traits->construct_max_vertex_2_object(); this->add_segment(ctr_min(curve), ctr_max(curve)); @@ -267,7 +294,7 @@ class Arr_2_basic_viewer_qt : public Basic_viewer_qt { /*! Draw an exact region. */ void draw_exact_region(Halfedge_const_handle curr) { - this->add_point_in_face(curr->source()->point()); + // this->add_point_in_face(curr->source()->point()); draw_exact_curve(curr->curve()); } @@ -300,9 +327,46 @@ class Arr_2_basic_viewer_qt : public Basic_viewer_qt { { draw_approximate_region(curr, traits.approximate_2_object()); } #endif + template + void draw_region_impl1 + (Halfedge_const_handle curr, + Arr_geodesic_arc_on_sphere_traits_2 const& traits, + int) { + // std::cout << "draw_region_impl1()\n"; + auto approx = traits.approximate_2_object(); + using Kernel = Kernel_; + using Traits = Arr_geodesic_arc_on_sphere_traits_2; + using Ak = typename Traits::Approximate_kernel; + using Ap = typename Traits::Approximate_point_2; + using Approx_point_3 = typename Ak::Point_3; + + std::vector polyline; + double error(0.01); + bool l2r = curr->direction() == ARR_LEFT_TO_RIGHT; + approx(curr->curve(), error, std::back_inserter(polyline), l2r); + if (polyline.empty()) return; + auto it = polyline.begin(); + auto x = it->dx(); + auto y = it->dy(); + auto z = it->dz(); + auto l = std::sqrt(x*x + y*y + z*z); + Approx_point_3 prev(x/l, y/l, z/l); + for (++it; it != polyline.end(); ++it) { + auto x = it->dx(); + auto y = it->dy(); + auto z = it->dz(); + auto l = std::sqrt(x*x + y*y + z*z); + Approx_point_3 next(x/l, y/l, z/l); + this->add_segment(prev, next); + prev = next; + // this->add_point_in_face(*prev); + } + } + /*! Draw a region. */ void draw_region(Ccb_halfedge_const_circulator circ) { + // std::cout << "draw_region()\n"; /* Check whether the traits has a member function called * approximate_2_object() and if so check whether the return type, namely * `Approximate_2` has an appropriate operator. @@ -321,8 +385,8 @@ class Arr_2_basic_viewer_qt : public Basic_viewer_qt { auto color = m_color_generator(circ->face()); this->face_begin(color); - const auto* traits = this->m_arr.geometry_traits(); - auto ext = find_smallest(circ); + const auto* traits = this->m_aos.geometry_traits(); + auto ext = find_smallest(circ, *traits); auto curr = ext; do { @@ -382,6 +446,37 @@ class Arr_2_basic_viewer_qt : public Basic_viewer_qt { { draw_approximate_curve(xcv, traits.approximate_2_object()); } #endif + template + void draw_curve_impl1 + (const X_monotone_curve& xcv, + Arr_geodesic_arc_on_sphere_traits_2 const& traits, + int) { + auto approx = traits.approximate_2_object(); + using Kernel = Kernel_; + using Traits = Arr_geodesic_arc_on_sphere_traits_2; + using Ak = typename Traits::Approximate_kernel; + using Ap = typename Traits::Approximate_point_2; + using Approx_point_3 = typename Ak::Point_3; + std::vector apoints; + double error(0.01); + approx(xcv, error, std::back_inserter(apoints)); + auto it = apoints.begin(); + auto x = it->dx(); + auto y = it->dy(); + auto z = it->dz(); + auto l = std::sqrt(x*x + y*y + z*z); + Approx_point_3 prev(x/l, y/l, z/l); + for (++it; it != apoints.end(); ++it) { + auto x = it->dx(); + auto y = it->dy(); + auto z = it->dz(); + auto l = std::sqrt(x*x + y*y + z*z); + Approx_point_3 next(x/l, y/l, z/l); + this->add_segment(prev, next); + prev = next; + } + } + /*! Draw a curve. */ template @@ -404,14 +499,14 @@ class Arr_2_basic_viewer_qt : public Basic_viewer_qt { #if 0 if constexpr (std::experimental::is_detected_v) { - const auto* traits = this->m_arr.geometry_traits(); + const auto* traits = this->m_aos.geometry_traits(); auto approx = traits->approximate_2_object(); draw_approximate_curve(curve, approx); return; } draw_exact_curve(curve); #else - const auto* traits = this->m_arr.geometry_traits(); + const auto* traits = this->m_aos.geometry_traits(); draw_curve_impl1(curve, *traits, 0); #endif } @@ -442,16 +537,35 @@ class Arr_2_basic_viewer_qt : public Basic_viewer_qt { { add_point(traits.approximate_2_object()(p)); } #endif + template + void draw_point_impl1 + (const Point& p, + Arr_geodesic_arc_on_sphere_traits_2 const& traits, + int) { + auto approx = traits.approximate_2_object(); + using Traits = Arr_geodesic_arc_on_sphere_traits_2; + using Ak = typename Traits::Approximate_kernel; + using Approx_point_3 = typename Ak::Point_3; + auto ap = approx(p); + auto x = ap.dx(); + auto y = ap.dy(); + auto z = ap.dz(); + auto l = std::sqrt(x*x + y*y + z*z); + Approx_point_3 p3(x/l, y/l, z/l); + add_point(p3); + } + /*! Draw a point. */ void draw_point(const Point& p) { - const auto* traits = m_arr.geometry_traits(); + const auto* traits = m_aos.geometry_traits(); draw_point_impl1(p, *traits, 0); } /*! Add a Connected Component of the Boundary. */ void add_ccb(Ccb_halfedge_const_circulator circ) { + // std::cout << "add_ccb()\n"; auto curr = circ; do { auto new_face = curr->twin()->face(); @@ -464,8 +578,9 @@ class Arr_2_basic_viewer_qt : public Basic_viewer_qt { /*! Add a face. */ void add_face(Face_const_handle face) { - using Inner_ccb_const_iterator = typename Arr::Inner_ccb_const_iterator; - using Outer_ccb_const_iterator = typename Arr::Outer_ccb_const_iterator; + // std::cout << "add_face()\n"; + using Inner_ccb_const_iterator = typename Aos::Inner_ccb_const_iterator; + using Outer_ccb_const_iterator = typename Aos::Outer_ccb_const_iterator; for (Inner_ccb_const_iterator it = face->inner_ccbs_begin(); it != face->inner_ccbs_end(); ++it) @@ -505,7 +620,7 @@ class Arr_2_basic_viewer_qt : public Basic_viewer_qt { double m_pixel_ratio = 1; //! The arrangement to draw. - const Arr& m_arr; + const Aos& m_aos; //! The color generator. Color_generator m_color_generator; @@ -514,32 +629,63 @@ class Arr_2_basic_viewer_qt : public Basic_viewer_qt { }; //! Basic viewer of a 2D arrangement. -template -class Arr_2_viewer_qt : public Arr_2_basic_viewer_qt { public: - using Arr = Arrangement_2_; + using Aos = ArrangementOnSurface_2; using Color_generator = ColorGenerator; - using Base = Arr_2_basic_viewer_qt; - using Point = typename Arr::Point_2; - using X_monotone_curve = typename Arr::X_monotone_curve_2; - using Halfedge_const_handle = typename Arr::Halfedge_const_handle; - using Face_const_handle = typename Arr::Face_const_handle; + using Base = Aos_2_basic_viewer_qt; + using Point = typename Aos::Point_2; + using X_monotone_curve = typename Aos::X_monotone_curve_2; + using Halfedge_const_handle = typename Aos::Halfedge_const_handle; + using Face_const_handle = typename Aos::Face_const_handle; using Ccb_halfedge_const_circulator = - typename Arr::Ccb_halfedge_const_circulator; + typename Aos::Ccb_halfedge_const_circulator; /// Construct the viewer. /// @param arr the arrangement to view /// @param title the title of the window - Arr_2_viewer_qt(QWidget* parent, const Arr& arr, + Aos_2_viewer_qt(QWidget* parent, const Aos& aos, Color_generator color_generator, - const char* title = "2D Arrangement Basic Viewer", + const char* title = "2D Arrangement on Surface Basic Viewer", bool draw_vertices = false) : - Base(parent, arr, color_generator, title, draw_vertices) + Base(parent, aos, color_generator, title, draw_vertices) {} }; +/*! Draw an arrangement on surface. + */ +template +void draw(const Arrangement_on_surface_2& aos, + const char* title = "2D Arrangement on Surface Basic Viewer", + bool draw_vertices = false) { +#if defined(CGAL_TEST_SUITE) + bool cgal_test_suite=true; +#else + bool cgal_test_suite = qEnvironmentVariableIsSet("CGAL_TEST_SUITE"); +#endif + + if (cgal_test_suite) return; + using Gt = GeometryTraits_2; + using Tt = TopologyTraits; + using Aos = CGAL::Arrangement_on_surface_2; + using Viewer = Aos_2_viewer_qt; + + CGAL::Qt::init_ogl_context(4,3); + + int argc = 1; + const char* argv[2] = {"t2_viewer", nullptr}; + QApplication app(argc, const_cast(argv)); + Default_color_generator color_generator; + Viewer mainwindow(app.activeWindow(), aos, color_generator, title, + draw_vertices); + mainwindow.add_elements(); + mainwindow.show(); + app.exec(); +} + /*! Draw an arrangement. */ template @@ -555,7 +701,7 @@ void draw(const Arrangement_2& arr, if (cgal_test_suite) return; using Gt = GeometryTraits_2; using Arr = CGAL::Arrangement_2; - using Viewer = Arr_2_viewer_qt; + using Viewer = Aos_2_viewer_qt; CGAL::Qt::init_ogl_context(4,3); @@ -589,7 +735,7 @@ void draw(const Arrangement_2& arr, using Color_generator = ColorGenerator; using Gt = GeometryTraits_2; using Arr = CGAL::Arrangement_2; - using Viewer = Arr_2_viewer_qt; + using Viewer = Aos_2_viewer_qt; CGAL::Qt::init_ogl_context(4,3); diff --git a/Arrangement_on_surface_2/test/Arrangement_on_surface_2/data/test_construction/geodesic_arcs_on_sphere/test40.txt b/Arrangement_on_surface_2/test/Arrangement_on_surface_2/data/test_construction/geodesic_arcs_on_sphere/test40.txt new file mode 100644 index 000000000000..adccf6ecad93 --- /dev/null +++ b/Arrangement_on_surface_2/test/Arrangement_on_surface_2/data/test_construction/geodesic_arcs_on_sphere/test40.txt @@ -0,0 +1,17 @@ +5 +0 -1 0 0 0 -1 0 +0 0 -1 0 0 0 -1 +0 0 0 -1 0 1 0 +0 0 1 0 -1 0 0 +0 0 0 -1 -1 0 0 +0 +4 5 3 +0 0 -1 +-1 0 0 +0 -1 0 +0 1 0 +0 0 0 -1 -1 0 0 1 +0 -1 0 0 0 -1 0 1 +0 0 -1 0 0 0 -1 1 +0 0 0 -1 0 1 0 1 +0 0 1 0 -1 0 0 1 diff --git a/Arrangement_on_surface_2/test/Arrangement_on_surface_2/test_construction.geodesic_arcs_on_sphere.cmd b/Arrangement_on_surface_2/test/Arrangement_on_surface_2/test_construction.geodesic_arcs_on_sphere.cmd index 547aa5188647..d995644a0954 100644 --- a/Arrangement_on_surface_2/test/Arrangement_on_surface_2/test_construction.geodesic_arcs_on_sphere.cmd +++ b/Arrangement_on_surface_2/test/Arrangement_on_surface_2/test_construction.geodesic_arcs_on_sphere.cmd @@ -37,3 +37,4 @@ data/test_construction/geodesic_arcs_on_sphere/test36.txt data/test_construction/geodesic_arcs_on_sphere/test37.txt data/test_construction/geodesic_arcs_on_sphere/test38.txt data/test_construction/geodesic_arcs_on_sphere/test39.txt +data/test_construction/geodesic_arcs_on_sphere/test40.txt diff --git a/Installation/CHANGES.md b/Installation/CHANGES.md index 01485063f191..8ee04d38577d 100644 --- a/Installation/CHANGES.md +++ b/Installation/CHANGES.md @@ -35,6 +35,10 @@ Release date: October 2023 - Removed the class templates `Gray_image_mesh_domain_3`, `Implicit_mesh_domain_3`, and `Labeled_image_mesh_domain_3` which are deprecated since CGAL-4.13. +### [2D Arrangements](https://doc.cgal.org/6.0/Manual/packages.html#PkgArrangementOnSurface2) +- Fixed a bug in the zone construction code applied to arrangements of geodesic arcs on a sphere, + when inserting an arc that lies on the identification curve. + ### [Tetrahedral Remeshing](https://doc.cgal.org/6.0/Manual/packages.html#PkgTetrahedralRemeshing) - **Breaking change**: The template parameters of `CGAL::Tetrahedral_remeshing::Remeshing_cell_base_3`