ceph/src/boost/libs/geometry/test/robustness/overlay/buffer/many_ring_buffer.cpp

   1 // Boost.Geometry (aka GGL, Generic Geometry Library)
   2 // Unit Test
   3
   4 // Copyright (c) 2012-2015 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_BUFFER_TEST_SVG_USE_ALTERNATE_BOX_FOR_INPUT
  11 #define BOOST_GEOMETRY_BUFFER_TEST_SVG_ALTERNATE_BOX "BOX(179 4, 180 5)"
  12
  13
  14 #include <test_buffer.hpp>
  15
  16 #include <boost/geometry/algorithms/difference.hpp>
  17 #include <boost/geometry/geometries/geometries.hpp>
  18
  19 #include <boost/random/linear_congruential.hpp>
  20 #include <boost/random/uniform_int.hpp>
  21 #include <boost/random/uniform_real.hpp>
  22 #include <boost/random/variate_generator.hpp>
  23
  24
  25
  26 const int point_count = 90; // points for a full circle
  27
  28 // Function to let buffer-distance depend on alpha, e.g.:
  29 inline double corrected_distance(double distance, double alpha)
  30 {
  31     return distance * 1.0 + 0.2 * sin(alpha * 6.0);
  32 }
  33
  34 class buffer_point_strategy_sample
  35 {
  36 public :
  37
  38     template
  39     <
  40         typename Point,
  41         typename OutputRange,
  42         typename DistanceStrategy
  43     >
  44     void apply(Point const& point,
  45             DistanceStrategy const& distance_strategy,
  46             OutputRange& output_range) const
  47     {
  48         double const distance = distance_strategy.apply(point, point,
  49                         bg::strategy::buffer::buffer_side_left);
  50
  51         double const angle_increment = 2.0 * M_PI / double(point_count);
  52
  53         double alpha = 0;
  54         for (std::size_t i = 0; i <= point_count; i++, alpha -= angle_increment)
  55         {
  56             double const cd = corrected_distance(distance, alpha);
  57
  58             typename boost::range_value<OutputRange>::type output_point;
  59             bg::set<0>(output_point, bg::get<0>(point) + cd * cos(alpha));
  60             bg::set<1>(output_point, bg::get<1>(point) + cd * sin(alpha));
  61             output_range.push_back(output_point);
  62         }
  63     }
  64 };
  65
  66 class buffer_join_strategy_sample
  67 {
  68 private :
  69     template
  70     <
  71         typename Point,
  72         typename DistanceType,
  73         typename RangeOut
  74     >
  75     inline void generate_points(Point const& vertex,
  76                 Point const& perp1, Point const& perp2,
  77                 DistanceType const& buffer_distance,
  78                 RangeOut& range_out) const
  79     {
  80         double dx1 = bg::get<0>(perp1) - bg::get<0>(vertex);
  81         double dy1 = bg::get<1>(perp1) - bg::get<1>(vertex);
  82         double dx2 = bg::get<0>(perp2) - bg::get<0>(vertex);
  83         double dy2 = bg::get<1>(perp2) - bg::get<1>(vertex);
  84
  85         // Assuming the corner is convex, angle2 < angle1
  86         double const angle1 = atan2(dy1, dx1);
  87         double angle2 = atan2(dy2, dx2);
  88
  89         while (angle2 > angle1)
  90         {
  91             angle2 -= 2 * M_PI;
  92         }
  93
  94         double const angle_increment = 2.0 * M_PI / double(point_count);
  95         double alpha = angle1 - angle_increment;
  96
  97         for (int i = 0; alpha >= angle2 && i < point_count; i++, alpha -= angle_increment)
  98         {
  99             double cd = corrected_distance(buffer_distance, alpha);
 100
 101             Point p;
 102             bg::set<0>(p, bg::get<0>(vertex) + cd * cos(alpha));
 103             bg::set<1>(p, bg::get<1>(vertex) + cd * sin(alpha));
 104             range_out.push_back(p);
 105         }
 106     }
 107
 108 public :
 109
 110     template <typename Point, typename DistanceType, typename RangeOut>
 111     inline bool apply(Point const& ip, Point const& vertex,
 112                 Point const& perp1, Point const& perp2,
 113                 DistanceType const& buffer_distance,
 114                 RangeOut& range_out) const
 115     {
 116         generate_points(vertex, perp1, perp2, buffer_distance, range_out);
 117         return true;
 118     }
 119
 120     template <typename NumericType>
 121     static inline NumericType max_distance(NumericType const& distance)
 122     {
 123         return distance;
 124     }
 125
 126 };
 127
 128 class buffer_side_sample
 129 {
 130 public :
 131     template
 132     <
 133         typename Point,
 134         typename OutputRange,
 135         typename DistanceStrategy
 136     >
 137     static inline void apply(
 138                 Point const& input_p1, Point const& input_p2,
 139                 bg::strategy::buffer::buffer_side_selector side,
 140                 DistanceStrategy const& distance,
 141                 OutputRange& output_range)
 142     {
 143         // Generate a block along (left or right of) the segment
 144
 145         double const dx = bg::get<0>(input_p2) - bg::get<0>(input_p1);
 146         double const dy = bg::get<1>(input_p2) - bg::get<1>(input_p1);
 147
 148         // For normalization [0,1] (=dot product d.d, sqrt)
 149         double const length = bg::math::sqrt(dx * dx + dy * dy);
 150
 151         if (bg::math::equals(length, 0))
 152         {
 153             return;
 154         }
 155
 156         // Generate the normalized perpendicular p, to the left (ccw)
 157         double const px = -dy / length;
 158         double const py = dx / length;
 159
 160         // Both vectors perpendicular to input p1 and input p2 have same angle
 161         double const alpha = atan2(py, px);
 162
 163         double const d = distance.apply(input_p1, input_p2, side);
 164
 165         double const cd = corrected_distance(d, alpha);
 166
 167         output_range.resize(2);
 168
 169         bg::set<0>(output_range.front(), bg::get<0>(input_p1) + px * cd);
 170         bg::set<1>(output_range.front(), bg::get<1>(input_p1) + py * cd);
 171         bg::set<0>(output_range.back(), bg::get<0>(input_p2) + px * cd);
 172         bg::set<1>(output_range.back(), bg::get<1>(input_p2) + py * cd);
 173     }
 174 };
 175
 176 #ifdef TEST_WITH_SVG
 177 template <typename Geometry1, typename Geometry2>
 178 void create_svg(std::string const& filename, Geometry1 const& original, Geometry2 const& buffer, std::string const& color)
 179 {
 180     typedef typename bg::point_type<Geometry1>::type point_type;
 181     std::ofstream svg(filename.c_str());
 182
 183     bg::svg_mapper<point_type> mapper(svg, 800, 800);
 184     mapper.add(buffer);
 185
 186     mapper.map(original, "fill-opacity:0.3;fill:rgb(255,0,0);stroke:rgb(0,0,0);stroke-width:1");
 187
 188     std::string style = "fill-opacity:0.3;fill:";
 189     style += color;
 190     style += ";stroke:rgb(0,0,0);stroke-width:1";
 191     mapper.map(buffer, style);
 192 }
 193 #endif
 194
 195
 196 void test_many_rings(int imn, int jmx, int count,
 197     double expected_area_exterior,
 198     double expected_area_interior)
 199 {
 200     typedef bg::model::point<double, 2, bg::cs::cartesian> point;
 201     typedef bg::model::polygon<point> polygon_type;
 202     typedef bg::model::multi_polygon<polygon_type> multi_polygon_type;
 203
 204     // Predefined strategies
 205     bg::strategy::buffer::distance_symmetric<double> distance_strategy(1.3);
 206     bg::strategy::buffer::end_flat end_strategy; // not effectively used
 207
 208     // Own strategies
 209     buffer_join_strategy_sample join_strategy;
 210     buffer_point_strategy_sample point_strategy;
 211     buffer_side_sample side_strategy;
 212
 213     // Declare output
 214
 215     bg::model::multi_point<point> mp;
 216
 217     // Use a bit of random disturbance in the otherwise too regular grid
 218     typedef boost::minstd_rand base_generator_type;
 219     base_generator_type generator(12345);
 220     boost::uniform_real<> random_range(0.0, 0.5);
 221     boost::variate_generator
 222     <
 223         base_generator_type&,
 224         boost::uniform_real<>
 225     > random(generator, random_range);
 226
 227     for (int i = 0; i < count; i++)
 228     {
 229         for (int j = 0; j < count; j++)
 230         {
 231             double x = i * 3.0 + random();
 232             double y = j * 3.0 + random();
 233             //if (i > 30 && j < 30)
 234             if (i > imn && j < jmx)
 235             {
 236                 point p(x, y);
 237                 mp.push_back(p);
 238             }
 239         }
 240     }
 241
 242     multi_polygon_type many_rings;
 243     // Create the buffer of a multi-point
 244     bg::buffer(mp, many_rings,
 245                 distance_strategy, side_strategy,
 246                 join_strategy, end_strategy, point_strategy);
 247
 248     bg::model::box<point> envelope;
 249     bg::envelope(many_rings, envelope);
 250     bg::buffer(envelope, envelope, 1.0);
 251
 252     multi_polygon_type many_interiors;
 253     bg::difference(envelope, many_rings, many_interiors);
 254
 255 #ifdef TEST_WITH_SVG
 256     create_svg("/tmp/many_interiors.svg", mp, many_interiors, "rgb(51,51,153)");
 257     create_svg("/tmp/buffer.svg", mp, many_rings, "rgb(51,51,153)");
 258 #endif
 259
 260     bg::strategy::buffer::join_round join_round(100);
 261     bg::strategy::buffer::end_flat end_flat;
 262
 263     {
 264         std::ostringstream out;
 265         out << "many_rings_" << count;
 266         out << "_" << imn << "_" << jmx;
 267         std::ostringstream wkt;
 268         wkt << std::setprecision(12) << bg::wkt(many_rings);
 269
 270         boost::timer t;
 271         test_one<multi_polygon_type, polygon_type>(out.str(), wkt.str(), join_round, end_flat, expected_area_exterior, 0.3);
 272         std::cout << "Exterior " << count << " " << t.elapsed() << std::endl;
 273     }
 274
 275     return;
 276     {
 277         std::ostringstream out;
 278         out << "many_interiors_" << count;
 279         std::ostringstream wkt;
 280         wkt << std::setprecision(12) << bg::wkt(many_interiors);
 281
 282         boost::timer t;
 283         test_one<multi_polygon_type, polygon_type>(out.str(), wkt.str(), join_round, end_flat, expected_area_interior, 0.3);
 284         std::cout << "Interior " << count << " " << t.elapsed() << std::endl;
 285     }
 286 }
 287
 288 int test_main(int, char* [])
 289 {
 290 //    test_many_rings(10, 795.70334, 806.7609);
 291 //    test_many_rings(30, 7136.7098, 6174.4496);
 292       test_many_rings(30, 30, 70, 38764.2721, 31910.3280);
 293 //    for (int i = 30; i < 60; i++)
 294 //    {
 295 //        for (int j = 5; j <= 30; j++)
 296 //        {
 297 //            test_many_rings(i, j, 70, 38764.2721, 31910.3280);
 298 //        }
 299 //    }
 300     return 0;
 301 }