remove bounding boxes from vector topology (Radim's TODO): update Vlib

git-svn-id: https://svn.osgeo.org/grass/grass/trunk@46901 15284696-431f-4ddb-bdfa-cd5b030d7da7

lib/vector/Vlib/area.c

@@ -345,20 +345,23 @@ int Vect_point_in_area(const struct Map_info *Map, int area, double x, double y)
     const struct Plus_head *Plus;
     struct P_area *Area;
     int poly;
+    struct bound_box box;
 
     Plus = &(Map->plus);
     Area = Plus->Area[area];
     if (Area == NULL)
 	return 0;
 
-    poly = Vect_point_in_area_outer_ring(x, y, Map, area);
+    Vect_get_area_box(Map, area, &box);
+    poly = Vect_point_in_area_outer_ring(x, y, Map, area, box);
     if (poly == 0)
 	return 0;		/* includes area boundary (poly == 2), OK? */
 
     /* check if in islands */
     for (i = 0; i < Area->n_isles; i++) {
 	isle = Area->isles[i];
-	poly = Vect_point_in_island(x, y, Map, isle);
+	Vect_get_isle_box(Map, isle, &box);
+	poly = Vect_point_in_island(x, y, Map, isle, box);
 	if (poly >= 1)
 	    return 0;		/* excludes island boundary (poly == 2), OK? */
     }

lib/vector/Vlib/box.c

@@ -17,6 +17,7 @@
 
 #include <stdlib.h>
 #include <grass/vector.h>
+#include <grass/glocale.h>
 
 /*!
    \brief Tests for point in box
@@ -206,6 +207,8 @@ int Vect_get_line_box(const struct Map_info *Map, int line, struct bound_box * B
 {
     const struct Plus_head *Plus;
     struct P_line *Line;
+    static struct line_pnts *Points = NULL;
+    static struct boxlist *list = NULL;
 
     Plus = &(Map->plus);
     Line = Plus->Line[line];
@@ -220,12 +223,56 @@ int Vect_get_line_box(const struct Map_info *Map, int line, struct bound_box * B
 	return 0;
     }
     else {
-	Box->N = Line->N;
-	Box->S = Line->S;
-	Box->E = Line->E;
-	Box->W = Line->W;
-	Box->T = Line->T;
-	Box->B = Line->B;
+	int type = Line->type;
+	
+	/* GV_POINTS: read line */
+	if (type & GV_POINTS) {
+	    if (Points == NULL)
+		Points = Vect_new_line_struct();
+
+	    Vect_read_line(Map, Points, NULL, line);
+	    dig_line_box(Points, Box);
+	}
+	/* all other: retrieve box from spatial index */
+	else {
+	    int node = 0;
+	    struct bound_box bbox;
+
+	    if (type == GV_LINE) {
+		struct P_topo_l *topo = (struct P_topo_l *)Line->topo;
+
+		node = topo->N1;
+	    }
+	    else if (type == GV_BOUNDARY) {
+		struct P_topo_b *topo = (struct P_topo_b *)Line->topo;
+
+		node = topo->N1;
+	    }
+	    
+	    if (list == NULL) {
+		list = Vect_new_boxlist();
+	    }
+	    Vect_reset_boxlist(list);
+	    
+	    bbox.N = Plus->Node[node]->y;
+	    bbox.S = Plus->Node[node]->y;
+	    bbox.E = Plus->Node[node]->x;
+	    bbox.W = Plus->Node[node]->x;
+	    bbox.T = Plus->Node[node]->z;
+	    bbox.B = Plus->Node[node]->z;
+
+	    dig_boxlist_add(list, line, bbox);
+	    
+	    if (dig_find_line_box(Plus, list) == 0)
+		G_fatal_error(_("Could not find line box"));
+
+	    Box->N = list->box[0].N;
+	    Box->S = list->box[0].S;
+	    Box->E = list->box[0].E;
+	    Box->W = list->box[0].W;
+	    Box->T = list->box[0].T;
+	    Box->B = list->box[0].B;
+	}
     }
 
     return 1;
@@ -245,6 +292,11 @@ int Vect_get_area_box(const struct Map_info *Map, int area, struct bound_box * B
 {
     const struct Plus_head *Plus;
     struct P_area *Area;
+    struct P_line *Line;
+    struct P_node *Node;
+    static struct boxlist *list = NULL;
+    struct bound_box bbox;
+    struct P_topo_b *topo;
 
     Plus = &(Map->plus);
     Area = Plus->Area[area];
@@ -256,15 +308,38 @@ int Vect_get_area_box(const struct Map_info *Map, int area, struct bound_box * B
 	Box->W = 0;
 	Box->T = 0;
 	Box->B = 0;
+
 	return 0;
     }
     else {
-	Box->N = Area->N;
-	Box->S = Area->S;
-	Box->E = Area->E;
-	Box->W = Area->W;
-	Box->T = Area->T;
-	Box->B = Area->B;
+
+	Line = Plus->Line[abs(Area->lines[0])];
+	topo = (struct P_topo_b *)Line->topo;
+	Node = Plus->Node[topo->N1];
+
+	if (list == NULL) {
+	    list = Vect_new_boxlist();
+	}
+	Vect_reset_boxlist(list);
+	
+	bbox.N = Node->y;
+	bbox.S = Node->y;
+	bbox.E = Node->x;
+	bbox.W = Node->x;
+	bbox.T = Node->z;
+	bbox.B = Node->z;
+
+	dig_boxlist_add(list, area, bbox);
+	
+	if (dig_find_area_box(Plus, list) == 0)
+	    G_fatal_error(_("Could not find area box"));
+
+	Box->N = list->box[0].N;
+	Box->S = list->box[0].S;
+	Box->E = list->box[0].E;
+	Box->W = list->box[0].W;
+	Box->T = list->box[0].T;
+	Box->B = list->box[0].B;
     }
 
     return 1;
@@ -284,6 +359,11 @@ int Vect_get_isle_box(const struct Map_info *Map, int isle, struct bound_box * B
 {
     const struct Plus_head *Plus;
     struct P_isle *Isle;
+    struct P_line *Line;
+    struct P_node *Node;
+    static struct boxlist *list = NULL;
+    struct bound_box bbox;
+    struct P_topo_b *topo;
 
     Plus = &(Map->plus);
     Isle = Plus->Isle[isle];
@@ -298,12 +378,34 @@ int Vect_get_isle_box(const struct Map_info *Map, int isle, struct bound_box * B
 	return 0;
     }
     else {
-	Box->N = Isle->N;
-	Box->S = Isle->S;
-	Box->E = Isle->E;
-	Box->W = Isle->W;
-	Box->T = Isle->T;
-	Box->B = Isle->B;
+
+	Line = Plus->Line[abs(Isle->lines[0])];
+	topo = (struct P_topo_b *)Line->topo;
+	Node = Plus->Node[topo->N1];
+
+	if (list == NULL) {
+	    list = Vect_new_boxlist();
+	}
+	Vect_reset_boxlist(list);
+	
+	bbox.N = Node->y;
+	bbox.S = Node->y;
+	bbox.E = Node->x;
+	bbox.W = Node->x;
+	bbox.T = Node->z;
+	bbox.B = Node->z;
+
+	dig_boxlist_add(list, isle, bbox);
+	
+	if (dig_find_isle_box(Plus, list) == 0)
+	    G_fatal_error(_("Could not find area box"));
+
+	Box->N = list->box[0].N;
+	Box->S = list->box[0].S;
+	Box->E = list->box[0].E;
+	Box->W = list->box[0].W;
+	Box->T = list->box[0].T;
+	Box->B = list->box[0].B;
     }
 
     return 1;

lib/vector/Vlib/break_lines.c

@@ -133,7 +133,7 @@ Vect_break_lines_list(struct Map_info *Map, struct ilist *List_break,
 	if (!(atype & type))
 	    continue;
 
-	Vect_get_line_box(Map, aline, &ABox);
+	Vect_line_box(APoints, &ABox);
 
 	/* Find which sides of the box are touched by intermediate (non-end) points of line */
 	if (!is3d) {
@@ -175,11 +175,13 @@ Vect_break_lines_list(struct Map_info *Map, struct ilist *List_break,
 	    }
 	    G_debug(3, "  j = %d bline = %d", j, bline);
 
+	    btype = Vect_read_line(Map, BPoints, BCats, bline);
+
 	    /* Check if thouch by end node only */
 	    if (!is3d) {
 		Vect_get_line_nodes(Map, aline, &anode1, &anode2);
 		Vect_get_line_nodes(Map, bline, &bnode1, &bnode2);
-		Vect_get_line_box(Map, bline, &BBox);
+		Vect_line_box(BPoints, &BBox);
 
 		if (anode1 == bnode1 || anode1 == bnode2)
 		    node = anode1;
@@ -190,22 +192,23 @@ Vect_break_lines_list(struct Map_info *Map, struct ilist *List_break,
 
 		if (node) {
 		    Vect_get_node_coor(Map, node, &nodex, &nodey, NULL);
-		    if ((node == anode1 && nodey == ABox.N && !touch1_n &&
-			 nodey == BBox.S) || (node == anode2 &&
-					      nodey == ABox.N && !touch2_n &&
-					      nodey == BBox.S) ||
-			(node == anode1 && nodey == ABox.S && !touch1_s &&
-			 nodey == BBox.N) || (node == anode2 &&
-					      nodey == ABox.S && !touch2_s &&
-					      nodey == BBox.N) ||
-			(node == anode1 && nodex == ABox.E && !touch1_e &&
-			 nodex == BBox.W) || (node == anode2 &&
-					      nodex == ABox.E && !touch2_e &&
-					      nodex == BBox.W) ||
-			(node == anode1 && nodex == ABox.W && !touch1_w &&
-			 nodex == BBox.E) || (node == anode2 &&
-					      nodex == ABox.W && !touch2_w &&
-					      nodex == BBox.E)) {
+		    if ((node == anode1 && nodey == ABox.N &&
+		         !touch1_n && nodey == BBox.S) ||
+		        (node == anode2 && nodey == ABox.N &&
+			 !touch2_n && nodey == BBox.S) ||
+			(node == anode1 && nodey == ABox.S &&
+			 !touch1_s && nodey == BBox.N) ||
+			(node == anode2 && nodey == ABox.S &&
+			 !touch2_s && nodey == BBox.N) ||
+			(node == anode1 && nodex == ABox.E &&
+			 !touch1_e && nodex == BBox.W) ||
+			(node == anode2 && nodex == ABox.E &&
+			 !touch2_e && nodex == BBox.W) ||
+			(node == anode1 && nodex == ABox.W &&
+			 !touch1_w && nodex == BBox.E) ||
+			(node == anode2 && nodex == ABox.W &&
+			 !touch2_w && nodex == BBox.E)) {
+
 			G_debug(3,
 				"lines %d and %d touching by end nodes only -> no intersection",
 				aline, bline);
@@ -214,8 +217,6 @@ Vect_break_lines_list(struct Map_info *Map, struct ilist *List_break,
 		}
 	    }
 
-	    btype = Vect_read_line(Map, BPoints, BCats, bline);
-
 	    AXLines = NULL;
 	    BXLines = NULL;
 	    Vect_line_intersection(APoints, BPoints, &AXLines, &BXLines,
@@ -374,6 +375,12 @@ Vect_break_lines_list(struct Map_info *Map, struct ilist *List_break,
 
     G_verbose_message(_("Intersections: %d"), nbreaks);
 
+    Vect_destroy_line_struct(APoints);
+    Vect_destroy_line_struct(BPoints);
+    Vect_destroy_line_struct(Points);
+    Vect_destroy_cats_struct(ACats);
+    Vect_destroy_cats_struct(BCats);
+    Vect_destroy_cats_struct(Cats);
     Vect_destroy_list(List);
 
     return nbreaks;

lib/vector/Vlib/break_polygons.c

@@ -73,7 +73,7 @@ typedef struct
 static int fpoint;
 
 /* Function called from RTreeSearch for point found */
-void srch(int id, int *arg)
+void srch(int id, struct Rect rect, int *arg)
 {
     fpoint = id;
 }

lib/vector/Vlib/build.c

@@ -180,16 +180,21 @@ int Vect_build_partial(struct Map_info *Map, int build)
 	    Line = Plus->Line[line];
 	    if (!Line)
 		continue;
-	    if (Line->type == GV_BOUNDARY &&
-		(Line->left == 0 || Line->right == 0)) {
-		G_debug(3, "line = %d left = %d right = %d", line, Line->left,
-			Line->right);
-		err_boundaries++;
+	    if (Line->type == GV_BOUNDARY) {
+		struct P_topo_b *topo = (struct P_topo_b *)Line->topo;
+
+		if (topo->left == 0 || topo->right == 0) {
+		    G_debug(3, "line = %d left = %d right = %d", line, 
+			    topo->left, topo->right);
+		    err_boundaries++;
+		}
 	    }
 	    if (Line->type == GV_CENTROID) {
-		if (Line->left == 0)
+		struct P_topo_c *topo = (struct P_topo_c *)Line->topo;
+
+		if (topo->area == 0)
 		    err_centr_out++;
-		else if (Line->left < 0)
+		else if (topo->area < 0)
 		    err_centr_dupl++;
 	    }
 	}
@@ -323,12 +328,45 @@ int Vect_topo_dump(const struct Map_info *Map, FILE *out)
 	    continue;
 	}
 	Line = plus->Line[i];
-	fprintf(out, "line = %d, type = %d, offset = %lu n1 = %d, n2 = %d, "
-		"left/area = %d, right = %d\n",
-		i, Line->type, (unsigned long)Line->offset, Line->N1,
-		Line->N2, Line->left, Line->right);
-	fprintf(out, "N,S,E,W,T,B: %f, %f, %f, %f, %f, %f\n", Line->N,
-		Line->S, Line->E, Line->W, Line->T, Line->B);
+	if (Line->type == GV_POINT) {
+	    fprintf(out, "line = %d, type = %d, offset = %lu\n",
+		    i, Line->type, (unsigned long)Line->offset);
+	}
+	else if (Line->type == GV_CENTROID) {
+	    struct P_topo_c *topo = (struct P_topo_c *)Line->topo;
+
+	    fprintf(out, "line = %d, type = %d, offset = %lu, area = %d\n",
+		    i, Line->type, (unsigned long)Line->offset, topo->area);
+	}
+	else if (Line->type == GV_LINE) {
+	    struct P_topo_l *topo = (struct P_topo_l *)Line->topo;
+
+	    fprintf(out, "line = %d, type = %d, offset = %lu, n1 = %d, n2 = %d\n",
+		    i, Line->type, (unsigned long)Line->offset, topo->N1,
+		    topo->N2);
+	}
+	else if (Line->type == GV_BOUNDARY) {
+	    struct P_topo_b *topo = (struct P_topo_b *)Line->topo;
+
+	    fprintf(out, "line = %d, type = %d, offset = %lu, n1 = %d, n2 = %d, "
+		    "left = %d, right = %d\n",
+		    i, Line->type, (unsigned long)Line->offset, topo->N1,
+		    topo->N2, topo->left, topo->right);
+	}
+	else if (Line->type == GV_FACE) {
+	    struct P_topo_f *topo = (struct P_topo_f *)Line->topo;
+
+	    fprintf(out, "line = %d, type = %d, offset = %lu, e1 = %d, e2 = %d, "
+		    "e3 = %d, left = %d, right = %d\n",
+		    i, Line->type, (unsigned long)Line->offset, topo->E[0],
+		    topo->E[1], topo->E[2], topo->left, topo->right);
+	}
+	else if (Line->type == GV_KERNEL) {
+	    struct P_topo_k *topo = (struct P_topo_k *)Line->topo;
+
+	    fprintf(out, "line = %d, type = %d, offset = %lu, volume = %d",
+		    i, Line->type, (unsigned long)Line->offset, topo->volume);
+	}
     }
 
     /* areas */
@@ -342,9 +380,6 @@ int Vect_topo_dump(const struct Map_info *Map, FILE *out)
 	fprintf(out, "area = %d, n_lines = %d, n_isles = %d centroid = %d\n",
 		i, Area->n_lines, Area->n_isles, Area->centroid);
 
-	fprintf(out, "N,S,E,W,T,B: %f, %f, %f, %f, %f, %f\n", Area->N,
-		Area->S, Area->E, Area->W, Area->T, Area->B);
-
 	for (j = 0; j < Area->n_lines; j++) {
 	    line = Area->lines[j];
 	    Line = plus->Line[abs(line)];
@@ -367,9 +402,6 @@ int Vect_topo_dump(const struct Map_info *Map, FILE *out)
 	fprintf(out, "isle = %d, n_lines = %d area = %d\n", i, Isle->n_lines,
 		Isle->area);
 
-	fprintf(out, "N,S,E,W,T,B: %f, %f, %f, %f, %f, %f\n", Isle->N,
-		Isle->S, Isle->E, Isle->W, Isle->T, Isle->B);
-
 	for (j = 0; j < Isle->n_lines; j++) {
 	    line = Isle->lines[j];
 	    Line = plus->Line[abs(line)];
@@ -413,95 +445,12 @@ int Vect_build_sidx(struct Map_info *Map)
  */
 int Vect_build_sidx_from_topo(struct Map_info *Map)
 {
-    int i, total, done;
-    struct Plus_head *plus;
-    struct bound_box box;
-    struct P_line *Line;
-    struct P_node *Node;
-    struct P_area *Area;
-    struct P_isle *Isle;
 
     G_debug(3, "Vect_build_sidx_from_topo()");
 
-    plus = &(Map->plus);
-
-    Vect_open_sidx(Map, 2);
-
-    total = plus->n_nodes + plus->n_lines + plus->n_areas + plus->n_isles;
-
-    /* Nodes */
-    for (i = 1; i <= plus->n_nodes; i++) {
-	G_percent(i, total, 3);
-
-	Node = plus->Node[i];
-	if (!Node)
-	    G_fatal_error(_("BUG (Vect_build_sidx_from_topo): node does not exist"));
-
-	dig_spidx_add_node(plus, i, Node->x, Node->y, Node->z);
-    }
-
-    /* Lines */
-    done = plus->n_nodes;
-    for (i = 1; i <= plus->n_lines; i++) {
-	G_percent(done + i, total, 3);
-
-	Line = plus->Line[i];
-	if (!Line)
-	    G_fatal_error(_("BUG (Vect_build_sidx_from_topo): line does not exist"));
-
-	box.N = Line->N;
-	box.S = Line->S;
-	box.E = Line->E;
-	box.W = Line->W;
-	box.T = Line->T;
-	box.B = Line->B;
-
-	dig_spidx_add_line(plus, i, &box);
-    }
-
-    /* Areas */
-    done += plus->n_lines;
-    for (i = 1; i <= plus->n_areas; i++) {
-	G_percent(done + i, total, 3);
-
-	Area = plus->Area[i];
-	if (!Area)
-	    G_fatal_error(_("BUG (Vect_build_sidx_from_topo): area does not exist"));
-
-	box.N = Area->N;
-	box.S = Area->S;
-	box.E = Area->E;
-	box.W = Area->W;
-	box.T = Area->T;
-	box.B = Area->B;
+    G_warning(_("Vect_build_sidx_from_topo() is no longer supported"));
 
-	dig_spidx_add_area(plus, i, &box);
-    }
-
-    /* Isles */
-    done += plus->n_areas;
-    for (i = 1; i <= plus->n_isles; i++) {
-	G_percent(done + i, total, 3);
-
-	Isle = plus->Isle[i];
-	if (!Isle)
-	    G_fatal_error(_("BUG (Vect_build_sidx_from_topo): isle does not exist"));
-
-	box.N = Isle->N;
-	box.S = Isle->S;
-	box.E = Isle->E;
-	box.W = Isle->W;
-	box.T = Isle->T;
-	box.B = Isle->B;
-
-	dig_spidx_add_isle(plus, i, &box);
-    }
-
-    Map->plus.Spidx_built = 1;
-
-    G_debug(3, "Spatial index was built");
-
-    return 0;
+    return 1;
 }
 
 /*!

lib/vector/Vlib/build_nat.c

@@ -40,6 +40,7 @@ int Vect_build_line_area(struct Map_info *Map, int iline, int side)
     struct Plus_head *plus;
     struct P_line *BLine;
     static struct line_pnts *Points, *APoints;
+    struct bound_box box;
     plus_t *lines;
     double area_size;
 
@@ -81,6 +82,7 @@ int Vect_build_line_area(struct Map_info *Map, int iline, int side)
 	Vect_append_points(APoints, Points, direction);
 	APoints->n_points--;	/* skip last point, avoids duplicates */
     }
+    dig_line_box(APoints, &box);
     APoints->n_points++;	/* close polygon */
 
     dig_find_area_poly(APoints, &area_size);
@@ -92,7 +94,7 @@ int Vect_build_line_area(struct Map_info *Map, int iline, int side)
 
     if (area_size > 0) {	/* CW: area */
 	/* add area structure to plus */
-	area = dig_add_area(plus, n_lines, lines);
+	area = dig_add_area(plus, n_lines, lines, &box);
 	if (area == -1) {	/* error */
 	    Vect_close(Map);
 	    G_fatal_error(_("Unable to add area (map closed, topo saved)"));
@@ -101,7 +103,7 @@ int Vect_build_line_area(struct Map_info *Map, int iline, int side)
 	return area;
     }
     else if (area_size < 0) {	/* CCW: island */
-	isle = dig_add_isle(plus, n_lines, lines);
+	isle = dig_add_isle(plus, n_lines, lines, &box);
 	if (isle == -1) {	/* error */
 	    Vect_close(Map);
 	    G_fatal_error(_("Unable to add isle (map closed, topo saved)"));
@@ -130,16 +132,17 @@ int Vect_build_line_area(struct Map_info *Map, int iline, int side)
  */
 int Vect_isle_find_area(struct Map_info *Map, int isle)
 {
-    int j, line, node, sel_area, first, area, poly;
+    int j, line, sel_area, first, area, poly;
     static int first_call = 1;
     const struct Plus_head *plus;
     struct P_line *Line;
     struct P_node *Node;
     struct P_isle *Isle;
     struct P_area *Area;
+    struct P_topo_b *topo;
     double size, cur_size;
     struct bound_box box, abox;
-    static struct ilist *List;
+    static struct boxlist *List;
     static struct line_pnts *APoints;
 
     /* Note: We should check all isle points (at least) because if topology is not clean
@@ -155,7 +158,7 @@ int Vect_isle_find_area(struct Map_info *Map, int isle)
     }
 
     if (first_call) {
-	List = Vect_new_list();
+	List = Vect_new_boxlist();
 	APoints = Vect_new_line_struct();
 	first_call = 0;
     }
@@ -163,8 +166,8 @@ int Vect_isle_find_area(struct Map_info *Map, int isle)
     Isle = plus->Isle[isle];
     line = abs(Isle->lines[0]);
     Line = plus->Line[line];
-    node = Line->N1;
-    Node = plus->Node[node];
+    topo = (struct P_topo_b *)Line->topo;
+    Node = plus->Node[topo->N1];
 
     /* select areas by box */
     box.E = Node->x;
@@ -173,7 +176,7 @@ int Vect_isle_find_area(struct Map_info *Map, int isle)
     box.S = Node->y;
     box.T = PORT_DOUBLE_MAX;
     box.B = -PORT_DOUBLE_MAX;
-    Vect_select_areas_by_box(Map, &box, List);
+    Vect_select_areas_by_box_with_box(Map, &box, List);
     G_debug(3, "%d areas overlap island boundary point", List->n_values);
 
     sel_area = 0;
@@ -181,7 +184,7 @@ int Vect_isle_find_area(struct Map_info *Map, int isle)
     first = 1;
     Vect_get_isle_box(Map, isle, &box);
     for (j = 0; j < List->n_values; j++) {
-	area = List->value[j];
+	area = List->id[j];
 	G_debug(3, "area = %d", area);
 
 	Area = plus->Area[area];
@@ -198,14 +201,16 @@ int Vect_isle_find_area(struct Map_info *Map, int isle)
 	 * box all overlapping areas selects all areas with box overlapping first node. 
 	 * Then reading coordinates for all those areas would take a long time -> check first 
 	 * if isle's box is completely within area box */
-	Vect_get_area_box(Map, area, &abox);
+
+	abox = List->box[j];
+
 	if (box.E > abox.E || box.W < abox.W || box.N > abox.N ||
 	    box.S < abox.S) {
 	    G_debug(3, "  isle not completely inside area box");
 	    continue;
 	}
 
-	poly = Vect_point_in_area_outer_ring(Node->x, Node->y, Map, area);
+	poly = Vect_point_in_area_outer_ring(Node->x, Node->y, Map, area, abox);
 	G_debug(3, "  poly = %d", poly);
 
 	if (poly == 1) {	/* point in area, but node is not part of area inside isle (would be poly == 2) */
@@ -330,7 +335,7 @@ int Vect_attach_isles(struct Map_info *Map, const struct bound_box * box)
 }
 
 /*!
-   \brief (Re)Attach centroids to areas in given bouding box
+   \brief (Re)Attach centroids to areas in given bounding box
 
    \param Map_info vector map
    \param box bounding box
@@ -342,8 +347,10 @@ int Vect_attach_centroids(struct Map_info *Map, const struct bound_box * box)
     int i, sel_area, centr;
     static int first = 1;
     static struct ilist *List;
+    static struct line_pnts *Points;
     struct P_area *Area;
     struct P_line *Line;
+    struct P_topo_c *topo;
     struct Plus_head *plus;
 
     G_debug(3, "Vect_attach_centroids ()");
@@ -352,6 +359,7 @@ int Vect_attach_centroids(struct Map_info *Map, const struct bound_box * box)
 
     if (first) {
 	List = Vect_new_list();
+	Points = Vect_new_line_struct();
 	first = 0;
     }
 
@@ -388,25 +396,27 @@ int Vect_attach_centroids(struct Map_info *Map, const struct bound_box * box)
 
 	centr = List->value[i];
 	Line = plus->Line[centr];
+	topo = (struct P_topo_c *)Line->topo;
 
 	/* only attach unregistered and duplicate centroids because 
 	 * 1) all properly attached centroids are properly attached, really! Don't touch.
 	 * 2) Vect_find_area() below does not always return the correct area
 	 * 3) it's faster
 	 */
-	if (Line->left > 0)
+	if (topo->area > 0)
 	    continue;
 
-	orig_area = Line->left;
+	orig_area = topo->area;
 
-	sel_area = Vect_find_area(Map, Line->E, Line->N);
+	Vect_read_line(Map, Points, NULL, centr);
+	sel_area = Vect_find_area(Map, Points->x[0], Points->y[0]);
 	G_debug(3, "  centroid %d is in area %d", centr, sel_area);
 	if (sel_area > 0) {
 	    Area = plus->Area[sel_area];
 	    if (Area->centroid == 0) {	/* first centroid */
 		G_debug(3, "  first centroid -> attach to area");
 		Area->centroid = centr;
-		Line->left = sel_area;
+		topo->area = sel_area;
 
 		if (sel_area != orig_area && plus->do_uplist)
 		    dig_line_add_updated(plus, centr);
@@ -415,7 +425,7 @@ int Vect_attach_centroids(struct Map_info *Map, const struct bound_box * box)
 		/* Note: it cannot happen that Area->centroid == centr, because the centroid
 		 * was not registered or a duplicate */
 		G_debug(3, "  duplicate centroid -> do not attach to area");
-		Line->left = -sel_area;
+		topo->area = -sel_area;
 
 		if (-sel_area != orig_area && plus->do_uplist)
 		    dig_line_add_updated(plus, centr);
@@ -438,15 +448,16 @@ int Vect_attach_centroids(struct Map_info *Map, const struct bound_box * box)
 int Vect_build_nat(struct Map_info *Map, int build)
 {
     struct Plus_head *plus;
-    int i, s, type, lineid;
+    int i, s, type, line;
     off_t offset;
-    int side, line, area;
+    int side, area;
     struct line_pnts *Points, *APoints;
     struct line_cats *Cats;
     struct P_line *Line;
     struct P_area *Area;
     struct bound_box box;
     struct ilist *List;
+    int print_counter = G_verbose() > G_verbose_min();
 
     G_debug(3, "Vect_build_nat() build = %d", build);
 
@@ -465,8 +476,10 @@ int Vect_build_nat(struct Map_info *Map, int build)
 
 	    for (line = 1; line <= nlines; line++) {
 		Line = plus->Line[line];
-		if (Line && Line->type == GV_CENTROID)
-		    Line->left = 0;
+		if (Line && Line->type == GV_CENTROID) {
+		    struct P_topo_c *topo = (struct P_topo_c *)Line->topo;
+		    topo->area = 0;
+		}
 	    }
 	    dig_free_plus_areas(plus);
 	    dig_spidx_free_areas(plus);
@@ -482,8 +495,9 @@ int Vect_build_nat(struct Map_info *Map, int build)
 	    for (line = 1; line <= nlines; line++) {
 		Line = plus->Line[line];
 		if (Line && Line->type == GV_BOUNDARY) {
-		    Line->left = 0;
-		    Line->right = 0;
+		    struct P_topo_b *topo = (struct P_topo_b *)Line->topo;
+		    topo->left = 0;
+		    topo->right = 0;
 		}
 	    }
 	    dig_free_plus_areas(plus);
@@ -508,7 +522,7 @@ int Vect_build_nat(struct Map_info *Map, int build)
     List = Vect_new_list();
 
     if (plus->built < GV_BUILD_BASE) {
-	int npoints, format;
+	register int npoints, format, c;
 
 	format = G_info_format();
 
@@ -521,7 +535,7 @@ int Vect_build_nat(struct Map_info *Map, int build)
 	/* register lines, create nodes */
 	Vect_rewind(Map);
 	G_message(_("Registering primitives..."));
-	i = 1;
+	i = 0;
 	npoints = 0;
 	while (1) {
 	    /* register line */
@@ -541,36 +555,34 @@ int Vect_build_nat(struct Map_info *Map, int build)
 	    offset = Map->head.last_offset;
 
 	    G_debug(3, "Register line: offset = %lu", (unsigned long)offset);
-	    lineid = dig_add_line(plus, type, Points, offset);
 	    dig_line_box(Points, &box);
-	    if (lineid == 1)
+	    line = dig_add_line(plus, type, Points, &box, offset);
+	    if (line == 1)
 		Vect_box_copy(&(plus->box), &box);
 	    else
 		Vect_box_extend(&(plus->box), &box);
 
 	    /* Add all categories to category index */
 	    if (build == GV_BUILD_ALL) {
-		int c;
-
 		for (c = 0; c < Cats->n_cats; c++) {
 		    dig_cidx_add_cat(plus, Cats->field[c], Cats->cat[c],
-				     lineid, type);
+				     line, type);
 		}
 		if (Cats->n_cats == 0)	/* add field 0, cat 0 */
-		    dig_cidx_add_cat(plus, 0, 0, lineid, type);
+		    dig_cidx_add_cat(plus, 0, 0, line, type);
 	    }
 
-	    if (G_verbose() > G_verbose_min() && i % 1000 == 0) {
+	    i++;
+	    if (i == 10000 && print_counter) {
 		if (format == G_INFO_FORMAT_PLAIN)
-		    fprintf(stderr, "%d..", i);
+		    fprintf(stderr, "%d..", plus->n_lines);
 		else
-		    fprintf(stderr, "%10d\b\b\b\b\b\b\b\b\b\b", i);
+		    fprintf(stderr, "%10d\b\b\b\b\b\b\b\b\b\b", plus->n_lines);
+		i = 0; 
 	    }
-
-	    i++;
 	}
 
-	if ((G_verbose() > G_verbose_min()) && format != G_INFO_FORMAT_PLAIN)
+	if ((print_counter) && format != G_INFO_FORMAT_PLAIN)
 	    fprintf(stderr, "\r");
 
 	G_message(_("%d primitives registered"), plus->n_lines);
@@ -586,14 +598,14 @@ int Vect_build_nat(struct Map_info *Map, int build)
 	/* Build areas */
 	/* Go through all bundaries and try to build area for both sides */
 	G_important_message(_("Building areas..."));
-	for (i = 1; i <= plus->n_lines; i++) {
-	    G_percent(i, plus->n_lines, 1);
+	for (line = 1; line <= plus->n_lines; line++) {
+	    G_percent(line, plus->n_lines, 1);
 
 	    /* build */
-	    if (plus->Line[i] == NULL) {
+	    if (plus->Line[line] == NULL) {
 		continue;
 	    }			/* dead line */
-	    Line = plus->Line[i];
+	    Line = plus->Line[line];
 	    if (Line->type != GV_BOUNDARY) {
 		continue;
 	    }
@@ -604,8 +616,8 @@ int Vect_build_nat(struct Map_info *Map, int build)
 		else
 		    side = GV_RIGHT;
 
-		G_debug(3, "Build area for line = %d, side = %d", i, side);
-		Vect_build_line_area(Map, i, side);
+		G_debug(3, "Build area for line = %d, side = %d", line, side);
+		Vect_build_line_area(Map, line, side);
 	    }
 	}
 	G_message(_("%d areas built"), plus->n_areas);
@@ -632,6 +644,7 @@ int Vect_build_nat(struct Map_info *Map, int build)
     /* Attach centroids to areas */
     if (plus->built < GV_BUILD_CENTROIDS) {
 	int nlines;
+	struct P_topo_c *topo;
 
 	G_important_message(_("Attaching centroids..."));
 
@@ -646,19 +659,21 @@ int Vect_build_nat(struct Map_info *Map, int build)
 	    if (Line->type != GV_CENTROID)
 		continue;
 
-	    area = Vect_find_area(Map, Line->E, Line->N);
+	    Vect_read_line(Map, Points, NULL, line);
+	    area = Vect_find_area(Map, Points->x[0], Points->y[0]);
 
 	    if (area > 0) {
 		G_debug(3, "Centroid (line=%d) in area %d", line, area);
 
 		Area = plus->Area[area];
+		topo = (struct P_topo_c *)Line->topo;
 
 		if (Area->centroid == 0) {	/* first */
 		    Area->centroid = line;
-		    Line->left = area;
+		    topo->area = area;
 		}
 		else {		/* duplicate */
-		    Line->left = -area;
+		    topo->area = -area;
 		}
 	    }
 	}
@@ -666,23 +681,23 @@ int Vect_build_nat(struct Map_info *Map, int build)
     }
 
     /* Add areas to category index */
-    for (area = 1; area <= plus->n_areas; area++) {
+    for (i = 1; i <= plus->n_areas; i++) {
 	int c;
 
-	if (plus->Area[area] == NULL)
+	if (plus->Area[i] == NULL)
 	    continue;
 
-	if (plus->Area[area]->centroid > 0) {
-	    Vect_read_line(Map, NULL, Cats, plus->Area[area]->centroid);
+	if (plus->Area[i]->centroid > 0) {
+	    Vect_read_line(Map, NULL, Cats, plus->Area[i]->centroid);
 
 	    for (c = 0; c < Cats->n_cats; c++) {
-		dig_cidx_add_cat(plus, Cats->field[c], Cats->cat[c], area,
+		dig_cidx_add_cat(plus, Cats->field[c], Cats->cat[c], i,
 				 GV_AREA);
 	    }
 	}
 
-	if (plus->Area[area]->centroid == 0 || Cats->n_cats == 0)	/* no centroid or no cats */
-	    dig_cidx_add_cat(plus, 0, 0, area, GV_AREA);
+	if (plus->Area[i]->centroid == 0 || Cats->n_cats == 0)	/* no centroid or no cats */
+	    dig_cidx_add_cat(plus, 0, 0, i, GV_AREA);
     }
 
     return 1;

lib/vector/Vlib/build_ogr.c

@@ -136,11 +136,11 @@ static int add_line(struct Map_info *Map, int type, struct line_pnts *Points,
 	offset = FID;		/* because centroids are read from topology, not from layer */
     }
     G_debug(4, "Register line: FID = %d offset = %ld", FID, offset);
-    line = dig_add_line(plus, type, Points, offset);
+    dig_line_box(Points, &box);
+    line = dig_add_line(plus, type, Points, &box, offset);
     G_debug(4, "Line registered with line = %d", line);
 
     /* Set box */
-    dig_line_box(Points, &box);
     if (line == 1)
 	Vect_box_copy(&(plus->box), &box);
     else
@@ -254,14 +254,12 @@ static int add_geometry(struct Map_info *Map, OGRGeometryH hGeom, int FID,
 	    else
 		lines[0] = -line;
 
-	    area = dig_add_area(plus, 1, lines);
-	    dig_area_set_box(plus, area, &box);
+	    area = dig_add_area(plus, 1, lines, &box);
 
 	    /* Each area is also isle */
 	    lines[0] = -lines[0];	/* island is counter clockwise */
 
-	    isle = dig_add_isle(plus, 1, lines);
-	    dig_isle_set_box(plus, isle, &box);
+	    isle = dig_add_isle(plus, 1, lines, &box);
 
 	    if (iPart == 0) {	/* outer ring */
 		outer_area = area;
@@ -287,16 +285,16 @@ static int add_geometry(struct Map_info *Map, OGRGeometryH hGeom, int FID,
 	}
 	else {
 	    struct P_area *Area;
+	    struct P_topo_c *topo;
 
 	    G_debug(4, "  Centroid: %f, %f", x, y);
 	    Vect_reset_line(Points[0]);
 	    Vect_append_point(Points[0], x, y, 0.0);
 	    line = add_line(Map, GV_CENTROID, Points[0], FID, parts);
-	    dig_line_box(Points[0], &box);
-	    dig_line_set_box(plus, line, &box);
 
 	    Line = plus->Line[line];
-	    Line->left = outer_area;
+	    topo = (struct P_topo_c *)Line->topo;
+	    topo->area = outer_area;
 
 	    /* register centroid to area */
 	    Area = plus->Area[outer_area];

lib/vector/Vlib/constraint.c

@@ -55,13 +55,13 @@ Vect_set_constraint_region(struct Map_info *Map,
 	return -1;
 
     Map->Constraint_region_flag = 1;
-    Map->Constraint_N = n;
-    Map->Constraint_S = s;
-    Map->Constraint_E = e;
-    Map->Constraint_W = w;
-    Map->Constraint_T = t;
-    Map->Constraint_B = b;
-    Map->proj = G_projection();
+    Map->Constraint_box.N = n;
+    Map->Constraint_box.S = s;
+    Map->Constraint_box.E = e;
+    Map->Constraint_box.W = w;
+    Map->Constraint_box.T = t;
+    Map->Constraint_box.B = b;
+    Map->head.proj = G_projection();
 
     return 0;
 }
@@ -76,12 +76,12 @@ Vect_set_constraint_region(struct Map_info *Map,
  */
 int Vect_get_constraint_box(const struct Map_info *Map, struct bound_box * Box)
 {
-    Box->N = Map->Constraint_N;
-    Box->S = Map->Constraint_S;
-    Box->E = Map->Constraint_E;
-    Box->W = Map->Constraint_W;
-    Box->T = Map->Constraint_T;
-    Box->B = Map->Constraint_B;
+    Box->N = Map->Constraint_box.N;
+    Box->S = Map->Constraint_box.S;
+    Box->E = Map->Constraint_box.E;
+    Box->W = Map->Constraint_box.W;
+    Box->T = Map->Constraint_box.T;
+    Box->B = Map->Constraint_box.B;
 
     return 0;
 }

lib/vector/Vlib/find.c

@@ -267,14 +267,14 @@ int Vect_find_area(struct Map_info *Map, double x, double y)
     int i, ret, area;
     static int first = 1;
     struct bound_box box;
-    static struct ilist *List;
+    static struct boxlist *List;
     static BOX_SIZE *size_list;
     static int alloc_size_list = 0;
 
     G_debug(3, "Vect_find_area() x = %f y = %f", x, y);
 
     if (first) {
-	List = Vect_new_list();
+	List = Vect_new_boxlist();
 	first = 0;
 	alloc_size_list = 10;
 	size_list = G_malloc(alloc_size_list * sizeof(BOX_SIZE));
@@ -287,7 +287,7 @@ int Vect_find_area(struct Map_info *Map, double x, double y)
     box.S = y;
     box.T = PORT_DOUBLE_MAX;
     box.B = -PORT_DOUBLE_MAX;
-    Vect_select_areas_by_box(Map, &box, List);
+    Vect_select_areas_by_box_with_box(Map, &box, List);
     G_debug(3, "  %d areas selected by box", List->n_values);
 
     /* sort areas by size, the smallest is likely to be the nearest */
@@ -297,16 +297,16 @@ int Vect_find_area(struct Map_info *Map, double x, double y)
     }
 
     for (i = 0; i < List->n_values; i++) {
-	size_list[i].i = area = List->value[i];
-	Vect_get_area_box(Map, area, &box);
+	size_list[i].i = List->id[i];
+	box = List->box[i];
 	size_list[i].size = (box.N - box.S) * (box.E - box.W);
     }
     
     if (List->n_values == 2) {
 	/* simple swap */
 	if (size_list[1].size < size_list[0].size) {
-	    size_list[0].i = List->value[1];
-	    size_list[1].i = List->value[0];
+	    size_list[0].i = List->id[1];
+	    size_list[1].i = List->id[0];
 	}
     }
     else if (List->n_values > 2)
@@ -339,13 +339,13 @@ int Vect_find_island(struct Map_info *Map, double x, double y)
     int i, ret, island, current, current_size, size;
     static int first = 1;
     struct bound_box box;
-    static struct ilist *List;
+    static struct boxlist *List;
     static struct line_pnts *Points;
 
     G_debug(3, "Vect_find_island() x = %f y = %f", x, y);
 
     if (first) {
-	List = Vect_new_list();
+	List = Vect_new_boxlist();
 	Points = Vect_new_line_struct();
 	first = 0;
     }
@@ -357,14 +357,14 @@ int Vect_find_island(struct Map_info *Map, double x, double y)
     box.S = y;
     box.T = PORT_DOUBLE_MAX;
     box.B = -PORT_DOUBLE_MAX;
-    Vect_select_isles_by_box(Map, &box, List);
+    Vect_select_isles_by_box_with_box(Map, &box, List);
     G_debug(3, "  %d islands selected by box", List->n_values);
 
     current_size = -1;
     current = 0;
     for (i = 0; i < List->n_values; i++) {
-	island = List->value[i];
-	ret = Vect_point_in_island(x, y, Map, island);
+	island = List->id[i];
+	ret = Vect_point_in_island(x, y, Map, island, List->box[i]);
 
 	if (ret >= 1) {		/* inside */
 	    if (current > 0) {	/* not first */

lib/vector/Vlib/header.c

@@ -434,8 +434,8 @@ const char *Vect_get_date(const struct Map_info *Map)
 */
 int Vect_set_person(struct Map_info *Map, const char *str)
 {
-    G_free(Map->head.your_name);
-    Map->head.your_name = G_store(str);
+    G_free(Map->head.user_name);
+    Map->head.user_name = G_store(str);
 
     return 0;
 }
@@ -449,7 +449,7 @@ int Vect_set_person(struct Map_info *Map, const char *str)
  */
 const char *Vect_get_person(const struct Map_info *Map)
 {
-    return Map->head.your_name;
+    return (Map->head.user_name);
 }
 
 /*!
@@ -545,8 +545,8 @@ int Vect_get_scale(const struct Map_info *Map)
  */
 int Vect_set_comment(struct Map_info *Map, const char *str)
 {
-    G_free(Map->head.line_3);
-    Map->head.line_3 = G_store(str);
+    G_free(Map->head.comment);
+    Map->head.comment = G_store(str);
     
     return 0;
 }
@@ -560,7 +560,7 @@ int Vect_set_comment(struct Map_info *Map, const char *str)
  */
 const char *Vect_get_comment(const struct Map_info *Map)
 {
-    return Map->head.line_3;
+    return (Map->head.comment);
 }
 
 /*!
@@ -606,7 +606,7 @@ int Vect_get_zone(const struct Map_info *Map)
  */
 int Vect_set_proj(struct Map_info *Map, int proj)
 {
-    Map->proj = proj;
+    Map->head.proj = proj;
     
     return 0;
 }
@@ -623,7 +623,7 @@ int Vect_set_proj(struct Map_info *Map, int proj)
 */
 int Vect_get_proj(const struct Map_info *Map)
 {
-    return Map->proj;
+    return (Map->head.proj);
 }
 
 /*!

lib/vector/Vlib/init_head.c

@@ -39,7 +39,7 @@ void Vect__init_head(struct Map_info *Map)
     Vect_set_date(Map, "");
 
     /* user name */
-    Map->head.your_name = NULL;
+    Map->head.user_name = NULL;
     sprintf(buf, "%s", G_whoami());
     Vect_set_person(Map, buf);
 
@@ -53,7 +53,7 @@ void Vect__init_head(struct Map_info *Map)
     Vect_set_map_date(Map, buf);
 
     /* comments */
-    Map->head.line_3 = NULL;
+    Map->head.comment = NULL;
     Vect_set_comment(Map, "");
 
     /* scale, threshold */

lib/vector/Vlib/intersect.c

@@ -81,8 +81,8 @@ static double dist2(double x1, double y1, double x2, double y2);
 #if 0
 static int ident(double x1, double y1, double x2, double y2, double thresh);
 #endif
-static int cross_seg(int id, int *arg);
-static int find_cross(int id, int *arg);
+static int cross_seg(int id, struct Rect rect, int *arg);
+static int find_cross(int id, struct Rect rect, int *arg);
 
 #define D  ((ax2-ax1)*(by1-by2) - (ay2-ay1)*(bx1-bx2))
 #define D1 ((bx1-ax1)*(by1-by2) - (by1-ay1)*(bx1-bx2))
@@ -462,7 +462,7 @@ typedef struct
 
 /* Current line in arrays is for some functions like cmp() set by: */
 static int current;
-static int second;		/* line whic is not current */
+static int second;		/* line which is not current */
 
 static int a_cross = 0;
 static int n_cross;
@@ -540,7 +540,7 @@ static int ident(double x1, double y1, double x2, double y2, double thresh)
 static struct line_pnts *APnts, *BPnts;
 
 /* break segments (called by rtree search) */
-static int cross_seg(int id, int *arg)
+static int cross_seg(int id, struct Rect rect, int *arg)
 {
     double x1, y1, z1, x2, y2, z2;
     int i, j, ret;
@@ -1059,6 +1059,8 @@ Vect_line_intersection(struct line_pnts *APoints,
 	    *BLines = XLines;
 	}
     }
+    
+    /* clean up */
 
     return 1;
 }
@@ -1068,7 +1070,7 @@ static struct line_pnts *APnts, *BPnts, *IPnts;
 static int cross_found;		/* set by find_cross() */
 
 /* break segments (called by rtree search) */
-static int find_cross(int id, int *arg)
+static int find_cross(int id, struct Rect rect, int *arg)
 {
     double x1, y1, z1, x2, y2, z2;
     int i, j, ret;

lib/vector/Vlib/legal_vname.c

@@ -90,6 +90,7 @@ int Vect_check_input_output_name(const char *input, const char *output,
 				 int error)
 {
     const char *mapset;
+    char nm[GNAME_MAX], ms[GMAPSET_MAX];
 
     if (Vect_legal_filename(output) == -1) {
 	if (error == GV_FATAL_EXIT) {
@@ -106,7 +107,14 @@ int Vect_check_input_output_name(const char *input, const char *output,
 	}
     }
 
-    mapset = G_find_vector2(input, "");
+    if (G_name_is_fully_qualified(input, nm, ms)) {
+	if (strcasecmp(ms, "ogr") != 0)
+	    mapset = G_find_vector2(input, "");
+	else
+	    mapset = ms;
+    }
+    else
+	mapset = G_find_vector2(input, "");
 
     if (mapset == NULL) {
 	if (error == GV_FATAL_EXIT) {
@@ -123,7 +131,6 @@ int Vect_check_input_output_name(const char *input, const char *output,
 
     if (strcmp(mapset, G_mapset()) == 0) {
 	const char *in;
-	char nm[GNAME_MAX], ms[GMAPSET_MAX];
 
 	if (G_name_is_fully_qualified(input, nm, ms)) {
 	    in = nm;

lib/vector/Vlib/level_two.c

@@ -23,9 +23,9 @@
 
    \return number of nodes
  */
-plus_t Vect_get_num_nodes(const struct Map_info *map)
+plus_t Vect_get_num_nodes(const struct Map_info *Map)
 {
-    return (map->plus.n_nodes);
+    return (Map->plus.n_nodes);
 }
 
 /*!
@@ -36,22 +36,22 @@ plus_t Vect_get_num_nodes(const struct Map_info *map)
 
    \return number of primitives
  */
-plus_t Vect_get_num_primitives(const struct Map_info *map, int type)
+plus_t Vect_get_num_primitives(const struct Map_info *Map, int type)
 {
     plus_t num = 0;
 
     if (type & GV_POINT)
-	num += map->plus.n_plines;
+	num += Map->plus.n_plines;
     if (type & GV_LINE)
-	num += map->plus.n_llines;
+	num += Map->plus.n_llines;
     if (type & GV_BOUNDARY)
-	num += map->plus.n_blines;
+	num += Map->plus.n_blines;
     if (type & GV_CENTROID)
-	num += map->plus.n_clines;
+	num += Map->plus.n_clines;
     if (type & GV_FACE)
-	num += map->plus.n_flines;
+	num += Map->plus.n_flines;
     if (type & GV_KERNEL)
-	num += map->plus.n_klines;
+	num += Map->plus.n_klines;
 
     return num;
 }
@@ -63,9 +63,9 @@ plus_t Vect_get_num_primitives(const struct Map_info *map, int type)
 
    \return number of features
  */
-plus_t Vect_get_num_lines(const struct Map_info *map)
+plus_t Vect_get_num_lines(const struct Map_info *Map)
 {
-    return (map->plus.n_lines);
+    return (Map->plus.n_lines);
 }
 
 /*!
@@ -75,9 +75,9 @@ plus_t Vect_get_num_lines(const struct Map_info *map)
 
    \return number of areas
  */
-plus_t Vect_get_num_areas(const struct Map_info *map)
+plus_t Vect_get_num_areas(const struct Map_info *Map)
 {
-    return (map->plus.n_areas);
+    return (Map->plus.n_areas);
 }
 
 /*!
@@ -87,9 +87,9 @@ plus_t Vect_get_num_areas(const struct Map_info *map)
 
    \return number of kernels
  */
-plus_t Vect_get_num_kernels(const struct Map_info *map)
+plus_t Vect_get_num_kernels(const struct Map_info *Map)
 {
-    return (map->plus.n_klines);
+    return (Map->plus.n_klines);
 }
 
 
@@ -100,9 +100,9 @@ plus_t Vect_get_num_kernels(const struct Map_info *map)
 
    \return number of faces
  */
-plus_t Vect_get_num_faces(const struct Map_info *map)
+plus_t Vect_get_num_faces(const struct Map_info *Map)
 {
-    return (map->plus.n_flines);
+    return (Map->plus.n_flines);
 }
 
 
@@ -113,9 +113,9 @@ plus_t Vect_get_num_faces(const struct Map_info *map)
 
    \return number of volumes
  */
-plus_t Vect_get_num_volumes(const struct Map_info *map)
+plus_t Vect_get_num_volumes(const struct Map_info *Map)
 {
-    return (map->plus.n_volumes);
+    return (Map->plus.n_volumes);
 }
 
 
@@ -126,9 +126,9 @@ plus_t Vect_get_num_volumes(const struct Map_info *map)
 
    \return number of islands
  */
-plus_t Vect_get_num_islands(const struct Map_info *map)
+plus_t Vect_get_num_islands(const struct Map_info *Map)
 {
-    return (map->plus.n_isles);
+    return (Map->plus.n_isles);
 }
 
 
@@ -139,9 +139,9 @@ plus_t Vect_get_num_islands(const struct Map_info *map)
 
    \return number of holes
  */
-plus_t Vect_get_num_holes(const struct Map_info *map)
+plus_t Vect_get_num_holes(const struct Map_info *Map)
 {
-    return (map->plus.n_holes);
+    return (Map->plus.n_holes);
 }
 
 
@@ -152,9 +152,10 @@ plus_t Vect_get_num_holes(const struct Map_info *map)
 
    \return number of dblinks
  */
-int Vect_get_num_dblinks(const struct Map_info *map)
+int Vect_get_num_dblinks(const struct Map_info *Map)
 {
-    return (map->dblnk->n_fields);
+    /* available on level 1 ? */
+    return (Map->dblnk->n_fields);
 }
 
 /*!
@@ -164,9 +165,9 @@ int Vect_get_num_dblinks(const struct Map_info *map)
 
    \return number of updated features
  */
-int Vect_get_num_updated_lines(const struct Map_info *map)
+int Vect_get_num_updated_lines(const struct Map_info *Map)
 {
-    return (map->plus.n_uplines);
+    return (Map->plus.n_uplines);
 }
 
 /*!
@@ -177,9 +178,9 @@ int Vect_get_num_updated_lines(const struct Map_info *map)
 
    \return updated line
  */
-int Vect_get_updated_line(const struct Map_info *map, int idx)
+int Vect_get_updated_line(const struct Map_info *Map, int idx)
 {
-    return (map->plus.uplines[idx]);
+    return (Map->plus.uplines[idx]);
 }
 
 /*!
@@ -189,9 +190,9 @@ int Vect_get_updated_line(const struct Map_info *map, int idx)
 
    \return number of updated nodes
  */
-int Vect_get_num_updated_nodes(const struct Map_info *map)
+int Vect_get_num_updated_nodes(const struct Map_info *Map)
 {
-    return (map->plus.n_upnodes);
+    return (Map->plus.n_upnodes);
 }
 
 /*!
@@ -202,9 +203,30 @@ int Vect_get_num_updated_nodes(const struct Map_info *map)
 
    \return updated node
  */
-int Vect_get_updated_node(const struct Map_info *map, int idx)
+int Vect_get_updated_node(const struct Map_info *Map, int idx)
 {
-    return (map->plus.upnodes[idx]);
+    return (Map->plus.upnodes[idx]);
+}
+
+/*!
+   \brief Get line type
+
+   \param map vector map
+   \param line line id
+
+   \return line type
+ */
+int
+Vect_get_line_type(const struct Map_info *Map, int line)
+{
+    if (Map->level < 2)
+	G_fatal_error(_("Vector map <%s> is not open on level >= 2"),
+		      Vect_get_full_name(Map));
+
+    if (!Vect_line_alive(Map, line))
+	return 0;
+	
+    return (Map->plus.Line[line]->type);
 }
 
 /*!
@@ -217,12 +239,12 @@ int Vect_get_updated_node(const struct Map_info *map, int idx)
    \return 0
  */
 int
-Vect_get_node_coor(const struct Map_info *map, int num, double *x, double *y,
+Vect_get_node_coor(const struct Map_info *Map, int num, double *x, double *y,
 		   double *z)
 {
     struct P_node *Node;
 
-    Node = map->plus.Node[num];
+    Node = Map->plus.Node[num];
     *x = Node->x;
     *y = Node->y;
 
@@ -243,16 +265,33 @@ Vect_get_node_coor(const struct Map_info *map, int num, double *x, double *y,
  */
 int Vect_get_line_nodes(const struct Map_info *Map, int line, int *n1, int *n2)
 {
+    char type;
 
     if (Map->level < 2)
 	G_fatal_error(_("Vector map <%s> is not open on level >= 2"),
 		      Vect_get_full_name(Map));
 
-    if (n1 != NULL)
-	*n1 = Map->plus.Line[line]->N1;
+    type = Vect_get_line_type(Map, line);
+
+    if (!(type & GV_LINES))
+	G_fatal_error(_("Nodes not available for line %d"), line);
+    
+    if (type == GV_LINE) {
+	struct P_topo_l *topo = (struct P_topo_l *)Map->plus.Line[line]->topo;
 
-    if (n2 != NULL)
-	*n2 = Map->plus.Line[line]->N2;
+	if (n1 != NULL)
+	    *n1 = topo->N1;
+	if (n2 != NULL)
+	    *n2 = topo->N2;
+    }
+    else if (type == GV_BOUNDARY) {
+	struct P_topo_b *topo = (struct P_topo_b *)Map->plus.Line[line]->topo;
+
+	if (n1 != NULL)
+	    *n1 = topo->N1;
+	if (n2 != NULL)
+	    *n2 = topo->N2;
+    }
 
     return 1;
 }
@@ -268,16 +307,24 @@ int Vect_get_line_nodes(const struct Map_info *Map, int line, int *n1, int *n2)
  */
 int Vect_get_line_areas(const struct Map_info *Map, int line, int *left, int *right)
 {
+    struct P_topo_b *topo;
 
     if (Map->level < 2)
 	G_fatal_error(_("Vector map <%s> is not open on level >= 2"),
 		      Vect_get_full_name(Map));
 
+    if (!Map->plus.Line[line]->topo)
+	G_fatal_error(_("Areas not available for line %d"), line);
+
+    if (Vect_get_line_type(Map, line) != GV_BOUNDARY)
+	G_fatal_error(_("Line %d is not a boundary"), line);
+
+    topo = (struct P_topo_b *)Map->plus.Line[line]->topo;
     if (left != NULL)
-	*left = Map->plus.Line[line]->left;
+	*left = topo->left;
 
     if (right != NULL)
-	*right = Map->plus.Line[line]->right;
+	*right = topo->right;
 
     return 1;
 }
@@ -349,9 +396,12 @@ float Vect_get_node_line_angle(const struct Map_info *Map, int node, int line)
  */
 int Vect_get_centroid_area(const struct Map_info *Map, int centroid)
 {
+    struct P_topo_c *topo;
+
     if (Map->level < 2)
 	G_fatal_error(_("Vector map <%s> is not open on level >= 2"),
 		      Vect_get_full_name(Map));
 
-    return (Map->plus.Line[centroid]->left);
+    topo = (struct P_topo_c *)Map->plus.Line[centroid]->topo;
+    return (topo->area);
 }

lib/vector/Vlib/list.c

@@ -207,3 +207,203 @@ int Vect_val_in_list(const struct ilist *list, int val)
 
     return 0;
 }
+
+/* box list routines */
+
+/**
+ * \brief Creates and initializes a struct boxlist.
+ *
+ * This structure is used as container for bounding boxes with id. The
+ * library routines handle all memory allocation.
+ *
+ * \return pointer to struct boxlist
+ * \return NULL on error
+ */
+struct boxlist *Vect_new_boxlist(void)
+{
+    struct boxlist *p;
+
+    p = (struct boxlist *)G_malloc(sizeof(struct boxlist));
+
+    if (p) {
+	p->id = NULL;
+	p->box = NULL;
+	p->n_values = 0;
+	p->alloc_values = 0;
+    }
+
+    return p;
+}
+
+/**
+ * \brief Reset boxlist structure.
+ *
+ * To make sure boxlist structure is clean to be re-used. List must have
+ * previously been created with Vect_new_boxlist().
+ *
+ * \param[in,out] list pointer to struct boxlist
+ * 
+ * \return 0
+ */
+int Vect_reset_boxlist(struct boxlist *list)
+{
+    list->n_values = 0;
+
+    return 0;
+}
+
+/**
+ * \brief Frees all memory associated with a struct boxlist, including
+ * the struct itself
+ *
+ * \param[in,out] list pointer to ilist structure
+ */
+void Vect_destroy_boxlist(struct boxlist *list)
+{
+    if (list) {			/* probably a moot test */
+	if (list->alloc_values) {
+	    G_free((void *)list->id);
+	    G_free((void *)list->box);
+	}
+	G_free((void *)list);
+    }
+    list = NULL;
+}
+
+/**
+ * \brief Append new item to the end of list if not yet present 
+ *
+ * \param[in,out] list pointer to ilist structure
+ * \param val new item to append to the end of list
+ *
+ * \return 0 on success
+ * \return 1 on error
+ */
+int Vect_boxlist_append(struct boxlist *list, int id, struct bound_box box)
+{
+    int i;
+    size_t size;
+
+    if (list == NULL)
+	return 1;
+
+    for (i = 0; i < list->n_values; i++) {
+	if (id == list->id[i])
+	    return 0;
+    }
+
+    if (list->n_values == list->alloc_values) {
+	size = (list->n_values + 1000) * sizeof(int);
+	list->id = (int *)G_realloc((void *)list->id, size);
+
+	size = (list->n_values + 1000) * sizeof(struct bound_box);
+	list->box = (struct bound_box *)G_realloc((void *)list->box, size);
+
+	list->alloc_values = list->n_values + 1000;
+    }
+
+    list->id[list->n_values] = id;
+    list->box[list->n_values] = box;
+    list->n_values++;
+
+    return 0;
+}
+
+/**
+ * \brief Append new items to the end of list if not yet present 
+ *
+ * \param[in,out] alist pointer to boxlist structure where items will be appended
+ * \param blist pointer to boxlist structure with new items
+ *
+ * \return 0 on success
+ * \return 1 on error
+ */
+int Vect_boxlist_append_boxlist(struct boxlist *alist, const struct boxlist *blist)
+{
+    int i;
+
+    if (alist == NULL || blist == NULL)
+	return 1;
+
+    for (i = 0; i < blist->n_values; i++)
+	Vect_boxlist_append(alist, blist->id[i], blist->box[i]);
+
+    return 0;
+}
+
+/**
+ * \brief Remove a given value (item) from list
+ *
+ * \param[in,out] list pointer to boxlist structure
+ * \param val to remove
+ *
+ * \return 0 on success
+ * \return 1 on error
+ */
+int Vect_boxlist_delete(struct boxlist *list, int id)
+{
+    int i, j;
+
+    if (list == NULL)
+	return 1;
+
+    for (i = 0; i < list->n_values; i++) {
+	if (id == list->id[i]) {
+	    for (j = i + 1; j < list->n_values; j++) {
+		list->id[j - 1] = list->id[j];
+		list->box[j - 1] = list->box[j];
+	    }
+
+	    list->n_values--;
+	    return 0;
+	}
+    }
+
+    return 0;
+}
+
+/**
+ * \brief Delete list from existing list 
+ *
+ * \param[in,out] alist pointer to original boxlist structure,
+ * \param blist pointer to boxlist structure with items to delete
+ *
+ * \return 0 on success
+ * \return 1 on error
+ */
+int Vect_boxlist_delete_boxlist(struct boxlist *alist, const struct boxlist *blist)
+{
+    int i;
+
+    if (alist == NULL || blist == NULL)
+	return 1;
+
+    for (i = 0; i < blist->n_values; i++)
+	Vect_boxlist_delete(alist, blist->id[i]);
+
+    return 0;
+}
+
+/**
+ * \brief Find a given item in the list
+ *
+ * \param list pointer to boxlist structure
+ * \param val value of item
+ *
+ * \return 1 if an item is found
+ * \return 0 no found item in the list
+*/
+int Vect_val_in_boxlist(const struct boxlist *list, int id)
+{
+    int i;
+
+    if (list == NULL)
+	return 0;
+
+    for (i = 0; i < list->n_values; i++) {
+	if (id == list->id[i])
+	    return 1;
+    }
+
+    return 0;
+}

lib/vector/Vlib/merge_lines.c

@@ -90,7 +90,15 @@ int Vect_merge_lines(struct Map_info *Map, int type, int *new_lines,
 
 	/* go backward as long as there is only one other line/boundary at the current node */
 	G_debug(3, "go backward");
-	next_node = Line->N1;
+	next_node = 0;
+	if (type == GV_LINE) {
+	    struct P_topo_l *topo = (struct P_topo_l *)Line->topo;
+	    next_node = topo->N1;
+	}
+	else if (type == GV_BOUNDARY) {
+	    struct P_topo_b *topo = (struct P_topo_b *)Line->topo;
+	    next_node = topo->N1;
+	}
 	first = -line;
 	while (1) {
 	    node_n_lines = Vect_get_node_n_lines(Map, next_node);
@@ -108,10 +116,26 @@ int Vect_merge_lines(struct Map_info *Map, int type, int *new_lines,
 		abs(next_line) != line) {
 		first = next_line;
 
-		if (first < 0)
-		    next_node = Plus->Line[-first]->N1;
-		else
-		    next_node = Plus->Line[first]->N2;
+		if (first < 0) {
+		    if (type == GV_LINE) {
+			struct P_topo_l *topo = (struct P_topo_l *)Plus->Line[-first]->topo;
+			next_node = topo->N1;
+		    }
+		    else if (type == GV_BOUNDARY) {
+			struct P_topo_b *topo = (struct P_topo_b *)Plus->Line[-first]->topo;
+			next_node = topo->N1;
+		    }
+		}
+		else {
+		    if (type == GV_LINE) {
+			struct P_topo_l *topo = (struct P_topo_l *)Plus->Line[first]->topo;
+			next_node = topo->N2;
+		    }
+		    else if (type == GV_BOUNDARY) {
+			struct P_topo_b *topo = (struct P_topo_b *)Plus->Line[first]->topo;
+			next_node = topo->N2;
+		    }
+		}
 	    }
 	    else
 		break;
@@ -123,10 +147,26 @@ int Vect_merge_lines(struct Map_info *Map, int type, int *new_lines,
 	/* reverse direction */
 	last = -first;
 
-	if (last < 0)
-	    next_node = Plus->Line[-last]->N1;
-	else
-	    next_node = Plus->Line[last]->N2;
+	if (last < 0) {
+	    if (type == GV_LINE) {
+		struct P_topo_l *topo = (struct P_topo_l *)Plus->Line[-last]->topo;
+		next_node = topo->N1;
+	    }
+	    else if (type == GV_BOUNDARY) {
+		struct P_topo_b *topo = (struct P_topo_b *)Plus->Line[-last]->topo;
+		next_node = topo->N1;
+	    }
+	}
+	else {
+	    if (type == GV_LINE) {
+		struct P_topo_l *topo = (struct P_topo_l *)Plus->Line[last]->topo;
+		next_node = topo->N2;
+	    }
+	    else if (type == GV_BOUNDARY) {
+		struct P_topo_b *topo = (struct P_topo_b *)Plus->Line[last]->topo;
+		next_node = topo->N2;
+	    }
+	}
 
 	Vect_reset_list(List);
 	while (1) {
@@ -147,10 +187,26 @@ int Vect_merge_lines(struct Map_info *Map, int type, int *new_lines,
 		abs(next_line) != abs(first)) {
 		last = next_line;
 
-		if (last < 0)
-		    next_node = Plus->Line[-last]->N1;
-		else
-		    next_node = Plus->Line[last]->N2;
+		if (last < 0) {
+		    if (type == GV_LINE) {
+			struct P_topo_l *topo = (struct P_topo_l *)Plus->Line[-last]->topo;
+			next_node = topo->N1;
+		    }
+		    else if (type == GV_BOUNDARY) {
+			struct P_topo_b *topo = (struct P_topo_b *)Plus->Line[-last]->topo;
+			next_node = topo->N1;
+		    }
+		}
+		else {
+		    if (type == GV_LINE) {
+			struct P_topo_l *topo = (struct P_topo_l *)Plus->Line[last]->topo;
+			next_node = topo->N2;
+		    }
+		    else if (type == GV_BOUNDARY) {
+			struct P_topo_b *topo = (struct P_topo_b *)Plus->Line[last]->topo;
+			next_node = topo->N2;
+		    }
+		}
 	    }
 	    else
 		break;

lib/vector/Vlib/net.c

@@ -223,11 +223,12 @@ Vect_net_build_graph(struct Map_info *Map,
     for (i = 1; i <= nlines; i++) {
 	G_percent(i, nlines, 1);	/* must be before any continue */
 	dofw = dobw = 1;
-	Vect_get_line_nodes(Map, i, &from, &to);
 	type = Vect_read_line(Map, Points, Cats, i);
 	if (!(type & ltype & (GV_LINE | GV_BOUNDARY)))
 	    continue;
 
+	Vect_get_line_nodes(Map, i, &from, &to);
+
 	if (afcol != NULL) {
 	    if (!(Vect_cat_get(Cats, afield, &cat))) {
 		G_debug(2,

lib/vector/Vlib/open.c

@@ -275,6 +275,8 @@ int Vect__open_old(struct Map_info *Map, const char *name, const char *mapset, c
 	G_warning(_("Unable to read header file of vector map <%s>"),
 		  Vect_get_full_name(Map));
     }
+
+    /* projection is not written to head but zone ??? */
     
     /* zone not set */
     if (Vect_get_zone(Map) == -1)

lib/vector/Vlib/overlap.c

@@ -31,11 +31,11 @@ V__map_overlap(struct Map_info *Map, double n, double s, double e, double w)
     struct Cell_head W;
 
     /* updated for Lat lon support 21 Jun 91 */
-    W.north = Map->Constraint_N;
-    W.south = Map->Constraint_S;
-    W.east = Map->Constraint_E;
-    W.west = Map->Constraint_W;
-    W.proj = Map->proj;
+    W.north = Map->Constraint_box.N;
+    W.south = Map->Constraint_box.S;
+    W.east = Map->Constraint_box.E;
+    W.west = Map->Constraint_box.W;
+    W.proj = Map->head.proj;
 
     return G_window_overlap(&W, n, s, e, w);
 }

lib/vector/Vlib/poly.c

@@ -681,6 +681,7 @@ int Vect_point_in_poly(double X, double Y, const struct line_pnts *Points)
    \param X,Y point coordinates
    \param Map vector map
    \param area area id
+   \param box area bounding box
 
    \return 0 - outside
    \return 1 - inside 
@@ -688,7 +689,7 @@ int Vect_point_in_poly(double X, double Y, const struct line_pnts *Points)
  */
 int
 Vect_point_in_area_outer_ring(double X, double Y, const struct Map_info *Map,
-			      int area)
+			      int area, struct bound_box box)
 {
     static int first = 1;
     int n_intersects, inter;
@@ -710,7 +711,7 @@ Vect_point_in_area_outer_ring(double X, double Y, const struct Map_info *Map,
     Area = Plus->Area[area];
 
     /* First it must be in box */
-    if (X < Area->W || X > Area->E || Y > Area->N || Y < Area->S)
+    if (X < box.W || X > box.E || Y > box.N || Y < box.S)
 	return 0;
 
     n_intersects = 0;
@@ -719,9 +720,11 @@ Vect_point_in_area_outer_ring(double X, double Y, const struct Map_info *Map,
 	G_debug(3, "  line[%d] = %d", i, line);
 
 	Line = Plus->Line[line];
-
+	
+	Vect_get_line_box(Map, line, &box);
+	
 	/* dont check lines that obviously do not intersect with test ray */
-	if ((Line->N < Y) || (Line->S > Y) || (Line->E < X))
+	if ((box.N < Y) || (box.S > Y) || (box.E < X))
 	    continue;
 
 	Vect_read_line(Map, Points, NULL, line);
@@ -747,12 +750,14 @@ Vect_point_in_area_outer_ring(double X, double Y, const struct Map_info *Map,
    \param X,Y point coordinates
    \param Map vector map
    \param isle isle id
+   \param box isle bounding box
 
    \return 0 - outside
    \return 1 - inside 
    \return 2 - on the boundary (exactly may be said only for vertex of vertical/horizontal line)
  */
-int Vect_point_in_island(double X, double Y, const struct Map_info *Map, int isle)
+int Vect_point_in_island(double X, double Y, const struct Map_info *Map,
+                         int isle, struct bound_box box)
 {
     static int first = 1;
     int n_intersects, inter;
@@ -772,7 +777,7 @@ int Vect_point_in_island(double X, double Y, const struct Map_info *Map, int isl
     Plus = &(Map->plus);
     Isle = Plus->Isle[isle];
 
-    if (X < Isle->W || X > Isle->E || Y > Isle->N || Y < Isle->S)
+    if (X < box.W || X > box.E || Y > box.N || Y < box.S)
 	return 0;
 
     n_intersects = 0;
@@ -781,8 +786,10 @@ int Vect_point_in_island(double X, double Y, const struct Map_info *Map, int isl
 
 	Line = Plus->Line[line];
 
+	Vect_get_line_box(Map, line, &box);
+	
 	/* dont check lines that obviously do not intersect with test ray */
-	if ((Line->N < Y) || (Line->S > Y) || (Line->E < X))
+	if ((box.N < Y) || (box.S > Y) || (box.E < X))
 	    continue;
 
 	Vect_read_line(Map, Points, NULL, line);

lib/vector/Vlib/read_nat.c

@@ -150,7 +150,7 @@ int
 V2_read_next_line_nat(struct Map_info *Map,
 		      struct line_pnts *line_p, struct line_cats *line_c)
 {
-    register int line;
+    register int line, ret;
     register struct P_line *Line;
     struct bound_box lbox, mbox;
 
@@ -177,15 +177,15 @@ V2_read_next_line_nat(struct Map_info *Map,
 	    continue;
 	}
 
+	ret = V2_read_line_nat(Map, line_p, line_c, Map->next_line++);
 	if (Map->Constraint_region_flag) {
-	    Vect_get_line_box(Map, line, &lbox);
+	    Vect_line_box(line_p, &lbox);
 	    if (!Vect_box_overlap(&lbox, &mbox)) {
-		Map->next_line++;
 		continue;
 	    }
 	}
 
-	return V2_read_line_nat(Map, line_p, line_c, Map->next_line++);
+	return ret;
     }
 
     /* NOTREACHED */ }
@@ -208,7 +208,7 @@ int
 Vect__Read_line_nat(struct Map_info *Map,
 		    struct line_pnts *p, struct line_cats *c, off_t offset)
 {
-    int i, dead = 0;
+    register int i, dead = 0;
     int n_points;
     off_t size;
     int n_cats, do_cats;

lib/vector/Vlib/read_ogr.c

@@ -462,15 +462,29 @@ int V2_read_line_ogr(struct Map_info *Map, struct line_pnts *line_p,
 		      _("Attempt to read dead line"), line);
 
     if (Line->type == GV_CENTROID) {
-	plus_t         node;
-	struct P_node *Node;
-	
 	G_debug(4, "Centroid");
-	node = Line->N1;
-	Node = Map->plus.Node[node];
-
+	
 	if (line_p != NULL) {
-	    Vect_append_point(line_p, Node->x, Node->y, 0.0);
+	    int i, found;
+	    struct bound_box box;
+	    struct boxlist list;
+	    struct P_topo_c *topo = (struct P_topo_c *)Line->topo;
+	    
+	    /* get area bbox */
+	    Vect_get_area_box(Map, topo->area, &box);
+	    /* search in spatial index for centroid with area bbox */
+	    dig_init_boxlist(&list);
+	    Vect_select_lines_by_box_with_box(Map, &box, Line->type, &list);
+	    
+	    found = 0;
+	    for (i = 0; i < list.n_values; i++) {
+		if (list.id[i] == line) {
+		    found = i;
+		    break;
+		}
+	    }
+
+	    Vect_append_point(line_p, list.box[found].E, list.box[found].N, 0.0);
 	}
 
 	if (line_c != NULL) {

lib/vector/Vlib/remove_duplicates.c

@@ -23,6 +23,7 @@
    Remove duplicate lines of given types from vector map. Duplicate
    lines may be optionally written to error map. Input map must be
    opened on level 2 for update. Categories are merged.
+   GV_BUILD_BASE is sufficient.
 
    \param[in,out] Map vector map where duplicate lines will be deleted
    \param type type of line to be delete

lib/vector/Vlib/select.c

@@ -106,7 +106,7 @@ void Vect_spatial_index_del_item(struct spatial_index * si, int id,
 
 /************************* SELECT BY BOX *********************************/
 /* This function is called by  RTreeSearch() to add selected item to the list */
-static int _add_item(int id, struct ilist *list)
+static int _add_item(int id, struct Rect rect, struct ilist *list)
 {
     dig_list_add(list, id);
     return 1;

lib/vector/Vlib/sindex.c

@@ -5,12 +5,12 @@
 
    Higher level functions for reading/writing/manipulating vectors.
 
-   (C) 2001-2009 by the GRASS Development Team
+   (C) 2001-2011 by the GRASS Development Team
 
    This program is free software under the GNU General Public License
    (>=v2).  Read the file COPYING that comes with GRASS for details.
 
-   \author Radim Blazek
+   \author Radim Blazek, Markus Metz
  */
 
 #include <stdlib.h>
@@ -72,6 +72,63 @@ Vect_select_lines_by_box(struct Map_info *Map, const struct bound_box * Box,
 }
 
 /*!
+   \brief Select lines with bounding boxes by box.
+
+   Select lines whose boxes overlap specified box!!!  It means that
+   selected line may or may not overlap the box.
+
+   \param Map vector map
+   \param Box bounding box
+   \param type line type
+   \param[out] list output list, must be initialized
+
+   \return number of lines
+ */
+int
+Vect_select_lines_by_box_with_box(struct Map_info *Map, const struct bound_box *Box,
+			 int type, struct boxlist *list)
+{
+    int i, line, nlines;
+    struct Plus_head *plus;
+    struct P_line *Line;
+    static struct boxlist *LocList = NULL;
+
+    G_debug(3, "Vect_select_lines_by_box()");
+    G_debug(3, "  Box(N,S,E,W,T,B): %e, %e, %e, %e, %e, %e", Box->N, Box->S,
+	    Box->E, Box->W, Box->T, Box->B);
+    plus = &(Map->plus);
+
+    if (!(plus->Spidx_built)) {
+	G_debug(3, "Building spatial index.");
+	Vect_build_sidx_from_topo(Map);
+    }
+
+    list->n_values = 0;
+    if (!LocList) {
+	LocList = (struct boxlist *)G_malloc(sizeof(struct boxlist));
+	dig_init_boxlist(LocList);
+    }
+
+    nlines = dig_select_lines_with_box(plus, Box, LocList);
+    G_debug(3, "  %d lines selected (all types)", nlines);
+
+    /* Remove lines of not requested types */
+    for (i = 0; i < nlines; i++) {
+	line = LocList->id[i];
+	if (plus->Line[line] == NULL)
+	    continue;		/* Should not happen */
+	Line = plus->Line[line];
+	if (!(Line->type & type))
+	    continue;
+	dig_boxlist_add(list, line, LocList->box[i]);
+    }
+
+    G_debug(3, "  %d lines of requested type", list->n_values);
+
+    return list->n_values;
+}
+
+/*!
    \brief Select areas by box.
 
    Select areas whose boxes overlap specified box!!!
@@ -124,6 +181,53 @@ Vect_select_areas_by_box(struct Map_info *Map, const struct bound_box * Box,
 
 
 /*!
+   \brief Select areas with bounding boxes by box.
+
+   Select areas whose boxes overlap specified box!!!
+   It means that selected area may or may not overlap the box.
+
+   \param Map vector map
+   \param Box bounding box
+   \param[out] output list, must be initialized
+
+   \return number of areas
+ */
+int
+Vect_select_areas_by_box_with_box(struct Map_info *Map, const struct bound_box * Box,
+			 struct boxlist *list)
+{
+    int i;
+    static int debug_level = -1;
+
+    if (debug_level == -1) {
+	const char *dstr = G__getenv("DEBUG");
+
+	if (dstr != NULL)
+	    debug_level = atoi(dstr);
+	else
+	    debug_level = 0;
+    }
+
+    G_debug(3, "Vect_select_areas_by_box()");
+    G_debug(3, "Box(N,S,E,W,T,B): %e, %e, %e, %e, %e, %e", Box->N, Box->S,
+	    Box->E, Box->W, Box->T, Box->B);
+
+    dig_select_areas_with_box(&(Map->plus), Box, list);
+    G_debug(3, "  %d areas selected", list->n_values);
+    /* avoid loop when not debugging */
+    if (debug_level > 2) {
+	for (i = 0; i < list->n_values; i++) {
+	    G_debug(3, "  area = %d pointer to area structure = %lx",
+		    list->id[i],
+		    (unsigned long)Map->plus.Area[list->id[i]]);
+	}
+    }
+
+    return list->n_values;
+}
+
+
+/*!
    \brief Select isles by box.
 
    Select isles whose boxes overlap specified box!!!
@@ -155,6 +259,32 @@ Vect_select_isles_by_box(struct Map_info *Map, const struct bound_box * Box,
 }
 
 /*!
+   \brief Select isles with bounding boxes by box.
+
+   Select isles whose boxes overlap specified box!!!
+   It means that selected isle may or may not overlap the box.
+
+   \param Map vector map
+   \param Box bounding box
+   \param[out] list output list, must be initialized
+
+   \return number of isles
+ */
+int
+Vect_select_isles_by_box_with_box(struct Map_info *Map, const struct bound_box * Box,
+			 struct boxlist *list)
+{
+    G_debug(3, "Vect_select_isles_by_box_with_box()");
+    G_debug(3, "Box(N,S,E,W,T,B): %e, %e, %e, %e, %e, %e", Box->N, Box->S,
+	    Box->E, Box->W, Box->T, Box->B);
+
+    dig_select_isles_with_box(&(Map->plus), Box, list);
+    G_debug(3, "  %d isles selected", list->n_values);
+
+    return list->n_values;
+}
+
+/*!
    \brief Select nodes by box.
 
    \param Map vector map

lib/vector/Vlib/snap.c

@@ -40,8 +40,8 @@ typedef struct
     double along;
 } NEW;
 
-/* This function is called by RTreeSearch() to add selected node/line/area/isle to thelist */
-int add_item(int id, struct ilist *list)
+/* This function is called by RTreeSearch() to add selected node/line/area/isle to the list */
+int add_item(int id, struct Rect rect, struct ilist *list)
 {
     dig_list_add(list, id);
     return 1;

lib/vector/Vlib/write_nat.c

@@ -152,10 +152,13 @@ static void V2__add_line_to_topo_nat(struct Map_info *Map, int line,
 
     if (plus->built >= GV_BUILD_AREAS) {
 	if (type == GV_BOUNDARY) {
+
 	    /* Delete neighbour areas/isles */
 	    first = 1;
 	    for (s = 0; s < 2; s++) {	/* for each node */
-		node = (s == 0 ? Line->N1 : Line->N2);
+		struct P_topo_b *topo = (struct P_topo_b *)Line->topo;
+
+		node = (s == 0 ? topo->N1 : topo->N2);
 		G_debug(3,
 			"  delete neighbour areas/isles: %s node = %d",
 			(s == 0 ? "first" : "second"), node);
@@ -182,10 +185,11 @@ static void V2__add_line_to_topo_nat(struct Map_info *Map, int line,
 
 		    if (next_line != 0) {	/* there is a boundary to the right */
 			NLine = plus->Line[abs(next_line)];
+			topo = (struct P_topo_b *)NLine->topo;
 			if (next_line > 0)	/* the boundary is connected by 1. node */
-			    area = NLine->right;	/* we are interested just in this side (close to our line) */
+			    area = topo->right;	/* we are interested just in this side (close to our line) */
 			else if (next_line < 0)	/* the boundary is connected by 2. node */
-			    area = NLine->left;
+			    area = topo->left;
 
 			G_debug(3, "  next_line = %d area = %d", next_line,
 				area);
@@ -265,16 +269,19 @@ static void V2__add_line_to_topo_nat(struct Map_info *Map, int line,
 
     /* Attach centroid */
     if (plus->built >= GV_BUILD_CENTROIDS) {
+	struct P_topo_c *topo;
+
 	if (type == GV_CENTROID) {
 	    sel_area = Vect_find_area(Map, points->x[0], points->y[0]);
 	    G_debug(3, "  new centroid %d is in area %d", line, sel_area);
 	    if (sel_area > 0) {
 		Area = plus->Area[sel_area];
 		Line = plus->Line[line];
+		topo = (struct P_topo_c *)Line->topo;
 		if (Area->centroid == 0) {	/* first centroid */
 		    G_debug(3, "  first centroid -> attach to area");
 		    Area->centroid = line;
-		    Line->left = sel_area;
+		    topo->area = sel_area;
 		    if (plus->update_cidx) {
 			V2__add_area_cats_to_cidx_nat(Map, sel_area);
 		    }
@@ -282,7 +289,7 @@ static void V2__add_line_to_topo_nat(struct Map_info *Map, int line,
 		else {		/* duplicate centroid */
 		    G_debug(3,
 			    "  duplicate centroid -> do not attach to area");
-		    Line->left = -sel_area;
+		    topo->area = -sel_area;
 		}
 	    }
 	}
@@ -353,10 +360,9 @@ off_t V2_write_line_nat(struct Map_info *Map, int type,
     plus = &(Map->plus);
     /* Add line */
     if (plus->built >= GV_BUILD_BASE) {
-	line = dig_add_line(plus, type, points, offset);
-	G_debug(3, "  line added to topo with id = %d", line);
 	dig_line_box(points, &box);
-	dig_line_set_box(plus, line, &box);
+	line = dig_add_line(plus, type, points, &box, offset);
+	G_debug(3, "  line added to topo with id = %d", line);
 	if (line == 1)
 	    Vect_box_copy(&(plus->box), &box);
 	else
@@ -475,10 +481,9 @@ off_t V2_rewrite_line_nat(struct Map_info *Map, int line, int type, off_t old_of
     plus = &(Map->plus);
     /* Add line */
     if (plus->built >= GV_BUILD_BASE) {
-	line = dig_add_line(plus, type, points, offset);
-	G_debug(3, "  line added to topo with id = %d", line);
 	dig_line_box(points, &box);
-	dig_line_set_box(plus, line, &box);
+	line = dig_add_line(plus, type, points, &box, offset);
+	G_debug(3, "  line added to topo with id = %d", line);
 	if (line == 1)
 	    Vect_box_copy(&(plus->box), &box);
 	else
@@ -654,6 +659,7 @@ int V2_delete_line_nat(struct Map_info *Map, off_t line)
     struct bound_box box, abox;
     int adjacent[4], n_adjacent;
     static struct line_cats *Cats = NULL;
+    static struct line_pnts *Points = NULL;
 
     G_debug(3, "V2_delete_line_nat(), line = %lu", (unsigned long) line);
 
@@ -672,10 +678,14 @@ int V2_delete_line_nat(struct Map_info *Map, off_t line)
     if (!Cats) {
 	Cats = Vect_new_cats_struct();
     }
+    if (!Points) {
+	Points = Vect_new_line_struct();
+    }
+
+    type = V2_read_line_nat(Map, Points, Cats, line);
 
     /* Update category index */
     if (plus->update_cidx) {
-	type = V2_read_line_nat(Map, NULL, Cats, line);
 
 	for (i = 0; i < Cats->n_cats; i++) {
 	    dig_cidx_del_cat(plus, Cats->field[i], Cats->cat[i], line, type);
@@ -691,6 +701,8 @@ int V2_delete_line_nat(struct Map_info *Map, off_t line)
 
     /* Update topology */
     if (plus->built >= GV_BUILD_AREAS && type == GV_BOUNDARY) {
+	struct P_topo_b *topo = (struct P_topo_b *)Line->topo;
+
 	/* Store adjacent boundaries at nodes (will be used to rebuild area/isle) */
 	/* Adjacent are stored: > 0 - we want right side; < 0 - we want left side */
 	n_adjacent = 0;
@@ -722,8 +734,8 @@ int V2_delete_line_nat(struct Map_info *Map, off_t line)
 
 	/* Delete area(s) and islands this line forms */
 	first = 1;
-	if (Line->left > 0) {	/* delete area */
-	    Vect_get_area_box(Map, Line->left, &box);
+	if (topo->left > 0) {	/* delete area */
+	    Vect_get_area_box(Map, topo->left, &box);
 	    if (first) {
 		Vect_box_copy(&abox, &box);
 		first = 0;
@@ -732,15 +744,15 @@ int V2_delete_line_nat(struct Map_info *Map, off_t line)
 		Vect_box_extend(&abox, &box);
 
 	    if (plus->update_cidx) {
-		V2__delete_area_cats_from_cidx_nat(Map, Line->left);
+		V2__delete_area_cats_from_cidx_nat(Map, topo->left);
 	    }
-	    dig_del_area(plus, Line->left);
+	    dig_del_area(plus, topo->left);
 	}
-	else if (Line->left < 0) {	/* delete isle */
-	    dig_del_isle(plus, -Line->left);
+	else if (topo->left < 0) {	/* delete isle */
+	    dig_del_isle(plus, -topo->left);
 	}
-	if (Line->right > 0) {	/* delete area */
-	    Vect_get_area_box(Map, Line->right, &box);
+	if (topo->right > 0) {	/* delete area */
+	    Vect_get_area_box(Map, topo->right, &box);
 	    if (first) {
 		Vect_box_copy(&abox, &box);
 		first = 0;
@@ -749,29 +761,31 @@ int V2_delete_line_nat(struct Map_info *Map, off_t line)
 		Vect_box_extend(&abox, &box);
 
 	    if (plus->update_cidx) {
-		V2__delete_area_cats_from_cidx_nat(Map, Line->right);
+		V2__delete_area_cats_from_cidx_nat(Map, topo->right);
 	    }
-	    dig_del_area(plus, Line->right);
+	    dig_del_area(plus, topo->right);
 	}
-	else if (Line->right < 0) {	/* delete isle */
-	    dig_del_isle(plus, -Line->right);
+	else if (topo->right < 0) {	/* delete isle */
+	    dig_del_isle(plus, -topo->right);
 	}
     }
 
     /* Delete reference from area */
     if (plus->built >= GV_BUILD_CENTROIDS && type == GV_CENTROID) {
-	if (Line->left > 0) {
-	    G_debug(3, "Remove centroid %d from area %d", (int) line, Line->left);
+	struct P_topo_c *topo = (struct P_topo_c *)Line->topo;
+
+	if (topo->area > 0) {
+	    G_debug(3, "Remove centroid %d from area %d", (int) line, topo->area);
 	    if (plus->update_cidx) {
-		V2__delete_area_cats_from_cidx_nat(Map, Line->left);
+		V2__delete_area_cats_from_cidx_nat(Map, topo->area);
 	    }
-	    Area = Map->plus.Area[Line->left];
+	    Area = Map->plus.Area[topo->area];
 	    Area->centroid = 0;
 	}
     }
 
     /* delete the line from topo */
-    dig_del_line(plus, line);
+    dig_del_line(plus, line, Points->x[0], Points->y[0], Points->z[0]);
 
     /* Rebuild areas/isles and attach centroids and isles */
     if (plus->built >= GV_BUILD_AREAS && type == GV_BOUNDARY) {
@@ -938,10 +952,9 @@ int V2_restore_line_nat(struct Map_info *Map, int line, off_t offset)
     
     /* restore the line from topo */		   
     if (plus->built >= GV_BUILD_BASE) {
-	dig_restore_line(plus, line, type, points, offset);
-	G_debug(3, "  line restored in topo with id = %d", line);
 	dig_line_box(points, &box);
-	dig_line_set_box(plus, line, &box);
+	dig_restore_line(plus, line, type, points, &box, offset);
+	G_debug(3, "  line restored in topo with id = %d", line);
 	Vect_box_extend(&(plus->box), &box);
     }
     

lib/vector/Vlib/write_ogr.c

@@ -174,10 +174,9 @@ off_t V2_write_line_ogr(struct Map_info *Map, int type,
     plus = &(Map->plus);
     /* Add line */
     if (plus->built >= GV_BUILD_BASE) {
-	line = dig_add_line(plus, type, points, offset);
-	G_debug(3, "  line added to topo with id = %d", line);
 	dig_line_box(points, &box);
-	dig_line_set_box(plus, line, &box);
+	line = dig_add_line(plus, type, points, &box, offset);
+	G_debug(3, "  line added to topo with id = %d", line);
 	if (line == 1)
 	    Vect_box_copy(&(plus->box), &box);
 	else
@@ -288,6 +287,7 @@ int V2_delete_line_ogr(struct Map_info *Map, off_t line)
     struct P_line *Line;
     struct Plus_head *plus;
     static struct line_cats *Cats = NULL;
+    static struct line_pnts *Points = NULL;
 
     G_debug(3, "V2_delete_line_nat(), line = %d", (int) line);
 
@@ -306,10 +306,14 @@ int V2_delete_line_ogr(struct Map_info *Map, off_t line)
     if (!Cats) {
 	Cats = Vect_new_cats_struct();
     }
+    if (!Points) {
+	Points = Vect_new_line_struct();
+    }
+
+    type = V2_read_line_nat(Map, Points, Cats, line);
 
     /* Update category index */
     if (plus->update_cidx) {
-	type = V2_read_line_ogr(Map, NULL, Cats, line);
 
 	for (i = 0; i < Cats->n_cats; i++) {
 	    dig_cidx_del_cat(plus, Cats->field[i], Cats->cat[i], line, type);
@@ -335,7 +339,7 @@ int V2_delete_line_ogr(struct Map_info *Map, off_t line)
     }
 
     /* delete the line from topo */
-    dig_del_line(plus, line);
+    dig_del_line(plus, line, Points->x[0], Points->y[0], Points->z[0]);
 
     /* Rebuild areas/isles and attach centroids and isles */
     if (plus->built >= GV_BUILD_AREAS && type == GV_BOUNDARY) {