1 #include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
2
3 #include <CGAL/Polygon_2.h>
4 #include <CGAL/Polygon_2_algorithms.h>
5 #include <CGAL/Straight_skeleton_builder_2.h>
6 #include <CGAL/Polygon_offset_builder_2.h>
7 #include <CGAL/compute_outer_frame_margin.h>
8 #include "print.h"
9
10 #include <boost/shared_ptr.hpp>
11
12 #include <vector>
13 #include <cassert>
14
15 //
16 // This example illustrates how to use the CGAL Straight Skeleton package
17 // to construct an offset contour on the outside of a polygon
18 //
19
20 // This is the recommended kernel
21 typedef CGAL::Exact_predicates_inexact_constructions_kernel Kernel;
22
23 typedef Kernel::Point_2 Point_2;
24 typedef CGAL::Polygon_2<Kernel> Contour;
25 typedef boost::shared_ptr<Contour> ContourPtr;
26 typedef std::vector<ContourPtr> ContourSequence ;
27
28 typedef CGAL::Straight_skeleton_2<Kernel> Ss;
29
30 typedef Ss::Halfedge_iterator Halfedge_iterator;
31 typedef Ss::Halfedge_handle Halfedge_handle;
32 typedef Ss::Vertex_handle Vertex_handle;
33
34 typedef CGAL::Straight_skeleton_builder_traits_2<Kernel> SsBuilderTraits;
35 typedef CGAL::Straight_skeleton_builder_2<SsBuilderTraits,Ss> SsBuilder;
36
37 typedef CGAL::Polygon_offset_builder_traits_2<Kernel> OffsetBuilderTraits;
38 typedef CGAL::Polygon_offset_builder_2<Ss,OffsetBuilderTraits,Contour> OffsetBuilder;
39
main()40 int main()
41 {
42 // A start-shaped polygon, oriented counter-clockwise as required for outer contours.
43 Point_2 pts[] = { Point_2(-1,-1)
44 , Point_2(0,-12)
45 , Point_2(1,-1)
46 , Point_2(12,0)
47 , Point_2(1,1)
48 , Point_2(0,12)
49 , Point_2(-1,1)
50 , Point_2(-12,0)
51 } ;
52
53 std::vector<Point_2> star(pts,pts+8);
54
55 assert(CGAL::orientation_2(pts,pts+8,Kernel()) == CGAL::COUNTERCLOCKWISE);
56
57 // We want an offset contour in the outside.
58 // Since the package doesn't support that operation directly, we use the following trick:
59 // (1) Place the polygon as a hole of a big outer frame.
60 // (2) Construct the skeleton on the interior of that frame (with the polygon as a hole)
61 // (3) Construc the offset contours
62 // (4) Identify the offset contour that corresponds to the frame and remove it from the result
63
64
65 double offset = 3 ; // The offset distance
66
67 // First we need to determine the proper separation between the polygon and the frame.
68 // We use this helper function provided in the package.
69 boost::optional<double> margin = CGAL::compute_outer_frame_margin(star.begin(),star.end(),offset);
70
71 // Proceed only if the margin was computed (an extremely sharp corner might cause overflow)
72 if ( margin )
73 {
74 // Get the bbox of the polygon
75 CGAL::Bbox_2 bbox = CGAL::bbox_2(star.begin(),star.end());
76
77 // Compute the boundaries of the frame
78 double fxmin = bbox.xmin() - *margin ;
79 double fxmax = bbox.xmax() + *margin ;
80 double fymin = bbox.ymin() - *margin ;
81 double fymax = bbox.ymax() + *margin ;
82
83 // Create the rectangular frame
84 Point_2 frame[4]= { Point_2(fxmin,fymin)
85 , Point_2(fxmax,fymin)
86 , Point_2(fxmax,fymax)
87 , Point_2(fxmin,fymax)
88 } ;
89
90 // Instantiate the skeleton builder
91 SsBuilder ssb ;
92
93 // Enter the frame
94 ssb.enter_contour(frame,frame+4);
95
96 // Enter the polygon as a hole of the frame (NOTE: as it is a hole we insert it in the opposite orientation)
97 ssb.enter_contour(star.rbegin(),star.rend());
98
99 // Construct the skeleton
100 boost::shared_ptr<Ss> ss = ssb.construct_skeleton();
101
102 // Proceed only if the skeleton was correctly constructed.
103 if ( ss )
104 {
105 print_straight_skeleton(*ss);
106
107 // Instantiate the container of offset contours
108 ContourSequence offset_contours ;
109
110 // Instantiate the offset builder with the skeleton
111 OffsetBuilder ob(*ss);
112
113 // Obtain the offset contours
114 ob.construct_offset_contours(offset, std::back_inserter(offset_contours));
115
116 // Locate the offset contour that corresponds to the frame
117 // That must be the outmost offset contour, which in turn must be the one
118 // with the largetst unsigned area.
119 ContourSequence::iterator f = offset_contours.end();
120 double lLargestArea = 0.0 ;
121 for (ContourSequence::iterator i = offset_contours.begin(); i != offset_contours.end(); ++ i )
122 {
123 double lArea = CGAL_NTS abs( (*i)->area() ) ; //Take abs() as Polygon_2::area() is signed.
124 if ( lArea > lLargestArea )
125 {
126 f = i ;
127 lLargestArea = lArea ;
128 }
129 }
130
131 // Remove the offset contour that corresponds to the frame.
132 offset_contours.erase(f);
133
134 print_polygons(offset_contours);
135 }
136 }
137
138 return 0;
139 }
140