Added segmentize()

git-svn-id: http://svn.osgeo.org/postgis/trunk@785 b70326c6-7e19-0410-871a-916f4a2858ee

diff --git a/lwgeom/MISSING_OBJECTS b/lwgeom/MISSING_OBJECTS
index 4a748c558861695f2bf132629bda4cdd750395e2..47011e7d968e68e25d471faa86d5b524056be627 100644 (file)
--- a/lwgeom/MISSING_OBJECTS
+++ b/lwgeom/MISSING_OBJECTS
@@ -7,7 +7,6 @@ FUNC: KEEPING FUNCTION: [asbinary(geometry, text)]
 FUNC: KEEPING FUNCTION: [line_interpolate_point(geometry, double precision)]
 
 --- UNDOCUMENTED
-FUNC: KEEPING FUNCTION: [segmentize(geometry, double precision)]
 FUNC: KEEPING FUNCTION: [optimistic_overlap(geometry, geometry, double precision)]
 
 --- OBSOLETED ?

diff --git a/lwgeom/lwgeom.h b/lwgeom/lwgeom.h
index 58ceac81964eee62aa2f0e4069b8ca33dc78a56c..a3c217163959fb44cf0043a5cfeabe9b52119d1d 100644 (file)
--- a/lwgeom/lwgeom.h
+++ b/lwgeom/lwgeom.h
@@ -769,3 +769,4 @@ double lwgeom_mindistance2d_recursive(char *lw1, char *lw2);
 void lwgeom_translate_recursive(char *serialized, double xoff, double yoff, double zoff);
 void lwgeom_translate_ptarray(POINTARRAY *pa, double xoff, double yoff, double zoff);
 int lwgeom_pt_inside_circle(POINT2D *p, double cx, double cy, double rad);
+POINTARRAY *segmentize2d_ptarray(POINTARRAY *ipa, double dist);

diff --git a/lwgeom/lwgeom_functions_analytic.c b/lwgeom/lwgeom_functions_analytic.c
index 430ddc84926a5be4d02d09812e63a515819c92b9..c22c18a39b0042fccd99182f5288faa17ca7d33d 100644 (file)
--- a/lwgeom/lwgeom_functions_analytic.c
+++ b/lwgeom/lwgeom_functions_analytic.c
@@ -330,3 +330,4 @@ Datum LWGEOM_simplify2d(PG_FUNCTION_ARGS)
 /***********************************************************************
  * --strk@keybit.net;
  ***********************************************************************/
+

diff --git a/lwgeom/lwgeom_functions_basic.c b/lwgeom/lwgeom_functions_basic.c
index 1244ffecc03ce9285958a0d918c35e3e11ef14e0..162e3a37d1ac49e221f52db19041ef1e80de4ed1 100644 (file)
--- a/lwgeom/lwgeom_functions_basic.c
+++ b/lwgeom/lwgeom_functions_basic.c
@@ -41,6 +41,7 @@ Datum LWGEOM_expand(PG_FUNCTION_ARGS);
 Datum LWGEOM_to_BOX(PG_FUNCTION_ARGS);
 Datum LWGEOM_envelope(PG_FUNCTION_ARGS);
 Datum LWGEOM_isempty(PG_FUNCTION_ARGS);
+Datum LWGEOM_segmentize2d(PG_FUNCTION_ARGS);
 
 // internal
 char * lwgeom_summary_recursive(char *serialized, int offset);
@@ -2401,3 +2402,115 @@ Datum centroid(PG_FUNCTION_ARGS)
        PG_RETURN_POINTER(result);
 }
 #endif // ! USE_GEOS
+
+// Returns a modified POINTARRAY so that no segment is 
+// longer then the given distance (computed using 2d).
+// Every input point is kept.
+// Z and M values for added points (if needed) are set to 0.
+POINTARRAY *
+segmentize2d_ptarray(POINTARRAY *ipa, double dist)
+{
+       double  segdist;
+       POINT4D *p1, *p2, *ip, *op;
+       POINT4D pbuf;
+       POINTARRAY *opa;
+       int maxpoints = ipa->npoints;
+       int ptsize = ipa->ndims * sizeof(double);
+       int ipoff=0; // input point offset
+
+       pbuf.x = pbuf.y = pbuf.z = pbuf.m = 0;
+
+       // Initial storage
+       opa = (POINTARRAY *)palloc(ptsize * maxpoints);
+       opa->ndims = ipa->ndims;
+       opa->npoints = 0;
+       opa->serialized_pointlist = (char *)palloc(maxpoints*ptsize);
+
+       // Add first point
+       opa->npoints++;
+       p1 = (POINT4D *)getPoint(ipa, ipoff);
+       op = (POINT4D *)getPoint(opa, opa->npoints-1);
+       memcpy(op, p1, ptsize); 
+       ipoff++;
+
+       while (ipoff<ipa->npoints)
+       {
+               p2 = (POINT4D *)getPoint(ipa, ipoff);
+
+               segdist = distance2d_pt_pt((POINT2D *)p1, (POINT2D *)p2);
+
+               if (segdist > dist) // add an intermediate point
+               {
+                       pbuf.x = p1->x + (p2->x-p1->x)/segdist * dist;
+                       pbuf.y = p1->y + (p2->y-p1->y)/segdist * dist;
+                       // might also compute z and m if available...
+                       ip = &pbuf;
+                       p1 = ip;
+               }
+               else // copy second point
+               {
+                       ip = p2;
+                       p1 = p2;
+                       ipoff++;
+               }
+
+               // Add point
+               if ( ++(opa->npoints) > maxpoints ) {
+                       maxpoints *= 1.5;
+                       opa->serialized_pointlist = (char *)repalloc(
+                               opa->serialized_pointlist,
+                               maxpoints*ptsize
+                       );
+               }
+               op = (POINT4D *)getPoint(opa, opa->npoints-1);
+               memcpy(op, ip, ptsize); 
+       }
+
+       return opa;
+}
+
+// Returns a modified [multi]polygon so that no ring segment is 
+// longer then the given distance (computed using 2d).
+// Every input point is kept.
+// Z and M values for added points (if needed) are set to 0.
+PG_FUNCTION_INFO_V1(LWGEOM_segmentize2d);
+Datum LWGEOM_segmentize2d(PG_FUNCTION_ARGS)
+{
+       LWGEOM *geom = (LWGEOM *)PG_DETOAST_DATUM(PG_GETARG_DATUM(0));
+       double dist = PG_GETARG_FLOAT8(1);
+       LWGEOM *result = NULL;
+       int type = lwgeom_getType(geom->type);
+       LWGEOM_EXPLODED *exp;
+       int i, j;
+       char *srl;
+
+       if ( (type != POLYGONTYPE) && (type != MULTIPOLYGONTYPE) )
+       {
+               elog(ERROR,"segmentize: 1st arg isnt a [multi-]polygon\n");
+               PG_RETURN_NULL();
+       }
+
+       exp = lwgeom_explode(SERIALIZED_FORM(geom));
+
+       for (i=0; i<exp->npolys; i++)
+       {
+               LWPOLY *poly = lwpoly_deserialize(exp->polys[i]);
+               for (j=0; j<poly->nrings; j++) 
+               {
+                       POINTARRAY *ring = poly->rings[j];
+                       poly->rings[j] = segmentize2d_ptarray(ring, dist);
+               } 
+               pfree_polygon(poly);
+               exp->polys[i] = lwpoly_serialize(poly);
+       } 
+
+       // Serialize exploded structure
+       srl = lwexploded_serialize(exp, lwgeom_hasBBOX(geom->type));
+       pfree_exploded(exp);
+
+       // Construct return value
+       result = LWGEOM_construct(srl, lwgeom_getSRID(geom),
+               lwgeom_hasBBOX(geom->type));
+
+       PG_RETURN_POINTER(result);
+}

diff --git a/lwgeom/lwpostgis.sql.in b/lwgeom/lwpostgis.sql.in
index 2313a3fe591df90fb5b8d4f06499d32f1465aa52..182e8c9934384e608ba413a5c047f33bc31f41f4 100644 (file)
--- a/lwgeom/lwpostgis.sql.in
+++ b/lwgeom/lwpostgis.sql.in
@@ -2503,6 +2503,11 @@ CREATEFUNCTION simplify(geometry, float8)
    AS '@MODULE_FILENAME@', 'LWGEOM_simplify2d'
    LANGUAGE 'C' WITH (isstrict);
 
+CREATEFUNCTION segmentize(geometry, float8)
+       RETURNS geometry
+       AS '@MODULE_FILENAME@', 'LWGEOM_segmentize2d'
+       LANGUAGE 'C' WITH (isstrict,iscachable);
+
 ---------------------------------------------------------------
 -- GEOS
 ---------------------------------------------------------------