]>
Commit | Line | Data |
---|---|---|
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 Adam Wulkiewicz, Lodz, Poland. | |
7 | ||
8 | // This file was modified by Oracle on 2014. | |
9 | // Modifications copyright (c) 2014 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_POINT_ON_SURFACE_HPP | |
18 | #define BOOST_GEOMETRY_ALGORITHMS_POINT_ON_SURFACE_HPP | |
19 | ||
20 | ||
21 | #include <cstddef> | |
22 | ||
23 | #include <numeric> | |
24 | ||
25 | #include <boost/concept_check.hpp> | |
26 | #include <boost/range.hpp> | |
27 | ||
28 | #include <boost/geometry/core/point_type.hpp> | |
29 | #include <boost/geometry/core/ring_type.hpp> | |
30 | ||
31 | #include <boost/geometry/geometries/concepts/check.hpp> | |
32 | ||
33 | #include <boost/geometry/algorithms/detail/extreme_points.hpp> | |
34 | ||
35 | #include <boost/geometry/strategies/cartesian/centroid_bashein_detmer.hpp> | |
36 | ||
37 | ||
38 | namespace boost { namespace geometry | |
39 | { | |
40 | ||
41 | ||
42 | #ifndef DOXYGEN_NO_DETAIL | |
43 | namespace detail { namespace point_on_surface | |
44 | { | |
45 | ||
46 | template <typename CoordinateType, int Dimension> | |
47 | struct specific_coordinate_first | |
48 | { | |
49 | CoordinateType const m_value_to_be_first; | |
50 | ||
51 | inline specific_coordinate_first(CoordinateType value_to_be_skipped) | |
52 | : m_value_to_be_first(value_to_be_skipped) | |
53 | {} | |
54 | ||
55 | template <typename Point> | |
56 | inline bool operator()(Point const& lhs, Point const& rhs) | |
57 | { | |
58 | CoordinateType const lh = geometry::get<Dimension>(lhs); | |
59 | CoordinateType const rh = geometry::get<Dimension>(rhs); | |
60 | ||
61 | // If both lhs and rhs equal m_value_to_be_first, | |
62 | // we should handle conform if lh < rh = FALSE | |
63 | // The first condition meets that, keep it first | |
64 | if (geometry::math::equals(rh, m_value_to_be_first)) | |
65 | { | |
66 | // Handle conform lh < -INF, which is always false | |
67 | return false; | |
68 | } | |
69 | ||
70 | if (geometry::math::equals(lh, m_value_to_be_first)) | |
71 | { | |
72 | // Handle conform -INF < rh, which is always true | |
73 | return true; | |
74 | } | |
75 | ||
76 | return lh < rh; | |
77 | } | |
78 | }; | |
79 | ||
80 | template <int Dimension, typename Collection, typename Value, typename Predicate> | |
81 | inline bool max_value(Collection const& collection, Value& the_max, Predicate const& predicate) | |
82 | { | |
83 | bool first = true; | |
84 | for (typename Collection::const_iterator it = collection.begin(); it != collection.end(); ++it) | |
85 | { | |
86 | if (! it->empty()) | |
87 | { | |
88 | Value the_value = geometry::get<Dimension>(*std::max_element(it->begin(), it->end(), predicate)); | |
89 | if (first || the_value > the_max) | |
90 | { | |
91 | the_max = the_value; | |
92 | first = false; | |
93 | } | |
94 | } | |
95 | } | |
96 | return ! first; | |
97 | } | |
98 | ||
99 | ||
100 | template <int Dimension, typename Value> | |
101 | struct select_below | |
102 | { | |
103 | Value m_value; | |
104 | inline select_below(Value const& v) | |
105 | : m_value(v) | |
106 | {} | |
107 | ||
108 | template <typename Intruder> | |
109 | inline bool operator()(Intruder const& intruder) const | |
110 | { | |
111 | if (intruder.empty()) | |
112 | { | |
113 | return true; | |
114 | } | |
115 | Value max = geometry::get<Dimension>(*std::max_element(intruder.begin(), intruder.end(), detail::extreme_points::compare<Dimension>())); | |
116 | return geometry::math::equals(max, m_value) || max < m_value; | |
117 | } | |
118 | }; | |
119 | ||
120 | template <int Dimension, typename Value> | |
121 | struct adapt_base | |
122 | { | |
123 | Value m_value; | |
124 | inline adapt_base(Value const& v) | |
125 | : m_value(v) | |
126 | {} | |
127 | ||
128 | template <typename Intruder> | |
129 | inline void operator()(Intruder& intruder) const | |
130 | { | |
131 | if (intruder.size() >= 3) | |
132 | { | |
133 | detail::extreme_points::move_along_vector<Dimension>(intruder, m_value); | |
134 | } | |
135 | } | |
136 | }; | |
137 | ||
138 | template <int Dimension, typename Value> | |
139 | struct min_of_intruder | |
140 | { | |
141 | template <typename Intruder> | |
142 | inline bool operator()(Intruder const& lhs, Intruder const& rhs) const | |
143 | { | |
144 | Value lhs_min = geometry::get<Dimension>(*std::min_element(lhs.begin(), lhs.end(), detail::extreme_points::compare<Dimension>())); | |
145 | Value rhs_min = geometry::get<Dimension>(*std::min_element(rhs.begin(), rhs.end(), detail::extreme_points::compare<Dimension>())); | |
146 | return lhs_min < rhs_min; | |
147 | } | |
148 | }; | |
149 | ||
150 | ||
151 | template <typename Point, typename P> | |
152 | inline void calculate_average(Point& point, std::vector<P> const& points) | |
153 | { | |
154 | typedef typename geometry::coordinate_type<Point>::type coordinate_type; | |
155 | typedef typename std::vector<P>::const_iterator iterator_type; | |
156 | typedef typename std::vector<P>::size_type size_type; | |
157 | ||
158 | coordinate_type x = 0; | |
159 | coordinate_type y = 0; | |
160 | ||
161 | iterator_type end = points.end(); | |
162 | for ( iterator_type it = points.begin() ; it != end ; ++it) | |
163 | { | |
164 | x += geometry::get<0>(*it); | |
165 | y += geometry::get<1>(*it); | |
166 | } | |
167 | ||
168 | size_type const count = points.size(); | |
169 | geometry::set<0>(point, x / count); | |
170 | geometry::set<1>(point, y / count); | |
171 | } | |
172 | ||
173 | ||
174 | template <int Dimension, typename Extremes, typename Intruders, typename CoordinateType> | |
175 | inline void replace_extremes_for_self_tangencies(Extremes& extremes, Intruders& intruders, CoordinateType const& max_intruder) | |
176 | { | |
177 | // This function handles self-tangencies. | |
178 | // Self-tangencies use, as usual, the major part of code... | |
179 | ||
180 | // ___ e | |
181 | // /|\ \ . | |
182 | // / | \ \ . | |
183 | // / | \ \ . | |
184 | // / | \ \ . | |
185 | // / /\ | \ \ . | |
186 | // i2 i1 | |
187 | ||
188 | // The picture above shows the extreme (outside, "e") and two intruders ("i1","i2") | |
189 | // Assume that "i1" is self-tangent with the extreme, in one point at the top | |
190 | // Now the "penultimate" value is searched, this is is the top of i2 | |
191 | // Then everything including and below (this is "i2" here) is removed | |
192 | // Then the base of "i1" and of "e" is adapted to this penultimate value | |
193 | // It then looks like: | |
194 | ||
195 | // b ___ e | |
196 | // /|\ \ . | |
197 | // / | \ \ . | |
198 | // / | \ \ . | |
199 | // / | \ \ . | |
200 | // a c i1 | |
201 | ||
202 | // Then intruders (here "i1" but there may be more) are sorted from left to right | |
203 | // Finally points "a","b" and "c" (in this order) are selected as a new triangle. | |
204 | // This triangle will have a centroid which is inside (assumed that intruders left segment | |
205 | // is not equal to extremes left segment, but that polygon would be invalid) | |
206 | ||
207 | // Find highest non-self tangent intrusion, if any | |
208 | CoordinateType penultimate_value; | |
209 | specific_coordinate_first<CoordinateType, Dimension> pu_compare(max_intruder); | |
210 | if (max_value<Dimension>(intruders, penultimate_value, pu_compare)) | |
211 | { | |
212 | // Throw away all intrusions <= this value, and of the kept one set this as base. | |
213 | select_below<Dimension, CoordinateType> predicate(penultimate_value); | |
214 | intruders.erase | |
215 | ( | |
216 | std::remove_if(boost::begin(intruders), boost::end(intruders), predicate), | |
217 | boost::end(intruders) | |
218 | ); | |
219 | adapt_base<Dimension, CoordinateType> fe_predicate(penultimate_value); | |
220 | // Sort from left to right (or bottom to top if Dimension=0) | |
221 | std::for_each(boost::begin(intruders), boost::end(intruders), fe_predicate); | |
222 | ||
223 | // Also adapt base of extremes | |
224 | detail::extreme_points::move_along_vector<Dimension>(extremes, penultimate_value); | |
225 | } | |
226 | // Then sort in 1-Dim. Take first to calc centroid. | |
227 | std::sort(boost::begin(intruders), boost::end(intruders), min_of_intruder<1 - Dimension, CoordinateType>()); | |
228 | ||
229 | Extremes triangle; | |
230 | triangle.reserve(3); | |
231 | ||
232 | // Make a triangle of first two points of extremes (the ramp, from left to right), and last point of first intruder (which goes from right to left) | |
233 | std::copy(extremes.begin(), extremes.begin() + 2, std::back_inserter(triangle)); | |
234 | triangle.push_back(intruders.front().back()); | |
235 | ||
236 | // (alternatively we could use the last two points of extremes, and first point of last intruder...): | |
237 | //// ALTERNATIVE: std::copy(extremes.rbegin(), extremes.rbegin() + 2, std::back_inserter(triangle)); | |
238 | //// ALTERNATIVE: triangle.push_back(intruders.back().front()); | |
239 | ||
240 | // Now replace extremes with this smaller subset, a triangle, such that centroid calculation will result in a point inside | |
241 | extremes = triangle; | |
242 | } | |
243 | ||
244 | template <int Dimension, typename Geometry, typename Point> | |
245 | inline bool calculate_point_on_surface(Geometry const& geometry, Point& point) | |
246 | { | |
247 | typedef typename geometry::point_type<Geometry>::type point_type; | |
248 | typedef typename geometry::coordinate_type<Geometry>::type coordinate_type; | |
249 | std::vector<point_type> extremes; | |
250 | ||
251 | typedef std::vector<std::vector<point_type> > intruders_type; | |
252 | intruders_type intruders; | |
253 | geometry::extreme_points<Dimension>(geometry, extremes, intruders); | |
254 | ||
255 | if (extremes.size() < 3) | |
256 | { | |
257 | return false; | |
258 | } | |
259 | ||
260 | // If there are intruders, find the max. | |
261 | if (! intruders.empty()) | |
262 | { | |
263 | coordinate_type max_intruder; | |
264 | detail::extreme_points::compare<Dimension> compare; | |
265 | if (max_value<Dimension>(intruders, max_intruder, compare)) | |
266 | { | |
267 | coordinate_type max_extreme = geometry::get<Dimension>(*std::max_element(extremes.begin(), extremes.end(), detail::extreme_points::compare<Dimension>())); | |
268 | if (max_extreme > max_intruder) | |
269 | { | |
270 | detail::extreme_points::move_along_vector<Dimension>(extremes, max_intruder); | |
271 | } | |
272 | else | |
273 | { | |
274 | replace_extremes_for_self_tangencies<Dimension>(extremes, intruders, max_intruder); | |
275 | } | |
276 | } | |
277 | } | |
278 | ||
279 | // Now calculate the average/centroid of the (possibly adapted) extremes | |
280 | calculate_average(point, extremes); | |
281 | ||
282 | return true; | |
283 | } | |
284 | ||
285 | }} // namespace detail::point_on_surface | |
286 | #endif // DOXYGEN_NO_DETAIL | |
287 | ||
288 | ||
289 | /*! | |
290 | \brief Assigns a Point guaranteed to lie on the surface of the Geometry | |
291 | \tparam Geometry geometry type. This also defines the type of the output point | |
292 | \param geometry Geometry to take point from | |
293 | \param point Point to assign | |
294 | */ | |
295 | template <typename Geometry, typename Point> | |
296 | inline void point_on_surface(Geometry const& geometry, Point & point) | |
297 | { | |
298 | concepts::check<Point>(); | |
299 | concepts::check<Geometry const>(); | |
300 | ||
301 | // First try in Y-direction (which should always succeed for valid polygons) | |
302 | if (! detail::point_on_surface::calculate_point_on_surface<1>(geometry, point)) | |
303 | { | |
304 | // For invalid polygons, we might try X-direction | |
305 | detail::point_on_surface::calculate_point_on_surface<0>(geometry, point); | |
306 | } | |
307 | } | |
308 | ||
309 | /*! | |
310 | \brief Returns point guaranteed to lie on the surface of the Geometry | |
311 | \tparam Geometry geometry type. This also defines the type of the output point | |
312 | \param geometry Geometry to take point from | |
313 | \return The Point guaranteed to lie on the surface of the Geometry | |
314 | */ | |
315 | template<typename Geometry> | |
316 | inline typename geometry::point_type<Geometry>::type | |
317 | return_point_on_surface(Geometry const& geometry) | |
318 | { | |
319 | typename geometry::point_type<Geometry>::type result; | |
320 | geometry::point_on_surface(geometry, result); | |
321 | return result; | |
322 | } | |
323 | ||
324 | }} // namespace boost::geometry | |
325 | ||
326 | ||
327 | #endif // BOOST_GEOMETRY_ALGORITHMS_POINT_ON_SURFACE_HPP |