libs/geometry/test/algorithms/overlay/traverse_ccw.cpp

   1 // Boost.Geometry (aka GGL, Generic Geometry Library)
   2 // Unit Test
   3
   4 // Copyright (c) 2010-2012 Barend Gehrels, Amsterdam, the Netherlands.
   5
   6 // Use, modification and distribution is subject to the Boost Software License,
   7 // Version 1.0. (See accompanying file LICENSE_1_0.txt or copy at
   8 // http://www.boost.org/LICENSE_1_0.txt)
   9
  10 #define BOOST_GEOMETRY_DEFINE_STREAM_OPERATOR_SEGMENT_RATIO
  11
  12 #include <fstream>
  13 #include <iostream>
  14 #include <iomanip>
  15
  16 #include <boost/foreach.hpp>
  17
  18 #include <geometry_test_common.hpp>
  19
  20 #include <boost/geometry/geometry.hpp>
  21 #include <boost/geometry/geometries/point_xy.hpp>
  22 #include <boost/geometry/geometries/polygon.hpp>
  23 #include <boost/geometry/io/wkt/wkt.hpp>
  24
  25 #if defined(TEST_WITH_SVG)
  26 #  include <boost/geometry/io/svg/svg_mapper.hpp>
  27 #endif
  28
  29 #include <boost/geometry/algorithms/detail/overlay/debug_turn_info.hpp>
  30 #include <boost/geometry/policies/robustness/get_rescale_policy.hpp>
  31
  32 #include <algorithms/overlay/overlay_cases.hpp>
  33
  34 template <typename Geometry>
  35 struct rev : boost::mpl::if_c<bg::point_order<Geometry>::value == bg::counterclockwise, boost::true_type, boost::false_type>::type
  36 {};
  37
  38 template <typename Geometry1, typename Geometry2>
  39 inline typename bg::coordinate_type<Geometry1>::type
  40 intersect(Geometry1 const& g1, Geometry2 const& g2, std::string const& name,
  41           bg::detail::overlay::operation_type op)
  42 {
  43     boost::ignore_unused(name);
  44
  45     typedef typename bg::strategy::side::services::default_strategy
  46     <
  47         typename bg::cs_tag<Geometry1>::type
  48     >::type side_strategy_type;
  49
  50
  51     typedef typename bg::point_type<Geometry1>::type point_type;
  52
  53     typedef typename bg::rescale_policy_type<point_type>::type
  54         rescale_policy_type;
  55
  56     rescale_policy_type rescale_policy
  57             = bg::get_rescale_policy<rescale_policy_type>(g1, g2);
  58
  59     typedef bg::detail::overlay::traversal_turn_info
  60     <
  61         point_type,
  62         typename bg::detail::segment_ratio_type<point_type, rescale_policy_type>::type
  63     > turn_info;
  64     std::vector<turn_info> turns;
  65
  66     bg::detail::get_turns::no_interrupt_policy policy;
  67     bg::get_turns
  68         <
  69             rev<Geometry1>::value,
  70             rev<Geometry2>::value,
  71             bg::detail::overlay::assign_null_policy
  72         >(g1, g2, rescale_policy, turns, policy);
  73
  74     bg::enrich_intersection_points
  75             <
  76                 rev<Geometry1>::value, rev<Geometry2>::value,
  77                 bg::overlay_intersection
  78             >(turns, bg::detail::overlay::operation_intersection,
  79         g1, g2, rescale_policy, side_strategy_type());
  80
  81     typedef bg::model::ring<typename bg::point_type<Geometry1>::type> ring_type;
  82     typedef std::deque<ring_type> out_vector;
  83     out_vector v;
  84
  85     bg::detail::overlay::traverse
  86         <
  87             rev<Geometry1>::value, rev<Geometry2>::value,
  88             Geometry1, Geometry2
  89         >::apply(g1, g2, op, rescale_policy, turns, v);
  90
  91     typename bg::coordinate_type<Geometry1>::type result = 0.0;
  92     BOOST_FOREACH(ring_type& ring, v)
  93     {
  94         result += bg::area(ring);
  95     }
  96
  97 #if defined(TEST_WITH_SVG)
  98     {
  99         std::ostringstream filename;
 100         filename
 101             << name << "_"
 102             << (op == bg::detail::overlay::operation_intersection ? "i" : "u")
 103             << "_" << (rev<Geometry1>::value ? "ccw" : "cw")
 104             << "_" << (rev<Geometry2>::value ? "ccw" : "cw")
 105             << ".svg";
 106
 107         std::ofstream svg(filename.str().c_str());
 108
 109         bg::svg_mapper<typename bg::point_type<Geometry1>::type> mapper(svg, 500, 500);
 110         mapper.add(g1);
 111         mapper.add(g2);
 112
 113         // Input shapes in green/blue
 114         mapper.map(g1, "fill-opacity:0.5;fill:rgb(153,204,0);"
 115                 "stroke:rgb(153,204,0);stroke-width:3");
 116         mapper.map(g2, "fill-opacity:0.3;fill:rgb(51,51,153);"
 117                 "stroke:rgb(51,51,153);stroke-width:3");
 118
 119         // Traversal rings in magenta/light yellow fill
 120         BOOST_FOREACH(ring_type const& ring, v)
 121         {
 122             mapper.map(ring, "fill-opacity:0.3;stroke-opacity:0.4;fill:rgb(255,255,0);"
 123                     "stroke:rgb(255,0,255);stroke-width:8");
 124         }
 125
 126
 127         // turn points in orange, + enrichment/traversal info
 128
 129         // Simple map to avoid two texts at same place (note that can still overlap!)
 130         std::map<std::pair<int, int>, int> offsets;
 131         int index = 0;
 132         int const lineheight = 10;
 133         int const margin = 5;
 134
 135         BOOST_FOREACH(turn_info const& turn, turns)
 136         {
 137             mapper.map(turn.point, "fill:rgb(255,128,0);"
 138                     "stroke:rgb(0,0,0);stroke-width:1", 3);
 139
 140             {
 141                 // Map characteristics
 142                 // Create a rounded off point
 143                 std::pair<int, int> p
 144                     = std::make_pair(
 145                         boost::numeric_cast<int>(0.5 + 10 * bg::get<0>(turn.point)),
 146                         boost::numeric_cast<int>(0.5 + 10 * bg::get<1>(turn.point))
 147                         );
 148                 std::string style =  "fill:rgb(0,0,0);font-family:Arial;font-size:10px";
 149
 150                 std::ostringstream out;
 151                 out << index
 152                     << ": " << bg::method_char(turn.method)
 153                     << std::endl
 154                     << "op: " << bg::operation_char(turn.operations[0].operation)
 155                     << " / " << bg::operation_char(turn.operations[1].operation)
 156                     << (turn.is_discarded() ? " (discarded) " : turn.blocked() ? " (blocked)" : "")
 157                     << std::endl;
 158
 159                 if (turn.operations[0].enriched.next_ip_index != -1)
 160                 {
 161                     out << "ip: " << turn.operations[0].enriched.next_ip_index;
 162                 }
 163                 else
 164                 {
 165                     out << "vx: " << turn.operations[0].enriched.travels_to_vertex_index;
 166                 }
 167                 out << " ";
 168                 if (turn.operations[1].enriched.next_ip_index != -1)
 169                 {
 170                     out << "ip: " << turn.operations[1].enriched.next_ip_index;
 171                 }
 172                 else
 173                 {
 174                     out << "vx: " << turn.operations[1].enriched.travels_to_vertex_index;
 175                 }
 176
 177                 out << std::endl;
 178
 179                 out
 180
 181                     << std::setprecision(3)
 182                     << "dist: " << turn.operations[0].fraction
 183                     << " / "  << turn.operations[1].fraction
 184                     << std::endl;
 185
 186                 offsets[p] += lineheight;
 187                 int offset = offsets[p];
 188                 offsets[p] += lineheight * 5;
 189                 mapper.text(turn.point, out.str(), style, margin, offset, lineheight);
 190             }
 191             index++;
 192         }
 193     }
 194 #endif
 195
 196     return result;
 197 }
 198
 199 template <typename Geometry1, typename Geometry2>
 200 inline typename bg::coordinate_type<Geometry1>::type intersect(std::string const& wkt1, std::string const& wkt2, std::string const& name,
 201             bg::detail::overlay::operation_type op)
 202 {
 203     Geometry1 geometry1;
 204     Geometry2 geometry2;
 205     bg::read_wkt(wkt1, geometry1);
 206     bg::read_wkt(wkt2, geometry2);
 207
 208     // Reverse if necessary: adapt to cw/ccw
 209     bg::correct(geometry1);
 210     bg::correct(geometry2);
 211
 212     return intersect(geometry1, geometry2, name, op);
 213 }
 214
 215 template <typename T>
 216 inline void test_polygon(std::string const& wkt1, std::string const& wkt2, std::string const& name)
 217 {
 218     typedef bg::model::d2::point_xy<T> point;
 219     typedef bg::model::polygon<point> clock;
 220     typedef bg::model::polygon<point, false> counter;
 221
 222     namespace ov = bg::detail::overlay;
 223     T area1 = intersect<clock, clock>(wkt1, wkt2, name, ov::operation_intersection);
 224     T area2 = intersect<counter, counter>(wkt1, wkt2, name, ov::operation_intersection);
 225     T area3 = intersect<clock, counter>(wkt1, wkt2, name, ov::operation_intersection);
 226     T area4 = intersect<counter, clock>(wkt1, wkt2, name, ov::operation_intersection);
 227     BOOST_CHECK_CLOSE(area1, area2, 0.001);
 228     BOOST_CHECK_CLOSE(area3, area4, 0.001);
 229     BOOST_CHECK_CLOSE(area1, area3, 0.001);
 230     BOOST_CHECK_CLOSE(area2, area4, 0.001);
 231
 232     area1 = intersect<clock, clock>(wkt1, wkt2, name, ov::operation_union);
 233     area2 = intersect<counter, counter>(wkt1, wkt2, name, ov::operation_union);
 234     area3 = intersect<clock, counter>(wkt1, wkt2, name, ov::operation_union);
 235     area4 = intersect<counter, clock>(wkt1, wkt2, name, ov::operation_union);
 236     BOOST_CHECK_CLOSE(area1, area2, 0.001);
 237     BOOST_CHECK_CLOSE(area3, area4, 0.001);
 238     BOOST_CHECK_CLOSE(area1, area3, 0.001);
 239     BOOST_CHECK_CLOSE(area2, area4, 0.001);
 240 }
 241
 242 template <typename T>
 243 inline void test_box_polygon(std::string const& wkt1, std::string const& wkt2, std::string const& name)
 244 {
 245     typedef bg::model::d2::point_xy<T> point;
 246     typedef bg::model::box<point> box;
 247     typedef bg::model::polygon<point> clock;
 248     typedef bg::model::polygon<point, false> counter;
 249
 250     namespace ov = bg::detail::overlay;
 251     T area1 = intersect<box, clock>(wkt1, wkt2, name + "_bp", ov::operation_intersection);
 252     T area2 = intersect<box, counter>(wkt1, wkt2, name + "_bp", ov::operation_intersection);
 253     T area3 = intersect<clock, box>(wkt2, wkt1, name + "_pb", ov::operation_intersection);
 254     T area4 = intersect<counter, box>(wkt2, wkt1, name + "_pb", ov::operation_intersection);
 255     BOOST_CHECK_CLOSE(area1, area2, 0.001);
 256     BOOST_CHECK_CLOSE(area3, area4, 0.001);
 257     BOOST_CHECK_CLOSE(area1, area3, 0.001);
 258     BOOST_CHECK_CLOSE(area2, area4, 0.001);
 259
 260     area1 = intersect<box, clock>(wkt1, wkt2, name + "_bp", ov::operation_union);
 261     area2 = intersect<box, counter>(wkt1, wkt2, name + "_bp", ov::operation_union);
 262     area3 = intersect<clock, box>(wkt2, wkt1, name + "_pb", ov::operation_union);
 263     area4 = intersect<counter, box>(wkt2, wkt1, name + "_pb", ov::operation_union);
 264     BOOST_CHECK_CLOSE(area1, area2, 0.001);
 265     BOOST_CHECK_CLOSE(area3, area4, 0.001);
 266     BOOST_CHECK_CLOSE(area1, area3, 0.001);
 267     BOOST_CHECK_CLOSE(area2, area4, 0.001);
 268 }
 269
 270 int test_main(int, char* [])
 271 {
 272     //bool const ig = true;
 273     test_polygon<double>(case_1[0], case_1[1], "c1");
 274     test_polygon<double>(case_2[0], case_2[1], "c2");
 275     test_polygon<double>(case_3[0], case_3[1], "c3");
 276     test_polygon<double>(case_4[0], case_4[1], "c4");
 277     test_polygon<double>(case_5[0], case_5[1], "c5");
 278     test_polygon<double>(case_6[0], case_6[1], "c6");
 279     test_polygon<double>(case_7[0], case_7[1], "c7");
 280     test_polygon<double>(case_8[0], case_8[1], "c8");
 281     test_polygon<double>(case_9[0], case_9[1], "c9" /*, ig */);
 282
 283
 284     test_polygon<double>(case_10[0], case_10[1], "c10");
 285     test_polygon<double>(case_11[0], case_11[1], "c11");
 286     test_polygon<double>(case_12[0], case_12[1], "c12");
 287     test_polygon<double>(case_13[0], case_13[1], "c13");
 288     test_polygon<double>(case_14[0], case_14[1], "c14");
 289     test_polygon<double>(case_15[0], case_15[1], "c15");
 290     test_polygon<double>(case_16[0], case_16[1], "c16");
 291     test_polygon<double>(case_17[0], case_17[1], "c17");
 292     test_polygon<double>(case_18[0], case_18[1], "c18");
 293     test_polygon<double>(case_19[0], case_19[1], "c19");
 294     test_polygon<double>(case_20[0], case_20[1], "c20");
 295     test_polygon<double>(case_21[0], case_21[1], "c21");
 296     test_polygon<double>(case_22[0], case_22[1], "c22" /*, ig */);
 297     test_polygon<double>(case_23[0], case_23[1], "c23");
 298     test_polygon<double>(case_24[0], case_24[1], "c24");
 299     test_polygon<double>(case_25[0], case_25[1], "c25" /*, ig */);
 300     test_polygon<double>(case_26[0], case_26[1], "c26" /*, ig */);
 301     test_polygon<double>(case_27[0], case_27[1], "c27");
 302     test_polygon<double>(case_28[0], case_28[1], "c28");
 303     test_polygon<double>(case_29[0], case_29[1], "c29");
 304     test_polygon<double>(case_30[0], case_30[1], "c30");
 305     test_polygon<double>(case_31[0], case_31[1], "c31" /*, ig */);
 306     test_polygon<double>(case_32[0], case_32[1], "c32" /*, ig */);
 307     test_polygon<double>(case_33[0], case_33[1], "c33" /*, ig */);
 308     test_polygon<double>(case_34[0], case_34[1], "c34");
 309     test_polygon<double>(case_35[0], case_35[1], "c35");
 310     test_polygon<double>(case_36[0], case_36[1], "c36" /*, ig */);
 311     test_polygon<double>(case_37[0], case_37[1], "c37" /*, ig */);
 312     test_polygon<double>(case_38[0], case_38[1], "c38" /*, ig */);
 313     test_polygon<double>(case_39[0], case_39[1], "c39");
 314     test_polygon<double>(case_40[0], case_40[1], "c40" /*, ig */);
 315     test_polygon<double>(case_41[0], case_41[1], "c41");
 316     test_polygon<double>(case_42[0], case_42[1], "c42");
 317     //test_polygon<double>(case_43[0], case_43[1], "c43", inv);
 318     test_polygon<double>(case_44[0], case_44[1], "c44");
 319     test_polygon<double>(case_45[0], case_45[1], "c45");
 320     //test_polygon<double>(case_46[0], case_46[1], "c46", inv);
 321     //test_polygon<double>(case_47[0], case_47[1], "c47" /*, ig */);
 322     //test_polygon<double>(case_48[0], case_48[1], "c48");
 323     test_polygon<double>(case_49[0], case_49[1], "c49");
 324     test_polygon<double>(case_50[0], case_50[1], "c50");
 325     test_polygon<double>(case_51[0], case_51[1], "c51");
 326     test_polygon<double>(case_52[0], case_52[1], "c52" /*, ig */);
 327     test_polygon<double>(case_53[0], case_53[1], "c53");
 328     // Invalid ones / overlaying intersection points / self tangencies
 329     //test_polygon<double>(case_54[0], case_54[1], "c54");
 330     //test_polygon<double>(case_55[0], case_55[1], "c55");
 331     //test_polygon<double>(case_56[0], case_56[1], "c56");
 332     //test_polygon<double>(case_57[0], case_57[1], "c57" /*, ig */);
 333     //test_polygon<double>(case_58[0], case_58[1], "c58");
 334     //test_polygon<double>(case_59[0], case_59[1], "c59");
 335
 336     test_polygon<double>(pie_16_4_12[0], pie_16_4_12[1], "pie_16_4_12");
 337     test_polygon<double>(pie_23_21_12_500[0], pie_23_21_12_500[1], "pie_23_21_12_500");
 338     test_polygon<double>(pie_23_23_3_2000[0], pie_23_23_3_2000[1], "pie_23_23_3_2000");
 339     test_polygon<double>(pie_23_16_16[0], pie_23_16_16[1], "pie_23_16_16");
 340     test_polygon<double>(pie_16_2_15_0[0], pie_16_2_15_0[1], "pie_16_2_15_0");
 341     test_polygon<double>(pie_4_13_15[0], pie_4_13_15[1], "pie_4_13_15");
 342     test_polygon<double>(pie_20_20_7_100[0], pie_20_20_7_100[1], "pie_20_20_7_100");
 343
 344     test_polygon<double>(hv_1[0], hv_1[1], "hv1");
 345     test_polygon<double>(hv_2[0], hv_2[1], "hv2");
 346     test_polygon<double>(hv_3[0], hv_3[1], "hv3");
 347     test_polygon<double>(hv_4[0], hv_4[1], "hv4");
 348     test_polygon<double>(hv_5[0], hv_5[1], "hv5");
 349     test_polygon<double>(hv_6[0], hv_6[1], "hv6");
 350     test_polygon<double>(hv_7[0], hv_7[1], "hv7");
 351     test_polygon<double>(dz_1[0], dz_1[1], "dz_1");
 352     test_polygon<double>(dz_2[0], dz_2[1], "dz_2");
 353     test_polygon<double>(dz_3[0], dz_3[1], "dz_3");
 354
 355     test_box_polygon<double>("POLYGON((1 1,4 4))", case_1[0], "bp1");
 356
 357     {
 358         static std::string example_box = "POLYGON((1.5 1.5, 4.5 2.5))";
 359         static std::string example_ring =
 360             "POLYGON((2 1.3,2.4 1.7,2.8 1.8,3.4 1.2,3.7 1.6,3.4 2,4.1 3,5.3 2.6,5.4 1.2,4.9 0.8,2.9 0.7,2 1.3))";
 361         test_box_polygon<double>(example_box, example_ring, "bp2");
 362     }
 363
 364     return 0;
 365 }