boost/geometry/algorithms/detail/overlay/sort_by_side.hpp

   1 // Boost.Geometry (aka GGL, Generic Geometry Library)
   2
   3 // Copyright (c) 2015 Barend Gehrels, Amsterdam, the Netherlands.
   4 // Copyright (c) 2017 Adam Wulkiewicz, Lodz, Poland.
   5
   6 // This file was modified by Oracle on 2017, 2019.
   7 // Modifications copyright (c) 2017, 2019 Oracle and/or its affiliates.
   8
   9 // Contributed and/or modified by Adam Wulkiewicz, on behalf of Oracle
  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_OVERLAY_SORT_BY_SIDE_HPP
  16 #define BOOST_GEOMETRY_ALGORITHMS_DETAIL_OVERLAY_SORT_BY_SIDE_HPP
  17
  18 #include <algorithm>
  19 #include <map>
  20 #include <vector>
  21
  22 #include <boost/geometry/algorithms/detail/overlay/copy_segment_point.hpp>
  23 #include <boost/geometry/algorithms/detail/overlay/get_ring.hpp>
  24 #include <boost/geometry/algorithms/detail/direction_code.hpp>
  25 #include <boost/geometry/algorithms/detail/overlay/turn_info.hpp>
  26
  27 #include <boost/geometry/util/condition.hpp>
  28
  29 namespace boost { namespace geometry
  30 {
  31
  32 #ifndef DOXYGEN_NO_DETAIL
  33 namespace detail { namespace overlay { namespace sort_by_side
  34 {
  35
  36 enum direction_type { dir_unknown = -1, dir_from = 0, dir_to = 1 };
  37
  38 typedef signed_size_type rank_type;
  39
  40
  41 // Point-wrapper, adding some properties
  42 template <typename Point>
  43 struct ranked_point
  44 {
  45     ranked_point()
  46         : rank(0)
  47         , turn_index(-1)
  48         , operation_index(-1)
  49         , direction(dir_unknown)
  50         , count_left(0)
  51         , count_right(0)
  52         , operation(operation_none)
  53     {}
  54
  55     ranked_point(Point const& p, signed_size_type ti, int oi,
  56                  direction_type d, operation_type op, segment_identifier const& si)
  57         : point(p)
  58         , rank(0)
  59         , zone(-1)
  60         , turn_index(ti)
  61         , operation_index(oi)
  62         , direction(d)
  63         , count_left(0)
  64         , count_right(0)
  65         , operation(op)
  66         , seg_id(si)
  67     {}
  68
  69     Point point;
  70     rank_type rank;
  71     signed_size_type zone; // index of closed zone, in uu turn there would be 2 zones
  72     signed_size_type turn_index;
  73     int operation_index; // 0,1
  74     direction_type direction;
  75     std::size_t count_left;
  76     std::size_t count_right;
  77     operation_type operation;
  78     segment_identifier seg_id;
  79 };
  80
  81 struct less_by_turn_index
  82 {
  83     template <typename T>
  84     inline bool operator()(const T& first, const T& second) const
  85     {
  86         return first.turn_index == second.turn_index
  87             ? first.index < second.index
  88             : first.turn_index < second.turn_index
  89             ;
  90     }
  91 };
  92
  93 struct less_by_index
  94 {
  95     template <typename T>
  96     inline bool operator()(const T& first, const T& second) const
  97     {
  98         // Length might be considered too
  99         // First order by from/to
 100         if (first.direction != second.direction)
 101         {
 102             return first.direction < second.direction;
 103         }
 104         // Then by turn index
 105         if (first.turn_index != second.turn_index)
 106         {
 107             return first.turn_index < second.turn_index;
 108         }
 109         // This can also be the same (for example in buffer), but seg_id is
 110         // never the same
 111         return first.seg_id < second.seg_id;
 112     }
 113 };
 114
 115 struct less_false
 116 {
 117     template <typename T>
 118     inline bool operator()(const T&, const T& ) const
 119     {
 120         return false;
 121     }
 122 };
 123
 124 template <typename Point, typename SideStrategy, typename LessOnSame, typename Compare>
 125 struct less_by_side
 126 {
 127     less_by_side(const Point& p1, const Point& p2, SideStrategy const& strategy)
 128         : m_origin(p1)
 129         , m_turn_point(p2)
 130         , m_strategy(strategy)
 131     {}
 132
 133     template <typename T>
 134     inline bool operator()(const T& first, const T& second) const
 135     {
 136         typedef typename SideStrategy::cs_tag cs_tag;
 137
 138         LessOnSame on_same;
 139         Compare compare;
 140
 141         int const side_first = m_strategy.apply(m_origin, m_turn_point, first.point);
 142         int const side_second = m_strategy.apply(m_origin, m_turn_point, second.point);
 143
 144         if (side_first == 0 && side_second == 0)
 145         {
 146             // Both collinear. They might point into different directions: <------*------>
 147             // If so, order the one going backwards as the very first.
 148
 149             int const first_code = direction_code<cs_tag>(m_origin, m_turn_point, first.point);
 150             int const second_code = direction_code<cs_tag>(m_origin, m_turn_point, second.point);
 151
 152             // Order by code, backwards first, then forward.
 153             return first_code != second_code
 154                 ? first_code < second_code
 155                 : on_same(first, second)
 156                 ;
 157         }
 158         else if (side_first == 0
 159                 && direction_code<cs_tag>(m_origin, m_turn_point, first.point) == -1)
 160         {
 161             // First collinear and going backwards.
 162             // Order as the very first, so return always true
 163             return true;
 164         }
 165         else if (side_second == 0
 166             && direction_code<cs_tag>(m_origin, m_turn_point, second.point) == -1)
 167         {
 168             // Second is collinear and going backwards
 169             // Order as very last, so return always false
 170             return false;
 171         }
 172
 173         // They are not both collinear
 174
 175         if (side_first != side_second)
 176         {
 177             return compare(side_first, side_second);
 178         }
 179
 180         // They are both left, both right, and/or both collinear (with each other and/or with p1,p2)
 181         // Check mutual side
 182         int const side_second_wrt_first = m_strategy.apply(m_turn_point, first.point, second.point);
 183
 184         if (side_second_wrt_first == 0)
 185         {
 186             return on_same(first, second);
 187         }
 188
 189         int const side_first_wrt_second = -side_second_wrt_first;
 190
 191         // Both are on same side, and not collinear
 192         // Union: return true if second is right w.r.t. first, so -1,
 193         // so other is 1. union has greater as compare functor
 194         // Intersection: v.v.
 195         return compare(side_first_wrt_second, side_second_wrt_first);
 196     }
 197
 198 private :
 199     Point const& m_origin;
 200     Point const& m_turn_point;
 201     SideStrategy const& m_strategy;
 202 };
 203
 204 // Sorts vectors in counter clockwise order (by default)
 205 template
 206 <
 207     bool Reverse1,
 208     bool Reverse2,
 209     overlay_type OverlayType,
 210     typename Point,
 211     typename SideStrategy,
 212     typename Compare
 213 >
 214 struct side_sorter
 215 {
 216     typedef ranked_point<Point> rp;
 217
 218 private :
 219     struct include_union
 220     {
 221         inline bool operator()(rp const& ranked_point) const
 222         {
 223             // New candidate if there are no polygons on left side,
 224             // but there are on right side
 225             return ranked_point.count_left == 0
 226                 && ranked_point.count_right > 0;
 227         }
 228     };
 229
 230     struct include_intersection
 231     {
 232         inline bool operator()(rp const& ranked_point) const
 233         {
 234             // New candidate if there are two polygons on right side,
 235             // and less on the left side
 236             return ranked_point.count_left < 2
 237                 && ranked_point.count_right >= 2;
 238         }
 239     };
 240
 241 public :
 242     side_sorter(SideStrategy const& strategy)
 243         : m_origin_count(0)
 244         , m_origin_segment_distance(0)
 245         , m_strategy(strategy)
 246     {}
 247
 248     void add_segment_from(signed_size_type turn_index, int op_index,
 249             Point const& point_from,
 250             operation_type op, segment_identifier const& si,
 251             bool is_origin)
 252     {
 253         m_ranked_points.push_back(rp(point_from, turn_index, op_index, dir_from, op, si));
 254         if (is_origin)
 255         {
 256             m_origin = point_from;
 257             m_origin_count++;
 258         }
 259     }
 260
 261     void add_segment_to(signed_size_type turn_index, int op_index,
 262             Point const& point_to,
 263             operation_type op, segment_identifier const& si)
 264     {
 265         m_ranked_points.push_back(rp(point_to, turn_index, op_index, dir_to, op, si));
 266     }
 267
 268     void add_segment(signed_size_type turn_index, int op_index,
 269             Point const& point_from, Point const& point_to,
 270             operation_type op, segment_identifier const& si,
 271             bool is_origin)
 272     {
 273         add_segment_from(turn_index, op_index, point_from, op, si, is_origin);
 274         add_segment_to(turn_index, op_index, point_to, op, si);
 275     }
 276
 277     template <typename Operation, typename Geometry1, typename Geometry2>
 278     Point add(Operation const& op, signed_size_type turn_index, int op_index,
 279             Geometry1 const& geometry1,
 280             Geometry2 const& geometry2,
 281             bool is_origin)
 282     {
 283         Point point1, point2, point3;
 284         geometry::copy_segment_points<Reverse1, Reverse2>(geometry1, geometry2,
 285                 op.seg_id, point1, point2, point3);
 286         Point const& point_to = op.fraction.is_one() ? point3 : point2;
 287         add_segment(turn_index, op_index, point1, point_to, op.operation, op.seg_id, is_origin);
 288         return point1;
 289     }
 290
 291     template <typename Operation, typename Geometry1, typename Geometry2>
 292     void add(Operation const& op, signed_size_type turn_index, int op_index,
 293             segment_identifier const& departure_seg_id,
 294             Geometry1 const& geometry1,
 295             Geometry2 const& geometry2,
 296             bool check_origin)
 297     {
 298         Point const point1 = add(op, turn_index, op_index, geometry1, geometry2, false);
 299
 300         if (check_origin)
 301         {
 302             bool const is_origin
 303                     = op.seg_id.source_index == departure_seg_id.source_index
 304                     && op.seg_id.ring_index == departure_seg_id.ring_index
 305                     && op.seg_id.multi_index == departure_seg_id.multi_index;
 306
 307             if (is_origin)
 308             {
 309                 signed_size_type const segment_distance = calculate_segment_distance(op, departure_seg_id, geometry1, geometry2);
 310                 if (m_origin_count == 0 ||
 311                         segment_distance < m_origin_segment_distance)
 312                 {
 313                     m_origin = point1;
 314                     m_origin_segment_distance = segment_distance;
 315                 }
 316                 m_origin_count++;
 317             }
 318         }
 319     }
 320
 321     template <typename Operation, typename Geometry1, typename Geometry2>
 322     static signed_size_type calculate_segment_distance(Operation const& op,
 323             segment_identifier const& departure_seg_id,
 324             Geometry1 const& geometry1,
 325             Geometry2 const& geometry2)
 326     {
 327         if (op.seg_id.segment_index >= departure_seg_id.segment_index)
 328         {
 329             // dep.seg_id=5, op.seg_id=7, distance=2, being segments 5,6
 330             return op.seg_id.segment_index - departure_seg_id.segment_index;
 331         }
 332         // Take wrap into account
 333         // Suppose point_count=10 (10 points, 9 segments), dep.seg_id=7, op.seg_id=2,
 334         // then distance=9-7+2=4, being segments 7,8,0,1
 335         std::size_t const segment_count
 336                     = op.seg_id.source_index == 0
 337                     ? segment_count_on_ring(geometry1, op.seg_id)
 338                     : segment_count_on_ring(geometry2, op.seg_id);
 339         return segment_count - departure_seg_id.segment_index + op.seg_id.segment_index;
 340     }
 341
 342     void apply(Point const& turn_point)
 343     {
 344         // We need three compare functors:
 345         // 1) to order clockwise (union) or counter clockwise (intersection)
 346         // 2) to order by side, resulting in unique ranks for all points
 347         // 3) to order by side, resulting in non-unique ranks
 348         //    to give colinear points
 349
 350         // Sort by side and assign rank
 351         less_by_side<Point, SideStrategy, less_by_index, Compare> less_unique(m_origin, turn_point, m_strategy);
 352         less_by_side<Point, SideStrategy, less_false, Compare> less_non_unique(m_origin, turn_point, m_strategy);
 353
 354         std::sort(m_ranked_points.begin(), m_ranked_points.end(), less_unique);
 355
 356         std::size_t colinear_rank = 0;
 357         for (std::size_t i = 0; i < m_ranked_points.size(); i++)
 358         {
 359             if (i > 0
 360                 && less_non_unique(m_ranked_points[i - 1], m_ranked_points[i]))
 361             {
 362                 // It is not collinear
 363                 colinear_rank++;
 364             }
 365
 366             m_ranked_points[i].rank = colinear_rank;
 367         }
 368     }
 369
 370     template <signed_size_type segment_identifier::*Member, typename Map>
 371     void find_open_generic(Map& handled, bool check)
 372     {
 373         for (std::size_t i = 0; i < m_ranked_points.size(); i++)
 374         {
 375             const rp& ranked = m_ranked_points[i];
 376             if (ranked.direction != dir_from)
 377             {
 378                 continue;
 379             }
 380
 381             signed_size_type const& index = ranked.seg_id.*Member;
 382             if (check && (index < 0 || index > 1))
 383             {
 384                 // Should not occur
 385                 continue;
 386             }
 387             if (! handled[index])
 388             {
 389                 find_polygons_for_source<Member>(index, i);
 390                 handled[index] = true;
 391             }
 392         }
 393     }
 394
 395     void find_open()
 396     {
 397         if (BOOST_GEOMETRY_CONDITION(OverlayType == overlay_buffer))
 398         {
 399             // For buffers, use piece index
 400             std::map<signed_size_type, bool> handled;
 401             find_open_generic
 402                 <
 403                     &segment_identifier::piece_index
 404                 >(handled, false);
 405         }
 406         else
 407         {
 408             // For other operations, by source (there should only source 0,1)
 409             bool handled[2] = {false, false};
 410             find_open_generic
 411                 <
 412                     &segment_identifier::source_index
 413                 >(handled, true);
 414         }
 415     }
 416
 417     void reverse()
 418     {
 419         if (m_ranked_points.empty())
 420         {
 421             return;
 422         }
 423
 424         std::size_t const last = 1 + m_ranked_points.back().rank;
 425
 426         // Move iterator after rank==0
 427         bool has_first = false;
 428         typename container_type::iterator it = m_ranked_points.begin() + 1;
 429         for (; it != m_ranked_points.end() && it->rank == 0; ++it)
 430         {
 431             has_first = true;
 432         }
 433
 434         if (has_first)
 435         {
 436             // Reverse first part (having rank == 0), if any,
 437             // but skip the very first row
 438             std::reverse(m_ranked_points.begin() + 1, it);
 439             for (typename container_type::iterator fit = m_ranked_points.begin();
 440                  fit != it; ++fit)
 441             {
 442                 BOOST_ASSERT(fit->rank == 0);
 443             }
 444         }
 445
 446         // Reverse the rest (main rank > 0)
 447         std::reverse(it, m_ranked_points.end());
 448         for (; it != m_ranked_points.end(); ++it)
 449         {
 450             BOOST_ASSERT(it->rank > 0);
 451             it->rank = last - it->rank;
 452         }
 453     }
 454
 455     bool has_origin() const
 456     {
 457         return m_origin_count > 0;
 458     }
 459
 460 //private :
 461
 462     typedef std::vector<rp> container_type;
 463     container_type m_ranked_points;
 464     Point m_origin;
 465     std::size_t m_origin_count;
 466     signed_size_type m_origin_segment_distance;
 467     SideStrategy m_strategy;
 468
 469 private :
 470
 471     //! Check how many open spaces there are
 472     template <typename Include>
 473     inline std::size_t open_count(Include const& include_functor) const
 474     {
 475         std::size_t result = 0;
 476         rank_type last_rank = 0;
 477         for (std::size_t i = 0; i < m_ranked_points.size(); i++)
 478         {
 479             rp const& ranked_point = m_ranked_points[i];
 480
 481             if (ranked_point.rank > last_rank
 482                 && ranked_point.direction == sort_by_side::dir_to
 483                 && include_functor(ranked_point))
 484             {
 485                 result++;
 486                 last_rank = ranked_point.rank;
 487             }
 488         }
 489         return result;
 490     }
 491
 492     std::size_t move(std::size_t index) const
 493     {
 494         std::size_t const result = index + 1;
 495         return result >= m_ranked_points.size() ? 0 : result;
 496     }
 497
 498     //! member is pointer to member (source_index or multi_index)
 499     template <signed_size_type segment_identifier::*Member>
 500     std::size_t move(signed_size_type member_index, std::size_t index) const
 501     {
 502         std::size_t result = move(index);
 503         while (m_ranked_points[result].seg_id.*Member != member_index)
 504         {
 505             result = move(result);
 506         }
 507         return result;
 508     }
 509
 510     void assign_ranks(rank_type min_rank, rank_type max_rank, int side_index)
 511     {
 512         for (std::size_t i = 0; i < m_ranked_points.size(); i++)
 513         {
 514             rp& ranked = m_ranked_points[i];
 515             // Suppose there are 8 ranks, if min=4,max=6: assign 4,5,6
 516             // if min=5,max=2: assign from 5,6,7,1,2
 517             bool const in_range
 518                 = max_rank >= min_rank
 519                 ? ranked.rank >= min_rank && ranked.rank <= max_rank
 520                 : ranked.rank >= min_rank || ranked.rank <= max_rank
 521                 ;
 522
 523             if (in_range)
 524             {
 525                 if (side_index == 1)
 526                 {
 527                     ranked.count_left++;
 528                 }
 529                 else if (side_index == 2)
 530                 {
 531                     ranked.count_right++;
 532                 }
 533             }
 534         }
 535     }
 536
 537     template <signed_size_type segment_identifier::*Member>
 538     void find_polygons_for_source(signed_size_type the_index,
 539                 std::size_t start_index)
 540     {
 541         bool in_polygon = true; // Because start_index is "from", arrives at the turn
 542         rp const& start_rp = m_ranked_points[start_index];
 543         rank_type last_from_rank = start_rp.rank;
 544         rank_type previous_rank = start_rp.rank;
 545
 546         for (std::size_t index = move<Member>(the_index, start_index);
 547              ;
 548              index = move<Member>(the_index, index))
 549         {
 550             rp& ranked = m_ranked_points[index];
 551
 552             if (ranked.rank != previous_rank && ! in_polygon)
 553             {
 554                 assign_ranks(last_from_rank, previous_rank - 1, 1);
 555                 assign_ranks(last_from_rank + 1, previous_rank, 2);
 556             }
 557
 558             if (index == start_index)
 559             {
 560                 return;
 561             }
 562
 563             if (ranked.direction == dir_from)
 564             {
 565                 last_from_rank = ranked.rank;
 566                 in_polygon = true;
 567             }
 568             else if (ranked.direction == dir_to)
 569             {
 570                 in_polygon = false;
 571             }
 572
 573             previous_rank = ranked.rank;
 574         }
 575     }
 576
 577     //! Find closed zones and assign it
 578     template <typename Include>
 579     std::size_t assign_zones(Include const& include_functor)
 580     {
 581         // Find a starting point (the first rank after an outgoing rank
 582         // with no polygons on the left side)
 583         rank_type start_rank = m_ranked_points.size() + 1;
 584         std::size_t start_index = 0;
 585         rank_type max_rank = 0;
 586         for (std::size_t i = 0; i < m_ranked_points.size(); i++)
 587         {
 588             rp const& ranked_point = m_ranked_points[i];
 589             if (ranked_point.rank > max_rank)
 590             {
 591                 max_rank = ranked_point.rank;
 592             }
 593             if (ranked_point.direction == sort_by_side::dir_to
 594                 && include_functor(ranked_point))
 595             {
 596                 start_rank = ranked_point.rank + 1;
 597             }
 598             if (ranked_point.rank == start_rank && start_index == 0)
 599             {
 600                 start_index = i;
 601             }
 602         }
 603
 604         // Assign the zones
 605         rank_type const undefined_rank = max_rank + 1;
 606         std::size_t zone_id = 0;
 607         rank_type last_rank = 0;
 608         rank_type rank_at_next_zone = undefined_rank;
 609         std::size_t index = start_index;
 610         for (std::size_t i = 0; i < m_ranked_points.size(); i++)
 611         {
 612             rp& ranked_point = m_ranked_points[index];
 613
 614             // Implement cyclic behavior
 615             index++;
 616             if (index == m_ranked_points.size())
 617             {
 618                 index = 0;
 619             }
 620
 621             if (ranked_point.rank != last_rank)
 622             {
 623                 if (ranked_point.rank == rank_at_next_zone)
 624                 {
 625                     zone_id++;
 626                     rank_at_next_zone = undefined_rank;
 627                 }
 628
 629                 if (ranked_point.direction == sort_by_side::dir_to
 630                     && include_functor(ranked_point))
 631                 {
 632                     rank_at_next_zone = ranked_point.rank + 1;
 633                     if (rank_at_next_zone > max_rank)
 634                     {
 635                         rank_at_next_zone = 0;
 636                     }
 637                 }
 638
 639                 last_rank = ranked_point.rank;
 640             }
 641
 642             ranked_point.zone = zone_id;
 643         }
 644         return zone_id;
 645     }
 646
 647 public :
 648     inline std::size_t open_count(operation_type for_operation) const
 649     {
 650         return for_operation == operation_union
 651             ? open_count(include_union())
 652             : open_count(include_intersection())
 653             ;
 654     }
 655
 656     inline std::size_t assign_zones(operation_type for_operation)
 657     {
 658         return for_operation == operation_union
 659             ? assign_zones(include_union())
 660             : assign_zones(include_intersection())
 661             ;
 662     }
 663
 664 };
 665
 666
 667 //! Metafunction to define side_order (clockwise, ccw) by operation_type
 668 template <operation_type OpType>
 669 struct side_compare {};
 670
 671 template <>
 672 struct side_compare<operation_union>
 673 {
 674     typedef std::greater<int> type;
 675 };
 676
 677 template <>
 678 struct side_compare<operation_intersection>
 679 {
 680     typedef std::less<int> type;
 681 };
 682
 683
 684 }}} // namespace detail::overlay::sort_by_side
 685 #endif //DOXYGEN_NO_DETAIL
 686
 687
 688 }} // namespace boost::geometry
 689
 690 #endif // BOOST_GEOMETRY_ALGORITHMS_DETAIL_OVERLAY_SORT_BY_SIDE_HPP