ceph/src/boost/libs/geometry/test/algorithms/densify.cpp

   1 // Boost.Geometry
   2 // Unit Test
   3
   4 // Copyright (c) 2017-2021, Oracle and/or its affiliates.
   5
   6 // Contributed and/or modified by Adam Wulkiewicz, on behalf of Oracle
   7
   8 // Licensed under the Boost Software License version 1.0.
   9 // http://www.boost.org/users/license.html
  10
  11
  12 #include <boost/variant/variant.hpp>
  13
  14 #include <geometry_test_common.hpp>
  15
  16 #include <boost/geometry/geometries/geometries.hpp>
  17
  18 #include <boost/geometry/algorithms/densify.hpp>
  19 #include <boost/geometry/algorithms/length.hpp>
  20 #include <boost/geometry/algorithms/num_points.hpp>
  21 #include <boost/geometry/algorithms/perimeter.hpp>
  22
  23 #include <boost/geometry/iterators/segment_iterator.hpp>
  24
  25 #include <boost/geometry/strategies/cartesian/densify.hpp>
  26 #include <boost/geometry/strategies/cartesian/distance_pythagoras.hpp>
  27 #include <boost/geometry/strategies/geographic/densify.hpp>
  28 #include <boost/geometry/strategies/geographic/distance.hpp>
  29 #include <boost/geometry/strategies/spherical/densify.hpp>
  30 #include <boost/geometry/strategies/spherical/distance_haversine.hpp>
  31
  32 #include <boost/geometry/io/wkt/wkt.hpp>
  33
  34
  35 struct check_lengths
  36 {
  37     template <typename G, typename S>
  38     void operator()(G const& g, G const& o, S const& s) const
  39     {
  40         double d1 = bg::length(g, s);
  41         double d2 = bg::length(o, s);
  42
  43         BOOST_CHECK_CLOSE(d1, d2, 0.0001);
  44     }
  45 };
  46
  47 struct check_perimeters
  48 {
  49     template <typename G, typename S>
  50     void operator()(G const& g, G const& o, S const& s) const
  51     {
  52         double d1 = bg::perimeter(g, s);
  53         double d2 = bg::perimeter(o, s);
  54
  55         BOOST_CHECK_CLOSE(d1, d2, 0.0001);
  56     }
  57 };
  58
  59 template <typename G, typename DistS>
  60 double inline shortest_length(G const& g, DistS const& dist_s)
  61 {
  62     double min_len = (std::numeric_limits<double>::max)();
  63     for (bg::segment_iterator<G const> it = bg::segments_begin(g);
  64             it != bg::segments_end(g); ++it)
  65     {
  66         double len = bg::length(*it, dist_s);
  67         min_len = (std::min)(min_len, len);
  68     }
  69     return min_len;
  70 }
  71
  72 template <typename G, typename DistS>
  73 double inline greatest_length(G const& o, DistS const& dist_s)
  74 {
  75     double max_len = 0.0;
  76     for (bg::segment_iterator<G const> it = bg::segments_begin(o);
  77             it != bg::segments_end(o); ++it)
  78     {
  79         double len = bg::length(*it, dist_s);
  80         max_len = (std::max)(max_len, len);
  81     }
  82     return max_len;
  83 }
  84
  85 template <typename G, typename CSTag = typename bg::cs_tag<G>::type>
  86 struct cs_data
  87 {};
  88
  89 template <typename G>
  90 struct cs_data<G, bg::cartesian_tag>
  91 {
  92     bg::strategy::densify::cartesian<> compl_s;
  93     bg::strategy::distance::pythagoras<> dist_s;
  94 };
  95
  96 template <typename G>
  97 struct cs_data<G, bg::spherical_equatorial_tag>
  98 {
  99     cs_data()
 100         : model(6378137.0)
 101         , compl_s(model)
 102         , dist_s(6378137.0)
 103     {}
 104
 105     bg::srs::sphere<double> model;
 106     bg::strategy::densify::spherical<> compl_s;
 107     bg::strategy::distance::haversine<double> dist_s;
 108 };
 109
 110 template <typename G>
 111 struct cs_data<G, bg::geographic_tag>
 112 {
 113     cs_data()
 114         : model(6378137.0, 6356752.3142451793)
 115         , compl_s(model)
 116         , dist_s(model)
 117     {}
 118
 119     bg::srs::spheroid<double> model;
 120     bg::strategy::densify::geographic<> compl_s;
 121     bg::strategy::distance::geographic<> dist_s;
 122 };
 123
 124 template <typename G, typename DistS, typename Check>
 125 inline void check_result(G const& g, G const& o, double max_distance,
 126                          DistS const& dist_s, Check const& check)
 127 {
 128     // geometry was indeed densified
 129     std::size_t g_count = bg::num_points(g);
 130     std::size_t o_count = bg::num_points(o);
 131     BOOST_CHECK(g_count < o_count);
 132
 133     // all segments have lengths smaller or equal to max_distance
 134     double gr_len = greatest_length(o, dist_s);
 135     // NOTE: Currently geographic strategies can generate segments that have
 136     //       lengths slightly greater than max_distance. In order to change
 137     //       this the generation of new points should e.g. be recursive with
 138     //       stop condition comparing the current distance calculated by
 139     //       inverse strategy.
 140     // NOTE: Closeness value tweaked for Andoyer
 141     bool is_close = (gr_len - max_distance) / (std::max)(gr_len, max_distance) < 0.0001;
 142     BOOST_CHECK(gr_len <= max_distance || is_close);
 143
 144     // the overall length or perimeter didn't change
 145     check(g, o, dist_s);
 146 }
 147
 148 template <typename G, typename Check>
 149 inline void test_geometry(std::string const& wkt, Check const& check)
 150 {
 151     cs_data<G> d;
 152
 153     G g;
 154     bg::read_wkt(wkt, g);
 155
 156     {
 157         bg::default_strategy def_s;
 158         double max_distance = shortest_length(g, def_s) / 3.0;
 159
 160         G o;
 161         bg::densify(g, o, max_distance);
 162
 163         check_result(g, o, max_distance, def_s, check);
 164
 165         using variant_t = boost::variant<G, typename bg::point_type<G>::type>;
 166         variant_t v = g, vo;
 167         bg::densify(v, vo, max_distance);
 168
 169         check(v, vo, def_s);
 170
 171         bg::model::geometry_collection<variant_t> gc{v}, gco;
 172         bg::densify(gc, gco, max_distance);
 173
 174         check(gc, gco, def_s);
 175     }
 176
 177     {
 178         double max_distance = shortest_length(g, d.dist_s) / 3.0;
 179
 180         G o;
 181         bg::densify(g, o, max_distance, d.compl_s);
 182
 183         check_result(g, o, max_distance, d.dist_s, check);
 184     }
 185 }
 186
 187 template <typename G>
 188 inline void test_linear(std::string const& wkt)
 189 {
 190     test_geometry<G>(wkt, check_lengths());
 191 }
 192
 193 template <typename G>
 194 inline void test_areal(std::string const& wkt)
 195 {
 196     test_geometry<G>(wkt, check_perimeters());
 197 }
 198
 199 template <typename P>
 200 void test_all()
 201 {
 202     typedef bg::model::linestring<P> ls_t;
 203     typedef bg::model::multi_linestring<ls_t> mls_t;
 204
 205     typedef bg::model::ring<P> ring_t;
 206     typedef bg::model::polygon<P> poly_t;
 207     typedef bg::model::multi_polygon<poly_t> mpoly_t;
 208
 209     typedef bg::model::ring<P, true, false> oring_t;
 210     typedef bg::model::polygon<P, true, false> opoly_t;
 211     typedef bg::model::multi_polygon<opoly_t> ompoly_t;
 212
 213     test_linear<ls_t>("LINESTRING(4 -4, 4 -1)");
 214     test_linear<ls_t>("LINESTRING(4 4, 4 1)");
 215     test_linear<ls_t>("LINESTRING(0 0, 180 0)");
 216     test_linear<ls_t>("LINESTRING(1 1, -179 -1)");
 217
 218     test_linear<ls_t>("LINESTRING(1 1, 2 2, 4 2)");
 219     test_linear<mls_t>("MULTILINESTRING((1 1, 2 2),(2 2, 4 2))");
 220
 221     test_areal<ring_t>("POLYGON((1 1, 1 2, 2 2, 1 1))");
 222     test_areal<poly_t>("POLYGON((1 1, 1 4, 4 4, 4 1, 1 1),(1 1, 2 2, 2 3, 1 1))");
 223     test_areal<mpoly_t>("MULTIPOLYGON(((1 1, 1 4, 4 4, 4 1, 1 1),(1 1, 2 2, 2 3, 1 1)),((4 4, 5 5, 5 4, 4 4)))");
 224
 225     test_areal<oring_t>("POLYGON((1 1, 1 2, 2 2))");
 226     test_areal<opoly_t>("POLYGON((1 1, 1 4, 4 4, 4 1),(1 1, 2 2, 2 3))");
 227     test_areal<ompoly_t>("MULTIPOLYGON(((1 1, 1 4, 4 4, 4 1),(1 1, 2 2, 2 3)),((4 4, 5 5, 5 4)))");
 228
 229     test_areal<ring_t>("POLYGON((0 0,0 40,40 40,40 0,0 0))");
 230     test_areal<oring_t>("POLYGON((0 0,0 40,40 40,40 0))");
 231 }
 232
 233 int test_main(int, char* [])
 234 {
 235     test_all< bg::model::point<double, 2, bg::cs::cartesian> >();
 236     test_all< bg::model::point<double, 2, bg::cs::spherical_equatorial<bg::degree> > >();
 237     test_all< bg::model::point<double, 2, bg::cs::geographic<bg::degree> > >();
 238
 239     return 0;
 240 }