diff --git a/Interpolation/include/CGAL/natural_neighbor_coordinates_2.h b/Interpolation/include/CGAL/natural_neighbor_coordinates_2.h index e625d0cb686..d5480fb0379 100644 --- a/Interpolation/include/CGAL/natural_neighbor_coordinates_2.h +++ b/Interpolation/include/CGAL/natural_neighbor_coordinates_2.h @@ -155,55 +155,7 @@ natural_neighbor_coordinates_vertex_2(const Dt& dt, std::list hole; - dt.get_boundary_of_conflicts(p, std::back_inserter(hole), fh); - - // The symbolic perturbation makes that the conflict zone - // is too big for the interpolation. As a consequence - // we would have points with lambda = 0 in the output - // Therefore we filter them out again - typedef typename std::list::iterator Edge_iterator; - Edge_iterator it = hole.begin(), nit = it; - ++nit; - - do { - Point_2 p0 = it->first->vertex(dt.cw(it->second))->point(); - Point_2 p1 = it->first->vertex(dt.ccw(it->second))->point(); - Point_2 p2 = nit->first->vertex(dt.ccw(nit->second))->point(); - if(dt.geom_traits().side_of_oriented_circle_2_object()(p, p0, p1, p2) - == ON_ORIENTED_BOUNDARY){ - Face_handle nh = it->first->neighbor(it->second); - int ni = nh->index(it->first->vertex(dt.ccw(it->second))); - Edge_iterator eit = hole.insert(it, Edge(nh, ni)); - hole.erase(it); - hole.erase(nit); - it = eit; - nit = it; - } else { - it = nit; - } - ++nit; - }while (nit != hole.end()); - - // now a similar test for the last and first edge - - nit = hole.begin(); - do { - Point_2 p0 = it->first->vertex(dt.cw(it->second))->point(); - Point_2 p1 = it->first->vertex(dt.ccw(it->second))->point(); - Point_2 p2 = nit->first->vertex(dt.ccw(nit->second))->point(); - if(dt.geom_traits().side_of_oriented_circle_2_object()(p, p0, p1, p2) - == ON_ORIENTED_BOUNDARY){ - Face_handle nh = it->first->neighbor(it->second); - int ni = nh->index(it->first->vertex(dt.ccw(it->second))); - Edge_iterator eit = hole.insert(it, Edge(nh, ni)); - hole.erase(it); - hole.erase(nit); - it = eit; - nit = hole.begin(); - } else { - ++it; - } - }while (it != hole.end()); + dt.get_boundary_of_conflicts(p, std::back_inserter(hole), fh, false); return natural_neighbor_coordinates_vertex_2 (dt, p, out, hole.begin(), hole.end()); diff --git a/Triangulation_2/include/CGAL/Delaunay_triangulation_2.h b/Triangulation_2/include/CGAL/Delaunay_triangulation_2.h index f59ce7e8280..d6ed6fb72fd 100644 --- a/Triangulation_2/include/CGAL/Delaunay_triangulation_2.h +++ b/Triangulation_2/include/CGAL/Delaunay_triangulation_2.h @@ -106,7 +106,7 @@ public: nearest_vertex(const Point& p, Face_handle f= Face_handle()) const; bool does_conflict(const Point &p, Face_handle fh) const;// deprecated - bool test_conflict(const Point &p, Face_handle fh) const; + bool test_conflict(const Point &p, Face_handle fh, bool strict = true) const; bool find_conflicts(const Point &p, //deprecated std::list& conflicts, Face_handle start= Face_handle() ) const; @@ -404,7 +404,8 @@ public: get_conflicts_and_boundary(const Point &p, OutputItFaces fit, OutputItBoundaryEdges eit, - Face_handle start = Face_handle()) const { + Face_handle start = Face_handle(), + bool strict = true) const { CGAL_triangulation_precondition( this->dimension() == 2); int li; Locate_type lt; @@ -419,9 +420,9 @@ public: *fit++ = fh; //put fh in OutputItFaces std::pair pit = std::make_pair(fit,eit); - pit = propagate_conflicts(p,fh,0,pit); - pit = propagate_conflicts(p,fh,1,pit); - pit = propagate_conflicts(p,fh,2,pit); + pit = propagate_conflicts(p,fh,0,pit, strict); + pit = propagate_conflicts(p,fh,1,pit, strict); + pit = propagate_conflicts(p,fh,2,pit, strict); return pit; } CGAL_triangulation_assertion(false); @@ -432,9 +433,10 @@ public: OutputItFaces get_conflicts (const Point &p, OutputItFaces fit, - Face_handle start= Face_handle()) const { + Face_handle start= Face_handle(), + bool strict = true) const { std::pair pp = - get_conflicts_and_boundary(p,fit,Emptyset_iterator(),start); + get_conflicts_and_boundary(p,fit,Emptyset_iterator(),start, strict); return pp.first; } @@ -442,9 +444,10 @@ public: OutputItBoundaryEdges get_boundary_of_conflicts(const Point &p, OutputItBoundaryEdges eit, - Face_handle start= Face_handle()) const { + Face_handle start= Face_handle(), + bool strict = true) const { std::pair pp = - get_conflicts_and_boundary(p,Emptyset_iterator(),eit,start); + get_conflicts_and_boundary(p,Emptyset_iterator(),eit,start,strict); return pp.second; } @@ -460,15 +463,16 @@ private: const Face_handle fh, const int i, std::pair - pit) const { + pit, + bool strict = true) const { Face_handle fn = fh->neighbor(i); - if (! test_conflict(p,fn)) { + if (! test_conflict(p,fn,strict)) { *(pit.second)++ = Edge(fn, fn->index(fh)); } else { *(pit.first)++ = fn; int j = fn->index(fh); - pit = propagate_conflicts(p,fn,ccw(j),pit); - pit = propagate_conflicts(p,fn,cw(j), pit); + pit = propagate_conflicts(p,fn,ccw(j),pit,strict); + pit = propagate_conflicts(p,fn,cw(j), pit,strict); } return pit; } @@ -478,7 +482,8 @@ private: non_recursive_propagate_conflicts ( const Point &p, const Face_handle fh, const int i, - std::pair pit) const + std::pair pit, + bool strict = true) const { std::stack > stack; stack.push( std::make_pair(fh,i) ); @@ -488,7 +493,7 @@ private: const int i=stack.top().second; stack.pop(); Face_handle fn = fh->neighbor(i); - if (! test_conflict(p,fn)) { + if (! test_conflict(p,fn,strict)) { *(pit.second)++ = Edge(fn, fn->index(fh)); } else { *(pit.first)++ = fn; @@ -507,19 +512,20 @@ private: const int i, std::pair pit, + bool strict = true, int depth=0) const { if (depth == 100) - return non_recursive_propagate_conflicts(p, fh, i, pit); + return non_recursive_propagate_conflicts(p, fh, i, pit, strict); Face_handle fn = fh->neighbor(i); - if (! test_conflict(p,fn)) { + if (! test_conflict(p,fn,strict)) { *(pit.second)++ = Edge(fn, fn->index(fh)); } else { *(pit.first)++ = fn; int j = fn->index(fh); - pit = propagate_conflicts(p,fn,ccw(j),pit,depth+1); - pit = propagate_conflicts(p,fn,cw(j), pit,depth+1); + pit = propagate_conflicts(p,fn,ccw(j),pit, strict, depth+1); + pit = propagate_conflicts(p,fn,cw(j), pit, strict, depth+1); } return pit; } @@ -618,8 +624,12 @@ protected: template < class Gt, class Tds > inline bool Delaunay_triangulation_2:: -test_conflict(const Point &p, Face_handle fh) const +test_conflict(const Point &p, Face_handle fh, bool strict) const { + if(! strict){ + Oriented_side os = side_of_oriented_circle(fh,p,false); + return os == ON_POSITIVE_SIDE; + } // 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