1 // Boost.Geometry (aka GGL, Generic Geometry Library)
2 
3 // Copyright (c) 2007-2012 Barend Gehrels, Amsterdam, the Netherlands.
4 
5 // This file was modified by Oracle on 2014.
6 // Modifications copyright (c) 2014 Oracle and/or its affiliates.
7 
8 // Contributed and/or modified by Adam Wulkiewicz, on behalf of Oracle
9 
10 // Parts of Boost.Geometry are redesigned from Geodan's Geographic Library
11 // (geolib/GGL), copyright (c) 1995-2010 Geodan, Amsterdam, the Netherlands.
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_STRATEGIES_AGNOSTIC_CONVEX_GRAHAM_ANDREW_HPP
18 #define BOOST_GEOMETRY_STRATEGIES_AGNOSTIC_CONVEX_GRAHAM_ANDREW_HPP
19 
20 
21 #include <cstddef>
22 #include <algorithm>
23 #include <vector>
24 
25 #include <boost/range.hpp>
26 
27 #include <boost/geometry/core/assert.hpp>
28 #include <boost/geometry/core/cs.hpp>
29 #include <boost/geometry/core/point_type.hpp>
30 #include <boost/geometry/strategies/convex_hull.hpp>
31 
32 #include <boost/geometry/views/detail/range_type.hpp>
33 
34 #include <boost/geometry/policies/compare.hpp>
35 
36 #include <boost/geometry/algorithms/detail/for_each_range.hpp>
37 #include <boost/geometry/views/reversible_view.hpp>
38 
39 
40 namespace boost { namespace geometry
41 {
42 
43 namespace strategy { namespace convex_hull
44 {
45 
46 #ifndef DOXYGEN_NO_DETAIL
47 namespace detail
48 {
49 
50 
51 template
52 <
53     typename InputRange,
54     typename RangeIterator,
55     typename StrategyLess,
56     typename StrategyGreater
57 >
58 struct get_extremes
59 {
60     typedef typename point_type<InputRange>::type point_type;
61 
62     point_type left, right;
63 
64     bool first;
65 
66     StrategyLess less;
67     StrategyGreater greater;
68 
get_extremesboost::geometry::strategy::convex_hull::detail::get_extremes69     inline get_extremes()
70         : first(true)
71     {}
72 
applyboost::geometry::strategy::convex_hull::detail::get_extremes73     inline void apply(InputRange const& range)
74     {
75         if (boost::size(range) == 0)
76         {
77             return;
78         }
79 
80         // First iterate through this range
81         // (this two-stage approach avoids many point copies,
82         //  because iterators are kept in memory. Because iterators are
83         //  not persistent (in MSVC) this approach is not applicable
84         //  for more ranges together)
85 
86         RangeIterator left_it = boost::begin(range);
87         RangeIterator right_it = boost::begin(range);
88 
89         for (RangeIterator it = boost::begin(range) + 1;
90             it != boost::end(range);
91             ++it)
92         {
93             if (less(*it, *left_it))
94             {
95                 left_it = it;
96             }
97 
98             if (greater(*it, *right_it))
99             {
100                 right_it = it;
101             }
102         }
103 
104         // Then compare with earlier
105         if (first)
106         {
107             // First time, assign left/right
108             left = *left_it;
109             right = *right_it;
110             first = false;
111         }
112         else
113         {
114             // Next time, check if this range was left/right from
115             // the extremes already collected
116             if (less(*left_it, left))
117             {
118                 left = *left_it;
119             }
120 
121             if (greater(*right_it, right))
122             {
123                 right = *right_it;
124             }
125         }
126     }
127 };
128 
129 
130 template
131 <
132     typename InputRange,
133     typename RangeIterator,
134     typename Container,
135     typename SideStrategy
136 >
137 struct assign_range
138 {
139     Container lower_points, upper_points;
140 
141     typedef typename point_type<InputRange>::type point_type;
142 
143     point_type const& most_left;
144     point_type const& most_right;
145 
assign_rangeboost::geometry::strategy::convex_hull::detail::assign_range146     inline assign_range(point_type const& left, point_type const& right)
147         : most_left(left)
148         , most_right(right)
149     {}
150 
applyboost::geometry::strategy::convex_hull::detail::assign_range151     inline void apply(InputRange const& range)
152     {
153         typedef SideStrategy side;
154 
155         // Put points in one of the two output sequences
156         for (RangeIterator it = boost::begin(range);
157             it != boost::end(range);
158             ++it)
159         {
160             // check if it is lying most_left or most_right from the line
161 
162             int dir = side::apply(most_left, most_right, *it);
163             switch(dir)
164             {
165                 case 1 : // left side
166                     upper_points.push_back(*it);
167                     break;
168                 case -1 : // right side
169                     lower_points.push_back(*it);
170                     break;
171 
172                 // 0: on line most_left-most_right,
173                 //    or most_left, or most_right,
174                 //    -> all never part of hull
175             }
176         }
177     }
178 };
179 
180 template <typename Range>
sort(Range & range)181 static inline void sort(Range& range)
182 {
183     typedef typename boost::range_value<Range>::type point_type;
184     typedef geometry::less<point_type> comparator;
185 
186     std::sort(boost::begin(range), boost::end(range), comparator());
187 }
188 
189 } // namespace detail
190 #endif // DOXYGEN_NO_DETAIL
191 
192 
193 /*!
194 \brief Graham scan strategy to calculate convex hull
195 \ingroup strategies
196 \note Completely reworked version inspired on the sources listed below
197 \see http://www.ddj.com/architect/201806315
198 \see http://marknelson.us/2007/08/22/convex
199  */
200 template <typename InputGeometry, typename OutputPoint>
201 class graham_andrew
202 {
203 public :
204     typedef OutputPoint point_type;
205     typedef InputGeometry geometry_type;
206 
207 private:
208 
209     typedef typename cs_tag<point_type>::type cs_tag;
210 
211     typedef typename std::vector<point_type> container_type;
212     typedef typename std::vector<point_type>::const_iterator iterator;
213     typedef typename std::vector<point_type>::const_reverse_iterator rev_iterator;
214 
215 
216     class partitions
217     {
218         friend class graham_andrew;
219 
220         container_type m_lower_hull;
221         container_type m_upper_hull;
222         container_type m_copied_input;
223     };
224 
225 
226 public:
227     typedef partitions state_type;
228 
229 
apply(InputGeometry const & geometry,partitions & state) const230     inline void apply(InputGeometry const& geometry, partitions& state) const
231     {
232         // First pass.
233         // Get min/max (in most cases left / right) points
234         // This makes use of the geometry::less/greater predicates
235 
236         // For the left boundary it is important that multiple points
237         // are sorted from bottom to top. Therefore the less predicate
238         // does not take the x-only template parameter (this fixes ticket #6019.
239         // For the right boundary it is not necessary (though also not harmful),
240         // because points are sorted from bottom to top in a later stage.
241         // For symmetry and to get often more balanced lower/upper halves
242         // we keep it.
243 
244         typedef typename geometry::detail::range_type<InputGeometry>::type range_type;
245 
246         typedef typename boost::range_iterator
247             <
248                 range_type const
249             >::type range_iterator;
250 
251         detail::get_extremes
252             <
253                 range_type,
254                 range_iterator,
255                 geometry::less<point_type>,
256                 geometry::greater<point_type>
257             > extremes;
258         geometry::detail::for_each_range(geometry, extremes);
259 
260         // Bounding left/right points
261         // Second pass, now that extremes are found, assign all points
262         // in either lower, either upper
263         detail::assign_range
264             <
265                 range_type,
266                 range_iterator,
267                 container_type,
268                 typename strategy::side::services::default_strategy<cs_tag>::type
269             > assigner(extremes.left, extremes.right);
270 
271         geometry::detail::for_each_range(geometry, assigner);
272 
273 
274         // Sort both collections, first on x(, then on y)
275         detail::sort(assigner.lower_points);
276         detail::sort(assigner.upper_points);
277 
278         //std::cout << boost::size(assigner.lower_points) << std::endl;
279         //std::cout << boost::size(assigner.upper_points) << std::endl;
280 
281         // And decide which point should be in the final hull
282         build_half_hull<-1>(assigner.lower_points, state.m_lower_hull,
283                 extremes.left, extremes.right);
284         build_half_hull<1>(assigner.upper_points, state.m_upper_hull,
285                 extremes.left, extremes.right);
286     }
287 
288 
289     template <typename OutputIterator>
result(partitions const & state,OutputIterator out,bool clockwise,bool closed) const290     inline void result(partitions const& state,
291                        OutputIterator out,
292                        bool clockwise,
293                        bool closed) const
294     {
295         if (clockwise)
296         {
297             output_ranges(state.m_upper_hull, state.m_lower_hull, out, closed);
298         }
299         else
300         {
301             output_ranges(state.m_lower_hull, state.m_upper_hull, out, closed);
302         }
303     }
304 
305 
306 private:
307 
308     template <int Factor>
build_half_hull(container_type const & input,container_type & output,point_type const & left,point_type const & right)309     static inline void build_half_hull(container_type const& input,
310             container_type& output,
311             point_type const& left, point_type const& right)
312     {
313         output.push_back(left);
314         for(iterator it = input.begin(); it != input.end(); ++it)
315         {
316             add_to_hull<Factor>(*it, output);
317         }
318         add_to_hull<Factor>(right, output);
319     }
320 
321 
322     template <int Factor>
add_to_hull(point_type const & p,container_type & output)323     static inline void add_to_hull(point_type const& p, container_type& output)
324     {
325         typedef typename strategy::side::services::default_strategy<cs_tag>::type side;
326 
327         output.push_back(p);
328         std::size_t output_size = output.size();
329         while (output_size >= 3)
330         {
331             rev_iterator rit = output.rbegin();
332             point_type const last = *rit++;
333             point_type const& last2 = *rit++;
334 
335             if (Factor * side::apply(*rit, last, last2) <= 0)
336             {
337                 // Remove last two points from stack, and add last again
338                 // This is much faster then erasing the one but last.
339                 output.pop_back();
340                 output.pop_back();
341                 output.push_back(last);
342                 output_size--;
343             }
344             else
345             {
346                 return;
347             }
348         }
349     }
350 
351 
352     template <typename OutputIterator>
output_ranges(container_type const & first,container_type const & second,OutputIterator out,bool closed)353     static inline void output_ranges(container_type const& first, container_type const& second,
354                                      OutputIterator out, bool closed)
355     {
356         std::copy(boost::begin(first), boost::end(first), out);
357 
358         BOOST_GEOMETRY_ASSERT(closed ? !boost::empty(second) : boost::size(second) > 1);
359         std::copy(++boost::rbegin(second), // skip the first Point
360                   closed ? boost::rend(second) : --boost::rend(second), // skip the last Point if open
361                   out);
362 
363         typedef typename boost::range_size<container_type>::type size_type;
364         size_type const count = boost::size(first) + boost::size(second) - 1;
365         // count describes a closed case but comparison with min size of closed
366         // gives the result compatible also with open
367         // here core_detail::closure::minimum_ring_size<closed> could be used
368         if (count < 4)
369         {
370             // there should be only one missing
371             *out++ = *boost::begin(first);
372         }
373     }
374 };
375 
376 }} // namespace strategy::convex_hull
377 
378 
379 #ifndef DOXYGEN_NO_STRATEGY_SPECIALIZATIONS
380 template <typename InputGeometry, typename OutputPoint>
381 struct strategy_convex_hull<InputGeometry, OutputPoint, cartesian_tag>
382 {
383     typedef strategy::convex_hull::graham_andrew<InputGeometry, OutputPoint> type;
384 };
385 #endif
386 
387 }} // namespace boost::geometry
388 
389 
390 #endif // BOOST_GEOMETRY_STRATEGIES_AGNOSTIC_CONVEX_GRAHAM_ANDREW_HPP
391