PostGIS  3.6.1dev-r@@SVN_REVISION@@
postgis/lwgeom_geos.c
Go to the documentation of this file.
1 /**********************************************************************
2  *
3  * PostGIS - Spatial Types for PostgreSQL
4  * http://postgis.net
5  *
6  * PostGIS is free software: you can redistribute it and/or modify
7  * it under the terms of the GNU General Public License as published by
8  * the Free Software Foundation, either version 2 of the License, or
9  * (at your option) any later version.
10  *
11  * PostGIS is distributed in the hope that it will be useful,
12  * but WITHOUT ANY WARRANTY; without even the implied warranty of
13  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14  * GNU General Public License for more details.
15  *
16  * You should have received a copy of the GNU General Public License
17  * along with PostGIS. If not, see <http://www.gnu.org/licenses/>.
18  *
19  **********************************************************************
20  *
21  * Copyright 2009-2014 Sandro Santilli <strk@kbt.io>
22  * Copyright 2008 Paul Ramsey <pramsey@cleverelephant.ca>
23  * Copyright 2001-2003 Refractions Research Inc.
24  *
25  **********************************************************************/
26 
27 
28 #include "../postgis_config.h"
29 
30 /* PostgreSQL */
31 #include "postgres.h"
32 #include "funcapi.h"
33 #include "utils/array.h"
34 #include "utils/builtins.h"
35 #include "utils/lsyscache.h"
36 #include "utils/numeric.h"
37 #include "access/htup_details.h"
38 
39 /* PostGIS */
40 #include "lwgeom_geos.h"
41 #include "liblwgeom.h"
42 #include "liblwgeom_internal.h"
43 #include "lwgeom_itree.h"
44 #include "lwgeom_geos_prepared.h"
45 #include "lwgeom_accum.h"
46 
47 
48 
49 
50 /*
51 ** Prototypes for SQL-bound functions
52 */
53 Datum isvalid(PG_FUNCTION_ARGS);
54 Datum isvalidreason(PG_FUNCTION_ARGS);
55 Datum isvaliddetail(PG_FUNCTION_ARGS);
56 Datum buffer(PG_FUNCTION_ARGS);
57 Datum ST_Intersection(PG_FUNCTION_ARGS);
58 Datum convexhull(PG_FUNCTION_ARGS);
59 Datum topologypreservesimplify(PG_FUNCTION_ARGS);
60 Datum ST_Difference(PG_FUNCTION_ARGS);
61 Datum boundary(PG_FUNCTION_ARGS);
62 Datum ST_SymDifference(PG_FUNCTION_ARGS);
63 Datum ST_Union(PG_FUNCTION_ARGS);
64 Datum issimple(PG_FUNCTION_ARGS);
65 Datum isring(PG_FUNCTION_ARGS);
66 Datum pointonsurface(PG_FUNCTION_ARGS);
67 Datum GEOSnoop(PG_FUNCTION_ARGS);
68 Datum postgis_geos_version(PG_FUNCTION_ARGS);
69 Datum postgis_geos_compiled_version(PG_FUNCTION_ARGS);
70 Datum centroid(PG_FUNCTION_ARGS);
71 Datum polygonize_garray(PG_FUNCTION_ARGS);
72 Datum clusterintersecting_garray(PG_FUNCTION_ARGS);
73 Datum cluster_within_distance_garray(PG_FUNCTION_ARGS);
74 Datum linemerge(PG_FUNCTION_ARGS);
75 Datum hausdorffdistance(PG_FUNCTION_ARGS);
76 Datum hausdorffdistancedensify(PG_FUNCTION_ARGS);
77 Datum ST_FrechetDistance(PG_FUNCTION_ARGS);
78 Datum ST_UnaryUnion(PG_FUNCTION_ARGS);
79 Datum ST_Equals(PG_FUNCTION_ARGS);
80 Datum ST_BuildArea(PG_FUNCTION_ARGS);
81 Datum ST_DelaunayTriangles(PG_FUNCTION_ARGS);
82 Datum ST_MaximumInscribedCircle(PG_FUNCTION_ARGS);
83 Datum ST_ConcaveHull(PG_FUNCTION_ARGS);
84 Datum ST_SimplifyPolygonHull(PG_FUNCTION_ARGS);
85 Datum pgis_union_geometry_array(PG_FUNCTION_ARGS);
86 
87 /*
88 ** Prototypes end
89 */
90 
92 Datum postgis_geos_version(PG_FUNCTION_ARGS)
93 {
94  const char *ver = lwgeom_geos_version();
95  text *result = cstring_to_text(ver);
96  PG_RETURN_POINTER(result);
97 }
98 
100 Datum postgis_geos_compiled_version(PG_FUNCTION_ARGS)
101 {
102  const char *ver = lwgeom_geos_compiled_version();
103  text *result = cstring_to_text(ver);
104  PG_RETURN_POINTER(result);
105 }
106 
114 Datum hausdorffdistance(PG_FUNCTION_ARGS)
115 {
116  GSERIALIZED *geom1;
117  GSERIALIZED *geom2;
118  GEOSGeometry *g1;
119  GEOSGeometry *g2;
120  double result;
121  int retcode;
122 
123  geom1 = PG_GETARG_GSERIALIZED_P(0);
124  geom2 = PG_GETARG_GSERIALIZED_P(1);
125 
126  if ( gserialized_is_empty(geom1) || gserialized_is_empty(geom2) )
127  PG_RETURN_NULL();
128 
129  initGEOS(lwpgnotice, lwgeom_geos_error);
130 
131  g1 = POSTGIS2GEOS(geom1);
132  if (!g1)
133  HANDLE_GEOS_ERROR("First argument geometry could not be converted to GEOS");
134 
135  g2 = POSTGIS2GEOS(geom2);
136  if (!g2)
137  {
138  GEOSGeom_destroy(g1);
139  HANDLE_GEOS_ERROR("Second argument geometry could not be converted to GEOS");
140  }
141 
142  retcode = GEOSHausdorffDistance(g1, g2, &result);
143  GEOSGeom_destroy(g1);
144  GEOSGeom_destroy(g2);
145 
146  if (retcode == 0) HANDLE_GEOS_ERROR("GEOSHausdorffDistance");
147 
148  PG_FREE_IF_COPY(geom1, 0);
149  PG_FREE_IF_COPY(geom2, 1);
150 
151  PG_RETURN_FLOAT8(result);
152 }
153 
163 Datum hausdorffdistancedensify(PG_FUNCTION_ARGS)
164 {
165  GSERIALIZED *geom1;
166  GSERIALIZED *geom2;
167  GEOSGeometry *g1;
168  GEOSGeometry *g2;
169  double densifyFrac;
170  double result;
171  int retcode;
172 
173  geom1 = PG_GETARG_GSERIALIZED_P(0);
174  geom2 = PG_GETARG_GSERIALIZED_P(1);
175  densifyFrac = PG_GETARG_FLOAT8(2);
176 
177  if ( gserialized_is_empty(geom1) || gserialized_is_empty(geom2) )
178  PG_RETURN_NULL();
179 
180  initGEOS(lwpgnotice, lwgeom_geos_error);
181 
182  g1 = POSTGIS2GEOS(geom1);
183  if (!g1)
184  HANDLE_GEOS_ERROR("First argument geometry could not be converted to GEOS");
185 
186  g2 = POSTGIS2GEOS(geom2);
187  if (!g2)
188  {
189  GEOSGeom_destroy(g1);
190  HANDLE_GEOS_ERROR("Second argument geometry could not be converted to GEOS");
191  }
192 
193  retcode = GEOSHausdorffDistanceDensify(g1, g2, densifyFrac, &result);
194  GEOSGeom_destroy(g1);
195  GEOSGeom_destroy(g2);
196 
197  if (retcode == 0) HANDLE_GEOS_ERROR("GEOSHausdorffDistanceDensify");
198 
199  PG_FREE_IF_COPY(geom1, 0);
200  PG_FREE_IF_COPY(geom2, 1);
201 
202  PG_RETURN_FLOAT8(result);
203 }
204 
213 Datum ST_FrechetDistance(PG_FUNCTION_ARGS)
214 {
215  GSERIALIZED *geom1;
216  GSERIALIZED *geom2;
217  GEOSGeometry *g1;
218  GEOSGeometry *g2;
219  double densifyFrac;
220  double result;
221  int retcode;
222 
223  geom1 = PG_GETARG_GSERIALIZED_P(0);
224  geom2 = PG_GETARG_GSERIALIZED_P(1);
225  densifyFrac = PG_GETARG_FLOAT8(2);
226 
227  if ( gserialized_is_empty(geom1) || gserialized_is_empty(geom2) )
228  PG_RETURN_NULL();
229 
230  initGEOS(lwpgnotice, lwgeom_geos_error);
231 
232  g1 = POSTGIS2GEOS(geom1);
233  if (!g1)
234  HANDLE_GEOS_ERROR("First argument geometry could not be converted to GEOS");
235 
236  g2 = POSTGIS2GEOS(geom2);
237  if (!g2)
238  {
239  GEOSGeom_destroy(g1);
240  HANDLE_GEOS_ERROR("Second argument geometry could not be converted to GEOS");
241  }
242 
243  if (densifyFrac <= 0.0)
244  {
245  retcode = GEOSFrechetDistance(g1, g2, &result);
246  }
247  else
248  {
249  retcode = GEOSFrechetDistanceDensify(g1, g2, densifyFrac, &result);
250  }
251 
252  GEOSGeom_destroy(g1);
253  GEOSGeom_destroy(g2);
254 
255  if (retcode == 0) HANDLE_GEOS_ERROR("GEOSFrechetDistance");
256 
257  PG_FREE_IF_COPY(geom1, 0);
258  PG_FREE_IF_COPY(geom2, 1);
259 
260  PG_RETURN_FLOAT8(result);
261 }
262 
263 
264 
266 Datum ST_MaximumInscribedCircle(PG_FUNCTION_ARGS)
267 {
268 #if POSTGIS_GEOS_VERSION < 30900
269 
270  lwpgerror("The GEOS version this PostGIS binary "
271  "was compiled against (%d) doesn't support "
272  "'GEOSMaximumInscribedCircle' function (3.9.0+ required)",
274  PG_RETURN_NULL();
275 
276 #else /* POSTGIS_GEOS_VERSION >= 30900 */
277  GSERIALIZED* geom;
278  GSERIALIZED* center;
279  GSERIALIZED* nearest;
280  TupleDesc resultTupleDesc;
281  HeapTuple resultTuple;
282  Datum result;
283  Datum result_values[3];
284  bool result_is_null[3];
285  double radius = 0.0;
286  int32 srid = SRID_UNKNOWN;
287  bool is3d;
288 
289  if (PG_ARGISNULL(0))
290  PG_RETURN_NULL();
291 
292  geom = PG_GETARG_GSERIALIZED_P(0);
293  srid = gserialized_get_srid(geom);
294  is3d = gserialized_has_z(geom);
295 
296  /* Empty geometry? Return POINT EMPTY with zero radius */
297  if (gserialized_is_empty(geom))
298  {
301  center = geometry_serialize(lwcenter);
302  nearest = geometry_serialize(lwnearest);
303  radius = 0.0;
304  }
305  else
306  {
307  GEOSGeometry *ginput, *gcircle, *gcenter, *gnearest;
308  double width, height, size, tolerance;
309  GBOX gbox;
310  int gtype;
311  LWGEOM *lwg;
312  lwg = lwgeom_from_gserialized(geom);
313  if (!lwgeom_isfinite(lwg))
314  {
315  lwpgerror("Geometry contains invalid coordinates");
316  PG_RETURN_NULL();
317  }
318  lwgeom_free(lwg);
319 
320  if (!gserialized_get_gbox_p(geom, &gbox))
321  PG_RETURN_NULL();
322 
323  width = gbox.xmax - gbox.xmin;
324  height = gbox.ymax - gbox.ymin;
325  size = width > height ? width : height;
326  tolerance = size / 1000.0;
327 
328  initGEOS(lwpgnotice, lwgeom_geos_error);
329 
330  ginput = POSTGIS2GEOS(geom);
331  if (!ginput)
332  HANDLE_GEOS_ERROR("Geometry could not be converted to GEOS");
333 
334  gtype = gserialized_get_type(geom);
335  if (gtype == POLYGONTYPE || gtype == MULTIPOLYGONTYPE)
336  {
337  gcircle = GEOSMaximumInscribedCircle(ginput, tolerance);
338  if (!gcircle)
339  {
340  lwpgerror("Error calculating GEOSMaximumInscribedCircle.");
341  GEOSGeom_destroy(ginput);
342  PG_RETURN_NULL();
343  }
344  }
345  else
346  {
347  gcircle = GEOSLargestEmptyCircle(ginput, NULL, tolerance);
348  if (!gcircle)
349  {
350  lwpgerror("Error calculating GEOSLargestEmptyCircle.");
351  GEOSGeom_destroy(ginput);
352  PG_RETURN_NULL();
353  }
354  }
355 
356  gcenter = GEOSGeomGetStartPoint(gcircle);
357  gnearest = GEOSGeomGetEndPoint(gcircle);
358  GEOSDistance(gcenter, gnearest, &radius);
359  GEOSSetSRID(gcenter, srid);
360  GEOSSetSRID(gnearest, srid);
361 
362  center = GEOS2POSTGIS(gcenter, is3d);
363  nearest = GEOS2POSTGIS(gnearest, is3d);
364  GEOSGeom_destroy(gcenter);
365  GEOSGeom_destroy(gnearest);
366  GEOSGeom_destroy(gcircle);
367  GEOSGeom_destroy(ginput);
368  }
369 
370  get_call_result_type(fcinfo, NULL, &resultTupleDesc);
371  BlessTupleDesc(resultTupleDesc);
372 
373  result_values[0] = PointerGetDatum(center);
374  result_is_null[0] = false;
375  result_values[1] = PointerGetDatum(nearest);
376  result_is_null[1] = false;
377  result_values[2] = Float8GetDatum(radius);
378  result_is_null[2] = false;
379  resultTuple = heap_form_tuple(resultTupleDesc, result_values, result_is_null);
380 
381  result = HeapTupleGetDatum(resultTuple);
382 
383  PG_RETURN_DATUM(result);
384 
385 #endif /* POSTGIS_GEOS_VERSION >= 30900 */
386 }
387 
388 
389 /* ST_LargestEmptyCircle(geom, boundary, tolerance) */
391 Datum ST_LargestEmptyCircle(PG_FUNCTION_ARGS)
392 {
393 #if POSTGIS_GEOS_VERSION < 30900
394 
395  lwpgerror("The GEOS version this PostGIS binary "
396  "was compiled against (%d) doesn't support "
397  "'GEOSMaximumInscribedCircle' function (3.9.0+ required)",
399  PG_RETURN_NULL();
400 
401 #else /* POSTGIS_GEOS_VERSION >= 30900 */
402  GSERIALIZED* geom;
404  GSERIALIZED* center;
405  GSERIALIZED* nearest;
406  TupleDesc resultTupleDesc;
407  HeapTuple resultTuple;
408  Datum result;
409  Datum result_values[3];
410  bool result_is_null[3];
411  double radius = 0.0, tolerance = 0.0;
412  int32 srid = SRID_UNKNOWN;
413  bool is3d = false, hasBoundary = false;
414 
415  if (PG_ARGISNULL(0))
416  PG_RETURN_NULL();
417 
418  geom = PG_GETARG_GSERIALIZED_P(0);
419  tolerance = PG_GETARG_FLOAT8(1);
420  boundary = PG_GETARG_GSERIALIZED_P(2);
421  srid = gserialized_get_srid(geom);
422  is3d = gserialized_has_z(geom);
423 
425  hasBoundary = true;
426 
427  /* Empty geometry? Return POINT EMPTY with zero radius */
428  if (gserialized_is_empty(geom))
429  {
432  center = geometry_serialize(lwcenter);
433  nearest = geometry_serialize(lwnearest);
434  radius = 0.0;
435  }
436  else
437  {
438  GEOSGeometry *ginput, *gcircle, *gcenter, *gnearest;
439  GEOSGeometry *gboundary = NULL;
440  double width, height, size;
441  GBOX gbox;
442  LWGEOM *lwg;
443  lwg = lwgeom_from_gserialized(geom);
444  if (!lwgeom_isfinite(lwg))
445  {
446  lwpgerror("Geometry contains invalid coordinates");
447  PG_RETURN_NULL();
448  }
449  lwgeom_free(lwg);
450 
451 
452  if (!gserialized_get_gbox_p(geom, &gbox))
453  PG_RETURN_NULL();
454 
455  if (tolerance <= 0)
456  {
457  width = gbox.xmax - gbox.xmin;
458  height = gbox.ymax - gbox.ymin;
459  size = width > height ? width : height;
460  tolerance = size / 1000.0;
461  }
462 
463  initGEOS(lwpgnotice, lwgeom_geos_error);
464 
465  ginput = POSTGIS2GEOS(geom);
466  if (!ginput)
467  HANDLE_GEOS_ERROR("Geometry could not be converted to GEOS");
468 
469  if (hasBoundary)
470  {
471  gboundary = POSTGIS2GEOS(boundary);
472  if (!gboundary)
473  HANDLE_GEOS_ERROR("Boundary could not be converted to GEOS");
474  }
475 
476  gcircle = GEOSLargestEmptyCircle(ginput, gboundary, tolerance);
477  if (!gcircle)
478  {
479  lwpgerror("Error calculating GEOSLargestEmptyCircle.");
480  GEOSGeom_destroy(ginput);
481  PG_RETURN_NULL();
482  }
483 
484  gcenter = GEOSGeomGetStartPoint(gcircle);
485  gnearest = GEOSGeomGetEndPoint(gcircle);
486  GEOSDistance(gcenter, gnearest, &radius);
487  GEOSSetSRID(gcenter, srid);
488  GEOSSetSRID(gnearest, srid);
489 
490  center = GEOS2POSTGIS(gcenter, is3d);
491  nearest = GEOS2POSTGIS(gnearest, is3d);
492  GEOSGeom_destroy(gcenter);
493  GEOSGeom_destroy(gnearest);
494  GEOSGeom_destroy(gcircle);
495  GEOSGeom_destroy(ginput);
496  if (gboundary) GEOSGeom_destroy(gboundary);
497  }
498 
499  get_call_result_type(fcinfo, NULL, &resultTupleDesc);
500  BlessTupleDesc(resultTupleDesc);
501 
502  result_values[0] = PointerGetDatum(center);
503  result_is_null[0] = false;
504  result_values[1] = PointerGetDatum(nearest);
505  result_is_null[1] = false;
506  result_values[2] = Float8GetDatum(radius);
507  result_is_null[2] = false;
508  resultTuple = heap_form_tuple(resultTupleDesc, result_values, result_is_null);
509 
510  result = HeapTupleGetDatum(resultTuple);
511 
512  PG_RETURN_DATUM(result);
513 
514 #endif /* POSTGIS_GEOS_VERSION >= 30900 */
515 }
516 
517 
518 
527 Datum pgis_union_geometry_array(PG_FUNCTION_ARGS)
528 {
529  ArrayType *array;
530 
531  ArrayIterator iterator;
532  Datum value;
533  bool isnull;
534 
535  int is3d = LW_FALSE, gotsrid = LW_FALSE;
536  int nelems = 0, geoms_size = 0, curgeom = 0, count = 0;
537 
538  GSERIALIZED *gser_out = NULL;
539 
540  GEOSGeometry *g = NULL;
541  GEOSGeometry *g_union = NULL;
542  GEOSGeometry **geoms = NULL;
543 
544  int32_t srid = SRID_UNKNOWN;
545 
546  int empty_type = 0;
547 
548  /* Null array, null geometry (should be empty?) */
549  if ( PG_ARGISNULL(0) )
550  PG_RETURN_NULL();
551 
552  array = PG_GETARG_ARRAYTYPE_P(0);
553  nelems = ArrayGetNItems(ARR_NDIM(array), ARR_DIMS(array));
554 
555  /* Empty array? Null return */
556  if ( nelems == 0 ) PG_RETURN_NULL();
557 
558  /* Quick scan for nulls */
559  iterator = array_create_iterator(array, 0, NULL);
560  while (array_iterate(iterator, &value, &isnull))
561  {
562  /* Skip null array items */
563  if (isnull) continue;
564  count++;
565  }
566  array_free_iterator(iterator);
567 
568 
569  /* All-nulls? Return null */
570  if ( count == 0 )
571  PG_RETURN_NULL();
572 
573  /* Ok, we really need GEOS now ;) */
574  initGEOS(lwpgnotice, lwgeom_geos_error);
575 
576  /* One geom, good geom? Return it */
577  if ( count == 1 && nelems == 1 )
578  {
579 #if __GNUC__ > 4 || (__GNUC__ == 4 && __GNUC_MINOR__ >= 6)
580 #pragma GCC diagnostic push
581 #pragma GCC diagnostic ignored "-Wsign-compare"
582 #endif
583  g = POSTGIS2GEOS((GSERIALIZED *)(ARR_DATA_PTR(array)));
584 #if __GNUC__ > 4 || (__GNUC__ == 4 && __GNUC_MINOR__ >= 6)
585 #pragma GCC diagnostic pop
586 #endif
587  srid = GEOSGetSRID(g);
588  g_union = GEOSUnaryUnion(g);
589  GEOSGeom_destroy(g);
590  if (!g_union) HANDLE_GEOS_ERROR("GEOSUnaryUnion");
591  GEOSSetSRID(g_union, srid);
592  gser_out = GEOS2POSTGIS(g_union, is3d);
593  GEOSGeom_destroy(g_union);
594  PG_RETURN_POINTER(gser_out);
595  }
596 
597  /*
598  ** Collect the non-empty inputs and stuff them into a GEOS collection
599  */
600  geoms_size = nelems;
601  geoms = palloc(sizeof(GEOSGeometry*) * geoms_size);
602 
603  /*
604  ** We need to convert the array of GSERIALIZED into a GEOS collection.
605  ** First make an array of GEOS geometries.
606  */
607  iterator = array_create_iterator(array, 0, NULL);
608  while (array_iterate(iterator, &value, &isnull))
609  {
610  GSERIALIZED *gser_in;
611 
612  /* Skip null array items */
613  if (isnull) continue;
614  gser_in = (GSERIALIZED *)DatumGetPointer(value);
615 
616  /* Check for SRID mismatch in array elements */
617  if ( gotsrid )
618  gserialized_error_if_srid_mismatch_reference(gser_in, srid, __func__);
619  else
620  {
621  /* Initialize SRID/dimensions info */
622  srid = gserialized_get_srid(gser_in);
623  is3d = gserialized_has_z(gser_in);
624  gotsrid = 1;
625  }
626 
627  /* Don't include empties in the union */
628  if ( gserialized_is_empty(gser_in) )
629  {
630  int gser_type = gserialized_get_type(gser_in);
631  if (gser_type > empty_type)
632  {
633  empty_type = gser_type;
634  POSTGIS_DEBUGF(4, "empty_type = %d gser_type = %d", empty_type, gser_type);
635  }
636  }
637  else
638  {
639  g = POSTGIS2GEOS(gser_in);
640 
641  /* Uh oh! Exception thrown at construction... */
642  if ( ! g )
643  {
645  "One of the geometries in the set "
646  "could not be converted to GEOS");
647  }
648 
649  /* Ensure we have enough space in our storage array */
650  if ( curgeom == geoms_size )
651  {
652  geoms_size *= 2;
653  geoms = repalloc( geoms, sizeof(GEOSGeometry*) * geoms_size );
654  }
655 
656  geoms[curgeom] = g;
657  curgeom++;
658  }
659 
660  }
661  array_free_iterator(iterator);
662 
663  /*
664  ** Take our GEOS geometries and turn them into a GEOS collection,
665  ** then pass that into cascaded union.
666  */
667  if (curgeom > 0)
668  {
669  g = GEOSGeom_createCollection(GEOS_GEOMETRYCOLLECTION, geoms, curgeom);
670  if (!g) HANDLE_GEOS_ERROR("Could not create GEOS COLLECTION from geometry array");
671 
672  g_union = GEOSUnaryUnion(g);
673  GEOSGeom_destroy(g);
674  if (!g_union) HANDLE_GEOS_ERROR("GEOSUnaryUnion");
675 
676  GEOSSetSRID(g_union, srid);
677  gser_out = GEOS2POSTGIS(g_union, is3d);
678  GEOSGeom_destroy(g_union);
679  }
680  /* No real geometries in our array, any empties? */
681  else
682  {
683  /* If it was only empties, we'll return the largest type number */
684  if ( empty_type > 0 )
685  {
686  PG_RETURN_POINTER(geometry_serialize(lwgeom_construct_empty(empty_type, srid, is3d, 0)));
687  }
688  /* Nothing but NULL, returns NULL */
689  else
690  {
691  PG_RETURN_NULL();
692  }
693  }
694 
695  if ( ! gser_out )
696  {
697  /* Union returned a NULL geometry */
698  PG_RETURN_NULL();
699  }
700 
701  PG_RETURN_POINTER(gser_out);
702 }
703 
704 
712 Datum ST_UnaryUnion(PG_FUNCTION_ARGS)
713 {
714  GSERIALIZED *geom1;
716  LWGEOM *lwgeom1, *lwresult ;
717  double prec = -1;
718 
719  geom1 = PG_GETARG_GSERIALIZED_P(0);
720  if (PG_NARGS() > 1 && ! PG_ARGISNULL(1))
721  prec = PG_GETARG_FLOAT8(1);
722 
723  lwgeom1 = lwgeom_from_gserialized(geom1) ;
724 
725  lwresult = lwgeom_unaryunion_prec(lwgeom1, prec);
726  result = geometry_serialize(lwresult) ;
727 
728  lwgeom_free(lwgeom1) ;
729  lwgeom_free(lwresult) ;
730 
731  PG_FREE_IF_COPY(geom1, 0);
732 
733  PG_RETURN_POINTER(result);
734 }
735 
737 Datum ST_Union(PG_FUNCTION_ARGS)
738 {
739  GSERIALIZED *geom1;
740  GSERIALIZED *geom2;
742  LWGEOM *lwgeom1, *lwgeom2, *lwresult;
743  double gridSize = -1;
744 
745  geom1 = PG_GETARG_GSERIALIZED_P(0);
746  geom2 = PG_GETARG_GSERIALIZED_P(1);
747  if (PG_NARGS() > 2 && ! PG_ARGISNULL(2))
748  gridSize = PG_GETARG_FLOAT8(2);
749 
750  lwgeom1 = lwgeom_from_gserialized(geom1);
751  lwgeom2 = lwgeom_from_gserialized(geom2);
752 
753  lwresult = lwgeom_union_prec(lwgeom1, lwgeom2, gridSize);
754  result = geometry_serialize(lwresult);
755 
756  lwgeom_free(lwgeom1);
757  lwgeom_free(lwgeom2);
758  lwgeom_free(lwresult);
759 
760  PG_FREE_IF_COPY(geom1, 0);
761  PG_FREE_IF_COPY(geom2, 1);
762 
763  PG_RETURN_POINTER(result);
764 }
765 
766 /* This is retained for backward ABI compatibility
767  * with PostGIS < 3.1.0 */
769 Datum symdifference(PG_FUNCTION_ARGS)
770 {
771  PG_RETURN_DATUM(DirectFunctionCall2(
773  PG_GETARG_DATUM(0), PG_GETARG_DATUM(1)
774  ));
775 }
776 
783 Datum ST_SymDifference(PG_FUNCTION_ARGS)
784 {
785  GSERIALIZED *geom1;
786  GSERIALIZED *geom2;
788  LWGEOM *lwgeom1, *lwgeom2, *lwresult ;
789  double prec = -1;
790 
791  geom1 = PG_GETARG_GSERIALIZED_P(0);
792  geom2 = PG_GETARG_GSERIALIZED_P(1);
793  if (PG_NARGS() > 2 && ! PG_ARGISNULL(2))
794  prec = PG_GETARG_FLOAT8(2);
795 
796  lwgeom1 = lwgeom_from_gserialized(geom1) ;
797  lwgeom2 = lwgeom_from_gserialized(geom2) ;
798 
799  lwresult = lwgeom_symdifference_prec(lwgeom1, lwgeom2, prec);
800  result = geometry_serialize(lwresult) ;
801 
802  lwgeom_free(lwgeom1) ;
803  lwgeom_free(lwgeom2) ;
804  lwgeom_free(lwresult) ;
805 
806  PG_FREE_IF_COPY(geom1, 0);
807  PG_FREE_IF_COPY(geom2, 1);
808 
809  PG_RETURN_POINTER(result);
810 }
811 
813 Datum convexhull(PG_FUNCTION_ARGS)
814 {
815  GSERIALIZED *geom1;
816  GEOSGeometry *g1, *g3;
818  LWGEOM *lwout;
819  int32_t srid;
820  GBOX bbox;
821 
822  geom1 = PG_GETARG_GSERIALIZED_P(0);
823 
824  /* Empty.ConvexHull() == Empty */
825  if ( gserialized_is_empty(geom1) )
826  PG_RETURN_POINTER(geom1);
827 
828  srid = gserialized_get_srid(geom1);
829 
830  initGEOS(lwpgnotice, lwgeom_geos_error);
831 
832  g1 = POSTGIS2GEOS(geom1);
833 
834  if (!g1)
835  HANDLE_GEOS_ERROR("First argument geometry could not be converted to GEOS");
836 
837  g3 = GEOSConvexHull(g1);
838  GEOSGeom_destroy(g1);
839 
840  if (!g3) HANDLE_GEOS_ERROR("GEOSConvexHull");
841 
842  POSTGIS_DEBUGF(3, "result: %s", GEOSGeomToWKT(g3));
843 
844  GEOSSetSRID(g3, srid);
845 
846  lwout = GEOS2LWGEOM(g3, gserialized_has_z(geom1));
847  GEOSGeom_destroy(g3);
848 
849  if (!lwout)
850  {
851  elog(ERROR,
852  "convexhull() failed to convert GEOS geometry to LWGEOM");
853  PG_RETURN_NULL(); /* never get here */
854  }
855 
856  /* Copy input bbox if any */
857  if ( gserialized_get_gbox_p(geom1, &bbox) )
858  {
859  /* Force the box to have the same dimensionality as the lwgeom */
860  bbox.flags = lwout->flags;
861  lwout->bbox = gbox_copy(&bbox);
862  }
863 
864  result = geometry_serialize(lwout);
865  lwgeom_free(lwout);
866 
867  if (!result)
868  {
869  elog(ERROR,"GEOS convexhull() threw an error (result postgis geometry formation)!");
870  PG_RETURN_NULL(); /* never get here */
871  }
872 
873  PG_FREE_IF_COPY(geom1, 0);
874  PG_RETURN_POINTER(result);
875 }
876 
877 
879 Datum ST_ConcaveHull(PG_FUNCTION_ARGS)
880 {
881 #if POSTGIS_GEOS_VERSION < 31100
882 
883  lwpgerror("The GEOS version this PostGIS binary "
884  "was compiled against (%d) doesn't support "
885  "'GEOSConcaveHull' function (3.11.0+ required)",
887  PG_RETURN_NULL();
888 
889 #else /* POSTGIS_GEOS_VERSION >= 31100 */
890  GSERIALIZED* geom = PG_GETARG_GSERIALIZED_P(0);
891  double ratio = PG_GETARG_FLOAT8(1);
892  bool allow_holes = PG_GETARG_BOOL(2);
893 
894  LWGEOM* lwgeom = lwgeom_from_gserialized(geom);
895  LWGEOM* lwresult = lwgeom_concavehull(lwgeom, ratio, allow_holes);
896  GSERIALIZED* result = geometry_serialize(lwresult);
897 
898  lwgeom_free(lwgeom);
899  lwgeom_free(lwresult);
900  PG_FREE_IF_COPY(geom, 0);
901  PG_RETURN_POINTER(result);
902 #endif
903 }
904 
905 
907 Datum ST_SimplifyPolygonHull(PG_FUNCTION_ARGS)
908 {
909 #if POSTGIS_GEOS_VERSION < 31100
910 
911  lwpgerror("The GEOS version this PostGIS binary "
912  "was compiled against (%d) doesn't support "
913  "'ST_SimplifyPolygonHull' function (3.11.0+ required)",
915  PG_RETURN_NULL();
916 
917 #else /* POSTGIS_GEOS_VERSION >= 31100 */
918  GSERIALIZED* geom = PG_GETARG_GSERIALIZED_P(0);
919  double vertex_fraction = PG_GETARG_FLOAT8(1);
920  uint32_t is_outer = PG_GETARG_BOOL(2);
921 
922  LWGEOM* lwgeom = lwgeom_from_gserialized(geom);
923  LWGEOM* lwresult = lwgeom_simplify_polygonal(lwgeom, vertex_fraction, is_outer);
924  GSERIALIZED* result = geometry_serialize(lwresult);
925 
926  lwgeom_free(lwgeom);
927  lwgeom_free(lwresult);
928  PG_FREE_IF_COPY(geom, 0);
929  PG_RETURN_POINTER(result);
930 #endif
931 }
932 
933 
935 Datum topologypreservesimplify(PG_FUNCTION_ARGS)
936 {
937  GSERIALIZED *gs1;
938  LWGEOM *lwg1;
939  double tolerance;
940  GEOSGeometry *g1, *g3;
942  uint32_t type;
943 
944  gs1 = PG_GETARG_GSERIALIZED_P(0);
945  tolerance = PG_GETARG_FLOAT8(1);
946  lwg1 = lwgeom_from_gserialized(gs1);
947 
948  /* Empty.Simplify() == Empty */
949  type = lwgeom_get_type(lwg1);
950  if (lwgeom_is_empty(lwg1) || type == TINTYPE || type == TRIANGLETYPE)
951  PG_RETURN_POINTER(gs1);
952 
953  if (!lwgeom_isfinite(lwg1))
954  {
955  lwpgerror("Geometry contains invalid coordinates");
956  PG_RETURN_NULL();
957  }
958 
959  initGEOS(lwpgnotice, lwgeom_geos_error);
960 
961  g1 = LWGEOM2GEOS(lwg1, LW_TRUE);
962  lwgeom_free(lwg1);
963  if (!g1)
964  HANDLE_GEOS_ERROR("First argument geometry could not be converted to GEOS");
965 
966  g3 = GEOSTopologyPreserveSimplify(g1,tolerance);
967  GEOSGeom_destroy(g1);
968 
969  if (!g3) HANDLE_GEOS_ERROR("GEOSTopologyPreserveSimplify");
970 
971  POSTGIS_DEBUGF(3, "result: %s", GEOSGeomToWKT(g3));
972 
973  GEOSSetSRID(g3, gserialized_get_srid(gs1));
974 
976  GEOSGeom_destroy(g3);
977 
978  if (!result)
979  {
980  elog(ERROR,"GEOS topologypreservesimplify() threw an error (result postgis geometry formation)!");
981  PG_RETURN_NULL(); /* never get here */
982  }
983 
984  PG_FREE_IF_COPY(gs1, 0);
985  PG_RETURN_POINTER(result);
986 }
987 
989 Datum buffer(PG_FUNCTION_ARGS)
990 {
991  GEOSBufferParams *bufferparams;
992  GEOSGeometry *g1, *g3 = NULL;
994  LWGEOM *lwg;
995  int quadsegs = 8; /* the default */
996  int singleside = 0; /* the default */
997  enum
998  {
999  ENDCAP_ROUND = 1,
1000  ENDCAP_FLAT = 2,
1001  ENDCAP_SQUARE = 3
1002  };
1003  enum
1004  {
1005  JOIN_ROUND = 1,
1006  JOIN_MITRE = 2,
1007  JOIN_BEVEL = 3
1008  };
1009  const double DEFAULT_MITRE_LIMIT = 5.0;
1010  const int DEFAULT_ENDCAP_STYLE = ENDCAP_ROUND;
1011  const int DEFAULT_JOIN_STYLE = JOIN_ROUND;
1012  double mitreLimit = DEFAULT_MITRE_LIMIT;
1013  int endCapStyle = DEFAULT_ENDCAP_STYLE;
1014  int joinStyle = DEFAULT_JOIN_STYLE;
1015 
1016  GSERIALIZED *geom1 = PG_GETARG_GSERIALIZED_P(0);
1017  double size = PG_GETARG_FLOAT8(1);
1018  text *params_text;
1019 
1020  if (PG_NARGS() > 2)
1021  {
1022  params_text = PG_GETARG_TEXT_P(2);
1023  }
1024  else
1025  {
1026  params_text = cstring_to_text("");
1027  }
1028 
1029  /* Empty.Buffer() == Empty[polygon] */
1030  if ( gserialized_is_empty(geom1) )
1031  {
1033  gserialized_get_srid(geom1),
1034  0, 0)); // buffer wouldn't give back z or m anyway
1035  PG_RETURN_POINTER(geometry_serialize(lwg));
1036  }
1037 
1038  lwg = lwgeom_from_gserialized(geom1);
1039 
1040  if (!lwgeom_isfinite(lwg))
1041  {
1042  lwpgerror("Geometry contains invalid coordinates");
1043  PG_RETURN_NULL();
1044  }
1045  lwgeom_free(lwg);
1046 
1047  initGEOS(lwpgnotice, lwgeom_geos_error);
1048 
1049  g1 = POSTGIS2GEOS(geom1);
1050  if (!g1)
1051  HANDLE_GEOS_ERROR("First argument geometry could not be converted to GEOS");
1052 
1053 
1054  if (VARSIZE_ANY_EXHDR(params_text) > 0)
1055  {
1056  char *param;
1057  char *params = text_to_cstring(params_text);
1058 
1059  for (param=params; ; param=NULL)
1060  {
1061  char *key, *val;
1062  param = strtok(param, " ");
1063  if (!param) break;
1064  POSTGIS_DEBUGF(3, "Param: %s", param);
1065 
1066  key = param;
1067  val = strchr(key, '=');
1068  if (!val || *(val + 1) == '\0')
1069  {
1070  lwpgerror("Missing value for buffer parameter %s", key);
1071  break;
1072  }
1073  *val = '\0';
1074  ++val;
1075 
1076  POSTGIS_DEBUGF(3, "Param: %s : %s", key, val);
1077 
1078  if ( !strcmp(key, "endcap") )
1079  {
1080  /* Supported end cap styles:
1081  * "round", "flat", "square"
1082  */
1083  if ( !strcmp(val, "round") )
1084  {
1085  endCapStyle = ENDCAP_ROUND;
1086  }
1087  else if ( !strcmp(val, "flat") ||
1088  !strcmp(val, "butt") )
1089  {
1090  endCapStyle = ENDCAP_FLAT;
1091  }
1092  else if ( !strcmp(val, "square") )
1093  {
1094  endCapStyle = ENDCAP_SQUARE;
1095  }
1096  else
1097  {
1098  lwpgerror("Invalid buffer end cap "
1099  "style: %s (accept: "
1100  "'round', 'flat', 'butt' "
1101  "or 'square'"
1102  ")", val);
1103  break;
1104  }
1105 
1106  }
1107  else if ( !strcmp(key, "join") )
1108  {
1109  if ( !strcmp(val, "round") )
1110  {
1111  joinStyle = JOIN_ROUND;
1112  }
1113  else if ( !strcmp(val, "mitre") ||
1114  !strcmp(val, "miter") )
1115  {
1116  joinStyle = JOIN_MITRE;
1117  }
1118  else if ( !strcmp(val, "bevel") )
1119  {
1120  joinStyle = JOIN_BEVEL;
1121  }
1122  else
1123  {
1124  lwpgerror("Invalid buffer end cap "
1125  "style: %s (accept: "
1126  "'round', 'mitre', 'miter' "
1127  " or 'bevel'"
1128  ")", val);
1129  break;
1130  }
1131  }
1132  else if ( !strcmp(key, "mitre_limit") ||
1133  !strcmp(key, "miter_limit") )
1134  {
1135  /* mitreLimit is a float */
1136  mitreLimit = atof(val);
1137  }
1138  else if ( !strcmp(key, "quad_segs") )
1139  {
1140  /* quadrant segments is an int */
1141  quadsegs = atoi(val);
1142  }
1143  else if ( !strcmp(key, "side") )
1144  {
1145  if ( !strcmp(val, "both") )
1146  {
1147  singleside = 0;
1148  }
1149  else if ( !strcmp(val, "left") )
1150  {
1151  singleside = 1;
1152  }
1153  else if ( !strcmp(val, "right") )
1154  {
1155  singleside = 1;
1156  size *= -1;
1157  }
1158  else
1159  {
1160  lwpgerror("Invalid side parameter: %s (accept: 'right', 'left', 'both')", val);
1161  break;
1162  }
1163  }
1164  else
1165  {
1166  lwpgerror(
1167  "Invalid buffer parameter: %s (accept: 'endcap', 'join', 'mitre_limit', 'miter_limit', 'quad_segs' and 'side')",
1168  key);
1169  break;
1170  }
1171  }
1172  pfree(params); /* was pstrduped */
1173  }
1174 
1175 
1176  POSTGIS_DEBUGF(3, "endCap:%d joinStyle:%d mitreLimit:%g",
1177  endCapStyle, joinStyle, mitreLimit);
1178 
1179  bufferparams = GEOSBufferParams_create();
1180  if (bufferparams)
1181  {
1182  if (GEOSBufferParams_setEndCapStyle(bufferparams, endCapStyle) &&
1183  GEOSBufferParams_setJoinStyle(bufferparams, joinStyle) &&
1184  GEOSBufferParams_setMitreLimit(bufferparams, mitreLimit) &&
1185  GEOSBufferParams_setQuadrantSegments(bufferparams, quadsegs) &&
1186  GEOSBufferParams_setSingleSided(bufferparams, singleside))
1187  {
1188  g3 = GEOSBufferWithParams(g1, bufferparams, size);
1189  }
1190  else
1191  {
1192  lwpgerror("Error setting buffer parameters.");
1193  }
1194  GEOSBufferParams_destroy(bufferparams);
1195  }
1196  else
1197  {
1198  lwpgerror("Error setting buffer parameters.");
1199  }
1200 
1201  GEOSGeom_destroy(g1);
1202 
1203  if (!g3) HANDLE_GEOS_ERROR("GEOSBuffer");
1204 
1205  POSTGIS_DEBUGF(3, "result: %s", GEOSGeomToWKT(g3));
1206 
1207  GEOSSetSRID(g3, gserialized_get_srid(geom1));
1208 
1209  result = GEOS2POSTGIS(g3, gserialized_has_z(geom1));
1210  GEOSGeom_destroy(g3);
1211 
1212  if (!result)
1213  {
1214  elog(ERROR,"GEOS buffer() threw an error (result postgis geometry formation)!");
1215  PG_RETURN_NULL(); /* never get here */
1216  }
1217 
1218  PG_FREE_IF_COPY(geom1, 0);
1219  PG_RETURN_POINTER(result);
1220 }
1221 
1222 /*
1223 * Generate a field of random points within the area of a
1224 * polygon or multipolygon. Throws an error for other geometry
1225 * types.
1226 */
1227 Datum ST_GeneratePoints(PG_FUNCTION_ARGS);
1229 Datum ST_GeneratePoints(PG_FUNCTION_ARGS)
1230 {
1231  GSERIALIZED *gser_input;
1232  GSERIALIZED *gser_result;
1233  LWGEOM *lwgeom_input;
1234  LWGEOM *lwgeom_result;
1235  int32 npoints;
1236  int32 seed = 0;
1237 
1238  gser_input = PG_GETARG_GSERIALIZED_P(0);
1239  npoints = PG_GETARG_INT32(1);
1240 
1241  if (npoints < 0)
1242  PG_RETURN_NULL();
1243 
1244  if (PG_NARGS() > 2 && ! PG_ARGISNULL(2))
1245  {
1246  seed = PG_GETARG_INT32(2);
1247  if (seed < 1)
1248  {
1249  lwpgerror("ST_GeneratePoints: seed must be greater than zero");
1250  PG_RETURN_NULL();
1251  }
1252  }
1253 
1254  /* Types get checked in the code, we'll keep things small out there */
1255  lwgeom_input = lwgeom_from_gserialized(gser_input);
1256  lwgeom_result = (LWGEOM*)lwgeom_to_points(lwgeom_input, npoints, seed);
1257  lwgeom_free(lwgeom_input);
1258  PG_FREE_IF_COPY(gser_input, 0);
1259 
1260  /* Return null as null */
1261  if (!lwgeom_result)
1262  PG_RETURN_NULL();
1263 
1264  /* Serialize and return */
1265  gser_result = geometry_serialize(lwgeom_result);
1266  lwgeom_free(lwgeom_result);
1267  PG_RETURN_POINTER(gser_result);
1268 }
1269 
1270 
1271 /*
1272 * Compute at offset curve to a line
1273 */
1274 Datum ST_OffsetCurve(PG_FUNCTION_ARGS);
1276 Datum ST_OffsetCurve(PG_FUNCTION_ARGS)
1277 {
1278  GSERIALIZED *gser_input;
1279  GSERIALIZED *gser_result;
1280  LWGEOM *lwgeom_input;
1281  LWGEOM *lwgeom_result;
1282  double size;
1283  int quadsegs = 8; /* the default */
1284  int nargs;
1285 
1286  enum
1287  {
1288  JOIN_ROUND = 1,
1289  JOIN_MITRE = 2,
1290  JOIN_BEVEL = 3
1291  };
1292 
1293  static const double DEFAULT_MITRE_LIMIT = 5.0;
1294  static const int DEFAULT_JOIN_STYLE = JOIN_ROUND;
1295  double mitreLimit = DEFAULT_MITRE_LIMIT;
1296  int joinStyle = DEFAULT_JOIN_STYLE;
1297  char *param = NULL;
1298  char *paramstr = NULL;
1299 
1300  /* Read SQL arguments */
1301  nargs = PG_NARGS();
1302  gser_input = PG_GETARG_GSERIALIZED_P(0);
1303  size = PG_GETARG_FLOAT8(1);
1304 
1305  /* For distance == 0, just return the input. */
1306  if (size == 0) PG_RETURN_POINTER(gser_input);
1307 
1308  /* Read the lwgeom, check for errors */
1309  lwgeom_input = lwgeom_from_gserialized(gser_input);
1310  if ( ! lwgeom_input )
1311  lwpgerror("ST_OffsetCurve: lwgeom_from_gserialized returned NULL");
1312 
1313  /* For empty inputs, just echo them back */
1314  if ( lwgeom_is_empty(lwgeom_input) )
1315  PG_RETURN_POINTER(gser_input);
1316 
1317  /* Process the optional arguments */
1318  if ( nargs > 2 )
1319  {
1320  text *wkttext = PG_GETARG_TEXT_P(2);
1321  paramstr = text_to_cstring(wkttext);
1322 
1323  POSTGIS_DEBUGF(3, "paramstr: %s", paramstr);
1324 
1325  for ( param=paramstr; ; param=NULL )
1326  {
1327  char *key, *val;
1328  param = strtok(param, " ");
1329  if (!param) break;
1330  POSTGIS_DEBUGF(3, "Param: %s", param);
1331 
1332  key = param;
1333  val = strchr(key, '=');
1334  if (!val || *(val + 1) == '\0')
1335  {
1336  lwpgerror("ST_OffsetCurve: Missing value for buffer parameter %s", key);
1337  break;
1338  }
1339  *val = '\0';
1340  ++val;
1341 
1342  POSTGIS_DEBUGF(3, "Param: %s : %s", key, val);
1343 
1344  if ( !strcmp(key, "join") )
1345  {
1346  if ( !strcmp(val, "round") )
1347  {
1348  joinStyle = JOIN_ROUND;
1349  }
1350  else if ( !(strcmp(val, "mitre") && strcmp(val, "miter")) )
1351  {
1352  joinStyle = JOIN_MITRE;
1353  }
1354  else if ( ! strcmp(val, "bevel") )
1355  {
1356  joinStyle = JOIN_BEVEL;
1357  }
1358  else
1359  {
1360  lwpgerror(
1361  "Invalid buffer end cap style: %s (accept: 'round', 'mitre', 'miter' or 'bevel')",
1362  val);
1363  break;
1364  }
1365  }
1366  else if ( !strcmp(key, "mitre_limit") ||
1367  !strcmp(key, "miter_limit") )
1368  {
1369  /* mitreLimit is a float */
1370  mitreLimit = atof(val);
1371  }
1372  else if ( !strcmp(key, "quad_segs") )
1373  {
1374  /* quadrant segments is an int */
1375  quadsegs = atoi(val);
1376  }
1377  else
1378  {
1379  lwpgerror(
1380  "Invalid buffer parameter: %s (accept: 'join', 'mitre_limit', 'miter_limit and 'quad_segs')",
1381  key);
1382  break;
1383  }
1384  }
1385  POSTGIS_DEBUGF(3, "joinStyle:%d mitreLimit:%g", joinStyle, mitreLimit);
1386  pfree(paramstr); /* alloc'ed in text_to_cstring */
1387  }
1388 
1389  lwgeom_result = lwgeom_offsetcurve(lwgeom_input, size, quadsegs, joinStyle, mitreLimit);
1390 
1391  if (!lwgeom_result)
1392  lwpgerror("ST_OffsetCurve: lwgeom_offsetcurve returned NULL");
1393 
1394  gser_result = geometry_serialize(lwgeom_result);
1395  lwgeom_free(lwgeom_input);
1396  lwgeom_free(lwgeom_result);
1397  PG_RETURN_POINTER(gser_result);
1398 }
1399 
1401 Datum ST_Intersection(PG_FUNCTION_ARGS)
1402 {
1403  GSERIALIZED *geom1;
1404  GSERIALIZED *geom2;
1406  LWGEOM *lwgeom1, *lwgeom2, *lwresult;
1407  double prec = -1;
1408 
1409  geom1 = PG_GETARG_GSERIALIZED_P(0);
1410  geom2 = PG_GETARG_GSERIALIZED_P(1);
1411  if (PG_NARGS() > 2 && ! PG_ARGISNULL(2))
1412  prec = PG_GETARG_FLOAT8(2);
1413 
1414  lwgeom1 = lwgeom_from_gserialized(geom1);
1415  lwgeom2 = lwgeom_from_gserialized(geom2);
1416 
1417  lwresult = lwgeom_intersection_prec(lwgeom1, lwgeom2, prec);
1418  result = geometry_serialize(lwresult);
1419 
1420  lwgeom_free(lwgeom1);
1421  lwgeom_free(lwgeom2);
1422  lwgeom_free(lwresult);
1423 
1424  PG_FREE_IF_COPY(geom1, 0);
1425  PG_FREE_IF_COPY(geom2, 1);
1426 
1427  PG_RETURN_POINTER(result);
1428 }
1429 
1431 Datum ST_Difference(PG_FUNCTION_ARGS)
1432 {
1433  GSERIALIZED *geom1;
1434  GSERIALIZED *geom2;
1436  LWGEOM *lwgeom1, *lwgeom2, *lwresult;
1437  double prec = -1;
1438 
1439  geom1 = PG_GETARG_GSERIALIZED_P(0);
1440  geom2 = PG_GETARG_GSERIALIZED_P(1);
1441  if (PG_NARGS() > 2 && ! PG_ARGISNULL(2))
1442  prec = PG_GETARG_FLOAT8(2);
1443 
1444  lwgeom1 = lwgeom_from_gserialized(geom1);
1445  lwgeom2 = lwgeom_from_gserialized(geom2);
1446 
1447  lwresult = lwgeom_difference_prec(lwgeom1, lwgeom2, prec);
1448  result = geometry_serialize(lwresult);
1449 
1450  lwgeom_free(lwgeom1);
1451  lwgeom_free(lwgeom2);
1452  lwgeom_free(lwresult);
1453 
1454  PG_FREE_IF_COPY(geom1, 0);
1455  PG_FREE_IF_COPY(geom2, 1);
1456 
1457  PG_RETURN_POINTER(result);
1458 }
1459 
1465 Datum pointonsurface(PG_FUNCTION_ARGS)
1466 {
1467  GSERIALIZED *geom, *result;
1468  LWGEOM *lwgeom, *lwresult;
1469 
1470  geom = PG_GETARG_GSERIALIZED_P(0);
1471 
1472  lwgeom = lwgeom_from_gserialized(geom);
1473  lwresult = lwgeom_pointonsurface(lwgeom);
1474  lwgeom_free(lwgeom);
1475  PG_FREE_IF_COPY(geom, 0);
1476 
1477  if (!lwresult) PG_RETURN_NULL();
1478 
1479  result = geometry_serialize(lwresult);
1480  lwgeom_free(lwresult);
1481  PG_RETURN_POINTER(result);
1482 }
1483 
1485 Datum centroid(PG_FUNCTION_ARGS)
1486 {
1487  GSERIALIZED *geom, *result;
1488  LWGEOM *lwgeom, *lwresult;
1489 
1490  geom = PG_GETARG_GSERIALIZED_P(0);
1491 
1492  lwgeom = lwgeom_from_gserialized(geom);
1493  lwresult = lwgeom_centroid(lwgeom);
1494  lwgeom_free(lwgeom);
1495  PG_FREE_IF_COPY(geom, 0);
1496 
1497  if (!lwresult) PG_RETURN_NULL();
1498 
1499  result = geometry_serialize(lwresult);
1500  lwgeom_free(lwresult);
1501  PG_RETURN_POINTER(result);
1502 }
1503 
1504 Datum ST_ReducePrecision(PG_FUNCTION_ARGS);
1506 Datum ST_ReducePrecision(PG_FUNCTION_ARGS)
1507 {
1508  GSERIALIZED *geom, *result;
1509  LWGEOM *lwgeom, *lwresult;
1510  double gridSize = PG_GETARG_FLOAT8(1);
1511  geom = PG_GETARG_GSERIALIZED_P(0);
1512 
1513  lwgeom = lwgeom_from_gserialized(geom);
1514  lwresult = lwgeom_reduceprecision(lwgeom, gridSize);
1515  lwgeom_free(lwgeom);
1516  PG_FREE_IF_COPY(geom, 0);
1517 
1518  if (!lwresult) PG_RETURN_NULL();
1519 
1520  result = geometry_serialize(lwresult);
1521  lwgeom_free(lwresult);
1522  PG_RETURN_POINTER(result);
1523 }
1524 
1525 Datum ST_ClipByBox2d(PG_FUNCTION_ARGS);
1527 Datum ST_ClipByBox2d(PG_FUNCTION_ARGS)
1528 {
1529  static const uint32_t geom_idx = 0;
1530  static const uint32_t box2d_idx = 1;
1532  LWGEOM *lwgeom1, *lwresult ;
1533  GBOX bbox1 = {0};
1534  GBOX *bbox2;
1535  uint8_t type;
1536  int32_t srid;
1537  lwflags_t flags;
1538 
1539  if (!gserialized_datum_get_internals_p(PG_GETARG_DATUM(geom_idx), &bbox1, &flags, &type, &srid))
1540  {
1541  /* empty clips to empty, no matter rect */
1542  PG_RETURN_DATUM(PG_GETARG_DATUM(geom_idx));
1543  }
1544 
1545  /* WARNING: this is really a BOX2DF, use only xmin and ymin fields */
1546  bbox2 = (GBOX *)PG_GETARG_POINTER(box2d_idx);
1547  bbox2->flags = 0;
1548 
1549  /* If bbox1 is strictly contained by bbox2, return input geometry */
1550  if (bbox2->xmin < bbox1.xmin && bbox2->xmax > bbox1.xmax &&
1551  bbox2->ymin < bbox1.ymin && bbox2->ymax > bbox1.ymax)
1552  {
1553  PG_RETURN_DATUM(PG_GETARG_DATUM(geom_idx));
1554  }
1555 
1556  /* If bbox1 outside of bbox2, return empty */
1557  if (!gbox_overlaps_2d(&bbox1, bbox2))
1558  {
1559  /* Get type and srid from datum */
1560  lwresult = lwgeom_construct_empty(type, srid, 0, 0);
1561  result = geometry_serialize(lwresult) ;
1562  lwgeom_free(lwresult) ;
1563  PG_RETURN_POINTER(result);
1564  }
1565 
1566  lwgeom1 = lwgeom_from_gserialized(PG_GETARG_GSERIALIZED_P(geom_idx));
1567  lwresult = lwgeom_clip_by_rect(lwgeom1, bbox2->xmin, bbox2->ymin,
1568  bbox2->xmax, bbox2->ymax);
1569 
1570  lwgeom_free(lwgeom1);
1571 
1572  if (!lwresult)
1573  PG_RETURN_NULL();
1574 
1575  result = geometry_serialize(lwresult) ;
1576  PG_RETURN_POINTER(result);
1577 }
1578 
1579 
1580 /*---------------------------------------------*/
1581 
1583 Datum isvalid(PG_FUNCTION_ARGS)
1584 {
1585  GSERIALIZED *geom1;
1586  LWGEOM *lwgeom;
1587  char result;
1588  GEOSGeom g1;
1589 
1590  geom1 = PG_GETARG_GSERIALIZED_P(0);
1591 
1592  /* Empty.IsValid() == TRUE */
1593  if ( gserialized_is_empty(geom1) )
1594  PG_RETURN_BOOL(true);
1595 
1596  initGEOS(lwpgnotice, lwgeom_geos_error);
1597 
1598  lwgeom = lwgeom_from_gserialized(geom1);
1599  if ( ! lwgeom )
1600  {
1601  lwpgerror("unable to deserialize input");
1602  }
1603  g1 = LWGEOM2GEOS(lwgeom, 0);
1604  lwgeom_free(lwgeom);
1605 
1606  if ( ! g1 )
1607  {
1608  PG_RETURN_BOOL(false);
1609  }
1610 
1611  result = GEOSisValid(g1);
1612  GEOSGeom_destroy(g1);
1613 
1614  if (result == 2)
1615  {
1616  elog(ERROR,"GEOS isvalid() threw an error!");
1617  PG_RETURN_NULL(); /*never get here */
1618  }
1619 
1620  PG_FREE_IF_COPY(geom1, 0);
1621  PG_RETURN_BOOL(result);
1622 }
1623 
1625 Datum isvalidreason(PG_FUNCTION_ARGS)
1626 {
1627  GSERIALIZED *geom = NULL;
1628  char *reason_str = NULL;
1629  text *result = NULL;
1630  const GEOSGeometry *g1 = NULL;
1631 
1632  geom = PG_GETARG_GSERIALIZED_P(0);
1633 
1634  initGEOS(lwpgnotice, lwgeom_geos_error);
1635 
1636  g1 = POSTGIS2GEOS(geom);
1637  if ( g1 )
1638  {
1639  reason_str = GEOSisValidReason(g1);
1640  GEOSGeom_destroy((GEOSGeometry *)g1);
1641  if (!reason_str) HANDLE_GEOS_ERROR("GEOSisValidReason");
1642  result = cstring_to_text(reason_str);
1643  GEOSFree(reason_str);
1644  }
1645  else
1646  {
1647  result = cstring_to_text(lwgeom_geos_errmsg);
1648  }
1649 
1650  PG_FREE_IF_COPY(geom, 0);
1651  PG_RETURN_POINTER(result);
1652 }
1653 
1655 Datum isvaliddetail(PG_FUNCTION_ARGS)
1656 {
1657  GSERIALIZED *geom = NULL;
1658  const GEOSGeometry *g1 = NULL;
1659  char *values[3]; /* valid bool, reason text, location geometry */
1660  char *geos_reason = NULL;
1661  char *reason = NULL;
1662  GEOSGeometry *geos_location = NULL;
1663  LWGEOM *location = NULL;
1664  char valid = 0;
1665  HeapTupleHeader result;
1666  TupleDesc tupdesc;
1667  HeapTuple tuple;
1668  AttInMetadata *attinmeta;
1669  int flags = 0;
1670 
1671  /*
1672  * Build a tuple description for a
1673  * valid_detail tuple
1674  */
1675  get_call_result_type(fcinfo, 0, &tupdesc);
1676  BlessTupleDesc(tupdesc);
1677 
1678  /*
1679  * generate attribute metadata needed later to produce
1680  * tuples from raw C strings
1681  */
1682  attinmeta = TupleDescGetAttInMetadata(tupdesc);
1683 
1684  geom = PG_GETARG_GSERIALIZED_P(0);
1685  flags = PG_GETARG_INT32(1);
1686 
1687  initGEOS(lwpgnotice, lwgeom_geos_error);
1688 
1689  g1 = POSTGIS2GEOS(geom);
1690 
1691  if ( g1 )
1692  {
1693  valid = GEOSisValidDetail(g1, flags, &geos_reason, &geos_location);
1694  GEOSGeom_destroy((GEOSGeometry *)g1);
1695  if ( geos_reason )
1696  {
1697  reason = pstrdup(geos_reason);
1698  GEOSFree(geos_reason);
1699  }
1700  if ( geos_location )
1701  {
1702  location = GEOS2LWGEOM(geos_location, GEOSHasZ(geos_location));
1703  GEOSGeom_destroy(geos_location);
1704  }
1705 
1706  if (valid == 2)
1707  {
1708  /* NOTE: should only happen on OOM or similar */
1709  lwpgerror("GEOS isvaliddetail() threw an exception!");
1710  PG_RETURN_NULL(); /* never gets here */
1711  }
1712  }
1713  else
1714  {
1715  /* TODO: check lwgeom_geos_errmsg for validity error */
1716  reason = pstrdup(lwgeom_geos_errmsg);
1717  }
1718 
1719  /* the boolean validity */
1720  values[0] = valid ? "t" : "f";
1721 
1722  /* the reason */
1723  values[1] = reason;
1724 
1725  /* the location */
1726  values[2] = location ? lwgeom_to_hexwkb_buffer(location, WKB_EXTENDED) : 0;
1727 
1728  tuple = BuildTupleFromCStrings(attinmeta, values);
1729  result = (HeapTupleHeader) palloc(tuple->t_len);
1730  memcpy(result, tuple->t_data, tuple->t_len);
1731  heap_freetuple(tuple);
1732 
1733  PG_RETURN_HEAPTUPLEHEADER(result);
1734 }
1735 
1736 
1738 Datum issimple(PG_FUNCTION_ARGS)
1739 {
1740  GSERIALIZED *geom;
1741  LWGEOM *lwgeom_in;
1742  int result;
1743 
1744  POSTGIS_DEBUG(2, "issimple called");
1745 
1746  geom = PG_GETARG_GSERIALIZED_P(0);
1747 
1748  if ( gserialized_is_empty(geom) )
1749  PG_RETURN_BOOL(true);
1750 
1751  lwgeom_in = lwgeom_from_gserialized(geom);
1752  result = lwgeom_is_simple(lwgeom_in);
1753  lwgeom_free(lwgeom_in) ;
1754  PG_FREE_IF_COPY(geom, 0);
1755 
1756  if (result == -1) {
1757  PG_RETURN_NULL(); /*never get here */
1758  }
1759 
1760  PG_RETURN_BOOL(result);
1761 }
1762 
1764 Datum isring(PG_FUNCTION_ARGS)
1765 {
1766  GSERIALIZED *geom;
1767  GEOSGeometry *g1;
1768  int result;
1769 
1770  geom = PG_GETARG_GSERIALIZED_P(0);
1771 
1772  /* Empty things can't close */
1773  if ( gserialized_is_empty(geom) )
1774  PG_RETURN_BOOL(false);
1775 
1776  initGEOS(lwpgnotice, lwgeom_geos_error);
1777 
1778  g1 = POSTGIS2GEOS(geom);
1779  if (!g1)
1780  HANDLE_GEOS_ERROR("First argument geometry could not be converted to GEOS");
1781 
1782  if ( GEOSGeomTypeId(g1) != GEOS_LINESTRING )
1783  {
1784  GEOSGeom_destroy(g1);
1785  elog(ERROR, "ST_IsRing() should only be called on a linear feature");
1786  }
1787 
1788  result = GEOSisRing(g1);
1789  GEOSGeom_destroy(g1);
1790 
1791  if (result == 2) HANDLE_GEOS_ERROR("GEOSisRing");
1792 
1793  PG_FREE_IF_COPY(geom, 0);
1794  PG_RETURN_BOOL(result);
1795 }
1796 
1797 GSERIALIZED *
1798 GEOS2POSTGIS(GEOSGeom geom, char want3d)
1799 {
1800  LWGEOM *lwgeom;
1802 
1803  lwgeom = GEOS2LWGEOM(geom, want3d);
1804  if ( ! lwgeom )
1805  {
1806  lwpgerror("%s: GEOS2LWGEOM returned NULL", __func__);
1807  return NULL;
1808  }
1809 
1810  POSTGIS_DEBUGF(4, "%s: GEOS2LWGEOM returned a %s", __func__, lwgeom_summary(lwgeom, 0));
1811 
1812  if (lwgeom_needs_bbox(lwgeom)) lwgeom_add_bbox(lwgeom);
1813 
1814  result = geometry_serialize(lwgeom);
1815  lwgeom_free(lwgeom);
1816 
1817  return result;
1818 }
1819 
1820 /*-----=POSTGIS2GEOS= */
1821 
1822 GEOSGeometry *
1823 POSTGIS2GEOS(const GSERIALIZED *pglwgeom)
1824 {
1825  GEOSGeometry *ret;
1826  LWGEOM *lwgeom = lwgeom_from_gserialized(pglwgeom);
1827  if ( ! lwgeom )
1828  {
1829  lwpgerror("POSTGIS2GEOS: unable to deserialize input");
1830  return NULL;
1831  }
1832  ret = LWGEOM2GEOS(lwgeom, 0);
1833  lwgeom_free(lwgeom);
1834 
1835  return ret;
1836 }
1837 
1838 static uint32_t
1839 array_nelems_not_null(ArrayType* array)
1840 {
1841  ArrayIterator iterator;
1842  Datum value;
1843  bool isnull;
1844  uint32_t nelems_not_null = 0;
1845  iterator = array_create_iterator(array, 0, NULL);
1846  while(array_iterate(iterator, &value, &isnull))
1847  {
1848  if (!isnull)
1849  nelems_not_null++;
1850  }
1851 
1852  array_free_iterator(iterator);
1853 
1854  return nelems_not_null;
1855 }
1856 
1857 /* ARRAY2LWGEOM: Converts the non-null elements of a Postgres array into a LWGEOM* array */
1858 LWGEOM** ARRAY2LWGEOM(ArrayType* array, uint32_t nelems, int* is3d, int* srid)
1859 {
1860  ArrayIterator iterator;
1861  Datum value;
1862  bool isnull;
1863  bool gotsrid = false;
1864  uint32_t i = 0;
1865 
1866  LWGEOM** lw_geoms = palloc(nelems * sizeof(LWGEOM*));
1867 
1868  iterator = array_create_iterator(array, 0, NULL);
1869 
1870  while (array_iterate(iterator, &value, &isnull))
1871  {
1872  GSERIALIZED *geom = (GSERIALIZED *)DatumGetPointer(value);
1873 
1874  if (isnull)
1875  continue;
1876 
1877  *is3d = *is3d || gserialized_has_z(geom);
1878 
1879  lw_geoms[i] = lwgeom_from_gserialized(geom);
1880  if (!lw_geoms[i]) /* error in creation */
1881  {
1882  lwpgerror("Geometry deserializing geometry");
1883  return NULL;
1884  }
1885  if (!gotsrid)
1886  {
1887  gotsrid = true;
1888  *srid = gserialized_get_srid(geom);
1889  }
1890  else
1891  gserialized_error_if_srid_mismatch_reference(geom, *srid, __func__);
1892 
1893  i++;
1894  }
1895 
1896  return lw_geoms;
1897 }
1898 
1899 /* ARRAY2GEOS: Converts the non-null elements of a Postgres array into a GEOSGeometry* array */
1900 GEOSGeometry** ARRAY2GEOS(ArrayType* array, uint32_t nelems, int* is3d, int* srid)
1901 {
1902  ArrayIterator iterator;
1903  Datum value;
1904  bool isnull;
1905  bool gotsrid = false;
1906  uint32_t i = 0;
1907 
1908  GEOSGeometry** geos_geoms = palloc(nelems * sizeof(GEOSGeometry*));
1909 
1910  iterator = array_create_iterator(array, 0, NULL);
1911 
1912  while(array_iterate(iterator, &value, &isnull))
1913  {
1914  GSERIALIZED *geom = (GSERIALIZED*) DatumGetPointer(value);
1915 
1916  if (isnull)
1917  continue;
1918 
1919  *is3d = *is3d || gserialized_has_z(geom);
1920 
1921  geos_geoms[i] = POSTGIS2GEOS(geom);
1922  if (!geos_geoms[i])
1923  {
1924  uint32_t j;
1925  lwpgerror("Geometry could not be converted to GEOS");
1926 
1927  for (j = 0; j < i; j++) {
1928  GEOSGeom_destroy(geos_geoms[j]);
1929  }
1930  return NULL;
1931  }
1932 
1933  if (!gotsrid)
1934  {
1935  *srid = gserialized_get_srid(geom);
1936  gotsrid = true;
1937  }
1938  else if (*srid != gserialized_get_srid(geom))
1939  {
1940  uint32_t j;
1941  for (j = 0; j <= i; j++) {
1942  GEOSGeom_destroy(geos_geoms[j]);
1943  }
1944  gserialized_error_if_srid_mismatch_reference(geom, *srid, __func__);
1945  return NULL;
1946  }
1947 
1948  i++;
1949  }
1950 
1951  array_free_iterator(iterator);
1952  return geos_geoms;
1953 }
1954 
1956 Datum GEOSnoop(PG_FUNCTION_ARGS)
1957 {
1958  GSERIALIZED *geom;
1959  GEOSGeometry *geosgeom;
1960  GSERIALIZED *lwgeom_result;
1961 
1962  initGEOS(lwpgnotice, lwgeom_geos_error);
1963 
1964  geom = PG_GETARG_GSERIALIZED_P(0);
1965  geosgeom = POSTGIS2GEOS(geom);
1966  if ( ! geosgeom ) PG_RETURN_NULL();
1967 
1968  lwgeom_result = GEOS2POSTGIS(geosgeom, gserialized_has_z(geom));
1969  GEOSGeom_destroy(geosgeom);
1970 
1971  PG_FREE_IF_COPY(geom, 0);
1972 
1973  PG_RETURN_POINTER(lwgeom_result);
1974 }
1975 
1977 Datum polygonize_garray(PG_FUNCTION_ARGS)
1978 {
1979  ArrayType *array;
1980  int is3d = 0;
1981  uint32 nelems, i;
1983  GEOSGeometry *geos_result;
1984  const GEOSGeometry **vgeoms;
1985  int32_t srid = SRID_UNKNOWN;
1986 
1987  if (PG_ARGISNULL(0))
1988  PG_RETURN_NULL();
1989 
1990  array = PG_GETARG_ARRAYTYPE_P(0);
1991  nelems = array_nelems_not_null(array);
1992 
1993  if (nelems == 0)
1994  PG_RETURN_NULL();
1995 
1996  POSTGIS_DEBUGF(3, "polygonize_garray: number of non-null elements: %d", nelems);
1997 
1998  /* Ok, we really need geos now ;) */
1999  initGEOS(lwpgnotice, lwgeom_geos_error);
2000 
2001  vgeoms = (const GEOSGeometry**) ARRAY2GEOS(array, nelems, &is3d, &srid);
2002 
2003  POSTGIS_DEBUG(3, "polygonize_garray: invoking GEOSpolygonize");
2004 
2005  geos_result = GEOSPolygonize(vgeoms, nelems);
2006 
2007  POSTGIS_DEBUG(3, "polygonize_garray: GEOSpolygonize returned");
2008 
2009  for (i=0; i<nelems; ++i) GEOSGeom_destroy((GEOSGeometry *)vgeoms[i]);
2010  pfree(vgeoms);
2011 
2012  if ( ! geos_result ) PG_RETURN_NULL();
2013 
2014  GEOSSetSRID(geos_result, srid);
2015  result = GEOS2POSTGIS(geos_result, is3d);
2016  GEOSGeom_destroy(geos_result);
2017  if (!result)
2018  {
2019  elog(ERROR, "%s returned an error", __func__);
2020  PG_RETURN_NULL(); /*never get here */
2021  }
2022 
2023  PG_RETURN_POINTER(result);
2024 }
2025 
2026 
2028 Datum clusterintersecting_garray(PG_FUNCTION_ARGS)
2029 {
2030  Datum* result_array_data;
2031  ArrayType *array, *result;
2032  int is3d = 0;
2033  uint32 nelems, nclusters, i;
2034  GEOSGeometry **geos_inputs, **geos_results;
2035  int32_t srid = SRID_UNKNOWN;
2036 
2037  /* Parameters used to construct a result array */
2038  int16 elmlen;
2039  bool elmbyval;
2040  char elmalign;
2041 
2042  /* Null array, null geometry (should be empty?) */
2043  if (PG_ARGISNULL(0))
2044  PG_RETURN_NULL();
2045 
2046  array = PG_GETARG_ARRAYTYPE_P(0);
2047  nelems = array_nelems_not_null(array);
2048 
2049  POSTGIS_DEBUGF(3, "clusterintersecting_garray: number of non-null elements: %d", nelems);
2050 
2051  if ( nelems == 0 ) PG_RETURN_NULL();
2052 
2053  /* TODO short-circuit for one element? */
2054 
2055  /* Ok, we really need geos now ;) */
2056  initGEOS(lwpgnotice, lwgeom_geos_error);
2057 
2058  geos_inputs = ARRAY2GEOS(array, nelems, &is3d, &srid);
2059  if(!geos_inputs)
2060  {
2061  PG_RETURN_NULL();
2062  }
2063 
2064  if (cluster_intersecting(geos_inputs, nelems, &geos_results, &nclusters) != LW_SUCCESS)
2065  {
2066  elog(ERROR, "clusterintersecting: Error performing clustering");
2067  PG_RETURN_NULL();
2068  }
2069  pfree(geos_inputs); /* don't need to destroy items because GeometryCollections have taken ownership */
2070 
2071  if (!geos_results) PG_RETURN_NULL();
2072 
2073  result_array_data = palloc(nclusters * sizeof(Datum));
2074  for (i=0; i<nclusters; ++i)
2075  {
2076  result_array_data[i] = PointerGetDatum(GEOS2POSTGIS(geos_results[i], is3d));
2077  GEOSGeom_destroy(geos_results[i]);
2078  }
2079  lwfree(geos_results);
2080 
2081  get_typlenbyvalalign(array->elemtype, &elmlen, &elmbyval, &elmalign);
2082  result = construct_array(result_array_data, nclusters, array->elemtype, elmlen, elmbyval, elmalign);
2083 
2084  if (!result)
2085  {
2086  elog(ERROR, "clusterintersecting: Error constructing return-array");
2087  PG_RETURN_NULL();
2088  }
2089 
2090  PG_RETURN_POINTER(result);
2091 }
2092 
2094 Datum cluster_within_distance_garray(PG_FUNCTION_ARGS)
2095 {
2096  Datum* result_array_data;
2097  ArrayType *array, *result;
2098  int is3d = 0;
2099  uint32 nelems, nclusters, i;
2100  LWGEOM** lw_inputs;
2101  LWGEOM** lw_results;
2102  double tolerance;
2103  int32_t srid = SRID_UNKNOWN;
2104 
2105  /* Parameters used to construct a result array */
2106  int16 elmlen;
2107  bool elmbyval;
2108  char elmalign;
2109 
2110  /* Null array, null geometry (should be empty?) */
2111  if (PG_ARGISNULL(0))
2112  PG_RETURN_NULL();
2113 
2114  array = PG_GETARG_ARRAYTYPE_P(0);
2115 
2116  tolerance = PG_GETARG_FLOAT8(1);
2117  if (tolerance < 0)
2118  {
2119  lwpgerror("Tolerance must be a positive number.");
2120  PG_RETURN_NULL();
2121  }
2122 
2123  nelems = array_nelems_not_null(array);
2124 
2125  POSTGIS_DEBUGF(3, "cluster_within_distance_garray: number of non-null elements: %d", nelems);
2126 
2127  if ( nelems == 0 ) PG_RETURN_NULL();
2128 
2129  /* TODO short-circuit for one element? */
2130 
2131  /* Ok, we really need geos now ;) */
2132  initGEOS(lwpgnotice, lwgeom_geos_error);
2133 
2134  lw_inputs = ARRAY2LWGEOM(array, nelems, &is3d, &srid);
2135  if (!lw_inputs)
2136  {
2137  PG_RETURN_NULL();
2138  }
2139 
2140  if (cluster_within_distance(lw_inputs, nelems, tolerance, &lw_results, &nclusters) != LW_SUCCESS)
2141  {
2142  elog(ERROR, "cluster_within: Error performing clustering");
2143  PG_RETURN_NULL();
2144  }
2145  pfree(lw_inputs); /* don't need to destroy items because GeometryCollections have taken ownership */
2146 
2147  if (!lw_results) PG_RETURN_NULL();
2148 
2149  result_array_data = palloc(nclusters * sizeof(Datum));
2150  for (i=0; i<nclusters; ++i)
2151  {
2152  result_array_data[i] = PointerGetDatum(geometry_serialize(lw_results[i]));
2153  lwgeom_free(lw_results[i]);
2154  }
2155  lwfree(lw_results);
2156 
2157  get_typlenbyvalalign(array->elemtype, &elmlen, &elmbyval, &elmalign);
2158  result = construct_array(result_array_data, nclusters, array->elemtype, elmlen, elmbyval, elmalign);
2159 
2160  if (!result)
2161  {
2162  elog(ERROR, "clusterwithin: Error constructing return-array");
2163  PG_RETURN_NULL();
2164  }
2165 
2166  PG_RETURN_POINTER(result);
2167 }
2168 
2170 Datum linemerge(PG_FUNCTION_ARGS)
2171 {
2172  GSERIALIZED *geom1;
2173  bool directed = false;
2175  LWGEOM *lwgeom1, *lwresult ;
2176 
2177  geom1 = PG_GETARG_GSERIALIZED_P(0);
2178 
2179  if ( PG_NARGS() > 1 )
2180  directed = PG_GETARG_BOOL(1);
2181 
2182  lwgeom1 = lwgeom_from_gserialized(geom1) ;
2183 
2184  lwresult = lwgeom_linemerge_directed(lwgeom1, directed ? LW_TRUE : LW_FALSE);
2185  result = geometry_serialize(lwresult) ;
2186 
2187  lwgeom_free(lwgeom1) ;
2188  lwgeom_free(lwresult) ;
2189 
2190  PG_FREE_IF_COPY(geom1, 0);
2191 
2192  PG_RETURN_POINTER(result);
2193 }
2194 
2195 /*
2196  * Take a geometry and return an areal geometry
2197  * (Polygon or MultiPolygon).
2198  * Actually a wrapper around GEOSpolygonize,
2199  * transforming the resulting collection into
2200  * a valid polygon Geometry.
2201  */
2203 Datum ST_BuildArea(PG_FUNCTION_ARGS)
2204 {
2206  GSERIALIZED *geom;
2207  LWGEOM *lwgeom_in, *lwgeom_out;
2208 
2209  geom = PG_GETARG_GSERIALIZED_P(0);
2210  lwgeom_in = lwgeom_from_gserialized(geom);
2211 
2212  lwgeom_out = lwgeom_buildarea(lwgeom_in);
2213  lwgeom_free(lwgeom_in) ;
2214 
2215  if ( ! lwgeom_out )
2216  {
2217  PG_FREE_IF_COPY(geom, 0);
2218  PG_RETURN_NULL();
2219  }
2220 
2221  result = geometry_serialize(lwgeom_out) ;
2222  lwgeom_free(lwgeom_out) ;
2223 
2224  PG_FREE_IF_COPY(geom, 0);
2225  PG_RETURN_POINTER(result);
2226 }
2227 
2228 /*
2229  * Take the vertices of a geometry and builds
2230  * Delaunay triangles around them.
2231  */
2233 Datum ST_DelaunayTriangles(PG_FUNCTION_ARGS)
2234 {
2236  GSERIALIZED *geom;
2237  LWGEOM *lwgeom_in, *lwgeom_out;
2238  double tolerance = 0.0;
2239  int flags = 0;
2240 
2241  geom = PG_GETARG_GSERIALIZED_P(0);
2242  tolerance = PG_GETARG_FLOAT8(1);
2243  flags = PG_GETARG_INT32(2);
2244 
2245  lwgeom_in = lwgeom_from_gserialized(geom);
2246  lwgeom_out = lwgeom_delaunay_triangulation(lwgeom_in, tolerance, flags);
2247  lwgeom_free(lwgeom_in) ;
2248 
2249  if ( ! lwgeom_out )
2250  {
2251  PG_FREE_IF_COPY(geom, 0);
2252  PG_RETURN_NULL();
2253  }
2254 
2255  result = geometry_serialize(lwgeom_out) ;
2256  lwgeom_free(lwgeom_out) ;
2257 
2258  PG_FREE_IF_COPY(geom, 0);
2259  PG_RETURN_POINTER(result);
2260 }
2261 
2262 /*
2263  * Take a polygon and build a constrained
2264  * triangulation that respect the edges of the
2265  * polygon.
2266  */
2268 Datum ST_TriangulatePolygon(PG_FUNCTION_ARGS)
2269 {
2270 #if POSTGIS_GEOS_VERSION < 31100
2271 
2272  lwpgerror("The GEOS version this PostGIS binary "
2273  "was compiled against (%d) doesn't support "
2274  "'GEOSConstrainedDelaunayTriangulation' function (3.11.0+ required)",
2276  PG_RETURN_NULL();
2277 
2278 #else /* POSTGIS_GEOS_VERSION >= 31100 */
2280  GSERIALIZED *geom = PG_GETARG_GSERIALIZED_P(0);
2281  LWGEOM *lwgeom_in = lwgeom_from_gserialized(geom);
2282  LWGEOM *lwgeom_out = lwgeom_triangulate_polygon(lwgeom_in);
2283  lwgeom_free(lwgeom_in);
2284 
2285  if (!lwgeom_out)
2286  {
2287  PG_FREE_IF_COPY(geom, 0);
2288  PG_RETURN_NULL();
2289  }
2290 
2291  result = geometry_serialize(lwgeom_out);
2292  lwgeom_free(lwgeom_out);
2293 
2294  PG_FREE_IF_COPY(geom, 0);
2295  PG_RETURN_POINTER(result);
2296 #endif
2297 }
2298 
2299 /*
2300  * ST_Snap
2301  *
2302  * Snap a geometry to another with a given tolerance
2303  */
2304 Datum ST_Snap(PG_FUNCTION_ARGS);
2306 Datum ST_Snap(PG_FUNCTION_ARGS)
2307 {
2308  GSERIALIZED *geom1, *geom2, *result;
2309  LWGEOM *lwgeom1, *lwgeom2, *lwresult;
2310  double tolerance;
2311 
2312  geom1 = PG_GETARG_GSERIALIZED_P(0);
2313  geom2 = PG_GETARG_GSERIALIZED_P(1);
2314  tolerance = PG_GETARG_FLOAT8(2);
2315 
2316  lwgeom1 = lwgeom_from_gserialized(geom1);
2317  lwgeom2 = lwgeom_from_gserialized(geom2);
2318 
2319  lwresult = lwgeom_snap(lwgeom1, lwgeom2, tolerance);
2320  lwgeom_free(lwgeom1);
2321  lwgeom_free(lwgeom2);
2322  PG_FREE_IF_COPY(geom1, 0);
2323  PG_FREE_IF_COPY(geom2, 1);
2324 
2325  result = geometry_serialize(lwresult);
2326  lwgeom_free(lwresult);
2327 
2328  PG_RETURN_POINTER(result);
2329 }
2330 
2331 /*
2332  * ST_Split
2333  *
2334  * Split polygon by line, line by line, line by point.
2335  * Returns at most components as a collection.
2336  * First element of the collection is always the part which
2337  * remains after the cut, while the second element is the
2338  * part which has been cut out. We arbitrarily take the part
2339  * on the *right* of cut lines as the part which has been cut out.
2340  * For a line cut by a point the part which remains is the one
2341  * from start of the line to the cut point.
2342  *
2343  *
2344  * Author: Sandro Santilli <strk@kbt.io>
2345  *
2346  * Work done for Faunalia (http://www.faunalia.it) with fundings
2347  * from Regione Toscana - Sistema Informativo per il Governo
2348  * del Territorio e dell'Ambiente (RT-SIGTA).
2349  *
2350  * Thanks to the PostGIS community for sharing poly/line ideas [1]
2351  *
2352  * [1] http://trac.osgeo.org/postgis/wiki/UsersWikiSplitPolygonWithLineString
2353  *
2354  */
2355 Datum ST_Split(PG_FUNCTION_ARGS);
2357 Datum ST_Split(PG_FUNCTION_ARGS)
2358 {
2359  GSERIALIZED *in, *blade_in, *out;
2360  LWGEOM *lwgeom_in, *lwblade_in, *lwgeom_out;
2361 
2362  in = PG_GETARG_GSERIALIZED_P(0);
2363  blade_in = PG_GETARG_GSERIALIZED_P(1);
2364  gserialized_error_if_srid_mismatch(in, blade_in, __func__);
2365 
2366  lwgeom_in = lwgeom_from_gserialized(in);
2367  lwblade_in = lwgeom_from_gserialized(blade_in);
2368 
2369  if (!lwgeom_isfinite(lwgeom_in))
2370  {
2371  lwpgerror("Input Geometry contains invalid coordinates");
2372  PG_RETURN_NULL();
2373  }
2374 
2375  if (!lwgeom_isfinite(lwblade_in))
2376  {
2377  lwpgerror("Blade Geometry contains invalid coordinates");
2378  PG_RETURN_NULL();
2379  }
2380 
2381 
2382  lwgeom_out = lwgeom_split(lwgeom_in, lwblade_in);
2383  lwgeom_free(lwgeom_in);
2384  lwgeom_free(lwblade_in);
2385 
2386  if ( ! lwgeom_out )
2387  {
2388  PG_FREE_IF_COPY(in, 0); /* possibly referenced by lwgeom_out */
2389  PG_FREE_IF_COPY(blade_in, 1);
2390  PG_RETURN_NULL();
2391  }
2392 
2393  out = geometry_serialize(lwgeom_out);
2394  lwgeom_free(lwgeom_out);
2395  PG_FREE_IF_COPY(in, 0); /* possibly referenced by lwgeom_out */
2396  PG_FREE_IF_COPY(blade_in, 1);
2397 
2398  PG_RETURN_POINTER(out);
2399 }
2400 
2401 /**********************************************************************
2402  *
2403  * ST_SharedPaths
2404  *
2405  * Return the set of paths shared between two linear geometries,
2406  * and their direction (same or opposite).
2407  *
2408  * Developed by Sandro Santilli (strk@kbt.io) for Faunalia
2409  * (http://www.faunalia.it) with funding from Regione Toscana - Sistema
2410  * Informativo per la Gestione del Territorio e dell' Ambiente
2411  * [RT-SIGTA]". For the project: "Sviluppo strumenti software per il
2412  * trattamento di dati geografici basati su QuantumGIS e Postgis (CIG
2413  * 0494241492)"
2414  *
2415  **********************************************************************/
2416 Datum ST_SharedPaths(PG_FUNCTION_ARGS);
2418 Datum ST_SharedPaths(PG_FUNCTION_ARGS)
2419 {
2420  GSERIALIZED *geom1, *geom2, *out;
2421  LWGEOM *g1, *g2, *lwgeom_out;
2422 
2423  geom1 = PG_GETARG_GSERIALIZED_P(0);
2424  geom2 = PG_GETARG_GSERIALIZED_P(1);
2425 
2426  g1 = lwgeom_from_gserialized(geom1);
2427  g2 = lwgeom_from_gserialized(geom2);
2428 
2429  lwgeom_out = lwgeom_sharedpaths(g1, g2);
2430  lwgeom_free(g1);
2431  lwgeom_free(g2);
2432 
2433  if ( ! lwgeom_out )
2434  {
2435  PG_FREE_IF_COPY(geom1, 0);
2436  PG_FREE_IF_COPY(geom2, 1);
2437  PG_RETURN_NULL();
2438  }
2439 
2440  out = geometry_serialize(lwgeom_out);
2441  lwgeom_free(lwgeom_out);
2442 
2443  PG_FREE_IF_COPY(geom1, 0);
2444  PG_FREE_IF_COPY(geom2, 1);
2445  PG_RETURN_POINTER(out);
2446 }
2447 
2448 
2449 /**********************************************************************
2450  *
2451  * ST_Node
2452  *
2453  * Fully node a set of lines using the least possible nodes while
2454  * preserving all of the input ones.
2455  *
2456  **********************************************************************/
2457 Datum ST_Node(PG_FUNCTION_ARGS);
2459 Datum ST_Node(PG_FUNCTION_ARGS)
2460 {
2461  GSERIALIZED *geom1, *out;
2462  LWGEOM *g1, *lwgeom_out;
2463 
2464  geom1 = PG_GETARG_GSERIALIZED_P(0);
2465 
2466  g1 = lwgeom_from_gserialized(geom1);
2467 
2468  lwgeom_out = lwgeom_node(g1);
2469  lwgeom_free(g1);
2470 
2471  if ( ! lwgeom_out )
2472  {
2473  PG_FREE_IF_COPY(geom1, 0);
2474  PG_RETURN_NULL();
2475  }
2476 
2477  out = geometry_serialize(lwgeom_out);
2478  lwgeom_free(lwgeom_out);
2479 
2480  PG_FREE_IF_COPY(geom1, 0);
2481  PG_RETURN_POINTER(out);
2482 }
2483 
2484 /******************************************
2485  *
2486  * ST_Voronoi
2487  *
2488  * Returns a Voronoi diagram constructed
2489  * from the points of the input geometry.
2490  *
2491  ******************************************/
2492 Datum ST_Voronoi(PG_FUNCTION_ARGS);
2494 Datum ST_Voronoi(PG_FUNCTION_ARGS)
2495 {
2496  GSERIALIZED* input;
2497  GSERIALIZED* clip;
2499  LWGEOM* lwgeom_input;
2500  LWGEOM* lwgeom_result;
2501  double tolerance;
2502  GBOX clip_envelope;
2503  int custom_clip_envelope;
2504  int return_polygons;
2505 
2506  /* Return NULL on NULL geometry */
2507  if (PG_ARGISNULL(0))
2508  PG_RETURN_NULL();
2509 
2510  /* Read our tolerance value */
2511  if (PG_ARGISNULL(2))
2512  {
2513  lwpgerror("Tolerance must be a positive number.");
2514  PG_RETURN_NULL();
2515  }
2516 
2517  tolerance = PG_GETARG_FLOAT8(2);
2518 
2519  if (tolerance < 0)
2520  {
2521  lwpgerror("Tolerance must be a positive number.");
2522  PG_RETURN_NULL();
2523  }
2524 
2525  /* Are we returning lines or polygons? */
2526  if (PG_ARGISNULL(3))
2527  {
2528  lwpgerror("return_polygons must be true or false.");
2529  PG_RETURN_NULL();
2530  }
2531  return_polygons = PG_GETARG_BOOL(3);
2532 
2533  /* Read our clipping envelope, if applicable. */
2534  custom_clip_envelope = !PG_ARGISNULL(1);
2535  if (custom_clip_envelope) {
2536  clip = PG_GETARG_GSERIALIZED_P(1);
2537  if (!gserialized_get_gbox_p(clip, &clip_envelope))
2538  {
2539  lwpgerror("Could not determine envelope of clipping geometry.");
2540  PG_FREE_IF_COPY(clip, 1);
2541  PG_RETURN_NULL();
2542  }
2543  PG_FREE_IF_COPY(clip, 1);
2544  }
2545 
2546  /* Read our input geometry */
2547  input = PG_GETARG_GSERIALIZED_P(0);
2548 
2549  lwgeom_input = lwgeom_from_gserialized(input);
2550 
2551  if(!lwgeom_input)
2552  {
2553  lwpgerror("Could not read input geometry.");
2554  PG_FREE_IF_COPY(input, 0);
2555  PG_RETURN_NULL();
2556  }
2557 
2558  lwgeom_result = lwgeom_voronoi_diagram(lwgeom_input, custom_clip_envelope ? &clip_envelope : NULL, tolerance, !return_polygons);
2559  lwgeom_free(lwgeom_input);
2560 
2561  if (!lwgeom_result)
2562  {
2563  lwpgerror("Error computing Voronoi diagram.");
2564  PG_FREE_IF_COPY(input, 0);
2565  PG_RETURN_NULL();
2566  }
2567 
2568  result = geometry_serialize(lwgeom_result);
2569  lwgeom_free(lwgeom_result);
2570 
2571  PG_FREE_IF_COPY(input, 0);
2572  PG_RETURN_POINTER(result);
2573 }
2574 
2575 /******************************************
2576  *
2577  * ST_MinimumClearance
2578  *
2579  * Returns the minimum clearance of a geometry.
2580  *
2581  ******************************************/
2582 Datum ST_MinimumClearance(PG_FUNCTION_ARGS);
2584 Datum ST_MinimumClearance(PG_FUNCTION_ARGS)
2585 {
2586  GSERIALIZED* input;
2587  GEOSGeometry* input_geos;
2588  int error;
2589  double result;
2590 
2591  initGEOS(lwpgnotice, lwgeom_geos_error);
2592 
2593  input = PG_GETARG_GSERIALIZED_P(0);
2594  input_geos = POSTGIS2GEOS(input);
2595  if (!input_geos)
2596  HANDLE_GEOS_ERROR("Geometry could not be converted to GEOS");
2597 
2598  error = GEOSMinimumClearance(input_geos, &result);
2599  GEOSGeom_destroy(input_geos);
2600  if (error) HANDLE_GEOS_ERROR("Error computing minimum clearance");
2601 
2602  PG_FREE_IF_COPY(input, 0);
2603  PG_RETURN_FLOAT8(result);
2604 }
2605 
2606 /******************************************
2607  *
2608  * ST_MinimumClearanceLine
2609  *
2610  * Returns the minimum clearance line of a geometry.
2611  *
2612  ******************************************/
2613 Datum ST_MinimumClearanceLine(PG_FUNCTION_ARGS);
2615 Datum ST_MinimumClearanceLine(PG_FUNCTION_ARGS)
2616 {
2617  GSERIALIZED* input;
2619  GEOSGeometry* input_geos;
2620  GEOSGeometry* result_geos;
2621  int32_t srid;
2622 
2623  initGEOS(lwpgnotice, lwgeom_geos_error);
2624 
2625  input = PG_GETARG_GSERIALIZED_P(0);
2626  srid = gserialized_get_srid(input);
2627  input_geos = POSTGIS2GEOS(input);
2628  if (!input_geos)
2629  HANDLE_GEOS_ERROR("Geometry could not be converted to GEOS");
2630 
2631  result_geos = GEOSMinimumClearanceLine(input_geos);
2632  GEOSGeom_destroy(input_geos);
2633  if (!result_geos)
2634  HANDLE_GEOS_ERROR("Error computing minimum clearance");
2635 
2636  GEOSSetSRID(result_geos, srid);
2637  result = GEOS2POSTGIS(result_geos, LW_FALSE);
2638  GEOSGeom_destroy(result_geos);
2639 
2640  PG_FREE_IF_COPY(input, 0);
2641  PG_RETURN_POINTER(result);
2642 }
2643 
2644 /******************************************
2645  *
2646  * ST_OrientedEnvelope
2647  *
2648  ******************************************/
2649 Datum ST_OrientedEnvelope(PG_FUNCTION_ARGS);
2651 Datum ST_OrientedEnvelope(PG_FUNCTION_ARGS)
2652 {
2653  GSERIALIZED* input;
2655  GEOSGeometry* input_geos;
2656  GEOSGeometry* result_geos;
2657  int32_t srid;
2658 
2659  initGEOS(lwpgnotice, lwgeom_geos_error);
2660 
2661  input = PG_GETARG_GSERIALIZED_P(0);
2662  srid = gserialized_get_srid(input);
2663  input_geos = POSTGIS2GEOS(input);
2664  if (!input_geos)
2665  HANDLE_GEOS_ERROR("Geometry could not be converted to GEOS");
2666 
2667  result_geos = GEOSMinimumRotatedRectangle(input_geos);
2668  GEOSGeom_destroy(input_geos);
2669  if (!result_geos)
2670  HANDLE_GEOS_ERROR("Error computing oriented envelope");
2671 
2672  GEOSSetSRID(result_geos, srid);
2673  result = GEOS2POSTGIS(result_geos, LW_FALSE);
2674  GEOSGeom_destroy(result_geos);
2675 
2676  PG_FREE_IF_COPY(input, 0);
2677  PG_RETURN_POINTER(result);
2678 }
2679 
2680 
2686 Datum LWGEOM_dfullywithin(PG_FUNCTION_ARGS);
2688 Datum LWGEOM_dfullywithin(PG_FUNCTION_ARGS)
2689 {
2690  GSERIALIZED *geom1 = PG_GETARG_GSERIALIZED_P(0);
2691  GSERIALIZED *geom2 = PG_GETARG_GSERIALIZED_P(1);
2692  LWGEOM *lwg1 = lwgeom_from_gserialized(geom1);
2693  LWGEOM *lwg2 = lwgeom_from_gserialized(geom2);
2694  double radius = PG_GETARG_FLOAT8(2);
2695  GEOSGeometry *buffer1 = NULL;
2696  GEOSGeometry *geos1 = NULL, *geos2 = NULL;
2697  char contained;
2698 
2699  if (radius < 0.0)
2700  {
2701  elog(ERROR, "Tolerance cannot be less than zero\n");
2702  PG_RETURN_NULL();
2703  }
2704 
2705  gserialized_error_if_srid_mismatch(geom1, geom2, __func__);
2706 
2707  if (lwgeom_is_empty(lwg1) || lwgeom_is_empty(lwg2))
2708  PG_RETURN_BOOL(false);
2709 
2710  if (!(lwgeom_isfinite(lwg1) && lwgeom_isfinite(lwg2)))
2711  PG_RETURN_BOOL(false);
2712 
2713  initGEOS(lwpgnotice, lwgeom_geos_error);
2714 
2715  geos1 = LWGEOM2GEOS(lwg1, true);
2716  geos2 = LWGEOM2GEOS(lwg2, true);
2717  lwgeom_free(lwg1);
2718  lwgeom_free(lwg2);
2719  if (!(geos1 && geos2))
2720  HANDLE_GEOS_ERROR("Geometry could not be converted to GEOS");
2721 
2722  buffer1 = GEOSBuffer(geos1, radius, 16);
2723  GEOSGeom_destroy(geos1);
2724  if (!(buffer1))
2725  HANDLE_GEOS_ERROR("Buffer operation failed");
2726 
2727  contained = GEOSCovers(buffer1, geos2);
2728  GEOSGeom_destroy(buffer1);
2729  GEOSGeom_destroy(geos2);
2730 
2731  if (contained == 2) HANDLE_GEOS_ERROR("GEOSContains");
2732 
2733  PG_FREE_IF_COPY(geom1, 0);
2734  PG_FREE_IF_COPY(geom2, 1);
2735 
2736  PG_RETURN_BOOL(contained == 1);
2737 }
char result[OUT_DOUBLE_BUFFER_SIZE]
Definition: cu_print.c:267
int gbox_overlaps_2d(const GBOX *g1, const GBOX *g2)
Return LW_TRUE if the GBOX overlaps on the 2d plane, LW_FALSE otherwise.
Definition: gbox.c:323
GBOX * gbox_copy(const GBOX *box)
Return a copy of the GBOX, based on dimensionality of flags.
Definition: gbox.c:438
void gserialized_error_if_srid_mismatch(const GSERIALIZED *g1, const GSERIALIZED *g2, const char *funcname)
Definition: gserialized.c:432
void gserialized_error_if_srid_mismatch_reference(const GSERIALIZED *g1, const int32_t srid2, const char *funcname)
Definition: gserialized.c:447
int32_t gserialized_get_srid(const GSERIALIZED *g)
Extract the SRID from the serialized form (it is packed into three bytes so this is a handy function)...
Definition: gserialized.c:155
int gserialized_get_gbox_p(const GSERIALIZED *g, GBOX *gbox)
Read the box from the GSERIALIZED or calculate it if necessary.
Definition: gserialized.c:94
LWGEOM * lwgeom_from_gserialized(const GSERIALIZED *g)
Allocate a new LWGEOM from a GSERIALIZED.
Definition: gserialized.c:268
int gserialized_is_empty(const GSERIALIZED *g)
Check if a GSERIALIZED is empty without deserializing first.
Definition: gserialized.c:181
int gserialized_has_z(const GSERIALIZED *g)
Check if a GSERIALIZED has a Z ordinate.
Definition: gserialized.c:203
uint32_t gserialized_get_type(const GSERIALIZED *g)
Extract the geometry type from the serialized form (it hides in the anonymous data area,...
Definition: gserialized.c:118
int gserialized_datum_get_internals_p(Datum gsdatum, GBOX *gbox, lwflags_t *flags, uint8_t *type, int32_t *srid)
Peak into a GSERIALIZED datum to find its bounding box and some other metadata.
char lwgeom_geos_errmsg[LWGEOM_GEOS_ERRMSG_MAXSIZE]
GEOSGeometry * LWGEOM2GEOS(const LWGEOM *lwgeom, uint8_t autofix)
void lwgeom_geos_error(const char *fmt,...)
void(*) LWGEOM GEOS2LWGEOM)(const GEOSGeometry *geom, uint8_t want3d)
int cluster_within_distance(LWGEOM **geoms, uint32_t num_geoms, double tolerance, LWGEOM ***clusterGeoms, uint32_t *num_clusters)
Takes an array of LWGEOM* and constructs an array of LWGEOM*, where each element in the constructed a...
int cluster_intersecting(GEOSGeometry **geoms, uint32_t num_geoms, GEOSGeometry ***clusterGeoms, uint32_t *num_clusters)
Takes an array of GEOSGeometry* and constructs an array of GEOSGeometry*, where each element in the c...
LWPOINT * lwpoint_construct_empty(int32_t srid, char hasz, char hasm)
Definition: lwpoint.c:151
LWGEOM * lwgeom_centroid(const LWGEOM *geom)
#define LW_FALSE
Definition: liblwgeom.h:94
LWGEOM * lwgeom_delaunay_triangulation(const LWGEOM *geom, double tolerance, int32_t edgeOnly)
Take vertices of a geometry and build a delaunay triangulation on them.
LWGEOM * lwgeom_node(const LWGEOM *lwgeom_in)
void lwgeom_free(LWGEOM *geom)
Definition: lwgeom.c:1218
LWGEOM * lwgeom_concavehull(const LWGEOM *geom, double ratio, uint32_t allow_holes)
Take a geometry and build the concave hull.
LWGEOM * lwgeom_split(const LWGEOM *lwgeom_in, const LWGEOM *blade_in)
LWGEOM * lwgeom_offsetcurve(const LWGEOM *geom, double size, int quadsegs, int joinStyle, double mitreLimit)
#define LW_SUCCESS
Definition: liblwgeom.h:97
uint16_t lwflags_t
Definition: liblwgeom.h:299
LWGEOM * lwpoly_as_lwgeom(const LWPOLY *obj)
Definition: lwgeom.c:329
int lwgeom_is_simple(const LWGEOM *lwgeom)
LWGEOM * lwgeom_snap(const LWGEOM *geom1, const LWGEOM *geom2, double tolerance)
Snap vertices and segments of a geometry to another using a given tolerance.
char * lwgeom_to_hexwkb_buffer(const LWGEOM *geom, uint8_t variant)
Definition: lwout_wkb.c:845
const char * lwgeom_geos_version(void)
Return GEOS version string (not to be freed)
int lwgeom_needs_bbox(const LWGEOM *geom)
Check whether or not a lwgeom is big enough to warrant a bounding box.
Definition: lwgeom.c:1271
int lwgeom_isfinite(const LWGEOM *lwgeom)
Check if a LWGEOM has any non-finite (NaN or Inf) coordinates.
Definition: lwgeom.c:2807
char * lwgeom_summary(const LWGEOM *lwgeom, int offset)
Definition: lwgeom_debug.c:166
LWGEOM * lwgeom_pointonsurface(const LWGEOM *geom)
LWGEOM * lwgeom_sharedpaths(const LWGEOM *geom1, const LWGEOM *geom2)
#define TINTYPE
Definition: liblwgeom.h:116
LWGEOM * lwgeom_simplify_polygonal(const LWGEOM *geom, double vertex_fraction, uint32_t is_outer)
Computes a boundary-respecting hull of a polygonal geometry, with hull shape determined by a target p...
LWGEOM * lwgeom_voronoi_diagram(const LWGEOM *g, const GBOX *env, double tolerance, int output_edges)
Take vertices of a geometry and build the Voronoi diagram.
#define MULTIPOLYGONTYPE
Definition: liblwgeom.h:107
LWGEOM * lwgeom_difference_prec(const LWGEOM *geom1, const LWGEOM *geom2, double gridSize)
void lwfree(void *mem)
Definition: lwutil.c:248
LWGEOM * lwgeom_union_prec(const LWGEOM *geom1, const LWGEOM *geom2, double gridSize)
#define POLYGONTYPE
Definition: liblwgeom.h:104
LWGEOM * lwgeom_unaryunion_prec(const LWGEOM *geom1, double gridSize)
const char * lwgeom_geos_compiled_version(void)
LWGEOM * lwgeom_buildarea(const LWGEOM *geom)
Take a geometry and return an areal geometry (Polygon or MultiPolygon).
#define WKB_EXTENDED
Definition: liblwgeom.h:2209
LWGEOM * lwgeom_symdifference_prec(const LWGEOM *geom1, const LWGEOM *geom2, double gridSize)
LWGEOM * lwgeom_intersection_prec(const LWGEOM *geom1, const LWGEOM *geom2, double gridSize)
LWGEOM * lwgeom_linemerge_directed(const LWGEOM *geom1, int directed)
#define TRIANGLETYPE
Definition: liblwgeom.h:115
LWGEOM * lwgeom_triangulate_polygon(const LWGEOM *geom)
Take vertices of a polygon and build a constrained triangulation that respects the boundary of the po...
LWGEOM * lwgeom_clip_by_rect(const LWGEOM *geom1, double x0, double y0, double x1, double y1)
LWGEOM * lwgeom_reduceprecision(const LWGEOM *geom, double gridSize)
LWMPOINT * lwgeom_to_points(const LWGEOM *lwgeom, uint32_t npoints, int32_t seed)
#define LW_TRUE
Return types for functions with status returns.
Definition: liblwgeom.h:93
#define SRID_UNKNOWN
Unknown SRID value.
Definition: liblwgeom.h:215
LWPOLY * lwpoly_construct_empty(int32_t srid, char hasz, char hasm)
Definition: lwpoly.c:161
LWGEOM * lwgeom_construct_empty(uint8_t type, int32_t srid, char hasz, char hasm)
Definition: lwgeom.c:2191
void lwgeom_add_bbox(LWGEOM *lwgeom)
Compute a bbox if not already computed.
Definition: lwgeom.c:695
This library is the generic geometry handling section of PostGIS.
static uint32_t lwgeom_get_type(const LWGEOM *geom)
Return LWTYPE number.
Definition: lwinline.h:141
static int lwgeom_is_empty(const LWGEOM *geom)
Return true or false depending on whether a geometry is an "empty" geometry (no vertices members)
Definition: lwinline.h:199
int value
Definition: genraster.py:62
int count
Definition: genraster.py:57
type
Definition: ovdump.py:42
Datum ST_ReducePrecision(PG_FUNCTION_ARGS)
Datum polygonize_garray(PG_FUNCTION_ARGS)
Datum ST_Intersection(PG_FUNCTION_ARGS)
Datum isvaliddetail(PG_FUNCTION_ARGS)
Datum isvalid(PG_FUNCTION_ARGS)
Datum ST_SimplifyPolygonHull(PG_FUNCTION_ARGS)
Datum ST_Snap(PG_FUNCTION_ARGS)
Datum ST_Union(PG_FUNCTION_ARGS)
Datum ST_BuildArea(PG_FUNCTION_ARGS)
Datum postgis_geos_compiled_version(PG_FUNCTION_ARGS)
Datum ST_Split(PG_FUNCTION_ARGS)
Datum ST_Node(PG_FUNCTION_ARGS)
Datum LWGEOM_dfullywithin(PG_FUNCTION_ARGS)
Returns boolean true if the second argument is fully contained in a buffer of the first argument.
Datum ST_GeneratePoints(PG_FUNCTION_ARGS)
Datum hausdorffdistancedensify(PG_FUNCTION_ARGS)
Datum ST_OffsetCurve(PG_FUNCTION_ARGS)
Datum cluster_within_distance_garray(PG_FUNCTION_ARGS)
Datum centroid(PG_FUNCTION_ARGS)
Datum ST_MinimumClearance(PG_FUNCTION_ARGS)
Datum GEOSnoop(PG_FUNCTION_ARGS)
Datum ST_FrechetDistance(PG_FUNCTION_ARGS)
Datum buffer(PG_FUNCTION_ARGS)
Datum convexhull(PG_FUNCTION_ARGS)
Datum ST_ClipByBox2d(PG_FUNCTION_ARGS)
static uint32_t array_nelems_not_null(ArrayType *array)
Datum ST_SymDifference(PG_FUNCTION_ARGS)
Datum issimple(PG_FUNCTION_ARGS)
Datum symdifference(PG_FUNCTION_ARGS)
LWGEOM ** ARRAY2LWGEOM(ArrayType *array, uint32_t nelems, int *is3d, int *srid)
Datum topologypreservesimplify(PG_FUNCTION_ARGS)
PG_FUNCTION_INFO_V1(postgis_geos_version)
Datum ST_MaximumInscribedCircle(PG_FUNCTION_ARGS)
GSERIALIZED * GEOS2POSTGIS(GEOSGeom geom, char want3d)
Datum ST_DelaunayTriangles(PG_FUNCTION_ARGS)
Datum ST_ConcaveHull(PG_FUNCTION_ARGS)
Datum boundary(PG_FUNCTION_ARGS)
Datum postgis_geos_version(PG_FUNCTION_ARGS)
GEOSGeometry * POSTGIS2GEOS(const GSERIALIZED *pglwgeom)
Datum ST_UnaryUnion(PG_FUNCTION_ARGS)
Datum ST_Equals(PG_FUNCTION_ARGS)
Datum ST_OrientedEnvelope(PG_FUNCTION_ARGS)
Datum ST_SharedPaths(PG_FUNCTION_ARGS)
Datum ST_MinimumClearanceLine(PG_FUNCTION_ARGS)
Datum ST_Voronoi(PG_FUNCTION_ARGS)
Datum clusterintersecting_garray(PG_FUNCTION_ARGS)
GEOSGeometry ** ARRAY2GEOS(ArrayType *array, uint32_t nelems, int *is3d, int *srid)
Datum ST_Difference(PG_FUNCTION_ARGS)
Datum pgis_union_geometry_array(PG_FUNCTION_ARGS)
Datum ST_TriangulatePolygon(PG_FUNCTION_ARGS)
Datum ST_LargestEmptyCircle(PG_FUNCTION_ARGS)
Datum hausdorffdistance(PG_FUNCTION_ARGS)
Datum isvalidreason(PG_FUNCTION_ARGS)
Datum isring(PG_FUNCTION_ARGS)
Datum linemerge(PG_FUNCTION_ARGS)
Datum pointonsurface(PG_FUNCTION_ARGS)
#define HANDLE_GEOS_ERROR(label)
unsigned int int32
Definition: shpopen.c:54
#define POSTGIS_GEOS_VERSION
Definition: sqldefines.h:11
double ymax
Definition: liblwgeom.h:357
double xmax
Definition: liblwgeom.h:355
double ymin
Definition: liblwgeom.h:356
double xmin
Definition: liblwgeom.h:354
lwflags_t flags
Definition: liblwgeom.h:353
GBOX * bbox
Definition: liblwgeom.h:458
lwflags_t flags
Definition: liblwgeom.h:461