1 /*!
2 \file diglib/spindex.c
3
4 \brief Vector library - spatial index (lower level functions)
5
6 Lower level functions for reading/writing/manipulating vectors.
7
8 (C) 2001-2009 by the GRASS Development Team
9
10 This program is free software under the GNU General Public License
11 (>=v2). Read the file COPYING that comes with GRASS for details.
12
13 \author Original author CERL, probably Dave Gerdes
14 \author Update to GRASS 5.7 Radim Blazek
15 \author Update to GRASS 7 Markus Metz
16 */
17
18 #include <stdlib.h>
19 #include <sys/types.h>
20 #include <sys/stat.h>
21 #include <fcntl.h>
22 #include <unistd.h>
23 #include <string.h>
24 #include <grass/vector.h>
25 #include <grass/glocale.h>
26
27 /*!
28 \brief Initit spatial index (nodes, lines, areas, isles)
29
30 \param Plus pointer to Plus_head structure
31
32 \return 1 OK
33 \return 0 on error
34 */
dig_spidx_init(struct Plus_head * Plus)35 int dig_spidx_init(struct Plus_head *Plus)
36 {
37 int ndims;
38
39 ndims = (Plus->with_z != 0) ? 3 : 2;
40 Plus->spidx_with_z = (Plus->with_z != 0);
41
42 G_debug(1, "dig_spidx_init(), %d dims", ndims);
43
44 if (Plus->Spidx_file) {
45 int fd;
46 char *filename;
47
48 filename = G_tempfile();
49 fd = open(filename, O_RDWR | O_CREAT | O_EXCL, 0600);
50 Plus->Node_spidx = RTreeCreateTree(fd, 0, ndims);
51 remove(filename);
52
53 filename = G_tempfile();
54 fd = open(filename, O_RDWR | O_CREAT | O_EXCL, 0600);
55 Plus->Line_spidx = RTreeCreateTree(fd, 0, ndims);
56 remove(filename);
57
58 filename = G_tempfile();
59 fd = open(filename, O_RDWR | O_CREAT | O_EXCL, 0600);
60 Plus->Area_spidx = RTreeCreateTree(fd, 0, ndims);
61 remove(filename);
62
63 filename = G_tempfile();
64 fd = open(filename, O_RDWR | O_CREAT | O_EXCL, 0600);
65 Plus->Isle_spidx = RTreeCreateTree(fd, 0, ndims);
66 remove(filename);
67
68 Plus->Face_spidx = NULL;
69 Plus->Volume_spidx = NULL;
70 Plus->Hole_spidx = NULL;
71
72 if (!Plus->Spidx_new) {
73 close(Plus->Node_spidx->fd);
74 close(Plus->Line_spidx->fd);
75 close(Plus->Area_spidx->fd);
76 close(Plus->Isle_spidx->fd);
77 }
78 }
79 else {
80 Plus->Node_spidx = RTreeCreateTree(-1, 0, ndims);
81 Plus->Line_spidx = RTreeCreateTree(-1, 0, ndims);
82 Plus->Area_spidx = RTreeCreateTree(-1, 0, ndims);
83 Plus->Isle_spidx = RTreeCreateTree(-1, 0, ndims);
84 Plus->Face_spidx = NULL;
85 Plus->Volume_spidx = NULL;
86 Plus->Hole_spidx = NULL;
87 }
88
89 Plus->Node_spidx_offset = 0L;
90 Plus->Line_spidx_offset = 0L;
91 Plus->Area_spidx_offset = 0L;
92 Plus->Isle_spidx_offset = 0L;
93 Plus->Face_spidx_offset = 0L;
94 Plus->Volume_spidx_offset = 0L;
95 Plus->Hole_spidx_offset = 0L;
96
97 Plus->Spidx_built = FALSE;
98
99 return 1;
100 }
101
102 /*!
103 \brief Free spatial index for nodes
104
105 \param Plus pointer to Plus_head structure
106 */
dig_spidx_free_nodes(struct Plus_head * Plus)107 void dig_spidx_free_nodes(struct Plus_head *Plus)
108 {
109 int ndims;
110
111 ndims = Plus->with_z ? 3 : 2;
112
113 /* Node spidx */
114 if (Plus->Node_spidx->fd > -1) {
115 int fd;
116 char *filename;
117
118 if (Plus->Spidx_new)
119 close(Plus->Node_spidx->fd);
120 RTreeDestroyTree(Plus->Node_spidx);
121 filename = G_tempfile();
122 fd = open(filename, O_RDWR | O_CREAT | O_EXCL, 0600);
123 Plus->Node_spidx = RTreeCreateTree(fd, 0, ndims);
124 remove(filename);
125 if (!Plus->Spidx_new)
126 close(Plus->Node_spidx->fd);
127 }
128 else {
129 RTreeDestroyTree(Plus->Node_spidx);
130 Plus->Node_spidx = RTreeCreateTree(-1, 0, ndims);
131 }
132 }
133
134 /*!
135 \brief Free spatial index for lines
136
137 \param Plus pointer to Plus_head structure
138 */
dig_spidx_free_lines(struct Plus_head * Plus)139 void dig_spidx_free_lines(struct Plus_head *Plus)
140 {
141 int ndims;
142
143 ndims = Plus->with_z ? 3 : 2;
144
145 /* Line spidx */
146 if (Plus->Line_spidx->fd > -1) {
147 int fd;
148 char *filename;
149
150 if (Plus->Spidx_new)
151 close(Plus->Line_spidx->fd);
152 RTreeDestroyTree(Plus->Line_spidx);
153 filename = G_tempfile();
154 fd = open(filename, O_RDWR | O_CREAT | O_EXCL, 0600);
155 Plus->Line_spidx = RTreeCreateTree(fd, 0, ndims);
156 remove(filename);
157 if (!Plus->Spidx_new)
158 close(Plus->Line_spidx->fd);
159 }
160 else {
161 RTreeDestroyTree(Plus->Line_spidx);
162 Plus->Line_spidx = RTreeCreateTree(-1, 0, ndims);
163 }
164 }
165
166 /*!
167 \brief Reset spatial index for areas
168
169 \param Plus pointer to Plus_head structure
170 */
dig_spidx_free_areas(struct Plus_head * Plus)171 void dig_spidx_free_areas(struct Plus_head *Plus)
172 {
173 int ndims;
174
175 ndims = Plus->with_z ? 3 : 2;
176
177 /* Area spidx */
178 if (Plus->Area_spidx->fd > -1) {
179 int fd;
180 char *filename;
181
182 if (Plus->Spidx_new)
183 close(Plus->Area_spidx->fd);
184 RTreeDestroyTree(Plus->Area_spidx);
185 filename = G_tempfile();
186 fd = open(filename, O_RDWR | O_CREAT | O_EXCL, 0600);
187 Plus->Area_spidx = RTreeCreateTree(fd, 0, ndims);
188 remove(filename);
189 if (!Plus->Spidx_new)
190 close(Plus->Area_spidx->fd);
191 }
192 else {
193 RTreeDestroyTree(Plus->Area_spidx);
194 Plus->Area_spidx = RTreeCreateTree(-1, 0, ndims);
195 }
196 }
197
198 /*!
199 \brief Reset spatial index for isles
200
201 \param Plus pointer to Plus_head structure
202 */
dig_spidx_free_isles(struct Plus_head * Plus)203 void dig_spidx_free_isles(struct Plus_head *Plus)
204 {
205 int ndims;
206
207 ndims = Plus->with_z ? 3 : 2;
208
209 /* Isle spidx */
210 if (Plus->Isle_spidx->fd > -1) {
211 int fd;
212 char *filename;
213
214 if (Plus->Spidx_new)
215 close(Plus->Isle_spidx->fd);
216 RTreeDestroyTree(Plus->Isle_spidx);
217 filename = G_tempfile();
218 fd = open(filename, O_RDWR | O_CREAT | O_EXCL, 0600);
219 Plus->Isle_spidx = RTreeCreateTree(fd, 0, ndims);
220 remove(filename);
221 if (!Plus->Spidx_new)
222 close(Plus->Isle_spidx->fd);
223 }
224 else {
225 RTreeDestroyTree(Plus->Isle_spidx);
226 Plus->Isle_spidx = RTreeCreateTree(-1, 0, ndims);
227 }
228 }
229
230 /*!
231 \brief Free spatial index (nodes, lines, areas, isles)
232
233 \param Plus pointer to Plus_head structure
234 */
dig_spidx_free(struct Plus_head * Plus)235 void dig_spidx_free(struct Plus_head *Plus)
236 {
237 /* close tmp files */
238 if (Plus->Spidx_new) {
239 /* Node spidx */
240 if (Plus->Node_spidx->fd > -1)
241 close(Plus->Node_spidx->fd);
242 /* Line spidx */
243 if (Plus->Spidx_new && Plus->Line_spidx->fd > -1)
244 close(Plus->Line_spidx->fd);
245 /* Area spidx */
246 if (Plus->Area_spidx->fd > -1)
247 close(Plus->Area_spidx->fd);
248 /* Isle spidx */
249 if (Plus->Isle_spidx->fd > -1)
250 close(Plus->Isle_spidx->fd);
251 }
252
253 /* destroy tree structures */
254 /* Node spidx */
255 if (Plus->Node_spidx)
256 RTreeDestroyTree(Plus->Node_spidx);
257 /* Line spidx */
258 if (Plus->Line_spidx)
259 RTreeDestroyTree(Plus->Line_spidx);
260 /* Area spidx */
261 if (Plus->Area_spidx)
262 RTreeDestroyTree(Plus->Area_spidx);
263 /* Isle spidx */
264 if (Plus->Isle_spidx)
265 RTreeDestroyTree(Plus->Isle_spidx);
266
267 /* 3D future : */
268 /* Face spidx */
269 /* Volume spidx */
270 /* Hole spidx */
271 }
272
273 /*!
274 \brief Add new node to spatial index
275
276 \param Plus pointer to Plus_head structure
277 \param node node id
278 \param x,y,z node coordinates
279
280 \return 1 OK
281 \return 0 on error
282 */
283 int
dig_spidx_add_node(struct Plus_head * Plus,int node,double x,double y,double z)284 dig_spidx_add_node(struct Plus_head *Plus, int node,
285 double x, double y, double z)
286 {
287 static struct RTree_Rect rect;
288 static int rect_init = 0;
289
290 if (!rect_init) {
291 /* always 6 sides for 3D */
292 rect.boundary = G_malloc(6 * sizeof(RectReal));
293 rect_init = 6;
294 }
295
296 G_debug(3, "dig_spidx_add_node(): node = %d, x,y,z = %f, %f, %f", node, x,
297 y, z);
298
299 rect.boundary[0] = x;
300 rect.boundary[1] = y;
301 rect.boundary[2] = z;
302 rect.boundary[3] = x;
303 rect.boundary[4] = y;
304 rect.boundary[5] = z;
305 RTreeInsertRect(&rect, node, Plus->Node_spidx);
306
307 return 1;
308 }
309
310 /*!
311 \brief Add new line to spatial index
312
313 \param Plus pointer to Plus_head structure
314 \param line line id
315 \param box bounding box
316
317 \return 0
318 */
dig_spidx_add_line(struct Plus_head * Plus,int line,const struct bound_box * box)319 int dig_spidx_add_line(struct Plus_head *Plus, int line,
320 const struct bound_box *box)
321 {
322 static struct RTree_Rect rect;
323 static int rect_init = 0;
324
325 if (!rect_init) {
326 /* always 6 sides for 3D */
327 rect.boundary = G_malloc(6 * sizeof(RectReal));
328 rect_init = 6;
329 }
330
331 G_debug(3, "dig_spidx_add_line(): line = %d", line);
332
333 rect.boundary[0] = box->W;
334 rect.boundary[1] = box->S;
335 rect.boundary[2] = box->B;
336 rect.boundary[3] = box->E;
337 rect.boundary[4] = box->N;
338 rect.boundary[5] = box->T;
339 RTreeInsertRect(&rect, line, Plus->Line_spidx);
340
341 return 0;
342 }
343
344 /*!
345 \brief Add new area to spatial index
346
347 \param Plus pointer to Plus_head structure
348 \param area area id
349 \param box bounding box
350
351 \return 0
352 */
dig_spidx_add_area(struct Plus_head * Plus,int area,const struct bound_box * box)353 int dig_spidx_add_area(struct Plus_head *Plus, int area,
354 const struct bound_box *box)
355 {
356 static struct RTree_Rect rect;
357 static int rect_init = 0;
358
359 if (!rect_init) {
360 /* always 6 sides for 3D */
361 rect.boundary = G_malloc(6 * sizeof(RectReal));
362 rect_init = 6;
363 }
364
365 G_debug(3, "dig_spidx_add_area(): area = %d", area);
366
367 rect.boundary[0] = box->W;
368 rect.boundary[1] = box->S;
369 rect.boundary[2] = box->B;
370 rect.boundary[3] = box->E;
371 rect.boundary[4] = box->N;
372 rect.boundary[5] = box->T;
373 RTreeInsertRect(&rect, area, Plus->Area_spidx);
374
375 return 0;
376 }
377
378 /*!
379 \brief Add new island to spatial index
380
381 \param Plus pointer to Plus_head structure
382 \param isle isle id
383 \param box bounding box
384
385 \return 0
386 */
387
dig_spidx_add_isle(struct Plus_head * Plus,int isle,const struct bound_box * box)388 int dig_spidx_add_isle(struct Plus_head *Plus, int isle,
389 const struct bound_box *box)
390 {
391 static struct RTree_Rect rect;
392 static int rect_init = 0;
393
394 if (!rect_init) {
395 /* always 6 sides for 3D */
396 rect.boundary = G_malloc(6 * sizeof(RectReal));
397 rect_init = 6;
398 }
399
400 G_debug(3, "dig_spidx_add_isle(): isle = %d", isle);
401
402 rect.boundary[0] = box->W;
403 rect.boundary[1] = box->S;
404 rect.boundary[2] = box->B;
405 rect.boundary[3] = box->E;
406 rect.boundary[4] = box->N;
407 rect.boundary[5] = box->T;
408 RTreeInsertRect(&rect, isle, Plus->Isle_spidx);
409
410 return 0;
411 }
412
413 /*!
414 \brief Delete node from spatial index
415
416 G_fatal_error() called on error.
417
418 \param Plus pointer to Plus_head structure
419 \param node node id
420
421 \return 0
422 */
dig_spidx_del_node(struct Plus_head * Plus,int node)423 int dig_spidx_del_node(struct Plus_head *Plus, int node)
424 {
425 int ret;
426 struct P_node *Node;
427 static struct RTree_Rect rect;
428 static int rect_init = 0;
429
430 if (!rect_init) {
431 /* always 6 sides for 3D */
432 rect.boundary = G_malloc(6 * sizeof(RectReal));
433 rect_init = 6;
434 }
435
436 G_debug(3, "dig_spidx_del_node(): node = %d", node);
437
438 Node = Plus->Node[node];
439
440 rect.boundary[0] = Node->x;
441 rect.boundary[1] = Node->y;
442 rect.boundary[2] = Node->z;
443 rect.boundary[3] = Node->x;
444 rect.boundary[4] = Node->y;
445 rect.boundary[5] = Node->z;
446
447 ret = RTreeDeleteRect(&rect, node, Plus->Node_spidx);
448
449 if (ret)
450 G_fatal_error(_("Unable to delete node %d from spatial index"), node);
451
452 return 0;
453 }
454
455 /*!
456 \brief Delete line from spatial index
457
458 G_fatal_error() called on error.
459
460 \param Plus pointer to Plus_head structure
461 \param line line id
462 \param x,y,z coordinates
463
464 \return 0
465 */
dig_spidx_del_line(struct Plus_head * Plus,int line,double x,double y,double z)466 int dig_spidx_del_line(struct Plus_head *Plus, int line,
467 double x, double y, double z)
468 {
469 int ret;
470 static struct RTree_Rect rect;
471 static int rect_init = 0;
472
473 if (!rect_init) {
474 /* always 6 sides for 3D */
475 rect.boundary = G_malloc(6 * sizeof(RectReal));
476 rect_init = 6;
477 }
478
479 G_debug(3, "dig_spidx_del_line(): line = %d", line);
480
481 rect.boundary[0] = x;
482 rect.boundary[1] = y;
483 rect.boundary[2] = z;
484 rect.boundary[3] = x;
485 rect.boundary[4] = y;
486 rect.boundary[5] = z;
487
488 ret = RTreeDeleteRect(&rect, line, Plus->Line_spidx);
489
490 G_debug(3, " ret = %d", ret);
491
492 if (ret)
493 G_fatal_error(_("Unable to delete line %d from spatial index"), line);
494
495 return 0;
496 }
497
498 /*!
499 \brief Delete area from spatial index
500
501 G_fatal_error() called on error.
502
503 \param Plus pointer to Plus_head structure
504 \param area area id
505
506 \return 0
507 */
dig_spidx_del_area(struct Plus_head * Plus,int area)508 int dig_spidx_del_area(struct Plus_head *Plus, int area)
509 {
510 int ret;
511 struct P_area *Area;
512 struct P_line *Line;
513 struct P_node *Node;
514 struct P_topo_b *topo;
515 static struct RTree_Rect rect;
516 static int rect_init = 0;
517
518 if (!rect_init) {
519 /* always 6 sides for 3D */
520 rect.boundary = G_malloc(6 * sizeof(RectReal));
521 rect_init = 6;
522 }
523
524 G_debug(3, "dig_spidx_del_area(): area = %d", area);
525
526 Area = Plus->Area[area];
527
528 if (Area == NULL) {
529 G_fatal_error(_("Attempt to delete sidx for dead area"));
530 }
531
532 Line = Plus->Line[abs(Area->lines[0])];
533 topo = (struct P_topo_b *)Line->topo;
534 Node = Plus->Node[topo->N1];
535
536 rect.boundary[0] = Node->x;
537 rect.boundary[1] = Node->y;
538 rect.boundary[2] = Node->z;
539 rect.boundary[3] = Node->x;
540 rect.boundary[4] = Node->y;
541 rect.boundary[5] = Node->z;
542
543 ret = RTreeDeleteRect(&rect, area, Plus->Area_spidx);
544
545 if (ret)
546 G_fatal_error(_("Unable to delete area %d from spatial index"), area);
547
548 return 0;
549 }
550
551 /*!
552 \brief Delete isle from spatial index
553
554 G_fatal_error() called on error.
555
556 \param Plus pointer to Plus_head structure
557 \param isle isle id
558
559 \return 0
560 */
dig_spidx_del_isle(struct Plus_head * Plus,int isle)561 int dig_spidx_del_isle(struct Plus_head *Plus, int isle)
562 {
563 int ret;
564 struct P_isle *Isle;
565 struct P_line *Line;
566 struct P_node *Node;
567 struct P_topo_b *topo;
568 static struct RTree_Rect rect;
569 static int rect_init = 0;
570
571 if (!rect_init) {
572 /* always 6 sides for 3D */
573 rect.boundary = G_malloc(6 * sizeof(RectReal));
574 rect_init = 6;
575 }
576
577 G_debug(3, "dig_spidx_del_isle(): isle = %d", isle);
578
579 Isle = Plus->Isle[isle];
580
581 Line = Plus->Line[abs(Isle->lines[0])];
582 topo = (struct P_topo_b *)Line->topo;
583 Node = Plus->Node[topo->N1];
584
585 rect.boundary[0] = Node->x;
586 rect.boundary[1] = Node->y;
587 rect.boundary[2] = Node->z;
588 rect.boundary[3] = Node->x;
589 rect.boundary[4] = Node->y;
590 rect.boundary[5] = Node->z;
591
592 ret = RTreeDeleteRect(&rect, isle, Plus->Isle_spidx);
593
594 if (ret)
595 G_fatal_error(_("Unable to delete isle %d from spatial index"), isle);
596
597 return 0;
598 }
599
600 /* This function is called by RTreeSearch() to add selected node/line/area/isle to the list */
_add_item(int id,const struct RTree_Rect * rect,struct ilist * list)601 static int _add_item(int id, const struct RTree_Rect *rect,
602 struct ilist *list)
603 {
604 G_ilist_add(list, id);
605 return 1;
606 }
607
608 /* This function is called by RTreeSearch() to add
609 * selected node/line/area/isle to the box list */
_add_item_with_box(int id,const struct RTree_Rect * rect,struct boxlist * list)610 static int _add_item_with_box(int id, const struct RTree_Rect *rect,
611 struct boxlist *list)
612 {
613 struct bound_box box;
614
615 box.W = rect->boundary[0];
616 box.S = rect->boundary[1];
617 box.B = rect->boundary[2];
618 box.E = rect->boundary[3];
619 box.N = rect->boundary[4];
620 box.T = rect->boundary[5];
621
622 dig_boxlist_add(list, id, &box);
623 return 1;
624 }
625
626 struct boxid
627 {
628 int id;
629 struct bound_box *box;
630 };
631
632 /* This function is called by RTreeSearch() to add
633 * selected node/line/area/isle to the box list */
_set_item_box(int id,const struct RTree_Rect * rect,struct boxid * box_id)634 static int _set_item_box(int id, const struct RTree_Rect *rect,
635 struct boxid *box_id)
636 {
637 if (id == box_id->id) {
638
639 box_id->box->W = rect->boundary[0];
640 box_id->box->S = rect->boundary[1];
641 box_id->box->B = rect->boundary[2];
642 box_id->box->E = rect->boundary[3];
643 box_id->box->N = rect->boundary[4];
644 box_id->box->T = rect->boundary[5];
645
646 return 0;
647 }
648
649 return 1;
650 }
651
652 /*!
653 \brief Select nodes by bbox
654
655 \param Plus pointer to Plus_head structure
656 \param box bounding box
657 \param list list of selected lines
658
659 \return number of selected nodes
660 \return -1 on error
661 */
662 int
dig_select_nodes(struct Plus_head * Plus,const struct bound_box * box,struct ilist * list)663 dig_select_nodes(struct Plus_head *Plus, const struct bound_box * box,
664 struct ilist *list)
665 {
666 static struct RTree_Rect rect;
667 static int rect_init = 0;
668
669 if (!rect_init) {
670 /* always 6 sides for 3D */
671 rect.boundary = G_malloc(6 * sizeof(RectReal));
672 rect_init = 6;
673 }
674
675 G_debug(3, "dig_select_nodes()");
676
677 list->n_values = 0;
678
679 rect.boundary[0] = box->W;
680 rect.boundary[1] = box->S;
681 rect.boundary[2] = box->B;
682 rect.boundary[3] = box->E;
683 rect.boundary[4] = box->N;
684 rect.boundary[5] = box->T;
685
686 if (Plus->Spidx_new)
687 RTreeSearch(Plus->Node_spidx, &rect, (void *)_add_item, list);
688 else
689 rtree_search(Plus->Node_spidx, &rect, (void *)_add_item, list, Plus);
690
691 return (list->n_values);
692 }
693
694 /* This function is called by RTreeSearch() for nodes to find the node id */
_add_node(int id,const struct RTree_Rect * rect,int * node)695 static int _add_node(int id, const struct RTree_Rect *rect, int *node)
696 {
697 *node = id;
698 return 0;
699 }
700
701 /*!
702 \brief Find one node by coordinates
703
704 \param Plus pointer to Plus_head structure
705 \param x,y,z coordinates
706
707 \return number of node
708 \return 0 not found
709 */
dig_find_node(struct Plus_head * Plus,double x,double y,double z)710 int dig_find_node(struct Plus_head *Plus, double x, double y, double z)
711 {
712 int node;
713 static struct RTree_Rect rect;
714 static int rect_init = 0;
715
716 if (!rect_init) {
717 /* always 6 sides for 3D */
718 rect.boundary = G_malloc(6 * sizeof(RectReal));
719 rect_init = 6;
720 }
721
722 G_debug(3, "dig_find_node()");
723
724 rect.boundary[0] = x;
725 rect.boundary[1] = y;
726 rect.boundary[2] = z;
727 rect.boundary[3] = x;
728 rect.boundary[4] = y;
729 rect.boundary[5] = z;
730
731 node = 0;
732 if (Plus->Spidx_new)
733 RTreeSearch(Plus->Node_spidx, &rect, (void *)_add_node, &node);
734 else
735 rtree_search(Plus->Node_spidx, &rect, (void *)_add_node, &node, Plus);
736
737 return node;
738 }
739
740 /*!
741 \brief Select lines with boxes by box
742
743 \param Plus pointer to Plus_head structure
744 \param box bounding box
745 \param list boxlist of selected lines
746
747 \return number of selected lines
748 \return 0 not found
749 */
dig_select_lines(struct Plus_head * Plus,const struct bound_box * box,struct boxlist * list)750 int dig_select_lines(struct Plus_head *Plus, const struct bound_box *box,
751 struct boxlist *list)
752 {
753 static struct RTree_Rect rect;
754 static int rect_init = 0;
755
756 if (!rect_init) {
757 /* always 6 sides for 3D */
758 rect.boundary = G_malloc(6 * sizeof(RectReal));
759 rect_init = 6;
760 }
761
762 G_debug(3, "dig_select_lines_with_box()");
763
764 list->n_values = 0;
765
766 rect.boundary[0] = box->W;
767 rect.boundary[1] = box->S;
768 rect.boundary[2] = box->B;
769 rect.boundary[3] = box->E;
770 rect.boundary[4] = box->N;
771 rect.boundary[5] = box->T;
772
773 if (Plus->Spidx_new)
774 RTreeSearch(Plus->Line_spidx, &rect, (void *)_add_item_with_box, list);
775 else
776 rtree_search(Plus->Line_spidx, &rect, (void *)_add_item_with_box, list, Plus);
777
778 return (list->n_values);
779 }
780
781 /*!
782 \brief Find box for line
783
784 \param Plus pointer to Plus_head structure
785 \param line line id
786 \param[out] box bounding box
787
788 \return > 0 bounding box for line found
789 \return 0 not found
790 */
dig_find_line_box(struct Plus_head * Plus,int line,struct bound_box * box)791 int dig_find_line_box(struct Plus_head *Plus, int line,
792 struct bound_box *box)
793 {
794 int ret, type;
795 struct P_line *Line;
796 struct boxid box_id;
797 static struct RTree_Rect rect;
798 static int rect_init = 0;
799
800 G_debug(3, "dig_find_line_box()");
801
802 if (!rect_init) {
803 /* always 6 sides for 3D */
804 rect.boundary = G_malloc(6 * sizeof(RectReal));
805 rect_init = 6;
806 }
807
808 Line = Plus->Line[line];
809 type = Line->type;
810
811 /* GV_LINES: retrieve box from spatial index */
812 if (type & GV_LINES) {
813 struct P_node *Node = NULL;
814
815 if (type == GV_LINE) {
816 struct P_topo_l *topo = (struct P_topo_l *)Line->topo;
817
818 Node = Plus->Node[topo->N1];
819 }
820 else if (type == GV_BOUNDARY) {
821 struct P_topo_b *topo = (struct P_topo_b *)Line->topo;
822
823 Node = Plus->Node[topo->N1];
824 }
825
826 rect.boundary[0] = Node->x;
827 rect.boundary[1] = Node->y;
828 rect.boundary[2] = Node->z;
829 rect.boundary[3] = Node->x;
830 rect.boundary[4] = Node->y;
831 rect.boundary[5] = Node->z;
832
833 box_id.id = line;
834 box_id.box = box;
835
836 if (Plus->Spidx_new)
837 ret = RTreeSearch(Plus->Line_spidx, &rect, (void *)_set_item_box, &box_id);
838 else
839 ret = rtree_search(Plus->Line_spidx, &rect, (void *)_set_item_box, &box_id, Plus);
840
841 return ret;
842 }
843
844 /* do not translate this error because
845 * 1. this error is not supposed to happen
846 * 2. the maintainer at which this message is directed prefers english */
847 G_fatal_error("Bug in vector lib: dig_find_line_box() may only be used for lines and boundaries.");
848
849 return 0;
850 }
851
852 /*!
853 \brief Select areas with boxes by box
854
855 \param Plus pointer to Plus_head structure
856 \param box bounding box
857 \param list boxlist of selected areas
858
859 \return number of selected areas
860 */
861 int
dig_select_areas(struct Plus_head * Plus,const struct bound_box * box,struct boxlist * list)862 dig_select_areas(struct Plus_head *Plus, const struct bound_box *box,
863 struct boxlist *list)
864 {
865 static struct RTree_Rect rect;
866 static int rect_init = 0;
867
868 if (!rect_init) {
869 /* always 6 sides for 3D */
870 rect.boundary = G_malloc(6 * sizeof(RectReal));
871 rect_init = 6;
872 }
873
874 G_debug(3, "dig_select_areas_with_box()");
875
876 list->n_values = 0;
877
878 rect.boundary[0] = box->W;
879 rect.boundary[1] = box->S;
880 rect.boundary[2] = box->B;
881 rect.boundary[3] = box->E;
882 rect.boundary[4] = box->N;
883 rect.boundary[5] = box->T;
884
885 if (Plus->Spidx_new)
886 RTreeSearch(Plus->Area_spidx, &rect, (void *)_add_item_with_box, list);
887 else
888 rtree_search(Plus->Area_spidx, &rect, (void *)_add_item_with_box, list, Plus);
889
890 return (list->n_values);
891 }
892
893 /*!
894 \brief Find bounding box for given area
895
896 \param Plus pointer to Plus_head structure
897 \param area area id
898 \param[out] box bounding box
899
900 \return > 0 bounding box for area found
901 \return 0 not found
902 */
dig_find_area_box(struct Plus_head * Plus,int area,struct bound_box * box)903 int dig_find_area_box(struct Plus_head *Plus, int area,
904 struct bound_box *box)
905 {
906 int ret;
907 struct boxid box_id;
908 struct P_area *Area;
909 struct P_line *Line;
910 struct P_node *Node;
911 struct P_topo_b *topo;
912 static struct RTree_Rect rect;
913 static int rect_init = 0;
914
915 G_debug(3, "dig_find_area_box()");
916
917 if (!rect_init) {
918 /* always 6 sides for 3D */
919 rect.boundary = G_malloc(6 * sizeof(RectReal));
920 rect_init = 6;
921 }
922
923 Area = Plus->Area[area];
924 Line = Plus->Line[abs(Area->lines[0])];
925 topo = (struct P_topo_b *)Line->topo;
926 Node = Plus->Node[topo->N1];
927
928 rect.boundary[0] = Node->x;
929 rect.boundary[1] = Node->y;
930 rect.boundary[2] = Node->z;
931 rect.boundary[3] = Node->x;
932 rect.boundary[4] = Node->y;
933 rect.boundary[5] = Node->z;
934
935 box_id.id = area;
936 box_id.box = box;
937
938 if (Plus->Spidx_new)
939 ret = RTreeSearch(Plus->Area_spidx, &rect, (void *)_set_item_box, &box_id);
940 else
941 ret = rtree_search(Plus->Area_spidx, &rect, (void *)_set_item_box, &box_id, Plus);
942
943 return ret;
944 }
945
946 /*!
947 \brief Select isles with boxes by box
948
949 \param Plus pointer to Plus_head structure
950 \param box bounding box
951 \param list boxlist of selected isles
952
953 \return number of selected isles
954 */
955 int
dig_select_isles(struct Plus_head * Plus,const struct bound_box * box,struct boxlist * list)956 dig_select_isles(struct Plus_head *Plus, const struct bound_box * box,
957 struct boxlist *list)
958 {
959 static struct RTree_Rect rect;
960 static int rect_init = 0;
961
962 if (!rect_init) {
963 /* always 6 sides for 3D */
964 rect.boundary = G_malloc(6 * sizeof(RectReal));
965 rect_init = 6;
966 }
967
968 G_debug(3, "dig_select_areas_with_box()");
969
970 list->n_values = 0;
971
972 rect.boundary[0] = box->W;
973 rect.boundary[1] = box->S;
974 rect.boundary[2] = box->B;
975 rect.boundary[3] = box->E;
976 rect.boundary[4] = box->N;
977 rect.boundary[5] = box->T;
978
979 if (Plus->Spidx_new)
980 RTreeSearch(Plus->Isle_spidx, &rect, (void *)_add_item_with_box, list);
981 else
982 rtree_search(Plus->Isle_spidx, &rect, (void *)_add_item_with_box, list, Plus);
983
984 return (list->n_values);
985 }
986
987 /*!
988 \brief Find box for isle
989
990 \param Plus pointer to Plus_head structure
991 \param isle isle id
992 \param[out] box bounding box
993
994 \return > 0 bounding box for isle found
995 \return 0 not found
996 */
dig_find_isle_box(struct Plus_head * Plus,int isle,struct bound_box * box)997 int dig_find_isle_box(struct Plus_head *Plus, int isle,
998 struct bound_box *box)
999 {
1000 int ret;
1001 struct boxid box_id;
1002 struct P_isle *Isle;
1003 struct P_line *Line;
1004 struct P_node *Node;
1005 struct P_topo_b *topo;
1006 static struct RTree_Rect rect;
1007 static int rect_init = 0;
1008
1009 G_debug(3, "dig_find_isle_box()");
1010
1011 if (!rect_init) {
1012 /* always 6 sides for 3D */
1013 rect.boundary = G_malloc(6 * sizeof(RectReal));
1014 rect_init = 6;
1015 }
1016
1017 Isle = Plus->Isle[isle];
1018 Line = Plus->Line[abs(Isle->lines[0])];
1019 topo = (struct P_topo_b *)Line->topo;
1020 Node = Plus->Node[topo->N1];
1021
1022 rect.boundary[0] = Node->x;
1023 rect.boundary[1] = Node->y;
1024 rect.boundary[2] = Node->z;
1025 rect.boundary[3] = Node->x;
1026 rect.boundary[4] = Node->y;
1027 rect.boundary[5] = Node->z;
1028
1029 box_id.id = isle;
1030 box_id.box = box;
1031
1032 if (Plus->Spidx_new)
1033 ret = RTreeSearch(Plus->Isle_spidx, &rect, (void *)_set_item_box, &box_id);
1034 else
1035 ret = rtree_search(Plus->Isle_spidx, &rect, (void *)_set_item_box, &box_id, Plus);
1036
1037 return ret;
1038 }
1039