From 14150ef95c2d653e3bb6d29f4d499c8a5dd7ebbb Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Mael=20Rouxel-Labb=C3=A9?= Date: Mon, 14 Jan 2019 11:15:28 +0100 Subject: [PATCH 1/2] Revert commit b307771 ("bug fix for natural_neighbor_coordinates_2... ... using a small feature in the conflict zone functions of Delaunay_triangulation_2") See concerns raised in the (not-actually-approved) small feature: https://cgal.geometryfactory.com/CGAL/Members/wiki/Features/Small_Features/DT2_strict_and_weak_conflict_zone --- .../CGAL/natural_neighbor_coordinates_2.h | 2 +- .../include/CGAL/Delaunay_triangulation_2.h | 55 ++++++++----------- 2 files changed, 23 insertions(+), 34 deletions(-) diff --git a/Interpolation/include/CGAL/natural_neighbor_coordinates_2.h b/Interpolation/include/CGAL/natural_neighbor_coordinates_2.h index 4b94385549f..0a18eb7d1a6 100644 --- a/Interpolation/include/CGAL/natural_neighbor_coordinates_2.h +++ b/Interpolation/include/CGAL/natural_neighbor_coordinates_2.h @@ -124,7 +124,7 @@ natural_neighbors_2(const Dt& dt, } std::list hole; - dt.get_boundary_of_conflicts(p, std::back_inserter(hole), fh, false); + dt.get_boundary_of_conflicts(p, std::back_inserter(hole), fh); return natural_neighbors_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 4e127cfc2ba..0ea93fcdd0f 100644 --- a/Triangulation_2/include/CGAL/Delaunay_triangulation_2.h +++ b/Triangulation_2/include/CGAL/Delaunay_triangulation_2.h @@ -114,7 +114,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, bool strict = true) const; + bool test_conflict(const Point &p, Face_handle fh) const; bool find_conflicts(const Point &p, //deprecated std::list& conflicts, Face_handle start= Face_handle()) const; @@ -411,8 +411,7 @@ public: get_conflicts_and_boundary(const Point &p, OutputItFaces fit, OutputItBoundaryEdges eit, - Face_handle start = Face_handle(), - bool strict = true) const + Face_handle start = Face_handle()) const { CGAL_triangulation_precondition(this->dimension() == 2); int li; @@ -427,9 +426,9 @@ public: case Triangulation::OUTSIDE_CONVEX_HULL: *fit++ = fh; //put fh in OutputItFaces std::pair pit = std::make_pair(fit,eit); - pit = propagate_conflicts(p,fh,0,pit, strict); - pit = propagate_conflicts(p,fh,1,pit, strict); - pit = propagate_conflicts(p,fh,2,pit, strict); + pit = propagate_conflicts(p,fh,0,pit); + pit = propagate_conflicts(p,fh,1,pit); + pit = propagate_conflicts(p,fh,2,pit); return pit; } CGAL_triangulation_assertion(false); @@ -440,11 +439,10 @@ public: OutputItFaces get_conflicts (const Point &p, OutputItFaces fit, - Face_handle start= Face_handle(), - bool strict = true) const + Face_handle start= Face_handle()) const { std::pair pp = - get_conflicts_and_boundary(p, fit, Emptyset_iterator(), start, strict); + get_conflicts_and_boundary(p, fit, Emptyset_iterator(), start); return pp.first; } @@ -452,11 +450,10 @@ public: OutputItBoundaryEdges get_boundary_of_conflicts(const Point &p, OutputItBoundaryEdges eit, - Face_handle start= Face_handle(), - bool strict = true) const + Face_handle start= Face_handle()) const { std::pair pp = - get_conflicts_and_boundary(p, Emptyset_iterator(), eit, start, strict); + get_conflicts_and_boundary(p, Emptyset_iterator(), eit, start); return pp.second; } @@ -467,18 +464,16 @@ private: propagate_conflicts (const Point &p, const Face_handle fh, const int i, - std::pair - pit, - bool strict = true) const + std::pair pit) const { Face_handle fn = fh->neighbor(i); - if(! test_conflict(p,fn,strict)) { + if(! test_conflict(p,fn)) { *(pit.second)++ = Edge(fn, fn->index(fh)); } else { *(pit.first)++ = fn; int j = fn->index(fh); - pit = propagate_conflicts(p,fn,ccw(j),pit,strict); - pit = propagate_conflicts(p,fn,cw(j), pit,strict); + pit = propagate_conflicts(p,fn,ccw(j),pit); + pit = propagate_conflicts(p,fn,cw(j), pit); } return pit; } @@ -488,8 +483,7 @@ private: non_recursive_propagate_conflicts(const Point &p, const Face_handle fh, const int i, - std::pair pit, - bool strict = true) const + std::pair pit) const { std::stack > stack; stack.push(std::make_pair(fh,i)); @@ -499,7 +493,8 @@ private: const int i=stack.top().second; stack.pop(); Face_handle fn = fh->neighbor(i); - if(! test_conflict(p,fn,strict)) { + if(! test_conflict(p,fn)) + { *(pit.second)++ = Edge(fn, fn->index(fh)); } else { *(pit.first)++ = fn; @@ -516,22 +511,20 @@ private: propagate_conflicts (const Point &p, const Face_handle fh, const int i, - std::pair - pit, - bool strict = true, + std::pair pit, int depth=0) const { if(depth == 100) - return non_recursive_propagate_conflicts(p, fh, i, pit, strict); + return non_recursive_propagate_conflicts(p, fh, i, pit); Face_handle fn = fh->neighbor(i); - if(! test_conflict(p,fn,strict)) { + if(! test_conflict(p,fn)) { *(pit.second)++ = Edge(fn, fn->index(fh)); } else { *(pit.first)++ = fn; int j = fn->index(fh); - pit = propagate_conflicts(p,fn,ccw(j),pit, strict, depth+1); - pit = propagate_conflicts(p,fn,cw(j), pit, strict, depth+1); + pit = propagate_conflicts(p,fn,ccw(j),pit, depth+1); + pit = propagate_conflicts(p,fn,cw(j), pit, depth+1); } return pit; } @@ -636,12 +629,8 @@ protected: template < class Gt, class Tds > inline bool Delaunay_triangulation_2:: -test_conflict(const Point &p, Face_handle fh, bool strict) const +test_conflict(const Point &p, Face_handle fh) 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 From 48185ac1533f327ffbd87ac8d3ccfc51610d763b Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Mael=20Rouxel-Labb=C3=A9?= Date: Mon, 14 Jan 2019 11:20:16 +0100 Subject: [PATCH 2/2] Filter out irrelevant (<= 0) natural and regular neighbor coordinates - coord == 0 is possible in theory, but not interesting to output - coord < 0 is not possible in theory, but since we do construction and use signed (polygon_area_2()) area computations, we could in theory get negative coordinates in output --- .../CGAL/natural_neighbor_coordinates_2.h | 34 +++++++++++++------ .../CGAL/regular_neighbor_coordinates_2.h | 14 +++++--- 2 files changed, 33 insertions(+), 15 deletions(-) diff --git a/Interpolation/include/CGAL/natural_neighbor_coordinates_2.h b/Interpolation/include/CGAL/natural_neighbor_coordinates_2.h index 0a18eb7d1a6..d2df19d72f0 100644 --- a/Interpolation/include/CGAL/natural_neighbor_coordinates_2.h +++ b/Interpolation/include/CGAL/natural_neighbor_coordinates_2.h @@ -99,21 +99,28 @@ natural_neighbors_2(const Dt& dt, Point_2 p1(v1->point()), p2(v2->point()); Coord_type coef1(0); - Coord_type coef2(0); Equal_x_2 equal_x_2; if(!equal_x_2(p1,p2)) - { coef1 = (p.x() - p2.x()) / (p1.x() - p2.x()); - coef2 = 1 - coef1; - *out++ = std::make_pair(v1,coef1); - *out++ = std::make_pair(v2,coef2); - } else { + else coef1 = (p.y() - p2.y()) / (p1.y() - p2.y()); - coef2 = 1 - coef1; - *out++ = std::make_pair(v1,coef1); - *out++ = std::make_pair(v2,coef2); + + if(coef1 == 0) + { + *out++ = std::make_pair(v2, Coord_type(1)); + return make_triple(out, Coord_type(1), true); } + Coord_type coef2 = 1 - coef1; + if(coef2 == 0) + { + *out++ = std::make_pair(v1, Coord_type(1)); + return make_triple(out, Coord_type(1), true); + } + + *out++ = std::make_pair(v1,coef1); + *out++ = std::make_pair(v2,coef2); + return make_triple(out, coef1+coef2, true); } @@ -142,6 +149,7 @@ natural_neighbors_2(const Dt& dt, EdgeIterator hole_begin, EdgeIterator hole_end) { CGAL_precondition(dt.dimension() == 2); + typedef typename Dt::Geom_traits Traits; typedef typename Traits::FT Coord_type; typedef typename Traits::Point_2 Point_2; @@ -187,8 +195,12 @@ natural_neighbors_2(const Dt& dt, p); area += polygon_area_2(vor.begin(), vor.end(), dt.geom_traits()); - *out++ = std::make_pair(current,area); - area_sum += area; + + if(area > 0) + { + *out++ = std::make_pair(current,area); + area_sum += area; + } //update prev and hit: prev = current; diff --git a/Interpolation/include/CGAL/regular_neighbor_coordinates_2.h b/Interpolation/include/CGAL/regular_neighbor_coordinates_2.h index bbb03227e5c..5a83dd079c1 100644 --- a/Interpolation/include/CGAL/regular_neighbor_coordinates_2.h +++ b/Interpolation/include/CGAL/regular_neighbor_coordinates_2.h @@ -131,9 +131,12 @@ regular_neighbor_coordinates_vertex_2(const Rt& rt, *vor_vertices++ = vor[2]; area += polygon_area_2(vor.begin(), vor.end(), rt.geom_traits()); - *out++= std::make_pair(current, area); - area_sum += area; + if(area > 0) + { + *out++= std::make_pair(current, area); + area_sum += area; + } //update prev and hit: prev = current; @@ -163,8 +166,11 @@ regular_neighbor_coordinates_vertex_2(const Rt& rt, ++fc; } - *out++ = std::make_pair((*hidden_vertices_begin), area); - area_sum += area; + if(area > 0) + { + *out++ = std::make_pair((*hidden_vertices_begin), area); + area_sum += area; + } } return make_triple(out, area_sum, true);