songsenand
/
cgal
mirror of https://github.com/CGAL/cgal
1
0
Fork 0
Code Issues Packages Projects Releases Wiki Activity

Use a symbolic perturbation so that in the grid case all diagonals have the same orientation

This commit is contained in:
Andreas Fabri 2010-02-25 15:29:19 +00:00
parent 4420a4ef71
commit 2ae72afd4d
6 changed files with 151 additions and 26 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

10
Triangulation_2/include/CGAL/Delaunay_triangulation_2.h
View File

@ -275,7 +275,7 @@ protected:
int mi = this->_tds.mirror_index(f, i);
Vertex_handle vm = this->_tds.mirror_vertex(f, i);
if(this->is_infinite(vm)) continue;
if(this->side_of_oriented_circle(f, vm->point()) == ON_POSITIVE_SIDE) {
if(this->side_of_oriented_circle(f, vm->point(),true) == ON_POSITIVE_SIDE) {
this->_tds.flip(f, i);
edges.push_back(Edge(f, i));
edges.push_back(Edge(f, this->cw(i)));
@ -296,7 +296,7 @@ test_conflict(const Point &p, Face_handle fh) const
// return true if P is inside the circumcircle of fh
// if fh is infinite, return true when p is in the positive
// halfspace or on the boundary and in the finite edge of fh
Oriented_side os = side_of_oriented_circle(fh,p);
Oriented_side os = side_of_oriented_circle(fh,p,true);
if (os == ON_POSITIVE_SIDE) return true;
if (os == ON_ORIENTED_BOUNDARY && is_infinite(fh)) {
@ -339,7 +339,7 @@ is_valid(bool verbose, int level) const
for(int i=0; i<3; i++) {
if ( ! is_infinite( this->mirror_vertex(it,i))) {
result = result && ON_POSITIVE_SIDE !=
side_of_oriented_circle( it, this->mirror_vertex(it,i)->point());
side_of_oriented_circle( it, this->mirror_vertex(it,i)->point(), false);
}
CGAL_triangulation_assertion( result );
}
@ -420,7 +420,7 @@ look_nearest_neighbor(const Point& p,
Vertex_handle& nn) const
{
Face_handle ni=f->neighbor(i);
if ( ON_POSITIVE_SIDE != side_of_oriented_circle(ni,p) ) return;
if ( ON_POSITIVE_SIDE != side_of_oriented_circle(ni,p,true) ) return;
typename Geom_traits::Compare_distance_2
compare_distance = this->geom_traits().compare_distance_2_object();
@ -576,7 +576,7 @@ propagating_flip(Face_handle& f,int i)
Face_handle n = f->neighbor(i);
if ( ON_POSITIVE_SIDE !=
side_of_oriented_circle(n, f->vertex(i)->point()) ) {
side_of_oriented_circle(n, f->vertex(i)->point(), true) ) {
return;
}
flip(f, i);

64
Triangulation_2/include/CGAL/Regular_triangulation_2.h
View File

@ -25,6 +25,8 @@
#include <CGAL/Regular_triangulation_vertex_base_2.h>
#include <CGAL/utility.h>
#include <boost/bind.hpp>
CGAL_BEGIN_NAMESPACE
template < typename K_ >
@ -192,14 +194,14 @@ public:
Oriented_side power_test(const Weighted_point &p,
const Weighted_point &q,
const Weighted_point &r,
const Weighted_point &s) const;
const Weighted_point &s, bool perturb) const;
Oriented_side power_test(const Weighted_point &p,
const Weighted_point &q,
const Weighted_point &r) const;
Oriented_side power_test(const Weighted_point &p,
const Weighted_point &r) const;
Oriented_side power_test(const Face_handle &f,
const Weighted_point &p) const;
const Weighted_point &p, bool perturb=false) const;
Oriented_side power_test(const Face_handle& f, int i,
const Weighted_point &p) const;
@ -615,7 +617,7 @@ operator=(const Self &tr)
template < class Gt, class Tds >
Oriented_side
Regular_triangulation_2<Gt,Tds>::
power_test(const Face_handle &f, const Weighted_point &p) const
power_test(const Face_handle &f, const Weighted_point &p, bool perturb) const
{
// p is supposed to be a finite point
// if f is a finite face,
@ -627,7 +629,7 @@ power_test(const Face_handle &f, const Weighted_point &p) const
if ( ! f->has_vertex(infinite_vertex(), i) )
return power_test(f->vertex(0)->point(),
f->vertex(1)->point(),
f->vertex(2)->point(),p);
f->vertex(2)->point(),p, perturb);
Orientation o = orientation(f->vertex(ccw(i))->point(),
f->vertex( cw(i))->point(),
@ -662,12 +664,48 @@ template < class Gt, class Tds >
inline
Oriented_side
Regular_triangulation_2<Gt,Tds>::
power_test(const Weighted_point &p,
const Weighted_point &q,
const Weighted_point &r,
const Weighted_point &s) const
power_test(const Weighted_point &p0,
const Weighted_point &p1,
const Weighted_point &p2,
const Weighted_point &p,
bool perturb) const
{
return geom_traits().power_test_2_object()(p,q,r,s);
CGAL_triangulation_precondition( orientation(p0, p1, p2) == POSITIVE );
using namespace boost;
Oriented_side os = geom_traits().power_test_2_object()(p0, p1, p2, p);
if ( (os != ON_ORIENTED_BOUNDARY) || (! perturb))
return os;
// We are now in a degenerate case => we do a symbolic perturbation.
// We sort the points lexicographically.
const Weighted_point * points[4] = {&p0, &p1, &p2, &p};
std::sort(points, points + 4,
bind(geom_traits().compare_xy_2_object(),
bind(Dereference<Weighted_point>(), _1),
bind(Dereference<Weighted_point>(), _2)) == SMALLER);
// We successively look whether the leading monomial, then 2nd monomial
// of the determinant has non null coefficient.
// 2 iterations are enough (cf paper)
for (int i=3; i>1; --i) {
if (points[i] == &p)
return ON_NEGATIVE_SIDE; // since p0 p1 p2 are non collinear
// and positively oriented
Orientation o;
if (points[i] == &p2 && (o = orientation(p0,p1,p)) != COLLINEAR )
return o;
if (points[i] == &p1 && (o = orientation(p0,p,p2)) != COLLINEAR )
return o;
if (points[i] == &p0 && (o = orientation(p,p1,p2)) != COLLINEAR )
return o;
}
CGAL_triangulation_assertion(false);
return ON_NEGATIVE_SIDE;
}
template < class Gt, class Tds >
@ -1090,7 +1128,7 @@ insert(const Weighted_point &p, Locate_type lt, Face_handle loc, int li)
{
CGAL_precondition (dimension() >= 1);
Oriented_side os = dimension() == 1 ? power_test (loc, li, p) :
power_test (loc, p);
power_test (loc, p, true);
if (os < 0) {
if (is_infinite (loc)) loc = loc->neighbor (li);
@ -1103,7 +1141,7 @@ insert(const Weighted_point &p, Locate_type lt, Face_handle loc, int li)
case Base::FACE:
{
CGAL_precondition (dimension() >= 2);
if (power_test (loc, p) < 0) {
if (power_test (loc, p, true) < 0) {
return hide_new_vertex(loc,p);
}
v = insert_in_face (p, loc);
@ -1514,7 +1552,7 @@ fill_hole_regular(std::list<Edge> & first_hole)
p2=p;
cut_after=hit;
}
else if (power_test(p0,p1,p2,p) ==
else if (power_test(p0,p1,p2,p,true) ==
ON_POSITIVE_SIDE)
{
v2=vv;
@ -1824,7 +1862,7 @@ stack_flip(Vertex_handle v, Faces_around_stack &faces_around)
// now dimension() == 2
//test the regularity of edge (f,i)
//if( power_test(n, v->point()) == ON_NEGATIVE_SIDE)
if( power_test(n, v->point()) != ON_POSITIVE_SIDE)
if( power_test(n, v->point(), true) != ON_POSITIVE_SIDE)
return;
if(is_infinite(f,i))

91
Triangulation_2/include/CGAL/Triangulation_2.h
View File

@ -91,6 +91,20 @@ public:
typedef typename Tds::Vertex_iterator All_vertices_iterator;
class Perturbation_order {
const Self *t;
public:
Perturbation_order(const Self *tr)
: t(tr) {}
bool operator()(const Point *p, const Point *q) const {
return t->compare_xy(*p, *q) == SMALLER;
}
};
friend class Perturbation_order;
// This class is used to generate the Finite_*_iterators.
class Infinite_tester
{
@ -363,14 +377,21 @@ public:
Oriented_side
oriented_side(Face_handle f, const Point &p) const;
Oriented_side
side_of_oriented_circle(const Point &p0, const Point &p1, const Point &p2,
const Point &p, bool perturb) const;
Oriented_side
side_of_oriented_circle(Face_handle f, const Point & p) const;
side_of_oriented_circle(Face_handle f, const Point & p, bool perturb = false) const;
bool
collinear_between(const Point& p, const Point& q, const Point& r)
const;
Comparison_result compare_x(const Point& p, const Point& q) const;
Comparison_result compare_xy(const Point& p, const Point& q) const;
Comparison_result compare_y(const Point& p, const Point& q) const;
bool xy_equal(const Point& p, const Point& q) const;
Orientation orientation(const Point& p,
@ -1459,7 +1480,8 @@ fill_hole_delaunay(std::list<Edge> & first_hole)
if (orientation(p0,p1,p) == COUNTERCLOCKWISE) {
if (is_infinite(v2)) { v2=vv; v3=vv; cut_after=hit;}
else{
if (in_circle(p0,p1,v3->point(),p) == ON_POSITIVE_SIDE){
//
if (this->side_of_oriented_circle(p0,p1,v3->point(),p,true) == ON_POSITIVE_SIDE){
v2=vv; v3=vv; cut_after=hit;}
}
}
@ -2432,17 +2454,67 @@ oriented_side(Face_handle f, const Point &p) const
p);
}
template <class Gt, class Tds >
Oriented_side
Triangulation_2<Gt, Tds>::
side_of_oriented_circle(const Point &p0, const Point &p1, const Point &p2,
const Point &p, bool perturb) const
{
CGAL_triangulation_precondition( orientation(p0, p1, p2) == POSITIVE );
Gt::Side_of_oriented_circle_2 pred = geom_traits().side_of_oriented_circle_2_object();
Oriented_side os =
pred(p0, p1, p2, p);
if ((os != ON_ORIENTED_BOUNDARY) || (! perturb))
return os;
// We are now in a degenerate case => we do a symbolic perturbation.
// We sort the points lexicographically.
const Point * points[4] = {&p0, &p1, &p2, &p};
std::sort(points, points+4, Perturbation_order(this) );
// We successively look whether the leading monomial, then 2nd monomial
// of the determinant has non null coefficient.
// 2 iterations are enough (cf paper)
for (int i=3; i>0; --i) {
if (points[i] == &p)
return ON_NEGATIVE_SIDE; // since p0 p1 p2 are non collinear
// and positively oriented
Orientation o;
if (points[i] == &p2 && (o = orientation(p0,p1,p)) != COLLINEAR )
return Oriented_side(o);
if (points[i] == &p1 && (o = orientation(p0,p,p2)) != COLLINEAR )
return Oriented_side(o);
if (points[i] == &p0 && (o = orientation(p,p1,p2)) != COLLINEAR )
return Oriented_side(o);
}
CGAL_triangulation_assertion(false);
return ON_NEGATIVE_SIDE;
}
template < class Gt, class Tds >
Oriented_side
Triangulation_2<Gt,Tds>::
side_of_oriented_circle(Face_handle f, const Point & p) const
side_of_oriented_circle(Face_handle f, const Point & p, bool perturb) const
{
if ( ! is_infinite(f) ) {
/*
typename Gt::Side_of_oriented_circle_2
in_circle = geom_traits().side_of_oriented_circle_2_object();
return in_circle(f->vertex(0)->point(),
f->vertex(1)->point(),
f->vertex(2)->point(),p);
*/
return this->side_of_oriented_circle(f->vertex(0)->point(),
f->vertex(1)->point(),
f->vertex(2)->point(),p, perturb);
}
int i = f->index(infinite_vertex());
@ -2488,6 +2560,19 @@ compare_x(const Point& p, const Point& q) const
return geom_traits().compare_x_2_object()(p,q);
}
template <class Gt, class Tds >
inline
Comparison_result
Triangulation_2<Gt, Tds>::
compare_xy(const Point& p, const Point& q) const
{
Comparison_result res = geom_traits().compare_x_2_object()(p,q);
if(res == EQUAL){
return geom_traits().compare_y_2_object()(p,q);
}
return res;
}
template <class Gt, class Tds >
inline
Comparison_result

2
Triangulation_2/test/Triangulation_2/include/CGAL/_test_cls_regular_triangulation_2.h
View File

@ -555,7 +555,7 @@ _test_cls_regular_triangulation_2( const Triangulation & )
assert(fc==fc2);
n=0;
do {fc2++ ; n = n+1;} while (fc2 != fc);
assert(n==3);
assert(n==4);
/*****************************/
/******** Miscellaneaous *****/

3
Triangulation_2/test/Triangulation_2/include/CGAL/_test_cls_triangulation_2.h
View File

@ -311,7 +311,8 @@ _test_cls_triangulation_2( const Triangul & )
assert(!T2_8.is_infinite(ff));
Face_handle f2 = ff->neighbor(li);
assert(!T2_8.is_infinite(f2));
T2_8.flip(ff,0);
int fli = ff->index(f2);
T2_8.flip(ff,fli);
assert( T2_8.is_valid() );
//make_hole star_hole

7
Triangulation_2/test/Triangulation_2/include/CGAL/_test_cls_triangulation_short_2.h
View File

@ -234,7 +234,8 @@ _test_cls_triangulation_short_2( const Triangul &)
assert(!T2_8.is_infinite(ff));
Face_handle f2 = ff->neighbor(li);
assert(!T2_8.is_infinite(f2));
T2_8.flip(ff,0);
int fli = ff->index(f2);
T2_8.flip(ff,fli);
assert( T2_8.is_valid() );
@ -412,7 +413,7 @@ _test_cls_triangulation_short_2( const Triangul &)
do {fc2++ ; n = n+1;} while (fc2 != fc);
assert(T2_8.number_of_vertices()>=2);
assert(T2_8.is_valid());
fc= T2_8.line_walk(Point(5,4,10),Point(5,5));
fc= T2_8.line_walk(Point(1,4,10),Point(20,5,10));
fc2=fc;
n=0;
assert(fc==fc2);
@ -431,7 +432,7 @@ _test_cls_triangulation_short_2( const Triangul &)
assert(fc==fc2);
n=0;
do {fc2++ ; n = n+1;} while (fc2 != fc);
assert(n==3);
assert(n==4 || n==3);
/*****************************/
/******** Miscellaneaous *****/