vlib(pg): store GRASS topology in DB (can be optionally disabled by the user)

          various fixes reading and writing topological elements from DB
          use offset to map nodes
          doxygen cosmetics


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

include/vect/dig_defines.h

@@ -213,9 +213,11 @@ typedef enum overlay_operator OVERLAY_OPERATOR;
 /*! \brief Maximum category value */
 #define GV_CAT_MAX   PORT_INT_MAX
 
-/*! \brief GRASS ASCII vector format */
+/*! \brief GRASS ASCII vector format - point format */
 #define GV_ASCII_FORMAT_POINT 0
+/*! \brief GRASS ASCII vector format - standard format */
 #define GV_ASCII_FORMAT_STD   1
+/*! \brief GRASS ASCII vector format - well-known-text format */
 #define GV_ASCII_FORMAT_WKT   2
 
 /*! \brief Simple feature types
@@ -250,6 +252,6 @@ typedef enum
 #define HEADSTR	50
 
 /*! \brief GRASS-PostGIS data provider - default fid column */
-#define GV_PG_FID_COLUMN       "fid"
+#define GV_PG_FID_COLUMN      "fid"
 /*! \brief GRASS-PostGIS data provider - default geometry column */
 #define GV_PG_GEOMETRY_COLUMN "geom"

include/vect/dig_externs.h

@@ -223,7 +223,7 @@ int rtree_search(struct RTree *, struct RTree_Rect *,
                  SearchHitCallback , void *, struct Plus_head *);
 
 /* struct_alloc.c */
-int dig_node_alloc_line(struct P_node *, int add);
+int dig_node_alloc_line(struct P_node *, int);
 int dig_alloc_nodes(struct Plus_head *, int);
 int dig_alloc_lines(struct Plus_head *, int);
 int dig_alloc_areas(struct Plus_head *, int);

include/vect/dig_structs.h

@@ -433,7 +433,7 @@ struct Format_info_offset
       5	   2	  third part (polygon)
       6	   0	  first ring of polygon
 
-      Topology:
+      GRASS topology:
       
       line idx
       -----------------
@@ -441,6 +441,9 @@ struct Format_info_offset
       2    2      line
       3    4      boundary
       4    1      centroid read from topo (idx == FID)
+
+      In PostGIS Topology mode the offset array is used for mapping
+      nodes.
       \endverbatim
     */
     int *array;
@@ -674,6 +677,13 @@ struct Format_info_pg
     */
     char    *toposchema_name;
     int      toposchema_id;
+    /*!
+      \brief Topology format
+
+      TRUE to store only Topo-Geo data in DB otherwise GRASS-like
+      topology is also maintained in DB
+    */
+    int      topo_geo_only;
 };
 
 /*!

lib/vector/Vlib/build_pg.c

@@ -28,6 +28,10 @@ static int build_topogeom_stmt(const struct Format_info_pg *, int, int, int, cha
 static int save_map_bbox(const struct Format_info_pg *, const struct bound_box*);
 static int create_topo_grass(const struct Format_info_pg *);
 static int has_topo_grass(const struct Format_info_pg *);
+static int write_nodes(const struct Plus_head *, const struct Format_info_pg *);
+static int write_areas(const struct Plus_head *, const struct Format_info_pg *);
+static int write_isles(const struct Plus_head *, const struct Format_info_pg *);
+static void build_stmt_id(const void *, int, int, const struct Plus_head *, char **, size_t *);
 #endif
 
 /*!
@@ -71,8 +75,8 @@ int Vect_build_pg(struct Map_info *Map, int build)
 
     /* TODO move this init to better place (Vect_open_ ?), because in
        theory build may be reused on level2 */
-    if (build >= plus->built && build > GV_BUILD_BASE) {
-        G_free((void *)pg_info->offset.array);
+    if (!pg_info->toposchema_name && build >= plus->built && build > GV_BUILD_BASE) {
+        G_free(pg_info->offset.array);
         G_zero(&(pg_info->offset), sizeof(struct Format_info_offset));
     }
 
@@ -162,6 +166,11 @@ int build_topo(struct Map_info *Map, int build)
     if (Vect__execute_pg(pg_info->conn, "BEGIN"))
         return 0;
     
+    /* write full node topo info to DB if requested */
+    if (!pg_info->topo_geo_only) {
+        write_nodes(plus, pg_info);
+    }
+
     /* update faces from GRASS Topology */
     if (build >= GV_BUILD_AREAS) {
         /* do clean up (1-3)
@@ -186,7 +195,7 @@ int build_topo(struct Map_info *Map, int build)
             return 0;
         }
 
-        /* 3) delete faces */        
+        /* 3) delete faces (areas/isles) */        
         sprintf(stmt, "DELETE FROM \"%s\".face WHERE "
                 "face_id != 0", pg_info->toposchema_name);
         G_debug(2, "SQL: %s", stmt);
@@ -194,7 +203,22 @@ int build_topo(struct Map_info *Map, int build)
             Vect__execute_pg(pg_info->conn, "ROLLBACK");
             return 0;
         }
-        
+        if (!pg_info->topo_geo_only) {
+            sprintf(stmt, "DELETE FROM \"%s\".%s", pg_info->toposchema_name, TOPO_TABLE_AREA);
+            G_debug(2, "SQL: %s", stmt);
+            if(Vect__execute_pg(pg_info->conn, stmt) == -1) {
+                Vect__execute_pg(pg_info->conn, "ROLLBACK");
+                return 0;
+            }
+            
+            sprintf(stmt, "DELETE FROM \"%s\".%s", pg_info->toposchema_name, TOPO_TABLE_ISLE);
+            G_debug(2, "SQL: %s", stmt);
+            if(Vect__execute_pg(pg_info->conn, stmt) == -1) {
+                Vect__execute_pg(pg_info->conn, "ROLLBACK");
+                return 0;
+            }
+        }
+
         /* 4) insert faces & update nodes (containing_face) based on
          * GRASS topology */
         G_message(_("Updating faces..."));
@@ -266,6 +290,11 @@ int build_topo(struct Map_info *Map, int build)
                 return 0;
             }
         }
+
+        /* write full area topo info to DB if requested */
+        if (!pg_info->topo_geo_only) {
+            write_areas(plus, pg_info);
+        }
     } /* build >= GV_BUILD_AREAS */
     
     if (build >= GV_BUILD_ATTACH_ISLES) {
@@ -275,6 +304,11 @@ int build_topo(struct Map_info *Map, int build)
             Isle = plus->Isle[isle];
             Vect__insert_face_pg(Map, -isle);
         }
+
+        /* write full isles topo info to DB if requested */
+        if (!pg_info->topo_geo_only) {
+            write_isles(plus, pg_info);
+        }
     } /* build >= GV_BUILD_ISLES */
 
     if (pg_info->feature_type == SF_POLYGON) {
@@ -479,4 +513,237 @@ int has_topo_grass(const struct Format_info_pg *pg_info)
 
     return has_topo;
 }
+
+/*!
+  \brief Insert node into 'node_grass' table
+  
+  Writes (see P_node struct):
+   - lines
+   - angles
+
+  Already stored in Topo-Geo:
+   - x,y,z (geom)
+   
+  \param plus pointer to Plus_head struct
+  \param pg_info pointer to Format_info_pg struct
+
+  \return 0 on success
+  \return -1 on error
+*/
+int write_nodes(const struct Plus_head *plus,
+                const struct Format_info_pg *pg_info)
+{
+    int i, node_id;
+    size_t stmt_lines_size, stmt_angles_size;
+    char *stmt_lines, *stmt_angles;
+    char stmt[DB_SQL_MAX];
+        
+    const struct P_node *Node;
+    const struct Format_info_offset *offset;
+
+    offset = &(pg_info->offset);
+    if (plus->n_nodes != offset->array_num)
+        return -1;
+    
+    stmt_lines = stmt_angles = NULL;
+    for (i = 1; i <= plus->n_nodes; i++) {
+        Node = plus->Node[i];
+        if (!Node)
+            continue; /* should not happen */
+        
+        node_id = offset->array[i-1]; 
+
+        /* 'lines' array */
+        build_stmt_id(Node->lines, Node->n_lines, TRUE, plus, &stmt_lines, &stmt_lines_size);
+        /* 'angle' array */
+        build_stmt_id(Node->angles, Node->n_lines, FALSE, NULL, &stmt_angles, &stmt_angles_size);
+        
+        /* build SQL statement to add new node into 'node_grass' */
+        sprintf(stmt, "INSERT INTO \"%s\".%s VALUES ("
+                "%d, '{%s}', '{%s}')", pg_info->toposchema_name, TOPO_TABLE_NODE,
+                node_id, stmt_lines, stmt_angles);
+        if (Vect__execute_pg(pg_info->conn, stmt) == -1) {
+            return -1;
+        }
+    }
+
+    G_free(stmt_lines);
+    G_free(stmt_angles);
+    
+    return 0;
+}
+
+/*!
+  \brief Insert area into 'area_grass' table
+
+  Writes (see P_area struct):
+   - lines
+   - centroid
+   - isles
+
+  \param plus pointer to Plus_head struct
+  \param pg_info pointer to Format_info_pg struct
+
+  \return 0 on success
+  \return -1 on error
+*/
+int write_areas(const struct Plus_head *plus,
+                const struct Format_info_pg *pg_info)
+{
+    int area, centroid;
+    size_t stmt_lines_size, stmt_isles_size;
+    char *stmt_lines, *stmt_isles;
+    char stmt[DB_SQL_MAX];
+        
+    const struct P_line *Line;
+    const struct P_area *Area;
+
+    stmt_lines = stmt_isles = NULL;
+    for (area = 1; area <= plus->n_areas; area++) {
+        Area = plus->Area[area];
+        if (!Area)
+            continue; /* should not happen */
+        
+        /* 'lines' array */
+        build_stmt_id(Area->lines, Area->n_lines, TRUE, NULL, &stmt_lines, &stmt_lines_size);
+        /* 'isles' array */
+        build_stmt_id(Area->isles, Area->n_isles, TRUE, NULL, &stmt_isles, &stmt_isles_size);
+        
+        Line = plus->Line[Area->centroid];
+        if (!Line)
+            return -1;
+        centroid = (int) Line->offset;
+        
+        /* build SQL statement to add new node into 'node_grass' */
+        sprintf(stmt, "INSERT INTO \"%s\".%s VALUES ("
+                "%d, '{%s}', %d, '{%s}')", pg_info->toposchema_name, TOPO_TABLE_AREA,
+                area, stmt_lines, centroid, stmt_isles);
+        if (Vect__execute_pg(pg_info->conn, stmt) == -1) {
+            return -1;
+        }
+    }
+
+    G_free(stmt_lines);
+    G_free(stmt_isles);
+    
+    return 0;
+}
+
+/*!
+  \brief Insert isle into 'isle_grass' table
+
+  Writes (see P_isle struct):
+   - lines
+   - area
+
+  \param plus pointer to Plus_head struct
+  \param pg_info pointer to Format_info_pg struct
+
+  \return 0 on success
+  \return -1 on error
+*/
+int write_isles(const struct Plus_head *plus,
+                const struct Format_info_pg *pg_info)
+{
+    int isle;
+    size_t stmt_lines_size;
+    char *stmt_lines;
+    char stmt[DB_SQL_MAX];
+        
+    const struct P_isle *Isle;
+
+    stmt_lines = NULL;
+    for (isle = 1; isle <= plus->n_isles; isle++) {
+        Isle = plus->Isle[isle];
+        if (!Isle)
+            continue; /* should not happen */
+        
+        /* 'lines' array */
+        build_stmt_id(Isle->lines, Isle->n_lines, TRUE, NULL, &stmt_lines, &stmt_lines_size);
+        
+        /* build SQL statement to add new node into 'node_grass' */
+        sprintf(stmt, "INSERT INTO \"%s\".%s VALUES ("
+                "%d, '{%s}', %d)", pg_info->toposchema_name, TOPO_TABLE_ISLE,
+                isle, stmt_lines, Isle->area);
+        if (Vect__execute_pg(pg_info->conn, stmt) == -1) {
+            return -1;
+        }
+    }
+
+    G_free(stmt_lines);
+    
+    return 0;
+}
+
+/*!
+  \brief Create PG-like array for int/float array
+
+  \param array array of items
+  \param nitems number of items in the array
+  \param is_int TRUE for array of integers otherwise floats
+  \param plus pointer to Plus_head struct
+  \param[in,out] output buffer (re-used)
+  \param[in,out] buffer size
+*/
+void build_stmt_id(const void *array, int nitems, int is_int, const struct Plus_head *plus,
+                   char **stmt, size_t *stmt_size)
+{
+    int i, ivalue;
+    int *iarray;
+    float *farray;
+    
+    size_t stmt_id_size;
+    char *stmt_id, buf_id[128];
+
+    struct P_line *Line;
+    
+    if (is_int)
+        iarray = (int *) array;
+    else
+        farray = (float *) array;
+    
+    if (!(*stmt)) {
+        stmt_id_size = DB_SQL_MAX;
+        stmt_id = (char *) G_malloc(stmt_id_size);
+    }
+    else {
+        stmt_id_size = *stmt_size;
+        stmt_id = *stmt;
+    }
+
+    /* reset array */
+    stmt_id[0] = '\0';
+    
+    for (i = 0; i < nitems; i++) {
+        /* realloc array if needed */
+        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);
+        }
+        
+        if (is_int) {
+            if (plus) {
+                Line = plus->Line[abs(iarray[i])];
+                ivalue = (int) Line->offset;
+                if (iarray[i] < 0)
+                    ivalue *= -1;
+            }
+            else {
+                ivalue = iarray[i];
+            }
+            sprintf(buf_id, "%d", ivalue);
+        }
+        else {
+            sprintf(buf_id, "%f", farray[i]);
+        }
+
+        if (i > 0)
+            strcat(stmt_id, ",");
+        strcat(stmt_id, buf_id);
+    }
+
+    *stmt = stmt_id;
+    *stmt_size = stmt_id_size;
+}
+
 #endif

lib/vector/Vlib/copy.c

@@ -37,7 +37,11 @@
 
 static int copy_lines_1(struct Map_info *, int, struct Map_info *);
 static int copy_lines_2(struct Map_info *, int, int, struct Map_info *);
+#if 0
 static int copy_nodes(const struct Map_info *, struct Map_info *);
+#endif
+static int copy_line_nodes(const struct Map_info *, int, int, struct line_pnts *,
+                           struct Map_info *);
 static int is_isle(const struct Map_info *, int);
 static int copy_areas(const struct Map_info *, int, struct Map_info *);
 
@@ -109,11 +113,12 @@ int Vect_copy_map_lines_field(struct Map_info *In, int field,
     ret = 0;
     if (Vect_level(In) >= 2) {
         /* -> copy features on level 2 */
+#if 0
         if (topo == TOPO_POSTGIS) {
             /* PostGIS topology - copy also nodes */
             copy_nodes(In, Out);
         }
-        
+#endif
         /* copy features */
         ret += copy_lines_2(In, field, topo, Out);
 
@@ -199,17 +204,20 @@ int copy_lines_1(struct Map_info *In, int field, struct Map_info *Out)
 */
 int copy_lines_2(struct Map_info *In, int field, int topo, struct Map_info *Out)
 {
-  int i, type, nlines, nskipped;
-    int ret, left, rite, centroid;
+    int i, type, nlines, nskipped;
+    int ret, left, rite, centroid, with_z;
 
-    struct line_pnts *Points, *CPoints;
+    struct line_pnts *Points, *CPoints, *NPoints;
     struct line_cats *Cats, *CCats;
 
     Points  = Vect_new_line_struct();
     CPoints = Vect_new_line_struct();
+    NPoints = Vect_new_line_struct();
     Cats    = Vect_new_cats_struct();
     CCats   = Vect_new_cats_struct();
     
+    with_z = Vect_is_3d(In);
+    
     ret = 0;
     nlines = Vect_get_num_lines(In);
     if (topo == TOPO_NONE) {
@@ -225,6 +233,7 @@ int copy_lines_2(struct Map_info *In, int field, int topo, struct Map_info *Out)
     else
         G_message(_("Copying features..."));    
     
+    Vect_append_point(NPoints, 0., 0., 0.);
     nskipped = 0;
     for (i = 1; i <= nlines; i++) {
         if (!Vect_line_alive(In, i))
@@ -288,6 +297,39 @@ int copy_lines_2(struct Map_info *In, int field, int topo, struct Map_info *Out)
             }
         }
         
+        /* copy also nodes connected to the line (PostGIS Topology
+         * mode only) */
+        if (topo == TOPO_POSTGIS && (type & GV_LINES)) {
+            int n1, n2;
+            
+            struct P_line *Line;
+            struct Format_info_offset *offset;
+            
+            offset = &(Out->fInfo.pg.offset);
+            
+            n1 = n2 = -1;
+            Line = In->plus.Line[i];
+            if (Line) {
+                if (type == GV_LINE) {
+                    struct P_topo_l *topo = (struct P_topo_l *)Line->topo;
+                    
+                    n1 = topo->N1;
+                    n2 = topo->N2;
+                }
+                else if (type == GV_BOUNDARY) {
+                    struct P_topo_b *topo = (struct P_topo_b *)Line->topo;
+                    
+                    n1 = topo->N1;
+                    n2 = topo->N2;
+                }
+            }
+                
+            if (n1 > 0 && (n1 > offset->array_num || offset->array[n1-1] == 0)) 
+                copy_line_nodes(In, n1, with_z, NPoints, Out);
+            if (n2 > 0 && (n2 > offset->array_num || offset->array[n2-1] == 0))
+                copy_line_nodes(In, n2, with_z, NPoints, Out);
+        }
+
         if (-1 == Vect_write_line(Out, type, Points, Cats)) {
             G_warning(_("Writing new feature failed"));
             return 1;
@@ -299,12 +341,14 @@ int copy_lines_2(struct Map_info *In, int field, int topo, struct Map_info *Out)
     
     Vect_destroy_line_struct(Points);
     Vect_destroy_line_struct(CPoints);
+    Vect_destroy_line_struct(NPoints);
     Vect_destroy_cats_struct(Cats);
     Vect_destroy_cats_struct(CCats);
 
     return ret;
 }
 
+#if 0
 /*!
   \brief Copy nodes as points (PostGIS Topology only)
 
@@ -317,40 +361,50 @@ int copy_lines_2(struct Map_info *In, int field, int topo, struct Map_info *Out)
 int copy_nodes(const struct Map_info *In, struct Map_info *Out)
 {
     int nnodes, node, with_z;
-    double x, y, z;
+
     struct line_pnts *Points;
-    
-    Points = Vect_new_line_struct();
+
+    Points  = Vect_new_line_struct();
     
     with_z = Vect_is_3d(In);
     
     nnodes = Vect_get_num_nodes(In);
-    if (nnodes > 0)
-        G_message(_("Exporting nodes..."));
+    G_message(_("Exporting nodes..."));
     Vect_append_point(Points, 0., 0., 0.);
     for (node = 1; node <= nnodes; node++) {
         G_debug(3, "Exporting GRASS node %d", node);
         
         G_percent(node, nnodes, 5);
-        Vect_get_node_coor(In, node, &x, &y, &z);
-        Points->x[0] = x;
-        Points->y[0] = y;
-        if (with_z)
-            Points->z[0] = z;
-        
+        copy_line_nodes(In, node, with_z, Points, Out);
+    }
+    
+    Vect_destroy_line_struct(Points);
+
+    return 0;
+}
+#endif
+
+int copy_line_nodes(const struct Map_info *In, int node, int with_z,
+                    struct line_pnts *Points, struct Map_info *Out)
+{
+    double x, y, z;
+ 
+    Vect_get_node_coor(In, node, &x, &y, &z);
+    Points->x[0] = x;
+    Points->y[0] = y;
+    if (with_z)
+        Points->z[0] = z;
+    
 #ifdef HAVE_POSTGRES
-        if (-1 == V2__write_node_pg(Out, Points)) {
-            G_warning(_("Writing node %d failed"), node);
-            return 1;
-        }
-#else
-        G_fatal_error(_("GRASS is not compiled with PostgreSQL support"));
+    if (-1 == V2__write_node_pg(Out, Points)) {
+        G_warning(_("Writing node %d failed"), node);
         return 1;
-#endif
     }
+#else
+    G_fatal_error(_("GRASS is not compiled with PostgreSQL support"));
+    return 1;
+#endif
 
-    Vect_destroy_line_struct(Points);
-    
     return 0;
 }
 

lib/vector/Vlib/open.c

@@ -1316,12 +1316,18 @@ int map_format(struct Map_info *Map)
             /* PostGIS topology enabled ? */
             p = G_find_key_value("topology", key_val);
             if (p && G_strcasecmp(p, "yes") == 0) {
-                /* define topology name
-                   this should be configurable by the user
-                */
-                G_asprintf(&(pg_info->toposchema_name), "topo_%s",
-                           pg_info->table_name);
+                /* define topology name */
+                p = G_find_key_value("toposchema_name", key_val);
+                if (p)
+                    pg_info->toposchema_name = G_store(p);
+                else
+                    G_asprintf(&(pg_info->toposchema_name), "topo_%s",
+                               pg_info->table_name);
+
+                G_debug(1, "PG: topology = yes, schema_name = %s",
+                        pg_info->toposchema_name);
             }
+            G_debug(1, "PG: topology = no");
             
             if (getenv("GRASS_VECTOR_EXTERNAL_IMMEDIATE")) {
                 /* vector features are written directly to PostGIS layer */

lib/vector/Vlib/open_pg.c

@@ -39,12 +39,20 @@ static void connect_db(struct Format_info_pg *);
 static int check_topo(struct Format_info_pg *, struct Plus_head *);
 static int parse_bbox(const char *, struct bound_box *);
 static struct P_node *read_p_node(struct Plus_head *, int, int,
-                                  const char *, struct Format_info_pg *);
+                                  const char *, const char *, const char *,
+                                  struct Format_info_pg *);
 static struct P_line *read_p_line(struct Plus_head *, int,
                                   const struct edge_data *,
                                   struct Format_info_cache *);
+static struct P_area *read_p_area(struct Plus_head *, int,
+                                  const char *, int, const char *);
+static struct P_isle *read_p_isle(struct Plus_head *, int,
+                                  const char *, int);
 static int load_plus_head(struct Format_info_pg *, struct Plus_head *);
 static void notice_processor(void *, const char *);
+static char **scan_array(const char *);
+static int remap_node(const struct Format_info_offset *, int);
+static int remap_line(const struct Plus_head*, int);
 #endif
 
 /*!
@@ -543,8 +551,14 @@ int check_topo(struct Format_info_pg *pg_info, struct Plus_head *plus)
     pg_info->toposchema_name = G_store(PQgetvalue(res, 0, 1));
     pg_info->topogeom_column = G_store(PQgetvalue(res, 0, 3));
 
-    G_debug(1, "PostGIS topology detected: schema = %s column = %s",
-            pg_info->toposchema_name, pg_info->topogeom_column);
+    /* check extra GRASS tables */
+    sprintf(stmt, "SELECT COUNT(*) FROM pg_tables WHERE schemaname = '%s' "
+            "AND tablename LIKE '%%_grass'", pg_info->toposchema_name);
+    if (Vect__execute_get_value_pg(pg_info->conn, stmt) != TOPO_TABLE_NUM)
+        pg_info->topo_geo_only = TRUE;
+    
+    G_debug(1, "PostGIS topology detected: schema = %s column = %s topo_geo_only = %d",
+            pg_info->toposchema_name, pg_info->topogeom_column, pg_info->topo_geo_only);
     
     /* check for 3D */
     if (strcmp(PQgetvalue(res, 0, 2), "t") == 0)
@@ -630,45 +644,60 @@ int parse_bbox(const char *value, struct bound_box *bbox)
   \param n index (starts at 1)
   \param id node id (table "node")
   \param wkb_data geometry data (wkb)
+  \param lines_data lines array or NULL
+  \param angles_data angles array or NULL
   \param pg_info pointer to Format_info_pg sttucture
 
   \return pointer to new P_node struct
   \return NULL on error
 */
 struct P_node *read_p_node(struct Plus_head *plus, int n,
-                           int id, const char *wkb_data,
-                           struct Format_info_pg *pg_info)
+                           int id, const char *wkb_data, const char *lines_data,
+                           const char *angles_data, struct Format_info_pg *pg_info)
 {
     int i, cnt;
-    char stmt[DB_SQL_MAX];
-    
+    char **lines, **angles;
+
     struct P_node *node;
     struct line_pnts *points;
     
     PGresult *res;
     
     /* get lines connected to the node */
-    sprintf(stmt,
-            "SELECT edge_id,'s' as node,"
-            "ST_Azimuth(ST_StartPoint(geom), ST_PointN(geom, 2)) AS angle"
-            " FROM \"%s\".edge WHERE start_node = %d UNION ALL "
-            "SELECT edge_id,'e' as node,"
-            "ST_Azimuth(ST_EndPoint(geom), ST_PointN(geom, ST_NumPoints(geom) - 1)) AS angle"
-            " FROM \"%s\".edge WHERE end_node = %d"
-            " ORDER BY angle DESC",
-            pg_info->toposchema_name, id,
-            pg_info->toposchema_name, id);
-    G_debug(2, "SQL: %s", stmt);
-    res = PQexec(pg_info->conn, stmt);
-    if (!res || PQresultStatus(res) != PGRES_TUPLES_OK) {
-        G_warning(_("Inconsistency in topology: unable to read node %d"), id);
-        if (res)
-            PQclear(res);
-        return NULL;
+    lines = angles = NULL;
+    if (!lines_data && !angles_data) { /* pg_info->topo_geo_only == TRUE */
+        char stmt[DB_SQL_MAX];
+        
+        sprintf(stmt,
+                "SELECT edge_id,'s' as node,"
+                "ST_Azimuth(ST_StartPoint(geom), ST_PointN(geom, 2)) AS angle"
+                " FROM \"%s\".edge WHERE start_node = %d UNION ALL "
+                "SELECT edge_id,'e' as node,"
+                "ST_Azimuth(ST_EndPoint(geom), ST_PointN(geom, ST_NumPoints(geom) - 1)) AS angle"
+                " FROM \"%s\".edge WHERE end_node = %d"
+                " ORDER BY angle DESC",
+                pg_info->toposchema_name, id,
+                pg_info->toposchema_name, id);
+        G_debug(2, "SQL: %s", stmt);
+        res = PQexec(pg_info->conn, stmt);
+        if (!res || PQresultStatus(res) != PGRES_TUPLES_OK) {
+            G_warning(_("Inconsistency in topology: unable to read node %d"), id);
+            if (res)
+                PQclear(res);
+            return NULL;
+        }
+        cnt = PQntuples(res);
     }
-    cnt = PQntuples(res);
-    
-    if (cnt == 0) { /* dead ??? */
+    else { /* pg_info->topo_geo_only != TRUE */
+        lines  = scan_array(lines_data);
+        angles = scan_array(angles_data);
+        
+        cnt = G_number_of_tokens(lines);
+        if (cnt != G_number_of_tokens(angles))
+            return NULL; /* 'lines' and 'angles' array must have the same size */
+    }
+
+    if (cnt == 0) { /* dead */
         plus->Node[n] = NULL;
         return NULL;
     }
@@ -681,22 +710,36 @@ struct P_node *read_p_node(struct Plus_head *plus, int n,
         return NULL;
 
     /* lines / angles */
-    for (i = 0; i < node->n_lines; i++) {
-        node->lines[i] = atoi(PQgetvalue(res, i, 0));
-        if (strcmp(PQgetvalue(res, i, 1), "s") != 0) {
-            /* end node */
-            node->lines[i] *= -1;
+    if (lines) {
+        for (i = 0; i < node->n_lines; i++) {
+            node->lines[i] = atoi(lines[i]);
+            node->angles[i] = atof(angles[i]);
+
+            G_debug(5, "\tline = %d angle = %f", node->lines[i],
+                    node->angles[i]);
         }
-        node->angles[i] = M_PI / 2 - atof(PQgetvalue(res, i, 2));
-        /* angles range <-PI; PI> */
-        if (node->angles[i] > M_PI)
-            node->angles[i] = node->angles[i] - 2 * M_PI;
-        if (node->angles[i] < -1.0 * M_PI)
-            node->angles[i] = node->angles[i] + 2 * M_PI;
-        G_debug(5, "\tline = %d angle = %f", node->lines[i],
-                node->angles[i]);
+        
+        G_free_tokens(lines);
+        G_free_tokens(angles);
+    }
+    else {
+        for (i = 0; i < node->n_lines; i++) {
+            node->lines[i] = atoi(PQgetvalue(res, i, 0));
+            if (strcmp(PQgetvalue(res, i, 1), "s") != 0) {
+                /* end node */
+                node->lines[i] *= -1;
+            }
+            node->angles[i] = M_PI / 2 - atof(PQgetvalue(res, i, 2));
+            /* angles range <-PI; PI> */
+            if (node->angles[i] > M_PI)
+                node->angles[i] = node->angles[i] - 2 * M_PI;
+            if (node->angles[i] < -1.0 * M_PI)
+                node->angles[i] = node->angles[i] + 2 * M_PI;
+            G_debug(5, "\tline = %d angle = %f", node->lines[i],
+                    node->angles[i]);
+        }
+        PQclear(res);
     }
-    PQclear(res);
     
     /* get node coordinates */
     if (SF_POINT != Vect__cache_feature_pg(wkb_data, FALSE, FALSE,
@@ -784,6 +827,9 @@ struct P_line *read_p_line(struct Plus_head *plus, int n,
     else {
         line->topo = dig_alloc_topo(line->type);
 
+        if ((line->type & GV_LINES) & (data->start_node < 0 || data->end_node < 0))
+            return NULL;
+        
         /* lines */
         if (line->type == GV_LINE) {
             struct P_topo_l *topo = (struct P_topo_l *)line->topo;
@@ -833,6 +879,117 @@ struct P_line *read_p_line(struct Plus_head *plus, int n,
 }
 
 /*!
+  \brief Read P_area structure
+  
+  \param plus pointer to Plus_head structure
+  \param n index (starts at 1)
+  \param lines_data lines array (see P_area struct)
+  \param centroid centroid id (see P_area struct)
+  \param isles_data lines array (see P_area struct)
+
+  \return pointer to P_area struct
+  \return NULL on error
+*/
+struct P_area *read_p_area(struct Plus_head *plus, int n,
+                           const char *lines_data, int centroid, const char *isles_data)
+{
+    int i;
+    int nlines, nisles;
+    char **lines, **isles;
+
+    struct P_area *area;
+
+    lines  = scan_array(lines_data);
+    nlines = G_number_of_tokens(lines);
+    isles  = scan_array(isles_data);
+    nisles = G_number_of_tokens(isles);
+
+    if (nlines < 1) {
+        G_warning(_("Area %d without boundary detected"), n);
+        return NULL;
+    }
+
+    G_debug(3, "read_p_area(): n = %d nlines = %d nisles = %d", n, nlines, nisles);
+    
+    /* allocate area */
+    area = dig_alloc_area();
+    dig_area_alloc_line(area, nlines);
+    dig_area_alloc_isle(area, nisles);
+
+    /* set lines */
+    area->n_lines = nlines;
+    for (i = 0; i < nlines; i++) {
+        area->lines[i] = atoi(lines[i]);
+    }
+
+    /* set isles */
+    area->n_isles = nisles;
+    for (i = 0; i < nisles; i++) {
+        area->isles[i] = atoi(isles[i]);
+    }
+
+    /* set centroid */
+    area->centroid = remap_line(plus, centroid);
+    
+    G_free_tokens(lines);
+    G_free_tokens(isles);
+
+    plus->Area[n] = area;
+
+    return area;
+}
+
+/*!
+  \brief Read P_isle structure
+  
+  \param plus pointer to Plus_head structure
+  \param n index (starts at 1)
+  \param lines_data lines array (see P_isle struct)
+  \param area area id (see P_isle struct)
+
+  \return pointer to P_isle struct
+  \return NULL on error
+*/
+struct P_isle *read_p_isle(struct Plus_head *plus, int n,
+                           const char *lines_data, int area)
+{
+    int i;
+    int nlines;
+    char **lines;
+
+    struct P_isle *isle;
+
+    lines  = scan_array(lines_data);
+    nlines = G_number_of_tokens(lines);
+
+    if (nlines < 1) {
+        G_warning(_("Isle %d without boundary detected"), n);
+        return NULL;
+    }
+
+    G_debug(3, "read_p_isle(): n = %d nlines = %d", n, nlines);
+    
+    /* allocate isle */
+    isle = dig_alloc_isle();
+    dig_isle_alloc_line(isle, nlines);
+    
+    /* set lines */
+    isle->n_lines = nlines;
+    for (i = 0; i < nlines; i++) {
+        isle->lines[i] = atoi(lines[i]);
+    }
+
+    /* set area */
+    isle->area = area;
+    
+    G_free_tokens(lines);
+
+    plus->Isle[n] = isle;
+
+    return isle;
+}
+
+/*!
   \brief Read topo from PostGIS topology schema -- header info only
 
   \param[in,out] plus pointer to Plus_head struct
@@ -890,6 +1047,19 @@ int load_plus_head(struct Format_info_pg *pg_info, struct Plus_head *plus)
             "AS node FROM \"%s\".edge GROUP BY end_node) AS foo",
             pg_info->toposchema_name, pg_info->toposchema_name);
     plus->n_nodes = Vect__execute_get_value_pg(pg_info->conn, stmt);
+    if (!pg_info->topo_geo_only) {
+        int n_nodes;
+        
+        /* check nodes consistency */
+        sprintf(stmt, "SELECT COUNT(*) FROM \"%s\".%s",
+                pg_info->toposchema_name, TOPO_TABLE_NODE);
+        n_nodes = Vect__execute_get_value_pg(pg_info->conn, stmt);
+        if (n_nodes != plus->n_nodes) {
+            G_warning(_("Different number of nodes detected (%d, %d)"),
+                      plus->n_nodes, n_nodes);
+            return -1;
+        }
+    }
     G_debug(3, "Vect_open_topo_pg(): n_nodes=%d", plus->n_nodes);
     
     /* lines (edges in PostGIS Topology model) */
@@ -905,6 +1075,19 @@ int load_plus_head(struct Format_info_pg *pg_info, struct Plus_head *plus)
             "SELECT COUNT(*) FROM \"%s\".face WHERE face_id > 0",
             pg_info->toposchema_name);
     plus->n_areas = Vect__execute_get_value_pg(pg_info->conn, stmt);
+    if (!pg_info->topo_geo_only) {
+        int n_areas;
+        
+        /* check areas consistency */
+        sprintf(stmt, "SELECT COUNT(*) FROM \"%s\".%s",
+                pg_info->toposchema_name, TOPO_TABLE_AREA);
+        n_areas = Vect__execute_get_value_pg(pg_info->conn, stmt);
+        if (n_areas != plus->n_areas) {
+            G_warning(_("Different number of areas detected (%d, %d)"),
+                      plus->n_areas, n_areas);
+            return -1;
+        }
+    }
     G_debug(3, "Vect_open_topo_pg(): n_areas=%d", plus->n_areas);
 
     /* isles (faces with face_id <=0 in PostGIS Topology model)
@@ -914,6 +1097,19 @@ int load_plus_head(struct Format_info_pg *pg_info, struct Plus_head *plus)
             "SELECT COUNT(*) FROM \"%s\".face WHERE face_id < 0",
             pg_info->toposchema_name);
     plus->n_isles = Vect__execute_get_value_pg(pg_info->conn, stmt);
+    if (!pg_info->topo_geo_only) {
+        int n_isles;
+        
+        /* check areas consistency */
+        sprintf(stmt, "SELECT COUNT(*) FROM \"%s\".%s",
+                pg_info->toposchema_name, TOPO_TABLE_ISLE);
+        n_isles = Vect__execute_get_value_pg(pg_info->conn, stmt);
+        if (n_isles != plus->n_isles) {
+            G_warning(_("Different number of areas detected (%d, %d)"),
+                      plus->n_isles, n_isles);
+            return -1;
+        }
+    }
     G_debug(3, "Vect_open_topo_pg(): n_isles=%d", plus->n_isles);
     
     /* number of features according the type */
@@ -982,6 +1178,7 @@ int Vect__load_plus_pg(struct Map_info *Map, int head_only)
     struct edge_data line_data;
     
     struct Format_info_pg *pg_info;
+    struct Format_info_offset *offset;
     struct Plus_head *plus;
     struct P_line *Line;
     struct P_area *Area;
@@ -992,6 +1189,7 @@ int Vect__load_plus_pg(struct Map_info *Map, int head_only)
   
     pg_info = &(Map->fInfo.pg);
     plus = &(Map->plus);
+    offset = &(pg_info->offset);
 
     if (load_plus_head(pg_info, plus) != 0)
         return -1;
@@ -1005,13 +1203,19 @@ int Vect__load_plus_pg(struct Map_info *Map, int head_only)
     /* read nodes (GRASS Topo)
        note: standalone nodes (ie. points/centroids) are ignored
     */
-    sprintf(stmt,
-            "SELECT node_id,geom FROM \"%s\".node WHERE node_id IN "
-            "(SELECT node FROM (SELECT start_node AS node FROM \"%s\".edge "
-            "GROUP BY start_node UNION ALL SELECT end_node AS node FROM "
-            "\"%s\".edge GROUP BY end_node) AS foo) ORDER BY node_id",
-            pg_info->toposchema_name, pg_info->toposchema_name,
-            pg_info->toposchema_name);
+    if (pg_info->topo_geo_only) 
+        sprintf(stmt,
+                "SELECT node_id,geom FROM \"%s\".node WHERE node_id IN "
+                "(SELECT node FROM (SELECT start_node AS node FROM \"%s\".edge "
+                "GROUP BY start_node UNION ALL SELECT end_node AS node FROM "
+                "\"%s\".edge GROUP BY end_node) AS foo) ORDER BY node_id",
+                pg_info->toposchema_name, pg_info->toposchema_name,
+                pg_info->toposchema_name);
+    else
+        sprintf(stmt, "SELECT node.node_id,geom,lines,angles FROM \"%s\".node AS node "
+                "join \"%s\".%s AS node_grass ON node.node_id = node_grass.node_id "
+                "ORDER BY node_id", pg_info->toposchema_name, pg_info->toposchema_name,
+                TOPO_TABLE_NODE);
     G_debug(2, "SQL: %s", stmt);
     res = PQexec(pg_info->conn, stmt);
     if (!res || PQresultStatus(res) != PGRES_TUPLES_OK ||
@@ -1026,11 +1230,18 @@ int Vect__load_plus_pg(struct Map_info *Map, int head_only)
 
     G_debug(3, "load_plus(): n_nodes = %d", plus->n_nodes);
     dig_alloc_nodes(plus, plus->n_nodes);
+    offset->array = (int *) G_malloc (sizeof(int) * plus->n_nodes);
+    offset->array_num = offset->array_alloc = plus->n_nodes;
     for (i = 0; i < plus->n_nodes; i++) {
         G_debug(5, "node: %d", i);
         id = atoi(PQgetvalue(res, i, 0));
         read_p_node(plus, i + 1, /* node index starts at 1 */
-                    id, (const char *) PQgetvalue(res, i, 1), pg_info);
+                    id, (const char *) PQgetvalue(res, i, 1),
+                    !pg_info->topo_geo_only ? (const char *) PQgetvalue(res, i, 2) : NULL,
+                    !pg_info->topo_geo_only ? (const char *) PQgetvalue(res, i, 3) : NULL,
+                    pg_info);
+        /* update offset */
+        offset->array[i] = id;
     }
     PQclear(res);
 
@@ -1045,14 +1256,20 @@ int Vect__load_plus_pg(struct Map_info *Map, int head_only)
     /* read PostGIS Topo standalone nodes (containing_face is null)
        -> points
     */
-    sprintf(stmt,
-            "SELECT node_id,geom FROM \"%s\".node WHERE containing_face "
-            "IS NULL AND node_id NOT IN "
-            "(SELECT node FROM (SELECT start_node AS node FROM \"%s\".edge "
-            "GROUP BY start_node UNION ALL SELECT end_node AS node FROM "
-            "\"%s\".edge GROUP BY end_node) AS foo) ORDER BY node_id",
-            pg_info->toposchema_name, pg_info->toposchema_name,
-            pg_info->toposchema_name);
+    if (pg_info->topo_geo_only)
+        sprintf(stmt,
+                "SELECT node_id,geom FROM \"%s\".node WHERE containing_face "
+                "IS NULL AND node_id NOT IN "
+                "(SELECT node FROM (SELECT start_node AS node FROM \"%s\".edge "
+                "GROUP BY start_node UNION ALL SELECT end_node AS node FROM "
+                "\"%s\".edge GROUP BY end_node) AS foo) ORDER BY node_id",
+                pg_info->toposchema_name, pg_info->toposchema_name,
+                pg_info->toposchema_name);
+    else
+        sprintf(stmt,
+                "SELECT node.node_id,geom FROM \"%s\".node AS node WHERE node_id NOT IN "
+                "(SELECT node_id FROM \"%s\".%s) AND containing_face IS NULL ORDER BY node_id",
+                pg_info->toposchema_name, pg_info->toposchema_name, TOPO_TABLE_NODE);
     G_debug(2, "SQL: %s", stmt);
     res = PQexec(pg_info->conn, stmt);
     if (!res || PQresultStatus(res) != PGRES_TUPLES_OK ||
@@ -1099,8 +1316,8 @@ int Vect__load_plus_pg(struct Map_info *Map, int head_only)
     ntuples = PQntuples(res);
     for (i = 0; i < ntuples; i++) {
         line_data.id         = atoi(PQgetvalue(res, i, 0));
-        line_data.start_node = atoi(PQgetvalue(res, i, 1));
-        line_data.end_node   = atoi(PQgetvalue(res, i, 2));
+        line_data.start_node = remap_node(offset, atoi(PQgetvalue(res, i, 1)));
+        line_data.end_node   = remap_node(offset, atoi(PQgetvalue(res, i, 2)));
         line_data.left_face  = atoi(PQgetvalue(res, i, 3));
         line_data.right_face = atoi(PQgetvalue(res, i, 4));
         line_data.wkb_geom   = (char *) PQgetvalue(res, i, 5);
@@ -1114,14 +1331,21 @@ int Vect__load_plus_pg(struct Map_info *Map, int head_only)
     /* read PostGIS Topo standalone nodes (containing_face is not null)
        -> centroids
     */
-    sprintf(stmt,
-            "SELECT node_id,geom,containing_face FROM \"%s\".node WHERE containing_face "
-            "IS NOT NULL AND node_id NOT IN "
-            "(SELECT node FROM (SELECT start_node AS node FROM \"%s\".edge "
-            "GROUP BY start_node UNION ALL SELECT end_node AS node FROM "
-            "\"%s\".edge GROUP BY end_node) AS foo) ORDER BY node_id",
-            pg_info->toposchema_name, pg_info->toposchema_name,
-            pg_info->toposchema_name);
+    if (pg_info->topo_geo_only)
+        sprintf(stmt,
+                "SELECT node_id,geom,containing_face FROM \"%s\".node WHERE containing_face "
+                "IS NOT NULL AND node_id NOT IN "
+                "(SELECT node FROM (SELECT start_node AS node FROM \"%s\".edge "
+                "GROUP BY start_node UNION ALL SELECT end_node AS node FROM "
+                "\"%s\".edge GROUP BY end_node) AS foo) ORDER BY node_id",
+                pg_info->toposchema_name, pg_info->toposchema_name,
+                pg_info->toposchema_name);
+    else
+        sprintf(stmt,
+                "SELECT node.node_id,geom,containing_face FROM \"%s\".node AS node WHERE "
+                "node_id NOT IN (SELECT node_id FROM \"%s\".%s) AND containing_face "
+                "IS NOT NULL ORDER BY node_id",
+                pg_info->toposchema_name, pg_info->toposchema_name, TOPO_TABLE_NODE);
     G_debug(2, "SQL: %s", stmt);
     res = PQexec(pg_info->conn, stmt);
     if (!res || PQresultStatus(res) != PGRES_TUPLES_OK ||
@@ -1145,29 +1369,89 @@ int Vect__load_plus_pg(struct Map_info *Map, int head_only)
         read_p_line(plus, id + i, &line_data, &(pg_info->cache));
     }
     PQclear(res);
+    plus->built = GV_BUILD_BASE;
+    
+    /* build areas */
+    if (pg_info->topo_geo_only) {
+        /* build areas for boundaries 
+           reset values -> build from scratch */
+        plus->n_areas = plus->n_isles = 0;
+        for (line = 1; line <= plus->n_lines; line++) {
+            Line = plus->Line[line]; /* centroids: Line is NULL */
+            if (!Line || Line->type != GV_BOUNDARY)
+                continue;
+            
+            for (i = 0; i < 2; i++) { /* for both sides build an area/isle */
+                side = i == 0 ? GV_LEFT : GV_RIGHT;
+                
+                G_debug(3, "Build area for line = %d, side = %d",
+                        id, side);
+                Vect_build_line_area(Map, line, side);
+            }
+        }
+    }
+    else {
+        /* read areas from 'area_grass' table */
+        sprintf(stmt,
+                "SELECT area_id,lines,centroid,isles FROM \"%s\".%s ORDER BY area_id",
+                pg_info->toposchema_name, TOPO_TABLE_AREA);
+        G_debug(2, "SQL: %s", stmt);
+        
+        res = PQexec(pg_info->conn, stmt);
+        if (!res || PQresultStatus(res) != PGRES_TUPLES_OK |
+            plus->n_areas != PQntuples(res)) {
+            if (res)
+                PQclear(res);
+            return -1;
+        }
 
-    /* build areas for boundaries
-       reset values -> build from scratch */
-    plus->n_areas = plus->n_isles = 0;
-    for (line = 1; line <= plus->n_lines; line++) {
-        Line = plus->Line[line]; /* centroids: Line is NULL */
-        if (!Line || Line->type != GV_BOUNDARY)
-            continue;
         
-        for (i = 0; i < 2; i++) { /* for both sides build an area/isle */
-            side = i == 0 ? GV_LEFT : GV_RIGHT;
-            
-            G_debug(3, "Build area for line = %d, side = %d",
-                    id, side);
-            Vect_build_line_area(Map, line, side);
+        dig_alloc_areas(plus, plus->n_areas);
+        G_zero(plus->Area, sizeof(struct P_area *) * (plus->n_areas + 1)); /* index starts at 1 */
+        G_debug(3, "Vect_open_topo_pg(): n_areas=%d", plus->n_areas);
+        
+        for (i = 0; i < plus->n_areas; i++) {
+            read_p_area(plus, i + 1, (char *)PQgetvalue(res, i, 1),
+                        atoi(PQgetvalue(res, i, 2)), (char *)PQgetvalue(res, i, 3));
         }
+        PQclear(res);
     }
     plus->built = GV_BUILD_AREAS;
 
-    /* TODO: attach isles */
+    /* attach isles */
+    if (pg_info->topo_geo_only) {
+        plus->n_isles = 0; /* reset isles */
+        G_warning(_("To be implented: isles not attached in Topo-Geo-only mode"));
+    }
+    else {
+        /* read isles from 'isle_grass' table */
+        sprintf(stmt,
+                "SELECT isle_id,lines,area FROM \"%s\".%s ORDER BY isle_id",
+                pg_info->toposchema_name, TOPO_TABLE_ISLE);
+        G_debug(2, "SQL: %s", stmt);
+        
+        res = PQexec(pg_info->conn, stmt);
+        if (!res || PQresultStatus(res) != PGRES_TUPLES_OK ||
+            plus->n_isles != PQntuples(res)) {
+            if (res)
+                PQclear(res);
+            return -1;
+        }
+        
+        dig_alloc_isles(plus, plus->n_isles);
+        G_zero(plus->Isle, sizeof(struct P_isle *) * (plus->n_isles + 1)); /* index starts at 1 */
+        G_debug(3, "Vect_open_topo_pg(): n_isles=%d", plus->n_isles);
+        
+        for (i = 0; i < plus->n_isles; i++) {
+            read_p_isle(plus, i + 1, (char *)PQgetvalue(res, i, 1),
+                        atoi(PQgetvalue(res, i, 2)));
+        }
+        PQclear(res);
+    }
     plus->built = GV_BUILD_ATTACH_ISLES;
     
-    /* attach centroids */
+    /* attach centroids (disabled: centroids are already attached) */
+#if 0
     if (plus->n_areas > 0) {
         struct P_topo_c *topo;
         
@@ -1184,6 +1468,7 @@ int Vect__load_plus_pg(struct Map_info *Map, int head_only)
             Area->centroid = Line->offset;
         }
     }
+#endif
     plus->built = GV_BUILD_CENTROIDS;
     
     /* done */
@@ -1206,4 +1491,87 @@ void notice_processor(void *arg, const char *message)
         fprintf(stderr, "%s", message);
     }
 }
+
+/*!
+  \brief Scan string array
+
+  Creates tokens based on string array, eg. '{1, 2, 3}' become
+  [1,2,3].
+
+  Allocated tokes should be freed by G_free_tokens().
+
+  \param sArray string array
+
+  \return tokens
+*/
+char **scan_array(const char *sarray)
+{
+    char *buf, **tokens;
+    int i, len;
+    
+    /* remove '{}' */
+    buf = (char *)G_malloc(sizeof(sarray) - 2);
+    len = strlen(sarray) - 1;
+    for (i = 1; i < len; i++)
+        buf[i-1] = sarray[i];
+    buf[len-1] = '\0';
+    
+    tokens = G_tokenize(buf, ",");
+    G_free(buf);
+    
+    return tokens;
+}
+
+/*!
+  \brief Get node id from offset
+
+  \todo speed up
+
+  \param offset pointer to Format_info_offset struct
+  \param node node to find
+
+  \return node id
+  \return -1 not found
+*/
+int remap_node(const struct Format_info_offset *offset, int node)
+{
+    int i;
+    
+    for (i = 0; i < offset->array_num; i++) {
+        if (offset->array[i] == node)
+            return i + 1; /* node id starts at 1 */
+    }
+
+    return -1;
+}
+
+/*!
+  \brief Get line id from offset
+
+  \todo speed up
+
+  \param plus pointer to Plus_head struct
+  \param line line to find
+
+  \return line id
+  \return -1 not found
+*/
+int remap_line(const struct Plus_head* plus, int line)
+{
+    int i;
+    
+    struct P_line *Line;
+    
+    for (i = 1; i <= plus->n_lines; i++) {
+        Line = plus->Line[i];
+        
+        if (!Line)
+            continue;
+        
+        if ((int) Line->offset == line)
+            return i;
+    }
+
+    return -1;
+}
 #endif

lib/vector/Vlib/pg_local_proto.h

@@ -13,6 +13,10 @@
 #define TOPO_ID     "topology_id"
 #define TOPO_TABLE  "grass"
 #define TOPO_BBOX   "bbox"
+#define TOPO_TABLE_NUM  3
+#define TOPO_TABLE_NODE "node_grass"
+#define TOPO_TABLE_AREA "area_grass"
+#define TOPO_TABLE_ISLE "isle_grass"
 
 #define SWAP32(x) \
         ((unsigned int)( \

lib/vector/Vlib/write_pg.c

@@ -458,6 +458,7 @@ int V2_delete_line_pg(struct Map_info *Map, int line)
    GV_BUILD_BASE. PostGIS Topology schema must be defined.
 
    \param Map pointer to Map_info structure
+   \param node node id (starts at 1)
    \param points pointer to line_pnts structure
    
    \return 0 on success
@@ -764,6 +765,8 @@ int check_schema(const struct Format_info_pg *pg_info)
   - create topology schema
   - add topology column to the feature table
   
+  \todo Add constraints for grass-like tables
+  
   \param pg_info pointer to Format_info_pg
 
   \return 0 on success
@@ -796,16 +799,24 @@ int create_topo_schema(struct Format_info_pg *pg_info, int with_z)
         fclose(fp);
 
         /* tolerance */
-        p = G_find_key_value("tolerance", key_val);
+        p = G_find_key_value("topo_tolerance", key_val);
         if (p)
             tolerance = atof(p);
-
+        G_debug(1, "PG: tolerance: %f", tolerance);
+        
         /* topogeom column */
-        p = G_find_key_value("topogeom_column", key_val);
+        p = G_find_key_value("topogeom_name", key_val);
         if (p)
             pg_info->topogeom_column = G_store(p);
         else
             pg_info->topogeom_column = G_store(TOPOGEOM_COLUMN);
+        G_debug(1, "PG: topogeom_column :%s", pg_info->topogeom_column);
+
+        /* topo-geo only (default: no) */
+        p = G_find_key_value("topo_geo_only", key_val);
+        if (p && G_strcasecmp(p, "yes") == 0)
+            pg_info->topo_geo_only = TRUE;
+        G_debug(1, "PG: topo_geo_only :%d", pg_info->topo_geo_only);
     }
 
     /* begin transaction (create topo schema) */
@@ -848,6 +859,53 @@ int create_topo_schema(struct Format_info_pg *pg_info, int with_z)
         return -1;
     }
 
+    /* create additional tables in topological schema to store
+       GRASS topology in DB */
+    if (!pg_info->topo_geo_only) {
+        /* (1) create 'node_grass' (see P_node struct)
+           
+           todo: add constraints for lines and angles
+        */
+        sprintf(stmt, "CREATE TABLE \"%s\".%s (node_id SERIAL PRIMARY KEY, "
+                "lines integer[], angles float[])", pg_info->toposchema_name, TOPO_TABLE_NODE);
+        if (Vect__execute_pg(pg_info->conn, stmt) == -1) {
+            Vect__execute_pg(pg_info->conn, "ROLLBACK");
+            return -1;
+        }
+
+        sprintf(stmt, "ALTER TABLE \"%s\".%s ADD CONSTRAINT node_exists "
+                "FOREIGN KEY (node_id) REFERENCES \"%s\".node (node_id)",
+                pg_info->toposchema_name, TOPO_TABLE_NODE, pg_info->toposchema_name);
+        if (Vect__execute_pg(pg_info->conn, stmt) == -1) {
+            Vect__execute_pg(pg_info->conn, "ROLLBACK");
+            return -1;
+        }
+
+        /* (2) create 'area_grass' (see P_area struct)
+           
+           todo: add constraints for lines, centtroid and isles
+        */
+        sprintf(stmt, "CREATE TABLE \"%s\".%s (area_id SERIAL PRIMARY KEY, "
+                "lines integer[], centroid integer, isles integer[])",
+                pg_info->toposchema_name, TOPO_TABLE_AREA);
+        if (Vect__execute_pg(pg_info->conn, stmt) == -1) {
+            Vect__execute_pg(pg_info->conn, "ROLLBACK");
+            return -1;
+        }
+
+        /* (3) create 'isle_grass' (see P_isle struct)
+           
+           todo: add constraints for lines and area
+        */
+        sprintf(stmt, "CREATE TABLE \"%s\".%s (isle_id SERIAL PRIMARY KEY, "
+                "lines integer[], area integer)",
+                pg_info->toposchema_name, TOPO_TABLE_ISLE);
+        if (Vect__execute_pg(pg_info->conn, stmt) == -1) {
+            Vect__execute_pg(pg_info->conn, "ROLLBACK");
+            return -1;
+        }
+    }
+
     /* close transaction (create topo schema) */
     if (Vect__execute_pg(pg_info->conn, "COMMIT") == -1) {
         return -1;
@@ -1148,13 +1206,14 @@ off_t write_line_sf(struct Map_info *Map, int type,
   \param is_node TRUE for nodes (written as points)
   
   \return feature id (build level >= GV_BUILD_BASE otherwise 0)
+  \return 0 for nodes
   \return -1 on error
 */
 off_t write_line_tp(struct Map_info *Map, int type, int is_node,
                     const struct line_pnts *points,
                     const struct line_cats *cats)
 {
-    int line, cat;
+    int line, cat, line_id;
 
     struct field_info *Fi;
     struct Format_info_pg *pg_info;
@@ -1242,11 +1301,30 @@ off_t write_line_tp(struct Map_info *Map, int type, int is_node,
        - feature table for simple features
        - feature table and topo schema for topological access
     */
-    if (-1 == write_feature(Map, line, type, &points, 1, cat, Fi)) {
+    line_id = write_feature(Map, line, type, &points, 1, cat, Fi);
+    if (line_id < 0) {
         Vect__execute_pg(pg_info->conn, "ROLLBACK");
         return -1;
     }
 
+    /* update offset array for nodes */    
+    if (is_node) {
+        int node;
+        
+        struct Format_info_offset *offset;
+
+        offset = &(pg_info->offset);
+        
+        node = abs(line);
+        if (node > offset->array_alloc) {
+            offset->array_alloc += 1000;
+            offset->array = (int *) G_realloc (offset->array, offset->array_alloc * sizeof(int));
+        }
+        
+        offset->array_num = node; 
+        offset->array[node-1] = (int) line_id; /* node id starts at 1 */
+    }
+
     /* update PostGIS-line topo */
     if (plus->built >= GV_BUILD_AREAS && type == GV_BOUNDARY)
         update_topo_face(Map, line); /* TODO: avoid extra statements */
@@ -1638,8 +1716,9 @@ char *line_to_wkb(struct Format_info_pg *pg_info,
    \param cat category number (-1 for no category)
    \param Fi pointer to field_info (attributes to copy, NULL for no attributes)
 
+   \return topo_id for PostGIS Topology 
+   \return 0 for simple features access
    \return -1 on error
-   \retirn 0 on success
  */
 int write_feature(struct Map_info *Map, int line, int type,
                   const struct line_pnts **points, int nparts,
@@ -1708,7 +1787,7 @@ int write_feature(struct Map_info *Map, int line, int type,
     }
     G_free(geom_data);
     
-    return 0;
+    return pg_info->toposchema_name ? topo_id : 0;
 }
 
 /*!
@@ -1954,6 +2033,12 @@ int insert_topo_element(struct Map_info *Map, int id, int type,
                 return -1;
             }
             Line = Map->plus.Line[topo_id];
+
+            if (Line->type & GV_POINTS) {
+                /* set topo_id for points */
+                topo_id = Vect_get_num_primitives(Map, GV_POINTS) + Vect_get_num_nodes(Map);
+            }
+
         }
         else if (id < 0) { /* node */
             topo_id = abs(id);
@@ -1965,6 +2050,10 @@ int insert_topo_element(struct Map_info *Map, int id, int type,
                 G_warning(_("Invalid node %d (%d)"), topo_id, Map->plus.n_nodes);
                 return -1;
             }
+
+            /* increment topo_id - also points and centroids are
+             * stored in 'node' table */
+            topo_id += Vect_get_num_primitives(Map, GV_POINTS);
         }
     }
 
@@ -1997,7 +2086,11 @@ int insert_topo_element(struct Map_info *Map, int id, int type,
                    pg_info->toposchema_name, geom_data);
 #else
         int n1, n2, nle, nre;
-
+        
+        struct Format_info_offset *offset;
+        
+        offset = &(pg_info->offset);
+        
         if (id == 0) {
             /* get edge_id */
             sprintf(stmt_id, "SELECT nextval('\"%s\".edge_data_edge_id_seq')",
@@ -2044,12 +2137,15 @@ int insert_topo_element(struct Map_info *Map, int id, int type,
         G_debug(3, "new edge: id=%d next_left_edge=%d next_right_edge=%d",
                 topo_id, nle, nre);
         
+        if (n1 > offset->array_num || n2 > offset->array_num) /* node id starts at 1 */
+            return -1;
+        
         /* build insert statement */
         G_asprintf(&stmt, "INSERT INTO \"%s\".edge_data (edge_id, start_node, end_node, "
                    "next_left_edge, abs_next_left_edge, next_right_edge, abs_next_right_edge, "
                    "left_face, right_face, geom) VALUES "
                    "(%d, %d, %d, %d, %d, %d, %d, 0, 0, '%s'::GEOMETRY)",
-                   pg_info->toposchema_name, topo_id, n1, n2,
+                   pg_info->toposchema_name, topo_id, offset->array[n1-1], offset->array[n2-1],
                    nle, abs(nle), nre, abs(nre), geom_data);
 #endif
         break;