v.out.vtk: port -n flag to trunk

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

vector/v.out.vtk/local_proto.h

@@ -23,7 +23,7 @@ extern double y_extent;
 
 /*Write the vtk output */
 int write_vtk(FILE * ascii, struct Map_info *, int layer, int *types,
-	      int typenum, int dp, double scale);
+	      int typenum, int dp, double scale, int numatts, int labels );
 /*Write the VTK header */
 int write_vtk_head(FILE * ascii, struct Map_info *Map);
 

vector/v.out.vtk/main.c

@@ -32,12 +32,12 @@ int main(int argc, char *argv[])
 {
     FILE *ascii;
     struct Option *input, *output, *type_opt, *dp_opt, *layer_opt, *scale;
-    struct Flag *coorcorr;
-    int *types = NULL, typenum = 0, dp, i;
+    struct Flag *coorcorr, *numatts, *labels;
+    int itype, *types = NULL, typenum = 0, dp, i;
     struct Map_info Map;
+    struct bound_box box;
     struct GModule *module;
-    int layer;
-    struct Cell_head region;
+    int layer, level;
     double zscale = 1.0, llscale = 1.0;
 
 
@@ -84,11 +84,23 @@ int main(int argc, char *argv[])
     layer_opt->answer = "1";
     layer_opt->description = _("Layer number");
 
-    coorcorr = G_define_flag();
-    coorcorr->key = 'c';
-    coorcorr->description =
-	_("Correct the coordinates to fit the VTK-OpenGL precision");
-
+    coorcorr = G_define_flag();                                            
+    coorcorr->key = 'c';                                                   
+    coorcorr->description = 
+    	_("Correct the coordinates to fit the VTK-OpenGL precision");
+    
+    numatts = G_define_flag();
+    numatts->key = 'n';
+    numatts->description = 
+    	_("Export numeric attribute table fields as VTK scalar variables");
+
+    labels = NULL; /* to avoid compiler warning about "unused variable"*/
+    /* not yet supported
+    labels = G_define_flag();
+    labels->key = 'l';
+    labels->description = _("Export text attribute table fields as VTK labels");
+    */
+    
     if (G_parser(argc, argv))
 	exit(EXIT_FAILURE);
 
@@ -130,47 +142,17 @@ int main(int argc, char *argv[])
 	}
 	i++;
     }
-
-    G_get_set_window(&region);
-
-    /*Correct the coordinates, so the precision of VTK is not hurt :( */
-    if (coorcorr->answer) {
-	/*Get the default region for coordiante correction */
-	G_get_default_window(&region);
-
-	/*Use the center of the current region as extent */
-	y_extent = (region.north + region.south) / 2;
-	x_extent = (region.west + region.east) / 2;
-    }
-    else {
-	x_extent = 0;
-	y_extent = 0;
-    }
-
+    itype = Vect_option_to_types(type_opt);
 
     /* read and compute the scale factor */
     sscanf(scale->answer, "%lf", &zscale);
     /*if LL projection, convert the elevation values to degrees */
-    if (region.proj == PROJECTION_LL) {
+    if (G_projection() == PROJECTION_LL) {
 	llscale = M_PI / (180) * 6378137;
 	zscale /= llscale;
 	printf("Scale %g\n", zscale);
     }
 
-    /*We need level 2 functions */
-    Vect_set_open_level(2);
-    Vect_open_old(&Map, input->answer, "");
-
-    if (output->answer) {
-	ascii = fopen(output->answer, "w");
-	if (ascii == NULL) {
-	    G_fatal_error(_("Unable to open file <%s>"), output->answer);
-	}
-    }
-    else {
-	ascii = stdout;
-    }
-
     /*The precision of the output */
     if (dp_opt->answer) {
 	if (sscanf(dp_opt->answer, "%d", &dp) != 1)
@@ -191,10 +173,77 @@ int main(int argc, char *argv[])
 	layer = 1;
     }
 
+    if (output->answer) {
+	ascii = fopen(output->answer, "w");
+	if (ascii == NULL) {
+	    G_fatal_error(_("Unable to open file <%s>"), output->answer);
+	}
+    }
+    else {
+	ascii = stdout;
+    }
+
+    /* Open input vector */
+    level = Vect_open_old(&Map, input->answer, "");
+    if (level < 2 && (itype & GV_AREA))
+	G_fatal_error(_("Export of areas requires topology. "
+	                "Please adjust '%s' option or rebuild topology."),
+			type_opt->key);
+
+    if (level == 2)
+	Vect_get_map_box(&Map, &box);
+    else {
+	int i, type, first = TRUE;
+	struct line_pnts *Points = Vect_new_line_struct();
+
+	Vect_rewind(&Map);
+	while ((type = Vect_read_next_line(&Map, Points, NULL)) > 0) {
+
+	    if (first) {
+		box.E = box.W = Points->x[0];
+		box.N = box.S = Points->y[0];
+		box.B = box.T = Points->z[0];
+		first = FALSE;
+	    }
+	    for (i = 1; i < Points->n_points; i++) {
+		if (Points->x[i] > box.E)
+		    box.E = Points->x[i];
+		else if (Points->x[i] < box.W)
+		    box.W = Points->x[i];
+
+		if (Points->y[i] > box.N)
+		    box.N = Points->y[i];
+		else if (Points->y[i] < box.S)
+		    box.S = Points->y[i];
+
+		if (Points->z[i] > box.T)
+		    box.T = Points->z[i];
+		else if (Points->z[i] < box.B)
+		    box.B = Points->z[i];
+	    }
+	}
+	Vect_destroy_line_struct(Points);
+    }
+
+    /*Correct the coordinates, so the precision of VTK is not hurt :( */
+    if (coorcorr->answer) {
+
+	/*Use the center of the vector's bbox as extent */
+	y_extent = (box.N + box.S) / 2;
+	x_extent = (box.W + box.E) / 2;
+    }
+    else {
+	x_extent = 0;
+	y_extent = 0;
+    }
+
     /*Write the header */
     write_vtk_head(ascii, &Map);
     /*Write the geometry and data */
-    write_vtk(ascii, &Map, layer, types, typenum, dp, zscale);
+    write_vtk(ascii, &Map, layer, types, typenum, dp, zscale, numatts->answer, 0 );
+    /* change to this, when labels get supported:
+    write_vtk(ascii, &Map, layer, types, typenum, dp, zscale, numatts->answer, labels->answer );
+    */
 
     if (ascii != NULL)
 	fclose(ascii);

vector/v.out.vtk/writeVTK.c

@@ -16,6 +16,7 @@
 
 #include <stdlib.h>
 #include <grass/vector.h>
+#include <grass/dbmi.h>
 #include <grass/gis.h>
 #include <grass/glocale.h>
 #include "writeVTK.h"
@@ -44,7 +45,7 @@ int write_vtk_points(FILE * ascii, struct Map_info *Map, VTKInfo * info,
     Points = Vect_new_line_struct();	/* init line_pnts struct */
     Cats = Vect_new_cats_struct();
 
-    G_message("Writing the coordinates");
+    G_message("Writing coordinates ...");
 
     /*For every available vector type */
     for (k = 0; k < typenum; k++) {
@@ -54,8 +55,19 @@ int write_vtk_points(FILE * ascii, struct Map_info *Map, VTKInfo * info,
 
 	    /*Get the number of the points to generate */
 	    info->typeinfo[types[k]]->pointoffset = pointoffset;
-	    info->typeinfo[types[k]]->numpoints =
-		Vect_get_num_primitives(Map, types[k]);
+
+	    /*count the number of line_nodes and lines */
+	    Vect_rewind(Map);
+	    while (1) {
+		if (-1 == (type = Vect_read_next_line(Map, Points, Cats)))
+		    break;
+		if (type == -2)	/* EOF */
+		    break;
+		if (type == types[k]) {
+		    info->typeinfo[types[k]]->numpoints++;
+		}
+	    }
+
 	    pointoffset += info->typeinfo[types[k]]->numpoints;
 
 	    info->typeinfo[types[k]]->numvertices =
@@ -102,7 +114,6 @@ int write_vtk_points(FILE * ascii, struct Map_info *Map, VTKInfo * info,
 	     * info->typeinfo[types[k]]->lineoffset);
 	     */
 	}
-
     }
 
     for (k = 0; k < typenum; k++) {
@@ -137,7 +148,6 @@ int write_vtk_points(FILE * ascii, struct Map_info *Map, VTKInfo * info,
 	     * info->typeinfo[types[k]]->polygonoffset);
 	     */
 	}
-
     }
 
     for (k = 0; k < typenum; k++) {
@@ -185,10 +195,8 @@ int write_vtk_points(FILE * ascii, struct Map_info *Map, VTKInfo * info,
     if (info->maxnumpoints == 0)
 	G_fatal_error(_("No coordinates to generate the output! Maybe an empty vector type chosen?"));
 
-
     /************************************************/
     /*Write the coordinates into the vtk ascii file */
-
     /************************************************/
 
     fprintf(ascii, "POINTS %i float\n", info->maxnumpoints);
@@ -217,7 +225,6 @@ int write_vtk_points(FILE * ascii, struct Map_info *Map, VTKInfo * info,
 		}
 		cur++;
 	    }
-
 	}
     }
 
@@ -238,7 +245,6 @@ int write_vtk_points(FILE * ascii, struct Map_info *Map, VTKInfo * info,
 		}
 		cur++;
 	    }
-
 	}
     }
 
@@ -259,7 +265,6 @@ int write_vtk_points(FILE * ascii, struct Map_info *Map, VTKInfo * info,
 		}
 		cur++;
 	    }
-
 	}
     }
 
@@ -275,7 +280,6 @@ int write_vtk_points(FILE * ascii, struct Map_info *Map, VTKInfo * info,
 		Vect_get_area_points(Map, i, Points);
 		write_point_coordinates(Points, dp, scale, ascii);
 	    }
-
 	}
     }
 
@@ -298,7 +302,7 @@ int write_vtk_cells(FILE * ascii, struct Map_info *Map, VTKInfo * info,
     int linekeyword = 1;
     int polykeyword = 1;
 
-    G_message("Writing vtk cells");
+    G_message("Writing vtk cells ...");
 
     Points = Vect_new_line_struct();	/* init line_pnts struct */
     Cats = Vect_new_cats_struct();
@@ -306,7 +310,6 @@ int write_vtk_cells(FILE * ascii, struct Map_info *Map, VTKInfo * info,
     /*For every available vector type */
     for (k = 0; k < typenum; k++) {
 
-
 	/*POINT KERNEL CENTROID */
 	if (types[k] == GV_POINT || types[k] == GV_KERNEL ||
 	    types[k] == GV_CENTROID) {
@@ -465,7 +468,7 @@ int write_vtk_cat_data(FILE * ascii, struct Map_info *Map, VTKInfo * info,
     Points = Vect_new_line_struct();	/* init line_pnts struct */
     Cats = Vect_new_cats_struct();
 
-    G_message("Writing category cell data");
+    G_message("Writing category cell data ...");
 
     if (numcelldata > 0) {
 	/*Write the pointdata */
@@ -542,25 +545,329 @@ int write_vtk_cat_data(FILE * ascii, struct Map_info *Map, VTKInfo * info,
 		}
 	    }
 	}
+	fprintf(ascii, "\n");
     }
 
-
     return 1;
 }
 
+
+/* 
+   Reads the attribute field "name" for current cat and returns the value as a string 
+   or NULL on error 
+
+   Memory for the result string is allocated by this function and must be free'd by
+   the caller.
+ */
+char *get_att(char *name, int cat, struct field_info *Fi, dbDriver * Driver,
+	      int ncol)
+{
+    char buf[2000];
+    int more;
+    dbTable *Table;
+    static dbString dbstring;
+    dbColumn *Column;
+    dbCursor cursor;
+    char *retval;
+    static int first = 1;
+
+    if (first) {
+	db_init_string(&dbstring);
+	first = 0;
+    }
+
+    sprintf(buf, "SELECT %s FROM %s WHERE %s = %d", name, Fi->table, Fi->key, cat);
+
+    db_set_string(&dbstring, buf);
+
+    if (db_open_select_cursor(Driver, &dbstring, &cursor, DB_SEQUENTIAL) !=
+	DB_OK) {
+	G_fatal_error(_("Cannot select attributes for cat = %d"), cat);
+    }
+
+    if (db_fetch(&cursor, DB_NEXT, &more) != DB_OK)
+	G_fatal_error(_("Unable to fetch data from table"));
+
+    Table = db_get_cursor_table(&cursor);
+
+    Column = db_get_table_column(Table, 0);
+    if (!strcmp(name, db_get_column_name(Column))) {
+	db_convert_column_value_to_string(Column, &dbstring);
+	retval =
+	    G_malloc(sizeof(char) *
+		     (strlen(db_get_string(&dbstring)) + 1));
+	retval = G_store(db_get_string(&dbstring));
+	db_close_cursor(&cursor);
+	return (retval);
+    }
+
+    db_close_cursor(&cursor);
+    return (NULL);
+}
+
+
+/* ************************************************************************* */
+/* This function writes numerical attribute table fields as VTK scalars **** */
+/* ************************************************************************* */
 int write_vtk_db_data(FILE * ascii, struct Map_info *Map, VTKInfo * info,
 		      int layer, int *types, int typenum, int dp)
 {
-    G_message("Writing database cell/point data");
+    int type, cat, i, k, centroid;
+    struct line_pnts *Points;
+    struct line_cats *Cats;
+
+    /*The keywords may only be written once! */
+    int numcelldata =
+	info->maxnumvertices + info->maxnumlines + info->maxnumpolygons;
+    /* attribute table info */
+    int ncol = 0, colsqltype, colctype, num_atts, cur_att, progress;
+    struct field_info *Fi = NULL;
+    dbDriver *Driver = NULL;
+    dbHandle handle;
+    dbTable *Table;
+    dbString dbstring;
+    dbColumn *Column;
+    char *valbuf;
+
+    if (layer < 1) {
+	G_warning(_("Cannot export attribute table fields for layer < 1. Skipping export"));
+	return 1;
+    }
+
+    /* attempt to open attribute table for selected layer */
+    db_init_string(&dbstring);
+    Fi = Vect_get_field(Map, layer);
+    if (Fi == NULL) {
+	G_fatal_error(_("No attribute table found"));
+    }
+    Driver = db_start_driver(Fi->driver);
+    if (Driver == NULL)
+	G_fatal_error(_("Unable to start driver <%s>"), Fi->driver);
+
+    db_init_handle(&handle);
+    db_set_handle(&handle, Fi->database, NULL);
+    if (db_open_database(Driver, &handle) != DB_OK)
+	G_fatal_error(_("Unable to open database <%s> by driver <%s>"),
+		      Fi->database, Fi->driver);
+
+    db_set_string(&dbstring, Fi->table);
+    if (db_describe_table(Driver, &dbstring, &Table) != DB_OK)
+	G_fatal_error(_("Unable to describe table <%s>"), Fi->table);
+
+    /* analyse field structure */
+    ncol = db_get_table_number_of_columns(Table);
+    num_atts = 0;
+    for (i = 0; i < ncol; i++) {
+	Column = db_get_table_column(Table, i);
+	colsqltype = db_get_column_sqltype(Column);
+	colctype = db_sqltype_to_Ctype(colsqltype);
+	if ((colctype == DB_C_TYPE_INT) || (colctype == DB_C_TYPE_DOUBLE)) {
+	    /* we don't want to export the category field twice */
+	    if (strcmp(db_get_column_name(Column), "cat")) {
+		num_atts++;
+		/* fprintf ( stderr, "%i: %s\n", num_atts, db_get_column_name(Column) ); */
+	    }
+	}
+    }
+    if (num_atts < 1) {
+	G_warning(_("No numerical attributes found. Skipping export"));
+	db_close_database(Driver);
+	db_shutdown_driver(Driver);
+	return 1;
+    }
+
+    Points = Vect_new_line_struct();	/* init line_pnts struct */
+    Cats = Vect_new_cats_struct();
+
+    G_message("Writing %i scalar variables as cell data ...", num_atts);
+
+    progress = 0;
+    for (cur_att = 0; cur_att < ncol; cur_att++) {
+
+	if (numcelldata > 0) {
+	    /*Write the pointdata */
+	    Column = db_get_table_column(Table, cur_att);
+	    colsqltype = db_get_column_sqltype(Column);
+	    colctype = db_sqltype_to_Ctype(colsqltype);
+
+	    if ((strcmp("cat", db_get_column_name(Column))) &&
+		((colctype == DB_C_TYPE_INT) ||
+		 (colctype == DB_C_TYPE_DOUBLE))) {
+
+		if (colctype == DB_C_TYPE_INT) {
+		    /* G_message("  Writing integer scalar %s", db_get_column_name(Column) ); */
+		    fprintf(ascii, "SCALARS %s int 1\n",
+			    db_get_column_name(Column));
+		}
+		if (colctype == DB_C_TYPE_DOUBLE) {
+		    /* *G_message("  Writing double scalar %s", db_get_column_name(Column) ); */
+		    fprintf(ascii, "SCALARS %s double 1\n",
+			    db_get_column_name(Column));
+		}
+
+		fprintf(ascii, "LOOKUP_TABLE default\n");
+		progress++;
+
+		/*For every available vector type */
+		for (k = 0; k < typenum; k++) {
+		    /*POINT KERNEL CENTROID */
+		    if (types[k] == GV_POINT || types[k] == GV_KERNEL ||
+			types[k] == GV_CENTROID) {
+
+			Vect_rewind(Map);
+
+			while (1) {
+			    if (-1 ==
+				(type =
+				 Vect_read_next_line(Map, NULL, Cats)))
+				break;
+			    if (type == -2)	/* EOF */
+				break;
+			    if (type == types[k]) {
+				Vect_cat_get(Cats, layer, &cat);
+				valbuf =
+				    get_att((char *)
+					    db_get_column_name(Column), cat,
+					    Fi, Driver, ncol);
+				if (valbuf == NULL) {
+				    db_close_database(Driver);
+				    db_shutdown_driver(Driver);
+				    G_fatal_error(_("Error reading value of attribute '%s'"),
+						  db_get_column_name(Column));
+				}
+				/* DEBUG 
+				   fprintf ( stderr, "%s (%i) = %s\n", db_get_column_name(Column), cat, valbuf );
+				 */
+				fprintf(ascii, " %s", valbuf);
+				G_free(valbuf);
+			    }
+			}
+		    }
+		}
+
+		for (k = 0; k < typenum; k++) {
+		    /*LINE BOUNDARY */
+		    if (types[k] == GV_LINE || types[k] == GV_BOUNDARY) {
+			Vect_rewind(Map);
+			while (1) {
+			    if (-1 ==
+				(type =
+				 Vect_read_next_line(Map, NULL, Cats)))
+				break;
+			    if (type == -2)	/* EOF */
+				break;
+			    if (type == types[k]) {
+				Vect_cat_get(Cats, layer, &cat);
+				valbuf =
+				    get_att((char *)
+					    db_get_column_name(Column), cat,
+					    Fi, Driver, ncol);
+				if (valbuf == NULL) {
+				    db_close_database(Driver);
+				    db_shutdown_driver(Driver);
+				    G_fatal_error(_("Error reading value of attribute '%s'"),
+						  db_get_column_name(Column));
+				}
+				/* DEBUG 
+				   fprintf ( stderr, "%s (%i) = %s\n", db_get_column_name(Column), cat, valbuf );
+				 */
+				fprintf(ascii, " %s", valbuf);
+				G_free(valbuf);
+			    }
+			}
+		    }
+		}
+
+		for (k = 0; k < typenum; k++) {
+		    /*FACE */
+		    if (types[k] == GV_FACE) {
+			Vect_rewind(Map);
+			while (1) {
+			    if (-1 ==
+				(type =
+				 Vect_read_next_line(Map, NULL, Cats)))
+				break;
+			    if (type == -2)	/* EOF */
+				break;
+			    if (type == types[k]) {
+				Vect_cat_get(Cats, layer, &cat);
+				valbuf =
+				    get_att((char *)
+					    db_get_column_name(Column), cat,
+					    Fi, Driver, ncol);
+				if (valbuf == NULL) {
+				    db_close_database(Driver);
+				    db_shutdown_driver(Driver);
+				    G_fatal_error(_("Error reading value of attribute '%s'"),
+						  db_get_column_name(Column));
+				}
+				/* DEBUG 
+				   fprintf ( stderr, "%s (%i) = %s\n", db_get_column_name(Column), cat, valbuf );
+				 */
+				fprintf(ascii, " %s", valbuf);
+				G_free(valbuf);
+			    }
+			}
+		    }
+		}
+
+		for (k = 0; k < typenum; k++) {
+		    /*AREA */
+		    if (types[k] == GV_AREA) {
+			Vect_rewind(Map);
+			for (i = 1;
+			     i <= info->typeinfo[types[k]]->numpolygons;
+			     i++) {
+			    centroid = Vect_get_area_centroid(Map, i);
+			    if (centroid > 0) {
+				Vect_read_line(Map, NULL, Cats, centroid);
+			    }
+			    Vect_cat_get(Cats, layer, &cat);
+			    valbuf =
+				get_att((char *)db_get_column_name(Column),
+					cat, Fi, Driver, ncol);
+			    if (valbuf == NULL) {
+				db_close_database(Driver);
+				db_shutdown_driver(Driver);
+				G_fatal_error(_("Error reading value of attribute '%s'"),
+					      db_get_column_name(Column));
+			    }
+			    /* DEBUG 
+			       fprintf ( stderr, "%s (%i) = %s\n", db_get_column_name(Column), cat, valbuf );
+			     */
+			    fprintf(ascii, " %s", valbuf);
+			    G_free(valbuf);
+			}
+		    }
+		}
+		fprintf(ascii, "\n");
+	    }			/* END (do for all scalars != cat */
+	}
+    }				/* END (step through all numerical attributes) */
+    fprintf(stdout, "\n");
+    fflush(stdout);
+
+    db_close_database(Driver);
+    db_shutdown_driver(Driver);
 
     return 1;
 }
 
+
+/* ************************************************************************* */
+/* This function writes attribute table fields as VTK labels            **** */
+/* ************************************************************************* */
+int write_vtk_db_labels(FILE * ascii, struct Map_info *Map, VTKInfo * info,
+			int layer, int *types, int typenum, int dp)
+{
+    return 1;
+}
+
 /* ************************************************************************* */
 /* This function writes the point coordinates and the geometric feature **** */
 /* ************************************************************************* */
 int write_vtk(FILE * ascii, struct Map_info *Map, int layer, int *types,
-	      int typenum, int dp, double scale)
+	      int typenum, int dp, double scale, int numatts, int labels)
 {
     VTKInfo *info;
     VTKTypeInfo **typeinfo;
@@ -569,7 +876,6 @@ int write_vtk(FILE * ascii, struct Map_info *Map, int layer, int *types,
 	GV_POINT + GV_KERNEL + GV_CENTROID + GV_LINE + GV_BOUNDARY + GV_FACE +
 	GV_AREA;
 
-
     /*Initiate the typeinfo structure for every supported type */
     typeinfo = (VTKTypeInfo **) calloc(infonum, sizeof(VTKTypeInfo *));
     for (i = 0; i < infonum; i++) {
@@ -604,9 +910,15 @@ int write_vtk(FILE * ascii, struct Map_info *Map, int layer, int *types,
     /*3. write the cat data */
     write_vtk_cat_data(ascii, Map, info, layer, types, typenum, dp);
 
-    /*4. write the DB data */
-    /* not yet implemented
-       write_vtk_db_data(ascii, Map, info, layer, types, typenum, dp);
+    /*4. write the DB data: numerical attributes */
+    if (numatts) {
+	write_vtk_db_data(ascii, Map, info, layer, types, typenum, dp);
+    }
+
+    /*5. Write labels (not yet supported) 
+       if ( labels ) {
+       write_vtk_db_labels(ascii, Map, info, layer, types, typenum, dp);
+       }
      */
 
     /*Release the memory */