detail/overlay/get_clusters.hpp

0001 // Boost.Geometry (aka GGL, Generic Geometry Library)
0002
0003 // Copyright (c) 2021 Barend Gehrels, Amsterdam, the Netherlands.
0004
0005 // This file was modified by Oracle on 2024.
0006 // Modifications copyright (c) 2024 Oracle and/or its affiliates.
0007 // Contributed and/or modified by Vissarion Fysikopoulos, on behalf of Oracle
0008
0009 // Use, modification and distribution is subject to the Boost Software License,
0010 // Version 1.0. (See accompanying file LICENSE_1_0.txt or copy at
0011 // http://www.boost.org/LICENSE_1_0.txt)
0012
0013 #ifndef BOOST_GEOMETRY_ALGORITHMS_DETAIL_OVERLAY_GET_CLUSTERS_HPP
0014 #define BOOST_GEOMETRY_ALGORITHMS_DETAIL_OVERLAY_GET_CLUSTERS_HPP
0015
0016 #include <algorithm>
0017 #include <map>
0018 #include <vector>
0019
0020 #include <boost/geometry/core/access.hpp>
0021 #include <boost/geometry/algorithms/detail/overlay/approximately_equals.hpp>
0022 #include <boost/geometry/algorithms/detail/overlay/cluster_info.hpp>
0023 #include <boost/geometry/algorithms/detail/overlay/get_ring.hpp>
0024 #include <boost/range/value_type.hpp>
0025 #include <boost/geometry/util/math.hpp>
0026
0027 namespace boost { namespace geometry
0028 {
0029
0030 #ifndef DOXYGEN_NO_DETAIL
0031 namespace detail { namespace overlay
0032 {
0033
0034 template <bool Integral = false>
0035 struct sweep_equal_policy
0036 {
0037
0038 public:
0039     // Returns true if point are considered equal
0040     template <typename P>
0041     static inline bool equals(P const& p1, P const& p2)
0042     {
0043         using coor_t = coordinate_type_t<P>;
0044         static auto const tolerance
0045             = common_approximately_equals_epsilon_multiplier<coor_t>::value();
0046         return approximately_equals(p1, p2, tolerance);
0047     }
0048
0049     template <typename T>
0050     static inline bool exceeds(T value)
0051     {
0052         static auto const tolerance
0053             = common_approximately_equals_epsilon_multiplier<T>::value();
0054         T const limit = T(1) / tolerance;
0055         return value > limit;
0056     }
0057 };
0058
0059 template <>
0060 struct sweep_equal_policy<true>
0061 {
0062     template <typename P>
0063     static inline bool equals(P const& p1, P const& p2)
0064     {
0065         return geometry::get<0>(p1) == geometry::get<0>(p2)
0066             && geometry::get<1>(p1) == geometry::get<1>(p2);
0067     }
0068
0069     template <typename T>
0070     static inline bool exceeds(T value)
0071     {
0072         return value > 0;
0073     }
0074 };
0075
0076 template <typename Point>
0077 struct turn_with_point
0078 {
0079     std::size_t turn_index;
0080     Point pnt;
0081 };
0082
0083 template <typename Cluster, typename Point>
0084 struct cluster_with_point
0085 {
0086     Cluster cluster;
0087     Point pnt;
0088 };
0089
0090 // Use a sweep algorithm to detect clusters
0091 template
0092 <
0093     typename Turns,
0094     typename Clusters
0095 >
0096 inline void get_clusters(Turns& turns, Clusters& clusters)
0097 {
0098     using turn_type = typename boost::range_value<Turns>::type;
0099     using cluster_type = typename Clusters::mapped_type;
0100     using point_type = typename turn_type::point_type;
0101
0102     sweep_equal_policy
0103         <
0104             std::is_integral<coordinate_type_t<point_type>>::value
0105         > equal_policy;
0106
0107     std::vector<turn_with_point<point_type>> points;
0108     std::size_t turn_index = 0;
0109     for (auto const& turn : turns)
0110     {
0111         if (! turn.discarded)
0112         {
0113             points.push_back({turn_index, turn.point});
0114         }
0115         turn_index++;
0116     }
0117
0118     // Sort the points from top to bottom
0119     std::sort(points.begin(), points.end(), [](auto const& e1, auto const& e2)
0120     {
0121        return geometry::get<1>(e1.pnt) > geometry::get<1>(e2.pnt);
0122     });
0123
0124     // The output vector will be sorted from bottom too
0125     std::vector<cluster_with_point<cluster_type, point_type>> clustered_points;
0126
0127     // Compare points with each other. Performance is O(n log(n)) because of the sorting.
0128     for (auto it1 = points.begin(); it1 != points.end(); ++it1)
0129     {
0130         // Inner loop, iterates until it exceeds coordinates in y-direction
0131         for (auto it2 = it1 + 1; it2 != points.end(); ++it2)
0132         {
0133             auto const d = geometry::get<1>(it1->pnt) - geometry::get<1>(it2->pnt);
0134             if (equal_policy.exceeds(d))
0135             {
0136                 // Points at this y-coordinate or below cannot be equal
0137                 break;
0138             }
0139             if (equal_policy.equals(it1->pnt, it2->pnt))
0140             {
0141                 std::size_t cindex = 0;
0142
0143                 // Most recent clusters (with this y-value) are at the bottom
0144                 // therefore we can stop as soon as the y-value is out of reach (TODO)
0145                 bool found = false;
0146                 for (auto cit = clustered_points.begin();
0147                      cit != clustered_points.end(); ++cit)
0148                 {
0149                     found = equal_policy.equals(cit->pnt, it1->pnt);
0150                     if (found)
0151                     {
0152                         break;
0153                     }
0154                     cindex++;
0155                 }
0156
0157                 // Add new cluster
0158                 if (! found)
0159                 {
0160                    cindex = clustered_points.size();
0161                    cluster_type newcluster;
0162                    clustered_points.push_back({newcluster, it1->pnt});
0163                 }
0164                 clustered_points[cindex].cluster.turn_indices.insert(it1->turn_index);
0165                 clustered_points[cindex].cluster.turn_indices.insert(it2->turn_index);
0166             }
0167         }
0168     }
0169
0170     // Convert to map
0171     signed_size_type cluster_id = 1;
0172     for (auto& trace : clustered_points)
0173     {
0174         clusters[cluster_id++] = trace.cluster;
0175     }
0176 }
0177
0178 }} // namespace detail::overlay
0179 #endif //DOXYGEN_NO_DETAIL
0180
0181
0182 }} // namespace boost::geometry
0183
0184 #endif // BOOST_GEOMETRY_ALGORITHMS_DETAIL_OVERLAY_GET_CLUSTERS_HPP