1 #ifndef CGAL_SDG_VERBOSE
2 #define CGAL_SDG_DEBUG(a)
3 #else
4 #define CGAL_SDG_DEBUG(a) { a }
5 #endif
6 
7 // Kernels
8 //#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
9 #include <CGAL/Exact_predicates_exact_constructions_kernel.h>
10 //#include <CGAL/Cartesian.h>
11 
12 #include <CGAL/CGAL_Ipelet_base.h>
13 #include <CGAL/Segment_Delaunay_graph_Linf_2.h>
14 #include <CGAL/Segment_Delaunay_graph_Linf_traits_2.h>
15 
16 namespace CGAL_svdlinf {
17 
18 //typedef CGAL::Cartesian<double>                               Kernel;
19   typedef CGAL::Exact_predicates_exact_constructions_kernel     Kernel;
20 //typedef CGAL::Exact_predicates_inexact_constructions_kernel   Kernel;
21   typedef CGAL::Segment_Delaunay_graph_Linf_traits_2<Kernel>    Gt;
22   typedef CGAL::Segment_Delaunay_graph_Linf_2<Gt>               SDG2;
23 
24   const unsigned int num_entries = 3;
25 
26   const std::string sublabel[] = {
27     "Segment VD Linf general",
28     "Segment skeleton Linf general",
29     "Help"
30   };
31 
32   const std::string helpmsg[] = {
33     "Draw the L_inf Voronoi diagram of segments in Linf",
34     "Draw the L_inf Voronoi skeleton of segments in Linf",
35   };
36 
37   class svdlinfIpelet
38     : public CGAL::Ipelet_base<Kernel,num_entries> {
39       public:
svdlinfIpelet()40         svdlinfIpelet()
41           :CGAL::Ipelet_base<Kernel,num_entries>
42              ("SVDLinf",sublabel,helpmsg){}
43         void protected_run(int);
44 
45     };
46   // --------------------------------------------------------------------
47 
protected_run(int fn)48   void svdlinfIpelet::protected_run(int fn)
49   {
50     SDG2   svd;
51 
52     if (fn == (num_entries-1)) {
53       show_help();
54       return;
55     }
56 
57     std::list<Point_2> pt_list;
58     std::list<Segment_2> sg_list;
59 
60     Iso_rectangle_2 bbox;
61 
62     // grab input
63 
64     bbox =
65       read_active_objects(
66           CGAL::dispatch_or_drop_output
67           <Point_2,Polygon_2,Segment_2>(
68           std::back_inserter(pt_list),
69           segment_grabber(std::back_inserter(sg_list)),
70           std::back_inserter(sg_list)
71           )
72           );
73 
74     // check input
75 
76     if (pt_list.empty() and sg_list.empty()) {
77       print_error_message(("Nothing selected"));
78       return;
79     }
80 
81     if ( (fn == 2) or (fn == 3) ) {
82       // check that segments are all axis-parallel
83       for (std::list<Segment_2>::iterator
84           sit  = sg_list.begin();
85           sit != sg_list.end();
86           ++sit)
87       {
88         if (not (sit->is_horizontal() or sit->is_vertical())) {
89           print_error_message(("Non axis-parallel segment"));
90           return;
91         }
92       }
93     }
94 
95 
96     Kernel::FT incr_len = 75;
97     // slightly increase the size of the bbox
98     bbox = Iso_rectangle_2(
99       (bbox.min)()+Kernel::Vector_2(-incr_len,-incr_len),
100       (bbox.max)()+Kernel::Vector_2(incr_len,incr_len));
101 
102     for (std::list<Segment_2>::iterator
103          sit  = sg_list.begin();
104          sit != sg_list.end();
105          ++sit)
106     {
107       CGAL_SDG_DEBUG( std::cout << "IPELET: inserting segment "
108           << *sit << std::endl ; );
109       svd.insert(sit->point(0),sit->point(1));
110     }
111 
112     CGAL_SDG_DEBUG( std::cout << "IPELET: inserting points"
113         << std::endl ; );
114 
115     svd.insert(pt_list.begin(),pt_list.end());
116 
117     if ( fn == 0 ) {
118       draw_dual_in_ipe(svd, bbox);
119     } else if ( fn == 1 ) {
120       draw_skeleton_in_ipe(svd, bbox);
121     }
122 
123   } // end of void svdlinfIpelet::protected_run(int fn)
124 
125 }
126 
127 CGAL_IPELET(CGAL_svdlinf::svdlinfIpelet)
128