[GRASS-SVN] r51019 - in grass/trunk: include/vect lib/vector/Vlib
svn_grass at osgeo.org
svn_grass at osgeo.org
Wed Mar 7 17:12:18 EST 2012
Author: martinl
Date: 2012-03-07 14:12:18 -0800 (Wed, 07 Mar 2012)
New Revision: 51019
Modified:
grass/trunk/include/vect/dig_structs.h
grass/trunk/lib/vector/Vlib/build_ogr.c
grass/trunk/lib/vector/Vlib/build_pg.c
grass/trunk/lib/vector/Vlib/build_sfa.c
grass/trunk/lib/vector/Vlib/field.c
grass/trunk/lib/vector/Vlib/pg_local_proto.h
grass/trunk/lib/vector/Vlib/poly.c
grass/trunk/lib/vector/Vlib/read.c
grass/trunk/lib/vector/Vlib/read_ogr.c
grass/trunk/lib/vector/Vlib/read_pg.c
Log:
vlib/pg: support reading geometry collections
Modified: grass/trunk/include/vect/dig_structs.h
===================================================================
--- grass/trunk/include/vect/dig_structs.h 2012-03-07 22:03:17 UTC (rev 51018)
+++ grass/trunk/include/vect/dig_structs.h 2012-03-07 22:12:18 UTC (rev 51019)
@@ -409,11 +409,42 @@
GRASS. This is not used for GV_CENTROID. Because one feature may
contain more elements (geometry collection also recursively),
offset for one line may be stored in more records. First record
- is FID, next records are part indexes if necessary. Example:
+ is FID, next records are part indexes if necessary.
+
+ Example 1:
+
5. ring in 3. polygon in 7. feature (multipolygon) of geometry
collection which has FID = 123 123 (feature 123: geometry
colletion) 6 (7. feature in geometry collection: multiPolygon) 2
(3. polygon) 4 (5. ring in the polygon)
+
+ Example 2: geometry collection FID '1' containing one point, one
+ linestring and one polygon
+
+ \verbatim
+ Offset:
+
+ idx offset note
+ ----------------
+ 0 1 FID
+ 1 0 first part (point)
+
+ 2 1 FID
+ 3 1 second part (linestring)
+
+ 4 1 FID
+ 5 2 third part (polygon)
+ 6 0 first ring of polygon
+
+ Topology:
+
+ line idx
+ -----------------
+ 1 0 point
+ 2 2 line
+ 3 4 boundary
+ 4 1 centroid read from topo (idx == FID)
+ \endverbatim
*/
int *array;
/*!
@@ -428,6 +459,39 @@
};
/*!
+ \brief Lines cache for reading feature (non-native formats)
+*/
+struct Format_info_cache {
+ /*!
+ \brief Lines array
+
+ Some features requires more allocated lines (eg. polygon
+ with more rings, multipoint, or geometrycollection)
+ */
+ struct line_pnts **lines;
+ /*!
+ \brief List of line types (GV_POINT, GV_LINE, ...)
+ */
+ int *lines_types;
+ /*!
+ \brief Number of allocated lines in cache
+ */
+ int lines_alloc;
+ /*!
+ \brief Number of lines which forms current feature
+ */
+ int lines_num;
+ /*!
+ \brief Next line to be read from cache
+ */
+ int lines_next;
+ /*!
+ \brief Feature id
+ */
+ long fid;
+};
+
+/*!
\brief Non-native format info (OGR)
\todo Structure size should not change depending on compilation I
@@ -484,34 +548,10 @@
char **layer_options;
/*!
- \brief Lines cache (per feature)
+ \brief Lines cache for reading feature
*/
- struct {
- /*!
- \brief Lines array
-
- Some features requires more allocated lines (eg. polygon
- with more rings, multipoint, or geometrycollection)
- */
- struct line_pnts **lines;
- /*!
- \brief List of line types
- */
- int *lines_types;
- /*!
- \brief Number of allocated lines in cache
- */
- int lines_alloc;
- /*!
- \brief Number of lines which forms current feature
- */
- int lines_num;
- /*!
- \brief Next line to be read from cache
- */
- int lines_next;
- } cache;
-
+ struct Format_info_cache cache;
+
/*!
\brief Cache to avoid repeated reading (level 2)
@@ -581,39 +621,11 @@
\brief Next line to be read for sequential access
*/
int next_line;
-
+
/*!
- \brief Lines cache (per feature)
+ \brief Lines cache for reading feature
*/
- struct {
- /*!
- \brief Lines array
-
- Some features requires more allocated lines (eg. polygon
- with more rings, multipoint, or geometrycollection)
- */
- struct line_pnts **lines;
- /*!
- \brief List of line types
- */
- int *lines_types;
- /*!
- \brief Number of allocated lines in cache
- */
- int lines_alloc;
- /*!
- \brief Number of lines which forms current feature
- */
- int lines_num;
- /*!
- \brief Next line to be read from cache
- */
- int lines_next;
- /*!
- \brief Feature id
- */
- long fid;
- } cache;
+ struct Format_info_cache cache;
/*!
\brief Offset list used for building pseudo-topology
Modified: grass/trunk/lib/vector/Vlib/build_ogr.c
===================================================================
--- grass/trunk/lib/vector/Vlib/build_ogr.c 2012-03-07 22:03:17 UTC (rev 51018)
+++ grass/trunk/lib/vector/Vlib/build_ogr.c 2012-03-07 22:12:18 UTC (rev 51019)
@@ -79,8 +79,8 @@
"PostgreSQL") == 0)
G_warning(_("Feature table <%s> has no primary key defined"),
ogr_info->layer_name);
- G_warning(_("Random read is not supported by OGR for this layer, "
- "unable to build topology"));
+ G_warning(_("Random read is not supported by OGR for this layer. "
+ "Unable to build topology."));
return 0;
}
@@ -159,6 +159,8 @@
offset->array_num, &fp))
return 0;
+ G_debug(3, "Vect_save_fidx(): offset_num = %d", offset->array_num);
+
fclose(fp.file);
return 1;
Modified: grass/trunk/lib/vector/Vlib/build_pg.c
===================================================================
--- grass/trunk/lib/vector/Vlib/build_pg.c 2012-03-07 22:03:17 UTC (rev 51018)
+++ grass/trunk/lib/vector/Vlib/build_pg.c 2012-03-07 22:12:18 UTC (rev 51019)
@@ -70,6 +70,14 @@
G_warning(_("No DB connection"));
return 0;
}
+
+ if (!pg_info->fid_column) {
+ G_warning(_("Feature table <%s> has no primary key defined"),
+ pg_info->table_name);
+ G_warning(_("Random read is not supported by OGR for this layer. "
+ "Unable to build topology."));
+ return 0;
+ }
G_message(_("Using external data format '%s' (feature type '%s')"),
Vect_get_finfo_format_info(Map),
Modified: grass/trunk/lib/vector/Vlib/build_sfa.c
===================================================================
--- grass/trunk/lib/vector/Vlib/build_sfa.c 2012-03-07 22:03:17 UTC (rev 51018)
+++ grass/trunk/lib/vector/Vlib/build_sfa.c 2012-03-07 22:12:18 UTC (rev 51019)
@@ -52,8 +52,9 @@
static int add_geometry_pg(struct Plus_head *,
struct Format_info_pg *,
- SF_FeatureType, int, int,
- struct geom_parts *);
+ struct feat_parts *, int,
+ int, int,
+ struct geom_parts *);
static void build_pg(struct Map_info *, int);
#endif
@@ -176,9 +177,10 @@
dig_cidx_add_cat(plus, 0, 0, line, type);
}
- if (type != GV_CENTROID) /* because centroids are read from topology, not from layer */
+ /* because centroids are read from topology, not from layer */
+ if (type != GV_CENTROID)
add_parts_to_offset(offset, parts);
-
+
return line;
}
@@ -188,15 +190,19 @@
*/
int add_geometry_pg(struct Plus_head *plus,
struct Format_info_pg *pg_info,
- SF_FeatureType ftype, int FID, int build,
- struct geom_parts *parts)
+ struct feat_parts *fparts, int ipart,
+ int FID, int build, struct geom_parts *parts)
{
- int line, iPart, area, isle, outer_area, ret;
+ int line, i, idx, area, isle, outer_area, ret;
int lines[1];
double area_size, x, y;
+ SF_FeatureType ftype;
struct bound_box box;
struct Format_info_offset *offset;
-
+ struct line_pnts *line_i;
+
+ ftype = fparts->ftype[ipart];
+
G_debug(4, "add_geometry_pg() FID = %d ftype = %d", FID, ftype);
offset = &(pg_info->offset);
@@ -206,28 +212,34 @@
switch(ftype) {
case SF_POINT:
G_debug(4, "Point");
- add_line(plus, offset, GV_POINT, pg_info->cache.lines[0], FID, parts);
+ line_i = pg_info->cache.lines[fparts->idx[ipart]];
+ add_line(plus, offset, GV_POINT, line_i,
+ FID, parts);
break;
case SF_LINESTRING:
G_debug(4, "LineString");
- add_line(plus, offset, GV_LINE, pg_info->cache.lines[0], FID, parts);
+ line_i = pg_info->cache.lines[fparts->idx[ipart]];
+ add_line(plus, offset, GV_LINE, line_i,
+ FID, parts);
break;
case SF_POLYGON:
G_debug(4, "Polygon");
/* register boundaries */
- for (iPart = 0; iPart < pg_info->cache.lines_num; iPart++) {
- add_part(parts, iPart);
- line = add_line(plus, offset, GV_BOUNDARY, pg_info->cache.lines[iPart],
- FID, parts);
+ idx = fparts->idx[ipart];
+ for (i = 0; i < fparts->nlines[ipart]; i++) {
+ line_i = pg_info->cache.lines[idx++];
+ add_part(parts, i);
+ line = add_line(plus, offset, GV_BOUNDARY,
+ line_i, FID, parts);
del_part(parts);
if (build < GV_BUILD_AREAS)
continue;
/* add area (each inner ring is also area) */
- dig_line_box(pg_info->cache.lines[iPart], &box);
- dig_find_area_poly(pg_info->cache.lines[iPart], &area_size);
+ dig_line_box(line_i, &box);
+ dig_find_area_poly(line_i, &area_size);
if (area_size > 0) /* area clockwise */
lines[0] = line;
@@ -244,7 +256,7 @@
if (build < GV_BUILD_ATTACH_ISLES)
continue;
- if (iPart == 0) { /* outer ring */
+ if (i == 0) { /* outer ring */
outer_area = area;
}
else { /* inner ring */
@@ -259,9 +271,9 @@
if (build >= GV_BUILD_CENTROIDS) {
/* create virtual centroid */
- ret = Vect_get_point_in_poly_isl((const struct line_pnts *) pg_info->cache.lines[0],
- (const struct line_pnts **) pg_info->cache.lines + 1,
- pg_info->cache.lines_num - 1, &x, &y);
+ ret = Vect_get_point_in_poly_isl((const struct line_pnts *) pg_info->cache.lines[fparts->idx[ipart]],
+ (const struct line_pnts **) pg_info->cache.lines[fparts->idx[ipart]] + 1,
+ fparts->nlines[ipart] - 1, &x, &y);
if (ret < -1) {
G_warning(_("Unable to calculate centroid for area %d"),
outer_area);
@@ -270,11 +282,12 @@
struct P_area *Area;
struct P_topo_c *topo;
struct P_line *Line;
-
+ struct line_pnts *line_c;
+
G_debug(4, " Centroid: %f, %f", x, y);
- Vect_reset_line(pg_info->cache.lines[0]);
- Vect_append_point(pg_info->cache.lines[0], x, y, 0.0);
- line = add_line(plus, offset, GV_CENTROID, pg_info->cache.lines[0], FID, parts);
+ line_c = Vect_new_line_struct();
+ Vect_append_point(line_c, x, y, 0.0);
+ line = add_line(plus, offset, GV_CENTROID, line_c, FID, parts);
Line = plus->Line[line];
topo = (struct P_topo_c *)Line->topo;
@@ -283,6 +296,7 @@
/* register centroid to area */
Area = plus->Area[outer_area];
Area->centroid = line;
+ Vect_destroy_line_struct(line_c);
}
}
break;
@@ -299,19 +313,21 @@
*/
void build_pg(struct Map_info *Map, int build)
{
- int iFeature, FID, nrecords;
+ int iFeature, ipart, fid, nrecords;
char stmt[DB_SQL_MAX];
- SF_FeatureType ftype;
+ char *wkb_data;
struct Format_info_pg *pg_info;
+ struct feat_parts fparts;
struct geom_parts parts;
pg_info = &(Map->fInfo.pg);
/* initialize data structures */
init_parts(&parts);
-
+ G_zero(&fparts, sizeof(struct feat_parts));
+
/* get records */
sprintf(stmt, "SELECT %s,%s FROM %s",
pg_info->fid_column,
@@ -330,23 +346,29 @@
G_debug(4, "build_pg(): nrecords = %d", nrecords);
for (iFeature = 0; iFeature < nrecords; iFeature++) {
/* get feature id */
- FID = atoi(PQgetvalue(pg_info->res, iFeature, 0));
+ fid = atoi(PQgetvalue(pg_info->res, iFeature, 0));
+ wkb_data = PQgetvalue(pg_info->res, iFeature, 1);
+
/* cache feature (lines) */
- ftype = cache_feature(PQgetvalue(pg_info->res, iFeature, 1),
- FALSE, pg_info);
-
- if (pg_info->cache.lines_num < 1) {
- G_warning(_("Feature %d without geometry ignored"), iFeature);
- continue;
- }
-
- G_debug(3, "Feature %d: cat=%d", iFeature, FID);
+ cache_feature(wkb_data,
+ FALSE, &(pg_info->cache), &fparts);
/* register topo */
reset_parts(&parts);
- add_part(&parts, FID);
- add_geometry_pg(&(Map->plus), pg_info,
- ftype, FID, build, &parts);
+ add_part(&parts, fid);
+ for (ipart = 0; ipart < fparts.n_parts; ipart++) {
+ if (fparts.nlines[ipart] < 1) {
+ G_warning(_("Feature %d without geometry skipped"), fid);
+ continue;
+ }
+
+ G_debug(4, "Feature: fid = %d part = %d", fid, ipart);
+
+ add_part(&parts, ipart);
+ add_geometry_pg(&(Map->plus), pg_info, &fparts, ipart,
+ fid, build, &parts);
+ del_part(&parts);
+ }
}
Map->plus.built = GV_BUILD_BASE;
@@ -594,14 +616,14 @@
hGeom = OGR_F_GetGeometryRef(hFeature);
if (hGeom == NULL) {
- G_warning(_("Feature %d without geometry ignored"), iFeature);
+ G_warning(_("Feature %d without geometry skipped"), iFeature);
OGR_F_Destroy(hFeature);
continue;
}
FID = (int) OGR_F_GetFID(hFeature);
if (FID == OGRNullFID) {
- G_warning(_("OGR feature %d without ID ignored"), iFeature);
+ G_warning(_("OGR feature %d without ID skipped"), iFeature);
OGR_F_Destroy(hFeature);
continue;
}
Modified: grass/trunk/lib/vector/Vlib/field.c
===================================================================
--- grass/trunk/lib/vector/Vlib/field.c 2012-03-07 22:03:17 UTC (rev 51018)
+++ grass/trunk/lib/vector/Vlib/field.c 2012-03-07 22:12:18 UTC (rev 51019)
@@ -817,7 +817,8 @@
pg_info = &(Map->fInfo.pg);
if (!pg_info->fid_column) {
- G_warning(_("No FID column defined"));
+ G_warning(_("Feature table <%s> has no primary key defined. "
+ "Unable to define DB links."), pg_info->table_name);
return -1;
}
G_debug(3, "Using FID column <%s>", pg_info->fid_column);
Modified: grass/trunk/lib/vector/Vlib/pg_local_proto.h
===================================================================
--- grass/trunk/lib/vector/Vlib/pg_local_proto.h 2012-03-07 22:03:17 UTC (rev 51018)
+++ grass/trunk/lib/vector/Vlib/pg_local_proto.h 2012-03-07 22:12:18 UTC (rev 51019)
@@ -4,10 +4,22 @@
#ifdef HAVE_POSTGRES
#include <libpq-fe.h>
+/* used for building pseudo-topology (requires some extra information
+ * about lines in cache) */
+struct feat_parts
+{
+ int a_parts; /* number of allocated items */
+ int n_parts; /* number of parts which forms given feature */
+ SF_FeatureType *ftype; /* simple feature type */
+ int *nlines; /* number of lines used in cache */
+ int *idx; /* index in cache where to start */
+};
+
/* functions used in *_pg.c files */
int execute(PGconn *, const char *);
-int cache_feature(const char *, int,
- struct Format_info_pg *);
+SF_FeatureType cache_feature(const char *, int,
+ struct Format_info_cache *,
+ struct feat_parts *);
#endif /* HAVE_POSTGRES */
Modified: grass/trunk/lib/vector/Vlib/poly.c
===================================================================
--- grass/trunk/lib/vector/Vlib/poly.c 2012-03-07 22:03:17 UTC (rev 51018)
+++ grass/trunk/lib/vector/Vlib/poly.c 2012-03-07 22:12:18 UTC (rev 51019)
@@ -398,22 +398,21 @@
\brief Get point inside polygon but outside the islands specifiled in IPoints.
Take a line and intersect it with the polygon and any islands.
- sort the list of X values from these intersections. This will
- be a list of segments alternating IN/OUT/IN/OUt of the polygon.
- Pick the largest IN segment and take the midpoint.
+ sort the list of X values from these intersections. This will be a
+ list of segments alternating IN/OUT/IN/OUT of the polygon. Pick the
+ largest IN segment and take the midpoint.
- \param Points polygon
- \param IPoints isles
+ \param Points polygon (boundary)
+ \param IPoints isles (list of isle boundaries)
\param n_isles number of isles
\param[out] att_x,att_y point coordinates
\return 0 on success
\return -1 on error
*/
-int
-Vect_get_point_in_poly_isl(const struct line_pnts *Points,
- const struct line_pnts **IPoints, int n_isles,
- double *att_x, double *att_y)
+int Vect_get_point_in_poly_isl(const struct line_pnts *Points,
+ const struct line_pnts **IPoints, int n_isles,
+ double *att_x, double *att_y)
{
static struct line_pnts *Intersects;
static int first_time = 1;
@@ -771,7 +770,6 @@
static struct line_pnts *Points;
struct bound_box lbox;
const struct Plus_head *Plus;
- struct P_line *Line;
struct P_isle *Isle;
G_debug(3, "Vect_point_in_island(): x = %f y = %f isle = %d", X, Y, isle);
@@ -791,8 +789,6 @@
for (i = 0; i < Isle->n_lines; i++) {
line = abs(Isle->lines[i]);
- Line = Plus->Line[line];
-
/* this is slow, but the fastest of all alternatives */
Vect_get_line_box(Map, line, &lbox);
Modified: grass/trunk/lib/vector/Vlib/read.c
===================================================================
--- grass/trunk/lib/vector/Vlib/read.c 2012-03-07 22:03:17 UTC (rev 51018)
+++ grass/trunk/lib/vector/Vlib/read.c 2012-03-07 22:12:18 UTC (rev 51019)
@@ -126,8 +126,8 @@
ret = (*Read_next_line_array[Map->format][Map->level]) (Map, line_p,
line_c);
if (ret == -1)
- G_fatal_error(_("Unable to read vector map <%s>"),
- Vect_get_full_name(Map));
+ G_fatal_error(_("Unable to read feature %d vector map <%s>"),
+ Map->next_line, Vect_get_full_name(Map));
return ret;
}
@@ -168,7 +168,7 @@
ret = (*Read_line_array[Map->format]) (Map, line_p, line_c, line);
if (ret == -1)
- G_fatal_error(_("Unable to read feature %d vector map <%s>"),
+ G_fatal_error(_("Unable to read feature %d from vector map <%s>"),
line, Vect_get_full_name(Map));
return ret;
Modified: grass/trunk/lib/vector/Vlib/read_ogr.c
===================================================================
--- grass/trunk/lib/vector/Vlib/read_ogr.c 2012-03-07 22:03:17 UTC (rev 51018)
+++ grass/trunk/lib/vector/Vlib/read_ogr.c 2012-03-07 22:12:18 UTC (rev 51019)
@@ -180,61 +180,61 @@
struct line_pnts *line_p, struct line_cats *line_c, off_t offset)
{
#ifdef HAVE_OGR
- long FID;
+ long fid;
int type;
OGRGeometryH hGeom;
struct Format_info_ogr *ogr_info;
ogr_info = &(Map->fInfo.ogr);
- G_debug(4, "V1_read_line_ogr(): offset = %lu offset_num = %lu",
+ G_debug(3, "V1_read_line_ogr(): offset = %lu offset_num = %lu",
(long) offset, (long) ogr_info->offset.array_num);
if (offset >= ogr_info->offset.array_num)
- return -2;
+ return -2; /* nothing to read */
if (line_p != NULL)
Vect_reset_line(line_p);
if (line_c != NULL)
Vect_reset_cats(line_c);
- FID = ogr_info->offset.array[offset];
- G_debug(4, " FID = %ld", FID);
+ fid = ogr_info->offset.array[offset];
+ G_debug(4, " fid = %ld", fid);
/* coordinates */
if (line_p != NULL) {
- /* Read feature to cache if necessary */
- if (ogr_info->feature_cache_id != FID) {
- G_debug(4, "Read feature (FID = %ld) to cache", FID);
+ /* read feature to cache if necessary */
+ if (ogr_info->feature_cache_id != fid) {
+ G_debug(4, "Read feature (fid = %ld) to cache", fid);
if (ogr_info->feature_cache) {
OGR_F_Destroy(ogr_info->feature_cache);
}
ogr_info->feature_cache =
- OGR_L_GetFeature(ogr_info->layer, FID);
+ OGR_L_GetFeature(ogr_info->layer, fid);
if (ogr_info->feature_cache == NULL) {
- G_warning(_("Unable to get feature geometry, FID %ld"),
- FID);
+ G_warning(_("Unable to get feature geometry, fid %ld"),
+ fid);
return -1;
}
- ogr_info->feature_cache_id = FID;
+ ogr_info->feature_cache_id = fid;
}
hGeom = OGR_F_GetGeometryRef(ogr_info->feature_cache);
if (hGeom == NULL) {
- G_warning(_("Unable to get feature geometry, FID %ld"),
- FID);
+ G_warning(_("Unable to get feature geometry, fid %ld"),
+ fid);
return -1;
}
type = read_line(Map, hGeom, offset + 1, line_p);
}
else {
- type = get_line_type(Map, FID);
+ type = get_line_type(Map, fid);
}
/* category */
if (line_c != NULL) {
- Vect_cat_set(line_c, 1, (int) FID);
+ Vect_cat_set(line_c, 1, (int) fid);
}
return type;
@@ -275,6 +275,7 @@
}
if (Line->type == GV_CENTROID) {
+ /* read centroid for topo */
if (line_p != NULL) {
int i, found;
struct bound_box box;
@@ -288,7 +289,7 @@
/* get area bbox */
Vect_get_area_box(Map, topo->area, &box);
/* search in spatial index for centroid with area bbox */
- dig_init_boxlist(&list, 1);
+ dig_init_boxlist(&list, TRUE);
Vect_select_lines_by_box(Map, &box, Line->type, &list);
found = 0;
@@ -304,7 +305,7 @@
}
if (line_c != NULL) {
- /* cat = FID and offset = FID for centroid */
+ /* cat = fid and offset = fid for centroid */
Vect_reset_cats(line_c);
Vect_cat_set(line_c, 1, (int) Line->offset);
}
@@ -605,7 +606,7 @@
\return feature type
\return -1 on error
*/
-int get_line_type(const struct Map_info *Map, long FID)
+int get_line_type(const struct Map_info *Map, long fid)
{
int eType;
@@ -614,11 +615,11 @@
OGRFeatureH hFeat;
OGRGeometryH hGeom;
- G_debug(4, "get_line_type() fid = %ld", FID);
+ G_debug(4, "get_line_type() fid = %ld", fid);
ogr_info = &(Map->fInfo.ogr);
- hFeat = OGR_L_GetFeature(ogr_info->layer, FID);
+ hFeat = OGR_L_GetFeature(ogr_info->layer, fid);
if (hFeat == NULL)
return -1;
Modified: grass/trunk/lib/vector/Vlib/read_pg.c
===================================================================
--- grass/trunk/lib/vector/Vlib/read_pg.c 2012-03-07 22:03:17 UTC (rev 51018)
+++ grass/trunk/lib/vector/Vlib/read_pg.c 2012-03-07 22:12:18 UTC (rev 51019)
@@ -60,12 +60,17 @@
static unsigned char *hex_to_wkb(const char *, int *);
static int point_from_wkb(const unsigned char *, int, int, int,
struct line_pnts *);
-static int line_from_wkb(const unsigned char *, int, int, int,
- struct line_pnts *, int);
+static int linestring_from_wkb(const unsigned char *, int, int, int,
+ struct line_pnts *, int);
static int polygon_from_wkb(const unsigned char *, int, int, int,
- struct Format_info_pg *);
+ struct Format_info_cache *);
+static int geometry_collection_from_wkb(const unsigned char *, int, int, int,
+ struct Format_info_cache *,
+ struct feat_parts *);
static int error_corrupted_data(const char *);
static int set_initial_query();
+static void reallocate_cache(struct Format_info_cache *, int);
+static void add_fpart(struct feat_parts *, SF_FeatureType, int, int);
#endif
/*!
@@ -170,7 +175,7 @@
dig_init_boxlist(&list, TRUE);
Vect_select_lines_by_box(Map, &box, Line->type, &list);
- found = 0;
+ found = -1;
for (i = 0; i < list.n_values; i++) {
if (list.id[i] == line) {
found = i;
@@ -178,8 +183,10 @@
}
}
- Vect_reset_line(line_p);
- Vect_append_point(line_p, list.box[found].E, list.box[found].N, 0.0);
+ if (found > -1) {
+ Vect_reset_line(line_p);
+ Vect_append_point(line_p, list.box[found].E, list.box[found].N, 0.0);
+ }
}
if (line_c != NULL) {
/* cat = FID and offset = FID for centroid */
@@ -235,7 +242,7 @@
{
#ifdef HAVE_POSTGRES
long fid;
- int i, type;
+ int i, ipart, type;
SF_FeatureType sf_type;
struct line_pnts *line_i;
@@ -246,7 +253,7 @@
(long) offset, (long) pg_info->offset.array_num);
if (offset >= pg_info->offset.array_num)
- return -2;
+ return -2; /* nothing to read */
if (line_p != NULL)
Vect_reset_line(line_p);
@@ -267,9 +274,12 @@
return (int) sf_type;
}
+ ipart = pg_info->offset.array[offset + 1];
+ G_debug(4, "read feature part: %d", ipart);
+
/* get data from cache */
- type = pg_info->cache.lines_types[pg_info->cache.lines_next];
- line_i = pg_info->cache.lines[pg_info->cache.lines_next];
+ type = pg_info->cache.lines_types[ipart];
+ line_i = pg_info->cache.lines[ipart];
for (i = 0; i < line_i->n_points; i++) {
Vect_append_point(line_p,
line_i->x[i], line_i->y[i], line_i->z[i]);
@@ -318,6 +328,7 @@
}
if (Line->type == GV_CENTROID) {
+ /* read centroid from topo */
if (line_p != NULL) {
int i, found;
struct bound_box box;
@@ -331,10 +342,10 @@
/* get area bbox */
Vect_get_area_box(Map, topo->area, &box);
/* search in spatial index for centroid with area bbox */
- dig_init_boxlist(&list, 1);
+ dig_init_boxlist(&list, TRUE);
Vect_select_lines_by_box(Map, &box, Line->type, &list);
- found = 0;
+ found = -1;
for (i = 0; i < list.n_values; i++) {
if (list.id[i] == line) {
found = i;
@@ -342,8 +353,17 @@
}
}
- Vect_append_point(line_p, list.box[found].E, list.box[found].N, 0.0);
+ if (found > -1) {
+ Vect_append_point(line_p, list.box[found].E, list.box[found].N, 0.0);
+ }
+ else {
+ G_warning(_("Unable to construct centroid for area %d. Skipped."),
+ topo->area);
+ }
}
+ else {
+ G_warning(_("Centroid %d: invalid area %d"), line, topo->area);
+ }
}
if (line_c != NULL) {
@@ -411,8 +431,9 @@
if ((int) sf_type < 0) /* -1 || - 2 */
return (int) sf_type;
- if (sf_type & (SF_UNKNOWN | SF_NONE)) {
+ if (sf_type == SF_UNKNOWN || sf_type == SF_NONE) {
G_warning(_("Feature without geometry. Skipped."));
+ pg_info->cache.lines_next = pg_info->cache.lines_num = 0;
continue;
}
@@ -547,7 +568,7 @@
}
data = (char *)PQgetvalue(pg_info->res, pg_info->next_line, 0);
- ftype = cache_feature(data, FALSE, pg_info);
+ ftype = cache_feature(data, FALSE, &(pg_info->cache), NULL);
if (fid < 0) {
pg_info->cache.fid = atoi(PQgetvalue(pg_info->res, pg_info->next_line, 1));
pg_info->next_line++;
@@ -632,14 +653,15 @@
\param data HEX data
\param skip_polygon skip polygons (level 1)
- \param[in,out] pg_info pointer to Format_info_pg struct
- (geometry is stored in lines cache)
+ \param[out] cache lines cache
+ \param[out] fparts used for building pseudo-topology (or NULL)
- \return 0 on success
- \return -1 on error
+ \return simple feature type
+ \return SF_UNKNOWN on error
*/
-int cache_feature(const char *data, int skip_polygon,
- struct Format_info_pg *pg_info)
+SF_FeatureType cache_feature(const char *data, int skip_polygon,
+ struct Format_info_cache *cache,
+ struct feat_parts *fparts)
{
int ret, byte_order, nbytes, is3D;
unsigned char *wkb_data;
@@ -652,7 +674,7 @@
if (nbytes < 5) {
G_warning(_("Invalid WKB content: %d bytes"), nbytes);
G_free(wkb_data);
- return -1;
+ return SF_UNKNOWN;
}
/* parsing M coordinate not supported */
@@ -665,7 +687,7 @@
G_warning(_("Reading EWKB with 4-dimensional coordinates (XYZM) "
"is not supported"));
G_free(wkb_data);
- return -1;
+ return SF_UNKNOWN;
}
/* PostGIS EWKB format includes an SRID, but this won't be
@@ -685,7 +707,7 @@
if (nbytes < 9 && nbytes != -1) {
G_free(wkb_data);
- return -1;
+ return SF_UNKNOWN;
}
/* Get the geometry feature type. For now we assume that geometry
@@ -699,45 +721,52 @@
ftype = (SF_FeatureType) wkb_data[4];
is3D = wkb_data[1] & 0x80 || wkb_data[3] & 0x80;
}
-
+ G_debug(5, "cache_feature(): sf_type = %d", ftype);
+
/* allocate space in lines cache - be minimalistic
more lines require eg. polygon with more rings, multi-features
or geometry collections
*/
- if (!pg_info->cache.lines) {
- pg_info->cache.lines_alloc = 1;
- pg_info->cache.lines = (struct line_pnts **) G_malloc(sizeof(struct line_pnts *));
-
- pg_info->cache.lines_types = (int *) G_malloc(sizeof(int));
- pg_info->cache.lines[0] = Vect_new_line_struct();
- pg_info->cache.lines_types[0] = -1;
+ if (!cache->lines) {
+ reallocate_cache(cache, 1);
}
- pg_info->cache.lines_num = 0;
+ cache->lines_num = 0;
ret = -1;
if (ftype == SF_POINT) {
- pg_info->cache.lines_num = 1;
- pg_info->cache.lines_types[0] = GV_POINT;
+ cache->lines_num = 1;
+ cache->lines_types[0] = GV_POINT;
ret = point_from_wkb(wkb_data, nbytes, byte_order,
- is3D, pg_info->cache.lines[0]);
+ is3D, cache->lines[0]);
+ add_fpart(fparts, ftype, 0, 1);
}
else if (ftype == SF_LINESTRING) {
- pg_info->cache.lines_num = 1;
- pg_info->cache.lines_types[0] = GV_LINE;
- ret = line_from_wkb(wkb_data, nbytes, byte_order,
- is3D, pg_info->cache.lines[0], FALSE);
+ cache->lines_num = 1;
+ cache->lines_types[0] = GV_LINE;
+ ret = linestring_from_wkb(wkb_data, nbytes, byte_order,
+ is3D, cache->lines[0], FALSE);
+ add_fpart(fparts, ftype, 0, 1);
}
- else if (ftype == SF_POLYGON && !skip_polygon)
+ else if (ftype == SF_POLYGON && !skip_polygon) {
ret = polygon_from_wkb(wkb_data, nbytes, byte_order,
- is3D, pg_info) > 0 ? 0 : -1;
+ is3D, cache);
+ add_fpart(fparts, ftype, 0, 1);
+ }
+ else if (ftype == SF_MULTIPOINT ||
+ ftype == SF_MULTILINESTRING ||
+ ftype == SF_MULTIPOLYGON ||
+ ftype == SF_GEOMETRYCOLLECTION) {
+ ret = geometry_collection_from_wkb(wkb_data, nbytes, byte_order,
+ is3D, cache, fparts);
+ }
else {
- G_warning(_("Unsupported geometry type %d"), ftype);
+ G_warning(_("Unsupported feature type %d"), ftype);
}
G_free(wkb_data);
- return ret;
+ return ret > 0 ? ftype : SF_UNKNOWN;
}
/*!
@@ -751,7 +780,7 @@
\param with_z WITH_Z for 3D data
\param[out] line_p point geometry (pointer to line_pnts struct)
- \return 0 on success
+ \return wkb size
\return -1 on error
*/
int point_from_wkb(const unsigned char *wkb_data, int nbytes, int byte_order,
@@ -788,7 +817,7 @@
Vect_append_point(line_p, x, y, z);
}
- return 0;
+ return 5 + 8 * (with_z == WITH_Z ? 3 : 2);
}
/*!
@@ -802,11 +831,11 @@
\param with_z WITH_Z for 3D data
\param[out] line_p line geometry (pointer to line_pnts struct)
- \return 0 on success
+ \return wkb size
\return -1 on error
*/
-int line_from_wkb(const unsigned char *wkb_data, int nbytes, int byte_order,
- int with_z, struct line_pnts *line_p, int is_ring)
+int linestring_from_wkb(const unsigned char *wkb_data, int nbytes, int byte_order,
+ int with_z, struct line_pnts *line_p, int is_ring)
{
int npoints, point_size, buff_min_size, offset;
int i;
@@ -862,7 +891,7 @@
Vect_append_point(line_p, x, y, z);
}
- return 0;
+ return (9 - offset) + (with_z == WITH_Z ? 3 : 2) * 8 * line_p->n_points;
}
/*!
@@ -876,13 +905,13 @@
\param with_z WITH_Z for 3D data
\param[out] line_p array of rings (pointer to line_pnts struct)
- \return number of rings
+ \return wkb size
\return -1 on error
*/
int polygon_from_wkb(const unsigned char *wkb_data, int nbytes, int byte_order,
- int with_z, struct Format_info_pg *pg_info)
+ int with_z, struct Format_info_cache *cache)
{
- int nrings, data_offset, i, nsize;
+ int nrings, data_offset, i, nsize, isize;
struct line_pnts *line_i;
if (nbytes < 9 && nbytes != -1)
@@ -898,21 +927,8 @@
}
/* reallocate space for islands if needed */
- if (nrings > pg_info->cache.lines_alloc) {
- pg_info->cache.lines_alloc += 20;
- pg_info->cache.lines = (struct line_pnts **) G_realloc(pg_info->cache.lines,
- pg_info->cache.lines_alloc *
- sizeof(struct line_pnts *));
- pg_info->cache.lines_types = (int *) G_realloc(pg_info->cache.lines_types,
- pg_info->cache.lines_alloc *
- sizeof(int));
-
- for (i = pg_info->cache.lines_alloc - 20; i < pg_info->cache.lines_alloc; i++) {
- pg_info->cache.lines[i] = Vect_new_line_struct();
- pg_info->cache.lines_types[i] = -1;
- }
- }
- pg_info->cache.lines_num = nrings;
+ reallocate_cache(cache, nrings);
+ cache->lines_num += nrings;
/* each ring has a minimum of 4 bytes (point count) */
if (nbytes != -1 && nbytes - 9 < nrings * 4) {
@@ -924,27 +940,131 @@
nbytes -= data_offset;
/* get the rings */
- nsize = 0;
+ nsize = 9;
for (i = 0; i < nrings; i++ ) {
- line_i = pg_info->cache.lines[i];
- pg_info->cache.lines_types[i] = GV_BOUNDARY;
+ line_i = cache->lines[cache->lines_next];
+ cache->lines_types[cache->lines_next++] = GV_BOUNDARY;
- line_from_wkb(wkb_data + data_offset, nbytes, byte_order,
- with_z, line_i, TRUE);
+ linestring_from_wkb(wkb_data + data_offset, nbytes, byte_order,
+ with_z, line_i, TRUE);
if (nbytes != -1) {
- if (with_z)
- nsize = 4 + 24 * line_i->n_points;
- else
- nsize = 4 + 16 * line_i->n_points;
-
+ isize = 4 + 8 * (with_z == WITH_Z ? 3 : 2) * line_i->n_points;
+ nbytes -= isize;
+ }
+
+ nsize += isize;
+ data_offset += isize;
+ }
+
+ return nsize;
+}
+
+/*!
+ \brief Read geometry collection for WKB data
+
+ See OGRGeometryCollection::importFromWkbInternal() from GDAL/OGR library
+
+ \param wkb_data WKB data
+ \param nbytes number of bytes (WKB data buffer)
+ \param byte_order byte order (ENDIAN_LITTLE, ENDIAN_BIG)
+ \param with_z WITH_Z for 3D data
+ \param ipart part to cache (starts at 0)
+ \param[out] cache lines cache
+ \param[in,out] fparts feature parts (required for building pseudo-topology)
+
+ \return number of parts
+ \return -1 on error
+*/
+int geometry_collection_from_wkb(const unsigned char *wkb_data, int nbytes, int byte_order,
+ int with_z, struct Format_info_cache *cache,
+ struct feat_parts *fparts)
+{
+ int ipart, nparts, data_offset, nsize;
+ unsigned char *wkb_subdata;
+ SF_FeatureType ftype;
+
+ if (nbytes < 9 && nbytes != -1)
+ return error_corrupted_data(NULL);
+
+ /* get the geometry count */
+ memcpy(&nparts, wkb_data + 5, 4 );
+ if (byte_order == ENDIAN_BIG) {
+ nparts = SWAP32(nparts);
+ }
+ if (nparts < 0 || nparts > INT_MAX / 9) {
+ return error_corrupted_data(NULL);
+ }
+ G_debug(5, "\t(geometry collections) parts: %d", nparts);
+
+ /* each geometry has a minimum of 9 bytes */
+ if (nbytes != -1 && nbytes - 9 < nparts * 9) {
+ return error_corrupted_data(_("Length of input WKB is too small"));
+ }
+
+ data_offset = 9;
+ if (nbytes != -1)
+ nbytes -= data_offset;
+
+ /* reallocate space for parts if needed */
+ reallocate_cache(cache, nparts);
+
+ /* get parts */
+ cache->lines_next = cache->lines_num = 0;
+ for (ipart = 0; ipart < nparts; ipart++) {
+ wkb_subdata = (unsigned char *)wkb_data + data_offset;
+ if (nbytes < 9 && nbytes != -1)
+ return error_corrupted_data(NULL);
+
+ if (byte_order == ENDIAN_LITTLE) {
+ ftype = (SF_FeatureType) wkb_subdata[1];
+ }
+ else {
+ ftype = (SF_FeatureType) wkb_subdata[4];
+ }
+
+ if (ftype == SF_POINT) {
+ cache->lines_types[cache->lines_next] = GV_POINT;
+ nsize = point_from_wkb(wkb_subdata, nbytes, byte_order, with_z,
+ cache->lines[cache->lines_next]);
+ cache->lines_num++;
+ add_fpart(fparts, ftype, cache->lines_next, 1);
+ cache->lines_next++;
+ }
+ else if (ftype == SF_LINESTRING) {
+ cache->lines_types[cache->lines_next] = GV_LINE;
+ nsize = linestring_from_wkb(wkb_subdata, nbytes, byte_order, with_z,
+ cache->lines[cache->lines_next],
+ FALSE);
+ cache->lines_num++;
+ add_fpart(fparts, ftype, cache->lines_next, 1);
+ cache->lines_next++;
+ }
+ else if (ftype == SF_POLYGON) {
+ int idx = cache->lines_next;
+ nsize = polygon_from_wkb(wkb_subdata, nbytes, byte_order,
+ with_z, cache);
+ add_fpart(fparts, ftype, idx, cache->lines_num - idx);
+ }
+ else if (ftype == SF_GEOMETRYCOLLECTION ||
+ ftype == SF_MULTIPOLYGON ||
+ ftype == SF_MULTILINESTRING ||
+ ftype == SF_MULTIPOLYGON) {
+ // geometry_collection_from_wkb();
+ }
+ else {
+ G_warning(_("Unsupported feature type %d"), ftype);
+ }
+
+ if (nbytes != -1) {
nbytes -= nsize;
}
data_offset += nsize;
}
+ cache->lines_next = 0;
- return nrings;
+ return nparts;
}
/*!
@@ -1025,4 +1145,67 @@
return 0;
}
+/*!
+ \brief Reallocate lines cache
+*/
+void reallocate_cache(struct Format_info_cache *cache, int num)
+{
+ int i;
+
+ if (cache->lines_alloc >= num)
+ return;
+
+ if (!cache->lines) {
+ /* most of features requires only one line cache */
+ cache->lines_alloc = 1;
+ }
+ else {
+ cache->lines_alloc += 20;
+ }
+
+ cache->lines = (struct line_pnts **) G_realloc(cache->lines,
+ cache->lines_alloc *
+ sizeof(struct line_pnts *));
+ cache->lines_types = (int *) G_realloc(cache->lines_types,
+ cache->lines_alloc *
+ sizeof(int));
+
+ if (cache->lines_alloc > 1) {
+ for (i = cache->lines_alloc - 20; i < cache->lines_alloc; i++) {
+ cache->lines[i] = Vect_new_line_struct();
+ cache->lines_types[i] = -1;
+ }
+ }
+ else {
+ cache->lines[0] = Vect_new_line_struct();
+ cache->lines_types[0] = -1;
+ }
+}
+
+void add_fpart(struct feat_parts *fparts, SF_FeatureType ftype,
+ int idx, int nlines)
+{
+ if (!fparts)
+ return;
+
+ if (fparts->a_parts == 0 || fparts->n_parts >= fparts->a_parts) {
+ if (fparts->a_parts == 0)
+ fparts->a_parts = 1;
+ else
+ fparts->a_parts += 20;
+
+ fparts->ftype = (SF_FeatureType *) G_realloc(fparts->ftype,
+ fparts->a_parts * sizeof(SF_FeatureType));
+ fparts->nlines = (int *) G_realloc(fparts->nlines,
+ fparts->a_parts * sizeof(int));
+ fparts->idx = (int *) G_realloc(fparts->idx,
+ fparts->a_parts * sizeof(int));
+ }
+
+ fparts->ftype[fparts->n_parts] = ftype;
+ fparts->idx[fparts->n_parts] = idx;
+ fparts->nlines[fparts->n_parts] = nlines;
+
+ fparts->n_parts++;
+}
#endif
More information about the grass-commit
mailing list