vlib(pg): speed up Vect_get_area|isle_points() for PostGIS Topology

access (work in progress)


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

lib/vector/Vlib/area.c

@@ -5,7 +5,7 @@
 
    Higher level functions for reading/writing/manipulating vectors.
 
-   (C) 2001-2009, 2011-2012 by the GRASS Development Team
+   (C) 2001-2009, 2011-2013 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.
@@ -41,14 +41,14 @@ int Vect_get_area_points(const struct Map_info *Map,
     struct P_area *Area;
 
     G_debug(3, "Vect_get_area_points(): area = %d", area);
-    BPoints->n_points = 0;
+    Vect_reset_line(BPoints);
 
     Plus = &(Map->plus);
     Area = Plus->Area[area];
 
     if (Area == NULL) {		/* dead area */
 	G_warning(_("Attempt to read points of nonexistent area"));
-	return -1;		/* error , because we should not read dead areas */
+	return -1;		/* error, because we should not read dead areas */
     }
 
     G_debug(3, "  n_lines = %d", Area->n_lines);
@@ -66,49 +66,35 @@ int Vect_get_area_points(const struct Map_info *Map,
    \return -1 on error
  */
 int Vect_get_isle_points(const struct Map_info *Map,
-			 int isle, struct line_pnts *BPoints)
+                         int isle, struct line_pnts *BPoints)
 {
-    int i, line, aline, dir;
     const struct Plus_head *Plus;
     struct P_isle *Isle;
-    static int first_time = 1;
-    static struct line_pnts *Points;
 
     G_debug(3, "Vect_get_isle_points(): isle = %d", isle);
-    BPoints->n_points = 0;
+    Vect_reset_line(BPoints);
 
     Plus = &(Map->plus);
     Isle = Plus->Isle[isle];
-
-    if (first_time == 1) {
-	Points = Vect_new_line_struct();
-	first_time = 0;
+    
+    if (Isle == NULL) {		/* dead isle */
+	G_warning(_("Attempt to read points of nonexistent isle"));
+	return -1;		/* error, because we should not read dead isles */
     }
 
     G_debug(3, "  n_lines = %d", Isle->n_lines);
-    for (i = 0; i < Isle->n_lines; i++) {
-	line = Isle->lines[i];
-	aline = abs(line);
-	G_debug(3, "  append line(%d) = %d", i, line);
 
-	if (0 > Vect_read_line(Map, Points, NULL, aline)) {
-	    G_fatal_error(_("Unable to read line %d"), aline);
-	}
-
-	G_debug(3, "  line n_points = %d", Points->n_points);
-
-	if (line > 0)
-	    dir = GV_FORWARD;
-	else
-	    dir = GV_BACKWARD;
-
-	Vect_append_points(BPoints, Points, dir);
-	if (i != (Isle->n_lines - 1))	/* all but not last */
-	    BPoints->n_points--;
-	G_debug(3, "  isle n_points = %d", BPoints->n_points);
+    if (Map->format == GV_FORMAT_POSTGIS &&
+        Map->fInfo.pg.toposchema_name) {
+#ifdef HAVE_POSTGRES
+        /* PostGIS Topology */
+        return Vect__get_area_points_pg(Map, Isle->lines, Isle->n_lines, BPoints);
+#else
+        G_fatal_error(_("GRASS is not compiled with PostgreSQL support"));
+#endif
     }
-
-    return (BPoints->n_points);
+    /* native format */
+    return Vect__get_area_points_nat(Map, Isle->lines, Isle->n_lines, BPoints);
 }
 
 /*!
@@ -452,6 +438,20 @@ int Vect_get_area_cat(const struct Map_info *Map, int area, int field)
     return -1;
 }
 
+/*!
+  \brief Get area boundary points (internal use only)
+  
+  For PostGIS Topology calls Vect__get_area_points_pg() otherwise
+  Vect__get_area_points_nat(),
+  
+  \param Map pointer to Map_info struct
+  \param lines array of boundary lines
+  \param n_lines number of lines in array
+  \param[out] APoints pointer to output line_pnts struct
+
+  \return number of points
+  \return -1 on error
+*/
 int Vect__get_area_points(const struct Map_info *Map, const plus_t *lines, int n_lines,
                           struct line_pnts *BPoints)
 {
@@ -467,3 +467,44 @@ int Vect__get_area_points(const struct Map_info *Map, const plus_t *lines, int n
     /* native format */
     return Vect__get_area_points_nat(Map, lines, n_lines, BPoints);
 }
+
+
+/*!
+  \brief Get area boundary points (native format)
+  
+  Used by Vect_build_line_area() and Vect_get_area_points().
+  
+  \param Map pointer to Map_info struct
+  \param lines array of boundary lines
+  \param n_lines number of lines in array
+  \param[out] APoints pointer to output line_pnts struct
+
+  \return number of points
+  \return -1 on error
+*/
+int Vect__get_area_points_nat(const struct Map_info *Map, const plus_t *lines, int n_lines,
+                              struct line_pnts *BPoints)
+{
+    int i, line, aline, dir;
+    static struct line_pnts *Points;
+    
+    if (!Points)
+        Points = Vect_new_line_struct();
+    
+    Vect_reset_line(BPoints);
+    for (i = 0; i < n_lines; i++) {
+        line = lines[i];
+        aline = abs(line);
+        G_debug(5, "  append line(%d) = %d", i, line);
+        
+        if (0 > Vect_read_line(Map, Points, NULL, aline))
+            return -1;
+        
+        dir = line > 0 ? GV_FORWARD : GV_BACKWARD;
+        Vect_append_points(BPoints, Points, dir);
+        BPoints->n_points--;    /* skip last point, avoids duplicates */
+    }
+    BPoints->n_points++;        /* close polygon */
+
+    return BPoints->n_points;
+}

lib/vector/Vlib/area_pg.c

@@ -0,0 +1,114 @@
+/*!
+   \file lib/vector/Vlib/area_pg.c
+
+   \brief Vector library - area-related functions (PostGIS Topology)
+
+   Higher level functions for reading/writing/manipulating vectors.
+
+   (C) 2013 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 Martin Landa <landa.martin gmail.com>
+*/
+
+#include <grass/vector.h>
+
+#ifdef HAVE_POSTGRES
+#include "pg_local_proto.h"
+
+static PGresult *build_stmt(const struct Plus_head *, const struct Format_info_pg *, const plus_t *, int);
+
+/*!
+  \brief Get area boundary points (PostGIS Topology)
+  
+  Used by Vect_build_line_area() and Vect_get_area_points().
+
+  \param Map pointer to Map_info struct
+  \param lines array of boundary lines
+  \param n_lines number of lines in array
+  \param[out] APoints pointer to output line_pnts struct
+
+  \return number of points
+  \return -1 on error
+*/
+int Vect__get_area_points_pg(const struct Map_info *Map, const plus_t *lines, int n_lines,
+                             struct line_pnts *APoints)
+{
+    int i, direction;
+    
+    struct Format_info_pg *pg_info;
+    
+    PGresult *res;
+    
+    pg_info = (struct Format_info_pg *)&(Map->fInfo.pg);
+   
+    Vect_reset_line(APoints);
+
+    res = build_stmt(&(Map->plus), pg_info, lines, n_lines);
+    if (!res)
+        return -1;
+    
+    for (i = 0; i < n_lines; i++) {
+        Vect__cache_feature_pg(PQgetvalue(res, i, 0), FALSE, FALSE,
+                               &(pg_info->cache), NULL); /* do caching in readable way */
+	direction = lines[i] > 0 ? GV_FORWARD : GV_BACKWARD;
+	Vect_append_points(APoints, pg_info->cache.lines[0], direction);
+	APoints->n_points--;	/* skip last point, avoids duplicates */
+    }
+    APoints->n_points++;	/* close polygon */
+
+    PQclear(res);
+    
+    return APoints->n_points;
+}
+
+PGresult *build_stmt(const struct Plus_head *plus, const struct Format_info_pg *pg_info,
+                     const plus_t *lines, int n_lines)
+{
+    int i, line;
+    size_t stmt_id_size;
+    char *stmt, *stmt_id, buf_id[128];
+
+    struct P_line *BLine;
+    
+    PGresult *res;
+    
+    stmt = NULL;
+    stmt_id_size = DB_SQL_MAX;
+    stmt_id = (char *) G_malloc(stmt_id_size);
+    stmt_id[0] = '\0';
+
+    for (i = 0; i < n_lines; i++) {
+        if (strlen(stmt_id) + 100 > stmt_id_size) {
+            stmt_id_size = strlen(stmt_id) + DB_SQL_MAX;
+            stmt_id = (char *) G_realloc(stmt_id, stmt_id_size);
+        }
+	line = abs(lines[i]);
+        BLine = plus->Line[line];
+        if (i > 0)
+            strcat(stmt_id, ",");
+        sprintf(buf_id, "%d", (int) BLine->offset);
+        strcat(stmt_id, buf_id);
+    }
+    G_asprintf(&stmt, "SELECT geom FROM \"%s\".edge_data WHERE edge_id IN (%s) "
+               "ORDER BY POSITION(edge_id::text in '%s')", pg_info->toposchema_name,
+               stmt_id, stmt_id);
+    G_free(stmt_id);
+    
+    G_debug(2, "SQL: %s", stmt);
+    res = PQexec(pg_info->conn, stmt);
+    G_free(stmt);
+
+    if (!res || PQresultStatus(res) != PGRES_TUPLES_OK ||
+        PQntuples(res) != n_lines) {
+        if (res)
+            PQclear(res);
+        
+        return NULL;
+    }
+
+    return res;
+}
+#endif

lib/vector/Vlib/build.c

@@ -128,7 +128,6 @@ int Vect_build_line_area(struct Map_info *Map, int iline, int side)
 	/* add area structure to plus */
 	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)"));
 	}
 	G_debug(3, "  -> area %d", area);
@@ -137,7 +136,6 @@ int Vect_build_line_area(struct Map_info *Map, int iline, int side)
     else if (area_size < 0) {	/* CCW: island */
 	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)"));
 	}
 	G_debug(3, "  -> isle %d", isle);

lib/vector/Vlib/build_nat.c

@@ -61,8 +61,8 @@ int Vect_build_nat(struct Map_info *Map, int build)
     Cats = Vect_new_cats_struct();
     
     if (plus->built < GV_BUILD_BASE) {
-      int npoints, c;
-
+        int npoints, c;
+        
 	/* 
 	 *  We shall go through all primitives in coor file and add
 	 *  new node for each end point to nodes structure if the node
@@ -233,47 +233,3 @@ int Vect_build_nat(struct Map_info *Map, int build)
     
     return 1;
 }
-
-/*!
-  \brief Get area boundary points (native format)
-  
-  Used by Vect_build_line_area().
-  
-  \param Map pointer to Map_info struct
-  \param lines array of boundary lines
-  \param n_lines number of lines in array
-  \param[out] APoints pointer to output line_pnts struct
-
-  \return number of points
-  \return -1 on error
-*/
-int Vect__get_area_points_nat(struct Map_info *Map, plus_t *lines, int n_lines,
-                              struct line_pnts *APoints)
-{
-    int j, line, direction;
-    struct Plus_head *plus;
-    struct P_line *BLine;
-
-    plus = &(Map->plus);
-    
-    if (!Points)
-        Points = Vect_new_line_struct();
-    
-    Vect_reset_line(APoints);
-    for (j = 0; j < n_lines; j++) {
-	line = abs(lines[j]);
-	BLine = plus->Line[line];
-	G_debug(5, "  line[%d] = %d, offset = %lu", j, line,
-		(unsigned long) BLine->offset);
-	if (0 > Vect_read_line(Map, Points, NULL, line))
-            return -1;
-        
-        direction = lines[j] > 0 ? GV_FORWARD : GV_BACKWARD;
-	Vect_append_points(APoints, Points, direction);
-	APoints->n_points--;	/* skip last point, avoids duplicates */
-    }
-    APoints->n_points++;	/* close polygon */
-
-    return APoints->n_points;
-}
-

lib/vector/Vlib/build_pg.c

@@ -476,82 +476,4 @@ int has_topo_grass(const struct Format_info_pg *pg_info)
 
     return has_topo;
 }
-
-/*!
-  \brief Get area boundary points (PostGIS Topology)
-  
-  Used by Vect_build_line_area().
-
-  \param Map pointer to Map_info struct
-  \param lines array of boundary lines
-  \param n_lines number of lines in array
-  \param[out] APoints pointer to output line_pnts struct
-
-  \return number of points
-  \return -1 on error
-*/
-int Vect__get_area_points_pg(const struct Map_info *Map, const plus_t *lines, int n_lines,
-                             struct line_pnts *APoints)
-{
-
-    int i, line, direction;
-    size_t stmt_id_size;
-    char *stmt, *stmt_id, buf_id[128];
-    
-    const struct Plus_head *plus;
-    struct Format_info_pg *pg_info;
-    struct P_line *BLine;
-    
-    PGresult *res;
-    
-    plus = &(Map->plus);
-    pg_info = (struct Format_info_pg *)&(Map->fInfo.pg);
-    
-    stmt = NULL;
-    stmt_id_size = DB_SQL_MAX;
-    stmt_id = (char *) G_malloc(stmt_id_size);
-    stmt_id[0] = '\0';
-    
-    Vect_reset_line(APoints);
-
-    for (i = 0; i < n_lines; i++) {
-        if (strlen(stmt_id) + 100 > stmt_id_size) {
-            stmt_id_size = strlen(stmt_id) + DB_SQL_MAX;
-            stmt_id = (char *) G_realloc(stmt_id, stmt_id_size);
-        }
-	line = abs(lines[i]);
-        BLine = plus->Line[line];
-        if (i > 0)
-            strcat(stmt_id, ",");
-        sprintf(buf_id, "%d", (int) BLine->offset);
-        strcat(stmt_id, buf_id);
-    }
-    G_asprintf(&stmt, "SELECT geom FROM \"%s\".edge_data WHERE edge_id IN (%s) "
-               "ORDER BY POSITION(edge_id::text in '%s')", pg_info->toposchema_name,
-               stmt_id, stmt_id);
-    G_free(stmt_id);
-    
-    G_debug(2, "SQL: %s", stmt);
-    res = PQexec(pg_info->conn, stmt);
-    G_free(stmt);
-
-    if (!res || PQresultStatus(res) != PGRES_TUPLES_OK ||
-        PQntuples(res) != n_lines) {
-        if (res)
-            PQclear(res);
-        
-        return -1;
-    }
-    
-    for (i = 0; i < n_lines; i++) {
-        Vect__cache_feature_pg(PQgetvalue(res, i, 0), FALSE, FALSE,
-                               &(pg_info->cache), NULL); /* do caching in readable way */
-	direction = lines[i] > 0 ? GV_FORWARD : GV_BACKWARD;
-	Vect_append_points(APoints, pg_info->cache.lines[0], direction);
-	APoints->n_points--;	/* skip last point, avoids duplicates */
-    }
-    APoints->n_points++;	/* close polygon */
-
-    return APoints->n_points;
-}
 #endif

lib/vector/Vlib/local_proto.h

@@ -8,8 +8,6 @@
 
 /* area.c */
 int Vect__get_area_points(const struct Map_info *, const plus_t *, int, struct line_pnts *);
-
-/* build_nat.c */
 int Vect__get_area_points_nat(const struct Map_info *, const plus_t *, int, struct line_pnts *);
 
 /* map.c */

lib/vector/Vlib/pg_local_proto.h

@@ -53,7 +53,7 @@ struct feat_parts
     int            *idx;     /* index in cache where to start */
 };
 
-/* build_pg.c */
+/* area_pg.c */
 int Vect__get_area_points_pg(const struct Map_info *, const plus_t *, int, struct line_pnts *);
 
 /* read_pg.c */