ceph/src/boost/libs/geometry/include/boost/geometry/algorithms/detail/equals/collect_vectors.hpp

   1 // Boost.Geometry (aka GGL, Generic Geometry Library)
   2
   3 // Copyright (c) 2007-2014 Barend Gehrels, Amsterdam, the Netherlands.
   4 // Copyright (c) 2008-2014 Bruno Lalande, Paris, France.
   5 // Copyright (c) 2009-2014 Mateusz Loskot, London, UK.
   6 // Copyright (c) 2014 Adam Wulkiewicz, Lodz, Poland.
   7
   8 // Parts of Boost.Geometry are redesigned from Geodan's Geographic Library
   9 // (geolib/GGL), copyright (c) 1995-2010 Geodan, Amsterdam, the Netherlands.
  10
  11 // Use, modification and distribution is subject to the Boost Software License,
  12 // Version 1.0. (See accompanying file LICENSE_1_0.txt or copy at
  13 // http://www.boost.org/LICENSE_1_0.txt)
  14
  15 #ifndef BOOST_GEOMETRY_ALGORITHMS_DETAIL_EQUALS_COLLECT_VECTORS_HPP
  16 #define BOOST_GEOMETRY_ALGORITHMS_DETAIL_EQUALS_COLLECT_VECTORS_HPP
  17
  18
  19 #include <boost/numeric/conversion/cast.hpp>
  20
  21 #include <boost/geometry/algorithms/detail/interior_iterator.hpp>
  22 #include <boost/geometry/algorithms/detail/normalize.hpp>
  23
  24 #include <boost/geometry/core/cs.hpp>
  25 #include <boost/geometry/core/interior_rings.hpp>
  26 #include <boost/geometry/core/tags.hpp>
  27
  28 #include <boost/geometry/formulas/spherical.hpp>
  29
  30 #include <boost/geometry/geometries/concepts/check.hpp>
  31
  32 #include <boost/geometry/util/math.hpp>
  33 #include <boost/geometry/util/range.hpp>
  34
  35
  36 namespace boost { namespace geometry
  37 {
  38
  39 // Since these vectors (though ray would be a better name) are used in the
  40 // implementation of equals() for Areal geometries the internal representation
  41 // should be consistent with the default side strategy for CS because currently
  42 // it's used in other relops.
  43
  44 template <
  45     typename T,
  46     typename Geometry,
  47     typename CSTag = typename cs_tag<Geometry>::type
  48 >
  49 struct collected_vector
  50 {
  51     typedef T type;
  52
  53     inline collected_vector()
  54     {}
  55
  56     inline collected_vector(T const& px, T const& py,
  57                             T const& pdx, T const& pdy)
  58         : x(px)
  59         , y(py)
  60         , dx(pdx)
  61         , dy(pdy)
  62         //, dx_0(dx)
  63         //, dy_0(dy)
  64     {}
  65
  66     template <typename Point>
  67     inline collected_vector(Point const& p1, Point const& p2)
  68         : x(get<0>(p1))
  69         , y(get<1>(p1))
  70         , dx(get<0>(p2) - x)
  71         , dy(get<1>(p2) - y)
  72         //, dx_0(dx)
  73         //, dy_0(dy)
  74     {}
  75
  76     bool normalize()
  77     {
  78         T magnitude = math::sqrt(
  79             boost::numeric_cast<T>(dx * dx + dy * dy));
  80
  81         // NOTE: shouldn't here math::equals() be called?
  82         if (magnitude > 0)
  83         {
  84             dx /= magnitude;
  85             dy /= magnitude;
  86             return true;
  87         }
  88
  89         return false;
  90     }
  91
  92     // For sorting
  93     inline bool operator<(collected_vector const& other) const
  94     {
  95         if (math::equals(x, other.x))
  96         {
  97             if (math::equals(y, other.y))
  98             {
  99                 if (math::equals(dx, other.dx))
 100                 {
 101                     return dy < other.dy;
 102                 }
 103                 return dx < other.dx;
 104             }
 105             return y < other.y;
 106         }
 107         return x < other.x;
 108     }
 109
 110     inline bool next_is_collinear(collected_vector const& other) const
 111     {
 112         return same_direction(other);
 113     }
 114
 115     // For std::equals
 116     inline bool operator==(collected_vector const& other) const
 117     {
 118         return math::equals(x, other.x)
 119             && math::equals(y, other.y)
 120             && same_direction(other);
 121     }
 122
 123 private:
 124     inline bool same_direction(collected_vector const& other) const
 125     {
 126         // For high precision arithmetic, we have to be
 127         // more relaxed then using ==
 128         // Because 2/sqrt( (0,0)<->(2,2) ) == 1/sqrt( (0,0)<->(1,1) )
 129         // is not always true (at least, it is not for ttmath)
 130         return math::equals_with_epsilon(dx, other.dx)
 131             && math::equals_with_epsilon(dy, other.dy);
 132     }
 133
 134     T x, y;
 135     T dx, dy;
 136     //T dx_0, dy_0;
 137 };
 138
 139 template <typename T, typename Geometry>
 140 struct collected_vector<T, Geometry, spherical_equatorial_tag>
 141 {
 142     typedef T type;
 143
 144     typedef typename coordinate_system<Geometry>::type cs_type;
 145     typedef model::point<T, 2, cs_type> point_type;
 146     typedef model::point<T, 3, cs::cartesian> vector_type;
 147
 148     collected_vector()
 149     {}
 150
 151     template <typename Point>
 152     collected_vector(Point const& p1, Point const& p2)
 153         : origin(get<0>(p1), get<1>(p1))
 154     {
 155         origin = detail::return_normalized<point_type>(origin);
 156
 157         using namespace geometry::formula;
 158         prev = sph_to_cart3d<vector_type>(p1);
 159         next = sph_to_cart3d<vector_type>(p2);
 160         direction = cross_product(prev, next);
 161     }
 162
 163     bool normalize()
 164     {
 165         T magnitude_sqr = dot_product(direction, direction);
 166
 167         // NOTE: shouldn't here math::equals() be called?
 168         if (magnitude_sqr > 0)
 169         {
 170             divide_value(direction, math::sqrt(magnitude_sqr));
 171             return true;
 172         }
 173
 174         return false;
 175     }
 176
 177     bool operator<(collected_vector const& other) const
 178     {
 179         if (math::equals(get<0>(origin), get<0>(other.origin)))
 180         {
 181             if (math::equals(get<1>(origin), get<1>(other.origin)))
 182             {
 183                 if (math::equals(get<0>(direction), get<0>(other.direction)))
 184                 {
 185                     if (math::equals(get<1>(direction), get<1>(other.direction)))
 186                     {
 187                         return get<2>(direction) < get<2>(other.direction);
 188                     }
 189                     return get<1>(direction) < get<1>(other.direction);
 190                 }
 191                 return get<0>(direction) < get<0>(other.direction);
 192             }
 193             return get<1>(origin) < get<1>(other.origin);
 194         }
 195         return get<0>(origin) < get<0>(other.origin);
 196     }
 197
 198     // For consistency with side and intersection strategies used by relops
 199     // IMPORTANT: this method should be called for previous vector
 200     // and next vector should be passed as parameter
 201     bool next_is_collinear(collected_vector const& other) const
 202     {
 203         return formula::sph_side_value(direction, other.next) == 0;
 204     }
 205
 206     // For std::equals
 207     bool operator==(collected_vector const& other) const
 208     {
 209         return math::equals(get<0>(origin), get<0>(other.origin))
 210             && math::equals(get<1>(origin), get<1>(other.origin))
 211             && is_collinear(other);
 212     }
 213
 214 private:
 215     // For consistency with side and intersection strategies used by relops
 216     bool is_collinear(collected_vector const& other) const
 217     {
 218         return formula::sph_side_value(direction, other.prev) == 0
 219             && formula::sph_side_value(direction, other.next) == 0;
 220     }
 221
 222     /*bool same_direction(collected_vector const& other) const
 223     {
 224         return math::equals_with_epsilon(get<0>(direction), get<0>(other.direction))
 225             && math::equals_with_epsilon(get<1>(direction), get<1>(other.direction))
 226             && math::equals_with_epsilon(get<2>(direction), get<2>(other.direction));
 227     }*/
 228
 229     point_type origin; // used for sorting and equality check
 230     vector_type direction; // used for sorting, only in operator<
 231     vector_type prev; // used for collinearity check, only in operator==
 232     vector_type next; // used for collinearity check
 233 };
 234
 235 template <typename T, typename Geometry>
 236 struct collected_vector<T, Geometry, spherical_polar_tag>
 237     : public collected_vector<T, Geometry, spherical_equatorial_tag>
 238 {
 239     typedef collected_vector<T, Geometry, spherical_equatorial_tag> base_type;
 240
 241     collected_vector() {}
 242
 243     template <typename Point>
 244     collected_vector(Point const& p1, Point const& p2)
 245         : base_type(to_equatorial(p1), to_equatorial(p2))
 246     {}
 247
 248 private:
 249     template <typename Point>
 250     Point polar_to_equatorial(Point const& p)
 251     {
 252         typedef typename coordinate_type<Point>::type coord_type;
 253
 254         typedef math::detail::constants_on_spheroid
 255             <
 256                 coord_type,
 257                 typename coordinate_system<Point>::type::units
 258             > constants;
 259
 260         coord_type const pi_2 = constants::half_period() / 2;
 261
 262         Point res = p;
 263         set<1>(res, pi_2 - get<1>(p));
 264         return res;
 265     }
 266 };
 267
 268 // This is consistent with the currently used default geographic side
 269 // and intersection strategies. Spherical strategies are used by default.
 270 // When default strategies are changed in the future this specialization
 271 // should be changed too.
 272 template <typename T, typename Geometry>
 273 struct collected_vector<T, Geometry, geographic_tag>
 274     : public collected_vector<T, Geometry, spherical_equatorial_tag>
 275 {
 276     typedef collected_vector<T, Geometry, spherical_equatorial_tag> base_type;
 277
 278     collected_vector() {}
 279
 280     template <typename Point>
 281     collected_vector(Point const& p1, Point const& p2)
 282         : base_type(p1, p2)
 283     {}
 284 };
 285
 286
 287 #ifndef DOXYGEN_NO_DETAIL
 288 namespace detail { namespace collect_vectors
 289 {
 290
 291
 292 template <typename Range, typename Collection>
 293 struct range_collect_vectors
 294 {
 295     typedef typename boost::range_value<Collection>::type item_type;
 296     typedef typename item_type::type calculation_type;
 297
 298     static inline void apply(Collection& collection, Range const& range)
 299     {
 300         if (boost::size(range) < 2)
 301         {
 302             return;
 303         }
 304
 305         typedef typename boost::range_size<Collection>::type collection_size_t;
 306         collection_size_t c_old_size = boost::size(collection);
 307
 308         typedef typename boost::range_iterator<Range const>::type iterator;
 309
 310         bool is_first = true;
 311         iterator it = boost::begin(range);
 312
 313         for (iterator prev = it++;
 314             it != boost::end(range);
 315             prev = it++)
 316         {
 317             typename boost::range_value<Collection>::type v(*prev, *it);
 318
 319             // Normalize the vector -> this results in points+direction
 320             // and is comparible between geometries
 321             // Avoid non-duplicate points (AND division by zero)
 322             if (v.normalize())
 323             {
 324                 // Avoid non-direction changing points
 325                 if (is_first || ! collection.back().next_is_collinear(v))
 326                 {
 327                     collection.push_back(v);
 328                 }
 329                 is_first = false;
 330             }
 331         }
 332
 333         // If first one has same direction as last one, remove first one
 334         collection_size_t collected_count = boost::size(collection) - c_old_size;
 335         if ( collected_count > 1 )
 336         {
 337             typedef typename boost::range_iterator<Collection>::type c_iterator;
 338             c_iterator first = range::pos(collection, c_old_size);
 339
 340             if (collection.back().next_is_collinear(*first) )
 341             {
 342                 //collection.erase(first);
 343                 // O(1) instead of O(N)
 344                 *first = collection.back();
 345                 collection.pop_back();
 346             }
 347         }
 348     }
 349 };
 350
 351
 352 // Default version (cartesian)
 353 template <typename Box, typename Collection, typename CSTag = typename cs_tag<Box>::type>
 354 struct box_collect_vectors
 355 {
 356     // Calculate on coordinate type, but if it is integer,
 357     // then use double
 358     typedef typename boost::range_value<Collection>::type item_type;
 359     typedef typename item_type::type calculation_type;
 360
 361     static inline void apply(Collection& collection, Box const& box)
 362     {
 363         typename point_type<Box>::type lower_left, lower_right,
 364             upper_left, upper_right;
 365         geometry::detail::assign_box_corners(box, lower_left, lower_right,
 366             upper_left, upper_right);
 367
 368         typedef typename boost::range_value<Collection>::type item;
 369
 370         collection.push_back(item(get<0>(lower_left), get<1>(lower_left), 0, 1));
 371         collection.push_back(item(get<0>(upper_left), get<1>(upper_left), 1, 0));
 372         collection.push_back(item(get<0>(upper_right), get<1>(upper_right), 0, -1));
 373         collection.push_back(item(get<0>(lower_right), get<1>(lower_right), -1, 0));
 374     }
 375 };
 376
 377 // NOTE: This is not fully correct because Box in spherical and geographic
 378 // cordinate systems cannot be seen as Polygon
 379 template <typename Box, typename Collection>
 380 struct box_collect_vectors<Box, Collection, spherical_equatorial_tag>
 381 {
 382     static inline void apply(Collection& collection, Box const& box)
 383     {
 384         typename point_type<Box>::type lower_left, lower_right,
 385                 upper_left, upper_right;
 386         geometry::detail::assign_box_corners(box, lower_left, lower_right,
 387                 upper_left, upper_right);
 388
 389         typedef typename boost::range_value<Collection>::type item;
 390
 391         collection.push_back(item(lower_left, upper_left));
 392         collection.push_back(item(upper_left, upper_right));
 393         collection.push_back(item(upper_right, lower_right));
 394         collection.push_back(item(lower_right, lower_left));
 395     }
 396 };
 397
 398 template <typename Box, typename Collection>
 399 struct box_collect_vectors<Box, Collection, spherical_polar_tag>
 400     : box_collect_vectors<Box, Collection, spherical_equatorial_tag>
 401 {};
 402
 403 template <typename Box, typename Collection>
 404 struct box_collect_vectors<Box, Collection, geographic_tag>
 405     : box_collect_vectors<Box, Collection, spherical_equatorial_tag>
 406 {};
 407
 408
 409 template <typename Polygon, typename Collection>
 410 struct polygon_collect_vectors
 411 {
 412     static inline void apply(Collection& collection, Polygon const& polygon)
 413     {
 414         typedef typename geometry::ring_type<Polygon>::type ring_type;
 415
 416         typedef range_collect_vectors<ring_type, Collection> per_range;
 417         per_range::apply(collection, exterior_ring(polygon));
 418
 419         typename interior_return_type<Polygon const>::type
 420             rings = interior_rings(polygon);
 421         for (typename detail::interior_iterator<Polygon const>::type
 422                 it = boost::begin(rings); it != boost::end(rings); ++it)
 423         {
 424             per_range::apply(collection, *it);
 425         }
 426     }
 427 };
 428
 429
 430 template <typename MultiGeometry, typename Collection, typename SinglePolicy>
 431 struct multi_collect_vectors
 432 {
 433     static inline void apply(Collection& collection, MultiGeometry const& multi)
 434     {
 435         for (typename boost::range_iterator<MultiGeometry const>::type
 436                 it = boost::begin(multi);
 437             it != boost::end(multi);
 438             ++it)
 439         {
 440             SinglePolicy::apply(collection, *it);
 441         }
 442     }
 443 };
 444
 445
 446 }} // namespace detail::collect_vectors
 447 #endif // DOXYGEN_NO_DETAIL
 448
 449
 450
 451 #ifndef DOXYGEN_NO_DISPATCH
 452 namespace dispatch
 453 {
 454
 455
 456 template
 457 <
 458     typename Tag,
 459     typename Collection,
 460     typename Geometry
 461 >
 462 struct collect_vectors
 463 {
 464     static inline void apply(Collection&, Geometry const&)
 465     {}
 466 };
 467
 468
 469 template <typename Collection, typename Box>
 470 struct collect_vectors<box_tag, Collection, Box>
 471     : detail::collect_vectors::box_collect_vectors<Box, Collection>
 472 {};
 473
 474
 475
 476 template <typename Collection, typename Ring>
 477 struct collect_vectors<ring_tag, Collection, Ring>
 478     : detail::collect_vectors::range_collect_vectors<Ring, Collection>
 479 {};
 480
 481
 482 template <typename Collection, typename LineString>
 483 struct collect_vectors<linestring_tag, Collection, LineString>
 484     : detail::collect_vectors::range_collect_vectors<LineString, Collection>
 485 {};
 486
 487
 488 template <typename Collection, typename Polygon>
 489 struct collect_vectors<polygon_tag, Collection, Polygon>
 490     : detail::collect_vectors::polygon_collect_vectors<Polygon, Collection>
 491 {};
 492
 493
 494 template <typename Collection, typename MultiPolygon>
 495 struct collect_vectors<multi_polygon_tag, Collection, MultiPolygon>
 496     : detail::collect_vectors::multi_collect_vectors
 497         <
 498             MultiPolygon,
 499             Collection,
 500             detail::collect_vectors::polygon_collect_vectors
 501             <
 502                 typename boost::range_value<MultiPolygon>::type,
 503                 Collection
 504             >
 505         >
 506 {};
 507
 508
 509
 510 } // namespace dispatch
 511 #endif
 512
 513
 514 /*!
 515     \ingroup collect_vectors
 516     \tparam Collection Collection type, should be e.g. std::vector<>
 517     \tparam Geometry geometry type
 518     \param collection the collection of vectors
 519     \param geometry the geometry to make collect_vectors
 520 */
 521 template <typename Collection, typename Geometry>
 522 inline void collect_vectors(Collection& collection, Geometry const& geometry)
 523 {
 524     concepts::check<Geometry const>();
 525
 526     dispatch::collect_vectors
 527         <
 528             typename tag<Geometry>::type,
 529             Collection,
 530             Geometry
 531         >::apply(collection, geometry);
 532 }
 533
 534
 535 }} // namespace boost::geometry
 536
 537
 538 #endif // BOOST_GEOMETRY_ALGORITHMS_DETAIL_EQUALS_COLLECT_VECTORS_HPP