Added polygonize interface (makepoly aggregate)

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

diff --git a/lwgeom/lwgeom_geos.c b/lwgeom/lwgeom_geos.c
index c6d656147ee2a69c898777983e728f18d5374e95..d5e4152e91078056fd5a7c977bb08ec284a1f0eb 100644 (file)
--- a/lwgeom/lwgeom_geos.c
+++ b/lwgeom/lwgeom_geos.c
@@ -34,6 +34,7 @@ Datum pointonsurface(PG_FUNCTION_ARGS);
 Datum GEOSnoop(PG_FUNCTION_ARGS);
 Datum postgis_geos_version(PG_FUNCTION_ARGS);
 Datum centroid(PG_FUNCTION_ARGS);
+Datum GEOS_makepoly_garray(PG_FUNCTION_ARGS);
 
 
 
@@ -48,6 +49,7 @@ Datum centroid(PG_FUNCTION_ARGS);
 #define DEBUG_POSTGIS2GEOS 1
 #define DEBUG_GEOS2POSTGIS 1
 #endif // DEBUG_CONVERTER
+#define DEBUG 1
 
 typedef  struct Geometry Geometry;
 
@@ -91,6 +93,7 @@ extern int      GEOSGetNumCoordinate(Geometry *g1);
 extern Geometry        *GEOSGetGeometryN(Geometry *g1, int n);
 extern Geometry        *GEOSGetExteriorRing(Geometry *g1);
 extern Geometry        *GEOSGetInteriorRingN(Geometry *g1,int n);
+extern Geometry *GEOSpolygonize(Geometry **geoms, unsigned int ngeoms);
 extern int     GEOSGetNumInteriorRings(Geometry *g1);
 extern int      GEOSGetSRID(Geometry *g1);
 extern int      GEOSGetNumGeometries(Geometry *g1);
@@ -2239,6 +2242,76 @@ Datum GEOSnoop(PG_FUNCTION_ARGS)
        PG_RETURN_POINTER(result);
 }
 
+PG_FUNCTION_INFO_V1(GEOS_makepoly_garray);
+Datum GEOS_makepoly_garray(PG_FUNCTION_ARGS)
+{
+       Datum datum;
+       ArrayType *array;
+       int is3d = 0;
+       unsigned int nelems, i;
+       PG_LWGEOM **geoms, *result, *pgis_geom;
+       Geometry *geos_result;
+       Geometry **vgeoms;
+#ifdef DEBUG
+       static int call=1;
+#endif
+
+#ifdef DEBUG
+       call++;
+#endif
+
+       datum = PG_GETARG_DATUM(0);
+
+       /* Null array, null geometry (should be empty?) */
+       if ( (Pointer *)datum == NULL ) PG_RETURN_NULL();
+
+       array = (ArrayType *) PG_DETOAST_DATUM(datum);
+
+       nelems = ArrayGetNItems(ARR_NDIM(array), ARR_DIMS(array));
+
+#ifdef DEBUG
+       elog(NOTICE, "GEOS_makepoly_garray: number of elements: %d", nelems);
+#endif
+
+       if ( nelems == 0 ) PG_RETURN_NULL();
+
+       geoms = (PG_LWGEOM **)ARR_DATA_PTR(array);
+
+       /* Ok, we really need geos now ;) */
+       initGEOS(MAXIMUM_ALIGNOF);
+
+       vgeoms = palloc(sizeof(Geometry *)*nelems);
+       for (i=0; i<nelems; i++)
+       {
+               //if ( TYPE_NDIMS(geoms[i]->type) > 2 ) is3d = 1;
+               vgeoms[i] = POSTGIS2GEOS(geoms[i]);
+       }
+
+#ifdef DEBUG
+       elog(NOTICE, "GEOS_makepoly_garray: invoking GEOSpolygonize");
+#endif
+
+       geos_result = GEOSpolygonize(vgeoms, nelems);
+#ifdef DEBUG
+       elog(NOTICE, "GEOS_makepoly_garray: GEOSpolygonize returned");
+#endif
+       //pfree(vgeoms);
+       if ( ! geos_result ) PG_RETURN_NULL();
+
+       result = GEOS2POSTGIS(geos_result, is3d);
+       GEOSdeleteGeometry(geos_result);
+       if ( result == NULL )
+       {
+               elog(ERROR, "GEOS2POSTGIS returned an error");
+               PG_RETURN_NULL(); //never get here
+       }
+
+       //compressType(result);
+
+       PG_RETURN_POINTER(result);
+
+}
+
 #else // ! USE_GEOS
 
 PG_FUNCTION_INFO_V1(postgis_geos_version);
@@ -2388,4 +2461,11 @@ Datum unite_garray(PG_FUNCTION_ARGS)
        PG_RETURN_NULL(); // never get here
 }
 
+PG_FUNCTION_INFO_V1(GEOS_makepoly_garray);
+Datum GEOS_makepoly_garray(PG_FUNCTION_ARGS)
+{
+       elog(ERROR,"GEOS_makepoly_garray:: operation not implemented - compile PostGIS with GEOS support");
+       PG_RETURN_NULL(); // never get here
+}
+
 #endif // ! USE_GEOS

diff --git a/lwgeom/lwgeom_geos_wrapper.cpp b/lwgeom/lwgeom_geos_wrapper.cpp
index 950e6efea27a6d8c142f7f6ae95f06e01cc7718d..7a7185761eb9029c8b40362f40d8261caa80e5eb 100644 (file)
--- a/lwgeom/lwgeom_geos_wrapper.cpp
+++ b/lwgeom/lwgeom_geos_wrapper.cpp
@@ -9,8 +9,12 @@
 #include "postgis_geos_version.h"
 #include "geos/geom.h"
 #include "geos/util.h"
+#if GEOS_FIRST_INTERFACE <= 3 
+#include "geos/opPolygonize.h"
+#endif // GEOS_FIRST_INTERFACE <= 3
 
 //#define DEBUG_POSTGIS2GEOS 1
+#define DEBUG 1
 
 using namespace geos;
 
@@ -102,6 +106,7 @@ extern "C" char GEOSrelateCrosses(Geometry *g1, Geometry*g2);
 extern "C" char GEOSrelateWithin(Geometry *g1, Geometry*g2);
 extern "C" char GEOSrelateContains(Geometry *g1, Geometry*g2);
 extern "C" char GEOSrelateOverlaps(Geometry *g1, Geometry*g2);
+extern "C" Geometry *GEOSpolygonize(Geometry **geoms, unsigned int ngeoms);
 extern "C" char *GEOSversion();
 extern "C" char *GEOSjtsport();
 
@@ -1578,6 +1583,65 @@ Geometry *GEOSGetCentroid(Geometry *g)
        }
 }
 
+#if GEOS_FIRST_INTERFACE <= 3 && GEOS_LAST_INTERFACE >= 3
+Geometry *GEOSpolygonize(Geometry **g, unsigned int ngeoms)
+{
+       unsigned int i;
+       Geometry *multipoly = NULL;
+
+       // construct vector
+       vector<Geometry *> *geoms = new vector<Geometry *>(ngeoms);
+       for (i=0; i<ngeoms; i++) (*geoms)[i] = g[i];
+
+#if DEBUG
+       NOTICE_MESSAGE("geometry vector constructed");
+#endif
+
+       try{
+               // Polygonize
+               Polygonizer plgnzr;
+               plgnzr.add(geoms);
+#if DEBUG
+       NOTICE_MESSAGE("geometry vector added to polygonizer");
+#endif
+
+               vector<Polygon *>*polys = plgnzr.getPolygons();
+
+#if DEBUG
+       NOTICE_MESSAGE("output polygons got");
+#endif
+
+               delete geoms;
+
+#if DEBUG
+       NOTICE_MESSAGE("geometry vector deleted");
+#endif
+
+               geoms = new vector<Geometry *>(polys->size());
+               for (i=0; i<polys->size(); i++) (*geoms)[i] = (*polys)[i];
+               multipoly = geomFactory->createMultiPolygon(geoms);
+       }
+       catch (GEOSException *ge)
+       {
+               NOTICE_MESSAGE((char *)ge->toString().c_str());
+               delete ge;
+               return NULL;
+       }
+       catch(...)
+       {
+               return NULL;
+       }
+
+       return multipoly;
+}
+#else // ! (GEOS_FIRST_INTERFACE <= 3 && GEOS_LAST_INTERFACE >= 3)
+Geometry *GEOSpolygonize(Geometry **g, unsigned int ngeoms)
+{
+       NOTICE_MESSAGE("GEOS library does not support required interface 3");
+       return NULL;
+}
+#endif // ! (GEOS_FIRST_INTERFACE <= 3 && GEOS_LAST_INTERFACE >= 3)
+
 
 int      GEOSGetSRID(Geometry *g1)
 {
@@ -1630,3 +1694,4 @@ GEOSjtsport()
 }
 
 
+