detail/overlay/enrich_intersection_points.hpp

0001 // Boost.Geometry (aka GGL, Generic Geometry Library)
0002
0003 // Copyright (c) 2007-2012 Barend Gehrels, Amsterdam, the Netherlands.
0004 // Copyright (c) 2017 Adam Wulkiewicz, Lodz, Poland.
0005
0006 // This file was modified by Oracle on 2017-2021.
0007 // Modifications copyright (c) 2017-2020 Oracle and/or its affiliates.
0008 // Contributed and/or modified by Adam Wulkiewicz, on behalf of Oracle
0009
0010 // Use, modification and distribution is subject to the Boost Software License,
0011 // Version 1.0. (See accompanying file LICENSE_1_0.txt or copy at
0012 // http://www.boost.org/LICENSE_1_0.txt)
0013
0014 #ifndef BOOST_GEOMETRY_ALGORITHMS_DETAIL_OVERLAY_ENRICH_HPP
0015 #define BOOST_GEOMETRY_ALGORITHMS_DETAIL_OVERLAY_ENRICH_HPP
0016
0017 #include <cstddef>
0018 #include <algorithm>
0019 #include <map>
0020 #include <set>
0021 #include <vector>
0022
0023 #ifdef BOOST_GEOMETRY_DEBUG_ENRICH
0024 #  include <iostream>
0025 #  include <boost/geometry/algorithms/detail/overlay/debug_turn_info.hpp>
0026 #  include <boost/geometry/io/wkt/wkt.hpp>
0027 #  if ! defined(BOOST_GEOMETRY_DEBUG_IDENTIFIER)
0028 #    define BOOST_GEOMETRY_DEBUG_IDENTIFIER
0029   #endif
0030 #endif
0031
0032 #include <boost/range/begin.hpp>
0033 #include <boost/range/end.hpp>
0034 #include <boost/range/value_type.hpp>
0035
0036 #include <boost/geometry/algorithms/detail/ring_identifier.hpp>
0037 #include <boost/geometry/algorithms/detail/overlay/discard_duplicate_turns.hpp>
0038 #include <boost/geometry/algorithms/detail/overlay/handle_colocations.hpp>
0039 #include <boost/geometry/algorithms/detail/overlay/handle_self_turns.hpp>
0040 #include <boost/geometry/algorithms/detail/overlay/is_self_turn.hpp>
0041 #include <boost/geometry/algorithms/detail/overlay/less_by_segment_ratio.hpp>
0042 #include <boost/geometry/algorithms/detail/overlay/overlay_type.hpp>
0043 #include <boost/geometry/policies/robustness/robust_type.hpp>
0044 #include <boost/geometry/util/for_each_with_index.hpp>
0045
0046 #ifdef BOOST_GEOMETRY_DEBUG_ENRICH
0047 #  include <boost/geometry/algorithms/detail/overlay/check_enrich.hpp>
0048 #endif
0049
0050
0051 namespace boost { namespace geometry
0052 {
0053
0054 #ifndef DOXYGEN_NO_DETAIL
0055 namespace detail { namespace overlay
0056 {
0057
0058 template <typename Turns>
0059 struct discarded_indexed_turn
0060 {
0061     discarded_indexed_turn(Turns const& turns)
0062         : m_turns(turns)
0063     {}
0064
0065     template <typename IndexedTurn>
0066     inline bool operator()(IndexedTurn const& indexed) const
0067     {
0068         return m_turns[indexed.turn_index].discarded;
0069     }
0070
0071     Turns const& m_turns;
0072 };
0073
0074 // Sorts IP-s of this ring on segment-identifier, and if on same segment,
0075 //  on distance.
0076 // Then assigns for each IP which is the next IP on this segment,
0077 // plus the vertex-index to travel to, plus the next IP
0078 // (might be on another segment)
0079 template
0080 <
0081     bool Reverse1, bool Reverse2,
0082     typename Operations,
0083     typename Turns,
0084     typename Geometry1, typename Geometry2,
0085     typename RobustPolicy,
0086     typename Strategy
0087 >
0088 inline void enrich_sort(Operations& operations,
0089             Turns const& turns,
0090             Geometry1 const& geometry1,
0091             Geometry2 const& geometry2,
0092             RobustPolicy const& robust_policy,
0093             Strategy const& strategy)
0094 {
0095     std::sort(std::begin(operations),
0096               std::end(operations),
0097               less_by_segment_ratio
0098                 <
0099                     Turns,
0100                     typename boost::range_value<Operations>::type,
0101                     Geometry1, Geometry2,
0102                     RobustPolicy,
0103                     Strategy,
0104                     Reverse1, Reverse2
0105                 >(turns, geometry1, geometry2, robust_policy, strategy));
0106 }
0107
0108
0109 template <typename Operations, typename Turns>
0110 inline void enrich_assign(Operations& operations, Turns& turns,
0111                           bool check_turns)
0112 {
0113     if (operations.empty())
0114     {
0115         return;
0116     }
0117
0118     // Assign travel-to-vertex/ip index for each turning point.
0119     // Iterator "next" is circular
0120
0121     geometry::ever_circling_range_iterator<Operations const> next(operations);
0122     ++next;
0123
0124     for (auto const& indexed : operations)
0125     {
0126         auto& turn = turns[indexed.turn_index];
0127         auto& op = turn.operations[indexed.operation_index];
0128
0129         if (check_turns && indexed.turn_index == next->turn_index)
0130         {
0131             // Normal behaviour: next points at next turn, increase next.
0132             // For dissolve this should not be done, turn_index is often
0133             // the same for two consecutive operations
0134             ++next;
0135         }
0136
0137         // Cluster behaviour: next should point after cluster, unless
0138         // their seg_ids are not the same
0139         // (For dissolve, this is still to be examined - TODO)
0140         while (turn.is_clustered()
0141                 && indexed.turn_index != next->turn_index
0142                 && turn.cluster_id == turns[next->turn_index].cluster_id
0143                 && op.seg_id == turns[next->turn_index].operations[next->operation_index].seg_id)
0144         {
0145             ++next;
0146         }
0147
0148         auto const& next_turn = turns[next->turn_index];
0149         auto const& next_op = next_turn.operations[next->operation_index];
0150
0151         op.enriched.travels_to_ip_index
0152                 = static_cast<signed_size_type>(next->turn_index);
0153         op.enriched.travels_to_vertex_index
0154                 = next->subject->seg_id.segment_index;
0155
0156         if (op.seg_id.segment_index == next_op.seg_id.segment_index
0157                 && op.fraction < next_op.fraction)
0158         {
0159             // Next turn is located further on same segment
0160             // assign next_ip_index
0161             // (this is one not circular therefore fraction is considered)
0162             op.enriched.next_ip_index = static_cast<signed_size_type>(next->turn_index);
0163         }
0164
0165         if (! check_turns)
0166         {
0167             ++next;
0168         }
0169     }
0170
0171     // DEBUG
0172 #ifdef BOOST_GEOMETRY_DEBUG_ENRICH
0173     for (auto const& indexed_op : operations)
0174     {
0175         auto const& op = turns[indexed_op.turn_index].operations[indexed_op.operation_index];
0176
0177         std::cout << indexed_op.turn_index
0178             << " cl=" << turns[indexed_op.turn_index].cluster_id
0179             << " meth=" << method_char(turns[indexed_op.turn_index].method)
0180             << " seg=" << op.seg_id
0181             << " dst=" << op.fraction // needs define
0182             << " op=" << operation_char(turns[indexed_op.turn_index].operations[0].operation)
0183             << operation_char(turns[indexed_op.turn_index].operations[1].operation)
0184             << " (" << operation_char(op.operation) << ")"
0185             << " nxt=" << op.enriched.next_ip_index
0186             << " / " << op.enriched.travels_to_ip_index
0187             << " [vx " << op.enriched.travels_to_vertex_index << "]"
0188             << (turns[indexed_op.turn_index].discarded ? " [discarded]" : "")
0189             << (op.enriched.startable ? "" : " [not startable]")
0190             << std::endl;
0191     }
0192 #endif
0193     // END DEBUG
0194
0195 }
0196
0197 template <typename Operations, typename Turns>
0198 inline void enrich_adapt(Operations& operations, Turns& turns)
0199 {
0200     // Operations is a vector of indexed_turn_operation<>
0201     // If it is empty, or contains one or two items, it makes no sense
0202     if (operations.size() < 3)
0203     {
0204         return;
0205     }
0206
0207     bool next_phase = false;
0208     std::size_t previous_index = operations.size() - 1;
0209
0210     for_each_with_index(operations, [&](std::size_t index, auto const& indexed)
0211     {
0212         auto& turn = turns[indexed.turn_index];
0213         auto& op = turn.operations[indexed.operation_index];
0214
0215         std::size_t const next_index = index + 1 < operations.size() ? index + 1 : 0;
0216         auto const& next_turn = turns[operations[next_index].turn_index];
0217         auto const& next_op = next_turn.operations[operations[next_index].operation_index];
0218
0219         if (op.seg_id.segment_index == next_op.seg_id.segment_index)
0220         {
0221             auto const& prev_turn = turns[operations[previous_index].turn_index];
0222             auto const& prev_op = prev_turn.operations[operations[previous_index].operation_index];
0223             if (op.seg_id.segment_index == prev_op.seg_id.segment_index)
0224             {
0225                 op.enriched.startable = false;
0226                 next_phase = true;
0227             }
0228         }
0229         previous_index = index;
0230     });
0231
0232     if (! next_phase)
0233     {
0234         return;
0235     }
0236
0237     // Discard turns which are both non-startable
0238     next_phase = false;
0239     for (auto& turn : turns)
0240     {
0241         if (! turn.operations[0].enriched.startable
0242             && ! turn.operations[1].enriched.startable)
0243         {
0244             turn.discarded = true;
0245             next_phase = true;
0246         }
0247     }
0248
0249     if (! next_phase)
0250     {
0251         return;
0252     }
0253
0254     // Remove discarded turns from operations to avoid having them as next turn
0255     discarded_indexed_turn<Turns> const predicate(turns);
0256     operations.erase(std::remove_if(std::begin(operations),
0257         std::end(operations), predicate), std::end(operations));
0258 }
0259
0260 struct enriched_map_default_include_policy
0261 {
0262     template <typename Operation>
0263     static inline bool include(Operation const& )
0264     {
0265         // By default include all operations
0266         return true;
0267     }
0268 };
0269
0270 // Add all (non discarded) operations on this ring
0271 // Blocked operations or uu on clusters (for intersection)
0272 // should be included, to block potential paths in clusters
0273 template <typename Turns, typename MappedVector, typename IncludePolicy>
0274 inline void create_map(Turns const& turns, MappedVector& mapped_vector,
0275                        IncludePolicy const& include_policy)
0276 {
0277     for_each_with_index(turns, [&](std::size_t index, auto const& turn)
0278     {
0279         if (! turn.discarded)
0280         {
0281             for_each_with_index(turn.operations, [&](std::size_t op_index, auto const& op)
0282             {
0283                 if (include_policy.include(op.operation))
0284                 {
0285                     ring_identifier const ring_id
0286                         (
0287                             op.seg_id.source_index,
0288                             op.seg_id.multi_index,
0289                             op.seg_id.ring_index
0290                         );
0291                     mapped_vector[ring_id].emplace_back
0292                         (
0293                             index, op_index, op, turn.operations[1 - op_index].seg_id
0294                         );
0295                 }
0296             });
0297         }
0298     });
0299 }
0300
0301 template <typename Point1, typename Point2>
0302 inline typename geometry::coordinate_type<Point1>::type
0303         distance_measure(Point1 const& a, Point2 const& b)
0304 {
0305     // TODO: use comparable distance for point-point instead - but that
0306     // causes currently cycling include problems
0307     using ctype = typename geometry::coordinate_type<Point1>::type;
0308     ctype const dx = get<0>(a) - get<0>(b);
0309     ctype const dy = get<1>(a) - get<1>(b);
0310     return dx * dx + dy * dy;
0311 }
0312
0313 template <typename Turns>
0314 inline void calculate_remaining_distance(Turns& turns)
0315 {
0316     for (auto& turn : turns)
0317     {
0318         auto& op0 = turn.operations[0];
0319         auto& op1 = turn.operations[1];
0320
0321         static decltype(op0.remaining_distance) const zero_distance = 0;
0322
0323         if (op0.remaining_distance != zero_distance
0324             || op1.remaining_distance != zero_distance)
0325         {
0326             continue;
0327         }
0328
0329         auto const to_index0 = op0.enriched.get_next_turn_index();
0330         auto const to_index1 = op1.enriched.get_next_turn_index();
0331         if (to_index0 >= 0
0332             && to_index1 >= 0
0333             && to_index0 != to_index1)
0334         {
0335             op0.remaining_distance = distance_measure(turn.point, turns[to_index0].point);
0336             op1.remaining_distance = distance_measure(turn.point, turns[to_index1].point);
0337         }
0338     }
0339 }
0340
0341
0342 }} // namespace detail::overlay
0343 #endif //DOXYGEN_NO_DETAIL
0344
0345
0346
0347 /*!
0348 \brief All intersection points are enriched with successor information
0349 \ingroup overlay
0350 \tparam Turns type of intersection container
0351             (e.g. vector of "intersection/turn point"'s)
0352 \tparam Clusters type of cluster container
0353 \tparam Geometry1 \tparam_geometry
0354 \tparam Geometry2 \tparam_geometry
0355 \tparam PointInGeometryStrategy point in geometry strategy type
0356 \param turns container containing intersection points
0357 \param clusters container containing clusters
0358 \param geometry1 \param_geometry
0359 \param geometry2 \param_geometry
0360 \param robust_policy policy to handle robustness issues
0361 \param strategy point in geometry strategy
0362  */
0363 template
0364 <
0365     bool Reverse1, bool Reverse2,
0366     overlay_type OverlayType,
0367     typename Turns,
0368     typename Clusters,
0369     typename Geometry1, typename Geometry2,
0370     typename RobustPolicy,
0371     typename IntersectionStrategy
0372 >
0373 inline void enrich_intersection_points(Turns& turns,
0374     Clusters& clusters,
0375     Geometry1 const& geometry1, Geometry2 const& geometry2,
0376     RobustPolicy const& robust_policy,
0377     IntersectionStrategy const& strategy)
0378 {
0379     constexpr detail::overlay::operation_type target_operation
0380             = detail::overlay::operation_from_overlay<OverlayType>::value;
0381     constexpr detail::overlay::operation_type opposite_operation
0382             = target_operation == detail::overlay::operation_union
0383             ? detail::overlay::operation_intersection
0384             : detail::overlay::operation_union;
0385     constexpr bool is_dissolve = OverlayType == overlay_dissolve;
0386
0387     using turn_type = typename boost::range_value<Turns>::type;
0388     using indexed_turn_operation = detail::overlay::indexed_turn_operation
0389         <
0390             typename turn_type::turn_operation_type
0391         > ;
0392
0393     using mapped_vector_type = std::map
0394         <
0395             ring_identifier,
0396             std::vector<indexed_turn_operation>
0397         >;
0398
0399     // From here on, turn indexes are used (in clusters, next_index, etc)
0400     // and turns may not be DELETED - they may only be flagged as discarded
0401     discard_duplicate_start_turns(turns, geometry1, geometry2);
0402
0403     bool has_cc = false;
0404     bool const has_colocations
0405         = detail::overlay::handle_colocations
0406             <
0407                 Reverse1, Reverse2, OverlayType, Geometry1, Geometry2
0408             >(turns, clusters, robust_policy);
0409
0410     // Discard turns not part of target overlay
0411     for (auto& turn : turns)
0412     {
0413         if (turn.both(detail::overlay::operation_none)
0414             || turn.both(opposite_operation)
0415             || turn.both(detail::overlay::operation_blocked)
0416             || (detail::overlay::is_self_turn<OverlayType>(turn)
0417                 && ! turn.is_clustered()
0418                 && ! turn.both(target_operation)))
0419         {
0420             // For all operations, discard xx and none/none
0421             // For intersections, remove uu to avoid the need to travel
0422             // a union (during intersection) in uu/cc clusters (e.g. #31,#32,#33)
0423             // The ux is necessary to indicate impossible paths
0424             // (especially if rescaling is removed)
0425
0426             // Similarly, for union, discard ii and ix
0427
0428             // For self-turns, only keep uu / ii
0429
0430             turn.discarded = true;
0431             turn.cluster_id = -1;
0432             continue;
0433         }
0434
0435         if (! turn.discarded
0436             && turn.both(detail::overlay::operation_continue))
0437         {
0438             has_cc = true;
0439         }
0440     }
0441
0442     if (! is_dissolve)
0443     {
0444         detail::overlay::discard_closed_turns
0445             <
0446                 OverlayType,
0447                 target_operation
0448             >::apply(turns, clusters, geometry1, geometry2,
0449                      strategy);
0450         detail::overlay::discard_open_turns
0451             <
0452                 OverlayType,
0453                 target_operation
0454             >::apply(turns, clusters, geometry1, geometry2,
0455                      strategy);
0456     }
0457
0458     // Create a map of vectors of indexed operation-types to be able
0459     // to sort intersection points PER RING
0460     mapped_vector_type mapped_vector;
0461
0462     detail::overlay::create_map(turns, mapped_vector,
0463                                 detail::overlay::enriched_map_default_include_policy());
0464
0465     for (auto& pair : mapped_vector)
0466     {
0467         detail::overlay::enrich_sort<Reverse1, Reverse2>(
0468                     pair.second, turns,
0469                     geometry1, geometry2,
0470                     robust_policy, strategy);
0471 #ifdef BOOST_GEOMETRY_DEBUG_ENRICH
0472         std::cout << "ENRICH-sort Ring " << pair.first << std::endl;
0473         for (auto const& op : pair.second)
0474         {
0475             std::cout << op.turn_index << " " << op.operation_index << std::endl;
0476         }
0477 #endif
0478     }
0479
0480     if (has_colocations)
0481     {
0482         // First gather cluster properties (using even clusters with
0483         // discarded turns - for open turns), then clean up clusters
0484         detail::overlay::gather_cluster_properties
0485             <
0486                 Reverse1,
0487                 Reverse2,
0488                 OverlayType
0489             >(clusters, turns, target_operation,
0490               geometry1, geometry2, strategy.side()); // TODO: pass strategy
0491
0492         detail::overlay::cleanup_clusters(turns, clusters);
0493     }
0494
0495     // After cleaning up clusters assign the next turns
0496
0497     for (auto& pair : mapped_vector)
0498     {
0499 #ifdef BOOST_GEOMETRY_DEBUG_ENRICH
0500     std::cout << "ENRICH-assign Ring " << pair.first << std::endl;
0501 #endif
0502         if (is_dissolve)
0503         {
0504             detail::overlay::enrich_adapt(pair.second, turns);
0505         }
0506
0507         detail::overlay::enrich_assign(pair.second, turns, ! is_dissolve);
0508     }
0509
0510     if (has_cc)
0511     {
0512         detail::overlay::calculate_remaining_distance(turns);
0513     }
0514
0515 #ifdef BOOST_GEOMETRY_DEBUG_ENRICH
0516     //detail::overlay::check_graph(turns, for_operation);
0517 #endif
0518
0519 }
0520
0521 }} // namespace boost::geometry
0522
0523 #endif // BOOST_GEOMETRY_ALGORITHMS_DETAIL_OVERLAY_ENRICH_HPP