diff --git a/Nef_3/include/CGAL/Nef_3/K3_tree.h b/Nef_3/include/CGAL/Nef_3/K3_tree.h index 8329767d8eb..6dd6ef1c84a 100644 --- a/Nef_3/include/CGAL/Nef_3/K3_tree.h +++ b/Nef_3/include/CGAL/Nef_3/K3_tree.h @@ -41,8 +41,7 @@ template class K3_tree { -template +template class Smaller_than { public: @@ -59,26 +58,13 @@ public: return false; } - bool operator()( const Object& o1, const Object& o2) { - Vertex v1,v2; - CGAL::assign(v1,o1); - CGAL::assign(v2,o2); - switch(coord) { - case 0: return CGAL::compare_x(v1->point(), v2->point()) == SMALLER; - case 1: return CGAL::compare_y(v1->point(), v2->point()) == SMALLER; - case 2: return CGAL::compare_z(v1->point(), v2->point()) == SMALLER; - default: CGAL_error(); - } - return false; - } private: Coordinate coord; }; -template - class Smaller_than, Object, Vertex, Coordinate> +template + class Smaller_than, Vertex, Coordinate> { public: Smaller_than(Coordinate c) : coord(c) { @@ -97,21 +83,6 @@ public: return false; } - bool operator()( const Object& o1, const Object& o2) { - Vertex v1,v2; - CGAL::assign(v1,o1); - CGAL::assign(v2,o2); - switch(coord) { - case 0: return CGAL::to_interval(v1->point().x()).second < - CGAL::to_interval(v2->point().x()).first; - case 1: return CGAL::to_interval(v1->point().y()).second < - CGAL::to_interval(v2->point().y()).first; - case 2: return CGAL::to_interval(v1->point().z()).second < - CGAL::to_interval(v2->point().z()).first; - default: CGAL_error(); - } - return false; - } private: Coordinate coord; }; @@ -126,13 +97,17 @@ public: typedef typename Traits::SNC_decorator SNC_decorator; typedef typename Traits::Infimaximal_box Infimaximal_box; typedef typename Traits::Vertex_handle Vertex_handle; +typedef typename Traits::Vertex_list Vertex_list; +typedef typename Vertex_list::iterator Vertex_iterator; +typedef typename Vertex_list::const_iterator Vertex_const_iterator; typedef typename Traits::Halfedge_handle Halfedge_handle; typedef typename Traits::Halffacet_handle Halffacet_handle; typedef typename Traits::Object_handle Object_handle; -typedef std::vector Object_list; -typedef typename Object_list::const_iterator Object_const_iterator; -typedef typename Object_list::iterator Object_iterator; -typedef typename Object_list::size_type size_type; +typedef typename Traits::Object_list Object_list; +typedef typename Traits::Halfedge_list Halfedge_list; +typedef typename Halfedge_list::const_iterator Halfedge_const_iterator; +typedef typename Traits::Halffacet_list Halffacet_list; +typedef typename Halffacet_list::const_iterator Halffacet_const_iterator; typedef typename Traits::Point_3 Point_3; typedef typename Traits::Segment_3 Segment_3; @@ -144,7 +119,6 @@ typedef typename Traits::Aff_transformation_3 Aff_transformation_3; typedef typename Traits::Bounding_box_3 Bounding_box_3; typedef typename Traits::Side_of_plane Side_of_plane; -typedef typename Traits::Objects_bbox Objects_bbox; typedef typename Traits::Kernel Kernel; typedef typename Kernel::RT RT; @@ -152,29 +126,36 @@ typedef typename Kernel::FT FT; typedef Smaller_than< Kernel, - Object_handle, Vertex_handle, - int> Smaller_; + int> Smaller; class Node { friend class K3_tree; public: typedef Node* Node_handle; - - Node( Node_handle p, Node_handle l, Node_handle r, Plane_3 pl, const Object_list& O) : - parent_node(p), left_node(l), right_node(r), splitting_plane(pl), - object_list(O) { - if(l == nullptr) - point_on_plane = Point_3(); - else - point_on_plane = pl.point(); + Node(const Vertex_list& V, const Halfedge_list& E, const Halffacet_list& F) : + left_node(nullptr), right_node(nullptr) + { + object_list.reserve(V.size()+E.size()+F.size()); + for(Vertex_const_iterator vi=V.begin(); vi!=V.end(); ++vi) + object_list.push_back(make_object(*vi)); + for(Halfedge_const_iterator ei=E.begin(); ei!=E.end(); ++ei) + object_list.push_back(make_object(*ei)); + for(Halffacet_const_iterator fi=F.begin(); fi!=F.end(); ++fi) + object_list.push_back(make_object(*fi)); } + + Node(Node_handle l, Node_handle r, const Plane_3& pl) : + left_node(l), right_node(r), splitting_plane(pl) + { + } + bool is_leaf() const { CGAL_assertion( (left_node != nullptr && right_node != nullptr) || (left_node == nullptr && right_node == nullptr)); return (left_node == nullptr && right_node == nullptr); } - Node_handle parent() const { return parent_node; } + Node_handle left() const { return left_node; } Node_handle right() const { return right_node; } const Plane_3& plane() const { return splitting_plane; } @@ -189,45 +170,42 @@ public: } } - template - void add_facet(Halffacet_handle f, Depth depth) { + void add_facet(Halffacet_handle f, int depth) { if(left_node == nullptr) { object_list.push_back(make_object(f)); return; } - Side_of_plane sop; - Oriented_side side = sop(splitting_plane.point(), f, depth); + Side_of_plane sop(splitting_plane.point(), depth); + Oriented_side side = sop(f); if( side == ON_NEGATIVE_SIDE || side == ON_ORIENTED_BOUNDARY) left_node->add_facet(f, depth+1); if( side == ON_POSITIVE_SIDE || side == ON_ORIENTED_BOUNDARY) right_node->add_facet(f, depth+1); } - template - void add_edge(Halfedge_handle e, Depth depth) { + void add_edge(Halfedge_handle e, int depth) { if(left_node == nullptr) { object_list.push_back(make_object(e)); return; } - Side_of_plane sop; - Oriented_side side = sop(splitting_plane.point(), e, depth); + Side_of_plane sop(splitting_plane.point(), depth); + Oriented_side side = sop(e); if( side == ON_NEGATIVE_SIDE || side == ON_ORIENTED_BOUNDARY) left_node->add_edge(e, depth+1); if( side == ON_POSITIVE_SIDE || side == ON_ORIENTED_BOUNDARY) right_node->add_edge(e, depth+1); } - template - void add_vertex(Vertex_handle v, Depth depth) { + void add_vertex(Vertex_handle v, int depth) { if(left_node == nullptr) { object_list.push_back(make_object(v)); return; } - Side_of_plane sop; - Oriented_side side = sop(splitting_plane.point(), v, depth); + Side_of_plane sop(splitting_plane.point(), depth); + Oriented_side side = sop(v); if( side == ON_NEGATIVE_SIDE || side == ON_ORIENTED_BOUNDARY) left_node->add_vertex(v, depth+1); if( side == ON_POSITIVE_SIDE || side == ON_ORIENTED_BOUNDARY) @@ -252,27 +230,11 @@ friend std::ostream& operator<< return os; } - /* -~Node() noexcept(!CGAL_ASSERTIONS_ENABLED) -{ - CGAL_NEF_TRACEN("~Node: deleting node..."); - CGAL_destructor_assertion_catch( - if( !is_leaf()) { - delete left_node; - delete right_node; - } - ); -} - */ private: - - - Node_handle parent_node; Node_handle left_node; Node_handle right_node; Plane_3 splitting_plane; - Point_3 point_on_plane; Object_list object_list; }; @@ -437,23 +399,15 @@ void divide_segment_by_plane( Segment_3 s, Plane_3 pl, // First of all, we need to find out wheather we are working over an extended kernel or on a standard kernel. As precondition we have that ray is oriented in the minus x axis direction. When having an extended kernel, the ray can be subtituted by a segment with the endpoint on the 'intersection' between the ray and the bounding infimaximal box. In the presence of a standard kernel, the intersection is computed with the bounding box with the vertices of the Nef polyhedron. Point_3 p(r.source()), q; Bounding_box_3 b = k.bounding_box; - int c = (CGAL::abs(vec[0]) > CGAL::abs(vec[1]) ? 0 : 1); - c = (CGAL::abs(vec[2]) > CGAL::abs(vec[c]) ? 2 : c); + typename Kernel::Non_zero_dimension_3 non_zero_dimension_3; + int c = non_zero_dimension_3(vec); Point_3 pt_on_minus_x_plane = vec[c] < 0 ? Point_3(FT(b.min_coord(0)), FT(b.min_coord(1)),FT(b.min_coord(2))) : Point_3(FT(b.max_coord(0)), FT(b.max_coord(1)),FT(b.max_coord(2))); // We compute the intersection between a plane with normal vector in // the minus x direction and located at the minimum point of the bounding box, and the input ray. When the ray does not intersect the bounding volume, there won't be any object hit, so it is safe to construct a segment that simply lay in the unbounded side of the bounding box. This approach is taken instead of somehow (efficiently) report that there was no hit object, in order to mantain a clear interface with the Iterator class. - Plane_3 pl_on_minus_x; - if(c==0) - pl_on_minus_x = Plane_3(pt_on_minus_x_plane, Vector_3( 1, 0, 0)); - else if(c==1) - pl_on_minus_x = Plane_3(pt_on_minus_x_plane, Vector_3( 0, 1, 0)); - else { - CGAL_assertion_msg(c==2, "wrong value"); - pl_on_minus_x = Plane_3(pt_on_minus_x_plane, Vector_3( 0, 0, 1)); - } + Plane_3 pl_on_minus_x = K3_tree::construct_splitting_plane(pt_on_minus_x_plane, c); Object o = traits.intersect_object()( pl_on_minus_x, r); if( !CGAL::assign( q, o) || pl_on_minus_x.has_on(p)) q = r.source() + vec; @@ -471,7 +425,7 @@ private: Node_handle root; - boost::container::deque nodes; + Node_range nodes; int max_depth; Bounding_box_3 bounding_box; @@ -488,39 +442,40 @@ public: typedef typename SNC_structure::Halffacet_iterator Halffacet_iterator; CGAL_assertion( W != nullptr); - Object_list objects; + + Vertex_list vertices; + vertices.reserve(W->number_of_vertices()); Vertex_iterator v; + CGAL_forall_vertices(v, *W) + vertices.push_back(v); + + Halfedge_list edges; + edges.reserve(W->number_of_halfedges()); Halfedge_iterator e; + CGAL_forall_edges(e, *W) + edges.push_back(e); + + Halffacet_list facets; + facets.reserve(W->number_of_halffacets()); Halffacet_iterator f; - CGAL_forall_vertices( v, *W) - objects.push_back(make_object(Vertex_handle(v))); - typename Object_list::difference_type v_end = objects.size(); - CGAL_forall_edges( e, *W) - objects.push_back(make_object(Halfedge_handle(e))); - CGAL_forall_facets( f, *W) { - objects.push_back(make_object(Halffacet_handle(f))); - } - Object_iterator oli=objects.begin()+v_end; - root = build_kdtree( objects, oli, 0); - } + CGAL_forall_facets(f, *W) + facets.push_back(f); - K3_tree(Object_list& objects, Object_iterator& v_end) { + CGAL_NEF_TRACEN("K3_tree(): n_vertices = " << vertices.size()); + std::frexp( (double) vertices.size(), &max_depth); -typename Object_list::difference_type n_vertices = std::distance(objects.begin(),v_end); - CGAL_NEF_TRACEN("K3_tree(): n_vertices = " << std::distance(objects.begin(),v_end)); - std::frexp( (double) n_vertices, &max_depth); - - // TODO: in the presence of a infimaximal bounding box, the bounding box does not have to be computed - Objects_bbox objects_bbox = traits.objects_bbox_object(); - bounding_box = objects_bbox(objects); - //CGAL_NEF_TRACEN("bounding box:"<point()); + //CGAL_NEF_TRACEN("bounding box:"< -Node_handle build_kdtree(Object_list& O, Object_iterator v_end, - Depth depth, Node_handle parent=nullptr, int non_efective_splits=0) { +Node_handle build_kdtree(Vertex_list& V, Halfedge_list& E, Halffacet_list& F, + int depth, int non_efective_splits=0) { CGAL_precondition( depth >= 0); - CGAL_NEF_TRACEN( "build_kdtree: "< max_depth) - nodes.push_back(Node( parent, nullptr, nullptr, Plane_3(), O)); + nodes.push_back(Node(V, E, F)); return &(nodes.back()); - } else { - CGAL_NEF_TRACEN("Sizes " << O1.size() << ", " << O2.size() << ", " << O.size()); - CGAL_assertion( O1.size() <= O.size() && O2.size() <= O.size()); - CGAL_assertion( O1.size() + O2.size() >= O.size()); - non_efective_split = ((O1.size() == O.size()) || (O2.size() == O.size())); } - if( non_efective_split) + auto O_size = V.size() + E.size() + F.size(); + auto O1_size = V1.size() + E1.size() + F1.size(); + auto O2_size = V2.size() + E2.size() + F2.size(); + CGAL_NEF_TRACEN("Sizes " << O1_size << ", " << O2_size << ", " << O_size); + CGAL_assertion( O1_size <= O_size && O2_size <= O_size); + CGAL_assertion( O1_size + O2_size >= O_size); + if((O1_size == O_size) || (O2_size == O_size)) non_efective_splits++; else non_efective_splits = 0; - if( non_efective_splits > 2) { + + if(non_efective_splits > 2) { CGAL_NEF_TRACEN("build_kdtree: non efective splits reached maximum"); - nodes.push_back(Node( parent, nullptr, nullptr, Plane_3(), O)); + nodes.push_back(Node(V, E, F)); return &(nodes.back()); } - nodes.push_back(Node( parent, nullptr, nullptr, partition_plane, Object_list())); - Node_handle node = &(nodes.back()); - node->left_node = build_kdtree( O1, O1.begin()+v_end1, depth + 1, node, non_efective_splits); - node->right_node = build_kdtree( O2, O2.begin()+v_end2, depth + 1, node, non_efective_splits); - return node; + + Node_handle left_node = build_kdtree(V1, E1, F1, depth + 1, non_efective_splits); + Node_handle right_node = build_kdtree(V2, E2, F2, depth + 1, non_efective_splits); + nodes.push_back(Node(left_node, right_node, construct_splitting_plane(point_on_plane, coord))); + return &(nodes.back()); } -template -bool can_set_be_divided(Object_iterator start, Object_iterator end, Depth depth) { +bool can_set_be_divided(int depth, typename Vertex_list::size_type size) { if(depth >= max_depth) return false; - if(std::distance(start,end)<2) + if(size <= 2) return false; return true; } -template -bool classify_objects(Object_iterator start, Object_iterator end, - Plane_3 partition_plane, Side_of_plane& sop, - OutputIterator o1, OutputIterator o2, Depth depth) { - typename Object_list::difference_type on_oriented_boundary = 0; - typename Object_list::const_iterator o; - - Point_3 point_on_plane(partition_plane.point()); - - for( o = start; o != end; ++o) { - Oriented_side side = sop( point_on_plane, *o, depth); - if( side == ON_NEGATIVE_SIDE || side == ON_ORIENTED_BOUNDARY) { - *o1 = *o; - ++o1; +template +static bool classify_objects(const List& l,Side_of_plane& sop, + List& l1, List& l2) { + typename List::size_type on_oriented_boundary = 0; + for(typename List::const_iterator i=l.begin(); i!=l.end(); ++i) { + Oriented_side side = sop(*i); + switch(side) { + case ON_NEGATIVE_SIDE: + l1.push_back(*i); + break; + case ON_POSITIVE_SIDE: + l2.push_back(*i); + break; + case ON_ORIENTED_BOUNDARY: + ++on_oriented_boundary; + l1.push_back(*i); + l2.push_back(*i); } - if( side == ON_POSITIVE_SIDE || side == ON_ORIENTED_BOUNDARY) { - *o2 = *o; - ++o2; - } - if( side == ON_ORIENTED_BOUNDARY) - ++on_oriented_boundary; } - return (on_oriented_boundary != std::distance(start,end)); + return (on_oriented_boundary != l.size()); } +static Point_3 find_median_point(Vertex_list& V, + Vertex_iterator& median, int coord) { + CGAL_assertion(V.size() > 1); -template -Plane_3 construct_splitting_plane(Object_iterator start, Object_iterator end, - Object_iterator& median, Depth depth) { - CGAL_precondition( depth >= 0); - typename Object_list::difference_type n=std::distance(start,end); - CGAL_assertion(n>1); + median = V.begin() + V.size()/2; + std::nth_element(V.begin(), median, V.end(), Smaller(coord)); - std::nth_element(start, start+n/2, end, - Smaller_(depth%3)); + return (*median)->point(); +} - Vertex_handle v; - median = start+n/2; - CGAL::assign(v,*median); - switch( depth % 3) { - case 0: return Plane_3( v->point(), Vector_3( 1, 0, 0)); break; - case 1: return Plane_3( v->point(), Vector_3( 0, 1, 0)); break; - case 2: return Plane_3( v->point(), Vector_3( 0, 0, 1)); break; +static Plane_3 construct_splitting_plane(const Point_3& pt, int coord) +{ + switch(coord) { + case 0: return Plane_3(1, 0, 0, -pt.x()); + case 1: return Plane_3(0, 1, 0, -pt.y()); + case 2: return Plane_3(0, 0, 1, -pt.z()); } CGAL_error_msg( "never reached"); diff --git a/Nef_3/include/CGAL/Nef_3/SNC_k3_tree_traits.h b/Nef_3/include/CGAL/Nef_3/SNC_k3_tree_traits.h index ff49f3e7b21..38201f2ce18 100644 --- a/Nef_3/include/CGAL/Nef_3/SNC_k3_tree_traits.h +++ b/Nef_3/include/CGAL/Nef_3/SNC_k3_tree_traits.h @@ -18,6 +18,7 @@ #include #include +#include #include #undef CGAL_NEF_DEBUG @@ -28,11 +29,11 @@ namespace CGAL { template -class ComparePoints { +class Compare_points { typedef typename Kernel::Point_3 Point_3; public: - ComparePoints(Coordinate c) : coord(c) { + Compare_points(Coordinate c) : coord(c) { CGAL_assertion( c >= 0 && c <=2); } CGAL::Comparison_result operator()(const Point_3& p1, const Point_3& p2) { @@ -56,15 +57,15 @@ private: template -class ComparePoints, Coordinate> { +class Compare_points, Coordinate> { typedef CGAL::Lazy_kernel Kernel; typedef typename Kernel::Point_3 Point_3; public: - ComparePoints(Coordinate c) : coord(c) { + Compare_points(Coordinate c) : coord(c) { CGAL_assertion( c >= 0 && c <=2); } - CGAL::Comparison_result operator()( const Point_3 p1, const Point_3 p2) { + CGAL::Comparison_result operator()( const Point_3& p1, const Point_3& p2) { switch(coord) { case 0: if(CGAL::to_interval(p1.x()).second < @@ -100,16 +101,13 @@ private: template class Side_of_plane { -public: - typedef typename SNC_decorator::SNC_structure SNC_structure; + typedef typename SNC_decorator::Decorator_traits Decorator_traits; typedef typename Decorator_traits::Vertex_handle Vertex_handle; typedef typename Decorator_traits::Halfedge_handle Halfedge_handle; typedef typename Decorator_traits::Halffacet_handle Halffacet_handle; - typedef typename SNC_structure::Object_handle Object_handle; - typedef typename Decorator_traits::Halffacet_cycle_iterator Halffacet_cycle_iterator; typedef typename Decorator_traits::SHalfedge_around_facet_circulator @@ -118,154 +116,34 @@ public: typedef typename SNC_decorator::Kernel Kernel; typedef typename Kernel::Point_3 Point_3; - typedef typename Kernel::Segment_3 Segment_3; - typedef typename Kernel::Plane_3 Plane_3; - typedef typename Kernel::Triangle_3 Triangle_3; - typedef typename Kernel::Vector_3 Vector_3; - typedef typename Kernel::RT RT; - - typedef ComparePoints ComparePoints_; - #ifdef CGAL_NEF_EXPLOIT_REFERENCE_COUNTING - Side_of_plane(bool rc = false) : reference_counted(rc) {} + typedef typename Kernel::RT RT; +#endif + typedef Compare_points Compare; + static constexpr Oriented_side unknown_side = static_cast(-2); + +public: +#ifdef CGAL_NEF_EXPLOIT_REFERENCE_COUNTING + Side_of_plane(const Point_3& p, int c, bool rc = false) : reference_counted(rc), coord(c), pop(p) {} #else - Side_of_plane() {} + Side_of_plane(const Point_3& p, int c) : OnSideMap(unknown_side), coord(c), pop(p) {} #endif - - template Oriented_side operator() - ( const Point_3& pop, Object_handle o, Depth depth); - template Oriented_side operator() - ( const Point_3& pop, Vertex_handle v, Depth depth); - template Oriented_side operator() - ( const Point_3& pop, Halfedge_handle e, Depth depth); - template Oriented_side operator() - ( const Point_3& pop, Halffacet_handle f, Depth depth); -#ifdef CGAL_NEF_EXPLOIT_REFERENCE_COUNTING - bool reference_counted; -#endif - Unique_hash_map OnSideMap; + Oriented_side operator()(Vertex_handle v); + Oriented_side operator()(Halfedge_handle e); + Oriented_side operator()(Halffacet_handle f); +private: + Unique_hash_map OnSideMap; #ifdef CGAL_NEF_EXPLOIT_REFERENCE_COUNTING Unique_hash_map OnSideMapRC; + bool reference_counted; #endif + int coord; + const Point_3& pop; }; template -class Objects_bbox { -public: - typedef typename SNC_decorator::SNC_structure SNC_structure; - typedef typename SNC_decorator::Decorator_traits Decorator_traits; - - typedef typename Decorator_traits::Vertex_handle Vertex_handle; - typedef typename Decorator_traits::Halfedge_handle Halfedge_handle; - typedef typename Decorator_traits::Halffacet_handle Halffacet_handle; - - typedef typename SNC_structure::Halffacet_triangle_handle Halffacet_triangle_handle; - typedef typename SNC_structure::Object_handle Object_handle; - typedef std::vector Object_list; - - typedef typename Decorator_traits::Halffacet_cycle_iterator - Halffacet_cycle_iterator; - typedef typename Decorator_traits::SHalfedge_around_facet_circulator - SHalfedge_around_facet_circulator; - typedef typename Decorator_traits::SHalfedge_handle SHalfedge_handle; - - typedef typename SNC_decorator::Kernel Kernel; - typedef typename Kernel::Plane_3 Plane_3; - typedef typename Kernel::Segment_3 Segment_3; - typedef typename Kernel::Point_3 Point_3; - typedef typename Kernel::Triangle_3 Triangle_3; - - typedef typename Kernel::RT RT; - typedef typename Kernel::FT FT; - typedef typename Kernel::Kernel_tag Kernel_tag; - typedef CGAL::Bounding_box_3 - Bounding_box_3; - // typedef CGAL::Bounding_box_3 - // ::value_type, Kernel> - // Bounding_box_3; - - Bounding_box_3 operator()( const Object_list& O) const { - Bounding_box_3 b; - typename Object_list::const_iterator o = O.begin(); - Vertex_handle v; - while(o != O.end() && !CGAL::assign(v, *o)) ++o; - if(o != O.end()) { - FT q[3]; - q[0] = v->point().x(); - q[1] = v->point().y(); - q[2] = v->point().z(); - Bounding_box_3 b(q); - for(++o; o != O.end(); ++o) { - if( CGAL::assign( v, *o)) { - b.extend(v->point()); - } - } - return b; - } - FT q[3]; - q[0] = q[1] = q[2] = 0; - return Bounding_box_3(q); - } - - /* - Bounding_box_3 operator()(Object_handle o) const { - Vertex_handle v; - Halfedge_handle e; - Halffacet_handle f; - if( CGAL::assign( v, o)) - return operator()(v); - else if( CGAL::assign( e, o)) - return operator()(e); - else if( CGAL::assign( f, o)) - return operator()(f); - else { - Halffacet_triangle_handle t; - if( CGAL::assign( t, o)) - return operator()(t); - else - CGAL_error_msg( "wrong handle"); - } - return Bounding_box_3(); // never reached - } - - Bounding_box_3 operator()(Vertex_handle v) const { - Bounding_box_3 b; - b.extend(v->point()); - return b; - } - - Bounding_box_3 operator()(Halfedge_handle e) const { - Bounding_box_3 b; - b.extend(e->source()->point()); - b.extend(e->twin()->source()->point()); - return b; - } - - Bounding_box_3 operator()(Halffacet_triangle_handle t) const { - Bounding_box_3 bbox; - Triangle_3 tr(t.get_triangle()); - for( int i = 0; i < 3; ++i) - bbox.extend(tr[i]); - return bbox; - } - - Bounding_box_3 operator()(Halffacet_handle f) const { - CGAL_assertion( f->facet_cycles_begin() != - Halffacet_cycle_iterator()); - Halffacet_cycle_iterator fc(f->facet_cycles_begin()); - SHalfedge_handle e; - CGAL_assertion(fc.is_shalfedge()); - e = SHalfedge_handle(fc); - SHalfedge_around_facet_circulator sc(e), send(sc); - CGAL_assertion( !is_empty_range( sc, send)); - Bounding_box_3 b; - CGAL_For_all( sc, send) - b.extend(sc->source()->source()->point()); - return b; - } - */ -}; +constexpr Oriented_side Side_of_plane::unknown_side; template class SNC_k3_tree_traits { @@ -283,6 +161,9 @@ public: typedef typename SNC_structure::Object_handle Object_handle; typedef std::vector Object_list; + typedef std::vector Vertex_list; + typedef std::vector Halfedge_list; + typedef std::vector Halffacet_list; typedef typename Kernel::Point_3 Point_3; typedef typename Kernel::Segment_3 Segment_3; @@ -294,55 +175,24 @@ public: typedef typename Kernel::RT RT; typedef typename Kernel::Kernel_tag Kernel_tag; - // typedef CGAL::Bounding_box_3 - // ::value_type, Kernel> - // Bounding_box_3; - typedef CGAL::Bounding_box_3 - Bounding_box_3; - + typedef CGAL::Bounding_box_3 Bounding_box_3; typedef typename Kernel::Intersect_3 Intersect; - typedef CGAL::Objects_bbox Objects_bbox; typedef CGAL::Side_of_plane Side_of_plane; Intersect intersect_object() const { return Intersect(); } - Objects_bbox objects_bbox_object() const { - return Objects_bbox(); - } }; template -template Oriented_side -Side_of_plane::operator() - (const Point_3& pop, Object_handle o, Depth depth) { - Vertex_handle v; - Halfedge_handle e; - Halffacet_handle f; - if( CGAL::assign( v, o)) - return (*this)(pop, v, depth); - else if( CGAL::assign( e, o)) - return (*this)(pop, e, depth); - else if( CGAL::assign( f, o)) - return (*this)(pop, f, depth); - else - CGAL_error_msg( "wrong handle"); - - return Oriented_side(); // never reached -} - -template -template -Oriented_side -Side_of_plane::operator() -( const Point_3& pop, Vertex_handle v, Depth depth) { - Comparison_result cr; +Side_of_plane::operator()(Vertex_handle v) { +Comparison_result cr; #ifdef CGAL_NEF_EXPLOIT_REFERENCE_COUNTING if(reference_counted) { if(!OnSideMapRC.is_defined(&(v->point().hw()))) - switch(depth%3) { + switch(coord) { case 0: cr = CGAL::compare_x(v->point(), pop); OnSideMapRC[&(v->point().hw())] = cr == LARGER ? ON_POSITIVE_SIDE : @@ -363,13 +213,14 @@ Side_of_plane::operator() return OnSideMapRC[&(v->point().hw())]; } else { #endif - ComparePoints_ compare(depth%3); - if(!OnSideMap.is_defined(v)) { - cr = compare(v->point(), pop); - OnSideMap[v] = cr == LARGER ? ON_POSITIVE_SIDE : - cr == SMALLER ? ON_NEGATIVE_SIDE : ON_ORIENTED_BOUNDARY; - } - return OnSideMap[v]; + Oriented_side& side = OnSideMap[v]; + if(side == unknown_side) { + Compare compare(coord); + cr = compare(v->point(), pop); + side = cr == LARGER ? ON_POSITIVE_SIDE : + cr == SMALLER ? ON_NEGATIVE_SIDE : ON_ORIENTED_BOUNDARY; + } + return side; #ifdef CGAL_NEF_EXPLOIT_REFERENCE_COUNTING } #endif @@ -383,57 +234,13 @@ Side_of_plane::operator() */ template -template Oriented_side -Side_of_plane::operator() -( const Point_3& pop, Halfedge_handle e, Depth depth) { +Side_of_plane::operator()(Halfedge_handle e) { Vertex_handle v = e->source(); Vertex_handle vt = e->twin()->source(); - /* - Comparison_result cr; - if(!OnSideMap.is_defined(v)) - switch(depth%3) { - case 0: - cr = CGAL::compare_x(v->point(), pop); - OnSideMap[v] = cr == LARGER ? ON_POSITIVE_SIDE : - cr == SMALLER ? ON_NEGATIVE_SIDE : ON_ORIENTED_BOUNDARY; - break; - case 1: - cr = CGAL::compare_y(v->point(), pop); - OnSideMap[v] = cr == LARGER ? ON_POSITIVE_SIDE : - cr == SMALLER ? ON_NEGATIVE_SIDE : ON_ORIENTED_BOUNDARY; - break; - case 2: - cr = CGAL::compare_z(v->point(), pop); - OnSideMap[v] = cr == LARGER ? ON_POSITIVE_SIDE : - cr == SMALLER ? ON_NEGATIVE_SIDE : ON_ORIENTED_BOUNDARY; - break; - default: CGAL_error_msg( "wrong value"); - } - if(!OnSideMap.is_defined(vt)) - switch(depth%3) { - case 0: - cr = CGAL::compare_x(vt->point(), pop); - OnSideMap[vt] = cr == LARGER ? ON_POSITIVE_SIDE : - cr == SMALLER ? ON_NEGATIVE_SIDE : ON_ORIENTED_BOUNDARY; - break; - case 1: - cr = CGAL::compare_y(vt->point(), pop); - OnSideMap[vt] = cr == LARGER ? ON_POSITIVE_SIDE : - cr == SMALLER ? ON_NEGATIVE_SIDE : ON_ORIENTED_BOUNDARY; - break; - case 2: - cr = CGAL::compare_z(vt->point(), pop); - OnSideMap[vt] = cr == LARGER ? ON_POSITIVE_SIDE : - cr == SMALLER ? ON_NEGATIVE_SIDE : ON_ORIENTED_BOUNDARY; - break; - default: CGAL_error_msg( "wrong value"); - } - Oriented_side src_side = OnSideMap[v]; - Oriented_side tgt_side = OnSideMap[vt]; -*/ - Oriented_side src_side = (*this) (pop, v, depth); - Oriented_side tgt_side = (*this) (pop, vt, depth); + + Oriented_side src_side = (*this) (v); + Oriented_side tgt_side = (*this) (vt); if( src_side == tgt_side) return src_side; if( src_side == ON_ORIENTED_BOUNDARY) @@ -457,37 +264,10 @@ Side_of_plane::operator() template -template Oriented_side -Side_of_plane::operator() - (const Point_3& pop, Halffacet_handle f, Depth depth) { +Side_of_plane::operator()(Halffacet_handle f) { CGAL_assertion( std::distance( f->facet_cycles_begin(), f->facet_cycles_end()) > 0); - /* -#ifdef CGAL_NEF3_FACET_WITH_BOX - switch(depth%3) { - case 0: - if(f->b.min_coord(0) > pop.x()) - return ON_POSITIVE_SIDE; - if(f->b.max_coord(0) < pop.x()) - return ON_NEGATIVE_SIDE; - break; - case 1: - if(f->b.min_coord(1) > pop.y()) - return ON_POSITIVE_SIDE; - if(f->b.max_coord(1) < pop.y()) - return ON_NEGATIVE_SIDE; - break; - case 2: - if(f->b.min_coord(2) > pop.z()) - return ON_POSITIVE_SIDE; - if(f->b.max_coord(2) < pop.z()) - return ON_NEGATIVE_SIDE; - break; - default: CGAL_error_msg( "wrong value"); - } - return ON_ORIENTED_BOUNDARY; -#else - */ + Halffacet_cycle_iterator fc(f->facet_cycles_begin()); SHalfedge_handle e; CGAL_assertion(fc.is_shalfedge()); @@ -499,7 +279,7 @@ Side_of_plane::operator() Vertex_handle v; do { v = sc->source()->center_vertex(); - facet_side = (*this) (pop, v, depth); + facet_side = (*this) (v); ++sc; } while( facet_side == ON_ORIENTED_BOUNDARY && sc != send); @@ -509,7 +289,7 @@ Side_of_plane::operator() Oriented_side point_side; while( sc != send) { v = sc->source()->center_vertex(); - point_side = (*this) (pop, v, depth); + point_side = (*this) (v); ++sc; if( point_side == ON_ORIENTED_BOUNDARY) continue; diff --git a/Nef_3/include/CGAL/Nef_3/SNC_point_locator.h b/Nef_3/include/CGAL/Nef_3/SNC_point_locator.h index 9a91e4259d8..0cb23307d66 100644 --- a/Nef_3/include/CGAL/Nef_3/SNC_point_locator.h +++ b/Nef_3/include/CGAL/Nef_3/SNC_point_locator.h @@ -190,37 +190,13 @@ public: virtual void initialize(SNC_structure* W) { -#ifdef CGAL_NEF_LIST_OF_TRIANGLES - this->set_snc(*W); - candidate_provider = new SNC_candidate_provider(W); -#else // CGAL_NEF_LIST_OF_TRIANGLES - CGAL_NEF_TIMER(ct_t.start()); - this->version_ = std::string("Point Locator by Spatial Subdivision (tm)"); - CGAL_NEF_CLOG(version()); - CGAL_assertion( W != nullptr); -// (Base) *this = SNC_decorator(*W); - this->set_snc(*W); - Object_list objects; - Vertex_iterator v; - Halfedge_iterator e; - Halffacet_iterator f; - CGAL_forall_vertices( v, *this->sncp()) - objects.push_back(make_object(Vertex_handle(v))); - typename Object_list::difference_type v_end = objects.size(); - CGAL_forall_edges( e, *this->sncp()) - objects.push_back(make_object(Halfedge_handle(e))); - CGAL_forall_facets( f, *this->sncp()) { - objects.push_back(make_object(Halffacet_handle(f))); - } + if(initialized) delete candidate_provider; - Object_list_iterator oli=objects.begin()+v_end; - candidate_provider = new SNC_candidate_provider(objects,oli); + this->set_snc(*W); + candidate_provider = new SNC_candidate_provider(W); - // CGAL_NEF_TRACEN(*candidate_provider); - CGAL_NEF_TIMER(ct_t.stop()); -#endif // CGAL_NEF_LIST_OF_TRIANGLES initialized = true; }