Finished the rsa_minkowski_sum2() function.

(Still needs debugging ...)
This commit is contained in:
Ron Wein 2006-04-17 13:13:14 +00:00
parent 6f2452b1e5
commit ea0c8b2927
5 changed files with 282 additions and 100 deletions

View File

@ -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);

View File

@ -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

View File

@ -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();

View File

@ -23,6 +23,9 @@
#include <CGAL/Polygon_2.h>
#include <CGAL/Gps_circle_segment_traits_2.h>
#include <CGAL/Gps_traits_2.h>
#include <CGAL/Minkowski_sum_2/Arr_labeled_traits_2.h>
#include <CGAL/Minkowski_sum_2/Union_of_curve_cycles_2.h>
#include <list>
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 <class ConicTraits_, class Container_, class DecompStrategy_>
template <class ConicTraits_, class Container_>
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_kernel, Container> Rat_polygon_2;
typedef Polygon_2<Rat_kernel, Container_> 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<Conic_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<Conic_traits_2> Labeled_traits_2;
typedef typename Labeled_traits_2::X_monotone_curve_2 Labeled_curve_2;
typedef std::list<Labeled_curve_2> Curves_list;
typedef Union_of_curve_cycles_2<Labeled_traits_2,
Sum_polygon_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 <class OutputIterator>
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 <class OutputIterator>
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);
}
};

View File

@ -20,9 +20,23 @@
#ifndef CGAL_RSA_MINKOWSKI_SUM_H
#define CGAL_RSA_MINKOWSKI_SUM_H
#include <CGAL/Minkowski_sum_2/RSA_Minkowski_sum_2.h>
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 <class ConicTraits, class Container>
typename Gps_traits_2<ConicTraits>::Polygon_with_holes_2
@ -32,31 +46,18 @@ rsa_minkowski_sum_2
const typename Gps_circle_segment_traits_2<typename
ConicTraits::Rat_kernel>::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<Rat_kernel, Container> 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<ConicTraits, Container> Base;
typedef Offset_by_convolution_2<Base> Exact_offset_2;
typedef typename Exact_offset_2::Offset_polygon_2 Offset_polygon_2;
typedef RSA_Minkowski_sum_2<ConicTraits, Container> 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_polygon_2> offset_holes;
RSA_sum_2 rsa_sum;
Sum_polygon_2 sum_bound;
std::list<Sum_polygon_2> 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<ConicTraits>::Polygon_with_holes_2
(offset_bound, offset_holes.begin(), offset_holes.end()));
(sum_bound, sum_holes.begin(), sum_holes.end()));
}
CGAL_END_NAMESPACE