GRASS Programmer's Manual
6.4.2(2012)
|
00001 00020 #include <stdlib.h> 00021 #include <grass/gis.h> 00022 #include <grass/Vect.h> 00023 00037 int 00038 Vect_select_lines_by_box(struct Map_info *Map, BOUND_BOX * Box, 00039 int type, struct ilist *list) 00040 { 00041 int i, line, nlines; 00042 struct Plus_head *plus; 00043 P_LINE *Line; 00044 static struct ilist *LocList = NULL; 00045 00046 G_debug(3, "Vect_select_lines_by_box()"); 00047 G_debug(3, " Box(N,S,E,W,T,B): %e, %e, %e, %e, %e, %e", Box->N, Box->S, 00048 Box->E, Box->W, Box->T, Box->B); 00049 plus = &(Map->plus); 00050 00051 if (!(plus->Spidx_built)) { 00052 G_debug(3, "Building spatial index."); 00053 Vect_build_sidx_from_topo(Map); 00054 } 00055 00056 list->n_values = 0; 00057 if (!LocList) 00058 LocList = Vect_new_list(); 00059 00060 nlines = dig_select_lines(plus, Box, LocList); 00061 G_debug(3, " %d lines selected (all types)", nlines); 00062 00063 /* Remove lines of not requested types */ 00064 for (i = 0; i < nlines; i++) { 00065 line = LocList->value[i]; 00066 if (plus->Line[line] == NULL) 00067 continue; /* Should not happen */ 00068 Line = plus->Line[line]; 00069 if (!(Line->type & type)) 00070 continue; 00071 dig_list_add(list, line); 00072 } 00073 00074 G_debug(3, " %d lines of requested type", list->n_values); 00075 00076 return list->n_values; 00077 } 00078 00091 int 00092 Vect_select_areas_by_box(struct Map_info *Map, BOUND_BOX * Box, 00093 struct ilist *list) 00094 { 00095 int i; 00096 00097 G_debug(3, "Vect_select_areas_by_box()"); 00098 G_debug(3, "Box(N,S,E,W,T,B): %e, %e, %e, %e, %e, %e", Box->N, Box->S, 00099 Box->E, Box->W, Box->T, Box->B); 00100 00101 if (!(Map->plus.Spidx_built)) { 00102 G_debug(3, "Building spatial index."); 00103 Vect_build_sidx_from_topo(Map); 00104 } 00105 00106 dig_select_areas(&(Map->plus), Box, list); 00107 G_debug(3, " %d areas selected", list->n_values); 00108 for (i = 0; i < list->n_values; i++) { 00109 G_debug(3, " area = %d pointer to area structure = %lx", 00110 list->value[i], 00111 (unsigned long)Map->plus.Area[list->value[i]]); 00112 00113 } 00114 return list->n_values; 00115 } 00116 00117 00130 int 00131 Vect_select_isles_by_box(struct Map_info *Map, BOUND_BOX * Box, 00132 struct ilist *list) 00133 { 00134 G_debug(3, "Vect_select_isles_by_box()"); 00135 G_debug(3, "Box(N,S,E,W,T,B): %e, %e, %e, %e, %e, %e", Box->N, Box->S, 00136 Box->E, Box->W, Box->T, Box->B); 00137 00138 if (!(Map->plus.Spidx_built)) { 00139 G_debug(3, "Building spatial index."); 00140 Vect_build_sidx_from_topo(Map); 00141 } 00142 00143 dig_select_isles(&(Map->plus), Box, list); 00144 G_debug(3, " %d isles selected", list->n_values); 00145 00146 return list->n_values; 00147 } 00148 00158 int 00159 Vect_select_nodes_by_box(struct Map_info *Map, BOUND_BOX * Box, 00160 struct ilist *list) 00161 { 00162 struct Plus_head *plus; 00163 00164 G_debug(3, "Vect_select_nodes_by_box()"); 00165 G_debug(3, "Box(N,S,E,W,T,B): %e, %e, %e, %e, %e, %e", Box->N, Box->S, 00166 Box->E, Box->W, Box->T, Box->B); 00167 00168 plus = &(Map->plus); 00169 00170 if (!(plus->Spidx_built)) { 00171 G_debug(3, "Building spatial index."); 00172 Vect_build_sidx_from_topo(Map); 00173 } 00174 00175 list->n_values = 0; 00176 00177 dig_select_nodes(plus, Box, list); 00178 G_debug(3, " %d nodes selected", list->n_values); 00179 00180 return list->n_values; 00181 } 00182 00197 int 00198 Vect_select_lines_by_polygon(struct Map_info *Map, struct line_pnts *Polygon, 00199 int nisles, struct line_pnts **Isles, int type, 00200 struct ilist *List) 00201 { 00202 int i; 00203 BOUND_BOX box; 00204 static struct line_pnts *LPoints = NULL; 00205 static struct ilist *LocList = NULL; 00206 00207 /* TODO: this function was not tested with isles */ 00208 G_debug(3, "Vect_select_lines_by_polygon() nisles = %d", nisles); 00209 00210 List->n_values = 0; 00211 if (!LPoints) 00212 LPoints = Vect_new_line_struct(); 00213 if (!LocList) 00214 LocList = Vect_new_list(); 00215 00216 /* Select first all lines by box */ 00217 dig_line_box(Polygon, &box); 00218 box.T = PORT_DOUBLE_MAX; 00219 box.B = -PORT_DOUBLE_MAX; 00220 Vect_select_lines_by_box(Map, &box, type, LocList); 00221 G_debug(3, " %d lines selected by box", LocList->n_values); 00222 00223 /* Check all lines if intersect the polygon */ 00224 for (i = 0; i < LocList->n_values; i++) { 00225 int j, line, intersect = 0; 00226 00227 line = LocList->value[i]; 00228 /* Read line points */ 00229 Vect_read_line(Map, LPoints, NULL, line); 00230 00231 /* Check if any of line vertices is within polygon */ 00232 for (j = 0; j < LPoints->n_points; j++) { 00233 if (Vect_point_in_poly(LPoints->x[j], LPoints->y[j], Polygon) >= 1) { /* inside polygon */ 00234 int k, inisle = 0; 00235 00236 for (k = 0; k < nisles; k++) { 00237 if (Vect_point_in_poly(LPoints->x[j], LPoints->y[j], Isles[k]) >= 1) { /* in isle */ 00238 inisle = 1; 00239 break; 00240 } 00241 } 00242 00243 if (!inisle) { /* inside polygon, outside isles -> select */ 00244 intersect = 1; 00245 break; 00246 } 00247 } 00248 } 00249 if (intersect) { 00250 dig_list_add(List, line); 00251 continue; 00252 } 00253 00254 /* Check intersections of the line with area/isles boundary */ 00255 /* Outer boundary */ 00256 if (Vect_line_check_intersection(LPoints, Polygon, 0)) { 00257 dig_list_add(List, line); 00258 continue; 00259 } 00260 00261 /* Islands */ 00262 for (j = 0; j < nisles; j++) { 00263 if (Vect_line_check_intersection(LPoints, Isles[j], 0)) { 00264 intersect = 1; 00265 break; 00266 } 00267 } 00268 if (intersect) { 00269 dig_list_add(List, line); 00270 } 00271 } 00272 00273 G_debug(4, " %d lines selected by polygon", List->n_values); 00274 00275 return List->n_values; 00276 } 00277 00278 00295 int 00296 Vect_select_areas_by_polygon(struct Map_info *Map, struct line_pnts *Polygon, 00297 int nisles, struct line_pnts **Isles, 00298 struct ilist *List) 00299 { 00300 int i, area; 00301 static struct ilist *BoundList = NULL; 00302 00303 /* TODO: this function was not tested with isles */ 00304 G_debug(3, "Vect_select_areas_by_polygon() nisles = %d", nisles); 00305 00306 List->n_values = 0; 00307 if (!BoundList) 00308 BoundList = Vect_new_list(); 00309 00310 /* Select boundaries by polygon */ 00311 Vect_select_lines_by_polygon(Map, Polygon, nisles, Isles, GV_BOUNDARY, 00312 BoundList); 00313 00314 /* Add areas on left/right side of selected boundaries */ 00315 for (i = 0; i < BoundList->n_values; i++) { 00316 int line, left, right; 00317 00318 line = BoundList->value[i]; 00319 00320 Vect_get_line_areas(Map, line, &left, &right); 00321 G_debug(4, "boundary = %d left = %d right = %d", line, left, right); 00322 00323 if (left > 0) { 00324 dig_list_add(List, left); 00325 } 00326 else if (left < 0) { /* island */ 00327 area = Vect_get_isle_area(Map, abs(left)); 00328 G_debug(4, " left island -> area = %d", area); 00329 if (area > 0) 00330 dig_list_add(List, area); 00331 } 00332 00333 if (right > 0) { 00334 dig_list_add(List, right); 00335 } 00336 else if (right < 0) { /* island */ 00337 area = Vect_get_isle_area(Map, abs(right)); 00338 G_debug(4, " right island -> area = %d", area); 00339 if (area > 0) 00340 dig_list_add(List, area); 00341 } 00342 } 00343 00344 /* But the Polygon may be completely inside the area (only one), in that case 00345 * we find the area by one polygon point and add it to the list */ 00346 area = Vect_find_area(Map, Polygon->x[0], Polygon->y[0]); 00347 if (area > 0) 00348 dig_list_add(List, area); 00349 00350 G_debug(3, " %d areas selected by polygon", List->n_values); 00351 00352 return List->n_values; 00353 }