1 // Boost.Geometry (aka GGL, Generic Geometry Library) // Robustness Test
2 
3 // Copyright (c) 2013-2015 Barend Gehrels, Amsterdam, the Netherlands.
4 
5 // Use, modification and distribution is subject to the Boost Software License,
6 // Version 1.0. (See accompanying file LICENSE_1_0.txt or copy at
7 // http://www.boost.org/LICENSE_1_0.txt)
8 
9 // Adapted from: the attachment of ticket 9081
10 
11 #define CHECK_SELF_INTERSECTIONS
12 #define LIST_WKT
13 
14  #include <iomanip>
15  #include <iostream>
16  #include <vector>
17  #include <boost/geometry.hpp>
18  #include <boost/geometry/geometries/point_xy.hpp>
19  #include <boost/geometry/geometries/polygon.hpp>
20  #include <boost/geometry/geometries/register/point.hpp>
21  #include <boost/geometry/geometries/register/ring.hpp>
22  #include <boost/geometry/io/wkt/wkt.hpp>
23  #include <boost/geometry/geometries/multi_polygon.hpp>
24 
25 #include <boost/foreach.hpp>
26 #include <boost/timer.hpp>
27 #include <boost/algorithm/string.hpp>
28 #include <boost/geometry/io/svg/svg_mapper.hpp>
29 #include <fstream>
30 
31 
32 typedef boost::geometry::model::d2::point_xy<double> pt;
33 typedef boost::geometry::model::polygon<pt> polygon;
34 typedef boost::geometry::model::segment<pt> segment;
35 typedef boost::geometry::model::multi_polygon<polygon> multi_polygon;
36 
37 template <typename Geometry>
debug_with_svg(int index,char method,Geometry const & a,Geometry const & b,std::string const & headera,std::string const & headerb)38 inline void debug_with_svg(int index, char method, Geometry const& a, Geometry const& b, std::string const& headera, std::string const& headerb)
39 {
40     multi_polygon output;
41     try
42     {
43         switch(method)
44         {
45             case 'i': boost::geometry::intersection(a, b, output); break;
46             case 'u': boost::geometry::union_(a, b, output); break;
47             case 'd': boost::geometry::difference(a, b, output); break;
48             case 'v': boost::geometry::difference(b, a, output); break;
49             default : return;
50         }
51     }
52     catch(...)
53     {}
54 
55     std::ostringstream filename;
56     filename << "ticket_9081_" << method << "_" << (1000000 + index) << ".svg";
57     std::ofstream svg(filename.str().c_str());
58 
59     boost::geometry::svg_mapper<pt> mapper(svg, 400, 400);
60     mapper.add(a);
61     mapper.add(b);
62 
63     mapper.map(a, "fill-opacity:0.5;fill:rgb(153,204,0);stroke:rgb(153,204,0);stroke-width:2");
64     mapper.map(b, "fill-opacity:0.3;fill:rgb(51,51,153);stroke:rgb(51,51,153);stroke-width:2");
65     BOOST_FOREACH(polygon const& g, output)
66     {
67         mapper.map(g, "opacity:0.8;fill:none;stroke:rgb(255,128,0);stroke-width:4;stroke-dasharray:1,7;stroke-linecap:round");
68     }
69 
70     std::ostringstream out;
71     out << headera << std::endl << headerb;
72     mapper.map(boost::geometry::return_centroid<pt>(a), "fill:rgb(152,204,0);stroke:rgb(153,204,0);stroke-width:0.1", 3);
73     mapper.map(boost::geometry::return_centroid<pt>(b), "fill:rgb(51,51,153);stroke:rgb(153,204,0);stroke-width:0.1", 3);
74     mapper.text(boost::geometry::return_centroid<pt>(a), headera, "fill:rgb(0,0,0);font-family:Arial;font-size:10px");
75     mapper.text(boost::geometry::return_centroid<pt>(b), headerb, "fill:rgb(0,0,0);font-family:Arial;font-size:10px");
76 }
77 
main()78 int main()
79 {
80     int num_orig = 50;
81     int num_rounds = 30000;
82     srand(1234);
83     std::cout << std::setprecision(16);
84     std::map<int, std::string> genesis;
85     int pj;
86 
87 
88     std::string wkt1, wkt2, operation;
89 
90     try
91     {
92 
93 
94     boost::timer t;
95     std::vector<multi_polygon> poly_list;
96 
97     for(int i=0;i<num_orig;i++)
98     {
99         multi_polygon mp;
100         polygon p;
101         for(int j=0;j<3;j++)
102         {
103             double x=(double)rand()/RAND_MAX;
104             double y=(double)rand()/RAND_MAX;
105             p.outer().push_back(pt(x,y));
106         }
107         boost::geometry::correct(p);
108         mp.push_back(p);
109         boost::geometry::detail::overlay::has_self_intersections(mp);
110 
111         std::ostringstream out;
112         out << "original " << poly_list.size();
113         genesis[poly_list.size()] = out.str();
114         poly_list.push_back(mp);
115 
116 #ifdef LIST_WKT
117         std::cout << "Original " << i << " " << boost::geometry::wkt(p) << std::endl;
118 #endif
119     }
120 
121 
122     for(int j=0;j<num_rounds;j++)
123     {
124         if (j % 100 == 0) { std::cout << " " << j; }
125         pj = j;
126         int a = rand() % poly_list.size();
127         int b = rand() % poly_list.size();
128 
129         debug_with_svg(j, 'i', poly_list[a], poly_list[b], genesis[a], genesis[b]);
130 
131         { std::ostringstream out; out << boost::geometry::wkt(poly_list[a]); wkt1 = out.str(); }
132         { std::ostringstream out; out << boost::geometry::wkt(poly_list[b]); wkt2 = out.str(); }
133 
134         multi_polygon mp_i, mp_u, mp_d, mp_e;
135         operation = "intersection";
136         boost::geometry::intersection(poly_list[a],poly_list[b],mp_i);
137         operation = "intersection";
138         boost::geometry::union_(poly_list[a],poly_list[b],mp_u);
139         operation = "difference";
140         boost::geometry::difference(poly_list[a],poly_list[b],mp_d);
141         boost::geometry::difference(poly_list[b],poly_list[a],mp_e);
142 
143 #ifdef LIST_WKT
144         std::cout << j << std::endl;
145         std::cout << "  Genesis a " << genesis[a] << std::endl;
146         std::cout << "  Genesis b " << genesis[b] << std::endl;
147         std::cout << "  Intersection " << boost::geometry::wkt(mp_i) << std::endl;
148         std::cout << "  Difference a " << boost::geometry::wkt(mp_d) << std::endl;
149         std::cout << "  Difference b " << boost::geometry::wkt(mp_e) << std::endl;
150 #endif
151 
152 #ifdef CHECK_SELF_INTERSECTIONS
153         try
154         {
155             boost::geometry::detail::overlay::has_self_intersections(mp_i);
156         }
157         catch(...)
158         {
159             std::cout << "FAILED TO INTERSECT " << j << std::endl;
160             std::cout << boost::geometry::wkt(poly_list[a]) << std::endl;
161             std::cout << boost::geometry::wkt(poly_list[b]) << std::endl;
162             std::cout << boost::geometry::wkt(mp_i) << std::endl;
163             try
164             {
165                 boost::geometry::detail::overlay::has_self_intersections(mp_i);
166             }
167             catch(...)
168             {
169             }
170             break;
171         }
172 
173         try
174         {
175             boost::geometry::detail::overlay::has_self_intersections(mp_d);
176         }
177         catch(...)
178         {
179             std::cout << "FAILED TO SUBTRACT " << j << std::endl;
180             std::cout << boost::geometry::wkt(poly_list[a]) << std::endl;
181             std::cout << boost::geometry::wkt(poly_list[b]) << std::endl;
182             std::cout << boost::geometry::wkt(mp_d) << std::endl;
183             break;
184         }
185         try
186         {
187             boost::geometry::detail::overlay::has_self_intersections(mp_e);
188         }
189         catch(...)
190         {
191             std::cout << "FAILED TO SUBTRACT " << j << std::endl;
192             std::cout << boost::geometry::wkt(poly_list[b]) << std::endl;
193             std::cout << boost::geometry::wkt(poly_list[a]) << std::endl;
194             std::cout << boost::geometry::wkt(mp_e) << std::endl;
195             break;
196         }
197 #endif
198 
199         if(boost::geometry::area(mp_i) > 0)
200         {
201             std::ostringstream out;
202             out << j << " intersection(" << genesis[a] << " , " << genesis[b] << ")";
203             genesis[poly_list.size()] = out.str();
204             poly_list.push_back(mp_i);
205         }
206         if(boost::geometry::area(mp_d) > 0)
207         {
208             std::ostringstream out;
209             out << j << " difference(" << genesis[a] << " - " << genesis[b] << ")";
210             genesis[poly_list.size()] = out.str();
211             poly_list.push_back(mp_d);
212         }
213         if(boost::geometry::area(mp_e) > 0)
214         {
215             std::ostringstream out;
216             out << j << " difference(" << genesis[b] << " - " << genesis[a] << ")";
217             genesis[poly_list.size()] = out.str();
218             poly_list.push_back(mp_e);
219         }
220     }
221 
222     std::cout << "FINISHED " << t.elapsed() << std::endl;
223 
224     }
225     catch(std::exception const& e)
226     {
227         std::cout << e.what()
228             << " in " << operation << " at " << pj << std::endl
229             << wkt1 << std::endl
230             << wkt2 << std::endl
231             << std::endl;
232     }
233     catch(...)
234     {
235         std::cout << "Other exception" << std::endl;
236     }
237 
238     return 0;
239 }
240