GRASS 8 Programmer's Manual 8.6.0dev(2026)-f6f2c534ea
Loading...
Searching...
No Matches
vector/Vlib/find.c
Go to the documentation of this file.
1/*!
2 * \file lib/vector/Vlib/find.c
3 *
4 * \brief Vector library - Find nearest vector feature.
5 *
6 * Higher 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 or Mike
14 * Higgins.
15 * \author Update to GRASS 5.7 Radim Blazek and David D. Gray.
16 */
17
18#include <stdlib.h>
19#include <math.h>
20#include <grass/vector.h>
21
22#ifndef HUGE_VAL
23#define HUGE_VAL 9999999999999.0
24#endif
25
26/* for qsort */
27
28typedef struct {
29 int i;
30 double size;
31 struct bound_box box;
32} BOX_SIZE;
33
34static int sort_by_size(const void *a, const void *b)
35{
36 BOX_SIZE *as = (BOX_SIZE *)a;
37 BOX_SIZE *bs = (BOX_SIZE *)b;
38
39 if (as->size < bs->size)
40 return -1;
41
42 return (as->size > bs->size);
43}
44
45/*!
46 * \brief Find the nearest node.
47 *
48 * \param Map vector map
49 * \param ux,uy,uz point coordinates
50 * \param maxdist max distance from the line
51 * \param with_z 3D (WITH_Z, WITHOUT_Z)
52 *
53 * \return number of nearest node
54 * \return 0 if not found
55 */
56int Vect_find_node(struct Map_info *Map, double ux, double uy, double uz,
57 double maxdist, int with_z)
58{
59 int i, nnodes, node;
60 struct bound_box box;
61 struct ilist *NList;
62 double x, y, z;
63 double cur_dist, dist;
64
65 G_debug(3, "Vect_find_node() for %f %f %f maxdist = %f", ux, uy, uz,
66 maxdist);
68
69 /* Select all nodes in box */
70 box.N = uy + maxdist;
71 box.S = uy - maxdist;
72 box.E = ux + maxdist;
73 box.W = ux - maxdist;
74 if (with_z) {
75 box.T = uz + maxdist;
76 box.B = uz - maxdist;
77 }
78 else {
79 box.T = HUGE_VAL;
80 box.B = -HUGE_VAL;
81 }
82
84 G_debug(3, " %d nodes in box", nnodes);
85
86 if (nnodes == 0)
87 return 0;
88
89 /* find nearest */
91 node = 0;
92 for (i = 0; i < nnodes; i++) {
93 Vect_get_node_coor(Map, NList->value[i], &x, &y, &z);
94 dist = Vect_points_distance(ux, uy, uz, x, y, z, with_z);
95 if (dist < cur_dist) {
96 cur_dist = dist;
97 node = i;
98 }
99 }
100 G_debug(3, " nearest node %d in distance %f", NList->value[node],
101 cur_dist);
102
103 /* Check if in max distance */
104 if (cur_dist <= maxdist)
105 return (NList->value[node]);
106 else
107 return 0;
108}
109
110/*!
111 * \brief Find the nearest line.
112 *
113 * \param map vector map
114 * \param ux,uy,uz points coordinates
115 * \param type feature type (GV_LINE, GV_POINT, GV_BOUNDARY or GV_CENTROID)
116 * if only want to search certain types of lines or -1 if search all lines
117 * \param maxdist max distance from the line
118 * \param with_z 3D (WITH_Z, WITHOUT_Z)
119 * \param exclude if > 0 number of line which should be excluded from selection.
120 * May be useful if we need line nearest to other one.
121 *
122 * \return number of nearest line
123 * \return 0 if not found
124 *
125 */
126int Vect_find_line(struct Map_info *map, double ux, double uy, double uz,
127 int type, double maxdist, int with_z, int exclude)
128{
129 int line;
130 struct ilist *exclude_list;
131
133
135
136 line = Vect_find_line_list(map, ux, uy, uz, type, maxdist, with_z,
138
140
141 return line;
142}
143
144/*!
145 * \brief Find the nearest line(s).
146 *
147 * \param map vector map
148 * \param ux,uy,uz points coordinates
149 * \param type feature type (GV_LINE, GV_POINT, GV_BOUNDARY or GV_CENTROID)
150 * if only want to search certain types of lines or -1 if search all lines
151 * \param maxdist max distance from the line
152 * \param with_z 3D (WITH_Z, WITHOUT_Z)
153 * \param exclude list of lines which should be excluded from selection
154 * \param found list of found lines (or NULL)
155 *
156 * \return number of nearest line
157 * \return 0 if not found
158 */
159int Vect_find_line_list(struct Map_info *map, double ux, double uy, double uz,
160 int type, double maxdist, int with_z,
161 const struct ilist *exclude, struct ilist *found)
162{
163 int choice;
164 double new_dist;
165 double cur_dist;
166 int gotone;
167 int i, line;
168 static struct line_pnts *Points;
169 static int first_time = 1;
170 struct bound_box box;
171 struct boxlist *List;
172
173 G_debug(3, "Vect_find_line_list() for %f %f %f type = %d maxdist = %f", ux,
174 uy, uz, type, maxdist);
175
176 if (first_time) {
177 Points = Vect_new_line_struct();
178 first_time = 0;
179 }
180
181 gotone = 0;
182 choice = 0;
184
185 box.N = uy + maxdist;
186 box.S = uy - maxdist;
187 box.E = ux + maxdist;
188 box.W = ux - maxdist;
189 if (with_z) {
190 box.T = uz + maxdist;
191 box.B = uz - maxdist;
192 }
193 else {
195 box.B = -PORT_DOUBLE_MAX;
196 }
197
199
200 if (found)
202
203 Vect_select_lines_by_box(map, &box, type, List);
204 for (i = 0; i < List->n_values; i++) {
205 line = List->id[i];
206 if (Vect_val_in_list(exclude, line)) {
207 G_debug(3, " line = %d exclude", line);
208 continue;
209 }
210
211 /* No more needed */
212 /*
213 Line = Plus->Line[line];
214 if ( Line == NULL ) continue;
215 if ( !(type & Line->type) ) continue;
216 */
217
218 Vect_read_line(map, Points, NULL, line);
219
220 Vect_line_distance(Points, ux, uy, uz, with_z, NULL, NULL, NULL,
221 &new_dist, NULL, NULL);
222 G_debug(3, " line = %d distance = %f", line, new_dist);
223
224 if (found && new_dist <= maxdist) {
225 Vect_list_append(found, line);
226 }
227
228 if ((++gotone == 1) || (new_dist <= cur_dist)) {
229 if (new_dist == cur_dist) {
230 /* TODO */
231 /* choice = dig_center_check (map->Line, choice, a, ux, uy); */
232 continue;
233 }
234
235 choice = line;
237 }
238 }
239
240 G_debug(3, "min distance found = %f", cur_dist);
241 if (cur_dist > maxdist)
242 choice = 0;
243
245
246 return (choice);
247}
248
249/*!
250 * \brief Find the nearest area
251 *
252 * \param Map vector map
253 * \param x,y point coordinates
254 *
255 * \return area number
256 * \return 0 if not found
257 */
258int Vect_find_area(struct Map_info *Map, double x, double y)
259{
260 int i, j, ret, area, isle;
261 struct bound_box box;
262 static struct boxlist *List = NULL;
263 static BOX_SIZE *size_list;
264 static int alloc_size_list = 0;
265 const struct Plus_head *Plus;
266 struct P_area *Area;
267
268 G_debug(3, "Vect_find_area() x = %f y = %f", x, y);
269
270 if (!List) {
272 alloc_size_list = 10;
273 size_list = G_malloc(alloc_size_list * sizeof(BOX_SIZE));
274 }
275
276 Plus = &(Map->plus);
277
278 /* select areas by box */
279 box.E = x;
280 box.W = x;
281 box.N = y;
282 box.S = y;
283 box.T = PORT_DOUBLE_MAX;
284 box.B = -PORT_DOUBLE_MAX;
286 G_debug(3, " %d areas selected by box", List->n_values);
287
288 /* sort areas by bbox size
289 * get the smallest area that contains the point
290 * using the bbox size is working because if 2 areas both contain
291 * the point, one of these areas must be inside the other area
292 * which means that the bbox of the outer area must be larger than
293 * the bbox of the inner area, and equal bbox sizes are not possible */
294
295 if (alloc_size_list < List->n_values) {
296 alloc_size_list = List->n_values;
297 size_list = G_realloc(size_list, alloc_size_list * sizeof(BOX_SIZE));
298 }
299
300 for (i = 0; i < List->n_values; i++) {
301 size_list[i].i = List->id[i];
302 box = List->box[i];
303 size_list[i].box = List->box[i];
304 size_list[i].size = (box.N - box.S) * (box.E - box.W);
305 }
306
307 if (List->n_values == 2) {
308 /* simple swap */
309 if (size_list[1].size < size_list[0].size) {
310 size_list[0].i = List->id[1];
311 size_list[1].i = List->id[0];
312 size_list[0].box = List->box[1];
313 size_list[1].box = List->box[0];
314 }
315 }
316 else if (List->n_values > 2)
317 qsort(size_list, List->n_values, sizeof(BOX_SIZE), sort_by_size);
318
319 for (i = 0; i < List->n_values; i++) {
320 area = size_list[i].i;
321 /* outer ring */
322 ret = Vect_point_in_area_outer_ring(x, y, Map, area, &size_list[i].box);
323
324 G_debug(3, " area = %d Vect_point_in_area_outer_ring() = %d", area,
325 ret);
326
327 if (ret >= 1) {
328 /* check if in islands */
329 Area = Plus->Area[area];
330 for (j = 0; j < Area->n_isles; j++) {
331 isle = Area->isles[j];
332 Vect_get_isle_box(Map, isle, &box);
333 ret = Vect_point_in_island(x, y, Map, isle, &box);
334
335 G_debug(3, " area = %d Vect_point_in_island() = %d", area,
336 ret);
337
338 if (ret >= 1) {
339 /* point is not in area
340 * point is also not in any inner area, those have
341 * been tested before (sorted list)
342 * -> area inside island could not be built */
343 return 0;
344 }
345 }
346 return (area);
347 }
348 }
349
350 return 0;
351}
352
353/*!
354 * \brief Find the nearest island
355 *
356 * \param Map vector map
357 * \param x,y points coordinates
358 *
359 * \return island number,
360 * \return 0 if not found
361 */
362int Vect_find_island(struct Map_info *Map, double x, double y)
363{
364 int i, ret, island, current, current_size, size;
365 static int first = 1;
366 struct bound_box box;
367 static struct boxlist *List;
368 static struct line_pnts *Points;
369
370 /* TODO: sync to Vect_find_area() */
371
372 G_debug(3, "Vect_find_island() x = %f y = %f", x, y);
373
374 if (first) {
376 Points = Vect_new_line_struct();
377 first = 0;
378 }
379
380 /* select islands by box */
381 box.E = x;
382 box.W = x;
383 box.N = y;
384 box.S = y;
385 box.T = PORT_DOUBLE_MAX;
386 box.B = -PORT_DOUBLE_MAX;
388 G_debug(3, " %d islands selected by box", List->n_values);
389
390 current_size = -1;
391 current = 0;
392 for (i = 0; i < List->n_values; i++) {
393 island = List->id[i];
394 ret = Vect_point_in_island(x, y, Map, island, &List->box[i]);
395
396 if (ret >= 1) { /* inside */
397 if (current > 0) { /* not first */
398 if (current_size == -1) { /* second */
400 Vect_get_isle_points(Map, current, Points);
401 current_size = G_area_of_polygon(Points->x, Points->y,
402 Points->n_points);
403 }
404
406 size =
407 G_area_of_polygon(Points->x, Points->y, Points->n_points);
408
409 if (size < current_size) {
410 current = island;
411 current_size = size;
412 }
413 }
414 else { /* first */
415 current = island;
416 }
417 }
418 }
419
420 return current;
421}
#define NULL
Definition ccmath.h:32
double G_area_of_polygon(const double *, const double *, int)
Area in square meters of polygon.
Definition gis/area.c:158
#define G_realloc(p, n)
Definition defs/gis.h:141
#define G_malloc(n)
Definition defs/gis.h:139
int G_debug(int, const char *,...) __attribute__((format(printf
int G_begin_polygon_area_calculations(void)
Begin polygon area calculations.
Definition gis/area.c:119
int Vect_get_node_coor(struct Map_info *, int, double *, double *, double *)
Get node coordinates.
Definition level_two.c:274
int Vect_select_isles_by_box(struct Map_info *, const struct bound_box *, struct boxlist *)
Select isles with bounding boxes by box.
Definition sindex.c:165
int Vect_point_in_island(double, double, struct Map_info *, int, struct bound_box *)
Determines if a point (X,Y) is inside an island.
Definition Vlib/poly.c:923
struct boxlist * Vect_new_boxlist(int)
Creates and initializes a struct boxlist.
void Vect_destroy_boxlist(struct boxlist *)
Frees all memory associated with a struct boxlist, including the struct itself.
int Vect_get_isle_points(struct Map_info *, int, struct line_pnts *)
Returns polygon array of points for given isle.
void Vect_destroy_list(struct ilist *)
Frees all memory associated with a struct ilist, including the struct itself.
int Vect_list_append(struct ilist *, int)
Append new item to the end of list if not yet present.
int Vect_read_line(struct Map_info *, struct line_pnts *, struct line_cats *, int)
Read vector feature (topological level required)
int Vect_line_distance(const struct line_pnts *, double, double, double, int, double *, double *, double *, double *, double *, double *)
Calculate distance of point to line.
Definition line.c:648
double Vect_points_distance(double, double, double, double, double, double, int)
Calculate distance of 2 points.
Definition line.c:866
int Vect_point_in_area_outer_ring(double, double, struct Map_info *, int, struct bound_box *)
Determines if a point (X,Y) is inside an area outer ring. Islands are not considered.
Definition Vlib/poly.c:854
int Vect_select_areas_by_box(struct Map_info *, const struct bound_box *, struct boxlist *)
Select areas with bounding boxes by box.
Definition sindex.c:121
struct ilist * Vect_new_list(void)
Creates and initializes a struct ilist.
int Vect_select_lines_by_box(struct Map_info *, const struct bound_box *, int, struct boxlist *)
Select lines with bounding boxes by box.
Definition sindex.c:32
int Vect_get_isle_box(struct Map_info *, int, struct bound_box *)
Get bounding box of isle.
int Vect_val_in_list(const struct ilist *, int)
Find a given item in the list.
struct line_pnts * Vect_new_line_struct(void)
Creates and initializes a line_pnts structure.
Definition line.c:45
int Vect_reset_list(struct ilist *)
Reset ilist structure.
int Vect_select_nodes_by_box(struct Map_info *, const struct bound_box *, struct ilist *)
Select nodes by box.
Definition sindex.c:187
#define PORT_DOUBLE_MAX
Limits for portable types.
Definition dig_defines.h:66
double b
Definition r_raster.c:39
Vector map info.
Area (topology) info.
plus_t n_isles
Number of islands inside.
plus_t * isles
1st generation interior islands
Basic topology-related info.
Bounding box.
Definition dig_structs.h:62
double W
West.
Definition dig_structs.h:78
double T
Top.
Definition dig_structs.h:82
double S
South.
Definition dig_structs.h:70
double N
North.
Definition dig_structs.h:66
double E
East.
Definition dig_structs.h:74
double B
Bottom.
Definition dig_structs.h:86
List of bounding boxes with id.
struct bound_box * box
Array of bounding boxes.
List of integers.
Definition gis.h:715
Feature geometry info - coordinates.
double * y
Array of Y coordinates.
double * x
Array of X coordinates.
int n_points
Number of points.
int Vect_find_area(struct Map_info *Map, double x, double y)
Find the nearest area.
int Vect_find_line(struct Map_info *map, double ux, double uy, double uz, int type, double maxdist, int with_z, int exclude)
Find the nearest line.
int Vect_find_line_list(struct Map_info *map, double ux, double uy, double uz, int type, double maxdist, int with_z, const struct ilist *exclude, struct ilist *found)
Find the nearest line(s).
int Vect_find_node(struct Map_info *Map, double ux, double uy, double uz, double maxdist, int with_z)
Find the nearest node.
int Vect_find_island(struct Map_info *Map, double x, double y)
Find the nearest island.
#define HUGE_VAL
#define x