1 // Boost.Geometry Index
2 //
3 // R-tree linear split algorithm implementation
4 //
5 // Copyright (c) 2008 Federico J. Fernandez.
6 // Copyright (c) 2011-2014 Adam Wulkiewicz, Lodz, Poland.
7 //
8 // Use, modification and distribution is subject to the Boost Software License,
9 // Version 1.0. (See accompanying file LICENSE_1_0.txt or copy at
10 // http://www.boost.org/LICENSE_1_0.txt)
11 
12 #ifndef BOOST_GEOMETRY_INDEX_DETAIL_RTREE_LINEAR_REDISTRIBUTE_ELEMENTS_HPP
13 #define BOOST_GEOMETRY_INDEX_DETAIL_RTREE_LINEAR_REDISTRIBUTE_ELEMENTS_HPP
14 
15 #include <boost/type_traits/is_unsigned.hpp>
16 
17 #include <boost/geometry/index/detail/algorithms/content.hpp>
18 #include <boost/geometry/index/detail/bounded_view.hpp>
19 
20 #include <boost/geometry/index/detail/rtree/node/node.hpp>
21 #include <boost/geometry/index/detail/rtree/visitors/insert.hpp>
22 #include <boost/geometry/index/detail/rtree/visitors/is_leaf.hpp>
23 
24 namespace boost { namespace geometry { namespace index {
25 
26 namespace detail { namespace rtree {
27 
28 namespace linear {
29 
30 template <typename R, typename T>
difference_dispatch(T const & from,T const & to,::boost::mpl::bool_<false> const &)31 inline R difference_dispatch(T const& from, T const& to, ::boost::mpl::bool_<false> const& /*is_unsigned*/)
32 {
33     return to - from;
34 }
35 
36 template <typename R, typename T>
difference_dispatch(T const & from,T const & to,::boost::mpl::bool_<true> const &)37 inline R difference_dispatch(T const& from, T const& to, ::boost::mpl::bool_<true> const& /*is_unsigned*/)
38 {
39     return from <= to ? R(to - from) : -R(from - to);
40 }
41 
42 template <typename R, typename T>
difference(T const & from,T const & to)43 inline R difference(T const& from, T const& to)
44 {
45     BOOST_MPL_ASSERT_MSG(!boost::is_unsigned<R>::value, RESULT_CANT_BE_UNSIGNED, (R));
46 
47     typedef ::boost::mpl::bool_<
48         boost::is_unsigned<T>::value
49     > is_unsigned;
50 
51     return difference_dispatch<R>(from, to, is_unsigned());
52 }
53 
54 // TODO: awulkiew
55 // In general, all aerial Indexables in the tree with box-like nodes will be analyzed as boxes
56 // because they must fit into larger box. Therefore the algorithm could be the same for Bounds type.
57 // E.g. if Bounds type is sphere, Indexables probably should be analyzed as spheres.
58 // 1. View could be provided to 'see' all Indexables as Bounds type.
59 //    Not ok in the case of big types like Ring, however it's possible that Rings won't be supported,
60 //    only simple types. Even then if we consider storing Box inside the Sphere we must calculate
61 //    the bounding sphere 2x for each box because there are 2 loops. For each calculation this means
62 //    4-2d or 8-3d expansions or -, / and sqrt().
63 // 2. Additional container could be used and reused if the Indexable type is other than the Bounds type.
64 
65 // IMPORTANT!
66 // Still probably the best way would be providing specialized algorithms for each Indexable-Bounds pair!
67 // Probably on pick_seeds algorithm level - For Bounds=Sphere seeds would be choosen differently
68 
69 // TODO: awulkiew
70 // there are loops inside find_greatest_normalized_separation::apply()
71 // iteration is done for each DimensionIndex.
72 // Separations and seeds for all DimensionIndex(es) could be calculated at once, stored, then the greatest would be choosen.
73 
74 // The following struct/method was adapted for the preliminary version of the R-tree. Then it was called:
75 // void find_normalized_separations(std::vector<Box> const& boxes, T& separation, unsigned int& first, unsigned int& second) const
76 
77 template <typename Elements, typename Parameters, typename Translator, typename Tag, size_t DimensionIndex>
78 struct find_greatest_normalized_separation
79 {
80     typedef typename Elements::value_type element_type;
81     typedef typename rtree::element_indexable_type<element_type, Translator>::type indexable_type;
82     typedef typename coordinate_type<indexable_type>::type coordinate_type;
83 
84     typedef typename boost::mpl::if_c<
85         boost::is_integral<coordinate_type>::value,
86         double,
87         coordinate_type
88     >::type separation_type;
89 
90     typedef typename geometry::point_type<indexable_type>::type point_type;
91     typedef geometry::model::box<point_type> bounds_type;
92     typedef index::detail::bounded_view<indexable_type, bounds_type> bounded_view_type;
93 
applyboost::geometry::index::detail::rtree::linear::find_greatest_normalized_separation94     static inline void apply(Elements const& elements,
95                              Parameters const& parameters,
96                              Translator const& translator,
97                              separation_type & separation,
98                              size_t & seed1,
99                              size_t & seed2)
100     {
101         const size_t elements_count = parameters.get_max_elements() + 1;
102         BOOST_GEOMETRY_INDEX_ASSERT(elements.size() == elements_count, "unexpected number of elements");
103         BOOST_GEOMETRY_INDEX_ASSERT(2 <= elements_count, "unexpected number of elements");
104 
105         // find the lowest low, highest high
106         bounded_view_type bounded_indexable_0(rtree::element_indexable(elements[0], translator));
107         coordinate_type lowest_low = geometry::get<min_corner, DimensionIndex>(bounded_indexable_0);
108         coordinate_type highest_high = geometry::get<max_corner, DimensionIndex>(bounded_indexable_0);
109 
110         // and the lowest high
111         coordinate_type lowest_high = highest_high;
112         size_t lowest_high_index = 0;
113         for ( size_t i = 1 ; i < elements_count ; ++i )
114         {
115             bounded_view_type bounded_indexable(rtree::element_indexable(elements[i], translator));
116             coordinate_type min_coord = geometry::get<min_corner, DimensionIndex>(bounded_indexable);
117             coordinate_type max_coord = geometry::get<max_corner, DimensionIndex>(bounded_indexable);
118 
119             if ( max_coord < lowest_high )
120             {
121                 lowest_high = max_coord;
122                 lowest_high_index = i;
123             }
124 
125             if ( min_coord < lowest_low )
126                 lowest_low = min_coord;
127 
128             if ( highest_high < max_coord )
129                 highest_high = max_coord;
130         }
131 
132         // find the highest low
133         size_t highest_low_index = lowest_high_index == 0 ? 1 : 0;
134         bounded_view_type bounded_indexable_hl(rtree::element_indexable(elements[highest_low_index], translator));
135         coordinate_type highest_low = geometry::get<min_corner, DimensionIndex>(bounded_indexable_hl);
136         for ( size_t i = highest_low_index ; i < elements_count ; ++i )
137         {
138             bounded_view_type bounded_indexable(rtree::element_indexable(elements[i], translator));
139             coordinate_type min_coord = geometry::get<min_corner, DimensionIndex>(bounded_indexable);
140             if ( highest_low < min_coord &&
141                  i != lowest_high_index )
142             {
143                 highest_low = min_coord;
144                 highest_low_index = i;
145             }
146         }
147 
148         coordinate_type const width = highest_high - lowest_low;
149 
150         // highest_low - lowest_high
151         separation = difference<separation_type>(lowest_high, highest_low);
152         // BOOST_GEOMETRY_INDEX_ASSERT(0 <= width);
153         if ( std::numeric_limits<coordinate_type>::epsilon() < width )
154             separation /= width;
155 
156         seed1 = highest_low_index;
157         seed2 = lowest_high_index;
158 
159         ::boost::ignore_unused_variable_warning(parameters);
160     }
161 };
162 
163 // Version for points doesn't calculate normalized separation since it would always be equal to 1
164 // It returns two seeds most distant to each other, separation is equal to distance
165 template <typename Elements, typename Parameters, typename Translator, size_t DimensionIndex>
166 struct find_greatest_normalized_separation<Elements, Parameters, Translator, point_tag, DimensionIndex>
167 {
168     typedef typename Elements::value_type element_type;
169     typedef typename rtree::element_indexable_type<element_type, Translator>::type indexable_type;
170     typedef typename coordinate_type<indexable_type>::type coordinate_type;
171 
172     typedef coordinate_type separation_type;
173 
applyboost::geometry::index::detail::rtree::linear::find_greatest_normalized_separation174     static inline void apply(Elements const& elements,
175                              Parameters const& parameters,
176                              Translator const& translator,
177                              separation_type & separation,
178                              size_t & seed1,
179                              size_t & seed2)
180     {
181         const size_t elements_count = parameters.get_max_elements() + 1;
182         BOOST_GEOMETRY_INDEX_ASSERT(elements.size() == elements_count, "unexpected number of elements");
183         BOOST_GEOMETRY_INDEX_ASSERT(2 <= elements_count, "unexpected number of elements");
184 
185         // find the lowest low, highest high
186         coordinate_type lowest = geometry::get<DimensionIndex>(rtree::element_indexable(elements[0], translator));
187         coordinate_type highest = geometry::get<DimensionIndex>(rtree::element_indexable(elements[0], translator));
188         size_t lowest_index = 0;
189         size_t highest_index = 0;
190         for ( size_t i = 1 ; i < elements_count ; ++i )
191         {
192             coordinate_type coord = geometry::get<DimensionIndex>(rtree::element_indexable(elements[i], translator));
193 
194             if ( coord < lowest )
195             {
196                 lowest = coord;
197                 lowest_index = i;
198             }
199 
200             if ( highest < coord )
201             {
202                 highest = coord;
203                 highest_index = i;
204             }
205         }
206 
207         separation = highest - lowest;
208         seed1 = lowest_index;
209         seed2 = highest_index;
210 
211         if ( lowest_index == highest_index )
212             seed2 = (lowest_index + 1) % elements_count; // % is just in case since if this is true lowest_index is 0
213 
214         ::boost::ignore_unused_variable_warning(parameters);
215     }
216 };
217 
218 template <typename Elements, typename Parameters, typename Translator, size_t Dimension>
219 struct pick_seeds_impl
220 {
221     BOOST_STATIC_ASSERT(0 < Dimension);
222 
223     typedef typename Elements::value_type element_type;
224     typedef typename rtree::element_indexable_type<element_type, Translator>::type indexable_type;
225 
226     typedef find_greatest_normalized_separation<
227         Elements, Parameters, Translator,
228         typename tag<indexable_type>::type, Dimension - 1
229     > find_norm_sep;
230 
231     typedef typename find_norm_sep::separation_type separation_type;
232 
applyboost::geometry::index::detail::rtree::linear::pick_seeds_impl233     static inline void apply(Elements const& elements,
234                              Parameters const& parameters,
235                              Translator const& tr,
236                              separation_type & separation,
237                              size_t & seed1,
238                              size_t & seed2)
239     {
240         pick_seeds_impl<Elements, Parameters, Translator, Dimension - 1>::apply(elements, parameters, tr, separation, seed1, seed2);
241 
242         separation_type current_separation;
243         size_t s1, s2;
244         find_norm_sep::apply(elements, parameters, tr, current_separation, s1, s2);
245 
246         // in the old implementation different operator was used: <= (y axis prefered)
247         if ( separation < current_separation )
248         {
249             separation = current_separation;
250             seed1 = s1;
251             seed2 = s2;
252         }
253     }
254 };
255 
256 template <typename Elements, typename Parameters, typename Translator>
257 struct pick_seeds_impl<Elements, Parameters, Translator, 1>
258 {
259     typedef typename Elements::value_type element_type;
260     typedef typename rtree::element_indexable_type<element_type, Translator>::type indexable_type;
261     typedef typename coordinate_type<indexable_type>::type coordinate_type;
262 
263     typedef find_greatest_normalized_separation<
264         Elements, Parameters, Translator,
265         typename tag<indexable_type>::type, 0
266     > find_norm_sep;
267 
268     typedef typename find_norm_sep::separation_type separation_type;
269 
applyboost::geometry::index::detail::rtree::linear::pick_seeds_impl270     static inline void apply(Elements const& elements,
271                              Parameters const& parameters,
272                              Translator const& tr,
273                              separation_type & separation,
274                              size_t & seed1,
275                              size_t & seed2)
276     {
277         find_norm_sep::apply(elements, parameters, tr, separation, seed1, seed2);
278     }
279 };
280 
281 // from void linear_pick_seeds(node_pointer const& n, unsigned int &seed1, unsigned int &seed2) const
282 
283 template <typename Elements, typename Parameters, typename Translator>
pick_seeds(Elements const & elements,Parameters const & parameters,Translator const & tr,size_t & seed1,size_t & seed2)284 inline void pick_seeds(Elements const& elements,
285                        Parameters const& parameters,
286                        Translator const& tr,
287                        size_t & seed1,
288                        size_t & seed2)
289 {
290     typedef typename Elements::value_type element_type;
291     typedef typename rtree::element_indexable_type<element_type, Translator>::type indexable_type;
292 
293     typedef pick_seeds_impl
294         <
295             Elements, Parameters, Translator,
296             geometry::dimension<indexable_type>::value
297         > impl;
298     typedef typename impl::separation_type separation_type;
299 
300     separation_type separation = 0;
301     impl::apply(elements, parameters, tr, separation, seed1, seed2);
302 }
303 
304 } // namespace linear
305 
306 // from void split_node(node_pointer const& n, node_pointer& n1, node_pointer& n2) const
307 
308 template <typename Value, typename Options, typename Translator, typename Box, typename Allocators>
309 struct redistribute_elements<Value, Options, Translator, Box, Allocators, linear_tag>
310 {
311     typedef typename Options::parameters_type parameters_type;
312 
313     typedef typename rtree::node<Value, parameters_type, Box, Allocators, typename Options::node_tag>::type node;
314     typedef typename rtree::internal_node<Value, parameters_type, Box, Allocators, typename Options::node_tag>::type internal_node;
315     typedef typename rtree::leaf<Value, parameters_type, Box, Allocators, typename Options::node_tag>::type leaf;
316 
317     template <typename Node>
applyboost::geometry::index::detail::rtree::redistribute_elements318     static inline void apply(Node & n,
319                              Node & second_node,
320                              Box & box1,
321                              Box & box2,
322                              parameters_type const& parameters,
323                              Translator const& translator,
324                              Allocators & allocators)
325     {
326         typedef typename rtree::elements_type<Node>::type elements_type;
327         typedef typename elements_type::value_type element_type;
328         typedef typename rtree::element_indexable_type<element_type, Translator>::type indexable_type;
329         typedef typename index::detail::default_content_result<Box>::type content_type;
330 
331         elements_type & elements1 = rtree::elements(n);
332         elements_type & elements2 = rtree::elements(second_node);
333         const size_t elements1_count = parameters.get_max_elements() + 1;
334 
335         BOOST_GEOMETRY_INDEX_ASSERT(elements1.size() == elements1_count, "unexpected number of elements");
336 
337         // copy original elements - use in-memory storage (std::allocator)
338         // TODO: move if noexcept
339         typedef typename rtree::container_from_elements_type<elements_type, element_type>::type
340             container_type;
341         container_type elements_copy(elements1.begin(), elements1.end());                                   // MAY THROW, STRONG (alloc, copy)
342 
343         // calculate initial seeds
344         size_t seed1 = 0;
345         size_t seed2 = 0;
346         linear::pick_seeds(elements_copy, parameters, translator, seed1, seed2);
347 
348         // prepare nodes' elements containers
349         elements1.clear();
350         BOOST_GEOMETRY_INDEX_ASSERT(elements2.empty(), "unexpected container state");
351 
352         BOOST_TRY
353         {
354             // add seeds
355             elements1.push_back(elements_copy[seed1]);                                                      // MAY THROW, STRONG (copy)
356             elements2.push_back(elements_copy[seed2]);                                                      // MAY THROW, STRONG (alloc, copy)
357 
358             // calculate boxes
359             detail::bounds(rtree::element_indexable(elements_copy[seed1], translator), box1);
360             detail::bounds(rtree::element_indexable(elements_copy[seed2], translator), box2);
361 
362             // initialize areas
363             content_type content1 = index::detail::content(box1);
364             content_type content2 = index::detail::content(box2);
365 
366             BOOST_GEOMETRY_INDEX_ASSERT(2 <= elements1_count, "unexpected elements number");
367             size_t remaining = elements1_count - 2;
368 
369             // redistribute the rest of the elements
370             for ( size_t i = 0 ; i < elements1_count ; ++i )
371             {
372                 if (i != seed1 && i != seed2)
373                 {
374                     element_type const& elem = elements_copy[i];
375                     indexable_type const& indexable = rtree::element_indexable(elem, translator);
376 
377                     // if there is small number of elements left and the number of elements in node is lesser than min_elems
378                     // just insert them to this node
379                     if ( elements1.size() + remaining <= parameters.get_min_elements() )
380                     {
381                         elements1.push_back(elem);                                                          // MAY THROW, STRONG (copy)
382                         geometry::expand(box1, indexable);
383                         content1 = index::detail::content(box1);
384                     }
385                     else if ( elements2.size() + remaining <= parameters.get_min_elements() )
386                     {
387                         elements2.push_back(elem);                                                          // MAY THROW, STRONG (alloc, copy)
388                         geometry::expand(box2, indexable);
389                         content2 = index::detail::content(box2);
390                     }
391                     // choose better node and insert element
392                     else
393                     {
394                         // calculate enlarged boxes and areas
395                         Box enlarged_box1(box1);
396                         Box enlarged_box2(box2);
397                         geometry::expand(enlarged_box1, indexable);
398                         geometry::expand(enlarged_box2, indexable);
399                         content_type enlarged_content1 = index::detail::content(enlarged_box1);
400                         content_type enlarged_content2 = index::detail::content(enlarged_box2);
401 
402                         content_type content_increase1 = enlarged_content1 - content1;
403                         content_type content_increase2 = enlarged_content2 - content2;
404 
405                         // choose group which box content have to be enlarged least or has smaller content or has fewer elements
406                         if ( content_increase1 < content_increase2 ||
407                                 ( content_increase1 == content_increase2 &&
408                                     ( content1 < content2 ||
409                                         ( content1 == content2 && elements1.size() <= elements2.size() ) ) ) )
410                         {
411                             elements1.push_back(elem);                                                      // MAY THROW, STRONG (copy)
412                             box1 = enlarged_box1;
413                             content1 = enlarged_content1;
414                         }
415                         else
416                         {
417                             elements2.push_back(elem);                                                      // MAY THROW, STRONG (alloc, copy)
418                             box2 = enlarged_box2;
419                             content2 = enlarged_content2;
420                         }
421                     }
422 
423                     BOOST_GEOMETRY_INDEX_ASSERT(0 < remaining, "unexpected value");
424                     --remaining;
425                 }
426             }
427         }
428         BOOST_CATCH(...)
429         {
430             elements1.clear();
431             elements2.clear();
432 
433             rtree::destroy_elements<Value, Options, Translator, Box, Allocators>::apply(elements_copy, allocators);
434             //elements_copy.clear();
435 
436             BOOST_RETHROW                                                                                     // RETHROW, BASIC
437         }
438         BOOST_CATCH_END
439     }
440 };
441 
442 }} // namespace detail::rtree
443 
444 }}} // namespace boost::geometry::index
445 
446 #endif // BOOST_GEOMETRY_INDEX_DETAIL_RTREE_LINEAR_REDISTRIBUTE_ELEMENTS_HPP
447