Branch data Line data Source code
1 : : //---------------------------------------------------------------------------
2 : : // Class Name: RTree
3 : : // Description: Rectangle tree. Multidimensional access method (efficient
4 : : // method to find ranges of boxes.
5 : : // The algorithm was taken from the following paper:
6 : : // Guttman, A., "R-Trees: A Dynamic Index Structure for
7 : : // Spatial Searching", Proceedings of the SIGMOD
8 : : // Conference, Boston, June 1984, p. 47-57.
9 : : // Creation Date: 3/29/02
10 : : // Owner: David R. White
11 : : //---------------------------------------------------------------------------
12 : :
13 : : //---------------------------------
14 : : //Include Files
15 : : //---------------------------------
16 : : #include "RTree.hpp"
17 : : #include "RTreeNode.hpp"
18 : : #include "CubitBox.hpp"
19 : : #include "CubitVector.hpp"
20 : : #include "DLIList.hpp"
21 : : #include "PriorityQueue.hpp"
22 : : //---------------------------
23 : : //Initialize Static Members
24 : : //---------------------------
25 : :
26 : : #ifdef INLINE_TEMPLATES
27 : : #define MY_INLINE inline
28 : : #else
29 : : #define MY_INLINE
30 : : #endif
31 : 190 : template <class Z> MY_INLINE RTree<Z>::RTree (double tol)
32 : : {
33 : 95 : myRoot = NULL;
34 : 95 : myTolerance = tol;
35 : 95 : maxChildren = 8;
36 : 95 : minChildren = 2;
37 : 95 : }
38 : : template <class Z> MY_INLINE RTree<Z>::RTree (double tol, int max_c, int min_c)
39 : : {
40 : : myRoot = NULL;
41 : : myTolerance = tol;
42 : : maxChildren = max_c;
43 : : minChildren = min_c;
44 : : }
45 : 95 : template <class Z> MY_INLINE RTree<Z>::~RTree()
46 : : {
47 [ + - ][ + - ]: 95 : if ( myRoot != NULL )
[ # # ][ + - ]
48 : : {
49 : : //Go through and get all the children in a list.
50 [ + - ][ + - ]: 95 : DLIList <RTreeNode<Z>*> to_delete;
[ # # ][ + - ]
51 [ + - ][ + - ]: 95 : to_list(to_delete, myRoot);
[ # # ][ + - ]
52 : : int ii;
53 [ + - ][ + + ]: 3607 : for(ii = to_delete.size(); ii > 0; ii-- )
[ + - ][ + + ]
[ # # ][ # # ]
[ + - ][ + + ]
54 [ + - ][ + - ]: 3512 : delete to_delete.pop();
[ + - ][ + - ]
[ + - ][ + - ]
[ # # ][ # # ]
[ # # ][ + - ]
[ + - ][ + - ]
55 [ + - ][ + - ]: 95 : delete myRoot;
[ + - ][ + - ]
[ + - ][ + - ]
[ # # ][ # # ]
[ # # ][ + - ]
[ + - ][ + - ]
56 : : }
57 [ - + ][ - + ]: 95 : }
[ # # ][ - + ]
58 : 715 : template <class Z> MY_INLINE void RTree<Z>::to_list(DLIList <RTreeNode<Z>*> &member_list,
59 : : RTreeNode<Z> *top)
60 : : {
61 : : //Get the children of the top into the list.
62 : : int ii;
63 : : RTreeNode <Z> *curr_node;
64 [ + - ][ + + ]: 4227 : for ( ii = 0; ii < top->num_children(); ii++ )
[ + - ][ + + ]
[ # # ][ # # ]
[ + - ][ + + ]
65 : : {
66 [ + - ][ + - ]: 3512 : curr_node = top->get_child(ii);
[ # # ][ + - ]
67 [ + - ][ + - ]: 3512 : member_list.append(curr_node);
[ # # ][ + - ]
68 : : //don't go below the bottom level...
69 [ + - ][ + + ]: 3512 : if ( curr_node->get_leaf_level() == 0 )
[ + - ][ + + ]
[ # # ][ # # ]
[ + - ][ + + ]
70 : 2892 : continue;
71 [ + - ][ + - ]: 620 : to_list(member_list, curr_node);
[ # # ][ + - ]
72 : : }
73 : 715 : return;
74 : : }
75 : :
76 : 2974 : template <class Z> MY_INLINE CubitStatus RTree<Z>::add(Z data)
77 : : {
78 : : CubitStatus stat;
79 : : RTreeNode<Z> *new_root;
80 : : RTreeNode<Z> *new_node = new RTreeNode<Z> (data, myTolerance, maxChildren,
81 [ + - ][ + - ]: 2974 : minChildren);
[ + - ][ + - ]
[ # # ][ # # ]
[ + - ][ + - ]
82 [ + + ][ + + ]: 2974 : if ( myRoot == NULL )
[ # # ][ + + ]
83 : : {
84 [ + - ][ + - ]: 95 : CubitBox b_box = data->bounding_box();
[ # # ][ + - ]
[ # # ]
85 [ + - ][ + - ]: 95 : myRoot = new RTreeNode <Z> (b_box, maxChildren, minChildren);
[ + - ][ + - ]
[ # # ][ # # ]
[ + - ][ + - ]
86 [ + - ][ + - ]: 95 : myRoot->set_leaf_level(LEAF_RNODE);
[ # # ][ + - ]
87 [ + - ][ + - ]: 95 : stat = myRoot->insert(new_node, new_root);
[ # # ][ + - ]
88 : : //this shouldn't change the root, or fail!
89 [ + - ][ - + ]: 95 : if ( stat != CUBIT_SUCCESS || new_root != NULL )
[ + - ][ - + ]
[ # # ][ # # ]
[ + - ][ - + ]
90 : : {
91 [ # # ][ # # ]: 0 : PRINT_ERROR("Insertion into RTree failed.\n");
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
92 : 0 : return CUBIT_FAILURE;
93 : : }
94 [ + - ][ + - ]: 95 : return CUBIT_SUCCESS;
[ # # ][ + - ]
95 : : }
96 : :
97 [ + - ][ + - ]: 2879 : stat = myRoot->insert(new_node, new_root);
[ # # ][ + - ]
98 [ - + ][ - + ]: 2879 : if ( stat != CUBIT_SUCCESS )
[ # # ][ - + ]
99 : : {
100 [ # # ][ # # ]: 0 : PRINT_ERROR("Insertion into RTree failed.\n");
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
101 : 0 : return CUBIT_FAILURE;
102 : : }
103 [ + + ][ + + ]: 2879 : if ( new_root != NULL )
[ # # ][ + + ]
104 : : {
105 : : //this is fine, it just means we are adding more
106 : : //so the root had to be split...
107 : 116 : myRoot = new_root;
108 : : }
109 : 2974 : return CUBIT_SUCCESS;
110 : : }
111 : 2892 : template <class Z> MY_INLINE CubitStatus RTree<Z>::find(const CubitBox &range_box,
112 : : DLIList <Z> &range_members )
113 : : {
114 : : //Find all of the members of the RTree that intersect this range_box.
115 [ - + ][ - + ]: 2892 : if ( myRoot == NULL )
[ # # ][ - + ]
116 : : {
117 : : // Nothing has been added to this Tree yet, so we are not going to find this
118 : : // object in it.
119 : 0 : return CUBIT_SUCCESS;
120 : : }
121 : 2892 : CubitStatus stat = recursive_find(myRoot, range_box, range_members);
122 [ - + - + : 2892 : if ( stat != CUBIT_SUCCESS )
# # - + ]
123 : 0 : return CUBIT_FAILURE;
124 : : else
125 : 2892 : return CUBIT_SUCCESS;
126 : : }
127 : 65988 : template <class Z> MY_INLINE CubitStatus RTree<Z>::recursive_find(RTreeNode<Z> *rect_tree,
128 : : const CubitBox &range_box,
129 : : DLIList <Z> &range_members )
130 : : {
131 [ + - ][ + - ]: 65988 : CubitBox rect_box = rect_tree->bounding_box();
[ + - ][ + - ]
[ # # ][ # # ]
[ + - ][ + - ]
132 [ + - ][ + + ]: 65988 : if ( !range_box.overlap(myTolerance, rect_box ) )
[ + - ][ + + ]
[ # # ][ # # ]
[ + - ][ + + ]
133 : 36154 : return CUBIT_SUCCESS;
134 : :
135 : : //Now see if this is a data member. If it is, append the data to the
136 : : //list.
137 [ + - ][ + + ]: 29834 : if (rect_tree->is_data() )
[ + - ][ + + ]
[ # # ][ # # ]
[ + - ][ + + ]
138 : : {
139 [ + - ][ + - ]: 17626 : range_members.append(rect_tree->get_data());
[ + - ][ + - ]
[ # # ][ # # ]
[ + - ][ + - ]
140 : 17626 : return CUBIT_SUCCESS;
141 : : }
142 : : //Now if this is anything else we need to keep iterating...
143 [ + - ][ + - ]: 12208 : int loop_size = rect_tree->num_children();
[ # # ][ + - ]
144 : : //We are doing a depth-first search of the tree. Not
145 : : //all branches will need to be followed since they won't
146 : : //all overlap...
147 : : int ii;
148 : : RTreeNode<Z> *curr_node;
149 : : CubitStatus stat;
150 [ + + ][ + + ]: 75304 : for ( ii = 0; ii < loop_size; ii++ )
[ # # ][ + + ]
151 : : {
152 [ + - ][ + - ]: 63096 : curr_node = rect_tree->get_child(ii);
[ # # ][ + - ]
153 [ - + ][ - + ]: 63096 : if ( curr_node == NULL )
[ # # ][ - + ]
154 : : {
155 [ # # ][ # # ]: 0 : PRINT_ERROR("Problems finding boxes in range.\n");
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
156 [ # # ][ # # ]: 0 : assert(curr_node != NULL);
[ # # ][ # # ]
157 : 0 : return CUBIT_FAILURE;
158 : : }
159 [ + - ][ + - ]: 63096 : stat = recursive_find(curr_node, range_box, range_members);
[ # # ][ + - ]
160 [ - + ][ - + ]: 63096 : if ( stat != CUBIT_SUCCESS )
[ # # ][ - + ]
161 : 0 : return stat;
162 : : }
163 : :
164 [ + - ][ + - ]: 65988 : return CUBIT_SUCCESS;
[ # # ][ + - ]
165 : : }
166 : 82 : template <class Z> MY_INLINE CubitBoolean RTree<Z>::remove( Z data )
167 : : {
168 : 82 : RTreeNode<Z> *new_root = NULL;
169 : 82 : CubitBoolean delete_root = CUBIT_FALSE;
170 [ # # ][ # # ]: 82 : CubitBoolean return_val = myRoot->remove( data, new_root, delete_root );
[ # # ][ + - ]
171 [ # # ][ # # ]: 82 : if ( new_root != NULL )
[ # # ][ - + ]
172 : : {
173 : : //Only if we are condensing the tree do we want to delete the root node.
174 : : //There are other reasons the root has changed (rebalance...), in which
175 : : //cases the root is now a child of the new root...
176 [ # # ][ # # ]: 0 : if ( delete_root )
[ # # ][ # # ]
177 [ # # ][ # # ]: 0 : delete myRoot;
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
178 : 0 : myRoot = new_root;
179 : : }
180 : 82 : return return_val;
181 : : }
182 : : //--------------------------------------------------------------------------
183 : : //Algorithm: min_dist_sq
184 : : //Description: Finds the minimum distance squared between the given
185 : : // point and the box. If the point is on or in the box, the
186 : : // min distance is zero.
187 : : //--------------------------------------------------------------------------
188 : : template <class Z> MY_INLINE
189 : 0 : double RTree<Z>::min_dist_sq(CubitVector &q,
190 : : CubitBox &b_box)
191 : : {
192 [ # # ][ # # ]: 0 : CubitVector b_min, b_max;
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
193 [ # # ][ # # ]: 0 : b_min = b_box.minimum();
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
194 [ # # ][ # # ]: 0 : b_max = b_box.maximum();
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
195 : : double dist;
196 [ # # ][ # # ]: 0 : CubitVector r;
[ # # ][ # # ]
197 : :
198 [ # # ][ # # ]: 0 : if ( q.x() < b_min.x() )
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
199 [ # # ][ # # ]: 0 : r.x(b_min.x());
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
200 [ # # ][ # # ]: 0 : else if ( q.x() > b_max.x() )
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
201 [ # # ][ # # ]: 0 : r.x(b_max.x());
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
202 : : else
203 [ # # ][ # # ]: 0 : r.x(q.x());
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
204 : :
205 [ # # ][ # # ]: 0 : if ( q.y() < b_min.y() )
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
206 [ # # ][ # # ]: 0 : r.y(b_min.y());
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
207 [ # # ][ # # ]: 0 : else if ( q.y() > b_max.y() )
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
208 [ # # ][ # # ]: 0 : r.y(b_max.y());
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
209 : : else
210 [ # # ][ # # ]: 0 : r.y(q.y());
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
211 : :
212 [ # # ][ # # ]: 0 : if ( q.z() < b_min.z() )
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
213 [ # # ][ # # ]: 0 : r.z(b_min.z());
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
214 [ # # ][ # # ]: 0 : else if ( q.z() > b_max.z() )
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
215 [ # # ][ # # ]: 0 : r.z(b_max.z());
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
216 : : else
217 [ # # ][ # # ]: 0 : r.z(q.z());
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
218 : :
219 [ # # ][ # # ]: 0 : dist = (q-r).length_squared();
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
220 : :
221 : 0 : return dist;
222 : : }
223 : : template <class Z> MY_INLINE
224 : 0 : bool RTree<Z>::less_than_func(RTreeNode<Z> *&node_a,
225 : : RTreeNode<Z> *&node_b)
226 : : {
227 [ # # ][ # # ]: 0 : if ( node_a->get_dist() < node_b->get_dist() )
[ # # ][ # # ]
228 : 0 : return true;
229 : : else
230 : 0 : return false;
231 : : }
232 : :
233 : : template <class Z> MY_INLINE
234 : 0 : CubitStatus RTree<Z>::k_nearest_neighbor(CubitVector &q,
235 : : int k,
236 : : double &closest_dist,
237 : : DLIList<Z> &nearest_neighbors,
238 : : typename RTree<Z>::DistSqFunc dist_sq_point_data)
239 : : {
240 : : //first create the priority queue.
241 [ # # ][ # # ]: 0 : PriorityQueue< RTreeNode<Z>*> near_queue(RTree<Z>::less_than_func);
[ # # ][ # # ]
242 : :
243 [ # # ][ # # ]: 0 : myRoot->set_dist(0.0);
[ # # ][ # # ]
244 [ # # ][ # # ]: 0 : near_queue.push(myRoot);
[ # # ][ # # ]
245 : : RTreeNode<Z> *element, *child_element;
246 : 0 : int num_found = 0;
247 : : int ii;
248 : : double data_dist, box_dist;
249 : : Z data;
250 : :
251 [ # # ][ # # ]: 0 : while( !near_queue.empty() )
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
252 : : {
253 [ # # ][ # # ]: 0 : element = near_queue.top();
[ # # ][ # # ]
254 [ # # ][ # # ]: 0 : near_queue.pop();
[ # # ][ # # ]
255 [ # # ][ # # ]: 0 : if ( element->is_data() )
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
256 : : {
257 [ # # ][ # # ]: 0 : data = element->get_data();
[ # # ][ # # ]
258 : : //calculate the exact distance.
259 [ # # ][ # # ]: 0 : data_dist = dist_sq_point_data(q, data);
[ # # ][ # # ]
260 : : //compare this distance with the next item's distance.
261 [ # # ][ # # ]: 0 : if ( element->dist_is_box() && !near_queue.empty() &&
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
262 [ # # ][ # # ]: 0 : data_dist > near_queue.top()->get_dist())
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
263 : : {
264 : : //If its bigger, add it back into the list
265 : : //with the updated distance.
266 [ # # ][ # # ]: 0 : element->set_dist(data_dist);
[ # # ][ # # ]
267 [ # # ][ # # ]: 0 : near_queue.push(element);
[ # # ][ # # ]
268 [ # # ][ # # ]: 0 : element->set_dist_is_box(0);
[ # # ][ # # ]
269 : : }
270 : : else
271 : : {
272 [ # # ][ # # ]: 0 : nearest_neighbors.append(element->get_data());
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
273 [ # # ][ # # ]: 0 : if ( num_found == 0 )
[ # # ][ # # ]
274 [ # # ][ # # ]: 0 : closest_dist = element->get_dist();
[ # # ][ # # ]
275 : 0 : num_found++;
276 [ # # ][ # # ]: 0 : if ( num_found == k )
[ # # ][ # # ]
277 : 0 : return CUBIT_SUCCESS;
278 : : }
279 : : }
280 : : else
281 : : {
282 [ # # ][ # # ]: 0 : for ( ii = 0; ii < element->num_children(); ii++ )
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
283 : : {
284 [ # # ][ # # ]: 0 : child_element = element->get_child(ii);
[ # # ][ # # ]
285 [ # # ][ # # ]: 0 : CubitBox bounding_box = child_element->bounding_box();
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
286 [ # # ][ # # ]: 0 : box_dist = min_dist_sq(q, bounding_box);
[ # # ][ # # ]
287 [ # # ][ # # ]: 0 : child_element->set_dist(box_dist);
[ # # ][ # # ]
288 [ # # ][ # # ]: 0 : near_queue.push(child_element);
[ # # ][ # # ]
289 : : }
290 : : }
291 : : }
292 [ # # ][ # # ]: 0 : return CUBIT_FAILURE;
[ # # ][ # # ]
293 : : }
294 : :
295 : : template <class Z> MY_INLINE
296 : : CubitStatus RTree<Z>::find( const CubitVector &ray_origin, const CubitVector &ray_direction,
297 : : DLIList <Z> &range_members)
298 : : {
299 : : //Find all of the members of the RTree that intersect this ray.
300 : : if ( myRoot == NULL )
301 : : {
302 : : // Nothing has been added to this Tree yet, so we are not going to find this
303 : : // object in it.
304 : : return CUBIT_SUCCESS;
305 : : }
306 : : CubitStatus stat = recursive_find(myRoot, ray_origin, ray_direction, range_members);
307 : : if ( stat != CUBIT_SUCCESS )
308 : : return CUBIT_FAILURE;
309 : : else
310 : : return CUBIT_SUCCESS;
311 : : }
312 : :
313 : : template <class Z> MY_INLINE
314 : : CubitStatus RTree<Z>::recursive_find(RTreeNode<Z> *rect_tree,
315 : : const CubitVector &ray_origin,
316 : : const CubitVector &ray_direction,
317 : : DLIList <Z> &range_members)
318 : : {
319 : : CubitBox rect_box = rect_tree->bounding_box();
320 : : //if ( !range_box.overlap(myTolerance, rect_box ) )
321 : : if ( !rect_box.intersect(ray_origin, ray_direction) )
322 : : return CUBIT_SUCCESS;
323 : :
324 : : //Now see if this is a data member. If it is, append the data to the
325 : : //list.
326 : : if (rect_tree->is_data() )
327 : : {
328 : : range_members.append(rect_tree->get_data());
329 : : return CUBIT_SUCCESS;
330 : : }
331 : : //Now if this is anything else we need to keep iterating...
332 : : int loop_size = rect_tree->num_children();
333 : : //We are doing a depth-first search of the tree. Not
334 : : //all branches will need to be followed since they won't
335 : : //all overlap...
336 : : int ii;
337 : : RTreeNode<Z> *curr_node;
338 : : CubitStatus stat;
339 : : for ( ii = 0; ii < loop_size; ii++ )
340 : : {
341 : : curr_node = rect_tree->get_child(ii);
342 : : if ( curr_node == NULL )
343 : : {
344 : : PRINT_ERROR("Problems finding boxes in range.\n");
345 : : assert(curr_node != NULL);
346 : : return CUBIT_FAILURE;
347 : : }
348 : : stat = recursive_find(curr_node, ray_origin, ray_direction, range_members);
349 : : if ( stat != CUBIT_SUCCESS )
350 : : return stat;
351 : : }
352 : :
353 : : return CUBIT_SUCCESS;
354 : : }
|
