1 /**********************************************************************
2  *
3  * PostGIS - Spatial Types for PostgreSQL
4  * http://postgis.net
5  *
6  * Copyright 2015 Daniel Baston
7  *
8  * This is free software; you can redistribute and/or modify it under
9  * the terms of the GNU General Public Licence. See the COPYING file.
10  *
11  **********************************************************************/
12 
13 #include <stdlib.h>
14 
15 #include "CUnit/Basic.h"
16 
17 #include "../lwgeom_log.h"
18 #include "../lwgeom_geos.h"
19 #include "cu_tester.h"
20 
assert_all_results_found(LWGEOM ** results,size_t num_outputs,LWGEOM ** expected,size_t num_expected_outputs)21 static void assert_all_results_found(LWGEOM** results, size_t num_outputs, LWGEOM** expected, size_t num_expected_outputs)
22 {
23 	size_t i, j;
24 
25 	char found_equal = 0;
26 	for (i = 0; i < num_outputs; i++)
27 	{
28 		for (j = 0; j < num_expected_outputs; j++)
29 		{
30 			if (lwgeom_same(results[i], expected[j]))
31 			{
32 				found_equal = 1;
33 				break;
34 			}
35 		}
36 
37 		CU_ASSERT_TRUE(found_equal);
38 	}
39 }
40 
LWGEOMARRAY2GEOS(LWGEOM ** lw_array,size_t num_geoms)41 static GEOSGeometry** LWGEOMARRAY2GEOS(LWGEOM** lw_array, size_t num_geoms)
42 {
43 	size_t i;
44 	GEOSGeometry** geos_geoms = lwalloc(num_geoms * sizeof(GEOSGeometry*));
45 
46 	for (i = 0; i < num_geoms; i++)
47 	{
48 		geos_geoms[i] = LWGEOM2GEOS(lw_array[i], 0);
49 	}
50 
51 	return geos_geoms;
52 }
53 
GEOSARRAY2LWGEOM(GEOSGeometry ** geos_array,size_t num_geoms)54 static LWGEOM** GEOSARRAY2LWGEOM(GEOSGeometry** geos_array, size_t num_geoms)
55 {
56 	size_t i;
57 	LWGEOM** lw_geoms = lwalloc(num_geoms * sizeof(LWGEOM*));
58 
59 	for (i = 0; i < num_geoms; i++)
60 	{
61 		lw_geoms[i] = GEOS2LWGEOM(geos_array[i], 0);
62 	}
63 
64 	return lw_geoms;
65 }
66 
WKTARRAY2LWGEOM(char ** wkt_array,size_t num_geoms)67 static LWGEOM** WKTARRAY2LWGEOM(char** wkt_array, size_t num_geoms)
68 {
69 	size_t i;
70 
71 	LWGEOM** lw_geoms = lwalloc(num_geoms * sizeof(LWGEOM*));
72 
73 	for (i = 0; i < num_geoms; i++)
74 	{
75 		lw_geoms[i] = lwgeom_from_wkt(wkt_array[i], LW_PARSER_CHECK_NONE);
76 	}
77 
78 	return lw_geoms;
79 }
80 
perform_cluster_intersecting_test(char ** wkt_inputs,uint32_t num_inputs,char ** wkt_outputs,uint32_t num_outputs)81 static void perform_cluster_intersecting_test(char** wkt_inputs, uint32_t num_inputs, char** wkt_outputs, uint32_t num_outputs)
82 {
83 	GEOSGeometry** geos_results;
84 	LWGEOM** lw_results;
85 	uint32_t num_clusters;
86 
87 	LWGEOM** expected_outputs = WKTARRAY2LWGEOM(wkt_outputs, num_outputs);
88 	LWGEOM** lw_inputs = WKTARRAY2LWGEOM(wkt_inputs, num_inputs);
89 	GEOSGeometry** geos_inputs = LWGEOMARRAY2GEOS(lw_inputs, num_inputs);
90 
91 	cluster_intersecting(geos_inputs, num_inputs, &geos_results, &num_clusters);
92 	CU_ASSERT_EQUAL(num_outputs, num_clusters);
93 
94 	lw_results = GEOSARRAY2LWGEOM(geos_results, num_clusters);
95 
96 	assert_all_results_found(lw_results, num_clusters, expected_outputs, num_outputs);
97 
98 	/* Cleanup */
99 	uint32_t i;
100 	for(i = 0; i < num_clusters; i++)
101 	{
102 		GEOSGeom_destroy(geos_results[i]);
103 	}
104 	lwfree(geos_inputs);
105 	lwfree(geos_results);
106 
107 	for(i = 0; i < num_outputs; i++)
108 	{
109 		lwgeom_free(expected_outputs[i]);
110 		lwgeom_free(lw_results[i]);
111 	}
112 	lwfree(expected_outputs);
113 	lwfree(lw_results);
114 
115 	for(i = 0; i < num_inputs; i++)
116 	{
117 		lwgeom_free(lw_inputs[i]);
118 	}
119 	lwfree(lw_inputs);
120 }
121 
perform_cluster_within_distance_test(double tolerance,char ** wkt_inputs,uint32_t num_inputs,char ** wkt_outputs,uint32_t num_outputs)122 static void perform_cluster_within_distance_test(double tolerance, char** wkt_inputs, uint32_t num_inputs, char** wkt_outputs, uint32_t num_outputs)
123 {
124 	LWGEOM** lw_results;
125 	uint32_t num_clusters;
126 
127 	LWGEOM** expected_outputs = WKTARRAY2LWGEOM(wkt_outputs, num_outputs);
128 	LWGEOM** lw_inputs = WKTARRAY2LWGEOM(wkt_inputs, num_inputs);
129 
130 	cluster_within_distance(lw_inputs, num_inputs, tolerance, &lw_results, &num_clusters);
131 
132 	CU_ASSERT_EQUAL(num_outputs, num_clusters);
133 
134 	assert_all_results_found(lw_results, num_clusters, expected_outputs, num_outputs);
135 
136 	/* Cleanup */
137 	uint32_t i;
138 	for(i = 0; i < num_outputs; i++)
139 	{
140 		lwgeom_free(expected_outputs[i]);
141 		lwgeom_free(lw_results[i]);
142 	}
143 	lwfree(lw_results);
144 	lwfree(expected_outputs);
145 	lwfree(lw_inputs);
146 }
147 
init_geos_cluster_suite(void)148 static int init_geos_cluster_suite(void)
149 {
150 	initGEOS(lwnotice, lwgeom_geos_error);
151 	return 0;
152 }
153 
clean_geos_cluster_suite(void)154 static int clean_geos_cluster_suite(void)
155 {
156 	finishGEOS();
157 	return 0;
158 }
159 
basic_test(void)160 static void basic_test(void)
161 {
162 	char* a = "LINESTRING (0 0, 1 1)";
163 	char* b = "LINESTRING (1 1, 2 2)";
164 	char* c = "LINESTRING (5 5, 6 6)";
165 
166 	char* wkt_inputs_a[] = {a, b, c};
167 	char* wkt_inputs_b[] = {b, c, a};
168 	char* wkt_inputs_c[] = {c, a, b};
169 
170 	char* expected_outputs_a[] = { "GEOMETRYCOLLECTION(LINESTRING (0 0, 1 1), LINESTRING (1 1, 2 2))",
171 	                               "GEOMETRYCOLLECTION(LINESTRING (5 5, 6 6))"
172 	                             };
173 
174 	char* expected_outputs_b[] = { "GEOMETRYCOLLECTION(LINESTRING (1 1, 2 2), LINESTRING (0 0, 1 1))",
175 	                               "GEOMETRYCOLLECTION(LINESTRING (5 5, 6 6))"
176 	                             };
177 
178 	char* expected_outputs_c[] = { "GEOMETRYCOLLECTION(LINESTRING (0 0, 1 1), LINESTRING (1 1, 2 2))",
179 	                               "GEOMETRYCOLLECTION(LINESTRING (5 5, 6 6))"
180 	                             };
181 
182 	perform_cluster_intersecting_test(wkt_inputs_a, 3, expected_outputs_a, 2);
183 	perform_cluster_intersecting_test(wkt_inputs_b, 3, expected_outputs_b, 2);
184 	perform_cluster_intersecting_test(wkt_inputs_c, 3, expected_outputs_c, 2);
185 
186 	perform_cluster_within_distance_test(0, wkt_inputs_a, 3, expected_outputs_a, 2);
187 	perform_cluster_within_distance_test(0, wkt_inputs_b, 3, expected_outputs_b, 2);
188 	perform_cluster_within_distance_test(0, wkt_inputs_c, 3, expected_outputs_c, 2);
189 }
190 
basic_distance_test(void)191 static void basic_distance_test(void)
192 {
193 	char* a = "LINESTRING (0 0, 1 1)";
194 	char* b = "LINESTRING (1 1, 2 2)";
195 	char* c = "LINESTRING (5 5, 6 6)";
196 
197 	char* wkt_inputs[] = {a, b, c};
198 
199 	char* expected_outputs_all[] = {"GEOMETRYCOLLECTION(LINESTRING(0 0, 1 1), LINESTRING(1 1, 2 2), LINESTRING(5 5, 6 6))"};
200 	char* expected_outputs_partial[] = {"GEOMETRYCOLLECTION(LINESTRING(0 0, 1 1), LINESTRING(1 1, 2 2))",
201 	                                    "GEOMETRYCOLLECTION(LINESTRING(5 5, 6 6))"
202 	                                   };
203 
204 	perform_cluster_within_distance_test(0, wkt_inputs, 3, expected_outputs_partial, 2);
205 	perform_cluster_within_distance_test(sqrt(18) - 0.0000001, wkt_inputs, 3, expected_outputs_partial, 2);
206 	perform_cluster_within_distance_test(sqrt(18) + 0.0000001, wkt_inputs, 3, expected_outputs_all, 1);
207 }
208 
nonsequential_test(void)209 static void nonsequential_test(void)
210 {
211 	char* wkt_inputs[] = { "LINESTRING (0 0, 1 1)",
212 	                       "LINESTRING (1 1, 2 2)",
213 	                       "LINESTRING (5 5, 6 6)",
214 	                       "LINESTRING (5 5, 4 4)",
215 	                       "LINESTRING (3 3, 2 2)",
216 	                       "LINESTRING (3 3, 4 4)"
217 	                     };
218 
219 	char* expected_outputs[] = { "GEOMETRYCOLLECTION (LINESTRING (0 0, 1 1), LINESTRING (1 1, 2 2), LINESTRING (5 5, 6 6), LINESTRING (5 5, 4 4), LINESTRING (3 3, 2 2), LINESTRING (3 3, 4 4))" };
220 
221 	perform_cluster_intersecting_test(wkt_inputs, 6, expected_outputs, 1);
222 	perform_cluster_within_distance_test(0, wkt_inputs, 6, expected_outputs, 1);
223 }
224 
single_input_test(void)225 static void single_input_test(void)
226 {
227 	char* wkt_inputs[] = { "POINT (0 0)" };
228 	char* expected_outputs[] = { "GEOMETRYCOLLECTION (POINT (0 0))" };
229 
230 	perform_cluster_intersecting_test(wkt_inputs, 1, expected_outputs, 1);
231 	perform_cluster_within_distance_test(1, wkt_inputs, 1, expected_outputs, 1);
232 }
233 
empty_inputs_test(void)234 static void empty_inputs_test(void)
235 {
236 	char* wkt_inputs[] = { "POLYGON EMPTY", "LINESTRING EMPTY"};
237 	char* expected_outputs[] = { "GEOMETRYCOLLECTION( LINESTRING EMPTY )", "GEOMETRYCOLLECTION( POLYGON EMPTY )" };
238 
239 	perform_cluster_intersecting_test(wkt_inputs, 2, expected_outputs, 2);
240 	perform_cluster_within_distance_test(1, wkt_inputs, 2, expected_outputs, 2);
241 }
242 
multipoint_test(void)243 static void multipoint_test(void)
244 {
245 	/* See #3433 */
246 	char* wkt_inputs_mp[] = { "MULTIPOINT ((0 0), (0 1))", "POINT (0 0)"};
247 	char* expected_outputs_mp[] = { "GEOMETRYCOLLECTION(MULTIPOINT ((0 0), (0 1)), POINT (0 0))"};
248 
249 	char* wkt_inputs_gc[] = { "GEOMETRYCOLLECTION (POINT (0 0), POINT (0 1))", "POINT (0 0)"};
250 	char* expected_outputs_gc[] = { "GEOMETRYCOLLECTION(GEOMETRYCOLLECTION (POINT (0 0), POINT (0 1)), POINT (0 0))"};
251 
252 	char* wkt_inputs_pt[] = { "POINT (3 3)", "POINT (3 3)"};
253 	char* expected_outputs_pt[] = { "GEOMETRYCOLLECTION(POINT (3 3), POINT (3 3))"};
254 
255 	perform_cluster_intersecting_test(wkt_inputs_mp, 2, expected_outputs_mp, 1);
256 	perform_cluster_intersecting_test(wkt_inputs_gc, 2, expected_outputs_gc, 1);
257 	perform_cluster_intersecting_test(wkt_inputs_pt, 2, expected_outputs_pt, 1);
258 }
259 
260 struct dbscan_test_info {
261 	double eps;
262 	uint32_t min_points;
263 	uint32_t num_geoms;
264 	char** wkt_inputs;
265 	uint32_t* expected_ids;
266 	int* expected_in_cluster;
267 };
268 
do_dbscan_test(struct dbscan_test_info test)269 static void do_dbscan_test(struct dbscan_test_info test)
270 {
271 	LWGEOM** geoms = WKTARRAY2LWGEOM(test.wkt_inputs, test.num_geoms);
272 	UNIONFIND* uf = UF_create(test.num_geoms);
273 	uint32_t* ids;
274 	char* in_a_cluster;
275 	uint32_t i;
276 
277 	union_dbscan(geoms, test.num_geoms, uf, test.eps, test.min_points, &in_a_cluster);
278 	ids = UF_get_collapsed_cluster_ids(uf, in_a_cluster);
279 
280 	for (i = 0; i < test.num_geoms; i++)
281 	{
282 		ASSERT_INT_EQUAL(in_a_cluster[i], test.expected_in_cluster[i]);
283 		if (in_a_cluster[i])
284 			ASSERT_INT_EQUAL(ids[i], test.expected_ids[i]);
285 	}
286 
287 	UF_destroy(uf);
288 	for (i = 0; i < test.num_geoms; i++)
289 	{
290 		lwgeom_free(geoms[i]);
291 	}
292 	lwfree(geoms);
293 	lwfree(in_a_cluster);
294 	lwfree(ids);
295 }
296 
dbscan_test(void)297 static void dbscan_test(void)
298 {
299 	struct dbscan_test_info test;
300 	char* wkt_inputs[] = { "POINT (0 0)", "POINT (-1 0)", "POINT (-1 -0.1)", "POINT (-1 0.1)",
301 		                   "POINT (1 0)",
302 						   "POINT (2 0)", "POINT (3  0)", "POINT ( 3 -0.1)", "POINT ( 3 0.1)" };
303 	/* Although POINT (1 0) and POINT (2 0) are within eps distance of each other,
304 	 * they do not connect the two clusters because POINT (1 0) is not a core point.
305 	 * See #3572
306 	 */
307 	test.eps = 1.01;
308 	test.min_points = 5;
309 	uint32_t expected_ids[] =   { 0, 0, 0, 0, 0, 1, 1, 1, 1, 1 };
310 	int expected_in_cluster[] = { 1, 1, 1, 1, 1, 1, 1, 1, 1, 1 };
311 	test.num_geoms = sizeof(wkt_inputs) / sizeof(char*);
312 
313 	test.expected_ids = expected_ids;
314 	test.wkt_inputs = wkt_inputs;
315 	test.expected_in_cluster = expected_in_cluster;
316 	do_dbscan_test(test);
317 }
318 
dbscan_test_3612a(void)319 static void dbscan_test_3612a(void)
320 {
321 	struct dbscan_test_info test;
322 	char* wkt_inputs[] = { "POINT (1 1)" };
323 
324 	test.eps = 0.0;
325 	test.min_points = 5;
326 	uint32_t expected_ids[] = { rand() };
327 	int expected_in_cluster[] = { 0 };
328 	test.num_geoms = sizeof(wkt_inputs) / sizeof(char*);
329 
330 	test.expected_ids = expected_ids;
331 	test.expected_in_cluster = expected_in_cluster;
332 	test.wkt_inputs = wkt_inputs;
333 	do_dbscan_test(test);
334 }
335 
dbscan_test_3612b(void)336 static void dbscan_test_3612b(void)
337 {
338 	struct dbscan_test_info test;
339 	char* wkt_inputs[] = { "POINT (1 1)" };
340 
341 	test.eps = 0.0;
342 	test.min_points = 1;
343 	uint32_t expected_ids[]   = { 0 };
344 	int expected_in_cluster[] = { 1 };
345 	test.num_geoms = sizeof(wkt_inputs) / sizeof(char*);
346 
347 	test.expected_ids = expected_ids;
348 	test.expected_in_cluster = expected_in_cluster;
349 	test.wkt_inputs = wkt_inputs;
350 	do_dbscan_test(test);
351 }
352 
dbscan_test_3612c(void)353 static void dbscan_test_3612c(void)
354 {
355 	struct dbscan_test_info test;
356 	char* wkt_inputs[] = { "POLYGONM((-71.1319 42.2503 1,-71.132 42.2502 3,-71.1323 42.2504 -2,-71.1322 42.2505 1,-71.1319 42.2503 0))",
357 						   "POLYGONM((-71.1319 42.2512 0,-71.1318 42.2511 20,-71.1317 42.2511 -20,-71.1317 42.251 5,-71.1317 42.2509 4,-71.132 42.2511 6,-71.1319 42.2512 30))" };
358 	test.eps = 20.1;
359 	test.min_points = 5;
360 	uint32_t expected_ids[]   = { rand(), rand() };
361 	int expected_in_cluster[] = { 0, 0 };
362 	test.num_geoms = sizeof(wkt_inputs) / sizeof(char*);
363 
364 	test.expected_ids = expected_ids;
365 	test.expected_in_cluster = expected_in_cluster;
366 	test.wkt_inputs = wkt_inputs;
367 	do_dbscan_test(test);
368 }
369 
370 void geos_cluster_suite_setup(void);
geos_cluster_suite_setup(void)371 void geos_cluster_suite_setup(void)
372 {
373 	CU_pSuite suite = CU_add_suite("clustering", init_geos_cluster_suite, clean_geos_cluster_suite);
374 	PG_ADD_TEST(suite, basic_test);
375 	PG_ADD_TEST(suite, nonsequential_test);
376 	PG_ADD_TEST(suite, basic_distance_test);
377 	PG_ADD_TEST(suite, single_input_test);
378 	PG_ADD_TEST(suite, empty_inputs_test);
379 	PG_ADD_TEST(suite, multipoint_test);
380 	PG_ADD_TEST(suite, dbscan_test);
381 	PG_ADD_TEST(suite, dbscan_test_3612a);
382 	PG_ADD_TEST(suite, dbscan_test_3612b);
383 	PG_ADD_TEST(suite, dbscan_test_3612c);
384 }
385