diff --git a/Polygon_mesh_processing/test/Polygon_mesh_processing/test_pmp_distance.cpp b/Polygon_mesh_processing/test/Polygon_mesh_processing/test_pmp_distance.cpp index 5ecd9b9530b4..1236f1b8d094 100644 --- a/Polygon_mesh_processing/test/Polygon_mesh_processing/test_pmp_distance.cpp +++ b/Polygon_mesh_processing/test/Polygon_mesh_processing/test_pmp_distance.cpp @@ -22,6 +22,7 @@ struct Custom_traits_Hausdorff FT operator+(FT)const{return FT();} FT operator*(FT)const{return FT();} bool operator<=(FT)const{return false;} + bool operator==(FT)const{return false;} bool operator>=(FT)const{return false;} bool operator!=(FT)const{return false;} bool operator<(FT)const{return false;} diff --git a/Spatial_searching/doc/Spatial_searching/CGAL/Kd_tree.h b/Spatial_searching/doc/Spatial_searching/CGAL/Kd_tree.h index b36cc8fb9733..c83d721efa8a 100644 --- a/Spatial_searching/doc/Spatial_searching/CGAL/Kd_tree.h +++ b/Spatial_searching/doc/Spatial_searching/CGAL/Kd_tree.h @@ -20,7 +20,7 @@ lot of cache misses, leading to non-optimal performance. This is the case for example when indices are stored inside the tree, or to a lesser extent when the points coordinates are stored in a dynamically allocated array (e.g., `Epick_d` with dynamic -dimension) — we says "to a lesser extent" because the points +dimension) — we say "to a lesser extent" because the points are re-created by the kd-tree in a cache-friendly order after its construction, so the coordinates are more likely to be stored in a near-optimal order on the heap. When `EnablePointsCache` is set to `Tag_true`, the points diff --git a/Spatial_searching/doc/Spatial_searching/Spatial_searching.txt b/Spatial_searching/doc/Spatial_searching/Spatial_searching.txt index afa3a3eb0ca8..9ccd82722ed4 100644 --- a/Spatial_searching/doc/Spatial_searching/Spatial_searching.txt +++ b/Spatial_searching/doc/Spatial_searching/Spatial_searching.txt @@ -127,7 +127,7 @@ which is the case for point queries. The following two classes implement the standard search strategy for orthogonal distances like the weighted Minkowski distance. The second one is a specialization for incremental neighbor -searching and distance browsing. Both require extendes nodes. +searching and distance browsing. Both require extended nodes. `Orthogonal_k_neighbor_search` @@ -175,9 +175,13 @@ because in general the query time will be less than using the default value. Instead of using the default splitting rule `Sliding_midpoint` described below, a user may, depending upon the data, select one from the following splitting rules, -which determine how a separating hyperplane is computed. Every splitter has +which determine how a separating hyperplane is computed. Some splitter have degenerated worst cases, which may lead to a linear tree and a stack overflow. Switching the splitting rule to one of a different kind will solve the problem. +The `Median_of_rectangle` and `Median_of_max_spread` are robust sliders that will +neither lead to a linear tree nor to a stack overflow. The `Median_of_rectangle` +splitter will detect if the data in a node is degenerated and applies the +`Median_of_max_spread` rule for that node to avoid a linear tree.
diff --git a/Spatial_searching/include/CGAL/Kd_tree.h b/Spatial_searching/include/CGAL/Kd_tree.h index 350b9f533d8e..a1d1e670d2b0 100644 --- a/Spatial_searching/include/CGAL/Kd_tree.h +++ b/Spatial_searching/include/CGAL/Kd_tree.h @@ -196,7 +196,7 @@ class Kd_tree { if (try_parallel_internal_node_creation (nh, c, c_low, tag)) return; - if (c_low.size() > split.bucket_size()) + if (c_low.size() > split.bucket_size() && !CGAL::is_zero(c_low.max_tight_spread())) { nh->lower_ch = new_internal_node(); create_internal_node (nh->lower_ch, c_low, tag); @@ -204,7 +204,7 @@ class Kd_tree { else nh->lower_ch = create_leaf_node(c_low); - if (c.size() > split.bucket_size()) + if (c.size() > split.bucket_size() && !CGAL::is_zero(c.max_tight_spread())) { nh->upper_ch = new_internal_node(); create_internal_node (nh->upper_ch, c, tag); @@ -341,7 +341,7 @@ class Kd_tree { Point_container c(dim_, data.begin(), data.end(),traits_); bbox = new Kd_tree_rectangle(c.bounding_box()); - if (c.size() <= split.bucket_size()){ + if (c.size() <= split.bucket_size() || CGAL::is_zero(c.max_tight_spread())){ tree_root = create_leaf_node(c); }else { tree_root = new_internal_node(); diff --git a/Spatial_searching/include/CGAL/Kd_tree_rectangle.h b/Spatial_searching/include/CGAL/Kd_tree_rectangle.h index ef2f13476d93..ba111699f079 100644 --- a/Spatial_searching/include/CGAL/Kd_tree_rectangle.h +++ b/Spatial_searching/include/CGAL/Kd_tree_rectangle.h @@ -111,11 +111,7 @@ namespace CGAL { explicit Kd_tree_rectangle(const Kd_tree_rectangle& r) - : max_span_coord_(r.max_span_coord_) - { - lower_ = r.lower_; - upper_ = r.upper_; - } + : lower_(r.lower_), upper_(r.upper_), max_span_coord_(r.max_span_coord_) {} template void update_from_point_pointers(PointPointerIter begin, diff --git a/Spatial_searching/include/CGAL/Point_container.h b/Spatial_searching/include/CGAL/Point_container.h index eb96024348e2..07ad4053079c 100644 --- a/Spatial_searching/include/CGAL/Point_container.h +++ b/Spatial_searching/include/CGAL/Point_container.h @@ -232,8 +232,9 @@ class Point_container { // building the container from a sequence of Point_d* Point_container(const int d, iterator begin, iterator end,const Traits& traits_) : - traits(traits_),m_b(begin), m_e(end), bbox(d, begin, end,traits.construct_cartesian_const_iterator_d_object()), tbox(bbox) + traits(traits_),m_b(begin), m_e(end), bbox(d, begin, end,traits.construct_cartesian_const_iterator_d_object()), tbox() { + tbox = bbox; built_coord = max_span_coord(); } @@ -419,7 +420,28 @@ class Point_container { typename Traits::Cartesian_const_iterator_d mpit = construct_it((*(*mid))); FT val1 = *(mpit+split_coord); + + // Avoid using the low coord value as it results in an empty split + if (val1 == tbox.min_coord(split_coord)) { + iterator it = std::min_element(mid, end(), [=](const Point_d* a, const Point_d* b) -> bool { + FT a_c = *(construct_it(*a) + split_coord); + FT b_c = *(construct_it(*b) + split_coord); + + if (a_c == val1) + return false; + + if (b_c == val1) + return true; + + return a_c < b_c; + }); + return *(construct_it(**it) + split_coord); + } + mid++; + if (mid == end()) + return val1; + mpit = construct_it((*(*mid))); FT val2 = *(mpit+split_coord); return (val1+val2)/FT(2); diff --git a/Spatial_searching/include/CGAL/Splitters.h b/Spatial_searching/include/CGAL/Splitters.h index 057c27aae43d..387144bcf8a3 100644 --- a/Spatial_searching/include/CGAL/Splitters.h +++ b/Spatial_searching/include/CGAL/Splitters.h @@ -19,6 +19,7 @@ #include +#include #include #include @@ -236,7 +237,9 @@ namespace CGAL { void operator() (Separator& sep, Container& c0, Container& c1) const { - sep = Separator(c0.max_span_coord(),FT(0)); + if (!CGAL::is_zero(c0.tight_bounding_box().max_coord(c0.max_span_coord()) - c0.tight_bounding_box().min_coord(c0.max_span_coord()))) + sep = Separator(c0.max_span_coord(),FT(0)); + else sep = Separator(c0.max_tight_span_coord(), FT(0)); sep.set_cutting_value(c0.median(sep.cutting_dimension())); c0.split(c1,sep,true); }