1 files changed, 301 insertions, 0 deletions
diff --git a/libs/geometry/test/robustness/overlay/buffer/many_ring_buffer.cpp b/libs/geometry/test/robustness/overlay/buffer/many_ring_buffer.cpp
new file mode 100644
index 000000000..82226a087
--- /dev/null
+++ b/libs/geometry/test/robustness/overlay/buffer/many_ring_buffer.cpp
@@ -0,0 +1,301 @@
+// Boost.Geometry (aka GGL, Generic Geometry Library)
+// Unit Test
+
+// Copyright (c) 2012-2014 Barend Gehrels, Amsterdam, the Netherlands.
+
+// Use, modification and distribution is subject to the Boost Software License,
+// Version 1.0. (See accompanying file LICENSE_1_0.txt or copy at
+// http://www.boost.org/LICENSE_1_0.txt)
+
+#define BOOST_GEOMETRY_BUFFER_TEST_SVG_USE_ALTERNATE_BOX_FOR_INPUT
+#define BOOST_GEOMETRY_BUFFER_TEST_SVG_ALTERNATE_BOX "BOX(179 4, 180 5)"
+
+
+#include <test_buffer.hpp>
+
+#include <boost/geometry/algorithms/difference.hpp>
+#include <boost/geometry/multi/geometries/multi_geometries.hpp>
+
+#include <boost/random/linear_congruential.hpp>
+#include <boost/random/uniform_int.hpp>
+#include <boost/random/uniform_real.hpp>
+#include <boost/random/variate_generator.hpp>
+
+
+
+const int point_count = 90; // points for a full circle
+
+// Function to let buffer-distance depend on alpha, e.g.:
+inline double corrected_distance(double distance, double alpha)
+{
+    return distance * 1.0 + 0.2 * sin(alpha * 6.0);
+}
+
+class buffer_point_strategy_sample
+{
+public :
+
+    template
+    <
+        typename Point,
+        typename OutputRange,
+        typename DistanceStrategy
+    >
+    void apply(Point const& point,
+            DistanceStrategy const& distance_strategy,
+            OutputRange& output_range) const
+    {
+        double const distance = distance_strategy.apply(point, point,
+                        bg::strategy::buffer::buffer_side_left);
+
+        double const angle_increment = 2.0 * M_PI / double(point_count);
+
+        double alpha = 0;
+        for (std::size_t i = 0; i <= point_count; i++, alpha -= angle_increment)
+        {
+            double const cd = corrected_distance(distance, alpha);
+
+            typename boost::range_value<OutputRange>::type output_point;
+            bg::set<0>(output_point, bg::get<0>(point) + cd * cos(alpha));
+            bg::set<1>(output_point, bg::get<1>(point) + cd * sin(alpha));
+            output_range.push_back(output_point);
+        }
+    }
+};
+
+class buffer_join_strategy_sample
+{
+private :
+    template
+    <
+        typename Point,
+        typename DistanceType,
+        typename RangeOut
+    >
+    inline void generate_points(Point const& vertex,
+                Point const& perp1, Point const& perp2,
+                DistanceType const& buffer_distance,
+                RangeOut& range_out) const
+    {
+        double dx1 = bg::get<0>(perp1) - bg::get<0>(vertex);
+        double dy1 = bg::get<1>(perp1) - bg::get<1>(vertex);
+        double dx2 = bg::get<0>(perp2) - bg::get<0>(vertex);
+        double dy2 = bg::get<1>(perp2) - bg::get<1>(vertex);
+
+        // Assuming the corner is convex, angle2 < angle1
+        double const angle1 = atan2(dy1, dx1);
+        double angle2 = atan2(dy2, dx2);
+
+        while (angle2 > angle1)
+        {
+            angle2 -= 2 * M_PI;
+        }
+
+        double const angle_increment = 2.0 * M_PI / double(point_count);
+        double alpha = angle1 - angle_increment;
+
+        for (int i = 0; alpha >= angle2 && i < point_count; i++, alpha -= angle_increment)
+        {
+            double cd = corrected_distance(buffer_distance, alpha);
+
+            Point p;
+            bg::set<0>(p, bg::get<0>(vertex) + cd * cos(alpha));
+            bg::set<1>(p, bg::get<1>(vertex) + cd * sin(alpha));
+            range_out.push_back(p);
+        }
+    }
+
+public :
+
+    template <typename Point, typename DistanceType, typename RangeOut>
+    inline bool apply(Point const& ip, Point const& vertex,
+                Point const& perp1, Point const& perp2,
+                DistanceType const& buffer_distance,
+                RangeOut& range_out) const
+    {
+        generate_points(vertex, perp1, perp2, buffer_distance, range_out);
+        return true;
+    }
+
+    template <typename NumericType>
+    static inline NumericType max_distance(NumericType const& distance)
+    {
+        return distance;
+    }
+
+};
+
+class buffer_side_sample
+{
+public :
+    template
+    <
+        typename Point,
+        typename OutputRange,
+        typename DistanceStrategy
+    >
+    static inline void apply(
+                Point const& input_p1, Point const& input_p2,
+                bg::strategy::buffer::buffer_side_selector side,
+                DistanceStrategy const& distance,
+                OutputRange& output_range)
+    {
+        // Generate a block along (left or right of) the segment
+
+        double const dx = bg::get<0>(input_p2) - bg::get<0>(input_p1);
+        double const dy = bg::get<1>(input_p2) - bg::get<1>(input_p1);
+
+        // For normalization [0,1] (=dot product d.d, sqrt)
+        double const length = bg::math::sqrt(dx * dx + dy * dy);
+
+        if (bg::math::equals(length, 0))
+        {
+            return;
+        }
+
+        // Generate the normalized perpendicular p, to the left (ccw)
+        double const px = -dy / length;
+        double const py = dx / length;
+
+        // Both vectors perpendicular to input p1 and input p2 have same angle
+        double const alpha = atan2(py, px);
+
+        double const d = distance.apply(input_p1, input_p2, side);
+
+        double const cd = corrected_distance(d, alpha);
+
+        output_range.resize(2);
+
+        bg::set<0>(output_range.front(), bg::get<0>(input_p1) + px * cd);
+        bg::set<1>(output_range.front(), bg::get<1>(input_p1) + py * cd);
+        bg::set<0>(output_range.back(), bg::get<0>(input_p2) + px * cd);
+        bg::set<1>(output_range.back(), bg::get<1>(input_p2) + py * cd);
+    }
+};
+
+#ifdef TEST_WITH_SVG
+template <typename Geometry1, typename Geometry2>
+void create_svg(std::string const& filename, Geometry1 const& original, Geometry2 const& buffer, std::string const& color)
+{
+    typedef typename bg::point_type<Geometry1>::type point_type;
+    std::ofstream svg(filename.c_str());
+
+    bg::svg_mapper<point_type> mapper(svg, 800, 800);
+    mapper.add(buffer);
+
+    mapper.map(original, "fill-opacity:0.3;fill:rgb(255,0,0);stroke:rgb(0,0,0);stroke-width:1");
+
+    std::string style = "fill-opacity:0.3;fill:";
+    style += color;
+    style += ";stroke:rgb(0,0,0);stroke-width:1";
+    mapper.map(buffer, style);
+}
+#endif
+
+
+void test_many_rings(int imn, int jmx, int count,
+    double expected_area_exterior,
+    double expected_area_interior)
+{
+    typedef bg::model::point<double, 2, bg::cs::cartesian> point;
+    typedef bg::model::polygon<point> polygon_type;
+    typedef bg::model::multi_polygon<polygon_type> multi_polygon_type;
+
+    // Predefined strategies
+    bg::strategy::buffer::distance_symmetric<double> distance_strategy(1.3);
+    bg::strategy::buffer::end_flat end_strategy; // not effectively used
+
+    // Own strategies
+    buffer_join_strategy_sample join_strategy;
+    buffer_point_strategy_sample point_strategy;
+    buffer_side_sample side_strategy;
+
+    // Declare output
+
+    bg::model::multi_point<point> mp;
+
+    // Use a bit of random disturbance in the otherwise too regular grid
+    typedef boost::minstd_rand base_generator_type;
+    base_generator_type generator(12345);
+    boost::uniform_real<> random_range(0.0, 0.5);
+    boost::variate_generator
+    <
+        base_generator_type&,
+        boost::uniform_real<>
+    > random(generator, random_range);
+
+    for (int i = 0; i < count; i++)
+    {
+        for (int j = 0; j < count; j++)
+        {
+            double x = i * 3.0 + random();
+            double y = j * 3.0 + random();
+            //if (i > 30 && j < 30)
+            if (i > imn && j < jmx)
+            {
+                point p(x, y);
+                mp.push_back(p);
+            }
+        }
+    }
+
+    multi_polygon_type many_rings;
+    // Create the buffer of a multi-point
+    bg::buffer(mp, many_rings,
+                distance_strategy, side_strategy,
+                join_strategy, end_strategy, point_strategy);
+
+    bg::model::box<point> envelope;
+    bg::envelope(many_rings, envelope);
+    bg::buffer(envelope, envelope, 1.0);
+
+    multi_polygon_type many_interiors;
+    bg::difference(envelope, many_rings, many_interiors);
+
+#ifdef TEST_WITH_SVG
+    create_svg("/tmp/many_interiors.svg", mp, many_interiors, "rgb(51,51,153)");
+    create_svg("/tmp/buffer.svg", mp, many_rings, "rgb(51,51,153)");
+#endif
+
+    bg::strategy::buffer::join_round join_round(100);
+    bg::strategy::buffer::end_flat end_flat;
+
+    {
+        std::ostringstream out;
+        out << "many_rings_" << count;
+        out << "_" << imn << "_" << jmx;
+        std::ostringstream wkt;
+        wkt << std::setprecision(12) << bg::wkt(many_rings);
+
+        boost::timer t;
+        test_one<multi_polygon_type, polygon_type>(out.str(), wkt.str(), join_round, end_flat, expected_area_exterior, 0.3);
+        std::cout << "Exterior " << count << " " << t.elapsed() << std::endl;
+    }
+
+    return;
+    {
+        std::ostringstream out;
+        out << "many_interiors_" << count;
+        std::ostringstream wkt;
+        wkt << std::setprecision(12) << bg::wkt(many_interiors);
+
+        boost::timer t;
+        test_one<multi_polygon_type, polygon_type>(out.str(), wkt.str(), join_round, end_flat, expected_area_interior, 0.3);
+        std::cout << "Interior " << count << " " << t.elapsed() << std::endl;
+    }
+}
+
+int test_main(int, char* [])
+{
+//    test_many_rings(10, 795.70334, 806.7609);
+//    test_many_rings(30, 7136.7098, 6174.4496);
+      test_many_rings(30, 30, 70, 38764.2721, 31910.3280);
+//    for (int i = 30; i < 60; i++)
+//    {
+//        for (int j = 5; j <= 30; j++)
+//        {
+//            test_many_rings(i, j, 70, 38764.2721, 31910.3280);
+//        }
+//    }
+    return 0;
+}