diff --git a/Minkowski_sum_2/include/CGAL/Minkowski_sum_2/Exact_offset_base_2.h b/Minkowski_sum_2/include/CGAL/Minkowski_sum_2/Exact_offset_base_2.h index 0534c42a086..b4a3ff8478b 100644 --- a/Minkowski_sum_2/include/CGAL/Minkowski_sum_2/Exact_offset_base_2.h +++ b/Minkowski_sum_2/include/CGAL/Minkowski_sum_2/Exact_offset_base_2.h @@ -219,7 +219,7 @@ protected: // The equation of the line connecting op1 and op2 is given by: // - // (y1 - y2)*x + (x2 - x1)*x + (r*len - y1*x2 - x1*y2) = 0 + // (y1 - y2)*x + (x2 - x1)*y + (r*len - y1*x2 - x1*y2) = 0 // a = nt_traits.convert (-delta_y); b = nt_traits.convert (delta_x); diff --git a/Minkowski_sum_2/include/CGAL/Minkowski_sum_2/Minkowski_sum_conv_2.h b/Minkowski_sum_2/include/CGAL/Minkowski_sum_2/Minkowski_sum_conv_2.h index 943e402d8dc..54d0900ddf7 100644 --- a/Minkowski_sum_2/include/CGAL/Minkowski_sum_2/Minkowski_sum_conv_2.h +++ b/Minkowski_sum_2/include/CGAL/Minkowski_sum_2/Minkowski_sum_conv_2.h @@ -362,20 +362,6 @@ public: } } - // RWRW: Output: - /* RWRW - DELETE THIS: - std::ofstream ofs ("conv_segments.txt"); - - ofs << conv_segments.size() << std::endl; - - typename Segments_list::const_iterator sit; - - for (sit = conv_segments.begin(); sit != conv_segments.end(); ++sit) - ofs << sit->source() << " " << sit->target() << std::endl; - - ofs.close(); - */ - #ifdef RWRW_STATS std::cout << "|P| = " << n1 << " (" << ref1 diff --git a/Minkowski_sum_2/include/CGAL/Minkowski_sum_2/Minkowski_sum_decomp_2.h b/Minkowski_sum_2/include/CGAL/Minkowski_sum_2/Minkowski_sum_decomp_2.h index 16f3ddfa329..d464a10ee72 100644 --- a/Minkowski_sum_2/include/CGAL/Minkowski_sum_2/Minkowski_sum_decomp_2.h +++ b/Minkowski_sum_2/include/CGAL/Minkowski_sum_2/Minkowski_sum_decomp_2.h @@ -157,21 +157,6 @@ public: } } - // RWRW: Output: - /* RWRW - DELTE THIS: - std::ofstream ofs ("decomp_segments.txt"); - - ofs << boundary_segments.size() << std::endl; - - typename Segments_list::const_iterator sit; - - for (sit = boundary_segments.begin(); - sit != boundary_segments.end(); ++sit) - ofs << sit->source() << " " << sit->target() << std::endl; - - ofs.close(); - */ - #ifdef RWRW_STATS _timer.stop(); diff --git a/Minkowski_sum_2/include/CGAL/Minkowski_sum_2/RSA_Minkowski_sum_2.h b/Minkowski_sum_2/include/CGAL/Minkowski_sum_2/RSA_Minkowski_sum_2.h index 897563d3bbe..5cd69a6f596 100644 --- a/Minkowski_sum_2/include/CGAL/Minkowski_sum_2/RSA_Minkowski_sum_2.h +++ b/Minkowski_sum_2/include/CGAL/Minkowski_sum_2/RSA_Minkowski_sum_2.h @@ -23,6 +23,9 @@ #include #include #include +#include +#include +#include CGAL_BEGIN_NAMESPACE @@ -31,7 +34,7 @@ CGAL_BEGIN_NAMESPACE * a polygon with circular arcs, which represents the rotational swept area * of some linear polygon. */ -template +template class RSA_Minkowski_sum_2 { private: @@ -42,14 +45,16 @@ private: typedef typename Rat_kernel::FT Rational; typedef typename Rat_kernel::Point_2 Rat_point_2; typedef typename Rat_kernel::Segment_2 Rat_segment_2; + typedef typename Rat_kernel::Circle_2 Rat_circle_2; typedef typename Rat_kernel::Vector_2 Rat_vector_2; typedef typename Rat_kernel::Direction_2 Rat_direction_2; - typedef Polygon_2 Rat_polygon_2; + typedef Polygon_2 Rat_polygon_2; // Linear kernel functors: typedef typename Rat_kernel::Equal_2 Equal_2; typedef typename Rat_kernel::Compare_xy_2 Compare_xy_2; typedef typename Rat_kernel::Construct_center_2 Construct_center_2; + typedef typename Rat_kernel::Construct_circle_2 Construct_circle_2; typedef typename Rat_kernel::Compute_squared_radius_2 Compute_sqr_radius_2; typedef typename Rat_kernel::Compute_squared_distance_2 @@ -75,6 +80,8 @@ private: typedef typename Conic_traits_2::Curve_2 Curve_2; typedef typename Conic_traits_2::X_monotone_curve_2 X_monotone_curve_2; + typedef typename Alg_kernel::Compare_xy_2 Alg_compare_xy_2; + typedef Gps_traits_2 Gps_traits_2; public: @@ -84,18 +91,24 @@ public: private: // Polygon-related types: - typedef typename Polygon_2::Vertex_const_circulator Rat_vertex_circulator; + typedef typename Rat_polygon_2::Vertex_const_circulator + Rat_vertex_circulator; typedef typename Circ_polygon_2::Curve_const_iterator Circ_edge_iterator; // Labeled traits: typedef Arr_labeled_traits_2 Labeled_traits_2; typedef typename Labeled_traits_2::X_monotone_curve_2 Labeled_curve_2; + typedef std::list Curves_list; + + typedef Union_of_curve_cycles_2 Union_2; // Data members: Nt_traits nt_traits; Equal_2 f_equal; Construct_center_2 f_center; + Construct_circle_2 f_circle; Compute_sqr_radius_2 f_sqr_radius; Compute_sqr_distance_2 f_sqr_distance; Construct_vector_2 f_vector; @@ -104,6 +117,7 @@ private: Ccw_in_between_2 f_ccw_in_between; Translate_point_2 f_add; Compare_xy_2 f_compare_xy; + Alg_compare_xy_2 f_compare_xy_alg; public: @@ -115,6 +129,7 @@ public: f_equal = ker.equal_2_object(); f_center = ker.construct_center_2_object(); + f_circle = ker.construct_circle_2_object(); f_sqr_radius = ker.compute_squared_radius_2_object(); f_sqr_distance = ker.compute_squared_distance_2_object(); f_vector = ker.construct_vector_2_object(); @@ -123,29 +138,75 @@ public: f_ccw_in_between = ker.counterclockwise_in_between_2_object(); f_add = ker.construct_translated_point_2_object(); f_compare_xy = ker.compare_xy_2_object(); + + // Obtain the algebraic-kernel functors. + Alg_kernel alg_ker; + + f_compare_xy_alg = alg_ker.compare_xy_2_object(); + } + + /*! + * Compute the Minkowski sum of a linear polygon and a polygon that contains + * circular arcs, which represents the rotational swept-area of a linear + * polygon. + * Note that as the input polygons may not be convex, the Minkowski sum may + * not be a simple polygon. The result is therefore represented as + * the outer boundary of the Minkowski sum (which is always a simple polygon) + * and a container of simple polygons, representing the holes inside this + * polygon. + * \param pgn1 The linear polygon. + * \param pgn2 The polygon with circular arcs. + * \param sum_bound Output: A polygon respresenting the outer boundary + * of the Minkowski sum. + * \param sum_holes Output: An output iterator for the holes in the sum, + * represented as simple polygons. + * \pre Both input polygons are simple. + * \pre The value-type of the output iterator is Sum_polygon_2. + * \return A past-the-end iterator for the holes in the sum. + */ + template + OutputIterator operator() (const Rat_polygon_2& pgn1, + const Circ_polygon_2& pgn2, + Sum_polygon_2& sum_bound, + OutputIterator sum_holes) const + { + CGAL_precondition (pgn1.is_simple()); + + // Compute the convolution of the two polygons. + Curves_list conv_cycle; + + _convolution_cycle (pgn1, pgn2, + 1, // Fictitious cycle ID. + conv_cycle); + + // Compute the union of the cycle(s) that represent the Minkowski sum. + Union_2 unite; + + sum_holes = unite (conv_cycle.begin(), conv_cycle.end(), + sum_bound, sum_holes); + + return (sum_holes); } protected: /*! - * Compute the curves that constitute the Minkowski sum of a convex linear - * polygon with another polygon with circular arcs. + * Compute the curves that constitute the convolution cycle(s) of a given + * linear polygon with another polygon with circular arcs. * \param pgn1 The linear polygon. * \param pgn2 The polygon with circular arcs. * \param cycle_id The index of the cycle. - * \param oi An output iterator for the offset curves. - * \pre The value type of the output iterator is Labeled_curve_2. - * \return A past-the-end iterator for the holes container. + * \param cycle Output: The curves that consitute the convolution cycle(s). */ - template - OutputIterator _sum_with_convex (const Rat_polygon_2& pgn1, - const Circ_polygon_2& pgn2, - unsigned int cycle_id, - OutputIterator oi) const + void _convolution_cycle (const Rat_polygon_2& pgn1, + const Circ_polygon_2& pgn2, + unsigned int cycle_id, + Curves_list& cycle) const { // Go over the vertices of pgn2 (at each iteration, the vertex we consider // is the target of prev2 and the source of curr2). - const bool forward1 = (pgn1.orientation() == COUNTERCLOCKWISE); + const bool forward1 = (pgn1.orientation() == + CGAL::COUNTERCLOCKWISE); Rat_vertex_circulator first1 = pgn1.vertices_circulator(); Rat_vertex_circulator curr1, next1; Rat_direction_2 dir_curr1; @@ -156,8 +217,9 @@ protected: Rat_direction_2 dir_prev2, dir_next2; bool equal_dirs; bool shift_edge; - Rat_point_2 ps, pt; unsigned int xcv_index = 0; + X_monotone_curve_2 shifted_xcv; + bool dir_right; prev2 = end2; --prev2; next2 = begin2; @@ -204,15 +266,12 @@ protected: if (shift_edge) { // Shift the current edge in pgn1 by the current vertex in pgn2. - ps = f_add (*curr1, f_vector(CGAL::ORIGIN, p2)); - pt = f_add (*next1, f_vector(CGAL::ORIGIN, p2)); - - res = f_compare_xy (ps, pt); - CGAL_assertion (res != EQUAL); + shifted_xcv = _shift_segment (*curr1, *next1, + p2, + dir_right); - Curve_2 seg (Rat_segment_2 (ps, pt)); - cycle.push_back (Labeled_curve_2 (X_monotone_curve_2 (seg), - X_curve_label ((res == SMALLER), + cycle.push_back (Labeled_curve_2 (shifted_xcv, + X_curve_label (dir_right, cycle_id, 2*xcv_index, false))); @@ -236,6 +295,8 @@ protected: Rat_direction_2 dir_s2, dir_t2; bool is_between_s2, is_between_t2; bool is_between_prev1, is_between_next1; + Alg_point_2 ps, pt; + bool shift_arc, shift_prev1, shift_next1; prev1 = pgn1.vertices_circulator(); if (forward1) @@ -272,18 +333,16 @@ protected: // Shift the current edge in pgn1 by the current vertex in pgn2. s2 = Rat_point_2 (curr2->source().x().alpha(), curr2->source().y().alpha()); - ps = f_add (s2, f_vector(CGAL::ORIGIN, *curr1)); - t2 = Rat_point_2 (curr2->target.x().alpha(), + t2 = Rat_point_2 (curr2->target().x().alpha(), curr2->target().y().alpha()); - pt = f_add (t2, f_vector(CGAL::ORIGIN, *curr1)); - res = f_compare_xy (ps, pt); - CGAL_assertion (res != EQUAL); + shifted_xcv = _shift_segment (s2, t2, + *curr1, + dir_right); - Curve_2 seg (Rat_segment_2 (ps, pt)); - cycle.push_back (Labeled_curve_2 (X_monotone_curve_2 (seg), - X_curve_label ((res == SMALLER), + cycle.push_back (Labeled_curve_2 (shifted_xcv, + X_curve_label (dir_right, cycle_id, 2*xcv_index, false))); @@ -298,7 +357,7 @@ protected: curr2->source().y().alpha()); dir_s2 = _direction (*curr2, s2); - t2 = Rat_point_2 (curr2->target.x().alpha(), + t2 = Rat_point_2 (curr2->target().x().alpha(), curr2->target().y().alpha()); dir_t2 = _direction (*curr2, s2); @@ -313,17 +372,35 @@ protected: f_equal (dir_t2, dir_next1); // Act according to the result. + shift_arc = shift_prev1 = shift_next1 = false; + if (is_between_s2 && is_between_t2) { - // RWRW: (case (b)) + // We should add the entire arc, shifted by curr1, to the cycle. + ps = _convert_point (curr2->source()); + pt = _convert_point (curr2->target()); + shift_arc = true; } else if (is_between_s2) { - // RWRW: (case (c)) + // We should split the arc and shift the first portion. We should + // also shift the segment after curr1 by the split point. + ps = _convert_point (curr2->source()); + pt = _point_on_circle (curr2->supporting_circle(), + *curr1, *next1); + shift_arc = true; + shift_next1 = true; } else if (is_between_t2) { - // RWRW: (case (d)) + // We should split the arc and shift the second portion. We should + // also shift the segment before curr1 by the split point. + ps = _point_on_circle (curr2->supporting_circle(), + *prev1, *curr1); + pt = _convert_point (curr2->target()); + + shift_arc = true; + shift_prev1 = true; } else { @@ -338,13 +415,66 @@ protected: if (is_between_prev1 && is_between_next1) { - // RWRW: (case (e)) + ps = _point_on_circle (curr2->supporting_circle(), + *prev1, *curr1); + pt = _point_on_circle (curr2->supporting_circle(), + *curr1, *next1); + + shift_arc = true; + shift_prev1 = true; + shift_next1 = true; } else { CGAL_assertion (! is_between_prev1 && ! is_between_next1); } } + + // If necessary, shifted the circular arc by curr1. + if (shift_arc) + { + shifted_xcv = _shift_circular_arc (curr2->supporting_circle(), + ps, pt, + *curr1, + dir_right); + + cycle.push_back (Labeled_curve_2 (shifted_xcv, + X_curve_label (dir_right, + cycle_id, + 2*xcv_index, + false))); + xcv_index++; + } + + // If necessary, shift the polygon edge before curr1 by ps. + if (shift_prev1) + { + shifted_xcv = _shift_segment (*prev1, *curr1, + ps, + dir_right); + + cycle.push_back (Labeled_curve_2 (shifted_xcv, + X_curve_label (dir_right, + cycle_id, + 2*xcv_index, + false))); + xcv_index++; + } + + // If necessary, shift the polygon edge after curr1 by pt. + if (shift_next1) + { + shifted_xcv = _shift_segment (*curr1, *next1, + pt, + dir_right); + + cycle.push_back (Labeled_curve_2 (shifted_xcv, + X_curve_label (dir_right, + cycle_id, + 2*xcv_index, + false))); + xcv_index++; + } } } @@ -354,7 +484,7 @@ protected: } while (curr1 != first1); - return (oi); + return; } /*! @@ -363,8 +493,8 @@ protected: * \param p The point (either its source or its target). * \return The direction. */ - Direction_2 _direction (const Circ_segment_2& cv, - const Rat_point_2& p) const + Rat_direction_2 _direction (const Circ_segment_2& cv, + const Rat_point_2& p) const { // Act according to the curve type. if (cv.is_linear()) @@ -381,10 +511,34 @@ protected: Rat_circle_2 circ = cv.supporting_circle(); Rat_vector_2 vec = f_vector (f_center (circ), p); - return (f_direction (f_perp_vec (vec, COUNTERCLOCKWISE))); + return (f_direction (f_perp_vector (vec, CGAL::COUNTERCLOCKWISE))); } } + /*! + * Convert a point with one-root coordinates to an algebraic point. + */ + Alg_point_2 _convert_point (const Circ_point_2& p) const + { + Algebraic x = nt_traits.convert (p.x().alpha()); + + if (! p.x().is_rational()) + { + x += nt_traits.convert (p.x().beta()) * + nt_traits.sqrt (nt_traits.convert (p.x().gamma())); + } + + Algebraic y = nt_traits.convert (p.y().alpha()); + + if (! p.y().is_rational()) + { + y += nt_traits.convert (p.y().beta()) * + nt_traits.sqrt (nt_traits.convert (p.y().gamma())); + } + + return (Alg_point_2 (x, y)); + } + /*! * Compute a point on the given circle, where the tangent to the circle * has the same direction as the given vector. @@ -429,12 +583,21 @@ protected: * \param ps The source point of the segment. * \param pt The target point of the segment. * \param q The shift point. - * \return The shifted segemnt, represented as an x-monotone conic curve. + * \param dir_right Output: Is the segment directed right. + * \return The shifted segment, represented as an x-monotone conic curve. */ X_monotone_curve_2 _shift_segment (const Rat_point_2& ps, const Rat_point_2& pt, - const Rat_point_2& q) const + const Rat_point_2& q, + bool& dir_right) const { + // Determine if the segment is directed left or right. + Comparison_result res = f_compare_xy (ps, pt); + + CGAL_assertion (res != CGAL::EQUAL); + dir_right = (res == CGAL::SMALLER); + + // Construct the shifted segment. Rat_point_2 shift_ps = f_add (ps, f_vector(CGAL::ORIGIN, q)); Rat_point_2 shift_pt = f_add (pt, f_vector(CGAL::ORIGIN, q)); Curve_2 seg (Rat_segment_2 (shift_ps, shift_pt)); @@ -447,27 +610,74 @@ protected: * \param ps The source point of the segment. * \param pt The target point of the segment. * \param q The shift point. - * \return The shifted segemnt, represented as an x-monotone conic curve. + * \param dir_right Output: Is the segment directed right. + * \return The shifted segment, represented as an x-monotone conic curve. */ X_monotone_curve_2 _shift_segment (const Rat_point_2& ps, const Rat_point_2& pt, - const Alg_point_2& q) const + const Alg_point_2& q, + bool& dir_right) const { + // Determine if the segment is directed left or right. + Comparison_result res = f_compare_xy (ps, pt); + + CGAL_assertion (res != CGAL::EQUAL); + dir_right = (res == CGAL::SMALLER); + + // Construct the shifted segment. Alg_point_2 shift_ps = Alg_point_2 (nt_traits.convert (ps.x()) + q.x(), nt_traits.convert (ps.y()) + q.y()); Alg_point_2 shift_pt = Alg_point_2 (nt_traits.convert (pt.x()) + q.x(), nt_traits.convert (pt.y()) + q.y()); - return (X_monotone_curve_2 (shift_ps, shift_pt)); + // The supprting line a*x + b*y + c = 0 connecting the two shifted points + // is given by: + const Algebraic a = nt_traits.convert (ps.y() - pt.y()); + const Algebraic b = nt_traits.convert (pt.x() - ps.x()); + const Algebraic c = nt_traits.convert (ps.x()*pt.y() - ps.y()*pt.x()) - + (a*q.x() + b*q.y()); + + return (X_monotone_curve_2 (a, b, c, + shift_ps, shift_pt)); } /*! + * Shift a circular arc, given by a rational supporting circle and two + * endpoints, by a point with rational coordinates. + * \param circ The supporting circle. + * \param ps The source point of the arc. + * \param pt The target point of the arc. + * \param q The shift point. + * \param dir_right Output: Is the segment directed right. + * \return The shifted arc. */ - X_monotone_curve_2 _shift_arc (const Rat_circle_2& pc, - const Alg_point_2& ps, - const Alg_point_2& pt, - const Rat_point_2& q) const + X_monotone_curve_2 _shift_circular_arc (const Rat_circle_2& circ, + const Alg_point_2& ps, + const Alg_point_2& pt, + const Rat_point_2& q, + bool& dir_right) const { + // Determine if the segment is directed left or right. + Comparison_result res = f_compare_xy_alg (ps, pt); + + CGAL_assertion (res != CGAL::EQUAL); + dir_right = (res == CGAL::SMALLER); + + // Shift the supporting circle. + Rat_point_2 shift_pc = f_add (f_center (circ), f_vector(CGAL::ORIGIN, q)); + Rat_circle_2 shift_circ = f_circle (shift_pc, + f_sqr_radius (circ)); + + // Shift the two endpoint and construct the circular arc. + Alg_point_2 shift_ps = Alg_point_2 (ps.x() + nt_traits.convert (q.x()), + ps.y() + nt_traits.convert (q.y())); + Alg_point_2 shift_pt = Alg_point_2 (pt.x() + nt_traits.convert (q.x()), + pt.y() + nt_traits.convert (q.y())); + Curve_2 arc (shift_circ, + CGAL::COUNTERCLOCKWISE, + shift_ps, shift_pt); + + return (arc); } }; diff --git a/Minkowski_sum_2/include/CGAL/rsa_minkowski_sum_2.h b/Minkowski_sum_2/include/CGAL/rsa_minkowski_sum_2.h index 9f26fae42f6..883dfdeae95 100644 --- a/Minkowski_sum_2/include/CGAL/rsa_minkowski_sum_2.h +++ b/Minkowski_sum_2/include/CGAL/rsa_minkowski_sum_2.h @@ -20,9 +20,23 @@ #ifndef CGAL_RSA_MINKOWSKI_SUM_H #define CGAL_RSA_MINKOWSKI_SUM_H +#include + CGAL_BEGIN_NAMESPACE /*! + * Compute the Minkowski sum of a linear polygon and a polygon that contains + * circular arcs, which represents the rotational swept-area of a linear + * polygon. + * Note that as the input polygons may not be convex, the Minkowski sum may + * not be a simple polygon. The result is therefore represented as + * polygon with holes whose arcs are conic curves (which are either line + * segments with irrational coefficients, or circular arcs with rational + * coefficients). + * \param pgn1 The linear polygon. + * \param pgn2 The polygon with circular arcs. + * \pre Both input polygons are simple. + * \return A polygon with holes that reprsents the sum. */ template typename Gps_traits_2::Polygon_with_holes_2 @@ -32,31 +46,18 @@ rsa_minkowski_sum_2 const typename Gps_circle_segment_traits_2::Polygon_2& pgn2) { - typedef ConicTraits Conic_traits_2; - typedef typename Conic_traits_2::Rat_kernel Rat_kernel; - typedef typename Rat_kernel::FT Rational; - typedef typename Rat_kernel::Point_2 Point_2; - typedef typename Rat_kernel:: - typedef Polygon_2 Polygon_2; - - typedef typename Conic_traits_2::Alg_kernel Alg_kernel; - typedef typename Alg_kernel::FT Algebraic; - typedef typename Alg_kernel::Point_2 Alg_point_2; - - typedef Exact_offset_base_2 Base; - typedef Offset_by_convolution_2 Exact_offset_2; - typedef typename Exact_offset_2::Offset_polygon_2 Offset_polygon_2; + typedef RSA_Minkowski_sum_2 RSA_sum_2; + typedef typename RSA_sum_2::Sum_polygon_2 Sum_polygon_2; - Base base; - Exact_offset_2 exact_offset (base); - Offset_polygon_2 offset_bound; - std::list offset_holes; + RSA_sum_2 rsa_sum; + Sum_polygon_2 sum_bound; + std::list sum_holes; - exact_offset (pgn, r, - offset_bound, std::back_inserter(offset_holes)); + rsa_sum (pgn1, pgn2, + sum_bound, std::back_inserter(sum_holes)); return (typename Gps_traits_2::Polygon_with_holes_2 - (offset_bound, offset_holes.begin(), offset_holes.end())); + (sum_bound, sum_holes.begin(), sum_holes.end())); } CGAL_END_NAMESPACE