ceph/src/boost/boost/geometry/algorithms/detail/extreme_points.hpp

   1 // Boost.Geometry (aka GGL, Generic Geometry Library)
   2
   3 // Copyright (c) 2007-2013 Barend Gehrels, Amsterdam, the Netherlands.
   4 // Copyright (c) 2008-2013 Bruno Lalande, Paris, France.
   5 // Copyright (c) 2009-2013 Mateusz Loskot, London, UK.
   6 // Copyright (c) 2013-2014 Adam Wulkiewicz, Lodz, Poland.
   7
   8 // This file was modified by Oracle on 2017.
   9 // Modifications copyright (c) 2017 Oracle and/or its affiliates.
  10
  11 // Contributed and/or modified by Adam Wulkiewicz, on behalf of Oracle
  12
  13 // Use, modification and distribution is subject to the Boost Software License,
  14 // Version 1.0. (See accompanying file LICENSE_1_0.txt or copy at
  15 // http://www.boost.org/LICENSE_1_0.txt)
  16
  17 #ifndef BOOST_GEOMETRY_ALGORITHMS_DETAIL_EXTREME_POINTS_HPP
  18 #define BOOST_GEOMETRY_ALGORITHMS_DETAIL_EXTREME_POINTS_HPP
  19
  20
  21 #include <cstddef>
  22
  23 #include <boost/range.hpp>
  24
  25 #include <boost/geometry/algorithms/detail/interior_iterator.hpp>
  26
  27 #include <boost/geometry/core/cs.hpp>
  28 #include <boost/geometry/core/point_type.hpp>
  29 #include <boost/geometry/core/ring_type.hpp>
  30 #include <boost/geometry/core/tags.hpp>
  31
  32 #include <boost/geometry/geometries/concepts/check.hpp>
  33 #include <boost/geometry/iterators/ever_circling_iterator.hpp>
  34
  35 #include <boost/geometry/algorithms/detail/assign_box_corners.hpp>
  36
  37 #include <boost/geometry/strategies/side.hpp>
  38
  39 #include <boost/geometry/util/math.hpp>
  40
  41
  42 namespace boost { namespace geometry
  43 {
  44
  45
  46 #ifndef DOXYGEN_NO_DETAIL
  47 namespace detail { namespace extreme_points
  48 {
  49
  50 template <std::size_t Dimension>
  51 struct compare
  52 {
  53     template <typename Point>
  54     inline bool operator()(Point const& lhs, Point const& rhs)
  55     {
  56         return geometry::get<Dimension>(lhs) < geometry::get<Dimension>(rhs);
  57     }
  58 };
  59
  60
  61 template <std::size_t Dimension, typename PointType, typename CoordinateType>
  62 inline void move_along_vector(PointType& point, PointType const& extreme, CoordinateType const& base_value)
  63 {
  64     // Moves a point along the vector (point, extreme) in the direction of the extreme  point
  65     // This adapts the possibly uneven legs of the triangle (or trapezium-like shape)
  66     //      _____extreme            _____
  67     //     /     \                 /     \                                    .
  68     //    /base   \         =>    /       \ point                             .
  69     //             \ point                                                    .
  70     //
  71     // For so-called intruders, it can be used to adapt both legs to the level of "base"
  72     // For the base, it can be used to adapt both legs to the level of the max-value of the intruders
  73     // If there are 2 or more extreme values, use the one close to 'point' to have a correct vector
  74
  75     CoordinateType const value = geometry::get<Dimension>(point);
  76     //if (geometry::math::equals(value, base_value))
  77     if (value >= base_value)
  78     {
  79         return;
  80     }
  81
  82     PointType vector = point;
  83     subtract_point(vector, extreme);
  84
  85     CoordinateType const diff = geometry::get<Dimension>(vector);
  86
  87     // diff should never be zero
  88     // because of the way our triangle/trapezium is build.
  89     // We just return if it would be the case.
  90     if (geometry::math::equals(diff, 0))
  91     {
  92         return;
  93     }
  94
  95     CoordinateType const base_diff = base_value - geometry::get<Dimension>(extreme);
  96
  97     multiply_value(vector, base_diff);
  98     divide_value(vector, diff);
  99
 100     // The real move:
 101     point = extreme;
 102     add_point(point, vector);
 103 }
 104
 105
 106 template <std::size_t Dimension, typename Range, typename CoordinateType>
 107 inline void move_along_vector(Range& range, CoordinateType const& base_value)
 108 {
 109     if (range.size() >= 3)
 110     {
 111         move_along_vector<Dimension>(range.front(), *(range.begin() + 1), base_value);
 112         move_along_vector<Dimension>(range.back(), *(range.rbegin() + 1), base_value);
 113     }
 114 }
 115
 116
 117 template<typename Ring, std::size_t Dimension>
 118 struct extreme_points_on_ring
 119 {
 120
 121     typedef typename geometry::coordinate_type<Ring>::type coordinate_type;
 122     typedef typename boost::range_iterator<Ring const>::type range_iterator;
 123     typedef typename geometry::point_type<Ring>::type point_type;
 124
 125     template <typename CirclingIterator, typename Points>
 126     static inline bool extend(CirclingIterator& it,
 127             std::size_t n,
 128             coordinate_type max_coordinate_value,
 129             Points& points, int direction)
 130     {
 131         std::size_t safe_index = 0;
 132         do
 133         {
 134             it += direction;
 135
 136             points.push_back(*it);
 137
 138             if (safe_index++ >= n)
 139             {
 140                 // E.g.: ring is completely horizontal or vertical (= invalid, but we don't want to have an infinite loop)
 141                 return false;
 142             }
 143         } while (geometry::math::equals(geometry::get<Dimension>(*it), max_coordinate_value));
 144
 145         return true;
 146     }
 147
 148     // Overload without adding to poinst
 149     template <typename CirclingIterator>
 150     static inline bool extend(CirclingIterator& it,
 151             std::size_t n,
 152             coordinate_type max_coordinate_value,
 153             int direction)
 154     {
 155         std::size_t safe_index = 0;
 156         do
 157         {
 158             it += direction;
 159
 160             if (safe_index++ >= n)
 161             {
 162                 // E.g.: ring is completely horizontal or vertical (= invalid, but we don't want to have an infinite loop)
 163                 return false;
 164             }
 165         } while (geometry::math::equals(geometry::get<Dimension>(*it), max_coordinate_value));
 166
 167         return true;
 168     }
 169
 170     template <typename CirclingIterator>
 171     static inline bool extent_both_sides(Ring const& ring,
 172             point_type extreme,
 173             CirclingIterator& left,
 174             CirclingIterator& right)
 175     {
 176         std::size_t const n = boost::size(ring);
 177         coordinate_type const max_coordinate_value = geometry::get<Dimension>(extreme);
 178
 179         if (! extend(left, n, max_coordinate_value, -1))
 180         {
 181             return false;
 182         }
 183         if (! extend(right, n, max_coordinate_value, +1))
 184         {
 185             return false;
 186         }
 187
 188         return true;
 189     }
 190
 191     template <typename Collection, typename CirclingIterator>
 192     static inline bool collect(Ring const& ring,
 193             point_type extreme,
 194             Collection& points,
 195             CirclingIterator& left,
 196             CirclingIterator& right)
 197     {
 198         std::size_t const n = boost::size(ring);
 199         coordinate_type const max_coordinate_value = geometry::get<Dimension>(extreme);
 200
 201         // Collects first left, which is reversed (if more than one point) then adds the top itself, then right
 202         if (! extend(left, n, max_coordinate_value, points, -1))
 203         {
 204             return false;
 205         }
 206         std::reverse(points.begin(), points.end());
 207         points.push_back(extreme);
 208         if (! extend(right, n, max_coordinate_value, points, +1))
 209         {
 210             return false;
 211         }
 212
 213         return true;
 214     }
 215
 216     template <typename Extremes, typename Intruders, typename CirclingIterator, typename SideStrategy>
 217     static inline void get_intruders(Ring const& ring, CirclingIterator left, CirclingIterator right,
 218             Extremes const& extremes,
 219             Intruders& intruders,
 220             SideStrategy const& strategy)
 221     {
 222         if (boost::size(extremes) < 3)
 223         {
 224             return;
 225         }
 226         coordinate_type const min_value = geometry::get<Dimension>(*std::min_element(boost::begin(extremes), boost::end(extremes), compare<Dimension>()));
 227
 228         // Also select left/right (if Dimension=1)
 229         coordinate_type const other_min = geometry::get<1 - Dimension>(*std::min_element(boost::begin(extremes), boost::end(extremes), compare<1 - Dimension>()));
 230         coordinate_type const other_max = geometry::get<1 - Dimension>(*std::max_element(boost::begin(extremes), boost::end(extremes), compare<1 - Dimension>()));
 231
 232         std::size_t defensive_check_index = 0; // in case we skip over left/right check, collect modifies right too
 233         std::size_t const n = boost::size(ring);
 234         while (left != right && defensive_check_index < n)
 235         {
 236             coordinate_type const coordinate = geometry::get<Dimension>(*right);
 237             coordinate_type const other_coordinate = geometry::get<1 - Dimension>(*right);
 238             if (coordinate > min_value && other_coordinate > other_min && other_coordinate < other_max)
 239             {
 240                 int const factor = geometry::point_order<Ring>::value == geometry::clockwise ? 1 : -1;
 241                 int const first_side = strategy.apply(*right, extremes.front(), *(extremes.begin() + 1)) * factor;
 242                 int const last_side = strategy.apply(*right, *(extremes.rbegin() + 1), extremes.back()) * factor;
 243
 244                 // If not lying left from any of the extemes side
 245                 if (first_side != 1 && last_side != 1)
 246                 {
 247                     //std::cout << "first " << first_side << " last " << last_side << std::endl;
 248
 249                     // we start at this intrusion until it is handled, and don't affect our initial left iterator
 250                     CirclingIterator left_intrusion_it = right;
 251                     typename boost::range_value<Intruders>::type intruder;
 252                     collect(ring, *right, intruder, left_intrusion_it, right);
 253
 254                     // Also moves these to base-level, makes sorting possible which can be done in case of self-tangencies
 255                     // (we might postpone this action, it is often not necessary. However it is not time-consuming)
 256                     move_along_vector<Dimension>(intruder, min_value);
 257                     intruders.push_back(intruder);
 258                     --right;
 259                 }
 260             }
 261             ++right;
 262             defensive_check_index++;
 263         }
 264     }
 265
 266     template <typename Extremes, typename Intruders, typename SideStrategy>
 267     static inline void get_intruders(Ring const& ring,
 268             Extremes const& extremes,
 269             Intruders& intruders,
 270             SideStrategy const& strategy)
 271     {
 272         std::size_t const n = boost::size(ring);
 273         if (n >= 3)
 274         {
 275             geometry::ever_circling_range_iterator<Ring const> left(ring);
 276             geometry::ever_circling_range_iterator<Ring const> right(ring);
 277             ++right;
 278
 279             get_intruders(ring, left, right, extremes, intruders, strategy);
 280         }
 281     }
 282
 283     template <typename Iterator, typename SideStrategy>
 284     static inline bool right_turn(Ring const& ring, Iterator it, SideStrategy const& strategy)
 285     {
 286         typename std::iterator_traits<Iterator>::difference_type const index
 287             = std::distance(boost::begin(ring), it);
 288         geometry::ever_circling_range_iterator<Ring const> left(ring);
 289         geometry::ever_circling_range_iterator<Ring const> right(ring);
 290         left += index;
 291         right += index;
 292
 293         if (! extent_both_sides(ring, *it, left, right))
 294         {
 295             return false;
 296         }
 297
 298         int const factor = geometry::point_order<Ring>::value == geometry::clockwise ? 1 : -1;
 299         int const first_side = strategy.apply(*(right - 1), *right, *left) * factor;
 300         int const last_side = strategy.apply(*left, *(left + 1), *right) * factor;
 301
 302 //std::cout << "Candidate at " << geometry::wkt(*it) << " first=" << first_side << " last=" << last_side << std::endl;
 303
 304         // Turn should not be left (actually, it should be right because extent removes horizontal/collinear cases)
 305         return first_side != 1 && last_side != 1;
 306     }
 307
 308
 309     // Gets the extreme segments (top point plus neighbouring points), plus intruders, if any, on the same ring
 310     template <typename Extremes, typename Intruders, typename SideStrategy>
 311     static inline bool apply(Ring const& ring,
 312                              Extremes& extremes,
 313                              Intruders& intruders,
 314                              SideStrategy const& strategy)
 315     {
 316         std::size_t const n = boost::size(ring);
 317         if (n < 3)
 318         {
 319             return false;
 320         }
 321
 322         // Get all maxima, usually one. In case of self-tangencies, or self-crossings,
 323         // the max might be is not valid. A valid max should make a right turn
 324         range_iterator max_it = boost::begin(ring);
 325         compare<Dimension> smaller;
 326         for (range_iterator it = max_it + 1; it != boost::end(ring); ++it)
 327         {
 328             if (smaller(*max_it, *it) && right_turn(ring, it, strategy))
 329             {
 330                 max_it = it;
 331             }
 332         }
 333
 334         if (max_it == boost::end(ring))
 335         {
 336             return false;
 337         }
 338
 339         typename std::iterator_traits<range_iterator>::difference_type const
 340             index = std::distance(boost::begin(ring), max_it);
 341 //std::cout << "Extreme point lies at " << index << " having " << geometry::wkt(*max_it) << std::endl;
 342
 343         geometry::ever_circling_range_iterator<Ring const> left(ring);
 344         geometry::ever_circling_range_iterator<Ring const> right(ring);
 345         left += index;
 346         right += index;
 347
 348         // Collect all points (often 3) in a temporary vector
 349         std::vector<point_type> points;
 350         points.reserve(3);
 351         if (! collect(ring, *max_it, points, left, right))
 352         {
 353             return false;
 354         }
 355
 356 //std::cout << "Built vector of " << points.size() << std::endl;
 357
 358         coordinate_type const front_value = geometry::get<Dimension>(points.front());
 359         coordinate_type const back_value = geometry::get<Dimension>(points.back());
 360         coordinate_type const base_value = (std::max)(front_value, back_value);
 361         if (front_value < back_value)
 362         {
 363             move_along_vector<Dimension>(points.front(), *(points.begin() + 1), base_value);
 364         }
 365         else
 366         {
 367             move_along_vector<Dimension>(points.back(), *(points.rbegin() + 1), base_value);
 368         }
 369
 370         std::copy(points.begin(), points.end(), std::back_inserter(extremes));
 371
 372         get_intruders(ring, left, right, extremes, intruders, strategy);
 373
 374         return true;
 375     }
 376 };
 377
 378
 379
 380
 381
 382 }} // namespace detail::extreme_points
 383 #endif // DOXYGEN_NO_DETAIL
 384
 385
 386 #ifndef DOXYGEN_NO_DISPATCH
 387 namespace dispatch
 388 {
 389
 390
 391 template
 392 <
 393     typename Geometry,
 394     std::size_t Dimension,
 395     typename GeometryTag = typename tag<Geometry>::type
 396 >
 397 struct extreme_points
 398 {};
 399
 400
 401 template<typename Ring, std::size_t Dimension>
 402 struct extreme_points<Ring, Dimension, ring_tag>
 403     : detail::extreme_points::extreme_points_on_ring<Ring, Dimension>
 404 {};
 405
 406
 407 template<typename Polygon, std::size_t Dimension>
 408 struct extreme_points<Polygon, Dimension, polygon_tag>
 409 {
 410     template <typename Extremes, typename Intruders, typename SideStrategy>
 411     static inline bool apply(Polygon const& polygon, Extremes& extremes, Intruders& intruders,
 412                              SideStrategy const& strategy)
 413     {
 414         typedef typename geometry::ring_type<Polygon>::type ring_type;
 415         typedef detail::extreme_points::extreme_points_on_ring
 416             <
 417                 ring_type, Dimension
 418             > ring_implementation;
 419
 420         if (! ring_implementation::apply(geometry::exterior_ring(polygon),
 421                                          extremes, intruders, strategy))
 422         {
 423             return false;
 424         }
 425
 426         // For a polygon, its interior rings can contain intruders
 427         typename interior_return_type<Polygon const>::type
 428             rings = interior_rings(polygon);
 429         for (typename detail::interior_iterator<Polygon const>::type
 430                 it = boost::begin(rings); it != boost::end(rings); ++it)
 431         {
 432             ring_implementation::get_intruders(*it, extremes,  intruders, strategy);
 433         }
 434
 435         return true;
 436     }
 437 };
 438
 439 template<typename Box>
 440 struct extreme_points<Box, 1, box_tag>
 441 {
 442     template <typename Extremes, typename Intruders, typename SideStrategy>
 443     static inline bool apply(Box const& box, Extremes& extremes, Intruders& ,
 444                              SideStrategy const& )
 445     {
 446         extremes.resize(4);
 447         geometry::detail::assign_box_corners_oriented<false>(box, extremes);
 448         // ll,ul,ur,lr, contains too exactly the right info
 449         return true;
 450     }
 451 };
 452
 453 template<typename Box>
 454 struct extreme_points<Box, 0, box_tag>
 455 {
 456     template <typename Extremes, typename Intruders, typename SideStrategy>
 457     static inline bool apply(Box const& box, Extremes& extremes, Intruders& ,
 458                              SideStrategy const& )
 459     {
 460         extremes.resize(4);
 461         geometry::detail::assign_box_corners_oriented<false>(box, extremes);
 462         // ll,ul,ur,lr, rotate one to start with UL and end with LL
 463         std::rotate(extremes.begin(), extremes.begin() + 1, extremes.end());
 464         return true;
 465     }
 466 };
 467
 468 template<typename MultiPolygon, std::size_t Dimension>
 469 struct extreme_points<MultiPolygon, Dimension, multi_polygon_tag>
 470 {
 471     template <typename Extremes, typename Intruders, typename SideStrategy>
 472     static inline bool apply(MultiPolygon const& multi, Extremes& extremes,
 473                              Intruders& intruders, SideStrategy const& strategy)
 474     {
 475         // Get one for the very first polygon, that is (for the moment) enough.
 476         // It is not guaranteed the "extreme" then, but for the current purpose
 477         // (point_on_surface) it can just be this point.
 478         if (boost::size(multi) >= 1)
 479         {
 480             return extreme_points
 481                 <
 482                     typename boost::range_value<MultiPolygon const>::type,
 483                     Dimension,
 484                     polygon_tag
 485                 >::apply(*boost::begin(multi), extremes, intruders, strategy);
 486         }
 487
 488         return false;
 489     }
 490 };
 491
 492 } // namespace dispatch
 493 #endif // DOXYGEN_NO_DISPATCH
 494
 495
 496 /*!
 497 \brief Returns extreme points (for Edge=1 in dimension 1, so the top,
 498        for Edge=0 in dimension 0, the right side)
 499 \note We could specify a strategy (less/greater) to get bottom/left side too. However, until now we don't need that.
 500  */
 501 template
 502 <
 503     std::size_t Edge,
 504     typename Geometry,
 505     typename Extremes,
 506     typename Intruders,
 507     typename SideStrategy
 508 >
 509 inline bool extreme_points(Geometry const& geometry,
 510                            Extremes& extremes,
 511                            Intruders& intruders,
 512                            SideStrategy const& strategy)
 513 {
 514     concepts::check<Geometry const>();
 515
 516     // Extremes is not required to follow a geometry concept (but it should support an output iterator),
 517     // but its elements should fulfil the point-concept
 518     concepts::check<typename boost::range_value<Extremes>::type>();
 519
 520     // Intruders should contain collections which value type is point-concept
 521     // Extremes might be anything (supporting an output iterator), but its elements should fulfil the point-concept
 522     concepts::check
 523         <
 524             typename boost::range_value
 525                 <
 526                     typename boost::range_value<Intruders>::type
 527                 >::type
 528             const
 529         >();
 530
 531     return dispatch::extreme_points
 532             <
 533                 Geometry,
 534                 Edge
 535             >::apply(geometry, extremes, intruders, strategy);
 536 }
 537
 538
 539
 540 }} // namespace boost::geometry
 541
 542
 543 #endif // BOOST_GEOMETRY_ALGORITHMS_DETAIL_EXTREME_POINTS_HPP