|
|
Boost-Commit : |
Subject: [Boost-commit] svn:boost r66453 - in sandbox/SOC/2010/sweepline: boost/sweepline boost/sweepline/detail libs/sweepline/test
From: sydorchuk.andriy_at_[hidden]
Date: 2010-11-08 09:38:21
Author: asydorchuk
Date: 2010-11-08 09:38:18 EST (Mon, 08 Nov 2010)
New Revision: 66453
URL: http://svn.boost.org/trac/boost/changeset/66453
Log:
Added epsilon robust comparator class.
Changed ppp circle event creation function.
Made simplification step epsilon robust for point input data sets.
Text files modified:
sandbox/SOC/2010/sweepline/boost/sweepline/detail/voronoi_formation.hpp | 351 ++++++++++++++++++++++++---------------
sandbox/SOC/2010/sweepline/boost/sweepline/voronoi_builder.hpp | 18 -
sandbox/SOC/2010/sweepline/boost/sweepline/voronoi_output.hpp | 4
sandbox/SOC/2010/sweepline/libs/sweepline/test/circle_event_creation_test.cpp | 6
sandbox/SOC/2010/sweepline/libs/sweepline/test/event_queue_test.cpp | 2
sandbox/SOC/2010/sweepline/libs/sweepline/test/event_types_test.cpp | 2
sandbox/SOC/2010/sweepline/libs/sweepline/test/sweepline_test.cpp | 9
7 files changed, 235 insertions(+), 157 deletions(-)
Modified: sandbox/SOC/2010/sweepline/boost/sweepline/detail/voronoi_formation.hpp
==============================================================================
--- sandbox/SOC/2010/sweepline/boost/sweepline/detail/voronoi_formation.hpp (original)
+++ sandbox/SOC/2010/sweepline/boost/sweepline/detail/voronoi_formation.hpp 2010-11-08 09:38:18 EST (Mon, 08 Nov 2010)
@@ -18,10 +18,6 @@
if (a >= 0) { res = static_cast<unsigned long long>(a); sign = true; } \
else { res = static_cast<unsigned long long>(-a); sign = false; }
-#define INT_PREDICATE_AVOID_CANCELLATION(val, first_expr, second_expr) \
- if ((val) >= 0) first_expr += (val); \
- else second_expr -= (val);
-
namespace boost {
namespace sweepline {
namespace detail {
@@ -31,25 +27,34 @@
///////////////////////////////////////////////////////////////////////////
template <typename T>
- struct beach_line_node;
+ class beach_line_node;
template <typename T>
- struct beach_line_node_data;
+ class beach_line_node_data;
template <typename BeachLineNode>
struct node_comparer;
+ template <typename T>
+ class epsilon_robust_comparator;
+
enum kOrientation {
RIGHT_ORIENTATION = -1,
COLINEAR = 0,
LEFT_ORIENTATION = 1,
};
+
+ enum kPredicateResult {
+ LESS = -1,
+ UNDEFINED = 0,
+ MORE = 1,
+ };
// Site event type.
// Occurs when sweepline sweeps over one of the initial sites.
// Contains index of a site among the other sorted sites.
template <typename T>
- struct site_event {
+ class site_event {
public:
typedef T coordinate_type;
typedef point_2d<T> Point2D;
@@ -243,7 +248,7 @@
// Circle event contains circle center, lowest x coordinate, event state and
// iterator to the corresponding bisectors.
template <typename T>
- struct circle_event {
+ class circle_event {
public:
typedef T coordinate_type;
typedef point_2d<T> Point2D;
@@ -254,29 +259,24 @@
circle_event() : is_active_(true) {}
- circle_event(coordinate_type c_x, coordinate_type c_y, coordinate_type lower_x) :
- center_(c_x, c_y), lower_x_(lower_x), is_active_(true) {}
-
- circle_event(const Point2D ¢er, coordinate_type lower_x) :
- center_(center), lower_x_(lower_x), is_active_(true) {}
-
- circle_event(const circle_event& c_event) {
- center_ = c_event.center_;
- lower_x_ = c_event.lower_x_;
- bisector_node_ = c_event.bisector_node_;
- is_active_ = c_event.is_active_;
- }
-
- void operator=(const circle_event& c_event) {
- center_ = c_event.center_;
- lower_x_ = c_event.lower_x_;
- bisector_node_ = c_event.bisector_node_;
- is_active_ = c_event.is_active_;
- }
+ circle_event(coordinate_type c_x, coordinate_type c_y,
+ coordinate_type lower_x) :
+ center_x_(c_x),
+ center_y_(c_y),
+ lower_x_(lower_x),
+ is_active_(true) {}
+
+ circle_event(coordinate_type c_x_pos, coordinate_type c_x_neg,
+ coordinate_type c_y_pos, coordinate_type c_y_neg,
+ coordinate_type lower_x_pos, coordinate_type lower_x_neg) :
+ center_x_(c_x_pos, c_x_neg),
+ center_y_(c_y_pos, c_y_neg),
+ lower_x_(lower_x_pos, lower_x_neg),
+ is_active_(true) {}
bool operator==(const circle_event &c_event) const {
- return (center_.y() == c_event.y()) &&
- (lower_x_ == c_event.lower_x_);
+ return (center_y_.compare(c_event.center_y_) == UNDEFINED &&
+ lower_x_.compare(c_event.lower_x_) == UNDEFINED);
}
bool operator!=(const circle_event &c_event) const {
@@ -284,9 +284,10 @@
}
bool operator<(const circle_event &c_event) const {
- if (lower_x_ != c_event.lower_x_)
- return lower_x_ < c_event.lower_x_;
- return center_.y() < c_event.y();
+ kPredicateResult pres = lower_x_.compare(c_event.lower_x_);
+ if (pres != UNDEFINED)
+ return (pres == LESS);
+ return (center_y_.compare(c_event.center_y_) == LESS);
}
bool operator<=(const circle_event &c_event) const {
@@ -305,29 +306,35 @@
// If circle point is less than site point return -1.
// If circle point is equal to site point return 0.
// If circle point is greater than site point return 1.
- int compare(const site_event<coordinate_type> &s_event) const {
- if (s_event.x() != lower_x_) {
- return (lower_x_ < s_event.x()) ? -1 : 1;
- }
- if (s_event.y() != center_.y())
- return (center_.y() < s_event.y()) ? -1 : 1;
- return 0;
+ kPredicateResult compare(const site_event<coordinate_type> &s_event) const {
+ kPredicateResult pres = lower_x_.compare(s_event.x());
+ if (pres != UNDEFINED)
+ return pres;
+ return center_y_.compare(s_event.y());
}
coordinate_type x() const {
- return center_.x();
+ return center_x_.get_dif();
}
coordinate_type y() const {
- return center_.y();
+ return center_y_.get_dif();
+ }
+
+ coordinate_type lower_x() const {
+ return lower_x_.get_dif();
}
- const Point2D &get_center() const {
- return center_;
+ Point2D get_center() const {
+ return make_point_2d(x(), y());
}
- const coordinate_type &get_lower_x() const {
- return lower_x_;
+ const epsilon_robust_comparator<T> &get_erc_x() const {
+ return center_x_;
+ }
+
+ const epsilon_robust_comparator<T> &get_erc_y() const {
+ return center_y_;
}
void set_bisector(beach_line_iterator iterator) {
@@ -347,8 +354,9 @@
}
private:
- Point2D center_;
- coordinate_type lower_x_;
+ epsilon_robust_comparator<T> center_x_;
+ epsilon_robust_comparator<T> center_y_;
+ epsilon_robust_comparator<T> lower_x_;
beach_line_iterator bisector_node_;
bool is_active_;
};
@@ -358,11 +366,6 @@
return circle_event<T>(c_x, c_y, lower_x);
}
- template <typename T>
- circle_event<T> make_circle_event(point_2d<T> center, T lower_x) {
- return circle_event<T>(center, lower_x);
- }
-
///////////////////////////////////////////////////////////////////////////
// VORONOI CIRCLE EVENTS QUEUE ////////////////////////////////////////////
///////////////////////////////////////////////////////////////////////////
@@ -432,15 +435,21 @@
///////////////////////////////////////////////////////////////////////////
// GEOMETRY PREDICATES ////////////////////////////////////////////////////
///////////////////////////////////////////////////////////////////////////
+ inline static void avoid_cancellation(double value, double &left_expr, double &right_expr) {
+ if (value >= 0)
+ left_expr += value;
+ else
+ right_expr -= value;
+ }
- bool abs_equal(double a, double b, double abs_error) {
+ inline static bool abs_equal(double a, double b, double abs_error) {
return fabs(a - b) < abs_error;
}
// If two floating-point numbers in the same format are ordered (x < y), then
// they are ordered the same way when their bits are reinterpreted as
// sign-magnitude integers.
- bool almost_equal(double a, double b, long long maxUlps) {
+ static bool almost_equal(double a, double b, int maxUlps) {
long long ll_a, ll_b;
// Reinterpret double bits as long long.
memcpy(&ll_a, &a, sizeof(double));
@@ -465,9 +474,9 @@
// TODO(asydorchuk): Make templates specification for integer coordinate type,
// as it is actually working with integer input.
template <typename T>
- kOrientation orientation_test(const point_2d<T> &point1,
- const point_2d<T> &point2,
- const point_2d<T> &point3) {
+ static kOrientation orientation_test(const point_2d<T> &point1,
+ const point_2d<T> &point2,
+ const point_2d<T> &point3) {
typedef long long ll;
typedef unsigned long long ull;
ull dif_x1, dif_x2, dif_y1, dif_y2;
@@ -511,7 +520,7 @@
}
template <typename T>
- kOrientation orientation_test(T dif_x1_, T dif_y1_, T dif_x2_, T dif_y2_) {
+ static kOrientation orientation_test(T dif_x1_, T dif_y1_, T dif_x2_, T dif_y2_) {
typedef unsigned long long ull;
ull dif_x1, dif_y1, dif_x2, dif_y2;
bool dif_x1_plus, dif_x2_plus, dif_y1_plus, dif_y2_plus;
@@ -547,14 +556,14 @@
}
template <typename T>
- kOrientation orientation_test(T value) {
+ static kOrientation orientation_test(T value) {
if (value == static_cast<T>(0.0))
return COLINEAR;
return (value < static_cast<T>(0.0)) ? RIGHT_ORIENTATION : LEFT_ORIENTATION;
}
template <typename T>
- T robust_cross_product(T a1_, T b1_, T a2_, T b2_) {
+ static T robust_cross_product(T a1_, T b1_, T a2_, T b2_) {
typedef unsigned long long ull;
ull a1, b1, a2, b2;
bool a1_plus, a2_plus, b1_plus, b2_plus;
@@ -591,10 +600,52 @@
}
}
- enum kPredicateResult {
- LESS = -1,
- UNDEFINED = 0,
- MORE = 1,
+ template <typename T>
+ class epsilon_robust_comparator {
+ public:
+ epsilon_robust_comparator() :
+ positive_sum_(0),
+ negative_sum_(0) {}
+
+ epsilon_robust_comparator(T value) :
+ positive_sum_((value>0)?value:0),
+ negative_sum_((value<0)?-value:0) {}
+
+ epsilon_robust_comparator(T pos, T neg) :
+ positive_sum_(pos),
+ negative_sum_(neg) {}
+
+ T get_dif() const {
+ return positive_sum_ - negative_sum_;
+ }
+
+ T get_positive_sum() const {
+ return positive_sum_;
+ }
+
+ T get_negative_sum() const {
+ return negative_sum_;
+ }
+
+ kPredicateResult compare(T value, int ulp = 0) const {
+ T lhs = positive_sum_ - ((value < 0) ? value : 0);
+ T rhs = negative_sum_ + ((value > 0) ? value : 0);
+ if (almost_equal(lhs, rhs, ulp))
+ return UNDEFINED;
+ return (lhs < rhs) ? LESS : MORE;
+ }
+
+ kPredicateResult compare(const epsilon_robust_comparator &rc, int ulp = 0) const {
+ T lhs = positive_sum_ + rc.get_negative_sum();
+ T rhs = negative_sum_ + rc.get_positive_sum();
+ if (almost_equal(lhs, rhs, ulp))
+ return UNDEFINED;
+ return (lhs < rhs) ? LESS : MORE;
+ }
+
+ private:
+ T positive_sum_;
+ T negative_sum_;
};
// Returns true if horizontal line going through new site intersects
@@ -609,9 +660,9 @@
// at first if x2(y*) > x1(y*) or:
// (x0-x2)*(x0-x1)*(x1-x2) + (x0-x2)*(y*-y1)^2 < (x0-x1)*(y*-y2)^2.
template <typename T>
- kPredicateResult fast_less_predicate(const point_2d<T> &left_point,
- const point_2d<T> &right_point,
- const point_2d<T> &new_point) {
+ static kPredicateResult fast_less_predicate(const point_2d<T> &left_point,
+ const point_2d<T> &right_point,
+ const point_2d<T> &new_point) {
double fast_a1 = static_cast<double>(new_point.x()) - static_cast<double>(left_point.x());
double fast_a2 = static_cast<double>(new_point.x()) - static_cast<double>(right_point.x());
double fast_b1 = static_cast<double>(new_point.y()) - static_cast<double>(left_point.y());
@@ -621,17 +672,16 @@
double fast_right_expr = fast_a2 * fast_b1 * fast_b1;
// Avoid cancellation.
- INT_PREDICATE_AVOID_CANCELLATION(fast_a1 * fast_a2 * fast_c,
- fast_left_expr, fast_right_expr);
+ avoid_cancellation(fast_a1 * fast_a2 * fast_c, fast_left_expr, fast_right_expr);
if (!almost_equal(fast_left_expr, fast_right_expr, 5))
return (fast_left_expr < fast_right_expr) ? LESS : MORE;
return UNDEFINED;
}
template <typename T>
- bool less_predicate(const point_2d<T> &left_point,
- const point_2d<T> &right_point,
- const point_2d<T> &new_point) {
+ static bool less_predicate(const point_2d<T> &left_point,
+ const point_2d<T> &right_point,
+ const point_2d<T> &new_point) {
kPredicateResult fast_res = fast_less_predicate(left_point, right_point, new_point);
if (fast_res != UNDEFINED)
return (fast_res == LESS);
@@ -688,8 +738,8 @@
}
template <typename T>
- kPredicateResult fast_less_predicate(point_2d<T> site_point, site_event<T> segment_site,
- point_2d<T> new_point, bool reverse_order) {
+ static kPredicateResult fast_less_predicate(point_2d<T> site_point, site_event<T> segment_site,
+ point_2d<T> new_point, bool reverse_order) {
typedef long long ll;
typedef unsigned long long ull;
if (orientation_test(segment_site.get_point0(true), segment_site.get_point1(true),
@@ -793,8 +843,9 @@
#ifdef USE_MULTI_PRECISION_LIBRARY
template<typename T>
- bool mpz_less_predicate(point_2d<T> segment_start, point_2d<T> segment_end,
- point_2d<T> site_point, point_2d<T> new_point, bool reverse_order) {
+ static bool mpz_less_predicate(point_2d<T> segment_start, point_2d<T> segment_end,
+ point_2d<T> site_point, point_2d<T> new_point,
+ bool reverse_order) {
mpz_class segment_start_x, segment_start_y, segment_end_x, segment_end_y,
site_point_x, site_point_y, new_point_x, new_point_y;
segment_start_x = static_cast<int>(segment_start.x());
@@ -831,8 +882,8 @@
// (point, segment) and (segment, point) are two distinct points, except
// case of vertical segment.
template <typename T>
- bool less_predicate(point_2d<T> site_point, site_event<T> segment_site,
- point_2d<T> new_point, bool reverse_order) {
+ static bool less_predicate(point_2d<T> site_point, site_event<T> segment_site,
+ point_2d<T> new_point, bool reverse_order) {
kPredicateResult fast_res = fast_less_predicate(
site_point, segment_site, new_point, reverse_order);
if (fast_res != UNDEFINED) {
@@ -858,8 +909,8 @@
double left_expr = (a_sqr * temp + (4.0 * dif_x * mul1) * b_sqr) * temp;
double right_expr = (4.0 * dif_x * dif_x) * ((mul1 * mul1) * b_sqr + (mul2 * mul2) * a_sqr);
double common_expr = (4.0 * dif_x * a) * (b * mul2);
- INT_PREDICATE_AVOID_CANCELLATION(common_expr * temp, right_expr, left_expr);
- INT_PREDICATE_AVOID_CANCELLATION(common_expr * (2.0 * dif_x * mul1), left_expr, right_expr);
+ avoid_cancellation(common_expr * temp, right_expr, left_expr);
+ avoid_cancellation(common_expr * (2.0 * dif_x * mul1), left_expr, right_expr);
if (!almost_equal(left_expr, right_expr, 18)) {
if (!reverse_order)
@@ -876,7 +927,9 @@
}
template <typename T>
- bool less_predicate(site_event<T> left_site, site_event<T> right_site, point_2d<T> new_point) {
+ static bool less_predicate(site_event<T> left_site,
+ site_event<T> right_site,
+ point_2d<T> new_point) {
if (left_site.get_site_index() == right_site.get_site_index()) {
return orientation_test(left_site.get_point0(), left_site.get_point1(),
new_point) == LEFT_ORIENTATION;
@@ -915,9 +968,9 @@
mul1 = 1.0 / (b1 - mul1);
}
}
- INT_PREDICATE_AVOID_CANCELLATION(a1 * mul1 * static_cast<double>(new_point.y()),
+ avoid_cancellation(a1 * mul1 * static_cast<double>(new_point.y()),
intersection_x1, intersection_x2);
- INT_PREDICATE_AVOID_CANCELLATION(-c1 * mul1, intersection_x1, intersection_x2);
+ avoid_cancellation(-c1 * mul1, intersection_x1, intersection_x2);
}
double a2 = static_cast<double>(segment2_end.x()) -
@@ -946,9 +999,9 @@
mul2 = 1.0 / (b2 - mul2);
}
}
- INT_PREDICATE_AVOID_CANCELLATION(a2 * static_cast<double>(new_point.y()) * mul2,
+ avoid_cancellation(a2 * static_cast<double>(new_point.y()) * mul2,
intersection_x2, intersection_x1);
- INT_PREDICATE_AVOID_CANCELLATION(-c2 * mul2, intersection_x2, intersection_x1);
+ avoid_cancellation(-c2 * mul2, intersection_x2, intersection_x1);
}
if (!almost_equal(intersection_x1, intersection_x2, 20)) {
@@ -962,6 +1015,10 @@
///////////////////////////////////////////////////////////////////////////
// CIRCLE EVENTS CREATION /////////////////////////////////////////////////
///////////////////////////////////////////////////////////////////////////
+ template <typename T>
+ inline static T get_sqr_distance(T dif_x, T dif_y) {
+ return dif_x * dif_x + dif_y * dif_y;
+ }
// Create circle event from three point sites.
// TODO (asydorchuk): make precision optimizations.
@@ -970,27 +1027,43 @@
const site_event<T> &site2,
const site_event<T> &site3,
circle_event<T> &c_event) {
- double dif_x1 = static_cast<double>(site1.x()) -
- static_cast<double>(site2.x());
- double dif_x2 = static_cast<double>(site2.x()) -
- static_cast<double>(site3.x());
- double dif_y1 = static_cast<double>(site1.y()) -
- static_cast<double>(site2.y());
- double dif_y2 = static_cast<double>(site2.y()) -
- static_cast<double>(site3.y());
+ double dif_x1 = static_cast<double>(site1.x()) - static_cast<double>(site2.x());
+ double dif_x2 = static_cast<double>(site2.x()) - static_cast<double>(site3.x());
+ double dif_y1 = static_cast<double>(site1.y()) - static_cast<double>(site2.y());
+ double dif_y2 = static_cast<double>(site2.y()) - static_cast<double>(site3.y());
double orientation = robust_cross_product(dif_x1, dif_y1, dif_x2, dif_y2);
if (orientation_test(orientation) != RIGHT_ORIENTATION)
return false;
- orientation *= 2.0;
- double b1 = dif_x1 * (site1.x() + site2.x()) + dif_y1 * (site1.y() + site2.y());
- double b2 = dif_x2 * (site2.x() + site3.x()) + dif_y2 * (site2.y() + site3.y());
-
- // Create new circle event.
- double c_x = (b1*dif_y2 - b2*dif_y1) / orientation;
- double c_y = (b2*dif_x1 - b1*dif_x2) / orientation;
- double radius = sqrt((c_x-site2.x())*(c_x-site2.x()) +
- (c_y-site2.y())*(c_y-site2.y()));
- c_event = make_circle_event<double>(c_x, c_y, c_x + radius);
+ double inv_orientation = 0.5 / orientation;
+ double sum_x1 = static_cast<double>(site1.x()) + static_cast<double>(site2.x());
+ double sum_x2 = static_cast<double>(site2.x()) + static_cast<double>(site3.x());
+ double sum_y1 = static_cast<double>(site1.y()) + static_cast<double>(site2.y());
+ double sum_y2 = static_cast<double>(site2.y()) + static_cast<double>(site3.y());
+ double dif_x3 = static_cast<double>(site1.x()) - static_cast<double>(site3.x());
+ double dif_y3 = static_cast<double>(site1.y()) - static_cast<double>(site3.y());
+ double c_x_pos, c_x_neg, c_y_pos, c_y_neg;
+ c_x_pos = c_x_neg = c_y_pos = c_y_neg = 0;
+ detail::avoid_cancellation(dif_x1 * sum_x1 * dif_y2, c_x_pos, c_x_neg);
+ detail::avoid_cancellation(dif_y1 * sum_y1 * dif_y2, c_x_pos, c_x_neg);
+ detail::avoid_cancellation(dif_x2 * sum_x2 * dif_y1, c_x_neg, c_x_pos);
+ detail::avoid_cancellation(dif_y2 * sum_y2 * dif_y1, c_x_neg, c_x_pos);
+ detail::avoid_cancellation(dif_x2 * sum_x2 * dif_x1, c_y_pos, c_y_neg);
+ detail::avoid_cancellation(dif_y2 * sum_y2 * dif_x1, c_y_pos, c_y_neg);
+ detail::avoid_cancellation(dif_x1 * sum_x1 * dif_x2, c_y_neg, c_y_pos);
+ detail::avoid_cancellation(dif_y1 * sum_y1 * dif_x2, c_y_neg, c_y_pos);
+ if (orientation < 0) {
+ (std::swap)(c_x_pos, c_x_neg);
+ (std::swap)(c_y_pos, c_y_neg);
+ inv_orientation = -inv_orientation;
+ }
+ c_x_pos *= inv_orientation;
+ c_x_neg *= inv_orientation;
+ c_y_pos *= inv_orientation;
+ c_y_neg *= inv_orientation;
+ double radius = sqrt(get_sqr_distance(dif_x1, dif_y1) * get_sqr_distance(dif_x2, dif_y2) *
+ get_sqr_distance(dif_x3, dif_y3)) * fabs(inv_orientation);
+ c_event = circle_event<double>(c_x_pos, c_x_neg, c_y_pos, c_y_neg,
+ c_x_pos + radius, c_x_neg);
return true;
}
@@ -1044,7 +1117,7 @@
(c_y-site2.y())*(c_y-site2.y()));
double lower_x = (site3.is_vertical() && site3.is_inverse()) ?
site3.get_point0().x() : c_x + radius;
- c_event = make_circle_event<double>(c_x, c_y, lower_x);
+ c_event = circle_event<double>(c_x, c_y, lower_x);
return true;
}
@@ -1109,7 +1182,7 @@
double c_x = 0.5 * (static_cast<double>(segm_start1.x() + segm_start2.x())) + a1 * t;
double c_y = 0.5 * (static_cast<double>(segm_start1.y() + segm_start2.y())) + b1 * t;
double radius = 0.5 * fabs(c) / sqrt(a);
- c_event = make_circle_event<double>(c_x, c_y, c_x + radius);
+ c_event = circle_event<double>(c_x, c_y, c_x + radius);
return true;
} else {
double c1 = robust_cross_product(b1, a1, segm_end1.y(), segm_end1.x());
@@ -1130,7 +1203,7 @@
double c_x = intersection_x + vec_x * t;
double c_y = intersection_y + vec_y * t;
double radius = fabs(denom * t);
- c_event = make_circle_event<double>(c_x, c_y, c_x + radius);
+ c_event = circle_event<double>(c_x, c_y, c_x + radius);
return true;
}
}
@@ -1170,7 +1243,7 @@
double c_x = (vec_a1 * vec_c2 - vec_c1 * vec_a2) / det;
double c_y = (vec_b1 * vec_c2 - vec_c1 * vec_b2) / det;
double radius = fabs(-b2 * c_x + a2 * c_y + c2) / sqr_sum2;
- c_event = make_circle_event<double>(c_x, c_y, c_x + radius);
+ c_event = circle_event<double>(c_x, c_y, c_x + radius);
return true;
}
@@ -1188,7 +1261,7 @@
// In general case two arcs will create two different bisectors. That's why
// the order of arcs is important to define unique bisector.
template <typename T>
- struct beach_line_node {
+ class beach_line_node {
public:
typedef T coordinate_type;
typedef point_2d<T> Point2D;
@@ -1311,25 +1384,34 @@
// Represents edge data sturcture (bisector segment), which is
// associated with beach line node key in the binary search tree.
template <typename T>
- struct beach_line_node_data {
- voronoi_edge<T> *edge;
- typename circle_events_queue<T>::circle_event_type_it circle_event_it;
- bool contains_circle_event;
-
+ class beach_line_node_data {
+ public:
explicit beach_line_node_data(voronoi_edge<T> *new_edge) :
- edge(new_edge),
- contains_circle_event(false) {}
+ edge_(new_edge),
+ contains_circle_event_(false) {}
void activate_circle_event(typename circle_events_queue<T>::circle_event_type_it it) {
- circle_event_it = it;
- contains_circle_event = true;
+ circle_event_it_ = it;
+ contains_circle_event_ = true;
}
void deactivate_circle_event() {
- if (contains_circle_event)
- circle_event_it->deactivate();
- contains_circle_event = false;
+ if (contains_circle_event_)
+ circle_event_it_->deactivate();
+ contains_circle_event_ = false;
+ }
+
+ voronoi_edge<T> *get_edge() const {
+ return edge_;
}
+
+ void set_edge(voronoi_edge<T> *new_edge) {
+ edge_ = new_edge;
+ }
+ private:
+ typename circle_events_queue<T>::circle_event_type_it circle_event_it_;
+ voronoi_edge<T> *edge_;
+ bool contains_circle_event_;
};
template <typename BeachLineNode>
@@ -1399,13 +1481,19 @@
typedef point_2d<T> Point2D;
voronoi_edge<coordinate_type> *incident_edge;
+ epsilon_robust_comparator<T> center_x;
+ epsilon_robust_comparator<T> center_y;
Point2D vertex;
int num_incident_edges;
typename std::list< voronoi_vertex<coordinate_type> >::iterator voronoi_record_it;
- voronoi_vertex(const Point2D &point, voronoi_edge<coordinate_type>* edge) :
+ voronoi_vertex(const epsilon_robust_comparator<T> &c_x,
+ const epsilon_robust_comparator<T> &c_y,
+ voronoi_edge<coordinate_type>* edge) :
incident_edge(edge),
- vertex(point),
+ center_x(c_x),
+ center_y(c_y),
+ vertex(c_x.get_dif(), c_y.get_dif()),
num_incident_edges(0) {}
};
@@ -1572,7 +1660,6 @@
const circle_event_type &circle,
edge_type *edge12,
edge_type *edge23) {
- //voronoi_rect_.update(circle.get_center());
// Update counters.
num_vertex_records_++;
num_edges_++;
@@ -1581,7 +1668,8 @@
//voronoi_rect_.update(circle.get_center());
// Add new voronoi vertex with point at center of the circle.
- vertex_records_.push_back(voronoi_vertex_type(circle.get_center(), edge12));
+ vertex_records_.push_back(voronoi_vertex_type(
+ circle.get_erc_x(), circle.get_erc_y(), edge12));
voronoi_vertex_type &new_vertex = vertex_records_.back();
new_vertex.num_incident_edges = 3;
new_vertex.voronoi_record_it = vertex_records_.end();
@@ -1696,14 +1784,10 @@
std::advance(edge_it, 2);
if (edge_it1->start_point && edge_it1->end_point &&
- (almost_equal(edge_it1->start_point->vertex.x(),
- edge_it1->end_point->vertex.x(), 1024) ||
- abs_equal(edge_it1->start_point->vertex.x(),
- edge_it1->end_point->vertex.x(), 1E-6)) &&
- (almost_equal(edge_it1->start_point->vertex.y(),
- edge_it1->end_point->vertex.y(), 1024) ||
- abs_equal(edge_it1->start_point->vertex.y(),
- edge_it1->end_point->vertex.y(), 1E-6))) {
+ edge_it1->start_point->center_x.compare(
+ edge_it1->end_point->center_x, 128) == UNDEFINED &&
+ edge_it1->start_point->center_y.compare(
+ edge_it1->end_point->center_y, 128) == UNDEFINED) {
// Decrease number of cell incident edges.
edge_it1->cell->num_incident_edges--;
edge_it1->twin->cell->num_incident_edges--;
@@ -1994,8 +2078,7 @@
} // boost
} // detail
-#undef INT_PREDICATE_COMPUTE_DIFFERENCE
#undef INT_PREDICATE_CONVERT_65_BIT
-#undef INT_PREDICATE_AVOID_CANCELLATION
+#undef INT_PREDICATE_COMPUTE_DIFFERENCE
#endif
Modified: sandbox/SOC/2010/sweepline/boost/sweepline/voronoi_builder.hpp
==============================================================================
--- sandbox/SOC/2010/sweepline/boost/sweepline/voronoi_builder.hpp (original)
+++ sandbox/SOC/2010/sweepline/boost/sweepline/voronoi_builder.hpp 2010-11-08 09:38:18 EST (Mon, 08 Nov 2010)
@@ -129,7 +129,7 @@
} else if (site_events_iterator_ == site_events_.end()) {
process_circle_event();
} else {
- if (circle_events_.top().compare(*site_events_iterator_) > 0) {
+ if (circle_events_.top().compare(*site_events_iterator_) == detail::MORE) {
process_site_event();
} else {
process_circle_event();
@@ -153,12 +153,6 @@
output_.clip(clipped_output);
}
- // Clip using defined rectangle.
- // TODO(asydorchuk): Define what exactly it means to clip some region of voronoi diagram.
- //void clip(const BRect<coordinate_type> &brect, ClippedOutput &clipped_output) {
- // output_.clip(brect, clipped_output);
- //}
-
protected:
typedef typename std::vector<site_event_type>::const_iterator site_events_iterator_type;
@@ -280,7 +274,7 @@
// Doesn't use special comparison function as it works fine only for
// the site events processing.
- void process_circle_event() {
+ void process_circle_event() {
const circle_event_type &circle_event = circle_events_.top();
// Retrieve the second bisector node that corresponds to the given
@@ -292,11 +286,11 @@
site_event_type site3 = it_first->first.get_right_site();
// Get second bisector;
- edge_type *bisector2 = it_first->second.edge;
+ edge_type *bisector2 = it_first->second.get_edge();
// Get first bisector;
it_first--;
- edge_type *bisector1 = it_first->second.edge;
+ edge_type *bisector1 = it_first->second.get_edge();
// Get the first site that creates given circle event.
site_event_type site1 = it_first->first.get_left_site();
@@ -315,8 +309,8 @@
if (site3.is_segment() && site3.get_point1(true) == site1.get_point0(true)) {
const_cast<Key &>(it_first->first).set_right_site_inverse();
}
- it_first->second.edge = output_.insert_new_edge(site1, site3, circle_event,
- bisector1, bisector2);
+ it_first->second.set_edge(output_.insert_new_edge(site1, site3, circle_event,
+ bisector1, bisector2));
beach_line_.erase(it_last);
it_last = it_first;
Modified: sandbox/SOC/2010/sweepline/boost/sweepline/voronoi_output.hpp
==============================================================================
--- sandbox/SOC/2010/sweepline/boost/sweepline/voronoi_output.hpp (original)
+++ sandbox/SOC/2010/sweepline/boost/sweepline/voronoi_output.hpp 2010-11-08 09:38:18 EST (Mon, 08 Nov 2010)
@@ -13,7 +13,7 @@
namespace boost {
namespace sweepline {
template <typename T>
- struct point_2d {
+ class point_2d {
public:
typedef T coordinate_type;
@@ -83,7 +83,7 @@
///////////////////////////////////////////////////////////////////////////
namespace detail {
template <typename T>
- struct site_event;
+ class site_event;
}
// Bounding rectangle data structure.
Modified: sandbox/SOC/2010/sweepline/libs/sweepline/test/circle_event_creation_test.cpp
==============================================================================
--- sandbox/SOC/2010/sweepline/libs/sweepline/test/circle_event_creation_test.cpp (original)
+++ sandbox/SOC/2010/sweepline/libs/sweepline/test/circle_event_creation_test.cpp 2010-11-08 09:38:18 EST (Mon, 08 Nov 2010)
@@ -15,9 +15,9 @@
#include <boost/test/test_case_template.hpp>
#define CHECK_CIRCLE_EVENT(c_e, c_x, c_y, l_x) \
- BOOST_CHECK_EQUAL(detail::almost_equal(c_e.get_center().x(), static_cast<T>(c_x), 10), true); \
- BOOST_CHECK_EQUAL(detail::almost_equal(c_e.get_center().y(), static_cast<T>(c_y), 10), true); \
- BOOST_CHECK_EQUAL(detail::almost_equal(c_e.get_lower_x(), static_cast<T>(l_x), 10), true)
+ BOOST_CHECK_EQUAL(detail::almost_equal(c_e.x(), static_cast<T>(c_x), 10), true); \
+ BOOST_CHECK_EQUAL(detail::almost_equal(c_e.y(), static_cast<T>(c_y), 10), true); \
+ BOOST_CHECK_EQUAL(detail::almost_equal(c_e.lower_x(), static_cast<T>(l_x), 10), true)
// Test for (point, point, point) case.
BOOST_AUTO_TEST_CASE_TEMPLATE(circle_event_creation_ppp_test1, T, test_types) {
Modified: sandbox/SOC/2010/sweepline/libs/sweepline/test/event_queue_test.cpp
==============================================================================
--- sandbox/SOC/2010/sweepline/libs/sweepline/test/event_queue_test.cpp (original)
+++ sandbox/SOC/2010/sweepline/libs/sweepline/test/event_queue_test.cpp 2010-11-08 09:38:18 EST (Mon, 08 Nov 2010)
@@ -17,7 +17,7 @@
#include <boost/test/test_case_template.hpp>
#define CHECK_TOP_ELEMENT_EQUALITY(TOP, X, Y) \
- BOOST_CHECK_EQUAL(TOP.get_lower_x() == static_cast<T>(X) && \
+ BOOST_CHECK_EQUAL(TOP.lower_x() == static_cast<T>(X) && \
TOP.y() == static_cast<T>(Y), true)
BOOST_AUTO_TEST_CASE_TEMPLATE(event_queue_test1, T, test_types) {
Modified: sandbox/SOC/2010/sweepline/libs/sweepline/test/event_types_test.cpp
==============================================================================
--- sandbox/SOC/2010/sweepline/libs/sweepline/test/event_types_test.cpp (original)
+++ sandbox/SOC/2010/sweepline/libs/sweepline/test/event_types_test.cpp 2010-11-08 09:38:18 EST (Mon, 08 Nov 2010)
@@ -168,7 +168,7 @@
BOOST_CHECK_EQUAL(circle1.x(), static_cast<T>(1));
BOOST_CHECK_EQUAL(circle1.y(), static_cast<T>(2));
- BOOST_CHECK_EQUAL(circle1.get_lower_x(), static_cast<T>(3));
+ BOOST_CHECK_EQUAL(circle1.lower_x(), static_cast<T>(3));
circle2 = make_circle_event<T>(static_cast<T>(1), static_cast<T>(2), static_cast<T>(3));
bool arr1[] = { false, false, true, true, true, false };
Modified: sandbox/SOC/2010/sweepline/libs/sweepline/test/sweepline_test.cpp
==============================================================================
--- sandbox/SOC/2010/sweepline/libs/sweepline/test/sweepline_test.cpp (original)
+++ sandbox/SOC/2010/sweepline/libs/sweepline/test/sweepline_test.cpp 2010-11-08 09:38:18 EST (Mon, 08 Nov 2010)
@@ -651,14 +651,15 @@
voronoi_output_clipped<T> test_output_small, test_output_large;
std::vector< point_2d<T> > points_small, points_large;
int num_points = 27;
- int half_num_points = num_points >> 1;
- int koef = std::numeric_limits<int>::max() / half_num_points;
+ int md = 10;
+ int half_md = md >> 1;
+ int koef = std::numeric_limits<int>::max() / half_md;
for (int i = 0; i < 10000; i++) {
points_small.clear();
points_large.clear();
for (int j = 0; j < num_points; j++) {
- T x_small = static_cast<T>(rand() % num_points - half_num_points);
- T y_small = static_cast<T>(rand() % num_points - half_num_points);
+ T x_small = static_cast<T>(rand() % md - half_md);
+ T y_small = static_cast<T>(rand() % md - half_md);
T x_large = x_small * koef;
T y_large = y_small * koef;
points_small.push_back(make_point_2d(x_small, y_small));
Boost-Commit list run by bdawes at acm.org, david.abrahams at rcn.com, gregod at cs.rpi.edu, cpdaniel at pacbell.net, john at johnmaddock.co.uk