1 // Copyright (c) 1997-2002  Max-Planck-Institute Saarbruecken (Germany).
2 // All rights reserved.
3 //
4 // This file is part of CGAL (www.cgal.org).
5 //
6 // $URL: https://github.com/CGAL/cgal/blob/v5.3/Nef_3/include/CGAL/Nef_3/SNC_external_structure.h $
7 // $Id: SNC_external_structure.h ce7d06d 2020-12-05T08:12:56+00:00 Giles Bathgate
8 // SPDX-License-Identifier: GPL-3.0-or-later OR LicenseRef-Commercial
9 //
10 //
11 // Author(s)     : Michael Seel       <seel@mpi-sb.mpg.de>
12 //                 Miguel Granados    <granados@mpi-sb.mpg.de>
13 //                 Susan Hert         <hert@mpi-sb.mpg.de>
14 //                 Lutz Kettner       <kettner@mpi-sb.mpg.de>
15 //                 Peter Hachenberger <hachenberger@mpi-sb.mpg.de>
16 #ifndef CGAL_SNC_EXTERNAL_STRUCTURE_H
17 #define CGAL_SNC_EXTERNAL_STRUCTURE_H
18 
19 #include <CGAL/license/Nef_3.h>
20 
21 
22 #include <CGAL/Nef_S2/Normalizing.h>
23 #include <CGAL/Nef_3/Pluecker_line_3.h>
24 #include <CGAL/Nef_3/SNC_decorator.h>
25 #include <CGAL/Nef_3/SNC_point_locator.h>
26 #include <CGAL/Nef_S2/SM_point_locator.h>
27 #include <CGAL/Nef_3/SNC_FM_decorator.h>
28 #include <CGAL/Nef_3/SNC_io_parser.h>
29 #include <CGAL/Nef_3/SNC_indexed_items.h>
30 #include <CGAL/Nef_3/SNC_simplify.h>
31 #include <map>
32 #include <list>
33 
34 #undef CGAL_NEF_DEBUG
35 #define CGAL_NEF_DEBUG 43
36 #include <CGAL/Nef_2/debug.h>
37 
38 #include <CGAL/use.h>
39 
40 namespace CGAL {
41 
42 struct int_lt {
operatorint_lt43   bool operator()(const int& i1, const int& i2) const { return i1<i2; }
44 };
45 template <typename Edge_handle>
46 struct Halfedge_key_lt4 {
47 
operatorHalfedge_key_lt448   bool operator()(const Edge_handle& e1, const Edge_handle& e2) const {
49     if(CGAL::sign(e1->point().x()) != 0) {
50       if(e1->source() != e2->source())
51         return CGAL::compare_x(e1->source()->point(), e2->source()->point()) < 0;
52       else
53         return e1->point().x() < 0;
54     }
55     if(CGAL::sign(e1->point().y()) != 0) {
56       if(e1->source() != e2->source())
57         return CGAL::compare_y(e1->source()->point(), e2->source()->point()) < 0;
58       else
59         return e1->point().y() < 0;
60     }
61     if(e1->source() != e2->source())
62       return CGAL::compare_z(e1->source()->point(), e2->source()->point()) < 0;
63     return e1->point().z() < 0;
64   }
65 };
66 
67 template <typename Edge_handle>
68 struct Halfedge_key_lt3 {
69 
operatorHalfedge_key_lt370   bool operator()(const Edge_handle& e1, const Edge_handle& e2) const {
71     if(e1->source() != e2->source())
72       return CGAL::lexicographically_xyz_smaller(e1->source()->point(), e2->source()->point());
73     if(CGAL::sign(e1->point().x()) != 0)
74       return e1->point().x() < 0;
75     if(CGAL::sign(e1->point().y()) != 0)
76       return e1->point().y() < 0;
77     return e1->point().z() < 0;
78   }
79 };
80 
81 template <typename Point, typename Edge>
82 struct Halfedge_key {
83   typedef Halfedge_key<Point,Edge> Self;
84   Point p; int i; Edge e;
Halfedge_keyHalfedge_key85   Halfedge_key(Point pi, int ii, Edge ei) :
86     p(pi), i(ii), e(ei) {}
Halfedge_keyHalfedge_key87   Halfedge_key(const Self& k) : p(k.p), i(k.i), e(k.e) {}
88   Self& operator=(const Self& k) { p=k.p; i=k.i; e=k.e; return *this; }
89   bool operator==(const Self& k) const { return p==k.p && i==k.i; }
90   bool operator!=(const Self& k) const { return !operator==(k); }
91 };
92 
93 template <typename Point, typename Edge, class Decorator>
94 struct Halfedge_key_lt {
95   typedef Halfedge_key<Point,Edge> Key;
96   typedef typename Point::R R;
97   typedef typename R::Vector_3 Vector;
98   typedef typename R::Direction_3 Direction;
operatorHalfedge_key_lt99   bool operator()( const Key& k1, const Key& k2) const {
100     if( k1.e->source() == k2.e->source())
101       return (k1.i < k2.i);
102     Direction l(k1.e->vector());
103     if( k1.i < 0) l = -l;
104     return (Direction( k2.p - k1.p) == l);
105   }
106 };
107 
108 template <typename Point, typename Edge>
109 std::ostream& operator<<(std::ostream& os,
110                          const Halfedge_key<Point,Edge>& k )
111 { os << k.p << " " << k.i; return os; }
112 
113 template <typename R>
sign_of(const CGAL::Plane_3<R> & h)114 int sign_of(const CGAL::Plane_3<R>& h)
115 { if ( h.c() != 0 ) return CGAL_NTS sign(h.c());
116   if ( h.b() != 0 ) return CGAL_NTS sign(-h.b());
117   return CGAL_NTS sign(h.a());
118 }
119 
120 struct Plane_lt {
121     template <typename R>
operatorPlane_lt122     bool operator()(const CGAL::Plane_3<R>& h1,
123                     const CGAL::Plane_3<R>& h2) const
124     {
125         typedef typename R::RT     RT;
126         typedef typename R::FT     FT;
127 
128         bool flag=false;
129         FT ratioa,ratiob,ratioc,ratiod;
130         FT a1,b1,c1,d1,a2,b2,c2,d2;
131         a1=h1.a();
132         a2=h2.a();
133         b1=h1.b();
134         b2=h2.b();
135         c1=h1.c();
136         c2=h2.c();
137         d1=h1.d();
138         d2=h2.d();
139         if(a2==0 || a1==0)
140         {
141             if(a2==a1) ratioa=1;
142             else flag=true;
143         }
144         else{ratioa=a1/a2;}
145         if(b2==0 || b1==0)
146         {
147             if(b2==b1) ratiob=ratioa;
148             else flag=true;
149         }
150         else{ ratiob=b1/b2;}
151         if(c2==0 || c1==0)
152         {
153             if(c2==c1) ratioc=ratioa;
154             else flag=true;
155         }
156         else{ ratioc=c1/c2;}
157         if(d2==0 || d1==0)
158         {
159             if(d2==d1) ratiod=ratioc;
160             else flag=true;
161         }
162         else{ ratiod=d1/d2;}
163         if(flag || !(ratioa==ratiob && ratioc==ratiod && ratioa==ratioc))
164         {
165             RT diff = h1.a()-h2.a();
166             if ( (diff) != 0 ) return CGAL_NTS sign(diff) < 0;
167             diff = h1.b()-h2.b();
168             if ( (diff) != 0 ) return CGAL_NTS sign(diff) < 0;
169             diff = h1.c()-h2.c();
170             if ( (diff) != 0 ) return CGAL_NTS sign(diff) < 0;
171             diff = h1.d()-h2.d();
172             if ( (diff) != 0 ) return CGAL_NTS sign(diff) < 0;
173         }
174         return 0;
175     }
176 };
177 
178 
179 struct Plane_RT_lt {
180   template <typename R>
operatorPlane_RT_lt181   bool operator()(const CGAL::Plane_3<R>& h1,
182                   const CGAL::Plane_3<R>& h2) const
183   { return (&(h1.a())) < (&(h2.a())); }
184 };
185 
186 // ----------------------------------------------------------------------------
187 // SNC_external_structure
188 // ----------------------------------------------------------------------------
189 
190 template <typename Items_, typename SNC_structure_>
191 class SNC_external_structure_base : public SNC_decorator<SNC_structure_>
192 {
193 public:
194   typedef SNC_structure_ SNC_structure;
195 
196   typedef typename SNC_structure::Infi_box                   Infi_box;
197   typedef typename Infi_box::Standard_kernel                 Standard_kernel;
198   typedef typename Standard_kernel::Point_3                  Standard_point_3;
199   typedef typename Infi_box::NT                              NT;
200 
201   typedef CGAL::SNC_decorator<SNC_structure>                 SNC_decorator;
202   typedef CGAL::SNC_point_locator<SNC_decorator>             SNC_point_locator;
203   typedef CGAL::SNC_FM_decorator<SNC_structure>              FM_decorator;
204   typedef CGAL::SNC_simplify<Items_, SNC_structure>          SNC_simplify;
205 
206   typedef typename SNC_structure::Sphere_map                 Sphere_map;
207   typedef CGAL::SM_decorator<Sphere_map>                     SM_decorator;
208   typedef CGAL::SM_const_decorator<Sphere_map>               SM_const_decorator;
209   typedef CGAL::SM_point_locator<SM_const_decorator>         SM_point_locator;
210 
211   typedef typename SNC_structure::Halfedge_iterator Halfedge_iterator;
212   typedef typename SNC_structure::Halffacet_iterator Halffacet_iterator;
213   typedef typename SNC_structure::Volume_iterator Volume_iterator;
214 
215   typedef typename SNC_structure::Vertex_handle Vertex_handle;
216   typedef typename SNC_structure::Halfedge_handle Halfedge_handle;
217   typedef typename SNC_structure::Halffacet_handle Halffacet_handle;
218   typedef typename SNC_structure::Volume_handle Volume_handle;
219 
220   typedef typename SNC_structure::Halffacet_const_handle Halffacet_const_handle;
221   typedef typename SNC_structure::SFace_const_handle SFace_const_handle;
222 
223   typedef typename SNC_structure::SHalfedge_iterator SHalfedge_iterator;
224   typedef typename SNC_structure::SFace_iterator SFace_iterator;
225   typedef typename SNC_structure::SHalfloop_iterator SHalfloop_iterator;
226 
227   typedef typename SNC_structure::SVertex_handle SVertex_handle;
228   typedef typename SNC_structure::SHalfedge_handle SHalfedge_handle;
229   typedef typename SNC_structure::SFace_handle SFace_handle;
230   typedef typename SNC_structure::SHalfloop_handle SHalfloop_handle;
231 
232   typedef typename SNC_structure::Object_handle Object_handle;
233 
234   typedef typename SNC_structure::SHalfedge_around_facet_circulator
235     SHalfedge_around_facet_circulator;
236 
237   typedef typename SNC_structure::Point_3 Point_3;
238   typedef typename SNC_structure::Direction_3 Direction_3;
239   typedef typename SNC_structure::Plane_3 Plane_3;
240   typedef typename SNC_structure::Ray_3 Ray_3;
241 
242   typedef typename SNC_structure::Sphere_point Sphere_point;
243   typedef typename SNC_structure::Sphere_circle Sphere_circle;
244 
245   typedef typename SM_decorator::SHalfedge_around_svertex_circulator
246     SHalfedge_around_svertex_circulator;
247 
248   SNC_point_locator* pl;
249 
250   typedef CGAL::Unique_hash_map<SFace_const_handle,unsigned int>
251                                                          Sface_shell_hash;
252   typedef CGAL::Unique_hash_map<Halffacet_const_handle,unsigned int>
253                                                          Face_shell_hash;
254   typedef CGAL::Unique_hash_map<SFace_const_handle,bool> SFace_visited_hash;
255   typedef CGAL::Unique_hash_map<SFace_const_handle,bool> Shell_closed_hash;
256 
257   using SNC_decorator::visit_shell_objects;
258   using SNC_decorator::link_as_inner_shell;
259   using SNC_decorator::link_as_outer_shell;
260   using SNC_decorator::link_as_prev_next_pair;
261   using SNC_decorator::get_visible_facet;
262   using SNC_decorator::adjacent_sface;
263   using SNC_decorator::make_twins;
264 
265   struct Shell_explorer {
266     const SNC_decorator& D;
267     Sface_shell_hash&  ShellSf;
268     Face_shell_hash&   ShellF;
269     //    Shell_closed_hash& Closed;
270     SFace_visited_hash& Done;
271     SFace_handle sf_min;
272     int n;
273 
Shell_explorerShell_explorer274     Shell_explorer(const SNC_decorator& Di, Sface_shell_hash& SSf,
275                    Face_shell_hash& SF, SFace_visited_hash& Vi)
276       : D(Di), ShellSf(SSf), ShellF(SF), Done(Vi), n(0) {}
277 
visitShell_explorer278     void visit(SFace_handle h) {
279       CGAL_NEF_TRACEN("visit sf "<<h->center_vertex()->point());
280       ShellSf[h]=n;
281       Done[h]=true;
282       if ( CGAL::lexicographically_xyz_smaller(h->center_vertex()->point(),
283                                                sf_min->center_vertex()->point()))
284         sf_min = h;
285     }
286 
visitShell_explorer287     void visit(Vertex_handle h) {
288       CGAL_USE(h);
289       CGAL_NEF_TRACEN("visit v  "<<h->point());
290     }
291 
visitShell_explorer292     void visit(Halfedge_handle h) {
293       CGAL_USE(h);
294       CGAL_NEF_TRACEN("visit he "<< h->source()->point());
295     }
296 
visitShell_explorer297     void visit(Halffacet_handle h) {
298       CGAL_NEF_TRACEN(h->plane());
299       ShellF[h]=n;
300     }
301 
visitShell_explorer302     void visit(SHalfedge_handle ) {}
visitShell_explorer303     void visit(SHalfloop_handle ) {}
304 
minimal_sfaceShell_explorer305     SFace_handle& minimal_sface() { return sf_min; }
306 
increment_shell_numberShell_explorer307     void increment_shell_number() {
308       CGAL_NEF_TRACEN("leaving shell "<<n);
309       ++n;
310     }
311   };
312 
313   SNC_external_structure_base( SNC_structure& W, SNC_point_locator* spl = nullptr)
SNC_decorator(W)314     : SNC_decorator(W), pl(spl) {}
315   /*{\Mcreate makes |\Mvar| a decorator of |W|.}*/
316 
317  public:
318   //#define CGAL_NEF_NO_HALFEDGE_KEYS
319 #ifdef CGAL_NEF_NO_HALFEDGE_KEYS
pair_up_halfedges()320   void pair_up_halfedges() const {
321   /*{\Mop pairs all halfedge stubs to create the edges in 3-space.}*/
322 
323   //  CGAL_NEF_SETDTHREAD(43*61);
324     CGAL_NEF_TRACEN(">>>>>pair_up_halfedges");
325     typedef Halfedge_key_lt3<Halfedge_handle>
326       Halfedge_key_lt;
327     typedef std::list<Halfedge_handle>  Halfedge_list;
328 
329     typedef typename Standard_kernel::Kernel_tag   Kernel_tag;
330     typedef CGAL::Pluecker_line_3<Kernel_tag,Standard_kernel> Pluecker_line_3;
331     typedef CGAL::Pluecker_line_lt        Pluecker_line_lt;
332     typedef std::map< Pluecker_line_3, Halfedge_list, Pluecker_line_lt>
333       Pluecker_line_map;
334 
335     Pluecker_line_map M;
336     Pluecker_line_map M2;
337     Pluecker_line_map M3;
338     Pluecker_line_map M4;
339 
340     NT eval(Infi_box::compute_evaluation_constant_for_halfedge_pairup(*this->sncp()));
341 
342     Halfedge_iterator e;
343     CGAL_forall_halfedges(e,*this->sncp()) {
344       //    progress++;
345       Point_3 p = e->source()->point();
346       Point_3 q = p + e->vector();
347       CGAL_NEF_TRACE(" segment("<<p<<", "<<q<<")"<<
348             " direction("<<e->vector()<<")");
349       Standard_point_3 sp = Infi_box::standard_point(p,eval);
350       Standard_point_3 sq = Infi_box::standard_point(q,eval);
351       Pluecker_line_3  l( sp, sq);
352 
353       int inverted;
354       l = categorize( l, inverted);
355 
356       if(Infi_box::is_edge_on_infibox(e))
357         if(Infi_box::is_type4(e))
358           M4[l].push_back(e);
359         else
360           if(Infi_box::is_type3(e))
361             M3[l].push_back(e);
362           else
363             M2[l].push_back(e);
364       else
365         M[l].push_back(e);
366 
367       // the following trace crashes when compiling with optimizations (-O{n})
368       //CGAL_NEF_TRACEN(Infi_box::standard_point(point(vertex(e)))+
369 
370       CGAL_NEF_TRACEN(" line("<<l<<")"<<" inverted="<<inverted);
371     }
372 
373     typename Pluecker_line_map::iterator it;
374 
375     CGAL_forall_iterators(it,M4) {
376       //    progress++;
377       it->second.sort(Halfedge_key_lt());
378       CGAL_NEF_TRACEN("search opposite  "<<it->first);
379       typename Halfedge_list::iterator itl;
380       CGAL_forall_iterators(itl,it->second) {
381         Halfedge_handle e1 = *itl;
382         ++itl;
383         CGAL_assertion(itl != it->second.end());
384         Halfedge_handle e2 = *itl;
385         while(normalized(e1->vector()) != normalized(-e2->vector())) {
386           ++itl;
387           make_twins(e1,*itl);
388           e1 = e2;
389           ++itl;
390           e2 = *itl;
391         }
392         CGAL_NEF_TRACEN("    " << e1->source()->point()
393                         << " -> " << e2->source()->point());
394         CGAL_NEF_TRACEN(e1->vector()<<" -> "<<-e2->vector());
395         make_twins(e1,e2);
396         CGAL_assertion(e1->mark()==e2->mark());
397 
398         // discard temporary sphere_point ?
399       }
400     }
401 
402     CGAL_forall_iterators(it,M3) {
403       //    progress++;
404       it->second.sort(Halfedge_key_lt());
405       CGAL_NEF_TRACEN("search opposite  "<<it->first);
406       typename Halfedge_list::iterator itl;
407       CGAL_forall_iterators(itl,it->second) {
408         Halfedge_handle e1 = *itl;
409         ++itl;
410         CGAL_assertion(itl != it->second.end());
411         Halfedge_handle e2 = *itl;
412         while(normalized(e1->vector()) != normalized(-e2->vector())) {
413           ++itl;
414           make_twins(e1,*itl);
415           e1 = e2;
416           ++itl;
417           e2 = *itl;
418         }
419         CGAL_NEF_TRACEN("    " << e1->source()->point()
420                         << " -> " << e2->source()->point());
421         CGAL_NEF_TRACEN(e1->vector()<<" -> "<<-e2->vector());
422         make_twins(e1,e2);
423         CGAL_assertion(e1->mark()==e2->mark());
424 
425         // discard temporary sphere_point ?
426       }
427     }
428 
429     CGAL_forall_iterators(it,M2) {
430       //    progress++;
431       it->second.sort(Halfedge_key_lt());
432       CGAL_NEF_TRACEN("search opposite  "<<it->first);
433       typename Halfedge_list::iterator itl;
434       CGAL_forall_iterators(itl,it->second) {
435         Halfedge_handle e1 = *itl;
436         ++itl;
437         CGAL_assertion(itl != it->second.end());
438         Halfedge_handle e2 = *itl;
439         while(normalized(e1->vector()) != normalized(-e2->vector())) {
440           ++itl;
441           this->make_twins(e1,*itl);
442           e1 = e2;
443           ++itl;
444           e2 = *itl;
445         }
446         CGAL_NEF_TRACEN("    " << e1->source()->point()
447                         << " -> " << e2->source()->point());
448         CGAL_NEF_TRACEN(e1->vector()<<" -> "<<-e2->vector());
449         this->make_twins(e1,e2);
450         CGAL_assertion(e1->mark()==e2->mark());
451 
452         // discard temporary sphere_point ?
453       }
454     }
455 
456     CGAL_forall_iterators(it,M) {
457       //    progress++;
458       it->second.sort(Halfedge_key_lt());
459       CGAL_NEF_TRACEN("search opposite  "<<it->first);
460       typename Halfedge_list::iterator itl;
461       CGAL_forall_iterators(itl,it->second) {
462         Halfedge_handle e1 = *itl;
463         ++itl;
464         CGAL_assertion(itl != it->second.end());
465         Halfedge_handle e2 = *itl;
466         while(normalized(e1->vector()) != normalized(-e2->vector())) {
467           ++itl;
468           this->make_twins(e1,*itl);
469           e1 = e2;
470           ++itl;
471           e2 = *itl;
472         }
473         CGAL_NEF_TRACEN("    " << e1->source()->point()
474                         << " -> " << e2->source()->point());
475         CGAL_NEF_TRACEN(e1->vector()<<" -> "<<-e2->vector());
476         this->make_twins(e1,e2);
477         CGAL_assertion(e1->mark()==e2->mark());
478 
479         // discard temporary sphere_point ?
480       }
481     }
482 
483   }
484 #else
485   void pair_up_halfedges() const {
486   /*{\Mop pairs all halfedge stubs to create the edges in 3-space.}*/
487 
488 //    CGAL_NEF_SETDTHREAD(43*61);
489     CGAL_NEF_TRACEN(">>>>>pair_up_halfedges");
490     typedef Halfedge_key< Point_3, Halfedge_handle>
491       Halfedge_key;
492     typedef Halfedge_key_lt< Point_3, Halfedge_handle, SNC_decorator>
493       Halfedge_key_lt;
494     typedef std::list<Halfedge_key>  Halfedge_list;
495 
496     typedef typename Standard_kernel::Kernel_tag   Kernel_tag;
497     typedef CGAL::Pluecker_line_3<Kernel_tag,Standard_kernel> Pluecker_line_3;
498     typedef CGAL::Pluecker_line_lt        Pluecker_line_lt;
499     typedef std::map< Pluecker_line_3, Halfedge_list, Pluecker_line_lt>
500       Pluecker_line_map;
501 
502     SNC_decorator D(*this);
503     Pluecker_line_map M;
504     Pluecker_line_map M2;
505     Pluecker_line_map M3;
506     Pluecker_line_map M4;
507 
508     NT eval(Infi_box::compute_evaluation_constant_for_halfedge_pairup(*this->sncp()));;
509 
510     Halfedge_iterator e;
511     CGAL_forall_halfedges(e,*this->sncp()) {
512       //    progress++;
513       Point_3 p = e->source()->point();
514       Point_3 q = p + e->vector();
515       CGAL_NEF_TRACE(" segment("<<p<<", "<<q<<")"<<
516             " direction("<<e->vector()<<")");
517       Standard_point_3 sp = Infi_box::standard_point(p,eval);
518       Standard_point_3 sq = Infi_box::standard_point(q,eval);
519       Pluecker_line_3  l( sp, sq);
520 
521       int inverted;
522       l = categorize( l, inverted);
523 
524       if(Infi_box::is_edge_on_infibox(e))
525         if(Infi_box::is_type4(e))
526           M4[l].push_back(Halfedge_key(p,inverted,e));
527         else
528           if(Infi_box::is_type3(e))
529             M3[l].push_back(Halfedge_key(p,inverted,e));
530           else
531             M2[l].push_back(Halfedge_key(p,inverted,e));
532       else
533         M[l].push_back(Halfedge_key(p,inverted,e));
534 
535       // the following trace crashes when compiling with optimizations (-O{n})
536       //CGAL_NEF_TRACEN(Infi_box::standard_point(point(vertex(e)))+
537 
538       CGAL_NEF_TRACEN(" line("<<l<<")"<<" inverted="<<inverted);
539     }
540 
541     typename Pluecker_line_map::iterator it;
542 
543     CGAL_forall_iterators(it,M4) {
544       //    progress++;
545       it->second.sort(Halfedge_key_lt());
546       CGAL_NEF_TRACEN("search opposite  "<<it->first);
547       typename Halfedge_list::iterator itl;
548       CGAL_forall_iterators(itl,it->second) {
549         Halfedge_handle e1 = itl->e;
550         ++itl;
551         CGAL_assertion(itl != it->second.end());
552         Halfedge_handle e2 = itl->e;
553         CGAL_NEF_TRACEN("    " << e1->source()->point()
554                         << " -> " << e2->source()->point());
555         CGAL_NEF_TRACEN(e1->vector()<<" -> "<<-e2->vector());
556         make_twins(e1,e2);
557         CGAL_assertion(e1->mark()==e2->mark());
558 
559         // discard temporary sphere_point ?
560       }
561     }
562 
563     CGAL_forall_iterators(it,M3) {
564       //    progress++;
565       it->second.sort(Halfedge_key_lt());
566       CGAL_NEF_TRACEN("search opposite  "<<it->first);
567       typename Halfedge_list::iterator itl;
568       CGAL_forall_iterators(itl,it->second) {
569         Halfedge_handle e1 = itl->e;
570         ++itl;
571         CGAL_assertion(itl != it->second.end());
572         Halfedge_handle e2 = itl->e;
573         CGAL_NEF_TRACEN("    " << e1->source()->point()
574                         << " -> " << e2->source()->point());
575         CGAL_NEF_TRACEN(e1->vector()<<" -> "<<-e2->vector());
576         make_twins(e1,e2);
577         CGAL_assertion(e1->mark()==e2->mark());
578 
579         // discard temporary sphere_point ?
580       }
581     }
582 
583     CGAL_forall_iterators(it,M2) {
584       //    progress++;
585       it->second.sort(Halfedge_key_lt());
586       CGAL_NEF_TRACEN("search opposite  "<<it->first);
587       typename Halfedge_list::iterator itl;
588       CGAL_forall_iterators(itl,it->second) {
589         Halfedge_handle e1 = itl->e;
590         ++itl;
591         CGAL_assertion(itl != it->second.end());
592         Halfedge_handle e2 = itl->e;
593         CGAL_NEF_TRACEN("    " << e1->source()->point()
594                         << " -> " << e2->source()->point());
595         CGAL_NEF_TRACEN(e1->vector()<<" -> "<<-e2->vector());
596         make_twins(e1,e2);
597         CGAL_assertion(e1->mark()==e2->mark());
598 
599         // discard temporary sphere_point ?
600       }
601     }
602 
603     CGAL_forall_iterators(it,M) {
604       //    progress++;
605       it->second.sort(Halfedge_key_lt());
606       CGAL_NEF_TRACEN("search opposite  "<<it->first);
607       typename Halfedge_list::iterator itl;
608       CGAL_forall_iterators(itl,it->second) {
609         Halfedge_handle e1 = itl->e;
610         ++itl;
611         CGAL_assertion(itl != it->second.end());
612         Halfedge_handle e2 = itl->e;
613         CGAL_NEF_TRACEN("    " << e1->source()->point()
614                         << " -> " << e2->source()->point());
615         CGAL_NEF_TRACEN(e1->vector()<<" -> "<< -e2->vector());
616         CGAL_assertion(e1->source()->point() != e2->source()->point());
617         CGAL_assertion(e1->mark()==e2->mark());
618         make_twins(e1,e2);
619 
620         // discard temporary sphere_point ?
621       }
622     }
623   }
624 #endif
625 
link_shalfedges_to_facet_cycles()626   void link_shalfedges_to_facet_cycles() const {
627   /*{\Mop creates all non-trivial facet cycles from sedges.
628   \precond |pair_up_halfedges()| was called before.}*/
629 
630   //  CGAL_NEF_SETDTHREAD(43*31);
631     CGAL_NEF_TRACEN(">>>>>link_shalfedges_to_facet_cycles");
632 
633 #ifdef CGAL_NEF_EXPLOIT_REFERENCE_COUNTING
634     Point_3 p1(1,2,7), p2(p1);
635     bool reference_counted = (&(p1.hx()) == &(p2.hx()));
636 #endif
637 
638     Halfedge_iterator e;
639     CGAL_forall_edges(e,*this->sncp()) {
640       //    progress++;
641       CGAL_NEF_TRACEN("");
642       CGAL_NEF_TRACEN(PH(e));
643       Halfedge_iterator et = e->twin();
644       SM_decorator D(&*e->source()), Dt(&*et->source());
645       CGAL_NEF_TRACEN(e->source()->point());
646       if ( D.is_isolated(e) ) continue;
647       SHalfedge_around_svertex_circulator ce(D.first_out_edge(e)),cee(ce);
648       SHalfedge_around_svertex_circulator cet(Dt.first_out_edge(et)),cete(cet);
649 
650 #ifdef CGAL_NEF_EXPLOIT_REFERENCE_COUNTING
651       if(reference_counted) {
652         CGAL_For_all(cet,cete)
653           if ( &(cet->circle().a()) == &(ce->circle().opposite().a()) &&
654                cet->source()->twin() == ce->source() )
655             break;
656       } else
657 #endif
658         CGAL_For_all(cet,cete)
659           if ( cet->circle() == ce->circle().opposite() &&
660                cet->source()->twin() == ce->source() )
661             break;
662 
663 #ifdef CGAL_USE_TRACE
664       if( cet->circle() != ce->circle().opposite() )
665         CGAL_NEF_TRACEN("assertion failed!");
666 
667       CGAL_NEF_TRACEN("vertices " << e->source()->point() <<
668       "    "      << et->source()->point());
669 
670 
671       SHalfedge_around_svertex_circulator sc(D.first_out_edge(e));
672       SHalfedge_around_svertex_circulator sct(Dt.first_out_edge(et));
673 
674       CGAL_NEF_TRACEN("");
675       CGAL_For_all(sc,cee)
676         CGAL_NEF_TRACEN("sseg@E addr="<<&*sc<<
677                         " src="<< sc->source()->point()<<
678                         " tgt="<< sc->target()->point()<<std::endl<<
679                         " circle=" << normalized(sc->circle()));
680       CGAL_NEF_TRACEN("");
681 
682       CGAL_For_all(sct,cete)
683       CGAL_NEF_TRACEN("sseg@ET addr="<<&*sct<<
684                       " src="<< sct->source()->point()<<
685                       " tgt="<<sct->target()->point()<<std::endl<<
686                       " circle=" << normalized(sct->circle()));
687       CGAL_NEF_TRACEN("");
688 #endif
689 
690       CGAL_assertion( normalized(cet->circle()) == normalized(ce->circle().opposite()) );
691       CGAL_assertion( cet->source()->twin() == ce->source());
692       CGAL_For_all(ce,cee) {
693         CGAL_NEF_TRACEN("circles " << normalized(cet->circle()) << "   " << normalized(ce->circle()) <<
694                         " sources " << cet->target()->point() <<
695                         "   " << ce->target()->point());
696         CGAL_assertion( normalized(cet->circle()) == normalized(ce->circle().opposite()));
697         CGAL_assertion( cet->source()->twin() == ce->source());
698         CGAL_assertion(ce->mark()==cet->mark());
699         link_as_prev_next_pair(cet->twin(),ce);
700         link_as_prev_next_pair(ce->twin(),cet);
701         --cet; // ce moves ccw, cet moves cw
702       }
703     }
704   }
705 
categorize_facet_cycles_and_create_facets()706   void categorize_facet_cycles_and_create_facets() const {
707   /*{\Mop collects all facet cycles incident to a facet and creates
708   the facets. \precond |link_shalfedges_to_facet_cycles()| was called
709   before.}*/
710 
711     //    CGAL_NEF_SETDTHREAD(43*31);
712     CGAL_NEF_TRACEN(">>>>>categorize_facet_cycles_and_create_facets");
713 
714     typedef std::list<Object_handle> Object_list;
715 #ifdef CGAL_NEF_EXPLOIT_REFERENCE_COUNTING
716     typedef std::map<Plane_3, Object_list, Plane_RT_lt>
717       Map_planes;
718 #else
719     typedef std::map<Plane_3, Object_list, Plane_lt>
720       Map_planes;
721 #endif
722 
723     Map_planes M;
724     SHalfedge_iterator e;
725     CGAL_forall_shalfedges(e,*this->sncp()) {
726       Sphere_circle c(e->circle());
727       Plane_3 h = c.plane_through(e->source()->source()->point());
728       CGAL_NEF_TRACEN("\n" << e->source()->twin()->source()->point() <<" - "
729                       << e->source()->source()->point() <<" - "<<
730                       e->twin()->source()->twin()->source()->point() <<
731                       " has plane " << h << " has circle " << e->circle() <<
732                       " has signum " << sign_of(h));
733       if ( sign_of(h)<0 ) continue;
734       M[normalized(h)].push_back(make_object(e->twin()));
735       CGAL_NEF_TRACEN(" normalized as " << normalized(h));
736       /*
737         Unique_hash_map<SHalfedge_handle, bool> Done(false);
738         SHalfedge_iterator ei;
739         CGAL_forall_sedges(ei,*this->sncp()) {
740         if(Done[ei]) continue;
741         Sphere_circle c(ei->circle());
742         Plane_3 h = c.plane_through(ei->source()->source()->point());
743         CGAL_NEF_TRACEN("\n" << ei->source()->twin()->source()->point() <<" - "
744         << ei->source()->source()->point() <<" - "<<
745                       ei->twin()->source()->twin()->source()->point() <<
746                       " has plane " << h << " has circle " << ei->circle() <<
747                       " has signum " << sign_of(h));
748       SHalfedge_handle e(ei);
749       if ( sign_of(h)<0 ) {
750         h = h.opposite();
751         e = e->twin();
752       }
753       SHalfedge_around_facet_circulator sfc(e), send(sfc);
754       CGAL_For_all(sfc, send) {
755         M[normalized(h)].push_back(make_object(e->twin()));
756         Done[sfc] = true;
757         Done[sfc->twin()] = true;
758         CGAL_NEF_TRACEN(" normalized as " << normalized(h));
759       }
760       */
761     }
762     SHalfloop_iterator l;
763     CGAL_forall_shalfloops(l,*this->sncp()) {
764       Sphere_circle c(l->circle());
765       Plane_3 h = c.plane_through(l->incident_sface()->center_vertex()->point());
766       if ( sign_of(h)<0 ) continue;
767       // CGAL_assertion( h == normalized(h));
768       M[normalized(h)].push_back(make_object(l->twin()));
769     }
770 
771 #ifdef CGAL_NEF3_TIMER_PLANE_SWEEPS
772   number_of_plane_sweeps=0;
773   timer_plane_sweeps.reset();
774 #endif
775 
776     typename Map_planes::iterator it;
777     CGAL_forall_iterators(it,M) {
778       //    progress2++;
779       //      CGAL_NEF_TRACEN("  plane "<<it->first<<"             "<<(it->first).point());
780       FM_decorator D(*this->sncp());
781       D.create_facet_objects(it->first,it->second.begin(),it->second.end());
782     }
783     //       CGAL_NEF_SETDTHREAD(1);
784   }
785 
create_volumes()786   void create_volumes() {
787   /*{\Mop collects all shells incident to a volume and creates the
788   volumes.  \precond |categorize_facet_cycles_and_creating_facets()| was
789   called before.}*/
790 
791 #ifdef CGAL_NEF3_TIMER_POINT_LOCATION
792     number_of_ray_shooting_queries=0;
793     timer_ray_shooting.reset();
794 #endif
795 
796     //    CGAL_NEF_SETDTHREAD(37*43*503*509);
797 
798     CGAL_NEF_TRACEN(">>>>>create_volumes");
799     Sface_shell_hash     ShellSf(0);
800     Face_shell_hash      ShellF(0);
801     SFace_visited_hash Done(false);
802     Shell_explorer V(*this,ShellSf,ShellF,Done);
803     std::vector<SFace_handle> MinimalSFace;
804     std::vector<SFace_handle> EntrySFace;
805     std::vector<bool> Closed;
806 
807     SFace_iterator f;
808     // First, we classify all the Shere Faces per Shell.  For each Shell we
809     //     determine its minimum lexicographyly vertex and we check wheter the
810     //     Shell encloses a region (closed surface) or not.
811     CGAL_forall_sfaces(f,*this->sncp()) {
812       //    progress++;
813       CGAL_NEF_TRACEN("sface in " << ShellSf[f]);
814       if ( Done[f] )
815         continue;
816       V.minimal_sface() = f;
817       visit_shell_objects(f,V);
818 
819       CGAL_NEF_TRACEN("minimal vertex " << V.minimal_sface()->center_vertex()->point());
820 
821       MinimalSFace.push_back(V.minimal_sface());
822       EntrySFace.push_back(f);
823       V.increment_shell_number();
824       CGAL_NEF_TRACEN("sface out " << ShellSf[f]);
825     }
826 
827     for(unsigned int i=0; i<EntrySFace.size(); ++i)
828       Closed.push_back(false);
829 
830     Halffacet_iterator hf;
831     CGAL_forall_facets(hf,*this)
832       if(ShellF[hf] != ShellF[hf->twin()]) {
833         Closed[ShellF[hf]] = true;
834         Closed[ShellF[hf->twin()]] = true;
835       }
836 
837     CGAL_assertion( pl != nullptr);
838 
839 #ifdef CGAL_NEF3_TIMER_INITIALIZE_KDTREE
840     CGAL::Timer timer_initialize_kdtree;
841     timer_initialize_kdtree.start();
842 #endif
843     pl->initialize(this->sncp()); // construct the point locator
844 #ifdef CGAL_NEF3_TIMER_INITIALIZE_KDTREE
845     timer_initialize_kdtree.stop();
846     if(cgal_nef3_timer_on)
847       std::cout << "Runtime_initialize_kdtree: "
848                 << timer_initialize_kdtree.time() << std::endl;
849 #endif
850 
851     //   then, we determine the Shells which correspond to Volumes via a ray
852     //   shootting in the direction (-1,0,0) over the Sphere_map of the minimal
853     //   vertex.  The Shell corresponds to a Volume if the object hit belongs
854     //   to another Shell.
855 
856     this->sncp()->new_volume(); // outermost volume (nirvana)
857     if(MinimalSFace.size() == 0) return;
858     Vertex_handle v_min = MinimalSFace[0]->center_vertex();
859     for( unsigned int i = 0; i < MinimalSFace.size(); ++i) {
860       //    progress2++;
861       Vertex_handle v = MinimalSFace[i]->center_vertex();
862       if(CGAL::lexicographically_xyz_smaller(v->point(),v_min->point()))
863         v_min=v;
864       CGAL_NEF_TRACEN( "Shell #" << i << " minimal vertex: " << v->point());
865       SM_point_locator D((Sphere_map*) &*v);
866       Object_handle o = D.locate(Sphere_point(-1,0,0));
867       SFace_const_handle sfc;
868       if( !CGAL::assign(sfc, o) || ShellSf[sfc] != i) {
869         CGAL_NEF_TRACEN("the shell encloses a volume");
870         CGAL_NEF_TRACEN("sface hit? "<<CGAL::assign(sfc,o));
871         if( CGAL::assign(sfc, o)) { CGAL_NEF_TRACEN("sface on another shell? "<<ShellSf[sfc]); }
872         SFace_handle f = MinimalSFace[i];
873         CGAL_assertion( ShellSf[f] == i );
874         if( Closed[i] ) {
875           CGAL_NEF_TRACEN("Shell #" << i << " is closed");
876           SM_decorator SD(&*v);
877           Volume_handle c = this->sncp()->new_volume();
878           c->mark() = f->mark(); // TODO test if line is redundant
879           link_as_inner_shell(f, c );
880           CGAL_NEF_TRACE( "Shell #" << i <<" linked as inner shell");
881           CGAL_NEF_TRACEN( "(sface" << (CGAL::assign(sfc,o)?"":" not") << " hit case)");
882         }
883       }
884     }
885 
886     // finaly, we go through all the Shells which do not correspond to a Volume
887     //     and we assign them to its enclosing Volume determined via a facet below
888     //     check.
889 
890     CGAL_forall_sfaces(f,*this->sncp()) {
891       //    progress3++;
892       if ( f->volume() != Volume_handle() )
893         continue;
894       CGAL_NEF_TRACEN( "Outer shell #" << ShellSf[f] << " volume?");
895       Volume_handle c = determine_volume( MinimalSFace[ShellSf[f]],
896                                           MinimalSFace, ShellSf );
897       c->mark() = f->mark();
898       link_as_outer_shell( f, c );
899     }
900   }
901 
get_facet_below(Vertex_handle vi,const std::vector<SFace_handle> & MinimalSFace,const Sface_shell_hash & Shell)902   Halffacet_handle get_facet_below( Vertex_handle vi,
903                                     const std::vector< SFace_handle>& MinimalSFace,
904                                     const Sface_shell_hash&  Shell) const {
905     // {\Mop determines the facet below a vertex |vi| via ray shooting. }
906 
907     Halffacet_handle f_below;
908     Point_3 p = vi->point();
909     if(!Infi_box::is_standard(p))
910       return Halffacet_handle();
911 
912     Ray_3 ray = Ray_3(p, Direction_3(-1,0,0));
913 #ifdef CGAL_NEF3_TIMER_POINT_LOCATION
914     number_of_ray_shooting_queries++;
915     timer_ray_shooting.start();
916 #endif
917     Object_handle o = pl->shoot(ray);
918 #ifdef CGAL_NEF3_TIMER_POINT_LOCATION
919     timer_ray_shooting.stop();
920 #endif
921     // The ray here has an special property since it is shooted from the lowest
922     // vertex in a shell, so it would be expected that the ray goes along the
923     // interior of a volume before it hits a 2-skeleton element.
924     // Unfortunatelly, it seems to be possible that several shells are incident
925     // to this lowest vertex, and in consequence, the ray could also go along
926     // an edge or a facet belonging to a different shell.
927     // This fact invalidates the precondition of the get_visible_facet method,
928     // (the ray must pierce the local view of the hit object in a sface).
929     // This situation has not been analyzed and has to be verified. Granados.
930     Vertex_handle v;
931     Halfedge_handle e;
932     Halffacet_handle f;
933     CGAL_NEF_TRACEN("get_facet_below");
934     if( CGAL::assign(v, o)) {
935       CGAL_NEF_TRACEN("facet below from from vertex...");
936       f_below = get_visible_facet(v, ray);
937       if( f_below == Halffacet_handle()) {
938         CGAL_assertion(v->sfaces_begin() == v->sfaces_last());
939         f_below = get_facet_below(MinimalSFace[Shell[v->sfaces_begin()]]->center_vertex(),
940                                   MinimalSFace,Shell);
941       }
942     }
943     else if( CGAL::assign(e, o)) {
944       CGAL_NEF_TRACEN("facet below from from edge...");
945       f_below = get_visible_facet(e, ray);
946       if( f_below == Halffacet_handle()) {
947         CGAL_assertion(e->source()->sfaces_begin() == e->source()->sfaces_last());
948         f_below = get_facet_below(MinimalSFace[Shell[e->source()->sfaces_begin()]]->center_vertex(),
949                                   MinimalSFace, Shell);
950       }
951     }
952     else if( CGAL::assign(f, o)) {
953       CGAL_NEF_TRACEN("facet below from from facet...");
954       f_below = get_visible_facet(f, ray);
955       CGAL_assertion( f_below != Halffacet_handle());
956     }
957     else { CGAL_NEF_TRACEN("no facet below found..."); }
958     return f_below;
959   }
960 
determine_volume(SFace_handle sf,const std::vector<SFace_handle> & MinimalSFace,const Sface_shell_hash & Shell)961   Volume_handle determine_volume( SFace_handle sf,
962                 const std::vector< SFace_handle>& MinimalSFace,
963                                   const Sface_shell_hash&  Shell ) const {
964     //{\Mop determines the volume |C| that a shell |S| pointed by |sf|
965     //  belongs to.  \precondition |S| separates the volume |C| from an enclosed
966     //  volume.}
967 
968     CGAL_NEF_TRACEN("determine volume");
969     Vertex_handle v_min = MinimalSFace[Shell[sf]]->center_vertex();
970 
971     Halffacet_handle f_below = get_facet_below(v_min, MinimalSFace, Shell);
972     if ( f_below == Halffacet_handle())
973       return SNC_decorator(*this).volumes_begin();
974     Volume_handle c = f_below->incident_volume();
975     if( c != Volume_handle()) {
976       CGAL_NEF_TRACE( "Volume " << &*c << " hit ");
977       CGAL_NEF_TRACEN("(Shell #" << Shell[adjacent_sface(f_below)] << ")");
978       return c;
979     }
980     SFace_handle sf_below = adjacent_sface(f_below);
981     CGAL_NEF_TRACE( "Shell not assigned to a volume hit ");
982     CGAL_NEF_TRACEN( "(Inner shell #" << Shell[sf_below] << ")");
983     c = determine_volume( sf_below, MinimalSFace, Shell);
984     link_as_inner_shell( sf_below, c);
985     return c;
986   }
987 
build_external_structure()988   void build_external_structure() {
989 
990 #ifdef CGAL_NEF3_TIMER_EXTERNAL_STRUCTURE
991     CGAL::Timer timer_external_structure;
992     timer_external_structure.start();
993 #endif
994 
995 #ifdef CGAL_NEF3_TIMER_PLUECKER
996     CGAL::Timer timer_pluecker;
997     timer_pluecker.start();
998 #endif
999     //    SNC_io_parser<SNC_structure> O0(std::cerr,*this->sncp());
1000     //    O0.print();
1001     pair_up_halfedges();
1002 #ifdef CGAL_NEF3_TIMER_PLUECKER
1003     timer_pluecker.stop();
1004      if(cgal_nef3_timer_on)
1005       std::cout << "Runtime_pluecker: "
1006                 << timer_pluecker.time() << std::endl;
1007 #endif
1008     link_shalfedges_to_facet_cycles();
1009     //    SNC_io_parser<SNC_structure> O0(std::cerr,*this->sncp());
1010     //    O0.print();
1011     categorize_facet_cycles_and_create_facets();
1012     create_volumes();
1013 
1014 #ifdef CGAL_NEF3_TIMER_EXTERNAL_STRUCTURE
1015     timer_external_structure.stop();
1016     if(cgal_nef3_timer_on)
1017       std::cout << "Runtime_external_structure: "
1018                 << timer_external_structure.time() << std::endl;
1019 #endif
1020   }
1021 
1022   template<typename Association>
build_after_binary_operation(Association)1023   void build_after_binary_operation(Association) {
1024 
1025 #ifdef CGAL_NEF3_TIMER_SIMPLIFICATION
1026     CGAL::Timer timer_simplification;
1027     timer_simplification.start();
1028 #endif
1029     SNC_simplify simp(*this->sncp());
1030     simp.vertex_simplification(SNC_simplify::NO_SNC);
1031 #ifdef CGAL_NEF3_TIMER_SIMPLIFICATION
1032     timer_simplification.stop();
1033     if(cgal_nef3_timer_on)
1034       std::cout << "Runtime_simplification: "
1035                 << timer_simplification.time() << std::endl;
1036 #endif
1037 
1038     CGAL_NEF_TRACEN("\nnumber of vertices (so far...) = "
1039                     << this->sncp()->number_of_vertices());
1040 
1041     CGAL_NEF_TRACEN("=> resultant vertices (after simplification): ");
1042 
1043     build_external_structure();
1044   }
1045 
clear_external_structure()1046   void clear_external_structure() {
1047     this->sncp()->clear_snc_boundary();
1048 
1049     while(this->sncp()->volumes_begin()!= this->sncp()->volumes_end())
1050       this->sncp()->delete_volume(this->sncp()->volumes_begin());
1051 
1052     while(this->sncp()->halffacets_begin()!= this->sncp()->halffacets_end())
1053       this->sncp()->delete_halffacet_pair(this->sncp()->halffacets_begin());
1054 
1055     SHalfedge_iterator se;
1056     CGAL_forall_shalfedges(se,*this)
1057       se->facet() = Halffacet_handle();
1058 
1059     SFace_iterator sf;
1060     CGAL_forall_sfaces(sf,*this)
1061       sf->volume() = Volume_handle();
1062   }
1063 };
1064 
1065 
1066 
1067 template <typename Items_, typename SNC_structure_>
1068 class SNC_external_structure : public SNC_external_structure_base<Items_, SNC_structure_>
1069 {
1070   typedef CGAL::SNC_decorator<SNC_structure_>                 SNC_decorator;
1071   typedef CGAL::SNC_point_locator<SNC_decorator>             SNC_point_locator;
1072 public:
1073 
1074   SNC_external_structure( SNC_structure_& W, SNC_point_locator* spl = nullptr)
1075     : SNC_external_structure_base<Items_, SNC_structure_>(W, spl)
1076   {}
1077 };
1078 
1079 
1080 
1081 
1082 template <typename SNC_structure_>
1083 class SNC_external_structure<SNC_indexed_items, SNC_structure_>
1084   : public SNC_external_structure_base<int, SNC_structure_> {
1085 
1086 public:
1087   typedef SNC_structure_                                SNC_structure;
1088   typedef SNC_external_structure_base<int, SNC_structure>    Base;
1089 
1090   typedef CGAL::SNC_decorator<SNC_structure>                 SNC_decorator;
1091   typedef CGAL::SNC_point_locator<SNC_decorator>             SNC_point_locator;
1092   typedef CGAL::SNC_FM_decorator<SNC_structure>              FM_decorator;
1093   typedef typename SNC_structure::Sphere_map                 Sphere_map;
1094   typedef CGAL::SM_decorator<Sphere_map>                     SM_decorator;
1095   typedef CGAL::SNC_simplify<SNC_indexed_items, SNC_structure> SNC_simplify;
1096 
1097   typedef typename SNC_structure::Halfedge_iterator Halfedge_iterator;
1098   typedef typename SNC_structure::SHalfedge_iterator SHalfedge_iterator;
1099   typedef typename SNC_structure::SHalfloop_iterator SHalfloop_iterator;
1100 
1101   typedef typename SNC_structure::Halfedge_handle Halfedge_handle;
1102   typedef typename SNC_structure::SHalfedge_handle SHalfedge_handle;
1103 
1104   typedef typename SNC_structure::Object_handle Object_handle;
1105 
1106   typedef typename SNC_structure::SHalfedge_around_facet_circulator
1107     SHalfedge_around_facet_circulator;
1108   typedef typename SM_decorator::SHalfedge_around_svertex_circulator
1109     SHalfedge_around_svertex_circulator;
1110 
1111   typedef typename SNC_structure::Plane_3 Plane_3;
1112 
1113   using Base::make_twins;
1114   using Base::link_as_prev_next_pair;
1115   using Base::link_as_inner_shell;
1116 
1117 
1118   SNC_external_structure( SNC_structure& W, SNC_point_locator* spl = nullptr)
Base(W,spl)1119     : Base(W, spl) {}
1120   /*{\Mcreate makes |\Mvar| a decorator of |W|.}*/
1121 
1122  public:
pair_up_halfedges()1123   void pair_up_halfedges() const {
1124     typedef Halfedge_key_lt4<Halfedge_handle>  Halfedge_key_lt;
1125     typedef std::list<Halfedge_handle>  Halfedge_list;
1126     typedef std::map<int, Halfedge_list, int_lt> index_map;
1127 
1128     CGAL_NEF_TRACEN("pair up by indexes");
1129 
1130     index_map i2he;
1131     Halfedge_iterator ei;
1132     CGAL_forall_halfedges(ei, *this->sncp())
1133       i2he[ei->get_index()].push_back(ei);
1134 
1135     typename index_map::iterator it;
1136     CGAL_forall_iterators(it,i2he) {
1137       CGAL_NEF_TRACEN("pair up " << it->first);
1138       it->second.sort(Halfedge_key_lt());
1139       typename Halfedge_list::iterator itl;
1140       CGAL_forall_iterators(itl,it->second) {
1141         Halfedge_handle e1 = *itl;
1142         CGAL_NEF_TRACEN(e1->source()->point() << ", " << e1->vector());
1143         ++itl;
1144         CGAL_assertion(itl != it->second.end());
1145         Halfedge_handle e2 = *itl;
1146         CGAL_NEF_TRACEN(" + " << e2->source()->point() << ", " << e2->vector());
1147         make_twins(e1,e2);
1148         //        SNC_io_parser<SNC_structure> O0(std::cerr,*this->sncp());
1149         //        O0.print();
1150         CGAL_assertion(e1->mark()==e2->mark());
1151       }
1152     }
1153   }
1154 
link_shalfedges_to_facet_cycles()1155   void link_shalfedges_to_facet_cycles() const {
1156   /*{\Mop creates all non-trivial facet cycles from sedges.
1157   \precond |pair_up_halfedges()| was called before.}*/
1158 
1159   //  CGAL_NEF_SETDTHREAD(43*31);
1160     CGAL_NEF_TRACEN(">>>>>link_shalfedges_to_facet_cycles");
1161 
1162     Halfedge_iterator e;
1163     CGAL_forall_edges(e,*this->sncp()) {
1164       //    progress++;
1165       CGAL_NEF_TRACEN("");
1166       CGAL_NEF_TRACEN(PH(e));
1167       Halfedge_iterator et = e->twin();
1168       SM_decorator D(&*e->source()), Dt(&*et->source());
1169       CGAL_NEF_TRACEN(e->source()->point());
1170       if ( D.is_isolated(e) ) continue;
1171       SHalfedge_around_svertex_circulator ce(D.first_out_edge(e)),cee(ce);
1172       SHalfedge_around_svertex_circulator cet(Dt.first_out_edge(et)),cete(cet);
1173 
1174       CGAL_For_all(cet,cete) {
1175         //        std::cerr << cet->get_index() << ", " << ce->twin()->get_index() << std::endl;
1176         if (cet->get_forward_index() == ce->twin()->get_backward_index())
1177             //            cet->source()->twin() == ce->source())
1178           break;
1179       }
1180 
1181       CGAL_NEF_TRACEN("vertices " << e->source()->point() <<
1182       "    "      << et->source()->point());
1183 
1184       SHalfedge_around_svertex_circulator sc(D.first_out_edge(e));
1185       SHalfedge_around_svertex_circulator sct(Dt.first_out_edge(et));
1186 
1187       CGAL_NEF_TRACEN("");
1188       CGAL_For_all(sc,cee)
1189         CGAL_NEF_TRACEN("sseg@E addr="<<&*sc<<
1190                         " src="<< sc->source()->point()<<
1191                         " tgt="<< sc->target()->point()<< std::endl <<
1192                         " circle=" << sc->circle()<< std::endl <<
1193                         " indexes=" << sc->get_forward_index() <<
1194                         "," << sc->get_backward_index() << std::endl <<
1195                         "         " << sc->twin()->get_forward_index() <<
1196                         "," << sc->twin()->get_backward_index());
1197 
1198       CGAL_NEF_TRACEN("");
1199 
1200       CGAL_For_all(sct,cete)
1201       CGAL_NEF_TRACEN("sseg@ET addr="<<&*sct<<
1202                       " src="<< sct->source()->point()<<
1203                       " tgt="<<sct->target()->point() <<std::endl<<
1204                       " circle=" << sct->circle() << std::endl<<
1205                       " indexes=" << sct->get_forward_index() <<
1206                       "," << sct->get_backward_index() << std::endl <<
1207                       "         " << sct->twin()->get_forward_index() <<
1208                       "," << sct->twin()->get_backward_index());
1209       CGAL_NEF_TRACEN("");
1210 
1211       CGAL_assertion( cet->get_index() == ce->twin()->get_index());
1212       CGAL_assertion( cet->twin()->get_index() == ce->get_index());
1213       CGAL_assertion( cet->source()->twin() == ce->source());
1214       CGAL_For_all(ce,cee) {
1215         CGAL_NEF_TRACEN("circles " << cet->circle() << "   " << ce->circle() <<
1216                         " sources " << cet->target()->point() <<
1217                         "   " << ce->target()->point());
1218         CGAL_assertion( cet->source()->twin() == ce->source());
1219         CGAL_assertion(ce->mark()==cet->mark());
1220         link_as_prev_next_pair(cet->twin(),ce);
1221         link_as_prev_next_pair(ce->twin(),cet);
1222         --cet; // ce moves ccw, cet moves cw
1223       }
1224     }
1225   }
1226 
categorize_facet_cycles_and_create_facets()1227   void categorize_facet_cycles_and_create_facets() const {
1228   /*{\Mop collects all facet cycles incident to a facet and creates
1229   the facets. \precond |link_shalfedges_to_facet_cycles()| was called
1230   before.}*/
1231 
1232     //    CGAL_NEF_SETDTHREAD(43*31);
1233     CGAL_NEF_TRACEN(">>>>>categorize_facet_cycles_and_create_facets");
1234 
1235     typedef std::list<Object_handle> Object_list;
1236     typedef std::map<int, Object_list>
1237       Map_planes;
1238 
1239     Map_planes M;
1240     SHalfedge_iterator e;
1241     CGAL_forall_shalfedges(e,*this->sncp()) {
1242       if(e->get_index() > e->twin()->get_index())
1243         continue;
1244       M[e->get_index()].push_back(make_object(e));
1245     }
1246     SHalfloop_iterator l;
1247     CGAL_forall_shalfloops(l,*this->sncp()) {
1248       if(l->get_index() > l->twin()->get_index())
1249         continue;
1250       M[l->get_index()].push_back(make_object(l));
1251     }
1252 
1253 #ifdef CGAL_NEF3_TIMER_PLANE_SWEEPS
1254   number_of_plane_sweeps=0;
1255   timer_plane_sweeps.reset();
1256 #endif
1257 
1258     typename Map_planes::iterator it;
1259     CGAL_forall_iterators(it,M) {
1260       CGAL_NEF_TRACEN("  plane "<< it->first);
1261       CGAL_NEF_TRACEN("  size "<< it->second.size());
1262       FM_decorator D(*this->sncp());
1263       Plane_3 h;
1264       Object_handle o(*it->second.begin());
1265       if(CGAL::assign(e, o))
1266         h = e->circle().opposite().plane_through(e->source()->source()->point());
1267       else if(CGAL::assign(l, o))
1268         h = l->circle().opposite().plane_through(l->incident_sface()->center_vertex()->point());
1269       else
1270         CGAL_error_msg( "wrong handle");
1271 
1272       D.create_facet_objects(h,it->second.begin(),it->second.end());
1273     }
1274     //       CGAL_NEF_SETDTHREAD(1);
1275   }
1276 
create_volumes()1277   void create_volumes() {
1278     Base::create_volumes();
1279   }
1280 
build_external_structure()1281   void build_external_structure() {
1282 //    CGAL_NEF_SETDTHREAD(503*509);
1283 
1284 #ifdef CGAL_NEF3_TIMER_EXTERNAL_STRUCTURE
1285     CGAL::Timer timer_external_structure;
1286     timer_external_structure.start();
1287 #endif
1288 
1289 #ifdef CGAL_NEF3_TIMER_PLUECKER
1290     CGAL::Timer timer_pluecker;
1291     timer_pluecker.start();
1292 #endif
1293     //    SNC_io_parser<SNC_structure> O0(std::cerr,*this->sncp());
1294     //    O0.print();
1295     pair_up_halfedges();
1296 #ifdef CGAL_NEF3_TIMER_PLUECKER
1297     timer_pluecker.stop();
1298      if(cgal_nef3_timer_on)
1299       std::cout << "Runtime_pluecker: "
1300                 << timer_pluecker.time() << std::endl;
1301 #endif
1302      //     SNC_io_parser<SNC_structure> O0(std::cerr,*this->sncp());
1303      //     O0.print();
1304     link_shalfedges_to_facet_cycles();
1305 
1306     std::map<int, int> hash;
1307     CGAL::Unique_hash_map<SHalfedge_handle, bool> done(false);
1308 
1309     SHalfedge_iterator sei;
1310     CGAL_forall_shalfedges(sei, *this->sncp()) {
1311       hash[sei->get_index()] = sei->get_index();
1312     }
1313 
1314     CGAL_forall_shalfedges(sei, *this->sncp()) {
1315       if(done[sei])
1316         continue;
1317       SHalfedge_around_facet_circulator circ(sei), end(circ);
1318       int index = circ->get_index();
1319       ++circ;
1320       CGAL_For_all(circ, end)
1321         if(circ->get_index() < index)
1322           index = circ->get_index();
1323       index = hash[index];
1324       CGAL_For_all(circ, end) {
1325         hash[circ->get_index()] = index;
1326         circ->set_index(index);
1327         done[circ] = true;
1328       }
1329     }
1330 
1331     SHalfloop_iterator sli;
1332     CGAL_forall_shalfloops(sli, *this->sncp())
1333       sli->set_index(hash[sli->get_index()]);
1334 
1335     categorize_facet_cycles_and_create_facets();
1336     create_volumes();
1337 
1338 #ifdef CGAL_NEF3_TIMER_EXTERNAL_STRUCTURE
1339     timer_external_structure.stop();
1340     if(cgal_nef3_timer_on)
1341       std::cout << "Runtime_external_structure: "
1342                 << timer_external_structure.time() << std::endl;
1343 #endif
1344   }
1345 
clear_external_structure()1346   void clear_external_structure() {
1347     Base::clear_external_structure();
1348   }
1349 
1350   template<typename Association>
build_after_binary_operation(Association & A)1351   void build_after_binary_operation(Association& A) {
1352 
1353     SHalfedge_iterator sei;
1354     CGAL_forall_shalfedges(sei, *this->sncp()) {
1355       CGAL_NEF_TRACEN("hash sedge " << sei->get_index()
1356                       << "->" << A.get_hash(sei->get_index()));
1357       sei->set_index(A.get_hash(sei->get_index()));
1358     }
1359 
1360     SHalfloop_iterator sli;
1361     CGAL_forall_shalfloops(sli, *this->sncp()) {
1362       CGAL_NEF_TRACEN("hash sloop " << sli->get_index()
1363                       << "->" << A.get_hash(sli->get_index()));
1364       sli->set_index(A.get_hash(sli->get_index()));
1365     }
1366 
1367     //    CGAL_NEF_SETDTHREAD(43);
1368     /*
1369     {    CGAL::SNC_io_parser<SNC_structure> O
1370       (std::cerr, *this->sncp(), false, true);
1371       O.print();}
1372     */
1373     pair_up_halfedges();
1374     /*
1375     {      CGAL::SNC_io_parser<SNC_structure> O
1376         (std::cerr, *this->sncp(), false, true);
1377       O.print();}
1378     */
1379     link_shalfedges_to_facet_cycles();
1380 
1381     SNC_simplify simp(*this->sncp());
1382     simp.vertex_simplificationI();
1383 
1384     //    std::map<int, int> hash;
1385     CGAL::Unique_hash_map<SHalfedge_handle, bool>
1386       done(false);
1387 
1388     /*
1389     SHalfedge_iterator sei;
1390     CGAL_forall_shalfedges(sei, *this->sncp()) {
1391       hash[sei->get_forward_index()] = sei->get_forward_index();
1392       hash[sei->get_backward_index()] = sei->get_backward_index();
1393     }
1394     */
1395 
1396     categorize_facet_cycles_and_create_facets();
1397     create_volumes();
1398   }
1399 };
1400 
1401 } //namespace CGAL
1402 #endif //CGAL_SNC_EXTERNAL_STRUCTURE_H
1403