From: Sandro Santilli <strk@keybit.net>
Date: Wed, 8 Sep 2004 14:19:25 +0000 (+0000)
Subject: Added segmentize()
X-Git-Tag: pgis_0_9_1~9
X-Git-Url: https://granicus.if.org/sourcecode?a=commitdiff_plain;h=f57cfd6a08e16cd228766a75eeae255b5f6407e1;p=postgis

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 4a748c558..47011e7d9 100644
--- 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 58ceac819..a3c217163 100644
--- 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 430ddc849..c22c18a39 100644
--- 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 1244ffecc..162e3a37d 100644
--- 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 2313a3fe5..182e8c993 100644
--- 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
 ---------------------------------------------------------------