// Copyright (c) 1997-2002 Max-Planck-Institute Saarbruecken (Germany). // All rights reserved. // // This file is part of CGAL (www.cgal.org). // // $URL: https://github.com/CGAL/cgal/blob/v5.2/Nef_3/include/CGAL/Nef_3/SNC_external_structure.h $ // $Id: SNC_external_structure.h 0779373 2020-03-26T13:31:46+01:00 Sébastien Loriot // SPDX-License-Identifier: GPL-3.0-or-later OR LicenseRef-Commercial // // // Author(s) : Michael Seel // Miguel Granados // Susan Hert // Lutz Kettner // Peter Hachenberger #ifndef CGAL_SNC_EXTERNAL_STRUCTURE_H #define CGAL_SNC_EXTERNAL_STRUCTURE_H #include #include #include #include #include #include #include #include #include #include #include #include #undef CGAL_NEF_DEBUG #define CGAL_NEF_DEBUG 43 #include #include namespace CGAL { struct int_lt { bool operator()(const int& i1, const int& i2) const { return i1 struct Halfedge_key_lt4 { bool operator()(const Edge_handle& e1, const Edge_handle& e2) const { if(CGAL::sign(e1->point().x()) != 0) { if(e1->source() != e2->source()) return CGAL::compare_x(e1->source()->point(), e2->source()->point()) < 0; else return e1->point().x() < 0; } if(CGAL::sign(e1->point().y()) != 0) { if(e1->source() != e2->source()) return CGAL::compare_y(e1->source()->point(), e2->source()->point()) < 0; else return e1->point().y() < 0; } if(e1->source() != e2->source()) return CGAL::compare_z(e1->source()->point(), e2->source()->point()) < 0; return e1->point().z() < 0; } }; template struct Halfedge_key_lt3 { bool operator()(const Edge_handle& e1, const Edge_handle& e2) const { if(e1->source() != e2->source()) return CGAL::lexicographically_xyz_smaller(e1->source()->point(), e2->source()->point()); if(CGAL::sign(e1->point().x()) != 0) return e1->point().x() < 0; if(CGAL::sign(e1->point().y()) != 0) return e1->point().y() < 0; return e1->point().z() < 0; } }; template struct Halfedge_key { typedef Halfedge_key Self; Point p; int i; Edge e; Halfedge_key(Point pi, int ii, Edge ei) : p(pi), i(ii), e(ei) {} Halfedge_key(const Self& k) : p(k.p), i(k.i), e(k.e) {} Self& operator=(const Self& k) { p=k.p; i=k.i; e=k.e; return *this; } bool operator==(const Self& k) const { return p==k.p && i==k.i; } bool operator!=(const Self& k) const { return !operator==(k); } }; template struct Halfedge_key_lt { typedef Halfedge_key Key; typedef typename Point::R R; typedef typename R::Vector_3 Vector; typedef typename R::Direction_3 Direction; bool operator()( const Key& k1, const Key& k2) const { if( k1.e->source() == k2.e->source()) return (k1.i < k2.i); Direction l(k1.e->vector()); if( k1.i < 0) l = -l; return (Direction( k2.p - k1.p) == l); } }; template std::ostream& operator<<(std::ostream& os, const Halfedge_key& k ) { os << k.p << " " << k.i; return os; } template int sign_of(const CGAL::Plane_3& h) { if ( h.c() != 0 ) return CGAL_NTS sign(h.c()); if ( h.b() != 0 ) return CGAL_NTS sign(-h.b()); return CGAL_NTS sign(h.a()); } struct Plane_lt { template bool operator()(const CGAL::Plane_3& h1, const CGAL::Plane_3& h2) const { typedef typename R::RT RT; typedef typename R::FT FT; bool flag=false; FT ratioa,ratiob,ratioc,ratiod; FT a1,b1,c1,d1,a2,b2,c2,d2; a1=h1.a(); a2=h2.a(); b1=h1.b(); b2=h2.b(); c1=h1.c(); c2=h2.c(); d1=h1.d(); d2=h2.d(); if(a2==0 || a1==0) { if(a2==a1) ratioa=1; else flag=true; } else{ratioa=a1/a2;} if(b2==0 || b1==0) { if(b2==b1) ratiob=ratioa; else flag=true; } else{ ratiob=b1/b2;} if(c2==0 || c1==0) { if(c2==c1) ratioc=ratioa; else flag=true; } else{ ratioc=c1/c2;} if(d2==0 || d1==0) { if(d2==d1) ratiod=ratioc; else flag=true; } else{ ratiod=d1/d2;} if(flag || !(ratioa==ratiob && ratioc==ratiod && ratioa==ratioc)) { RT diff = h1.a()-h2.a(); if ( (diff) != 0 ) return CGAL_NTS sign(diff) < 0; diff = h1.b()-h2.b(); if ( (diff) != 0 ) return CGAL_NTS sign(diff) < 0; diff = h1.c()-h2.c(); if ( (diff) != 0 ) return CGAL_NTS sign(diff) < 0; diff = h1.d()-h2.d(); if ( (diff) != 0 ) return CGAL_NTS sign(diff) < 0; } return 0; } }; struct Plane_RT_lt { template bool operator()(const CGAL::Plane_3& h1, const CGAL::Plane_3& h2) const { return (&(h1.a())) < (&(h2.a())); } }; // ---------------------------------------------------------------------------- // SNC_external_structure // ---------------------------------------------------------------------------- template class SNC_external_structure_base : public SNC_decorator { public: typedef SNC_structure_ SNC_structure; typedef typename SNC_structure::Infi_box Infi_box; typedef typename Infi_box::Standard_kernel Standard_kernel; typedef typename Standard_kernel::Point_3 Standard_point_3; typedef typename Infi_box::NT NT; typedef CGAL::SNC_decorator SNC_decorator; typedef CGAL::SNC_point_locator SNC_point_locator; typedef CGAL::SNC_FM_decorator FM_decorator; typedef CGAL::SNC_simplify SNC_simplify; typedef typename SNC_structure::Sphere_map Sphere_map; typedef CGAL::SM_decorator SM_decorator; typedef CGAL::SM_const_decorator SM_const_decorator; typedef CGAL::SM_point_locator SM_point_locator; typedef typename SNC_structure::Halfedge_iterator Halfedge_iterator; typedef typename SNC_structure::Halffacet_iterator Halffacet_iterator; typedef typename SNC_structure::Volume_iterator Volume_iterator; typedef typename SNC_structure::Vertex_handle Vertex_handle; typedef typename SNC_structure::Halfedge_handle Halfedge_handle; typedef typename SNC_structure::Halffacet_handle Halffacet_handle; typedef typename SNC_structure::Volume_handle Volume_handle; typedef typename SNC_structure::Halffacet_const_handle Halffacet_const_handle; typedef typename SNC_structure::SFace_const_handle SFace_const_handle; typedef typename SNC_structure::SHalfedge_iterator SHalfedge_iterator; typedef typename SNC_structure::SFace_iterator SFace_iterator; typedef typename SNC_structure::SHalfloop_iterator SHalfloop_iterator; typedef typename SNC_structure::SVertex_handle SVertex_handle; typedef typename SNC_structure::SHalfedge_handle SHalfedge_handle; typedef typename SNC_structure::SFace_handle SFace_handle; typedef typename SNC_structure::SHalfloop_handle SHalfloop_handle; typedef typename SNC_structure::Object_handle Object_handle; typedef typename SNC_structure::SHalfedge_around_facet_circulator SHalfedge_around_facet_circulator; typedef typename SNC_structure::Point_3 Point_3; typedef typename SNC_structure::Direction_3 Direction_3; typedef typename SNC_structure::Plane_3 Plane_3; typedef typename SNC_structure::Ray_3 Ray_3; typedef typename SNC_structure::Sphere_point Sphere_point; typedef typename SNC_structure::Sphere_circle Sphere_circle; typedef typename SM_decorator::SHalfedge_around_svertex_circulator SHalfedge_around_svertex_circulator; SNC_point_locator* pl; typedef CGAL::Unique_hash_map Sface_shell_hash; typedef CGAL::Unique_hash_map Face_shell_hash; typedef CGAL::Unique_hash_map SFace_visited_hash; typedef CGAL::Unique_hash_map Shell_closed_hash; using SNC_decorator::visit_shell_objects; using SNC_decorator::link_as_inner_shell; using SNC_decorator::link_as_outer_shell; using SNC_decorator::link_as_prev_next_pair; using SNC_decorator::get_visible_facet; using SNC_decorator::adjacent_sface; using SNC_decorator::make_twins; struct Shell_explorer { const SNC_decorator& D; Sface_shell_hash& ShellSf; Face_shell_hash& ShellF; // Shell_closed_hash& Closed; SFace_visited_hash& Done; SFace_handle sf_min; int n; Shell_explorer(const SNC_decorator& Di, Sface_shell_hash& SSf, Face_shell_hash& SF, SFace_visited_hash& Vi) : D(Di), ShellSf(SSf), ShellF(SF), Done(Vi), n(0) {} void visit(SFace_handle h) { CGAL_NEF_TRACEN("visit sf "<center_vertex()->point()); ShellSf[h]=n; Done[h]=true; if ( CGAL::lexicographically_xyz_smaller(h->center_vertex()->point(), sf_min->center_vertex()->point())) sf_min = h; } void visit(Vertex_handle h) { CGAL_USE(h); CGAL_NEF_TRACEN("visit v "<point()); } void visit(Halfedge_handle h) { CGAL_USE(h); CGAL_NEF_TRACEN("visit he "<< h->source()->point()); } void visit(Halffacet_handle h) { CGAL_NEF_TRACEN(h->plane()); ShellF[h]=n; } void visit(SHalfedge_handle ) {} void visit(SHalfloop_handle ) {} SFace_handle& minimal_sface() { return sf_min; } void increment_shell_number() { CGAL_NEF_TRACEN("leaving shell "<>>>>pair_up_halfedges"); typedef Halfedge_key_lt3 Halfedge_key_lt; typedef std::list Halfedge_list; typedef typename Standard_kernel::Kernel_tag Kernel_tag; typedef CGAL::Pluecker_line_3 Pluecker_line_3; typedef CGAL::Pluecker_line_lt Pluecker_line_lt; typedef std::map< Pluecker_line_3, Halfedge_list, Pluecker_line_lt> Pluecker_line_map; Pluecker_line_map M; Pluecker_line_map M2; Pluecker_line_map M3; Pluecker_line_map M4; NT eval(Infi_box::compute_evaluation_constant_for_halfedge_pairup(*this->sncp())); Halfedge_iterator e; CGAL_forall_halfedges(e,*this->sncp()) { // progress++; Point_3 p = e->source()->point(); Point_3 q = p + e->vector(); CGAL_NEF_TRACE(" segment("<vector()<<")"); Standard_point_3 sp = Infi_box::standard_point(p,eval); Standard_point_3 sq = Infi_box::standard_point(q,eval); Pluecker_line_3 l( sp, sq); int inverted; l = categorize( l, inverted); if(Infi_box::is_edge_on_infibox(e)) if(Infi_box::is_type4(e)) M4[l].push_back(e); else if(Infi_box::is_type3(e)) M3[l].push_back(e); else M2[l].push_back(e); else M[l].push_back(e); // the following trace crashes when compiling with optimizations (-O{n}) //CGAL_NEF_TRACEN(Infi_box::standard_point(point(vertex(e)))+ CGAL_NEF_TRACEN(" line("<first); typename Halfedge_list::iterator itl; CGAL_forall_iterators(itl,it->second) { Halfedge_handle e1 = *itl; ++itl; CGAL_assertion(itl != it->second.end()); Halfedge_handle e2 = *itl; while(normalized(e1->vector()) != normalized(-e2->vector())) { ++itl; make_twins(e1,*itl); e1 = e2; ++itl; e2 = *itl; } CGAL_NEF_TRACEN(" " << e1->source()->point() << " -> " << e2->source()->point()); CGAL_NEF_TRACEN(e1->vector()<<" -> "<<-e2->vector()); make_twins(e1,e2); CGAL_assertion(e1->mark()==e2->mark()); // discard temporary sphere_point ? } } CGAL_forall_iterators(it,M3) { // progress++; it->second.sort(Halfedge_key_lt()); CGAL_NEF_TRACEN("search opposite "<first); typename Halfedge_list::iterator itl; CGAL_forall_iterators(itl,it->second) { Halfedge_handle e1 = *itl; ++itl; CGAL_assertion(itl != it->second.end()); Halfedge_handle e2 = *itl; while(normalized(e1->vector()) != normalized(-e2->vector())) { ++itl; make_twins(e1,*itl); e1 = e2; ++itl; e2 = *itl; } CGAL_NEF_TRACEN(" " << e1->source()->point() << " -> " << e2->source()->point()); CGAL_NEF_TRACEN(e1->vector()<<" -> "<<-e2->vector()); make_twins(e1,e2); CGAL_assertion(e1->mark()==e2->mark()); // discard temporary sphere_point ? } } CGAL_forall_iterators(it,M2) { // progress++; it->second.sort(Halfedge_key_lt()); CGAL_NEF_TRACEN("search opposite "<first); typename Halfedge_list::iterator itl; CGAL_forall_iterators(itl,it->second) { Halfedge_handle e1 = *itl; ++itl; CGAL_assertion(itl != it->second.end()); Halfedge_handle e2 = *itl; while(normalized(e1->vector()) != normalized(-e2->vector())) { ++itl; this->make_twins(e1,*itl); e1 = e2; ++itl; e2 = *itl; } CGAL_NEF_TRACEN(" " << e1->source()->point() << " -> " << e2->source()->point()); CGAL_NEF_TRACEN(e1->vector()<<" -> "<<-e2->vector()); this->make_twins(e1,e2); CGAL_assertion(e1->mark()==e2->mark()); // discard temporary sphere_point ? } } CGAL_forall_iterators(it,M) { // progress++; it->second.sort(Halfedge_key_lt()); CGAL_NEF_TRACEN("search opposite "<first); typename Halfedge_list::iterator itl; CGAL_forall_iterators(itl,it->second) { Halfedge_handle e1 = *itl; ++itl; CGAL_assertion(itl != it->second.end()); Halfedge_handle e2 = *itl; while(normalized(e1->vector()) != normalized(-e2->vector())) { ++itl; this->make_twins(e1,*itl); e1 = e2; ++itl; e2 = *itl; } CGAL_NEF_TRACEN(" " << e1->source()->point() << " -> " << e2->source()->point()); CGAL_NEF_TRACEN(e1->vector()<<" -> "<<-e2->vector()); this->make_twins(e1,e2); CGAL_assertion(e1->mark()==e2->mark()); // discard temporary sphere_point ? } } } #else void pair_up_halfedges() const { /*{\Mop pairs all halfedge stubs to create the edges in 3-space.}*/ // CGAL_NEF_SETDTHREAD(43*61); CGAL_NEF_TRACEN(">>>>>pair_up_halfedges"); typedef Halfedge_key< Point_3, Halfedge_handle> Halfedge_key; typedef Halfedge_key_lt< Point_3, Halfedge_handle, SNC_decorator> Halfedge_key_lt; typedef std::list Halfedge_list; typedef typename Standard_kernel::Kernel_tag Kernel_tag; typedef CGAL::Pluecker_line_3 Pluecker_line_3; typedef CGAL::Pluecker_line_lt Pluecker_line_lt; typedef std::map< Pluecker_line_3, Halfedge_list, Pluecker_line_lt> Pluecker_line_map; SNC_decorator D(*this); Pluecker_line_map M; Pluecker_line_map M2; Pluecker_line_map M3; Pluecker_line_map M4; NT eval(Infi_box::compute_evaluation_constant_for_halfedge_pairup(*this->sncp()));; Halfedge_iterator e; CGAL_forall_halfedges(e,*this->sncp()) { // progress++; Point_3 p = e->source()->point(); Point_3 q = p + e->vector(); CGAL_NEF_TRACE(" segment("<vector()<<")"); Standard_point_3 sp = Infi_box::standard_point(p,eval); Standard_point_3 sq = Infi_box::standard_point(q,eval); Pluecker_line_3 l( sp, sq); int inverted; l = categorize( l, inverted); if(Infi_box::is_edge_on_infibox(e)) if(Infi_box::is_type4(e)) M4[l].push_back(Halfedge_key(p,inverted,e)); else if(Infi_box::is_type3(e)) M3[l].push_back(Halfedge_key(p,inverted,e)); else M2[l].push_back(Halfedge_key(p,inverted,e)); else M[l].push_back(Halfedge_key(p,inverted,e)); // the following trace crashes when compiling with optimizations (-O{n}) //CGAL_NEF_TRACEN(Infi_box::standard_point(point(vertex(e)))+ CGAL_NEF_TRACEN(" line("<first); typename Halfedge_list::iterator itl; CGAL_forall_iterators(itl,it->second) { Halfedge_handle e1 = itl->e; ++itl; CGAL_assertion(itl != it->second.end()); Halfedge_handle e2 = itl->e; CGAL_NEF_TRACEN(" " << e1->source()->point() << " -> " << e2->source()->point()); CGAL_NEF_TRACEN(e1->vector()<<" -> "<<-e2->vector()); make_twins(e1,e2); CGAL_assertion(e1->mark()==e2->mark()); // discard temporary sphere_point ? } } CGAL_forall_iterators(it,M3) { // progress++; it->second.sort(Halfedge_key_lt()); CGAL_NEF_TRACEN("search opposite "<first); typename Halfedge_list::iterator itl; CGAL_forall_iterators(itl,it->second) { Halfedge_handle e1 = itl->e; ++itl; CGAL_assertion(itl != it->second.end()); Halfedge_handle e2 = itl->e; CGAL_NEF_TRACEN(" " << e1->source()->point() << " -> " << e2->source()->point()); CGAL_NEF_TRACEN(e1->vector()<<" -> "<<-e2->vector()); make_twins(e1,e2); CGAL_assertion(e1->mark()==e2->mark()); // discard temporary sphere_point ? } } CGAL_forall_iterators(it,M2) { // progress++; it->second.sort(Halfedge_key_lt()); CGAL_NEF_TRACEN("search opposite "<first); typename Halfedge_list::iterator itl; CGAL_forall_iterators(itl,it->second) { Halfedge_handle e1 = itl->e; ++itl; CGAL_assertion(itl != it->second.end()); Halfedge_handle e2 = itl->e; CGAL_NEF_TRACEN(" " << e1->source()->point() << " -> " << e2->source()->point()); CGAL_NEF_TRACEN(e1->vector()<<" -> "<<-e2->vector()); make_twins(e1,e2); CGAL_assertion(e1->mark()==e2->mark()); // discard temporary sphere_point ? } } CGAL_forall_iterators(it,M) { // progress++; it->second.sort(Halfedge_key_lt()); CGAL_NEF_TRACEN("search opposite "<first); typename Halfedge_list::iterator itl; CGAL_forall_iterators(itl,it->second) { Halfedge_handle e1 = itl->e; ++itl; CGAL_assertion(itl != it->second.end()); Halfedge_handle e2 = itl->e; CGAL_NEF_TRACEN(" " << e1->source()->point() << " -> " << e2->source()->point()); CGAL_NEF_TRACEN(e1->vector()<<" -> "<< -e2->vector()); CGAL_assertion(e1->source()->point() != e2->source()->point()); CGAL_assertion(e1->mark()==e2->mark()); make_twins(e1,e2); // discard temporary sphere_point ? } } } #endif void link_shalfedges_to_facet_cycles() const { /*{\Mop creates all non-trivial facet cycles from sedges. \precond |pair_up_halfedges()| was called before.}*/ // CGAL_NEF_SETDTHREAD(43*31); CGAL_NEF_TRACEN(">>>>>link_shalfedges_to_facet_cycles"); #ifdef CGAL_NEF_EXPLOIT_REFERENCE_COUNTING Point_3 p1(1,2,7), p2(p1); bool reference_counted = (&(p1.hx()) == &(p2.hx())); #endif Halfedge_iterator e; CGAL_forall_edges(e,*this->sncp()) { // progress++; CGAL_NEF_TRACEN(""); CGAL_NEF_TRACEN(PH(e)); Halfedge_iterator et = e->twin(); SM_decorator D(&*e->source()), Dt(&*et->source()); CGAL_NEF_TRACEN(e->source()->point()); if ( D.is_isolated(e) ) continue; SHalfedge_around_svertex_circulator ce(D.first_out_edge(e)),cee(ce); SHalfedge_around_svertex_circulator cet(Dt.first_out_edge(et)),cete(cet); #ifdef CGAL_NEF_EXPLOIT_REFERENCE_COUNTING if(reference_counted) { CGAL_For_all(cet,cete) if ( &(cet->circle().a()) == &(ce->circle().opposite().a()) && cet->source()->twin() == ce->source() ) break; } else #endif CGAL_For_all(cet,cete) if ( cet->circle() == ce->circle().opposite() && cet->source()->twin() == ce->source() ) break; #ifndef NDEBUG if( cet->circle() != ce->circle().opposite() ) CGAL_NEF_TRACEN("assertion failed!"); CGAL_NEF_TRACEN("vertices " << e->source()->point() << " " << et->source()->point()); SHalfedge_around_svertex_circulator sc(D.first_out_edge(e)); SHalfedge_around_svertex_circulator sct(Dt.first_out_edge(et)); CGAL_NEF_TRACEN(""); CGAL_For_all(sc,cee) CGAL_NEF_TRACEN("sseg@E addr="<<&*sc<< " src="<< sc->source()->point()<< " tgt="<< sc->target()->point()<circle()) == normalized(ce->circle().opposite()) ); CGAL_assertion( cet->source()->twin() == ce->source()); CGAL_For_all(ce,cee) { CGAL_NEF_TRACEN("circles " << normalized(cet->circle()) << " " << normalized(ce->circle()) << " sources " << cet->target()->point() << " " << ce->target()->point()); CGAL_assertion( normalized(cet->circle()) == normalized(ce->circle().opposite())); CGAL_assertion( cet->source()->twin() == ce->source()); CGAL_assertion(ce->mark()==cet->mark()); link_as_prev_next_pair(cet->twin(),ce); link_as_prev_next_pair(ce->twin(),cet); --cet; // ce moves ccw, cet moves cw } } } void categorize_facet_cycles_and_create_facets() const { /*{\Mop collects all facet cycles incident to a facet and creates the facets. \precond |link_shalfedges_to_facet_cycles()| was called before.}*/ // CGAL_NEF_SETDTHREAD(43*31); CGAL_NEF_TRACEN(">>>>>categorize_facet_cycles_and_create_facets"); typedef std::list Object_list; #ifdef CGAL_NEF_EXPLOIT_REFERENCE_COUNTING typedef std::map Map_planes; #else typedef std::map Map_planes; #endif Map_planes M; SHalfedge_iterator e; CGAL_forall_shalfedges(e,*this->sncp()) { Sphere_circle c(e->circle()); Plane_3 h = c.plane_through(e->source()->source()->point()); CGAL_NEF_TRACEN("\n" << e->source()->twin()->source()->point() <<" - " << e->source()->source()->point() <<" - "<< e->twin()->source()->twin()->source()->point() << " has plane " << h << " has circle " << e->circle() << " has signum " << sign_of(h)); if ( sign_of(h)<0 ) continue; M[normalized(h)].push_back(make_object(e->twin())); CGAL_NEF_TRACEN(" normalized as " << normalized(h)); /* Unique_hash_map Done(false); SHalfedge_iterator ei; CGAL_forall_sedges(ei,*this->sncp()) { if(Done[ei]) continue; Sphere_circle c(ei->circle()); Plane_3 h = c.plane_through(ei->source()->source()->point()); CGAL_NEF_TRACEN("\n" << ei->source()->twin()->source()->point() <<" - " << ei->source()->source()->point() <<" - "<< ei->twin()->source()->twin()->source()->point() << " has plane " << h << " has circle " << ei->circle() << " has signum " << sign_of(h)); SHalfedge_handle e(ei); if ( sign_of(h)<0 ) { h = h.opposite(); e = e->twin(); } SHalfedge_around_facet_circulator sfc(e), send(sfc); CGAL_For_all(sfc, send) { M[normalized(h)].push_back(make_object(e->twin())); Done[sfc] = true; Done[sfc->twin()] = true; CGAL_NEF_TRACEN(" normalized as " << normalized(h)); } */ } SHalfloop_iterator l; CGAL_forall_shalfloops(l,*this->sncp()) { Sphere_circle c(l->circle()); Plane_3 h = c.plane_through(l->incident_sface()->center_vertex()->point()); if ( sign_of(h)<0 ) continue; // CGAL_assertion( h == normalized(h)); M[normalized(h)].push_back(make_object(l->twin())); } #ifdef CGAL_NEF3_TIMER_PLANE_SWEEPS number_of_plane_sweeps=0; timer_plane_sweeps.reset(); #endif typename Map_planes::iterator it; CGAL_forall_iterators(it,M) { // progress2++; // CGAL_NEF_TRACEN(" plane "<first<<" "<<(it->first).point()); FM_decorator D(*this->sncp()); D.create_facet_objects(it->first,it->second.begin(),it->second.end()); } // CGAL_NEF_SETDTHREAD(1); } void create_volumes() { /*{\Mop collects all shells incident to a volume and creates the volumes. \precond |categorize_facet_cycles_and_creating_facets()| was called before.}*/ #ifdef CGAL_NEF3_TIMER_POINT_LOCATION number_of_ray_shooting_queries=0; timer_ray_shooting.reset(); #endif // CGAL_NEF_SETDTHREAD(37*43*503*509); CGAL_NEF_TRACEN(">>>>>create_volumes"); Sface_shell_hash ShellSf(0); Face_shell_hash ShellF(0); SFace_visited_hash Done(false); Shell_explorer V(*this,ShellSf,ShellF,Done); std::vector MinimalSFace; std::vector EntrySFace; std::vector Closed; SFace_iterator f; // First, we classify all the Shere Faces per Shell. For each Shell we // determine its minimum lexicographyly vertex and we check wheter the // Shell encloses a region (closed surface) or not. CGAL_forall_sfaces(f,*this->sncp()) { // progress++; CGAL_NEF_TRACEN("sface in " << ShellSf[f]); if ( Done[f] ) continue; V.minimal_sface() = f; visit_shell_objects(f,V); CGAL_NEF_TRACEN("minimal vertex " << V.minimal_sface()->center_vertex()->point()); MinimalSFace.push_back(V.minimal_sface()); EntrySFace.push_back(f); V.increment_shell_number(); CGAL_NEF_TRACEN("sface out " << ShellSf[f]); } for(unsigned int i=0; itwin()]) { Closed[ShellF[hf]] = true; Closed[ShellF[hf->twin()]] = true; } CGAL_assertion( pl != nullptr); #ifdef CGAL_NEF3_TIMER_INITIALIZE_KDTREE CGAL::Timer timer_initialize_kdtree; timer_initialize_kdtree.start(); #endif pl->initialize(this->sncp()); // construct the point locator #ifdef CGAL_NEF3_TIMER_INITIALIZE_KDTREE timer_initialize_kdtree.stop(); if(cgal_nef3_timer_on) std::cout << "Runtime_initialize_kdtree: " << timer_initialize_kdtree.time() << std::endl; #endif // then, we determine the Shells which correspond to Volumes via a ray // shootting in the direction (-1,0,0) over the Sphere_map of the minimal // vertex. The Shell corresponds to a Volume if the object hit belongs // to another Shell. this->sncp()->new_volume(); // outermost volume (nirvana) if(MinimalSFace.size() == 0) return; Vertex_handle v_min = MinimalSFace[0]->center_vertex(); for( unsigned int i = 0; i < MinimalSFace.size(); ++i) { // progress2++; Vertex_handle v = MinimalSFace[i]->center_vertex(); if(CGAL::lexicographically_xyz_smaller(v->point(),v_min->point())) v_min=v; CGAL_NEF_TRACEN( "Shell #" << i << " minimal vertex: " << v->point()); SM_point_locator D((Sphere_map*) &*v); Object_handle o = D.locate(Sphere_point(-1,0,0)); SFace_const_handle sfc; if( !CGAL::assign(sfc, o) || ShellSf[sfc] != i) { CGAL_NEF_TRACEN("the shell encloses a volume"); CGAL_NEF_TRACEN("sface hit? "<sncp()->new_volume(); c->mark() = f->mark(); // TODO test if line is redundant link_as_inner_shell(f, c ); CGAL_NEF_TRACE( "Shell #" << i <<" linked as inner shell"); CGAL_NEF_TRACEN( "(sface" << (CGAL::assign(sfc,o)?"":" not") << " hit case)"); } } } // finaly, we go through all the Shells which do not correspond to a Volume // and we assign them to its enclosing Volume determined via a facet below // check. CGAL_forall_sfaces(f,*this->sncp()) { // progress3++; if ( f->volume() != Volume_handle() ) continue; CGAL_NEF_TRACEN( "Outer shell #" << ShellSf[f] << " volume?"); Volume_handle c = determine_volume( MinimalSFace[ShellSf[f]], MinimalSFace, ShellSf ); c->mark() = f->mark(); link_as_outer_shell( f, c ); } } Halffacet_handle get_facet_below( Vertex_handle vi, const std::vector< SFace_handle>& MinimalSFace, const Sface_shell_hash& Shell) const { // {\Mop determines the facet below a vertex |vi| via ray shooting. } Halffacet_handle f_below; Point_3 p = vi->point(); if(!Infi_box::is_standard(p)) return Halffacet_handle(); Ray_3 ray = Ray_3(p, Direction_3(-1,0,0)); #ifdef CGAL_NEF3_TIMER_POINT_LOCATION number_of_ray_shooting_queries++; timer_ray_shooting.start(); #endif Object_handle o = pl->shoot(ray); #ifdef CGAL_NEF3_TIMER_POINT_LOCATION timer_ray_shooting.stop(); #endif // The ray here has an special property since it is shooted from the lowest // vertex in a shell, so it would be expected that the ray goes along the // interior of a volume before it hits a 2-skeleton element. // Unfortunatelly, it seems to be possible that several shells are incident // to this lowest vertex, and in consequence, the ray could also go along // an edge or a facet belonging to a different shell. // This fact invalidates the precondition of the get_visible_facet method, // (the ray must pierce the local view of the hit object in a sface). // This situation has not been analyzed and has to be verified. Granados. Vertex_handle v; Halfedge_handle e; Halffacet_handle f; CGAL_NEF_TRACEN("get_facet_below"); if( CGAL::assign(v, o)) { CGAL_NEF_TRACEN("facet below from from vertex..."); f_below = get_visible_facet(v, ray); if( f_below == Halffacet_handle()) { CGAL_assertion(v->sfaces_begin() == v->sfaces_last()); f_below = get_facet_below(MinimalSFace[Shell[v->sfaces_begin()]]->center_vertex(), MinimalSFace,Shell); } } else if( CGAL::assign(e, o)) { CGAL_NEF_TRACEN("facet below from from edge..."); f_below = get_visible_facet(e, ray); if( f_below == Halffacet_handle()) { CGAL_assertion(e->source()->sfaces_begin() == e->source()->sfaces_last()); f_below = get_facet_below(MinimalSFace[Shell[e->source()->sfaces_begin()]]->center_vertex(), MinimalSFace, Shell); } } else if( CGAL::assign(f, o)) { CGAL_NEF_TRACEN("facet below from from facet..."); f_below = get_visible_facet(f, ray); CGAL_assertion( f_below != Halffacet_handle()); } else { CGAL_NEF_TRACEN("no facet below found..."); } return f_below; } Volume_handle determine_volume( SFace_handle sf, const std::vector< SFace_handle>& MinimalSFace, const Sface_shell_hash& Shell ) const { //{\Mop determines the volume |C| that a shell |S| pointed by |sf| // belongs to. \precondition |S| separates the volume |C| from an enclosed // volume.} CGAL_NEF_TRACEN("determine volume"); Vertex_handle v_min = MinimalSFace[Shell[sf]]->center_vertex(); Halffacet_handle f_below = get_facet_below(v_min, MinimalSFace, Shell); if ( f_below == Halffacet_handle()) return SNC_decorator(*this).volumes_begin(); Volume_handle c = f_below->incident_volume(); if( c != Volume_handle()) { CGAL_NEF_TRACE( "Volume " << &*c << " hit "); CGAL_NEF_TRACEN("(Shell #" << Shell[adjacent_sface(f_below)] << ")"); return c; } SFace_handle sf_below = adjacent_sface(f_below); CGAL_NEF_TRACE( "Shell not assigned to a volume hit "); CGAL_NEF_TRACEN( "(Inner shell #" << Shell[sf_below] << ")"); c = determine_volume( sf_below, MinimalSFace, Shell); link_as_inner_shell( sf_below, c); return c; } void build_external_structure() { #ifdef CGAL_NEF3_TIMER_EXTERNAL_STRUCTURE CGAL::Timer timer_external_structure; timer_external_structure.start(); #endif #ifdef CGAL_NEF3_TIMER_PLUECKER CGAL::Timer timer_pluecker; timer_pluecker.start(); #endif // SNC_io_parser O0(std::cerr,*this->sncp()); // O0.print(); pair_up_halfedges(); #ifdef CGAL_NEF3_TIMER_PLUECKER timer_pluecker.stop(); if(cgal_nef3_timer_on) std::cout << "Runtime_pluecker: " << timer_pluecker.time() << std::endl; #endif link_shalfedges_to_facet_cycles(); // SNC_io_parser O0(std::cerr,*this->sncp()); // O0.print(); categorize_facet_cycles_and_create_facets(); create_volumes(); #ifdef CGAL_NEF3_TIMER_EXTERNAL_STRUCTURE timer_external_structure.stop(); if(cgal_nef3_timer_on) std::cout << "Runtime_external_structure: " << timer_external_structure.time() << std::endl; #endif } template void build_after_binary_operation(Association) { #ifdef CGAL_NEF3_TIMER_SIMPLIFICATION CGAL::Timer timer_simplification; timer_simplification.start(); #endif SNC_simplify simp(*this->sncp()); simp.vertex_simplification(SNC_simplify::NO_SNC); #ifdef CGAL_NEF3_TIMER_SIMPLIFICATION timer_simplification.stop(); if(cgal_nef3_timer_on) std::cout << "Runtime_simplification: " << timer_simplification.time() << std::endl; #endif CGAL_NEF_TRACEN("\nnumber of vertices (so far...) = " << this->sncp()->number_of_vertices()); CGAL_NEF_TRACEN("=> resultant vertices (after simplification): "); build_external_structure(); } void clear_external_structure() { this->sncp()->clear_snc_boundary(); while(this->sncp()->volumes_begin()!= this->sncp()->volumes_end()) this->sncp()->delete_volume(this->sncp()->volumes_begin()); while(this->sncp()->halffacets_begin()!= this->sncp()->halffacets_end()) this->sncp()->delete_halffacet_pair(this->sncp()->halffacets_begin()); SHalfedge_iterator se; CGAL_forall_shalfedges(se,*this) se->facet() = Halffacet_handle(); SFace_iterator sf; CGAL_forall_sfaces(sf,*this) sf->volume() = Volume_handle(); } }; template class SNC_external_structure : public SNC_external_structure_base { typedef CGAL::SNC_decorator SNC_decorator; typedef CGAL::SNC_point_locator SNC_point_locator; public: SNC_external_structure( SNC_structure_& W, SNC_point_locator* spl = nullptr) : SNC_external_structure_base(W, spl) {} }; template class SNC_external_structure : public SNC_external_structure_base { public: typedef SNC_structure_ SNC_structure; typedef SNC_external_structure_base Base; typedef CGAL::SNC_decorator SNC_decorator; typedef CGAL::SNC_point_locator SNC_point_locator; typedef CGAL::SNC_FM_decorator FM_decorator; typedef typename SNC_structure::Sphere_map Sphere_map; typedef CGAL::SM_decorator SM_decorator; typedef CGAL::SNC_simplify SNC_simplify; typedef typename SNC_structure::Halfedge_iterator Halfedge_iterator; typedef typename SNC_structure::SHalfedge_iterator SHalfedge_iterator; typedef typename SNC_structure::SHalfloop_iterator SHalfloop_iterator; typedef typename SNC_structure::Halfedge_handle Halfedge_handle; typedef typename SNC_structure::SHalfedge_handle SHalfedge_handle; typedef typename SNC_structure::Object_handle Object_handle; typedef typename SNC_structure::SHalfedge_around_facet_circulator SHalfedge_around_facet_circulator; typedef typename SM_decorator::SHalfedge_around_svertex_circulator SHalfedge_around_svertex_circulator; typedef typename SNC_structure::Plane_3 Plane_3; using Base::make_twins; using Base::link_as_prev_next_pair; using Base::link_as_inner_shell; SNC_external_structure( SNC_structure& W, SNC_point_locator* spl = nullptr) : Base(W, spl) {} /*{\Mcreate makes |\Mvar| a decorator of |W|.}*/ public: void pair_up_halfedges() const { typedef Halfedge_key_lt4 Halfedge_key_lt; typedef std::list Halfedge_list; typedef std::map index_map; CGAL_NEF_TRACEN("pair up by indexes"); index_map i2he; Halfedge_iterator ei; CGAL_forall_halfedges(ei, *this->sncp()) i2he[ei->get_index()].push_back(ei); typename index_map::iterator it; CGAL_forall_iterators(it,i2he) { CGAL_NEF_TRACEN("pair up " << it->first); it->second.sort(Halfedge_key_lt()); typename Halfedge_list::iterator itl; CGAL_forall_iterators(itl,it->second) { Halfedge_handle e1 = *itl; CGAL_NEF_TRACEN(e1->source()->point() << ", " << e1->vector()); ++itl; CGAL_assertion(itl != it->second.end()); Halfedge_handle e2 = *itl; CGAL_NEF_TRACEN(" + " << e2->source()->point() << ", " << e2->vector()); make_twins(e1,e2); // SNC_io_parser O0(std::cerr,*this->sncp()); // O0.print(); CGAL_assertion(e1->mark()==e2->mark()); } } } void link_shalfedges_to_facet_cycles() const { /*{\Mop creates all non-trivial facet cycles from sedges. \precond |pair_up_halfedges()| was called before.}*/ // CGAL_NEF_SETDTHREAD(43*31); CGAL_NEF_TRACEN(">>>>>link_shalfedges_to_facet_cycles"); Halfedge_iterator e; CGAL_forall_edges(e,*this->sncp()) { // progress++; CGAL_NEF_TRACEN(""); CGAL_NEF_TRACEN(PH(e)); Halfedge_iterator et = e->twin(); SM_decorator D(&*e->source()), Dt(&*et->source()); CGAL_NEF_TRACEN(e->source()->point()); if ( D.is_isolated(e) ) continue; SHalfedge_around_svertex_circulator ce(D.first_out_edge(e)),cee(ce); SHalfedge_around_svertex_circulator cet(Dt.first_out_edge(et)),cete(cet); CGAL_For_all(cet,cete) { // std::cerr << cet->get_index() << ", " << ce->twin()->get_index() << std::endl; if (cet->get_forward_index() == ce->twin()->get_backward_index()) // cet->source()->twin() == ce->source()) break; } CGAL_NEF_TRACEN("vertices " << e->source()->point() << " " << et->source()->point()); SHalfedge_around_svertex_circulator sc(D.first_out_edge(e)); SHalfedge_around_svertex_circulator sct(Dt.first_out_edge(et)); CGAL_NEF_TRACEN(""); CGAL_For_all(sc,cee) CGAL_NEF_TRACEN("sseg@E addr="<<&*sc<< " src="<< sc->source()->point()<< " tgt="<< sc->target()->point()<< std::endl << " circle=" << sc->circle()<< std::endl << " indexes=" << sc->get_forward_index() << "," << sc->get_backward_index() << std::endl << " " << sc->twin()->get_forward_index() << "," << sc->twin()->get_backward_index()); CGAL_NEF_TRACEN(""); CGAL_For_all(sct,cete) CGAL_NEF_TRACEN("sseg@ET addr="<<&*sct<< " src="<< sct->source()->point()<< " tgt="<target()->point() <get_backward_index() << std::endl << " " << sct->twin()->get_forward_index() << "," << sct->twin()->get_backward_index()); CGAL_NEF_TRACEN(""); CGAL_assertion( cet->get_index() == ce->twin()->get_index()); CGAL_assertion( cet->twin()->get_index() == ce->get_index()); CGAL_assertion( cet->source()->twin() == ce->source()); CGAL_For_all(ce,cee) { CGAL_NEF_TRACEN("circles " << cet->circle() << " " << ce->circle() << " sources " << cet->target()->point() << " " << ce->target()->point()); CGAL_assertion( cet->source()->twin() == ce->source()); CGAL_assertion(ce->mark()==cet->mark()); link_as_prev_next_pair(cet->twin(),ce); link_as_prev_next_pair(ce->twin(),cet); --cet; // ce moves ccw, cet moves cw } } } void categorize_facet_cycles_and_create_facets() const { /*{\Mop collects all facet cycles incident to a facet and creates the facets. \precond |link_shalfedges_to_facet_cycles()| was called before.}*/ // CGAL_NEF_SETDTHREAD(43*31); CGAL_NEF_TRACEN(">>>>>categorize_facet_cycles_and_create_facets"); typedef std::list Object_list; typedef std::map Map_planes; Map_planes M; SHalfedge_iterator e; CGAL_forall_shalfedges(e,*this->sncp()) { if(e->get_index() > e->twin()->get_index()) continue; M[e->get_index()].push_back(make_object(e)); } SHalfloop_iterator l; CGAL_forall_shalfloops(l,*this->sncp()) { if(l->get_index() > l->twin()->get_index()) continue; M[l->get_index()].push_back(make_object(l)); } #ifdef CGAL_NEF3_TIMER_PLANE_SWEEPS number_of_plane_sweeps=0; timer_plane_sweeps.reset(); #endif typename Map_planes::iterator it; CGAL_forall_iterators(it,M) { CGAL_NEF_TRACEN(" plane "<< it->first); CGAL_NEF_TRACEN(" size "<< it->second.size()); FM_decorator D(*this->sncp()); Plane_3 h; Object_handle o(*it->second.begin()); if(CGAL::assign(e, o)) h = e->circle().opposite().plane_through(e->source()->source()->point()); else if(CGAL::assign(l, o)) h = l->circle().opposite().plane_through(l->incident_sface()->center_vertex()->point()); else CGAL_error_msg( "wrong handle"); D.create_facet_objects(h,it->second.begin(),it->second.end()); } // CGAL_NEF_SETDTHREAD(1); } void create_volumes() { Base::create_volumes(); } void build_external_structure() { // CGAL_NEF_SETDTHREAD(503*509); #ifdef CGAL_NEF3_TIMER_EXTERNAL_STRUCTURE CGAL::Timer timer_external_structure; timer_external_structure.start(); #endif #ifdef CGAL_NEF3_TIMER_PLUECKER CGAL::Timer timer_pluecker; timer_pluecker.start(); #endif // SNC_io_parser O0(std::cerr,*this->sncp()); // O0.print(); pair_up_halfedges(); #ifdef CGAL_NEF3_TIMER_PLUECKER timer_pluecker.stop(); if(cgal_nef3_timer_on) std::cout << "Runtime_pluecker: " << timer_pluecker.time() << std::endl; #endif // SNC_io_parser O0(std::cerr,*this->sncp()); // O0.print(); link_shalfedges_to_facet_cycles(); std::map hash; CGAL::Unique_hash_map done(false); SHalfedge_iterator sei; CGAL_forall_shalfedges(sei, *this->sncp()) { hash[sei->get_index()] = sei->get_index(); } CGAL_forall_shalfedges(sei, *this->sncp()) { if(done[sei]) continue; SHalfedge_around_facet_circulator circ(sei), end(circ); int index = circ->get_index(); ++circ; CGAL_For_all(circ, end) if(circ->get_index() < index) index = circ->get_index(); index = hash[index]; CGAL_For_all(circ, end) { hash[circ->get_index()] = index; circ->set_index(index); done[circ] = true; } } SHalfloop_iterator sli; CGAL_forall_shalfloops(sli, *this->sncp()) sli->set_index(hash[sli->get_index()]); categorize_facet_cycles_and_create_facets(); create_volumes(); #ifdef CGAL_NEF3_TIMER_EXTERNAL_STRUCTURE timer_external_structure.stop(); if(cgal_nef3_timer_on) std::cout << "Runtime_external_structure: " << timer_external_structure.time() << std::endl; #endif } void clear_external_structure() { Base::clear_external_structure(); } template void build_after_binary_operation(Association& A) { SHalfedge_iterator sei; CGAL_forall_shalfedges(sei, *this->sncp()) { CGAL_NEF_TRACEN("hash sedge " << sei->get_index() << "->" << A.get_hash(sei->get_index())); sei->set_index(A.get_hash(sei->get_index())); } SHalfloop_iterator sli; CGAL_forall_shalfloops(sli, *this->sncp()) { CGAL_NEF_TRACEN("hash sloop " << sli->get_index() << "->" << A.get_hash(sli->get_index())); sli->set_index(A.get_hash(sli->get_index())); } // CGAL_NEF_SETDTHREAD(43); /* { CGAL::SNC_io_parser O (std::cerr, *this->sncp(), false, true); O.print();} */ pair_up_halfedges(); /* { CGAL::SNC_io_parser O (std::cerr, *this->sncp(), false, true); O.print();} */ link_shalfedges_to_facet_cycles(); SNC_simplify simp(*this->sncp()); simp.vertex_simplificationI(); // std::map hash; CGAL::Unique_hash_map done(false); /* SHalfedge_iterator sei; CGAL_forall_shalfedges(sei, *this->sncp()) { hash[sei->get_forward_index()] = sei->get_forward_index(); hash[sei->get_backward_index()] = sei->get_backward_index(); } */ categorize_facet_cycles_and_create_facets(); create_volumes(); } }; } //namespace CGAL #endif //CGAL_SNC_EXTERNAL_STRUCTURE_H