diff --git a/capi/geos_c.cpp b/capi/geos_c.cpp index 26a8cd0db3..bbcdb7c7e6 100644 --- a/capi/geos_c.cpp +++ b/capi/geos_c.cpp @@ -343,29 +343,34 @@ extern "C" { return GEOSNearestPoints_r(handle, g1, g2); } - GEOSGeometry* - GEOSClusterDBSCAN(Geometry* g, double eps, unsigned minPoints) { - return GEOSClusterDBSCAN_r(handle, g, eps, minPoints); + unsigned* + GEOSClusterDBSCAN(const GEOSGeometry* g, double eps, unsigned minPoints, unsigned* numClusters) + { + return GEOSClusterDBSCAN_r(handle, g, eps, minPoints, numClusters); } - GEOSGeometry* - GEOSClusterGeometryDistance(GEOSGeometry* g, double d) { - return GEOSClusterGeometryDistance_r(handle, g, d); + unsigned* + GEOSClusterGeometryDistance(const GEOSGeometry* g, double d, unsigned* numClusters) + { + return GEOSClusterGeometryDistance_r(handle, g, d, numClusters); } - GEOSGeometry* - GEOSClusterGeometryIntersects(GEOSGeometry* g) { - return GEOSClusterGeometryIntersects_r(handle, g); + unsigned* + GEOSClusterGeometryIntersects(const GEOSGeometry* g, unsigned* numClusters) + { + return GEOSClusterGeometryIntersects_r(handle, g, numClusters); } - GEOSGeometry* - GEOSClusterEnvelopeDistance(GEOSGeometry* g, double d) { - return GEOSClusterEnvelopeDistance_r(handle, g, d); + unsigned* + GEOSClusterEnvelopeDistance(const GEOSGeometry* g, double d, unsigned* numClusters) + { + return GEOSClusterEnvelopeDistance_r(handle, g, d, numClusters); } - GEOSGeometry* - GEOSClusterEnvelopeIntersects(GEOSGeometry* g) { - return GEOSClusterEnvelopeIntersects_r(handle, g); + unsigned* + GEOSClusterEnvelopeIntersects(const GEOSGeometry* g, unsigned* numClusters) + { + return GEOSClusterEnvelopeIntersects_r(handle, g, numClusters); } Geometry* diff --git a/capi/geos_c.h.in b/capi/geos_c.h.in index 6e15e994fe..69ef1f8c52 100644 --- a/capi/geos_c.h.in +++ b/capi/geos_c.h.in @@ -1822,34 +1822,38 @@ extern GEOSGeometry GEOS_DLL *GEOSGeom_transformXY_r( void* userdata); /** \see GEOSClusterDBSCAN */ -extern GEOSGeometry GEOS_DLL *GEOSClusterDBSCAN_r( - GEOSContextHandle_t handle, - GEOSGeometry* g, - double eps, - unsigned minPoints); - +extern unsigned GEOS_DLL* GEOSClusterDBSCAN_r( + GEOSContextHandle_t handle, + const GEOSGeometry* g, + double eps, + unsigned minPoints, + unsigned* numClusters); /** \see GEOSClusterGeometryDistance */ -extern GEOSGeometry GEOS_DLL *GEOSClusterGeometryDistance_r( - GEOSContextHandle_t handle, - GEOSGeometry* g, - double d); +extern unsigned GEOS_DLL* GEOSClusterGeometryDistance_r( + GEOSContextHandle_t handle, + const GEOSGeometry* g, + double d, + unsigned* numClusters); /** \see GEOSClusterGeometryIntersects */ -extern GEOSGeometry GEOS_DLL *GEOSClusterGeometryIntersects_r( - GEOSContextHandle_t handle, - GEOSGeometry* g); +extern unsigned GEOS_DLL* GEOSClusterGeometryIntersects_r( + GEOSContextHandle_t handle, + const GEOSGeometry* g, + unsigned* numClusters); /** \see GEOSClusterEnvelopeDistance */ -extern GEOSGeometry GEOS_DLL *GEOSClusterEnvelopeDistance_r( - GEOSContextHandle_t handle, - GEOSGeometry* g, - double d); +extern unsigned GEOS_DLL* GEOSClusterEnvelopeDistance_r( + GEOSContextHandle_t handle, + const GEOSGeometry* g, + double d, + unsigned* numClusters); /** \see GEOSClusterEnvelopeIntersects */ -extern GEOSGeometry GEOS_DLL *GEOSClusterEnvelopeIntersects_r( - GEOSContextHandle_t handle, - GEOSGeometry* g); +extern unsigned GEOS_DLL* GEOSClusterEnvelopeIntersects_r( + GEOSContextHandle_t handle, + const GEOSGeometry* g, + unsigned* numClusters); /* ========= Algorithms ========= */ @@ -3764,10 +3768,13 @@ extern GEOSGeometry GEOS_DLL *GEOSSharedPaths( const GEOSGeometry* g2); ///@} -/* ========== Buffer related functions ========== */ +/* ========== Clustering functions ========== */ /** @name Clustering * Functions for clustering geometries */ + +static const unsigned GEOS_CLUSTER_NONE = 4294967295; + ///@{ /** * @brief GEOSClusterDBSCAN @@ -3777,9 +3784,11 @@ extern GEOSGeometry GEOS_DLL *GEOSSharedPaths( * @param g a collection of geometries to be clustered * @param eps distance parameter for clustering * @param minPoints density parameter for clustering - * @return GeometryCollection whose components represent individual clusters + * @param numClusters the number of clusters formed + * @return an array of cluster identifiers */ -extern GEOSGeometry GEOS_DLL *GEOSClusterDBSCAN(GEOSGeometry* g, double eps, unsigned minPoints); +extern unsigned GEOS_DLL* GEOSClusterDBSCAN(const GEOSGeometry* g, double eps, unsigned minPoints, + unsigned* numClusters); /** * @brief GEOSClusterGeometryDistance @@ -3787,10 +3796,11 @@ extern GEOSGeometry GEOS_DLL *GEOSClusterDBSCAN(GEOSGeometry* g, double eps, uns * Cluster geometries according to a distance threshold * * @param g a collection of geometries to be clustered - * @param d minimim distance between geometries in the same cluster - * @return GeometryCollection whose components represent individual clusters + * @param d minimum distance between geometries in the same cluster + * @param numClusters the number of clusters formed + * @return an array of cluster identifiers */ -extern GEOSGeometry GEOS_DLL *GEOSClusterGeometryDistance(GEOSGeometry* g, double d); +extern unsigned GEOS_DLL* GEOSClusterGeometryDistance(const GEOSGeometry* g, double d, unsigned* numClusters); /** * @brief GEOSClusterGeometryIntersects @@ -3798,9 +3808,10 @@ extern GEOSGeometry GEOS_DLL *GEOSClusterGeometryDistance(GEOSGeometry* g, doubl * Cluster geometries that intersect * * @param g a collection of geometries to be clustered - * @return GeometryCollection whose components represent individual clusters + * @param numClusters the number of clusters formed + * @return an array of cluster identifiers */ -extern GEOSGeometry GEOS_DLL *GEOSClusterGeometryIntersects(GEOSGeometry* g); +extern unsigned GEOS_DLL* GEOSClusterGeometryIntersects(const GEOSGeometry* g, unsigned* numClusters); /** * @brief GEOSClusterEnvelopeDistance @@ -3808,10 +3819,11 @@ extern GEOSGeometry GEOS_DLL *GEOSClusterGeometryIntersects(GEOSGeometry* g); * Cluster geometries according to an envelope distance threshold * * @param g a collection of geometries to be clustered - * @param d minimim envelope distance between geometries in the same cluster - * @return GeometryCollection whose components represent individual clusters + * @param d minimum envelope distance between geometries in the same cluster + * @param numClusters the number of clusters formed + * @return an array of cluster identifiers */ -extern GEOSGeometry GEOS_DLL *GEOSClusterEnvelopeDistance(GEOSGeometry* g, double d); +extern unsigned GEOS_DLL* GEOSClusterEnvelopeDistance(const GEOSGeometry* g, double d, unsigned* numClusters); /** * @brief GEOSClusterEnvelopeIntersects @@ -3819,9 +3831,10 @@ extern GEOSGeometry GEOS_DLL *GEOSClusterEnvelopeDistance(GEOSGeometry* g, doubl * Cluster geometries whose envelopes intersect * * @param g - * @return GeometryCollection whose components represent individual clusters + * @param numClusters the number of clusters formed + * @return an array of cluster identifiers */ -extern GEOSGeometry GEOS_DLL *GEOSClusterEnvelopeIntersects(GEOSGeometry* g); +extern unsigned GEOS_DLL* GEOSClusterEnvelopeIntersects(const GEOSGeometry* g, unsigned* numClusters); ///@} /* ========== Buffer related functions ========== */ @@ -5018,9 +5031,9 @@ extern char GEOS_DLL GEOSEqualsIdentical( /** * Calculate the [DE9IM](https://en.wikipedia.org/wiki/DE-9IM) string for a geometry pair * and compare against a DE9IM pattern to check for -* consistency. +* consistency. * If the result matches the pattern return true. -* The pattern is a 9-character string +* The pattern is a 9-character string * containing symbols in the set "012TF*". * "012F" match the corresponding dimension symbol; * "T" matches any non-empty dimension; "*" matches any dimension. diff --git a/capi/geos_ts_c.cpp b/capi/geos_ts_c.cpp index b4cf71fc2b..cee250ac0d 100644 --- a/capi/geos_ts_c.cpp +++ b/capi/geos_ts_c.cpp @@ -953,53 +953,74 @@ extern "C" { }); } - Geometry* - GEOSClusterDBSCAN_r(GEOSContextHandle_t extHandle, Geometry* g, double eps, unsigned minPoints) + static unsigned* capi_clusters(const Geometry* g, + geos::operation::cluster::AbstractClusterFinder& finder, + unsigned* numClusters) + { + std::vector input{g->getNumGeometries()}; + for (std::size_t i = 0; i < input.size(); i++) { + input[i] = g->getGeometryN(i); + } + + auto result = finder.cluster(input); + + if (numClusters) { + *numClusters = static_cast(result.getNumClusters()); + } + + auto cluster_ids = result.getClusterIds(GEOS_CLUSTER_NONE); + unsigned* ids = (unsigned*) malloc(cluster_ids.size() * sizeof(unsigned)); + for (std::size_t i = 0; i < cluster_ids.size(); i++) { + ids[i] = static_cast(cluster_ids[i]); + } + + return ids; + } + + unsigned* + GEOSClusterDBSCAN_r(GEOSContextHandle_t extHandle, const Geometry* g, double eps, unsigned minPoints, + unsigned* numClusters) { return execute(extHandle, [&]() { - std::unique_ptr in(g); geos::operation::cluster::DBSCANClusterFinder finder(eps, minPoints); - return finder.clusterToCollection(std::move(in)).release(); + + return capi_clusters(g, finder, numClusters); }); } - Geometry* - GEOSClusterGeometryIntersects_r(GEOSContextHandle_t extHandle, Geometry* g) + unsigned* + GEOSClusterGeometryIntersects_r(GEOSContextHandle_t extHandle, const Geometry* g, unsigned* numClusters) { return execute(extHandle, [&]() { - std::unique_ptr in(g); geos::operation::cluster::GeometryIntersectsClusterFinder finder; - return finder.clusterToCollection(std::move(in)).release(); + return capi_clusters(g, finder, numClusters); }); } - Geometry* - GEOSClusterEnvelopeIntersects_r(GEOSContextHandle_t extHandle, Geometry* g) + unsigned* + GEOSClusterEnvelopeIntersects_r(GEOSContextHandle_t extHandle, const Geometry* g, unsigned* numClusters) { return execute(extHandle, [&]() { - std::unique_ptr in(g); geos::operation::cluster::EnvelopeIntersectsClusterFinder finder; - return finder.clusterToCollection(std::move(in)).release(); + return capi_clusters(g, finder, numClusters); }); } - Geometry* - GEOSClusterEnvelopeDistance_r(GEOSContextHandle_t extHandle, Geometry* g, double d) + unsigned* + GEOSClusterEnvelopeDistance_r(GEOSContextHandle_t extHandle, const Geometry* g, double d, unsigned* numClusters) { return execute(extHandle, [&]() { - std::unique_ptr in(g); geos::operation::cluster::EnvelopeDistanceClusterFinder finder(d); - return finder.clusterToCollection(std::move(in)).release(); + return capi_clusters(g, finder, numClusters); }); } - Geometry* - GEOSClusterGeometryDistance_r(GEOSContextHandle_t extHandle, Geometry* g, double d) + unsigned* + GEOSClusterGeometryDistance_r(GEOSContextHandle_t extHandle, const Geometry* g, double d, unsigned* numClusters) { return execute(extHandle, [&]() { - std::unique_ptr in(g); geos::operation::cluster::GeometryDistanceClusterFinder finder(d); - return finder.clusterToCollection(std::move(in)).release(); + return capi_clusters(g, finder, numClusters); }); } diff --git a/tests/unit/capi/GEOSClusterTest.cpp b/tests/unit/capi/GEOSClusterTest.cpp index ed87cdc2b1..9a8ea015db 100644 --- a/tests/unit/capi/GEOSClusterTest.cpp +++ b/tests/unit/capi/GEOSClusterTest.cpp @@ -14,19 +14,46 @@ namespace tut { // Common data used in test cases. struct test_capicluster_data : public capitest::utility { - - test_capicluster_data() { - m_reader = GEOSWKTReader_create(); - } - - ~test_capicluster_data() { - GEOSWKTReader_destroy(m_reader); + static GEOSGeometry* construct_clusters(const GEOSGeometry* input, + unsigned numClusters, + unsigned* clusterIds) + { + GEOSGeometry* input_cloned = GEOSGeom_clone(input); + unsigned ngeoms; + GEOSGeometry** components = GEOSGeom_releaseCollection(input_cloned, &ngeoms); + GEOSGeom_destroy(input_cloned); + + // collect the indices associated with each cluster + std::vector> cluster_component_ids(numClusters); + for (unsigned i = 0; i < ngeoms; i++) { + cluster_component_ids[clusterIds[i]].push_back(i); + } + + std::vector cluster_geoms(numClusters); + + // assemble a GeometryCollection for each cluster + for (unsigned cluster_id = 0; cluster_id < numClusters; cluster_id++) { + const std::vector component_ids = cluster_component_ids[cluster_id]; + + std::vector cluster_component_geoms(component_ids.size()); + for (std::size_t i = 0; i < component_ids.size(); i++) { + cluster_component_geoms[i] = components[component_ids[i]]; + } + + cluster_geoms[cluster_id] = GEOSGeom_createCollection(GEOS_GEOMETRYCOLLECTION, + cluster_component_geoms.data(), + static_cast(cluster_component_geoms.size())); + } + + GEOSFree(components); + + // combine the clusters into a single nested GeometryCollection + return GEOSGeom_createCollection(GEOS_GEOMETRYCOLLECTION, + cluster_geoms.data(), + static_cast(cluster_geoms.size())); } - - GEOSWKTReader* m_reader; }; - typedef test_group group; typedef group::object object; @@ -39,61 +66,104 @@ template<> template<> void object::test<1> () { - GEOSGeometry* input = GEOSWKTReader_read(m_reader, - "GEOMETRYCOLLECTION (" - "POINT (0 1)," - "LINESTRING (0 0, 0 0.1)," - "LINESTRING (0 0, 1.0 1.0)," - "POINT (0.9 1.0)," - "POINT (0 7))"); - - GEOSGeometry* result; - result = GEOSClusterEnvelopeIntersects(GEOSGeom_clone(input)); - ensure_equals("two clusters by envelope intersection", GEOSGetNumGeometries(result), 2); - GEOSGeom_destroy(result); - - result = GEOSClusterEnvelopeDistance(GEOSGeom_clone(input), 6); - ensure_equals("one cluster by envelope distance", GEOSGetNumGeometries(result), 1); - GEOSGeom_destroy(result); - - result = GEOSClusterGeometryIntersects(GEOSGeom_clone(input)); - ensure_equals("four clusters by geometry intersection", GEOSGetNumGeometries(result), 4); - GEOSGeom_destroy(result); - - result = GEOSClusterGeometryDistance(GEOSGeom_clone(input), 0.2); - ensure_equals("three clusters by distance", GEOSGetNumGeometries(result), 3); - GEOSGeom_destroy(result); - - GEOSGeom_destroy(input); + input_ = fromWKT( + "GEOMETRYCOLLECTION (" + "POINT (0 1)," + "LINESTRING (0 0, 0 0.1)," + "LINESTRING (0 0, 1.0 1.0)," + "POINT (0.9 1.0)," + "POINT (0 7))"); + + { + unsigned numClusters = 123; + unsigned* cluster_ids = GEOSClusterEnvelopeIntersects(input_, &numClusters); + ensure_equals("two clusters by envelope intersection", numClusters, 2u); + + GEOSGeometry* geom_result = construct_clusters(input_, numClusters, cluster_ids); + GEOSGeometry* geom_expected = fromWKT("GEOMETRYCOLLECTION (" + " GEOMETRYCOLLECTION (" + " POINT (0 1)," + " LINESTRING (0 0, 0 0.1)," + " LINESTRING (0 0, 1.0 1.0)," + " POINT (0.9 1.0))," + " GEOMETRYCOLLECTION (" + " POINT (0 7)))"); + + ensure_geometry_equals_identical(geom_expected, geom_result); + GEOSGeom_destroy(geom_result); + GEOSGeom_destroy(geom_expected); + + GEOSFree(cluster_ids); + } + + { + unsigned numClusters = 123; + unsigned* cluster_ids = GEOSClusterEnvelopeDistance(input_, 6, &numClusters); + ensure_equals("one cluster by envelope distance", numClusters, 1u); + + GEOSFree(cluster_ids); + } + + { + unsigned numClusters = 123; + unsigned* cluster_ids = GEOSClusterGeometryIntersects(input_, &numClusters); + ensure_equals("four clusters by geometry intersection", numClusters, 4u); + + GEOSFree(cluster_ids); + } + + { + unsigned numClusters = 123; + unsigned* cluster_ids = GEOSClusterGeometryDistance(input_, 0.2, &numClusters); + ensure_equals("three clusters by distance", numClusters, 3u); + + GEOSFree(cluster_ids); + } } template<> template<> void object::test<2>() { - GEOSGeometry* input = GEOSWKTReader_read(m_reader, - "GEOMETRYCOLLECTION (" - "POINT (0 0)," - "POINT (-1 0)," - "POINT (-1 -0.1)," - "POINT (-1 0.1)," - "POINT (1 0)," - "POINT (2 0)," - "POINT (3 0)," - "POINT ( 3 -0.1)," - "POINT ( 3 0.1)" - ")"); - - GEOSGeometry* result; - result = GEOSClusterDBSCAN(GEOSGeom_clone(input), 1.01, 5); - ensure_equals("two clusters with minPoints = 5", GEOSGetNumGeometries(result), 2); - GEOSGeom_destroy(result); - - result = GEOSClusterDBSCAN(GEOSGeom_clone(input), 1.01, 2); - ensure_equals("one cluster with minPoints = 2", GEOSGetNumGeometries(result), 1); - GEOSGeom_destroy(result); - - GEOSGeom_destroy(input); + input_ = fromWKT( + "GEOMETRYCOLLECTION (" + "POINT (0 0)," + "POINT (-1 0)," + "POINT (-1 -0.1)," + "POINT (-1 0.1)," + "POINT (1 0)," + "POINT (2 0)," + "POINT (3 0)," + "POINT ( 3 -0.1)," + "POINT ( 3 0.1)" + ")"); + + { + unsigned numClusters = 123; + unsigned* cluster_ids = GEOSClusterDBSCAN(input_, 1.01, 5, &numClusters); + ensure_equals("two clusters with minPoints = 5", numClusters, 2u); + + GEOSFree(cluster_ids); + } + + { + unsigned numClusters = 123; + unsigned* cluster_ids = GEOSClusterDBSCAN(input_, 1.01, 2, &numClusters); + ensure_equals("one cluster with minPoints = 2", numClusters, 1u); + + GEOSFree(cluster_ids); + } + + { + unsigned numClusters = 123; + unsigned* cluster_ids = GEOSClusterDBSCAN(input_, 1.01, 20, &numClusters); + ensure_equals("no clusters with minPoints = 20", numClusters, 0u); + + for (int i = 0; i < GEOSGetNumGeometries(input_); i++) { + ensure_equals(cluster_ids[i], GEOS_CLUSTER_NONE); + } + GEOSFree(cluster_ids); + } }