Changeset: 24420eca351a for MonetDB
URL: http://dev.monetdb.org/hg/MonetDB?cmd=changeset;node=24420eca351a
Added Files:
geom/monetdb5/geom_clean.c
Branch: sfcgal
Log Message:
First steps for MakeValid
diffs (truncated from 737 to 300 lines):
diff --git a/geom/monetdb5/geom_clean.c b/geom/monetdb5/geom_clean.c
new file mode 100644
--- /dev/null
+++ b/geom/monetdb5/geom_clean.c
@@ -0,0 +1,732 @@
+#include "geom.h"
+
+static str
+geos_verify(const GEOSGeometry **res, const GEOSGeometry *geom)
+{
+ int geometryType = GEOSGeomTypeId(geosGeometry) + 1;
+ str msg = MAL_SUCCEED;
+
+ switch (geom->type)
+ {
+ case POINTTYPE:
+ case MULTIPOINTTYPE:
+ *res = geom;
+ break;
+ case LINETYPE:
+ msg = geos_line_verify(res, geom);
+ break;
+ case POLYGONTYPE:
+ msg = geos_polygon_verify(res, geom);
+ break;
+ case MULTILINETYPE:
+ case MULTIPOLYGONTYPE:
+ case COLLECTIONTYPE:
+ msg = geos_collection_verify(res, geom);
+ break;
+ case CIRCSTRINGTYPE:
+ case COMPOUNDTYPE:
+ case CURVEPOLYTYPE:
+ case MULTISURFACETYPE:
+ case MULTICURVETYPE:
+ default:
+ *res = NULL;
+ msg = createException(MAL, "geos_verify", "Unknown geometry type");
+ break;
+ }
+
+ return msg;
+}
+
+static str
+geos_ring_verify(const GEOSGeometry **res, const GEOSGeometry *geosGeometry)
+{
+ int closed = 0;
+ str err = MAL_SUCCEED;
+ unsigned int pointsNum;
+ int srid = 0;
+
+ /*Check if the ring is closed*/
+ if ((closed = GEOSisClosed(geosGeometry)) == 2)
+ throw(MAL, "geos_ring_verify", "GEOSisClosed failed");
+
+ /*If it is not closed, close it*/
+ if (closed != 1) {
+ const GEOSCoordSequence *gcs_new = NULL, *gcs_old;
+ srid = GEOSGetSRID(geosGeometry);
+ double x, y, z;
+ int lineDim = 0;
+
+ if ((err = numPointsLineString(&pointsNum, geosGeometry)) !=
MAL_SUCCEED) {
+ *res = NULL;
+ msg = createException(MAL, "geos_ring_verify",
"numPointsLineString failed:%s", err;
+ GDKfree(err);
+ return msg;
+ } else {
+ assert(pointsNum < 4);
+ //get the coordinate sequences of the LineString
+ if (!(gcs_old = GEOSGeom_getCoordSeq(line))) {
+ *res = NULL;
+ throw(MAL, "geom.AddPoint", "GEOSGeom_getCoordSeq failed");
+ }
+
+ //Get dimension of the LineString
+ if (GEOSCoordSeq_getDimensions(gcs_old, &lineDim) == 0) {
+ *res = NULL;
+ throw(MAL, "geom.AddPoint", "GEOSGeom_getDimensions failed");
+ }
+
+ if ((gcs_new = GEOSCoordSeq_create(4, lineDim)) == NULL) {
+ *res = NULL;
+ throw(MAL, "geom.AddPoint", "GEOSCoordSeq_create failed");
+ }
+
+ if (!GEOSCoordSeq_getX(gcs_old, 0, &x))
+ throw(MAL, "geom.AddPoint", "GEOSCoordSeq_getX failed");
+
+ if (!GEOSCoordSeq_getY(gcs_old, 0, &y))
+ throw(MAL, "geom.AddPoint", "GEOSCoordSeq_getY failed");
+
+
+ if (lineDim > 2) {
+ if (!GEOSCoordSeq_getZ(gcs_old, 0, &z))
+ throw(MAL, "geom.AddPoint", "GEOSCoordSeq_getZ failed");
+ }
+
+ for (j = pointsNum-1; j < 4; j++) {
+ if (!GEOSCoordSeq_setX(gcs_new, j, x))
+ throw(MAL, "geom.AddPoint", "GEOSCoordSeq_setX failed");
+
+ if (!GEOSCoordSeq_setY(gcs_new, j, y))
+ throw(MAL, "geom.AddPoint", "GEOSCoordSeq_setY failed");
+
+ if (lineDim > 2) {
+ if (!GEOSCoordSeq_setZ(gcs_new, i, z))
+ throw(MAL, "geom.AddPoint", "GEOSCoordSeq_setZ
failed");
+ }
+ }
+ //Create new LineString
+ if (!(res = GEOSGeom_createLineString(gcs_new))) {
+ throw(MAL, "geom.AddPoint", "GEOSGeom_createLineString
failed");
+ }
+
+ GEOSSetSRID(res, srid);
+ }
+ } else {
+ /*We might have to clone it*/
+ *res = geosGeometry;
+ }
+
+ return MAL_SUCCEED;
+}
+
+str
+geos_poly_verify(const GEOSGeometry **res, const GEOSGeometry *geosGeometry)
+{
+ const GEOSGeometry *extRing = NULL, *extRres = NULL, **intRings = NULL;
+ int numInteriorRings = 0, i = 0, j = 0;
+ POINTARRAY **new_rings;
+ int i = 0, numIntRings = 0;
+ bit untouched = 1;
+
+ /* get the exterior ring of the polygon */
+ extRing = GEOSGetExteriorRing(geosGeometry);
+ if (extRing == NULL) {
+ *res = NULL;
+ throw(MAL, "geos_poly_verify", "GEOSGetExteriorRing failed");
+ }
+
+ /*verify exterior ring*/
+ if ( (geos_ring_verify(&extRes, extRing)) != MAL_SUCCEED) {
+ } else if (extRes != extRing)
+ untouched = 0;
+
+ numIntRings = GEOSGetNumInteriorRings(geosGeometry);
+ if (numInteriorRings == -1) {
+ *res = NULL;
+ throw(MAL, "geos_poly_verify", "GEOSGetNumInteriorRings
failed.");
+ } else if(numIntRings) {
+
+ if ( (intRings = (const GEOSGeometry **) GDKzalloc(sizeof(const
GEOSGeometry *)*numInteriorRings)) == NULL) {
+ *res = NULL;
+ throw(MAL, "geos_poly_verify", MAL_MALLOC_FAIL);
+ }
+
+ for (i = 0; i < numInteriorRings; i++) {
+ const GEOSGeometry *intRing = NULL;
+ if ((intRing = GEOSGetInteriorRingN(geosGeometry, i)) == NULL) {
+ msg = createException(MAL, "geos_poly_verify",
"GEOSGetInteriorRingN failed.");
+ break;
+ }
+ if ( (msg = geos_ring_verify(intRings[i], intRing)) !=
MAL_SUCCEED) {
+ break;
+ } else if (intRings[i] != intRing)
+ untouched = 0;
+ }
+ if (msg != MAL_SUCCEED) {
+ *res = NULL;
+ /*TODO: You should destroy if you are cloning*/
+ for (j = 0; j < i; j++)
+ GEOSGeom_destroy(intRings[j]);
+ GDKfree(intRings);
+ return msg;
+ }
+ }
+
+ /*Create a new geometry*/
+ if (!untouched) {
+ if ( (*res = GEOSGeom_createPolygon(exteriorRingGeometry, intRings,
numInteriorRings)) == NULL) {
+ GEOSGeom_destroy(exteriorRingGeometry);
+ if (numInteriorRings) {
+ /*TODO: You should destroy if you are cloning*/
+ for (i = 0; i < numInteriorRings; i++)
+ GEOSGeom_destroy(intRings[i]);
+ GDKfree(intRings);
+ }
+ return createException(MAL, "geos_poly_verify",
"GEOSGeom_createPolygon failed");
+ }
+ } else {
+ /*TODO: Maybe we should clone it*/
+ *res = geosGeometry;
+ }
+
+ if (intRings)
+ GDKfree(intRings);
+
+ return MAL_SUCCEED;
+}
+
+/*Line has 0 or more than one point. If missing duplicate points*/
+str
+geos_line_verify(const GEOSGeometry **res, const GEOSGeometry *geosGeometry)
+{
+ unsigned int *pointsN;
+ GEOSGeom *ret;
+ str err = MAL_SUCCEED;
+ int srid = 0;
+
+
+ //get the points in the exterior ring
+ if ((err = numPointsLineString(pointsN, lineGeometry)) != MAL_SUCCEED) {
+ msg = createException(MAL,"geos_line_verify", "numPointsLineString
failed:%s", err);
+ GDKfree(err);
+ *res = NULL;
+ return err;
+ }
+
+ if (pointsN == 1) {
+ const GEOSCoordSequence *gcs_new = NULL, *gcs_old;
+ srid = GEOSGetSRID(geosGeometry);
+ double x, y, z;
+ int lineDim = 0;
+
+ if (!(gcs_old = GEOSGeom_getCoordSeq(line))) {
+ *res = NULL;
+ throw(MAL, "geos_line_verify", "GEOSGeom_getCoordSeq failed");
+ }
+
+ //Get dimension of the LineString
+ if (GEOSCoordSeq_getDimensions(gcs_old, &lineDim) == 0) {
+ *res = NULL;
+ throw(MAL, "geos_line_verify", "GEOSGeom_getDimensions failed");
+ }
+
+ if ((gcs_new = GEOSCoordSeq_create(2, lineDim)) == NULL) {
+ *res = NULL;
+ throw(MAL, "geos_line_verify", "GEOSCoordSeq_create failed");
+ }
+
+ if (!GEOSCoordSeq_getX(gcs_old, 0, &x))
+ throw(MAL, "geos_line_verify", "GEOSCoordSeq_getX failed");
+
+ if (!GEOSCoordSeq_getY(gcs_old, 0, &y))
+ throw(MAL, "geos_line_verify", "GEOSCoordSeq_getY failed");
+
+ if (lineDim > 2) {
+ if (!GEOSCoordSeq_getZ(gcs_old, 0, &z))
+ throw(MAL, "geos_line_verify", "GEOSCoordSeq_getZ failed");
+ }
+
+ if (!GEOSCoordSeq_setX(gcs_new, 1, x))
+ throw(MAL, "geos_line_verify", "GEOSCoordSeq_setX failed");
+
+ if (!GEOSCoordSeq_setY(gcs_new, 1, y))
+ throw(MAL, "geos_line_verify", "GEOSCoordSeq_setY failed");
+
+ if (lineDim > 2) {
+ if (!GEOSCoordSeq_setZ(gcs_new, 1, z))
+ throw(MAL, "geos_line_verify", "GEOSCoordSeq_setZ failed");
+ }
+
+ //Create new LineString
+ if (!(res = GEOSGeom_createLineString(gcs_new))) {
+ throw(MAL, "geos_line_verify", "GEOSGeom_createLineString failed");
+ }
+
+ GEOSSetSRID(res, srid);
+ }
+ else {
+ /*TODO: Maybe here we should clone it*/
+ *res = line;
+ }
+
+ return MAL_SUCCEED;
+}
+
+str
+geos_collection_verify(const GEOSGeometry **res, const GEOSGeometry
*geosGeometry)
+{
+ GEOSGeom **newGeoms = NULL;
+ uint32_t i, numGeoms=0;
+
+ numGeoms = GEOSGetNumGeometries(geosGeometry);
+ if (numGeoms < 0)
+ throw(MAL, "geos_collection_verify", "GEOSGetNumGeometries failed");
+
+ if ( (newGeoms = (const GEOSGeometry **) GDKmalloc(sizeof(const
GEOSGeometry *)*numGeoms)) == NULL) {
+ throw(MAL, "geos_collection_verify", MAL_MALLOC_FAIL);
+ }
+
+ for (i = 0; i < numGeoms; i++) {
+ const GEOSGeometry *geom = GEOSGetGeometryN(geosGeometry, i);
+ if (!geom)
+ throw(MAL, "geos_collection_verify", "GEOSGetGeometryN failed");
+
+ if ((err = geos_geom_verify(newGeoms[i], geom)) != MAL_SUCCEED) {
+ msg = createException(MAL, "geos_collection_verify",
"geos_geom_verify failed:%s",err);
_______________________________________________
checkin-list mailing list
checkin-list@monetdb.org
https://www.monetdb.org/mailman/listinfo/checkin-list
