[GRASS-SVN] r54134 - grass/branches/develbranch_6/vector/v.overlay

svn_grass at osgeo.org svn_grass at osgeo.org
Sat Dec 1 12:03:35 PST 2012


Author: mmetz
Date: 2012-12-01 12:03:35 -0800 (Sat, 01 Dec 2012)
New Revision: 54134

Modified:
   grass/branches/develbranch_6/vector/v.overlay/area_area.c
   grass/branches/develbranch_6/vector/v.overlay/description.html
   grass/branches/develbranch_6/vector/v.overlay/line_area.c
   grass/branches/develbranch_6/vector/v.overlay/local.h
   grass/branches/develbranch_6/vector/v.overlay/main.c
Log:
v.overlay: use Vect_snap_line()

Modified: grass/branches/develbranch_6/vector/v.overlay/area_area.c
===================================================================
--- grass/branches/develbranch_6/vector/v.overlay/area_area.c	2012-12-01 19:59:23 UTC (rev 54133)
+++ grass/branches/develbranch_6/vector/v.overlay/area_area.c	2012-12-01 20:03:35 UTC (rev 54134)
@@ -15,9 +15,18 @@
 #include <grass/glocale.h>
 #include "local.h"
 
+/* for ilist qsort'ing and bsearch'ing */
+static int cmp_int(const void *a, const void *b)
+{
+    int ai = *(int *)a;
+    int bi = *(int *)b;
+    
+    return (ai < bi ? -1 : (ai > bi));
+}
+
 int area_area(struct Map_info *In, int *field, struct Map_info *Out,
 	      struct field_info *Fi, dbDriver * driver, int operator,
-	      int *ofield, ATTRIBUTES * attr)
+	      int *ofield, ATTRIBUTES * attr, struct ilist *BList, double snap)
 {
     int ret, input, line, nlines, area, nareas;
     int in_area, in_centr, out_cat;
@@ -28,16 +37,78 @@
     char buf[1000];
     dbString stmt;
     int nmodif;
+    int verbose;
 
+    verbose = G_verbose();
+
     Points = Vect_new_line_struct();
     Cats = Vect_new_cats_struct();
 
-    /* Vect_clean_small_angles_at_nodes() can change the geometry so that new intersections
+    /* optional snap */
+    if (snap > 0) {
+	int i, j, snapped_lines;
+	struct bound_box box;
+	struct ilist *list = Vect_new_list();
+	struct ilist *reflist = Vect_new_list();
+	
+	G_message(_("Snapping boundaries with %g ..."), snap);
+
+	/* snap boundaries in B to boundaries in A
+	 * not modifying boundaries in A */
+
+	if (BList->n_values > 1)
+	    qsort(BList->value, BList->n_values, sizeof(int), cmp_int);
+
+	snapped_lines = 0;
+	nlines = BList->n_values;
+	for (i = 0; i < nlines; i++) {
+	    line = BList->value[i];
+	    Vect_read_line(Out, Points, Cats, line);
+	    /* select lines by box */
+	    Vect_get_line_box(Out, line, &box);
+	    box.E += snap;
+	    box.W -= snap;
+	    box.N += snap;
+	    box.S -= snap;
+	    box.T = 0.0;
+	    box.B = 0.0;
+	    Vect_select_lines_by_box(Out, &box, GV_BOUNDARY, list);
+	    
+	    if (list->n_values > 0) {
+		Vect_reset_list(reflist);
+		for (j = 0; j < list->n_values; j++) {
+		    int aline = list->value[j];
+
+		    if (!bsearch(&aline, BList->value, BList->n_values,
+			sizeof(int), cmp_int)) {
+			Vect_list_append(reflist, aline);
+		    }
+		}
+		
+		/* snap bline to lines */
+		if (Vect_snap_line(Out, reflist, Points, snap, NULL, NULL)) {
+		    /* rewrite bline*/
+		    Vect_delete_line(Out, line);
+		    ret = Vect_write_line(Out, GV_BOUNDARY, Points, Cats);
+		    Vect_list_append(BList, ret);
+		    snapped_lines++;
+		    G_debug(3, "line %d snapped", line);
+		}
+	    }
+	}
+	Vect_destroy_list(list);
+	Vect_destroy_list(reflist);
+
+	G_verbose_message(_("%d boundaries snapped"), snapped_lines);
+    }
+
+    /* same procedure like for v.in.ogr
+     * Vect_clean_small_angles_at_nodes() can change the geometry so that new intersections
      * are created. We must call Vect_break_lines(), Vect_remove_duplicates()
      * and Vect_clean_small_angles_at_nodes() until no more small dangles are found */
     do {
 	G_message(_("Breaking lines..."));
-	Vect_break_lines(Out, GV_LINE | GV_BOUNDARY, NULL);
+	Vect_break_lines_list(Out, NULL, BList, GV_BOUNDARY, NULL);
 
 	/* Probably not necessary for LINE x AREA */
 	G_message(_("Removing duplicates..."));
@@ -51,8 +122,40 @@
     /* ?: May be result of Vect_break_lines() + Vect_remove_duplicates() any dangle or bridge?
      * In that case, calls to Vect_remove_dangles() and Vect_remove_bridges() would be also necessary */
 
+    G_set_verbose(0);
+    Vect_build_partial(Out, GV_BUILD_AREAS);
+    G_set_verbose(verbose);
+    nlines = Vect_get_num_lines(Out);
+    ret = 0;
+    for (line = 1; line <= nlines; line++) {
+	if (!Vect_line_alive(Out, line))
+	    continue;
+	if (Vect_read_line(Out, NULL, NULL, line) == GV_BOUNDARY) {
+	    int left, rite;
+	    
+	    Vect_get_line_areas(Out, line, &left, &rite);
+	    
+	    if (left == 0 || rite == 0) {
+		ret = 1;
+		break;
+	    }
+	}
+    }
+    if (ret) {
+	Vect_remove_dangles(Out, GV_BOUNDARY, -1, NULL);
+	Vect_remove_bridges(Out, NULL);
+    }
+
+    G_set_verbose(0);
+    Vect_build_partial(Out, GV_BUILD_NONE);
+    Vect_build_partial(Out, GV_BUILD_BASE);
+    G_set_verbose(verbose);
+    G_message(_("Merging lines..."));
+    Vect_merge_lines(Out, GV_BOUNDARY, NULL, NULL);
+
     /* Attach islands */
     G_message(_("Attaching islands..."));
+    Vect_build_partial(Out, GV_BUILD_NONE);
     Vect_build_partial(Out, GV_BUILD_ATTACH_ISLES);
 
 
@@ -81,6 +184,8 @@
 	for (area = 1; area <= nareas; area++) {
 	    Centr[area].cat[input] = Vect_new_cats_struct();
 
+	    G_percent(area, nareas, 1);
+
 	    in_area =
 		Vect_find_area(&(In[input]), Centr[area].x, Centr[area].y);
 	    if (in_area > 0) {
@@ -94,7 +199,7 @@
 			if (Cats->field[i] == field[input]) {
 			    ATTR *at;
 
-			    Vect_cat_set(Centr[area].cat[input], field[input],
+			    Vect_cat_set(Centr[area].cat[input], ofield[input + 1],
 					 Cats->cat[i]);
 
 			    /* Mark as used */
@@ -107,7 +212,6 @@
 		    }
 		}
 	    }
-	    G_percent(area, nareas, 1);
 	}
     }
 
@@ -118,6 +222,8 @@
     for (area = 1; area <= nareas; area++) {
 	int i;
 
+	G_percent(area, nareas, 1);
+
 	/* check the condition */
 	switch (operator) {
 	case OP_AND:
@@ -152,112 +258,114 @@
 
 	Vect_append_point(Points, Centr[area].x, Centr[area].y, 0.0);
 
-	/* Add new cats for all combinations of input cats (-1 in cycle for null) */
-	for (i = -1; i < Centr[area].cat[0]->n_cats; i++) {
-	    int j;
+	if (ofield[0] > 0) {
+	    /* Add new cats for all combinations of input cats (-1 in cycle for null) */
+	    for (i = -1; i < Centr[area].cat[0]->n_cats; i++) {
+		int j;
 
-	    if (i == -1 && Centr[area].cat[0]->n_cats > 0)
-		continue;	/* no need to make null */
-
-	    for (j = -1; j < Centr[area].cat[1]->n_cats; j++) {
-		if (j == -1 && Centr[area].cat[1]->n_cats > 0)
+		if (i == -1 && Centr[area].cat[0]->n_cats > 0)
 		    continue;	/* no need to make null */
 
-		if (ofield[0] > 0)
-		    Vect_cat_set(Cats, ofield[0], out_cat);
+		for (j = -1; j < Centr[area].cat[1]->n_cats; j++) {
+		    if (j == -1 && Centr[area].cat[1]->n_cats > 0)
+			continue;	/* no need to make null */
 
-		/* attributes */
-		if (driver) {
-		    ATTR *at;
+		    if (ofield[0] > 0)
+			Vect_cat_set(Cats, ofield[0], out_cat);
 
-		    sprintf(buf, "insert into %s values ( %d", Fi->table,
-			    out_cat);
-		    db_set_string(&stmt, buf);
+		    /* attributes */
+		    if (driver) {
+			ATTR *at;
 
-		    /* cata */
-		    if (i >= 0) {
-			if (attr[0].columns) {
-			    at = find_attr(&(attr[0]),
-					   Centr[area].cat[0]->cat[i]);
-			    if (!at)
-				G_fatal_error(_("Attribute not found"));
+			sprintf(buf, "insert into %s values ( %d", Fi->table,
+				out_cat);
+			db_set_string(&stmt, buf);
 
-			    if (at->values)
-				db_append_string(&stmt, at->values);
-			    else
-				db_append_string(&stmt, attr[0].null_values);
+			/* cata */
+			if (i >= 0) {
+			    if (attr[0].columns) {
+				at = find_attr(&(attr[0]),
+					       Centr[area].cat[0]->cat[i]);
+				if (!at)
+				    G_fatal_error(_("Attribute not found"));
+
+				if (at->values)
+				    db_append_string(&stmt, at->values);
+				else
+				    db_append_string(&stmt, attr[0].null_values);
+			    }
+			    else {
+				sprintf(buf, ", %d", Centr[area].cat[0]->cat[i]);
+				db_append_string(&stmt, buf);
+			    }
 			}
 			else {
-			    sprintf(buf, ", %d", Centr[area].cat[0]->cat[i]);
-			    db_append_string(&stmt, buf);
+			    if (attr[0].columns) {
+				db_append_string(&stmt, attr[0].null_values);
+			    }
+			    else {
+				sprintf(buf, ", null");
+				db_append_string(&stmt, buf);
+			    }
 			}
-		    }
-		    else {
-			if (attr[0].columns) {
-			    db_append_string(&stmt, attr[0].null_values);
-			}
-			else {
-			    sprintf(buf, ", null");
-			    db_append_string(&stmt, buf);
-			}
-		    }
 
-		    /* catb */
-		    if (j >= 0) {
-			if (attr[1].columns) {
-			    at = find_attr(&(attr[1]),
-					   Centr[area].cat[1]->cat[j]);
-			    if (!at)
-				G_fatal_error(_("Attribute not found"));
+			/* catb */
+			if (j >= 0) {
+			    if (attr[1].columns) {
+				at = find_attr(&(attr[1]),
+					       Centr[area].cat[1]->cat[j]);
+				if (!at)
+				    G_fatal_error(_("Attribute not found"));
 
-			    if (at->values)
-				db_append_string(&stmt, at->values);
-			    else
-				db_append_string(&stmt, attr[1].null_values);
+				if (at->values)
+				    db_append_string(&stmt, at->values);
+				else
+				    db_append_string(&stmt, attr[1].null_values);
+			    }
+			    else {
+				sprintf(buf, ", %d", Centr[area].cat[1]->cat[j]);
+				db_append_string(&stmt, buf);
+			    }
 			}
 			else {
-			    sprintf(buf, ", %d", Centr[area].cat[1]->cat[j]);
-			    db_append_string(&stmt, buf);
+			    if (attr[1].columns) {
+				db_append_string(&stmt, attr[1].null_values);
+			    }
+			    else {
+				sprintf(buf, ", null");
+				db_append_string(&stmt, buf);
+			    }
 			}
-		    }
-		    else {
-			if (attr[1].columns) {
-			    db_append_string(&stmt, attr[1].null_values);
-			}
-			else {
-			    sprintf(buf, ", null");
-			    db_append_string(&stmt, buf);
-			}
-		    }
 
-		    db_append_string(&stmt, " )");
+			db_append_string(&stmt, " )");
 
-		    G_debug(3, db_get_string(&stmt));
+			G_debug(3, db_get_string(&stmt));
 
-		    if (db_execute_immediate(driver, &stmt) != DB_OK)
-			G_warning(_("Unable to insert new record: '%s'"),
-				  db_get_string(&stmt));
+			if (db_execute_immediate(driver, &stmt) != DB_OK)
+			    G_warning(_("Unable to insert new record: '%s'"),
+				      db_get_string(&stmt));
+		    }
+		    out_cat++;
 		}
-		out_cat++;
 	    }
 	}
 
-	/* Add all cats from imput vectors */
-	if (ofield[1] > 0) {
+	/* Add all cats from input vectors */
+	if (ofield[1] > 0 && field[0] > 0) {
 	    for (i = 0; i < Centr[area].cat[0]->n_cats; i++) {
-		Vect_cat_set(Cats, ofield[1], Centr[area].cat[0]->cat[i]);
+		if (Centr[area].cat[0]->field[i] == field[0])
+		    Vect_cat_set(Cats, ofield[1], Centr[area].cat[0]->cat[i]);
 	    }
 	}
 
-	if (ofield[2] > 0) {
+	if (ofield[2] > 0 && field[1] > 0 && ofield[1] != ofield[2]) {
 	    for (i = 0; i < Centr[area].cat[1]->n_cats; i++) {
-		Vect_cat_set(Cats, ofield[2], Centr[area].cat[1]->cat[i]);
+		if (Centr[area].cat[1]->field[i] == field[1])
+		    Vect_cat_set(Cats, ofield[2], Centr[area].cat[1]->cat[i]);
 	    }
 	}
 
 	Vect_write_line(Out, GV_CENTROID, Points, Cats);
-
-	G_percent(area, nareas, 1);
     }
 
     /* Build topology and remove boundaries with area without centroid on both sides */
@@ -306,6 +414,10 @@
     }
 
     /* Delete boundaries */
+    G_set_verbose(0);
+    Vect_build_partial(Out, GV_BUILD_NONE);
+    Vect_build_partial(Out, GV_BUILD_BASE);
+    G_set_verbose(verbose);
     for (line = 1; line <= nlines; line++) {
 	if (Del[line])
 	    Vect_delete_line(Out, line);

Modified: grass/branches/develbranch_6/vector/v.overlay/description.html
===================================================================
--- grass/branches/develbranch_6/vector/v.overlay/description.html	2012-12-01 19:59:23 UTC (rev 54133)
+++ grass/branches/develbranch_6/vector/v.overlay/description.html	2012-12-01 20:03:35 UTC (rev 54134)
@@ -1,51 +1,85 @@
 <h2>DESCRIPTION</h2>
 
-<em>v.overlay</em> allows the user to overlay two vector area maps.
+<em>v.overlay</em> allows the user to overlay two vector maps. Features 
+in <em>ainput</em> can be lines or areas and are cut with areas in 
+<em>binput</em>. Simple clipping can be performed with the <em>and</em> 
+oerator.
+<p>
+If areas in <em>ainput</em> are overlaid with areas in <em>binput</em>, 
+it is sometimes necessary to snap areas of <em>binput</em> to those of 
+<em>ainput</em>, otherwise areas can go missing or many sliver areas 
+can be created. Snapping is enabled by default and can be disabled by 
+setting the <em>snap</em> option to a negative value. Recommended values 
+are between 0.00000001 and 0.0001. Using larger values for snapping can 
+have undesired side-effects, but may sometimes be necessary to get a 
+clean output (see example below). In general, it is recommended to start 
+with a small snapping threshold, gradually increasing the threshold until 
+the result is reasonably clean. Snapping modifies only boundaries in 
+binput, which are snapped to boundaries in ainput. Boundaries in ainput 
+are not modified.
 <!-- This is outdated 
 There are 3 links attached to features in output map, 
 <ul>
 <li><b>field 1</b>: link to the new table, new table has 3 columns
-    </ul>
+    <ul>
         <li><b>cat</b> - key column linking rows to features
         <li><b>cata</b> - category of <i>afield</i> from <i>ainput</i>
         <li><b>catb</b> - category of <i>bfield</i> from <i>binput</i>
-    <ul>
+    </ul>
 <li><b>field 2</b>: category of <i>afield</i> from <i>ainput</i>
 <li><b>field 3</b>: category of <i>bfield</i> from <i>binput</i>
 </ul>
 -->
-The resulting output map has a merged attribute-table. The origin column-names
-have a prefix (<em>a_</em> and <em>b_</em>) which results from the ainput- and binput-map.
+<p>
+If the first number of the <em>olayer</em> option is > 0, then the 
+resulting output map has a merged attribute table in the given layer 
+number. The original column names have a prefix (<em>a_</em> and 
+<em>b_</em>) corresponding to <em>ainput</em> and <em>binput</em> map.
+<p>
+If the second number of the <em>olayer</em> option is > 0, then the 
+categories of <em>ainput</em> in layer <em>alayer</em> are transfered to
+the output layer with the second number.
+<p>
+If the third number of the <em>olayer</em> option is > 0, then the 
+categories of <em>binput</em> in layer <em>blayer</em> are transfered to
+the output layer with the third number.
 
 <h2>NOTES</h2>
-Currently only areas are supported for the operators <em>or</em> and <em>xor</em>! See also <a href="v.select.html">v.select</a>.
+Currently only areas in <em>ainput</em> are supported for the operators 
+<em>or</em> and <em>xor</em>! See also <a href="v.select.html">v.select</a>.
 
-The operator defines what kind of operation will be done. Features are written to output,
-if the result of an operation 'ainput operator binput' is true.
+The operator defines what kind of operation will be done. Features are 
+written to output, if the result of an operation 'ainput operator binput' 
+is true.
 <p>
-Attributes of the tables from ainput and binput are joined into a new table
-linked to the output maps new cat-column. 
+If the first number of the <em>olayer</em> option is > 0, then attributes 
+of the tables from ainput and binput are joined into a new table linked 
+to the output map with a new cat column. 
+<p>
+If the second number of the <em>olayer</em> option is > 0, then the 
+attribute table of ainput is copied to the output map. 
+<p>
+If the third number of the <em>olayer</em> option is > 0, then the 
+attribute table of binput is copied to the output map. 
 
 <!-- This is outdated
-<p>
-<div class="code"><pre>
+<p><div class="code"><pre>
 v.db.connect map=outputmap table=ainput.dbf field=2
 v.db.connect map=outputmap table=binput.dbf field=3
 
 </pre></div>
 
-<p>
-<b>Attention:</b> Removing the output map will also delete all tables linked to
+<p><b>Attention:</b> Removing the output map will also delete all tables linked to
 it! Therefore it is advisable to copy tables from ainput and binput first and
 connect the copied tables to the output map.-->
 
-<h2>EXAMPLE POLYGON TO POLYGON UNION</h2>
+<h2>EXAMPLES</h2>
+
+<h4>Polygons overlaid with ploygons</h4>
 <div class="code"><pre>
 v.overlay ainput=lake binput=province output=lakeXprovince operator=or
 </pre></div>
 
-<h2>EXAMPLE POLYGON TO POLYGON UNION</h2>
-
 Polygon union of urban area and Census 2000 areas (North Carolina dataset):
 
 <div class="code"><pre>
@@ -89,24 +123,28 @@
 <img src="v_overlay_census_wake2000.png" alt="GRASS v.overlay: polygon to polygon union (input 2)" border=1>
 <img src="v_overlay_urban_census2000.png" alt="GRASS v.overlay: polygon to polygon union (result)" border=1>
 <br>
-<i>v.overlay: Polygon union of urban area and Census 2000 areas (North Carolina dataset)</i>
+<i>v.overlay: Polygon union (right) of urban area (left) and Census 2000 (middle) areas (North Carolina dataset)</i>
 </center>
 
+As can be seen by the resulting large number of centroids on boundaries, 
+the urban areas do not match exactly the Census 2000 areas. In this case 
+a clean result can be obtained by snapping with a threshold of 0.1 (1 cm).
 
-<h2>EXAMPLE LINE TO POLYGON CLIPPING</h2>
+<h4>Lines overlaid with polygons</h4>
 
 Using the North Carolina sample dataset, we clip the roads map to the area
-of city of Raleigh:
+of city of Raleigh, preserving road attributes in layer 1:
 <div class="code"><pre>
 g.region vect=zipcodes_wake
 
-# extract Raleigh city
+# extract Raleigh city:
 v.extract in=zipcodes_wake out=raleigh \
             where="ZIPNAME = 'RALEIGH'"
 
 # clip road network to city polygon:
 v.overlay ainput=roadsmajor atype=line binput=raleigh \
-            out=roadsmajor_raleigh operator=and
+            out=roadsmajor_raleigh operator=and \
+	    olayer=0,1,0
 </pre></div>
 
 <center>
@@ -124,10 +162,12 @@
 <em>
 <a href="v.db.connect.html">v.db.connect</a>,
 <a href="v.select.html">v.select</a>,
-<a href="g.copy.html">g.copy</a></em>
+<a href="g.copy.html">g.copy</a>
+</em>
 
 <h2>AUTHORS</h2>
 
-Radim Blazek, ITC-Irst, Trento, Italy
+Radim Blazek, ITC-Irst, Trento, Italy<br>
+Markus Metz
 
 <p><i>Last changed: $Date$</i>

Modified: grass/branches/develbranch_6/vector/v.overlay/line_area.c
===================================================================
--- grass/branches/develbranch_6/vector/v.overlay/line_area.c	2012-12-01 19:59:23 UTC (rev 54133)
+++ grass/branches/develbranch_6/vector/v.overlay/line_area.c	2012-12-01 20:03:35 UTC (rev 54134)
@@ -22,7 +22,7 @@
 	       struct line_cats *Cats)
 {
     int i, area, centr;
-    struct line_cats *CCats = NULL;
+    static struct line_cats *CCats = NULL;
 
     Vect_reset_cats(Cats);
     area = Vect_find_area(Map, x, y);
@@ -51,7 +51,7 @@
 
 int line_area(struct Map_info *In, int *field, struct Map_info *Out,
 	      struct field_info *Fi, dbDriver * driver, int operator,
-	      int *ofield, ATTRIBUTES * attr)
+	      int *ofield, ATTRIBUTES * attr, struct ilist *BList)
 {
     int line, nlines, ncat;
     struct line_pnts *Points;
@@ -67,11 +67,10 @@
     db_init_string(&stmt);
 
     G_message(_("Breaking lines..."));
-    Vect_break_lines(Out, GV_LINE | GV_BOUNDARY, NULL);
+    Vect_break_lines_list(Out, NULL, BList, GV_LINE | GV_BOUNDARY, NULL);
+    G_message(_("Merging lines..."));
+    Vect_merge_lines(Out, GV_LINE, NULL, NULL);
 
-    /* Basic topology needed only */
-    Vect_build_partial(Out, GV_BUILD_BASE);
-
     nlines = Vect_get_num_lines(Out);
 
     /* Warning!: cleaning process (break) creates new vertices which are usually slightly 
@@ -81,6 +80,7 @@
      */
 
     /* Check if the line is inside or outside binput area */
+    G_message(_("Selecting lines..."));
     ncat = 1;
     for (line = 1; line <= nlines; line++) {
 	int ltype;
@@ -131,107 +131,111 @@
 
 	    Vect_reset_cats(OCats);
 
-	    /* rewrite with all combinations of acat - bcat (-1 in cycle for null) */
-	    for (i = -1; i < Cats->n_cats; i++) {	/* line cats */
-		int j;
+	    if (ofield[0]  > 0) {
+		/* rewrite with all combinations of acat - bcat (-1 in cycle for null) */
+		for (i = -1; i < Cats->n_cats; i++) {	/* line cats */
+		    int j;
 
-		if (i == -1 && Cats->n_cats > 0)
-		    continue;	/* no need to make null */
-
-		for (j = -1; j < ACats->n_cats; j++) {
-		    if (j == -1 && ACats->n_cats > 0)
+		    if (i == -1 && Cats->n_cats > 0)
 			continue;	/* no need to make null */
 
-		    if (ofield[0] > 0)
-			Vect_cat_set(OCats, ofield[0], ncat);
+		    for (j = -1; j < ACats->n_cats; j++) {
+			if (j == -1 && ACats->n_cats > 0)
+			    continue;	/* no need to make null */
 
-		    /* Attributes */
-		    if (driver) {
-			ATTR *at;
+			if (ofield[0] > 0)
+			    Vect_cat_set(OCats, ofield[0], ncat);
 
-			sprintf(buf, "insert into %s values ( %d", Fi->table,
-				ncat);
-			db_set_string(&stmt, buf);
+			/* Attributes */
+			if (driver) {
+			    ATTR *at;
 
-			/* cata */
-			if (i >= 0) {
-			    if (attr[0].columns) {
-				at = find_attr(&(attr[0]), Cats->cat[i]);
-				if (!at)
-				    G_fatal_error(_("Attribute not found"));
+			    sprintf(buf, "insert into %s values ( %d", Fi->table,
+				    ncat);
+			    db_set_string(&stmt, buf);
 
-				if (at->values)
-				    db_append_string(&stmt, at->values);
-				else
-				    db_append_string(&stmt,
-						     attr[0].null_values);
+			    /* cata */
+			    if (i >= 0) {
+				if (attr[0].columns) {
+				    at = find_attr(&(attr[0]), Cats->cat[i]);
+				    if (!at)
+					G_fatal_error(_("Attribute not found"));
+
+				    if (at->values)
+					db_append_string(&stmt, at->values);
+				    else
+					db_append_string(&stmt,
+							 attr[0].null_values);
+				}
+				else {
+				    sprintf(buf, ", %d", Cats->cat[i]);
+				    db_append_string(&stmt, buf);
+				}
 			    }
 			    else {
-				sprintf(buf, ", %d", Cats->cat[i]);
-				db_append_string(&stmt, buf);
+				if (attr[0].columns) {
+				    db_append_string(&stmt, attr[0].null_values);
+				}
+				else {
+				    sprintf(buf, ", null");
+				    db_append_string(&stmt, buf);
+				}
 			    }
-			}
-			else {
-			    if (attr[0].columns) {
-				db_append_string(&stmt, attr[0].null_values);
-			    }
-			    else {
-				sprintf(buf, ", null");
-				db_append_string(&stmt, buf);
-			    }
-			}
 
-			/* catb */
-			if (j >= 0) {
-			    if (attr[1].columns) {
-				at = find_attr(&(attr[1]), ACats->cat[j]);
-				if (!at)
-				    G_fatal_error(_("Attribute not found"));
+			    /* catb */
+			    if (j >= 0) {
+				if (attr[1].columns) {
+				    at = find_attr(&(attr[1]), ACats->cat[j]);
+				    if (!at)
+					G_fatal_error(_("Attribute not found"));
 
-				if (at->values)
-				    db_append_string(&stmt, at->values);
-				else
-				    db_append_string(&stmt,
-						     attr[1].null_values);
+				    if (at->values)
+					db_append_string(&stmt, at->values);
+				    else
+					db_append_string(&stmt,
+							 attr[1].null_values);
+				}
+				else {
+				    sprintf(buf, ", %d", ACats->cat[j]);
+				    db_append_string(&stmt, buf);
+				}
 			    }
 			    else {
-				sprintf(buf, ", %d", ACats->cat[j]);
-				db_append_string(&stmt, buf);
+				if (attr[1].columns) {
+				    db_append_string(&stmt, attr[1].null_values);
+				}
+				else {
+				    sprintf(buf, ", null");
+				    db_append_string(&stmt, buf);
+				}
 			    }
-			}
-			else {
-			    if (attr[1].columns) {
-				db_append_string(&stmt, attr[1].null_values);
-			    }
-			    else {
-				sprintf(buf, ", null");
-				db_append_string(&stmt, buf);
-			    }
-			}
 
-			db_append_string(&stmt, " )");
+			    db_append_string(&stmt, " )");
 
-			G_debug(3, db_get_string(&stmt));
+			    G_debug(3, db_get_string(&stmt));
 
-			if (db_execute_immediate(driver, &stmt) != DB_OK)
-			    G_warning(_("Unable to insert new record: '%s'"),
-				      db_get_string(&stmt));
+			    if (db_execute_immediate(driver, &stmt) != DB_OK)
+				G_warning(_("Unable to insert new record: '%s'"),
+					  db_get_string(&stmt));
+			}
+
+			ncat++;
 		    }
-
-		    ncat++;
 		}
 	    }
 
-	    /* Add all cats from imput vectors */
-	    if (ofield[1] > 0) {
+	    /* Add cats from input vectors */
+	    if (ofield[1] > 0 && field[0] > 0) {
 		for (i = 0; i < Cats->n_cats; i++) {
-		    Vect_cat_set(OCats, ofield[1], Cats->cat[i]);
+		    if (Cats->field[i] == field[0])
+			Vect_cat_set(OCats, ofield[1], Cats->cat[i]);
 		}
 	    }
 
-	    if (ofield[2] > 0) {
+	    if (ofield[2] > 0 && field[1] > 0 && ofield[1] != ofield[2]) {
 		for (i = 0; i < ACats->n_cats; i++) {
-		    Vect_cat_set(OCats, ofield[2], ACats->cat[i]);
+		    if (Cats->field[i] == field[1])
+			Vect_cat_set(OCats, ofield[2], ACats->cat[i]);
 		}
 	    }
 

Modified: grass/branches/develbranch_6/vector/v.overlay/local.h
===================================================================
--- grass/branches/develbranch_6/vector/v.overlay/local.h	2012-12-01 19:59:23 UTC (rev 54133)
+++ grass/branches/develbranch_6/vector/v.overlay/local.h	2012-12-01 20:03:35 UTC (rev 54134)
@@ -32,7 +32,7 @@
 
 int area_area(struct Map_info *In, int *field, struct Map_info *Out,
 	      struct field_info *Fi, dbDriver * driver, int operator,
-	      int *ofield, ATTRIBUTES * attr);
+	      int *ofield, ATTRIBUTES * attr, struct ilist *BList, double snap_thresh);
 int line_area(struct Map_info *In, int *field, struct Map_info *Out,
 	      struct field_info *Fi, dbDriver * driver, int operator,
-	      int *ofield, ATTRIBUTES * attr);
+	      int *ofield, ATTRIBUTES * attr, struct ilist *BList);

Modified: grass/branches/develbranch_6/vector/v.overlay/main.c
===================================================================
--- grass/branches/develbranch_6/vector/v.overlay/main.c	2012-12-01 19:59:23 UTC (rev 54133)
+++ grass/branches/develbranch_6/vector/v.overlay/main.c	2012-12-01 20:03:35 UTC (rev 54134)
@@ -8,6 +8,7 @@
  *               Jachym Cepicky <jachym les-ejk.cz>,
  *               Markus Neteler <neteler itc.it>,
  *               Paul Kelly <paul-grass stjohnspoint.co.uk>
+ *               Markus Metz
  * PURPOSE:      
  * COPYRIGHT:    (C) 2003-2008 by the GRASS Development Team
  *
@@ -16,6 +17,7 @@
  *               for details.
  *
  *****************************************************************************/
+
 #include <stdlib.h>
 #include <string.h>
 #include <stdio.h>
@@ -29,15 +31,18 @@
 {
     int i, input, line, nlines, operator;
     int type[2], field[2], ofield[3];
+    double snap_thresh;
     char *mapset[2];
     char *pre[2];
     struct GModule *module;
     struct Option *in_opt[2], *out_opt, *type_opt[2], *field_opt[2],
-	*ofield_opt, *operator_opt;
+	*ofield_opt, *operator_opt, *snap_opt;
     struct Flag *table_flag;
     struct Map_info In[2], Out;
-    struct line_pnts *Points;
+    struct line_pnts *Points, *Points2;
     struct line_cats *Cats;
+    struct ilist *BList;
+    int verbose;
 
     struct field_info *Fi = NULL;
     char buf[1000];
@@ -56,7 +61,7 @@
     module->description = _("Overlays two vector maps.");
 
     in_opt[0] = G_define_standard_option(G_OPT_V_INPUT);
-    in_opt[0]->description = _("Name of input vector map (A)");
+    in_opt[0]->label = _("Name of input vector map (A)");
     in_opt[0]->key = "ainput";
 
     type_opt[0] = G_define_standard_option(G_OPT_V_TYPE);
@@ -70,7 +75,7 @@
     field_opt[0]->key = "alayer";
 
     in_opt[1] = G_define_standard_option(G_OPT_V_INPUT);
-    in_opt[1]->description = _("Name of input vector map (B)");
+    in_opt[1]->label = _("Name of input vector map (B)");
     in_opt[1]->key = "binput";
 
     type_opt[1] = G_define_standard_option(G_OPT_V_TYPE);
@@ -115,6 +120,13 @@
     ofield_opt->description = _("If 0 or not given, "
 				"the category is not written");
 
+    snap_opt = G_define_option();
+    snap_opt->key = "snap";
+    snap_opt->label = _("Snapping threshold for boundaries");
+    snap_opt->description = _("Disable snapping with snap <= 0");
+    snap_opt->type = TYPE_DOUBLE;
+    snap_opt->answer = "1e-8";
+
     table_flag = G_define_flag();
     table_flag->key = 't';
     table_flag->description = _("Do not create attribute table");
@@ -126,6 +138,8 @@
 	type[input] = Vect_option_to_types(type_opt[input]);
 	field[input] = atoi(field_opt[input]->answer);
     }
+    if (type[0] & GV_AREA)
+	type[0] = GV_AREA;
 
     ofield[0] = ofield[1] = ofield[2] = 0;
     i = 0;
@@ -136,15 +150,14 @@
 
     if (operator_opt->answer[0] == 'a')
 	operator = OP_AND;
-
     else if (operator_opt->answer[0] == 'o')
 	operator = OP_OR;
-
     else if (operator_opt->answer[0] == 'n')
 	operator = OP_NOT;
-
     else if (operator_opt->answer[0] == 'x')
 	operator = OP_XOR;
+    else
+	G_fatal_error(_("Unknown operator '%s'"), operator_opt->answer);
 
     /* OP_OR, OP_XOR is not supported for lines,
        mostly because I'am not sure if they make enouhg sense */
@@ -157,7 +170,10 @@
     Vect_check_input_output_name(in_opt[1]->answer, out_opt->answer,
 				 GV_FATAL_EXIT);
 
+    snap_thresh = atof(snap_opt->answer);
+
     Points = Vect_new_line_struct();
+    Points2 = Vect_new_line_struct();
     Cats = Vect_new_cats_struct();
 
     /* Open output */
@@ -193,8 +209,13 @@
     }
 
     /* Copy lines to output */
+    BList = Vect_new_list();
+    verbose = G_verbose();
+    G_set_verbose(0);
+    Vect_build_partial(&Out, GV_BUILD_BASE);
+    G_set_verbose(verbose);
     for (input = 0; input < 2; input++) {
-	int ncats, index, nlines_out;
+	int ncats, index, nlines_out, newline;
 
 	G_message(_("Copying vector objects from vector map <%s>..."),
 		  in_opt[input]->answer);
@@ -213,6 +234,7 @@
 	nlines_out = 0;
 	for (line = 1; line <= nlines; line++) {
 	    int ltype;
+	    int vertices = 100; /* max number of vertices per line */
 
 	    G_percent(line, nlines, 1);	/* must be before any continue */
 
@@ -226,8 +248,44 @@
 		if (!(ltype & type[input]))
 		    continue;
 	    }
+	    
+	    /* lines and boundaries must have at least 2 distinct vertices */
+	    Vect_line_prune(Points);
+	    if (Points->n_points < 2)
+		continue;
 
-	    Vect_write_line(&Out, ltype, Points, Cats);
+	    /* TODO: figure out a reasonable threshold */
+	    if (Points->n_points > vertices) {
+		int start = 0;	/* number of coordinates written */
+		
+		vertices = Points->n_points / (Points->n_points / vertices + 1);
+		
+		/* split */
+		while (start < Points->n_points - 1) {
+		    int v = 0;
+
+		    Vect_reset_line(Points2);
+		    for (i = 0; i < vertices; i++) {
+			v = start + i;
+			if (v == Points->n_points)
+			    break;
+
+			Vect_append_point(Points2, Points->x[v], Points->y[v],
+					  Points->z[v]);
+		    }
+
+		    newline = Vect_write_line(&Out, ltype, Points2, Cats);
+		    if (input == 1)
+			Vect_list_append(BList, newline);
+
+		    start = v;
+		}
+	    }
+	    else {
+		newline = Vect_write_line(&Out, ltype, Points, Cats);
+		if (input == 1)
+		    Vect_list_append(BList, newline);
+	    }
 	    nlines_out++;
 	}
 	if (nlines_out == 0) {
@@ -328,6 +386,9 @@
 		    sprintf(buf, "varchar(%d)", db_get_column_length(Column));
 		    db_append_string(&col_defs, buf);
 		    break;
+		case DB_SQL_TYPE_TEXT:
+		    db_append_string(&col_defs, "varchar(250)");
+		    break;
 		case DB_SQL_TYPE_SMALLINT:
 		case DB_SQL_TYPE_INTEGER:
 		    db_append_string(&col_defs, "integer");
@@ -349,8 +410,8 @@
 		    db_append_string(&col_defs, "datetime");
 		    break;
 		default:
-		    G_warning(_("Unknown column type '%s'"),
-			      db_get_column_name(Column));
+		    G_warning(_("Unknown column type '%s' of column '%s'"),
+			      db_sqltype_name(sqltype), db_get_column_name(Column));
 		    sprintf(buf, "varchar(250)");
 		}
 	    }
@@ -358,7 +419,7 @@
 	    attr[input].columns = G_store(db_get_string(&col_defs));
 
 	    while (1) {
-		int cat;
+		int cat = -1;
 		ATTR *at;
 
 		if (db_fetch(&cursor, DB_NEXT, &more) != DB_OK)
@@ -414,8 +475,8 @@
 			}
 			break;
 		    default:
-			G_warning(_("Unknown column type '%s' values lost"),
-				  db_get_column_name(Column));
+			G_warning(_("Unknown column type '%s' of column '%s', values lost"),
+			      db_sqltype_name(sqltype), db_get_column_name(Column));
 			db_append_string(&sql, "null");
 		    }
 		}
@@ -479,16 +540,12 @@
 			    Fi->database, Fi->driver);
     }
 
-    G_message(_("Building partial topology..."));
-    /* do not print output, because befor cleaning it is nonsense */
-    Vect_build_partial(&Out, GV_BUILD_BASE);
-
     /* AREA x AREA */
     if (type[0] == GV_AREA) {
-	area_area(In, field, &Out, Fi, driver, operator, ofield, attr);
+	area_area(In, field, &Out, Fi, driver, operator, ofield, attr, BList, snap_thresh);
     }
     else {			/* LINE x AREA */
-	line_area(In, field, &Out, Fi, driver, operator, ofield, attr);
+	line_area(In, field, &Out, Fi, driver, operator, ofield, attr, BList);
     }
 
     G_message(_("Rebuilding topology..."));
@@ -501,6 +558,23 @@
 	db_commit_transaction(driver);
 	db_close_database_shutdown_driver(driver);
     }
+    if (ofield[0] < 1 && !table_flag->answer) {
+	int otype;
+	
+	if (type[0] == GV_AREA)
+	    otype = GV_CENTROID;
+	else
+	    otype = GV_LINE;
+	
+	/* copy attributes from ainput */
+	if (ofield[1] > 0 && field[0] > 0) {
+	    Vect_copy_table(&In[0], &Out, field[0], ofield[1], NULL, otype);
+	}
+	/* copy attributes from binput */
+	if (ofield[2] > 0 && field[1] > 0 && ofield[1] != ofield[2]) {
+	    Vect_copy_table(&In[1], &Out, field[1], ofield[2], NULL, otype);
+	}
+    }
 
     Vect_close(&(In[0]));
     Vect_close(&(In[1]));



More information about the grass-commit mailing list