PostGIS  2.5.7dev-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 
38 #include "access/htup_details.h"
39 
40 
41 /* PostGIS */
42 #include "lwgeom_functions_analytic.h" /* for point_in_polygon */
43 #include "lwgeom_geos.h"
44 #include "liblwgeom.h"
45 #include "lwgeom_rtree.h"
46 #include "lwgeom_geos_prepared.h"
47 
48 #include "float.h" /* for DBL_DIG */
49 
50 
51 /* Return NULL on GEOS error
52  *
53  * Prints error message only if it was not for interruption, in which
54  * case we let PostgreSQL deal with the error.
55  */
56 #define HANDLE_GEOS_ERROR(label) \
57  { \
58  if (strstr(lwgeom_geos_errmsg, "InterruptedException")) \
59  ereport(ERROR, \
60  (errcode(ERRCODE_QUERY_CANCELED), errmsg("canceling statement due to user request"))); \
61  else \
62  lwpgerror("%s: %s", (label), lwgeom_geos_errmsg); \
63  PG_RETURN_NULL(); \
64  }
65 
66 /*
67 ** Prototypes for SQL-bound functions
68 */
69 Datum relate_full(PG_FUNCTION_ARGS);
70 Datum relate_pattern(PG_FUNCTION_ARGS);
71 Datum disjoint(PG_FUNCTION_ARGS);
72 Datum touches(PG_FUNCTION_ARGS);
73 Datum geos_intersects(PG_FUNCTION_ARGS);
74 Datum crosses(PG_FUNCTION_ARGS);
75 Datum contains(PG_FUNCTION_ARGS);
76 Datum containsproperly(PG_FUNCTION_ARGS);
77 Datum covers(PG_FUNCTION_ARGS);
78 Datum overlaps(PG_FUNCTION_ARGS);
79 Datum isvalid(PG_FUNCTION_ARGS);
80 Datum isvalidreason(PG_FUNCTION_ARGS);
81 Datum isvaliddetail(PG_FUNCTION_ARGS);
82 Datum buffer(PG_FUNCTION_ARGS);
83 Datum geos_intersection(PG_FUNCTION_ARGS);
84 Datum convexhull(PG_FUNCTION_ARGS);
85 Datum topologypreservesimplify(PG_FUNCTION_ARGS);
86 Datum geos_difference(PG_FUNCTION_ARGS);
87 Datum boundary(PG_FUNCTION_ARGS);
88 Datum symdifference(PG_FUNCTION_ARGS);
89 Datum geos_geomunion(PG_FUNCTION_ARGS);
90 Datum issimple(PG_FUNCTION_ARGS);
91 Datum isring(PG_FUNCTION_ARGS);
92 Datum pointonsurface(PG_FUNCTION_ARGS);
93 Datum GEOSnoop(PG_FUNCTION_ARGS);
94 Datum postgis_geos_version(PG_FUNCTION_ARGS);
95 Datum centroid(PG_FUNCTION_ARGS);
96 Datum polygonize_garray(PG_FUNCTION_ARGS);
97 Datum clusterintersecting_garray(PG_FUNCTION_ARGS);
98 Datum cluster_within_distance_garray(PG_FUNCTION_ARGS);
99 Datum linemerge(PG_FUNCTION_ARGS);
100 Datum coveredby(PG_FUNCTION_ARGS);
101 Datum hausdorffdistance(PG_FUNCTION_ARGS);
102 Datum hausdorffdistancedensify(PG_FUNCTION_ARGS);
103 Datum ST_FrechetDistance(PG_FUNCTION_ARGS);
104 Datum ST_UnaryUnion(PG_FUNCTION_ARGS);
105 Datum ST_Equals(PG_FUNCTION_ARGS);
106 Datum ST_BuildArea(PG_FUNCTION_ARGS);
107 Datum ST_DelaunayTriangles(PG_FUNCTION_ARGS);
108 
109 Datum pgis_union_geometry_array(PG_FUNCTION_ARGS);
110 
111 /*
112 ** Prototypes end
113 */
114 
115 
116 PG_FUNCTION_INFO_V1(postgis_geos_version);
117 Datum postgis_geos_version(PG_FUNCTION_ARGS)
118 {
119  const char *ver = lwgeom_geos_version();
120  text *result = cstring_to_text(ver);
121  PG_RETURN_POINTER(result);
122 }
123 
124 static char
125 is_poly(const GSERIALIZED* g)
126 {
127  int type = gserialized_get_type(g);
128  return type == POLYGONTYPE || type == MULTIPOLYGONTYPE;
129 }
130 
131 static char
132 is_point(const GSERIALIZED* g)
133 {
134  int type = gserialized_get_type(g);
135  return type == POINTTYPE || type == MULTIPOINTTYPE;
136 }
137 
138 /* Utility function that checks a LWPOINT and a GSERIALIZED poly against a cache.
139  * Serialized poly may be a multipart.
140  */
141 static int
142 pip_short_circuit(RTREE_POLY_CACHE* poly_cache, LWPOINT* point, GSERIALIZED* gpoly)
143 {
144  int result;
145 
146  if ( poly_cache && poly_cache->ringIndices )
147  {
148  result = point_in_multipolygon_rtree(poly_cache->ringIndices, poly_cache->polyCount, poly_cache->ringCounts, point);
149  }
150  else
151  {
152  LWGEOM* poly = lwgeom_from_gserialized(gpoly);
153  if ( lwgeom_get_type(poly) == POLYGONTYPE )
154  {
155  result = point_in_polygon(lwgeom_as_lwpoly(poly), point);
156  }
157  else
158  {
159  result = point_in_multipolygon(lwgeom_as_lwmpoly(poly), point);
160  }
161  lwgeom_free(poly);
162  }
163 
164  return result;
165 }
166 
174 PG_FUNCTION_INFO_V1(hausdorffdistance);
175 Datum hausdorffdistance(PG_FUNCTION_ARGS)
176 {
177  GSERIALIZED *geom1;
178  GSERIALIZED *geom2;
179  GEOSGeometry *g1;
180  GEOSGeometry *g2;
181  double result;
182  int retcode;
183 
184  geom1 = PG_GETARG_GSERIALIZED_P(0);
185  geom2 = PG_GETARG_GSERIALIZED_P(1);
186 
187  if ( gserialized_is_empty(geom1) || gserialized_is_empty(geom2) )
188  PG_RETURN_NULL();
189 
190  initGEOS(lwpgnotice, lwgeom_geos_error);
191 
192  g1 = POSTGIS2GEOS(geom1);
193  if (!g1)
194  HANDLE_GEOS_ERROR("First argument geometry could not be converted to GEOS");
195 
196  g2 = POSTGIS2GEOS(geom2);
197  if (!g2)
198  {
199  GEOSGeom_destroy(g1);
200  HANDLE_GEOS_ERROR("Second argument geometry could not be converted to GEOS");
201  }
202 
203  retcode = GEOSHausdorffDistance(g1, g2, &result);
204  GEOSGeom_destroy(g1);
205  GEOSGeom_destroy(g2);
206 
207  if (retcode == 0) HANDLE_GEOS_ERROR("GEOSHausdorffDistance");
208 
209  PG_FREE_IF_COPY(geom1, 0);
210  PG_FREE_IF_COPY(geom2, 1);
211 
212  PG_RETURN_FLOAT8(result);
213 }
214 
223 PG_FUNCTION_INFO_V1(hausdorffdistancedensify);
224 Datum hausdorffdistancedensify(PG_FUNCTION_ARGS)
225 {
226  GSERIALIZED *geom1;
227  GSERIALIZED *geom2;
228  GEOSGeometry *g1;
229  GEOSGeometry *g2;
230  double densifyFrac;
231  double result;
232  int retcode;
233 
234  geom1 = PG_GETARG_GSERIALIZED_P(0);
235  geom2 = PG_GETARG_GSERIALIZED_P(1);
236  densifyFrac = PG_GETARG_FLOAT8(2);
237 
238  if ( gserialized_is_empty(geom1) || gserialized_is_empty(geom2) )
239  PG_RETURN_NULL();
240 
241  initGEOS(lwpgnotice, lwgeom_geos_error);
242 
243  g1 = POSTGIS2GEOS(geom1);
244  if (!g1)
245  HANDLE_GEOS_ERROR("First argument geometry could not be converted to GEOS");
246 
247  g2 = POSTGIS2GEOS(geom2);
248  if (!g2)
249  {
250  GEOSGeom_destroy(g1);
251  HANDLE_GEOS_ERROR("Second argument geometry could not be converted to GEOS");
252  }
253 
254  retcode = GEOSHausdorffDistanceDensify(g1, g2, densifyFrac, &result);
255  GEOSGeom_destroy(g1);
256  GEOSGeom_destroy(g2);
257 
258  if (retcode == 0) HANDLE_GEOS_ERROR("GEOSHausdorffDistanceDensify");
259 
260  PG_FREE_IF_COPY(geom1, 0);
261  PG_FREE_IF_COPY(geom2, 1);
262 
263  PG_RETURN_FLOAT8(result);
264 }
265 
273 PG_FUNCTION_INFO_V1(ST_FrechetDistance);
274 Datum ST_FrechetDistance(PG_FUNCTION_ARGS)
275 {
276 #if POSTGIS_GEOS_VERSION < 37
277 
278  lwpgerror("The GEOS version this PostGIS binary "
279  "was compiled against (%d) doesn't support "
280  "'GEOSFechetDistance' function (3.7.0+ required)",
281  POSTGIS_GEOS_VERSION);
282  PG_RETURN_NULL();
283 
284 #else /* POSTGIS_GEOS_VERSION >= 37 */
285  GSERIALIZED *geom1;
286  GSERIALIZED *geom2;
287  GEOSGeometry *g1;
288  GEOSGeometry *g2;
289  double densifyFrac;
290  double result;
291  int retcode;
292 
293  geom1 = PG_GETARG_GSERIALIZED_P(0);
294  geom2 = PG_GETARG_GSERIALIZED_P(1);
295  densifyFrac = PG_GETARG_FLOAT8(2);
296 
297  if ( gserialized_is_empty(geom1) || gserialized_is_empty(geom2) )
298  PG_RETURN_NULL();
299 
300  initGEOS(lwpgnotice, lwgeom_geos_error);
301 
302  g1 = POSTGIS2GEOS(geom1);
303  if (!g1)
304  HANDLE_GEOS_ERROR("First argument geometry could not be converted to GEOS");
305 
306  g2 = POSTGIS2GEOS(geom2);
307  if (!g2)
308  {
309  GEOSGeom_destroy(g1);
310  HANDLE_GEOS_ERROR("Second argument geometry could not be converted to GEOS");
311  }
312 
313  if (densifyFrac <= 0.0)
314  {
315  retcode = GEOSFrechetDistance(g1, g2, &result);
316  }
317  else
318  {
319  retcode = GEOSFrechetDistanceDensify(g1, g2, densifyFrac, &result);
320  }
321 
322  GEOSGeom_destroy(g1);
323  GEOSGeom_destroy(g2);
324 
325  if (retcode == 0) HANDLE_GEOS_ERROR("GEOSFrechetDistance");
326 
327  PG_FREE_IF_COPY(geom1, 0);
328  PG_FREE_IF_COPY(geom2, 1);
329 
330  PG_RETURN_FLOAT8(result);
331 
332 #endif /* POSTGIS_GEOS_VERSION >= 37 */
333 }
334 
342 PG_FUNCTION_INFO_V1(pgis_union_geometry_array);
343 Datum pgis_union_geometry_array(PG_FUNCTION_ARGS)
344 {
345  ArrayType *array;
346 
347  ArrayIterator iterator;
348  Datum value;
349  bool isnull;
350 
351  int is3d = LW_FALSE, gotsrid = LW_FALSE;
352  int nelems = 0, geoms_size = 0, curgeom = 0, count = 0;
353 
354  GSERIALIZED *gser_out = NULL;
355 
356  GEOSGeometry *g = NULL;
357  GEOSGeometry *g_union = NULL;
358  GEOSGeometry **geoms = NULL;
359 
360  int srid = SRID_UNKNOWN;
361 
362  int empty_type = 0;
363 
364  /* Null array, null geometry (should be empty?) */
365  if ( PG_ARGISNULL(0) )
366  PG_RETURN_NULL();
367 
368  array = PG_GETARG_ARRAYTYPE_P(0);
369  nelems = ArrayGetNItems(ARR_NDIM(array), ARR_DIMS(array));
370 
371  /* Empty array? Null return */
372  if ( nelems == 0 ) PG_RETURN_NULL();
373 
374  /* Quick scan for nulls */
375 #if POSTGIS_PGSQL_VERSION >= 95
376  iterator = array_create_iterator(array, 0, NULL);
377 #else
378  iterator = array_create_iterator(array, 0);
379 #endif
380  while( array_iterate(iterator, &value, &isnull) )
381  {
382  /* Skip null array items */
383  if (isnull) continue;
384  count++;
385  }
386  array_free_iterator(iterator);
387 
388 
389  /* All-nulls? Return null */
390  if ( count == 0 )
391  PG_RETURN_NULL();
392 
393  /* One geom, good geom? Return it */
394  if ( count == 1 && nelems == 1 )
395  {
396 #pragma GCC diagnostic push
397 #pragma GCC diagnostic ignored "-Wsign-compare"
398  PG_RETURN_POINTER((GSERIALIZED *)(ARR_DATA_PTR(array)));
399 #pragma GCC diagnostic pop
400  }
401 
402  /* Ok, we really need GEOS now ;) */
403  initGEOS(lwpgnotice, lwgeom_geos_error);
404 
405  /*
406  ** Collect the non-empty inputs and stuff them into a GEOS collection
407  */
408  geoms_size = nelems;
409  geoms = palloc(sizeof(GEOSGeometry*) * geoms_size);
410 
411  /*
412  ** We need to convert the array of GSERIALIZED into a GEOS collection.
413  ** First make an array of GEOS geometries.
414  */
415 #if POSTGIS_PGSQL_VERSION >= 95
416  iterator = array_create_iterator(array, 0, NULL);
417 #else
418  iterator = array_create_iterator(array, 0);
419 #endif
420  while( array_iterate(iterator, &value, &isnull) )
421  {
422  GSERIALIZED *gser_in;
423 
424  /* Skip null array items */
425  if (isnull) continue;
426  gser_in = (GSERIALIZED *)DatumGetPointer(value);
427 
428  /* Check for SRID mismatch in array elements */
429  if ( gotsrid )
430  {
431  error_if_srid_mismatch(srid, gserialized_get_srid(gser_in));
432  }
433  else
434  {
435  /* Initialize SRID/dimensions info */
436  srid = gserialized_get_srid(gser_in);
437  is3d = gserialized_has_z(gser_in);
438  gotsrid = 1;
439  }
440 
441  /* Don't include empties in the union */
442  if ( gserialized_is_empty(gser_in) )
443  {
444  int gser_type = gserialized_get_type(gser_in);
445  if (gser_type > empty_type)
446  {
447  empty_type = gser_type;
448  POSTGIS_DEBUGF(4, "empty_type = %d gser_type = %d", empty_type, gser_type);
449  }
450  }
451  else
452  {
453  g = POSTGIS2GEOS(gser_in);
454 
455  /* Uh oh! Exception thrown at construction... */
456  if ( ! g )
457  {
458  HANDLE_GEOS_ERROR(
459  "One of the geometries in the set "
460  "could not be converted to GEOS");
461  }
462 
463  /* Ensure we have enough space in our storage array */
464  if ( curgeom == geoms_size )
465  {
466  geoms_size *= 2;
467  geoms = repalloc( geoms, sizeof(GEOSGeometry*) * geoms_size );
468  }
469 
470  geoms[curgeom] = g;
471  curgeom++;
472  }
473 
474  }
475  array_free_iterator(iterator);
476 
477  /*
478  ** Take our GEOS geometries and turn them into a GEOS collection,
479  ** then pass that into cascaded union.
480  */
481  if (curgeom > 0)
482  {
483  g = GEOSGeom_createCollection(GEOS_GEOMETRYCOLLECTION, geoms, curgeom);
484  if (!g) HANDLE_GEOS_ERROR("Could not create GEOS COLLECTION from geometry array");
485 
486  g_union = GEOSUnaryUnion(g);
487  GEOSGeom_destroy(g);
488  if (!g_union) HANDLE_GEOS_ERROR("GEOSUnaryUnion");
489 
490  GEOSSetSRID(g_union, srid);
491  gser_out = GEOS2POSTGIS(g_union, is3d);
492  GEOSGeom_destroy(g_union);
493  }
494  /* No real geometries in our array, any empties? */
495  else
496  {
497  /* If it was only empties, we'll return the largest type number */
498  if ( empty_type > 0 )
499  {
500  PG_RETURN_POINTER(geometry_serialize(lwgeom_construct_empty(empty_type, srid, is3d, 0)));
501  }
502  /* Nothing but NULL, returns NULL */
503  else
504  {
505  PG_RETURN_NULL();
506  }
507  }
508 
509  if ( ! gser_out )
510  {
511  /* Union returned a NULL geometry */
512  PG_RETURN_NULL();
513  }
514 
515  PG_RETURN_POINTER(gser_out);
516 }
517 
524 PG_FUNCTION_INFO_V1(ST_UnaryUnion);
525 Datum ST_UnaryUnion(PG_FUNCTION_ARGS)
526 {
527  GSERIALIZED *geom1;
528  GSERIALIZED *result;
529  LWGEOM *lwgeom1, *lwresult ;
530 
531  geom1 = PG_GETARG_GSERIALIZED_P(0);
532 
533  lwgeom1 = lwgeom_from_gserialized(geom1) ;
534 
535  lwresult = lwgeom_unaryunion(lwgeom1);
536  result = geometry_serialize(lwresult) ;
537 
538  lwgeom_free(lwgeom1) ;
539  lwgeom_free(lwresult) ;
540 
541  PG_FREE_IF_COPY(geom1, 0);
542 
543  PG_RETURN_POINTER(result);
544 }
545 
553 PG_FUNCTION_INFO_V1(geos_geomunion);
554 Datum geos_geomunion(PG_FUNCTION_ARGS)
555 {
556  GSERIALIZED *geom1;
557  GSERIALIZED *geom2;
558  GSERIALIZED *result;
559  LWGEOM *lwgeom1, *lwgeom2, *lwresult ;
560 
561  geom1 = PG_GETARG_GSERIALIZED_P(0);
562  geom2 = PG_GETARG_GSERIALIZED_P(1);
563 
564  lwgeom1 = lwgeom_from_gserialized(geom1) ;
565  lwgeom2 = lwgeom_from_gserialized(geom2) ;
566 
567  lwresult = lwgeom_union(lwgeom1, lwgeom2) ;
568  result = geometry_serialize(lwresult) ;
569 
570  lwgeom_free(lwgeom1) ;
571  lwgeom_free(lwgeom2) ;
572  lwgeom_free(lwresult) ;
573 
574  PG_FREE_IF_COPY(geom1, 0);
575  PG_FREE_IF_COPY(geom2, 1);
576 
577  PG_RETURN_POINTER(result);
578 }
579 
585 PG_FUNCTION_INFO_V1(symdifference);
586 Datum symdifference(PG_FUNCTION_ARGS)
587 {
588  GSERIALIZED *geom1;
589  GSERIALIZED *geom2;
590  GSERIALIZED *result;
591  LWGEOM *lwgeom1, *lwgeom2, *lwresult ;
592 
593  geom1 = PG_GETARG_GSERIALIZED_P(0);
594  geom2 = PG_GETARG_GSERIALIZED_P(1);
595 
596  lwgeom1 = lwgeom_from_gserialized(geom1) ;
597  lwgeom2 = lwgeom_from_gserialized(geom2) ;
598 
599  lwresult = lwgeom_symdifference(lwgeom1, lwgeom2) ;
600  result = geometry_serialize(lwresult) ;
601 
602  lwgeom_free(lwgeom1) ;
603  lwgeom_free(lwgeom2) ;
604  lwgeom_free(lwresult) ;
605 
606  PG_FREE_IF_COPY(geom1, 0);
607  PG_FREE_IF_COPY(geom2, 1);
608 
609  PG_RETURN_POINTER(result);
610 }
611 
612 
613 PG_FUNCTION_INFO_V1(boundary);
614 Datum boundary(PG_FUNCTION_ARGS)
615 {
616  GSERIALIZED *geom1;
617  GEOSGeometry *g1, *g3;
618  GSERIALIZED *result;
619  LWGEOM *lwgeom;
620  int srid;
621 
622  geom1 = PG_GETARG_GSERIALIZED_P(0);
623 
624  /* Empty.Boundary() == Empty */
625  if ( gserialized_is_empty(geom1) )
626  PG_RETURN_POINTER(geom1);
627 
628  srid = gserialized_get_srid(geom1);
629 
630  lwgeom = lwgeom_from_gserialized(geom1);
631  if ( ! lwgeom ) {
632  lwpgerror("POSTGIS2GEOS: unable to deserialize input");
633  PG_RETURN_NULL();
634  }
635 
636  /* GEOS doesn't do triangle type, so we special case that here */
637  if (lwgeom->type == TRIANGLETYPE)
638  {
639  lwgeom->type = LINETYPE;
640  result = geometry_serialize(lwgeom);
641  lwgeom_free(lwgeom);
642  PG_RETURN_POINTER(result);
643  }
644 
645  initGEOS(lwpgnotice, lwgeom_geos_error);
646 
647  g1 = LWGEOM2GEOS(lwgeom, 0);
648  lwgeom_free(lwgeom);
649 
650  if (!g1)
651  HANDLE_GEOS_ERROR("First argument geometry could not be converted to GEOS");
652 
653  g3 = GEOSBoundary(g1);
654 
655  if (!g3)
656  {
657  GEOSGeom_destroy(g1);
658  HANDLE_GEOS_ERROR("GEOSBoundary");
659  }
660 
661  POSTGIS_DEBUGF(3, "result: %s", GEOSGeomToWKT(g3));
662 
663  GEOSSetSRID(g3, srid);
664 
665  result = GEOS2POSTGIS(g3, gserialized_has_z(geom1));
666 
667  if (!result)
668  {
669  GEOSGeom_destroy(g1);
670  GEOSGeom_destroy(g3);
671  elog(NOTICE,
672  "GEOS2POSTGIS threw an error (result postgis geometry "
673  "formation)!");
674  PG_RETURN_NULL(); /* never get here */
675  }
676 
677  GEOSGeom_destroy(g1);
678  GEOSGeom_destroy(g3);
679 
680  PG_FREE_IF_COPY(geom1, 0);
681 
682  PG_RETURN_POINTER(result);
683 }
684 
685 PG_FUNCTION_INFO_V1(convexhull);
686 Datum convexhull(PG_FUNCTION_ARGS)
687 {
688  GSERIALIZED *geom1;
689  GEOSGeometry *g1, *g3;
690  GSERIALIZED *result;
691  LWGEOM *lwout;
692  int srid;
693  GBOX bbox;
694 
695  geom1 = PG_GETARG_GSERIALIZED_P(0);
696 
697  /* Empty.ConvexHull() == Empty */
698  if ( gserialized_is_empty(geom1) )
699  PG_RETURN_POINTER(geom1);
700 
701  srid = gserialized_get_srid(geom1);
702 
703  initGEOS(lwpgnotice, lwgeom_geos_error);
704 
705  g1 = POSTGIS2GEOS(geom1);
706 
707  if (!g1)
708  HANDLE_GEOS_ERROR("First argument geometry could not be converted to GEOS");
709 
710  g3 = GEOSConvexHull(g1);
711  GEOSGeom_destroy(g1);
712 
713  if (!g3) HANDLE_GEOS_ERROR("GEOSConvexHull");
714 
715  POSTGIS_DEBUGF(3, "result: %s", GEOSGeomToWKT(g3));
716 
717  GEOSSetSRID(g3, srid);
718 
719  lwout = GEOS2LWGEOM(g3, gserialized_has_z(geom1));
720  GEOSGeom_destroy(g3);
721 
722  if (!lwout)
723  {
724  elog(ERROR,
725  "convexhull() failed to convert GEOS geometry to LWGEOM");
726  PG_RETURN_NULL(); /* never get here */
727  }
728 
729  /* Copy input bbox if any */
730  if ( gserialized_get_gbox_p(geom1, &bbox) )
731  {
732  /* Force the box to have the same dimensionality as the lwgeom */
733  bbox.flags = lwout->flags;
734  lwout->bbox = gbox_copy(&bbox);
735  }
736 
737  result = geometry_serialize(lwout);
738  lwgeom_free(lwout);
739 
740  if (!result)
741  {
742  elog(ERROR,"GEOS convexhull() threw an error (result postgis geometry formation)!");
743  PG_RETURN_NULL(); /* never get here */
744  }
745 
746  PG_FREE_IF_COPY(geom1, 0);
747  PG_RETURN_POINTER(result);
748 }
749 
750 PG_FUNCTION_INFO_V1(topologypreservesimplify);
751 Datum topologypreservesimplify(PG_FUNCTION_ARGS)
752 {
753  GSERIALIZED *geom1;
754  double tolerance;
755  GEOSGeometry *g1, *g3;
756  GSERIALIZED *result;
757  uint32_t type;
758 
759  geom1 = PG_GETARG_GSERIALIZED_P(0);
760  tolerance = PG_GETARG_FLOAT8(1);
761 
762  /* Empty.Simplify() == Empty */
763  type = gserialized_get_type(geom1);
764  if ( gserialized_is_empty(geom1) || type == TINTYPE || type == TRIANGLETYPE )
765  PG_RETURN_POINTER(geom1);
766 
767  initGEOS(lwpgnotice, lwgeom_geos_error);
768 
769  g1 = POSTGIS2GEOS(geom1);
770  if (!g1)
771  HANDLE_GEOS_ERROR("First argument geometry could not be converted to GEOS");
772 
773  g3 = GEOSTopologyPreserveSimplify(g1,tolerance);
774  GEOSGeom_destroy(g1);
775 
776  if (!g3) HANDLE_GEOS_ERROR("GEOSTopologyPreserveSimplify");
777 
778  POSTGIS_DEBUGF(3, "result: %s", GEOSGeomToWKT(g3));
779 
780  GEOSSetSRID(g3, gserialized_get_srid(geom1));
781 
782  result = GEOS2POSTGIS(g3, gserialized_has_z(geom1));
783  GEOSGeom_destroy(g3);
784 
785  if (!result)
786  {
787  elog(ERROR,"GEOS topologypreservesimplify() threw an error (result postgis geometry formation)!");
788  PG_RETURN_NULL(); /* never get here */
789  }
790 
791  PG_FREE_IF_COPY(geom1, 0);
792  PG_RETURN_POINTER(result);
793 }
794 
795 PG_FUNCTION_INFO_V1(buffer);
796 Datum buffer(PG_FUNCTION_ARGS)
797 {
798  GSERIALIZED *geom1;
799  double size;
800  GEOSBufferParams *bufferparams;
801  GEOSGeometry *g1, *g3 = NULL;
802  GSERIALIZED *result;
803  int quadsegs = 8; /* the default */
804  int nargs;
805  int singleside = 0; /* the default */
806  enum
807  {
808  ENDCAP_ROUND = 1,
809  ENDCAP_FLAT = 2,
810  ENDCAP_SQUARE = 3
811  };
812  enum
813  {
814  JOIN_ROUND = 1,
815  JOIN_MITRE = 2,
816  JOIN_BEVEL = 3
817  };
818  static const double DEFAULT_MITRE_LIMIT = 5.0;
819  static const int DEFAULT_ENDCAP_STYLE = ENDCAP_ROUND;
820  static const int DEFAULT_JOIN_STYLE = JOIN_ROUND;
821 
822  double mitreLimit = DEFAULT_MITRE_LIMIT;
823  int endCapStyle = DEFAULT_ENDCAP_STYLE;
824  int joinStyle = DEFAULT_JOIN_STYLE;
825  char *param;
826  char *params = NULL;
827  LWGEOM *lwg;
828 
829  geom1 = PG_GETARG_GSERIALIZED_P(0);
830  size = PG_GETARG_FLOAT8(1);
831 
832  /* Empty.Buffer() == Empty[polygon] */
833  if ( gserialized_is_empty(geom1) )
834  {
835  lwg = lwpoly_as_lwgeom(lwpoly_construct_empty(
836  gserialized_get_srid(geom1),
837  0, 0)); // buffer wouldn't give back z or m anyway
838  PG_RETURN_POINTER(geometry_serialize(lwg));
839  }
840 
841  nargs = PG_NARGS();
842 
843  initGEOS(lwpgnotice, lwgeom_geos_error);
844 
845  g1 = POSTGIS2GEOS(geom1);
846  if (!g1)
847  HANDLE_GEOS_ERROR("First argument geometry could not be converted to GEOS");
848 
849  if (nargs > 2)
850  {
851  /* We strdup `cause we're going to modify it */
852  params = pstrdup(PG_GETARG_CSTRING(2));
853 
854  POSTGIS_DEBUGF(3, "Params: %s", params);
855 
856  for (param=params; ; param=NULL)
857  {
858  char *key, *val;
859  param = strtok(param, " ");
860  if (!param) break;
861  POSTGIS_DEBUGF(3, "Param: %s", param);
862 
863  key = param;
864  val = strchr(key, '=');
865  if (!val || *(val + 1) == '\0')
866  {
867  lwpgerror("Missing value for buffer parameter %s", key);
868  break;
869  }
870  *val = '\0';
871  ++val;
872 
873  POSTGIS_DEBUGF(3, "Param: %s : %s", key, val);
874 
875  if ( !strcmp(key, "endcap") )
876  {
877  /* Supported end cap styles:
878  * "round", "flat", "square"
879  */
880  if ( !strcmp(val, "round") )
881  {
882  endCapStyle = ENDCAP_ROUND;
883  }
884  else if ( !strcmp(val, "flat") ||
885  !strcmp(val, "butt") )
886  {
887  endCapStyle = ENDCAP_FLAT;
888  }
889  else if ( !strcmp(val, "square") )
890  {
891  endCapStyle = ENDCAP_SQUARE;
892  }
893  else
894  {
895  lwpgerror("Invalid buffer end cap "
896  "style: %s (accept: "
897  "'round', 'flat', 'butt' "
898  "or 'square'"
899  ")", val);
900  break;
901  }
902 
903  }
904  else if ( !strcmp(key, "join") )
905  {
906  if ( !strcmp(val, "round") )
907  {
908  joinStyle = JOIN_ROUND;
909  }
910  else if ( !strcmp(val, "mitre") ||
911  !strcmp(val, "miter") )
912  {
913  joinStyle = JOIN_MITRE;
914  }
915  else if ( !strcmp(val, "bevel") )
916  {
917  joinStyle = JOIN_BEVEL;
918  }
919  else
920  {
921  lwpgerror("Invalid buffer end cap "
922  "style: %s (accept: "
923  "'round', 'mitre', 'miter' "
924  " or 'bevel'"
925  ")", val);
926  break;
927  }
928  }
929  else if ( !strcmp(key, "mitre_limit") ||
930  !strcmp(key, "miter_limit") )
931  {
932  /* mitreLimit is a float */
933  mitreLimit = atof(val);
934  }
935  else if ( !strcmp(key, "quad_segs") )
936  {
937  /* quadrant segments is an int */
938  quadsegs = atoi(val);
939  }
940  else if ( !strcmp(key, "side") ||
941  !strcmp(key, "side") )
942  {
943  if ( !strcmp(val, "both") )
944  {
945  singleside = 0;
946  }
947  else if ( !strcmp(val, "left") )
948  {
949  singleside = 1;
950  }
951  else if ( !strcmp(val, "right") )
952  {
953  singleside = 1;
954  size *= -1;
955  }
956  else
957  {
958  lwpgerror("Invalid side parameter: %s (accept: 'right', 'left', 'both')", val);
959  break;
960  }
961  }
962  else
963  {
964  lwpgerror(
965  "Invalid buffer parameter: %s (accept: 'endcap', 'join', 'mitre_limit', 'miter_limit', 'quad_segs' and 'side')",
966  key);
967  break;
968  }
969  }
970 
971  pfree(params); /* was pstrduped */
972 
973  POSTGIS_DEBUGF(3, "endCap:%d joinStyle:%d mitreLimit:%g",
974  endCapStyle, joinStyle, mitreLimit);
975 
976  }
977 
978  bufferparams = GEOSBufferParams_create();
979  if (bufferparams)
980  {
981  if (GEOSBufferParams_setEndCapStyle(bufferparams, endCapStyle) &&
982  GEOSBufferParams_setJoinStyle(bufferparams, joinStyle) &&
983  GEOSBufferParams_setMitreLimit(bufferparams, mitreLimit) &&
984  GEOSBufferParams_setQuadrantSegments(bufferparams, quadsegs) &&
985  GEOSBufferParams_setSingleSided(bufferparams, singleside))
986  {
987  g3 = GEOSBufferWithParams(g1, bufferparams, size);
988  }
989  else
990  {
991  lwpgerror("Error setting buffer parameters.");
992  }
993  GEOSBufferParams_destroy(bufferparams);
994  }
995  else
996  {
997  lwpgerror("Error setting buffer parameters.");
998  }
999 
1000  GEOSGeom_destroy(g1);
1001 
1002  if (!g3) HANDLE_GEOS_ERROR("GEOSBuffer");
1003 
1004  POSTGIS_DEBUGF(3, "result: %s", GEOSGeomToWKT(g3));
1005 
1006  GEOSSetSRID(g3, gserialized_get_srid(geom1));
1007 
1008  result = GEOS2POSTGIS(g3, gserialized_has_z(geom1));
1009  GEOSGeom_destroy(g3);
1010 
1011  if (!result)
1012  {
1013  elog(ERROR,"GEOS buffer() threw an error (result postgis geometry formation)!");
1014  PG_RETURN_NULL(); /* never get here */
1015  }
1016 
1017  PG_FREE_IF_COPY(geom1, 0);
1018  PG_RETURN_POINTER(result);
1019 }
1020 
1021 /*
1022 * Generate a field of random points within the area of a
1023 * polygon or multipolygon. Throws an error for other geometry
1024 * types.
1025 */
1026 Datum ST_GeneratePoints(PG_FUNCTION_ARGS);
1027 PG_FUNCTION_INFO_V1(ST_GeneratePoints);
1028 Datum ST_GeneratePoints(PG_FUNCTION_ARGS)
1029 {
1030  GSERIALIZED *gser_input;
1031  GSERIALIZED *gser_result;
1032  LWGEOM *lwgeom_input;
1033  LWGEOM *lwgeom_result;
1034  int32 npoints;
1035 
1036  gser_input = PG_GETARG_GSERIALIZED_P(0);
1037  npoints = DatumGetInt32(DirectFunctionCall1(numeric_int4, PG_GETARG_DATUM(1)));
1038 
1039  if (npoints < 0)
1040  PG_RETURN_NULL();
1041 
1042  /* Types get checked in the code, we'll keep things small out there */
1043  lwgeom_input = lwgeom_from_gserialized(gser_input);
1044  lwgeom_result = (LWGEOM*)lwgeom_to_points(lwgeom_input, npoints);
1045  lwgeom_free(lwgeom_input);
1046  PG_FREE_IF_COPY(gser_input, 0);
1047 
1048  /* Return null as null */
1049  if (!lwgeom_result)
1050  PG_RETURN_NULL();
1051 
1052  /* Serialize and return */
1053  gser_result = gserialized_from_lwgeom(lwgeom_result, 0);
1054  lwgeom_free(lwgeom_result);
1055  PG_RETURN_POINTER(gser_result);
1056 }
1057 
1058 
1059 /*
1060 * Compute at offset curve to a line
1061 */
1062 Datum ST_OffsetCurve(PG_FUNCTION_ARGS);
1063 PG_FUNCTION_INFO_V1(ST_OffsetCurve);
1064 Datum ST_OffsetCurve(PG_FUNCTION_ARGS)
1065 {
1066  GSERIALIZED *gser_input;
1067  GSERIALIZED *gser_result;
1068  LWGEOM *lwgeom_input;
1069  LWGEOM *lwgeom_result;
1070  double size;
1071  int quadsegs = 8; /* the default */
1072  int nargs;
1073 
1074  enum
1075  {
1076  JOIN_ROUND = 1,
1077  JOIN_MITRE = 2,
1078  JOIN_BEVEL = 3
1079  };
1080 
1081  static const double DEFAULT_MITRE_LIMIT = 5.0;
1082  static const int DEFAULT_JOIN_STYLE = JOIN_ROUND;
1083  double mitreLimit = DEFAULT_MITRE_LIMIT;
1084  int joinStyle = DEFAULT_JOIN_STYLE;
1085  char *param = NULL;
1086  char *paramstr = NULL;
1087 
1088  /* Read SQL arguments */
1089  nargs = PG_NARGS();
1090  gser_input = PG_GETARG_GSERIALIZED_P(0);
1091  size = PG_GETARG_FLOAT8(1);
1092 
1093  /* For distance == 0, just return the input. */
1094  if (size == 0) PG_RETURN_POINTER(gser_input);
1095 
1096  /* Read the lwgeom, check for errors */
1097  lwgeom_input = lwgeom_from_gserialized(gser_input);
1098  if ( ! lwgeom_input )
1099  lwpgerror("ST_OffsetCurve: lwgeom_from_gserialized returned NULL");
1100 
1101  /* For empty inputs, just echo them back */
1102  if ( lwgeom_is_empty(lwgeom_input) )
1103  PG_RETURN_POINTER(gser_input);
1104 
1105  /* Process the optional arguments */
1106  if ( nargs > 2 )
1107  {
1108  text *wkttext = PG_GETARG_TEXT_P(2);
1109  paramstr = text_to_cstring(wkttext);
1110 
1111  POSTGIS_DEBUGF(3, "paramstr: %s", paramstr);
1112 
1113  for ( param=paramstr; ; param=NULL )
1114  {
1115  char *key, *val;
1116  param = strtok(param, " ");
1117  if (!param) break;
1118  POSTGIS_DEBUGF(3, "Param: %s", param);
1119 
1120  key = param;
1121  val = strchr(key, '=');
1122  if (!val || *(val + 1) == '\0')
1123  {
1124  lwpgerror("ST_OffsetCurve: Missing value for buffer parameter %s", key);
1125  break;
1126  }
1127  *val = '\0';
1128  ++val;
1129 
1130  POSTGIS_DEBUGF(3, "Param: %s : %s", key, val);
1131 
1132  if ( !strcmp(key, "join") )
1133  {
1134  if ( !strcmp(val, "round") )
1135  {
1136  joinStyle = JOIN_ROUND;
1137  }
1138  else if ( !(strcmp(val, "mitre") && strcmp(val, "miter")) )
1139  {
1140  joinStyle = JOIN_MITRE;
1141  }
1142  else if ( ! strcmp(val, "bevel") )
1143  {
1144  joinStyle = JOIN_BEVEL;
1145  }
1146  else
1147  {
1148  lwpgerror(
1149  "Invalid buffer end cap style: %s (accept: 'round', 'mitre', 'miter' or 'bevel')",
1150  val);
1151  break;
1152  }
1153  }
1154  else if ( !strcmp(key, "mitre_limit") ||
1155  !strcmp(key, "miter_limit") )
1156  {
1157  /* mitreLimit is a float */
1158  mitreLimit = atof(val);
1159  }
1160  else if ( !strcmp(key, "quad_segs") )
1161  {
1162  /* quadrant segments is an int */
1163  quadsegs = atoi(val);
1164  }
1165  else
1166  {
1167  lwpgerror(
1168  "Invalid buffer parameter: %s (accept: 'join', 'mitre_limit', 'miter_limit and 'quad_segs')",
1169  key);
1170  break;
1171  }
1172  }
1173  POSTGIS_DEBUGF(3, "joinStyle:%d mitreLimit:%g", joinStyle, mitreLimit);
1174  pfree(paramstr); /* alloc'ed in text_to_cstring */
1175  }
1176 
1177  lwgeom_result = lwgeom_offsetcurve(lwgeom_input, size, quadsegs, joinStyle, mitreLimit);
1178 
1179  if (!lwgeom_result)
1180  lwpgerror("ST_OffsetCurve: lwgeom_offsetcurve returned NULL");
1181 
1182  gser_result = gserialized_from_lwgeom(lwgeom_result, 0);
1183  lwgeom_free(lwgeom_input);
1184  lwgeom_free(lwgeom_result);
1185  PG_RETURN_POINTER(gser_result);
1186 }
1187 
1188 
1189 PG_FUNCTION_INFO_V1(geos_intersection);
1190 Datum geos_intersection(PG_FUNCTION_ARGS)
1191 {
1192  GSERIALIZED *geom1;
1193  GSERIALIZED *geom2;
1194  GSERIALIZED *result;
1195  LWGEOM *lwgeom1, *lwgeom2, *lwresult ;
1196 
1197  geom1 = PG_GETARG_GSERIALIZED_P(0);
1198  geom2 = PG_GETARG_GSERIALIZED_P(1);
1199 
1200  lwgeom1 = lwgeom_from_gserialized(geom1) ;
1201  lwgeom2 = lwgeom_from_gserialized(geom2) ;
1202 
1203  lwresult = lwgeom_intersection(lwgeom1, lwgeom2) ;
1204  result = geometry_serialize(lwresult) ;
1205 
1206  lwgeom_free(lwgeom1) ;
1207  lwgeom_free(lwgeom2) ;
1208  lwgeom_free(lwresult) ;
1209 
1210  PG_FREE_IF_COPY(geom1, 0);
1211  PG_FREE_IF_COPY(geom2, 1);
1212 
1213  PG_RETURN_POINTER(result);
1214 }
1215 
1221 PG_FUNCTION_INFO_V1(geos_difference);
1222 Datum geos_difference(PG_FUNCTION_ARGS)
1223 {
1224  GSERIALIZED *geom1;
1225  GSERIALIZED *geom2;
1226  GSERIALIZED *result;
1227  LWGEOM *lwgeom1, *lwgeom2, *lwresult ;
1228 
1229  geom1 = PG_GETARG_GSERIALIZED_P(0);
1230  geom2 = PG_GETARG_GSERIALIZED_P(1);
1231 
1232  lwgeom1 = lwgeom_from_gserialized(geom1) ;
1233  lwgeom2 = lwgeom_from_gserialized(geom2) ;
1234 
1235  lwresult = lwgeom_difference(lwgeom1, lwgeom2) ;
1236  result = geometry_serialize(lwresult) ;
1237 
1238  lwgeom_free(lwgeom1) ;
1239  lwgeom_free(lwgeom2) ;
1240  lwgeom_free(lwresult) ;
1241 
1242  PG_FREE_IF_COPY(geom1, 0);
1243  PG_FREE_IF_COPY(geom2, 1);
1244 
1245  PG_RETURN_POINTER(result);
1246 }
1247 
1252 PG_FUNCTION_INFO_V1(pointonsurface);
1253 Datum pointonsurface(PG_FUNCTION_ARGS)
1254 {
1255  GSERIALIZED *geom, *result;
1256  LWGEOM *lwgeom, *lwresult;
1257 
1258  geom = PG_GETARG_GSERIALIZED_P(0);
1259 
1260  lwgeom = lwgeom_from_gserialized(geom);
1261  lwresult = lwgeom_pointonsurface(lwgeom);
1262  lwgeom_free(lwgeom);
1263  PG_FREE_IF_COPY(geom, 0);
1264 
1265  if (!lwresult) PG_RETURN_NULL();
1266 
1267  result = geometry_serialize(lwresult);
1268  lwgeom_free(lwresult);
1269  PG_RETURN_POINTER(result);
1270 }
1271 
1272 PG_FUNCTION_INFO_V1(centroid);
1273 Datum centroid(PG_FUNCTION_ARGS)
1274 {
1275  GSERIALIZED *geom, *result;
1276  LWGEOM *lwgeom, *lwresult;
1277 
1278  geom = PG_GETARG_GSERIALIZED_P(0);
1279 
1280  lwgeom = lwgeom_from_gserialized(geom);
1281  lwresult = lwgeom_centroid(lwgeom);
1282  lwgeom_free(lwgeom);
1283  PG_FREE_IF_COPY(geom, 0);
1284 
1285  if (!lwresult) PG_RETURN_NULL();
1286 
1287  result = geometry_serialize(lwresult);
1288  lwgeom_free(lwresult);
1289  PG_RETURN_POINTER(result);
1290 }
1291 
1292 Datum ST_ClipByBox2d(PG_FUNCTION_ARGS);
1293 PG_FUNCTION_INFO_V1(ST_ClipByBox2d);
1294 Datum ST_ClipByBox2d(PG_FUNCTION_ARGS)
1295 {
1296  GSERIALIZED *geom1;
1297  GSERIALIZED *result;
1298  LWGEOM *lwgeom1, *lwresult ;
1299  const GBOX *bbox1;
1300  GBOX *bbox2;
1301 
1302  geom1 = PG_GETARG_GSERIALIZED_P(0);
1303  lwgeom1 = lwgeom_from_gserialized(geom1) ;
1304 
1305  bbox1 = lwgeom_get_bbox(lwgeom1);
1306  if ( ! bbox1 )
1307  {
1308  /* empty clips to empty, no matter rect */
1309  lwgeom_free(lwgeom1);
1310  PG_RETURN_POINTER(geom1);
1311  }
1312 
1313  /* WARNING: this is really a BOX2DF, use only xmin and ymin fields */
1314  bbox2 = (GBOX *)PG_GETARG_POINTER(1);
1315  bbox2->flags = 0;
1316 
1317  /* If bbox1 outside of bbox2, return empty */
1318  if ( ! gbox_overlaps_2d(bbox1, bbox2) )
1319  {
1320  lwresult = lwgeom_construct_empty(lwgeom1->type, lwgeom1->srid, 0, 0);
1321  lwgeom_free(lwgeom1);
1322  PG_FREE_IF_COPY(geom1, 0);
1323  result = geometry_serialize(lwresult) ;
1324  lwgeom_free(lwresult) ;
1325  PG_RETURN_POINTER(result);
1326  }
1327 
1328  /* if bbox1 is covered by bbox2, return lwgeom1 */
1329  if ( gbox_contains_2d(bbox2, bbox1) )
1330  {
1331  lwgeom_free(lwgeom1);
1332  PG_RETURN_POINTER(geom1);
1333  }
1334 
1335  lwresult = lwgeom_clip_by_rect(lwgeom1, bbox2->xmin, bbox2->ymin,
1336  bbox2->xmax, bbox2->ymax);
1337 
1338  lwgeom_free(lwgeom1);
1339  PG_FREE_IF_COPY(geom1, 0);
1340 
1341  if (!lwresult) PG_RETURN_NULL();
1342 
1343  result = geometry_serialize(lwresult) ;
1344  lwgeom_free(lwresult) ;
1345  PG_RETURN_POINTER(result);
1346 }
1347 
1348 
1349 /*---------------------------------------------*/
1350 
1357 void errorIfGeometryCollection(GSERIALIZED *g1, GSERIALIZED *g2)
1358 {
1359  int t1 = gserialized_get_type(g1);
1360  int t2 = gserialized_get_type(g2);
1361 
1362  char *hintmsg;
1363  char *hintwkt;
1364  size_t hintsz;
1365  LWGEOM *lwgeom;
1366 
1367  if ( t1 == COLLECTIONTYPE)
1368  {
1369  lwgeom = lwgeom_from_gserialized(g1);
1370  hintwkt = lwgeom_to_wkt(lwgeom, WKT_SFSQL, DBL_DIG, &hintsz);
1371  lwgeom_free(lwgeom);
1372  hintmsg = lwmessage_truncate(hintwkt, 0, hintsz-1, 80, 1);
1373  ereport(ERROR,
1374  (errmsg("Relate Operation called with a LWGEOMCOLLECTION type. This is unsupported."),
1375  errhint("Change argument 1: '%s'", hintmsg))
1376  );
1377  pfree(hintwkt);
1378  pfree(hintmsg);
1379  }
1380  else if (t2 == COLLECTIONTYPE)
1381  {
1382  lwgeom = lwgeom_from_gserialized(g2);
1383  hintwkt = lwgeom_to_wkt(lwgeom, WKT_SFSQL, DBL_DIG, &hintsz);
1384  hintmsg = lwmessage_truncate(hintwkt, 0, hintsz-1, 80, 1);
1385  lwgeom_free(lwgeom);
1386  ereport(ERROR,
1387  (errmsg("Relate Operation called with a LWGEOMCOLLECTION type. This is unsupported."),
1388  errhint("Change argument 2: '%s'", hintmsg))
1389  );
1390  pfree(hintwkt);
1391  pfree(hintmsg);
1392  }
1393 }
1394 
1395 PG_FUNCTION_INFO_V1(isvalid);
1396 Datum isvalid(PG_FUNCTION_ARGS)
1397 {
1398  GSERIALIZED *geom1;
1399  LWGEOM *lwgeom;
1400  char result;
1401  GEOSGeom g1;
1402 
1403  geom1 = PG_GETARG_GSERIALIZED_P(0);
1404 
1405  /* Empty.IsValid() == TRUE */
1406  if ( gserialized_is_empty(geom1) )
1407  PG_RETURN_BOOL(true);
1408 
1409  initGEOS(lwpgnotice, lwgeom_geos_error);
1410 
1411  lwgeom = lwgeom_from_gserialized(geom1);
1412  if ( ! lwgeom )
1413  {
1414  lwpgerror("unable to deserialize input");
1415  }
1416  g1 = LWGEOM2GEOS(lwgeom, 0);
1417  lwgeom_free(lwgeom);
1418 
1419  if ( ! g1 )
1420  {
1421  /* should we drop the following
1422  * notice now that we have ST_isValidReason ?
1423  */
1424  lwpgnotice("%s", lwgeom_geos_errmsg);
1425  PG_RETURN_BOOL(false);
1426  }
1427 
1428  result = GEOSisValid(g1);
1429  GEOSGeom_destroy(g1);
1430 
1431  if (result == 2)
1432  {
1433  elog(ERROR,"GEOS isvalid() threw an error!");
1434  PG_RETURN_NULL(); /*never get here */
1435  }
1436 
1437  PG_FREE_IF_COPY(geom1, 0);
1438  PG_RETURN_BOOL(result);
1439 }
1440 
1441 PG_FUNCTION_INFO_V1(isvalidreason);
1442 Datum isvalidreason(PG_FUNCTION_ARGS)
1443 {
1444  GSERIALIZED *geom = NULL;
1445  char *reason_str = NULL;
1446  text *result = NULL;
1447  const GEOSGeometry *g1 = NULL;
1448 
1449  geom = PG_GETARG_GSERIALIZED_P(0);
1450 
1451  initGEOS(lwpgnotice, lwgeom_geos_error);
1452 
1453  g1 = POSTGIS2GEOS(geom);
1454  if ( g1 )
1455  {
1456  reason_str = GEOSisValidReason(g1);
1457  GEOSGeom_destroy((GEOSGeometry *)g1);
1458  if (!reason_str) HANDLE_GEOS_ERROR("GEOSisValidReason");
1459  result = cstring_to_text(reason_str);
1460  GEOSFree(reason_str);
1461  }
1462  else
1463  {
1464  result = cstring_to_text(lwgeom_geos_errmsg);
1465  }
1466 
1467  PG_FREE_IF_COPY(geom, 0);
1468  PG_RETURN_POINTER(result);
1469 }
1470 
1471 PG_FUNCTION_INFO_V1(isvaliddetail);
1472 Datum isvaliddetail(PG_FUNCTION_ARGS)
1473 {
1474  GSERIALIZED *geom = NULL;
1475  const GEOSGeometry *g1 = NULL;
1476  char *values[3]; /* valid bool, reason text, location geometry */
1477  char *geos_reason = NULL;
1478  char *reason = NULL;
1479  GEOSGeometry *geos_location = NULL;
1480  LWGEOM *location = NULL;
1481  char valid = 0;
1482  HeapTupleHeader result;
1483  TupleDesc tupdesc;
1484  HeapTuple tuple;
1485  AttInMetadata *attinmeta;
1486  int flags = 0;
1487 
1488  /*
1489  * Build a tuple description for a
1490  * valid_detail tuple
1491  */
1492  get_call_result_type(fcinfo, 0, &tupdesc);
1493  BlessTupleDesc(tupdesc);
1494 
1495  /*
1496  * generate attribute metadata needed later to produce
1497  * tuples from raw C strings
1498  */
1499  attinmeta = TupleDescGetAttInMetadata(tupdesc);
1500 
1501  geom = PG_GETARG_GSERIALIZED_P(0);
1502 
1503  if ( PG_NARGS() > 1 && ! PG_ARGISNULL(1) )
1504  {
1505  flags = PG_GETARG_INT32(1);
1506  }
1507 
1508  initGEOS(lwpgnotice, lwgeom_geos_error);
1509 
1510  g1 = POSTGIS2GEOS(geom);
1511 
1512  if ( g1 )
1513  {
1514  valid = GEOSisValidDetail(g1, flags, &geos_reason, &geos_location);
1515  GEOSGeom_destroy((GEOSGeometry *)g1);
1516  if ( geos_reason )
1517  {
1518  reason = pstrdup(geos_reason);
1519  GEOSFree(geos_reason);
1520  }
1521  if ( geos_location )
1522  {
1523  location = GEOS2LWGEOM(geos_location, GEOSHasZ(geos_location));
1524  GEOSGeom_destroy(geos_location);
1525  }
1526 
1527  if (valid == 2)
1528  {
1529  /* NOTE: should only happen on OOM or similar */
1530  lwpgerror("GEOS isvaliddetail() threw an exception!");
1531  PG_RETURN_NULL(); /* never gets here */
1532  }
1533  }
1534  else
1535  {
1536  /* TODO: check lwgeom_geos_errmsg for validity error */
1537  reason = pstrdup(lwgeom_geos_errmsg);
1538  }
1539 
1540  /* the boolean validity */
1541  values[0] = valid ? "t" : "f";
1542 
1543  /* the reason */
1544  values[1] = reason;
1545 
1546  /* the location */
1547  values[2] = location ? lwgeom_to_hexwkb(location, WKB_EXTENDED, 0) : 0;
1548 
1549  tuple = BuildTupleFromCStrings(attinmeta, values);
1550  result = (HeapTupleHeader) palloc(tuple->t_len);
1551  memcpy(result, tuple->t_data, tuple->t_len);
1552  heap_freetuple(tuple);
1553 
1554  PG_RETURN_HEAPTUPLEHEADER(result);
1555 }
1556 
1564 PG_FUNCTION_INFO_V1(overlaps);
1565 Datum overlaps(PG_FUNCTION_ARGS)
1566 {
1567  GSERIALIZED *geom1;
1568  GSERIALIZED *geom2;
1569  GEOSGeometry *g1, *g2;
1570  char result;
1571  GBOX box1, box2;
1572 
1573  geom1 = PG_GETARG_GSERIALIZED_P(0);
1574  geom2 = PG_GETARG_GSERIALIZED_P(1);
1575 
1576  errorIfGeometryCollection(geom1,geom2);
1577  error_if_srid_mismatch(gserialized_get_srid(geom1), gserialized_get_srid(geom2));
1578 
1579  /* A.Overlaps(Empty) == FALSE */
1580  if ( gserialized_is_empty(geom1) || gserialized_is_empty(geom2) )
1581  PG_RETURN_BOOL(false);
1582 
1583  /*
1584  * short-circuit 1: if geom2 bounding box does not overlap
1585  * geom1 bounding box we can return FALSE.
1586  */
1587  if ( gserialized_get_gbox_p(geom1, &box1) &&
1588  gserialized_get_gbox_p(geom2, &box2) )
1589  {
1590  if ( ! gbox_overlaps_2d(&box1, &box2) )
1591  {
1592  PG_RETURN_BOOL(false);
1593  }
1594  }
1595 
1596  initGEOS(lwpgnotice, lwgeom_geos_error);
1597 
1598  g1 = POSTGIS2GEOS(geom1);
1599  if (!g1)
1600  HANDLE_GEOS_ERROR("First argument geometry could not be converted to GEOS");
1601 
1602  g2 = POSTGIS2GEOS(geom2);
1603 
1604  if (!g2)
1605  {
1606  GEOSGeom_destroy(g1);
1607  HANDLE_GEOS_ERROR("Second argument geometry could not be converted to GEOS");
1608  }
1609 
1610  result = GEOSOverlaps(g1,g2);
1611 
1612  GEOSGeom_destroy(g1);
1613  GEOSGeom_destroy(g2);
1614  if (result == 2) HANDLE_GEOS_ERROR("GEOSOverlaps");
1615 
1616  PG_FREE_IF_COPY(geom1, 0);
1617  PG_FREE_IF_COPY(geom2, 1);
1618 
1619  PG_RETURN_BOOL(result);
1620 }
1621 
1622 
1623 PG_FUNCTION_INFO_V1(contains);
1624 Datum contains(PG_FUNCTION_ARGS)
1625 {
1626  GSERIALIZED *geom1;
1627  GSERIALIZED *geom2;
1628  GEOSGeometry *g1, *g2;
1629  GBOX box1, box2;
1630  int result;
1631  PrepGeomCache *prep_cache;
1632 
1633  geom1 = PG_GETARG_GSERIALIZED_P(0);
1634  geom2 = PG_GETARG_GSERIALIZED_P(1);
1635 
1636  errorIfGeometryCollection(geom1,geom2);
1637  error_if_srid_mismatch(gserialized_get_srid(geom1), gserialized_get_srid(geom2));
1638 
1639  /* A.Contains(Empty) == FALSE */
1640  if ( gserialized_is_empty(geom1) || gserialized_is_empty(geom2) )
1641  PG_RETURN_BOOL(false);
1642 
1643  POSTGIS_DEBUG(3, "contains called.");
1644 
1645  /*
1646  ** short-circuit 1: if geom2 bounding box is not completely inside
1647  ** geom1 bounding box we can return FALSE.
1648  */
1649  if ( gserialized_get_gbox_p(geom1, &box1) &&
1650  gserialized_get_gbox_p(geom2, &box2) )
1651  {
1652  if (!gbox_contains_2d(&box1, &box2)) PG_RETURN_BOOL(false);
1653  }
1654 
1655  /*
1656  ** short-circuit 2: if geom2 is a point and geom1 is a polygon
1657  ** call the point-in-polygon function.
1658  */
1659  if (is_poly(geom1) && is_point(geom2))
1660  {
1661  GSERIALIZED* gpoly = is_poly(geom1) ? geom1 : geom2;
1662  GSERIALIZED* gpoint = is_point(geom1) ? geom1 : geom2;
1663  RTREE_POLY_CACHE* cache = GetRtreeCache(fcinfo, gpoly);
1664  int retval;
1665 
1666  POSTGIS_DEBUG(3, "Point in Polygon test requested...short-circuiting.");
1667  if (gserialized_get_type(gpoint) == POINTTYPE)
1668  {
1669  LWGEOM* point = lwgeom_from_gserialized(gpoint);
1670  int pip_result = pip_short_circuit(cache, lwgeom_as_lwpoint(point), gpoly);
1671  lwgeom_free(point);
1672 
1673  retval = (pip_result == 1); /* completely inside */
1674  }
1675  else if (gserialized_get_type(gpoint) == MULTIPOINTTYPE)
1676  {
1677  LWMPOINT* mpoint = lwgeom_as_lwmpoint(lwgeom_from_gserialized(gpoint));
1678  uint32_t i;
1679  int found_completely_inside = LW_FALSE;
1680 
1681  retval = LW_TRUE;
1682  for (i = 0; i < mpoint->ngeoms; i++)
1683  {
1684  /* We need to find at least one point that's completely inside the
1685  * polygons (pip_result == 1). As long as we have one point that's
1686  * completely inside, we can have as many as we want on the boundary
1687  * itself. (pip_result == 0)
1688  */
1689  int pip_result = pip_short_circuit(cache, mpoint->geoms[i], gpoly);
1690  if (pip_result == 1)
1691  found_completely_inside = LW_TRUE;
1692 
1693  if (pip_result == -1) /* completely outside */
1694  {
1695  retval = LW_FALSE;
1696  break;
1697  }
1698  }
1699 
1700  retval = retval && found_completely_inside;
1701  lwmpoint_free(mpoint);
1702  }
1703  else
1704  {
1705  /* Never get here */
1706  elog(ERROR,"Type isn't point or multipoint!");
1707  PG_RETURN_NULL();
1708  }
1709 
1710  PG_FREE_IF_COPY(geom1, 0);
1711  PG_FREE_IF_COPY(geom2, 1);
1712  PG_RETURN_BOOL(retval);
1713  }
1714  else
1715  {
1716  POSTGIS_DEBUGF(3, "Contains: type1: %d, type2: %d", gserialized_get_type(geom1), gserialized_get_type(geom2));
1717  }
1718 
1719  initGEOS(lwpgnotice, lwgeom_geos_error);
1720 
1721  prep_cache = GetPrepGeomCache( fcinfo, geom1, 0 );
1722 
1723  if ( prep_cache && prep_cache->prepared_geom && prep_cache->gcache.argnum == 1 )
1724  {
1725  g1 = POSTGIS2GEOS(geom2);
1726  if (!g1)
1727  HANDLE_GEOS_ERROR("Geometry could not be converted to GEOS");
1728 
1729  POSTGIS_DEBUG(4, "containsPrepared: cache is live, running preparedcontains");
1730  result = GEOSPreparedContains( prep_cache->prepared_geom, g1);
1731  GEOSGeom_destroy(g1);
1732  }
1733  else
1734  {
1735  g1 = POSTGIS2GEOS(geom1);
1736  if (!g1) HANDLE_GEOS_ERROR("First argument geometry could not be converted to GEOS");
1737  g2 = POSTGIS2GEOS(geom2);
1738  if (!g2)
1739  {
1740  HANDLE_GEOS_ERROR("Second argument geometry could not be converted to GEOS");
1741  GEOSGeom_destroy(g1);
1742  }
1743  POSTGIS_DEBUG(4, "containsPrepared: cache is not ready, running standard contains");
1744  result = GEOSContains( g1, g2);
1745  GEOSGeom_destroy(g1);
1746  GEOSGeom_destroy(g2);
1747  }
1748 
1749  if (result == 2) HANDLE_GEOS_ERROR("GEOSContains");
1750 
1751  PG_FREE_IF_COPY(geom1, 0);
1752  PG_FREE_IF_COPY(geom2, 1);
1753 
1754  PG_RETURN_BOOL(result);
1755 
1756 }
1757 
1758 PG_FUNCTION_INFO_V1(containsproperly);
1759 Datum containsproperly(PG_FUNCTION_ARGS)
1760 {
1761  GSERIALIZED * geom1;
1762  GSERIALIZED * geom2;
1763  char result;
1764  GBOX box1, box2;
1765  PrepGeomCache * prep_cache;
1766 
1767  geom1 = PG_GETARG_GSERIALIZED_P(0);
1768  geom2 = PG_GETARG_GSERIALIZED_P(1);
1769 
1770  errorIfGeometryCollection(geom1,geom2);
1771  error_if_srid_mismatch(gserialized_get_srid(geom1), gserialized_get_srid(geom2));
1772 
1773  /* A.ContainsProperly(Empty) == FALSE */
1774  if ( gserialized_is_empty(geom1) || gserialized_is_empty(geom2) )
1775  PG_RETURN_BOOL(false);
1776 
1777  /*
1778  * short-circuit: if geom2 bounding box is not completely inside
1779  * geom1 bounding box we can return FALSE.
1780  */
1781  if ( gserialized_get_gbox_p(geom1, &box1) &&
1782  gserialized_get_gbox_p(geom2, &box2) )
1783  {
1784  if ( ! gbox_contains_2d(&box1, &box2) )
1785  PG_RETURN_BOOL(false);
1786  }
1787 
1788  initGEOS(lwpgnotice, lwgeom_geos_error);
1789 
1790  prep_cache = GetPrepGeomCache( fcinfo, geom1, 0 );
1791 
1792  if ( prep_cache && prep_cache->prepared_geom && prep_cache->gcache.argnum == 1 )
1793  {
1794  GEOSGeometry *g = POSTGIS2GEOS(geom2);
1795  if (!g) HANDLE_GEOS_ERROR("First argument geometry could not be converted to GEOS");
1796  result = GEOSPreparedContainsProperly( prep_cache->prepared_geom, g);
1797  GEOSGeom_destroy(g);
1798  }
1799  else
1800  {
1801  GEOSGeometry *g2;
1802  GEOSGeometry *g1;
1803 
1804  g1 = POSTGIS2GEOS(geom1);
1805  if (!g1) HANDLE_GEOS_ERROR("First argument geometry could not be converted to GEOS");
1806  g2 = POSTGIS2GEOS(geom2);
1807  if (!g2)
1808  {
1809  GEOSGeom_destroy(g1);
1810  HANDLE_GEOS_ERROR("Second argument geometry could not be converted to GEOS");
1811  }
1812  result = GEOSRelatePattern( g1, g2, "T**FF*FF*" );
1813  GEOSGeom_destroy(g1);
1814  GEOSGeom_destroy(g2);
1815  }
1816 
1817  if (result == 2) HANDLE_GEOS_ERROR("GEOSContains");
1818 
1819  PG_FREE_IF_COPY(geom1, 0);
1820  PG_FREE_IF_COPY(geom2, 1);
1821 
1822  PG_RETURN_BOOL(result);
1823 }
1824 
1825 /*
1826  * Described at
1827  * http://lin-ear-th-inking.blogspot.com/2007/06/subtleties-of-ogc-covers-spatial.html
1828  */
1829 PG_FUNCTION_INFO_V1(covers);
1830 Datum covers(PG_FUNCTION_ARGS)
1831 {
1832  GSERIALIZED *geom1;
1833  GSERIALIZED *geom2;
1834  int result;
1835  GBOX box1, box2;
1836  PrepGeomCache *prep_cache;
1837 
1838  geom1 = PG_GETARG_GSERIALIZED_P(0);
1839  geom2 = PG_GETARG_GSERIALIZED_P(1);
1840 
1841  /* A.Covers(Empty) == FALSE */
1842  if ( gserialized_is_empty(geom1) || gserialized_is_empty(geom2) )
1843  PG_RETURN_BOOL(false);
1844 
1845  errorIfGeometryCollection(geom1,geom2);
1846  error_if_srid_mismatch(gserialized_get_srid(geom1), gserialized_get_srid(geom2));
1847 
1848  /*
1849  * short-circuit 1: if geom2 bounding box is not completely inside
1850  * geom1 bounding box we can return FALSE.
1851  */
1852  if ( gserialized_get_gbox_p(geom1, &box1) &&
1853  gserialized_get_gbox_p(geom2, &box2) )
1854  {
1855  if ( ! gbox_contains_2d(&box1, &box2) )
1856  {
1857  PG_RETURN_BOOL(false);
1858  }
1859  }
1860  /*
1861  * short-circuit 2: if geom2 is a point and geom1 is a polygon
1862  * call the point-in-polygon function.
1863  */
1864  if (is_poly(geom1) && is_point(geom2))
1865  {
1866  GSERIALIZED* gpoly = is_poly(geom1) ? geom1 : geom2;
1867  GSERIALIZED* gpoint = is_point(geom1) ? geom1 : geom2;
1868  RTREE_POLY_CACHE* cache = GetRtreeCache(fcinfo, gpoly);
1869  int retval;
1870 
1871  POSTGIS_DEBUG(3, "Point in Polygon test requested...short-circuiting.");
1872  if (gserialized_get_type(gpoint) == POINTTYPE)
1873  {
1874  LWGEOM* point = lwgeom_from_gserialized(gpoint);
1875  int pip_result = pip_short_circuit(cache, lwgeom_as_lwpoint(point), gpoly);
1876  lwgeom_free(point);
1877 
1878  retval = (pip_result != -1); /* not outside */
1879  }
1880  else if (gserialized_get_type(gpoint) == MULTIPOINTTYPE)
1881  {
1882  LWMPOINT* mpoint = lwgeom_as_lwmpoint(lwgeom_from_gserialized(gpoint));
1883  uint32_t i;
1884 
1885  retval = LW_TRUE;
1886  for (i = 0; i < mpoint->ngeoms; i++)
1887  {
1888  int pip_result = pip_short_circuit(cache, mpoint->geoms[i], gpoly);
1889  if (pip_result == -1)
1890  {
1891  retval = LW_FALSE;
1892  break;
1893  }
1894  }
1895 
1896  lwmpoint_free(mpoint);
1897  }
1898  else
1899  {
1900  /* Never get here */
1901  elog(ERROR,"Type isn't point or multipoint!");
1902  PG_RETURN_NULL();
1903  }
1904 
1905  PG_FREE_IF_COPY(geom1, 0);
1906  PG_FREE_IF_COPY(geom2, 1);
1907  PG_RETURN_BOOL(retval);
1908  }
1909  else
1910  {
1911  POSTGIS_DEBUGF(3, "Covers: type1: %d, type2: %d", gserialized_get_type(geom1), gserialized_get_type(geom2));
1912  }
1913 
1914  initGEOS(lwpgnotice, lwgeom_geos_error);
1915 
1916  prep_cache = GetPrepGeomCache( fcinfo, geom1, 0 );
1917 
1918  if ( prep_cache && prep_cache->prepared_geom && prep_cache->gcache.argnum == 1 )
1919  {
1920  GEOSGeometry *g1 = POSTGIS2GEOS(geom2);
1921  if (!g1) HANDLE_GEOS_ERROR("First argument geometry could not be converted to GEOS");
1922  result = GEOSPreparedCovers( prep_cache->prepared_geom, g1);
1923  GEOSGeom_destroy(g1);
1924  }
1925  else
1926  {
1927  GEOSGeometry *g1;
1928  GEOSGeometry *g2;
1929 
1930  g1 = POSTGIS2GEOS(geom1);
1931  if (!g1) HANDLE_GEOS_ERROR("First argument geometry could not be converted to GEOS");
1932  g2 = POSTGIS2GEOS(geom2);
1933  if (!g2)
1934  {
1935  GEOSGeom_destroy(g1);
1936  HANDLE_GEOS_ERROR("Second argument geometry could not be converted to GEOS");
1937  }
1938  result = GEOSRelatePattern( g1, g2, "******FF*" );
1939  GEOSGeom_destroy(g1);
1940  GEOSGeom_destroy(g2);
1941  }
1942 
1943  if (result == 2) HANDLE_GEOS_ERROR("GEOSCovers");
1944 
1945  PG_FREE_IF_COPY(geom1, 0);
1946  PG_FREE_IF_COPY(geom2, 1);
1947 
1948  PG_RETURN_BOOL(result);
1949 
1950 }
1951 
1952 
1960 /*
1961  * Described at:
1962  * http://lin-ear-th-inking.blogspot.com/2007/06/subtleties-of-ogc-covers-spatial.html
1963  */
1964 PG_FUNCTION_INFO_V1(coveredby);
1965 Datum coveredby(PG_FUNCTION_ARGS)
1966 {
1967  GSERIALIZED *geom1;
1968  GSERIALIZED *geom2;
1969  GEOSGeometry *g1, *g2;
1970  int result;
1971  GBOX box1, box2;
1972  char *patt = "**F**F***";
1973 
1974  geom1 = PG_GETARG_GSERIALIZED_P(0);
1975  geom2 = PG_GETARG_GSERIALIZED_P(1);
1976 
1977  errorIfGeometryCollection(geom1,geom2);
1978  error_if_srid_mismatch(gserialized_get_srid(geom1), gserialized_get_srid(geom2));
1979 
1980  /* A.CoveredBy(Empty) == FALSE */
1981  if ( gserialized_is_empty(geom1) || gserialized_is_empty(geom2) )
1982  PG_RETURN_BOOL(false);
1983 
1984  /*
1985  * short-circuit 1: if geom1 bounding box is not completely inside
1986  * geom2 bounding box we can return FALSE.
1987  */
1988  if ( gserialized_get_gbox_p(geom1, &box1) &&
1989  gserialized_get_gbox_p(geom2, &box2) )
1990  {
1991  if ( ! gbox_contains_2d(&box2, &box1) )
1992  {
1993  PG_RETURN_BOOL(false);
1994  }
1995 
1996  POSTGIS_DEBUG(3, "bounding box short-circuit missed.");
1997  }
1998  /*
1999  * short-circuit 2: if geom1 is a point and geom2 is a polygon
2000  * call the point-in-polygon function.
2001  */
2002  if (is_point(geom1) && is_poly(geom2))
2003  {
2004  GSERIALIZED* gpoly = is_poly(geom1) ? geom1 : geom2;
2005  GSERIALIZED* gpoint = is_point(geom1) ? geom1 : geom2;
2006  RTREE_POLY_CACHE* cache = GetRtreeCache(fcinfo, gpoly);
2007  int retval;
2008 
2009  POSTGIS_DEBUG(3, "Point in Polygon test requested...short-circuiting.");
2010  if (gserialized_get_type(gpoint) == POINTTYPE)
2011  {
2012  LWGEOM* point = lwgeom_from_gserialized(gpoint);
2013  int pip_result = pip_short_circuit(cache, lwgeom_as_lwpoint(point), gpoly);
2014  lwgeom_free(point);
2015 
2016  retval = (pip_result != -1); /* not outside */
2017  }
2018  else if (gserialized_get_type(gpoint) == MULTIPOINTTYPE)
2019  {
2020  LWMPOINT* mpoint = lwgeom_as_lwmpoint(lwgeom_from_gserialized(gpoint));
2021  uint32_t i;
2022 
2023  retval = LW_TRUE;
2024  for (i = 0; i < mpoint->ngeoms; i++)
2025  {
2026  int pip_result = pip_short_circuit(cache, mpoint->geoms[i], gpoly);
2027  if (pip_result == -1)
2028  {
2029  retval = LW_FALSE;
2030  break;
2031  }
2032  }
2033 
2034  lwmpoint_free(mpoint);
2035  }
2036  else
2037  {
2038  /* Never get here */
2039  elog(ERROR,"Type isn't point or multipoint!");
2040  PG_RETURN_NULL();
2041  }
2042 
2043  PG_FREE_IF_COPY(geom1, 0);
2044  PG_FREE_IF_COPY(geom2, 1);
2045  PG_RETURN_BOOL(retval);
2046  }
2047  else
2048  {
2049  POSTGIS_DEBUGF(3, "CoveredBy: type1: %d, type2: %d", gserialized_get_type(geom1), gserialized_get_type(geom2));
2050  }
2051 
2052  initGEOS(lwpgnotice, lwgeom_geos_error);
2053 
2054  g1 = POSTGIS2GEOS(geom1);
2055 
2056  if (!g1)
2057  HANDLE_GEOS_ERROR("First argument geometry could not be converted to GEOS");
2058 
2059  g2 = POSTGIS2GEOS(geom2);
2060 
2061  if (!g2)
2062  {
2063  GEOSGeom_destroy(g1);
2064  HANDLE_GEOS_ERROR("Second argument geometry could not be converted to GEOS");
2065  }
2066 
2067  result = GEOSRelatePattern(g1,g2,patt);
2068 
2069  GEOSGeom_destroy(g1);
2070  GEOSGeom_destroy(g2);
2071 
2072  if (result == 2) HANDLE_GEOS_ERROR("GEOSCoveredBy");
2073 
2074  PG_FREE_IF_COPY(geom1, 0);
2075  PG_FREE_IF_COPY(geom2, 1);
2076 
2077  PG_RETURN_BOOL(result);
2078 }
2079 
2080 PG_FUNCTION_INFO_V1(crosses);
2081 Datum crosses(PG_FUNCTION_ARGS)
2082 {
2083  GSERIALIZED *geom1;
2084  GSERIALIZED *geom2;
2085  GEOSGeometry *g1, *g2;
2086  int result;
2087  GBOX box1, box2;
2088 
2089  geom1 = PG_GETARG_GSERIALIZED_P(0);
2090  geom2 = PG_GETARG_GSERIALIZED_P(1);
2091 
2092  errorIfGeometryCollection(geom1,geom2);
2093  error_if_srid_mismatch(gserialized_get_srid(geom1), gserialized_get_srid(geom2));
2094 
2095  /* A.Crosses(Empty) == FALSE */
2096  if ( gserialized_is_empty(geom1) || gserialized_is_empty(geom2) )
2097  PG_RETURN_BOOL(false);
2098 
2099  /*
2100  * short-circuit 1: if geom2 bounding box does not overlap
2101  * geom1 bounding box we can return FALSE.
2102  */
2103  if ( gserialized_get_gbox_p(geom1, &box1) &&
2104  gserialized_get_gbox_p(geom2, &box2) )
2105  {
2106  if ( gbox_overlaps_2d(&box1, &box2) == LW_FALSE )
2107  {
2108  PG_RETURN_BOOL(false);
2109  }
2110  }
2111 
2112  initGEOS(lwpgnotice, lwgeom_geos_error);
2113 
2114  g1 = POSTGIS2GEOS(geom1);
2115  if (!g1)
2116  HANDLE_GEOS_ERROR("First argument geometry could not be converted to GEOS");
2117 
2118  g2 = POSTGIS2GEOS(geom2);
2119  if (!g2)
2120  {
2121  GEOSGeom_destroy(g1);
2122  HANDLE_GEOS_ERROR("Second argument geometry could not be converted to GEOS");
2123  }
2124 
2125  result = GEOSCrosses(g1,g2);
2126 
2127  GEOSGeom_destroy(g1);
2128  GEOSGeom_destroy(g2);
2129 
2130  if (result == 2) HANDLE_GEOS_ERROR("GEOSCrosses");
2131 
2132  PG_FREE_IF_COPY(geom1, 0);
2133  PG_FREE_IF_COPY(geom2, 1);
2134 
2135  PG_RETURN_BOOL(result);
2136 }
2137 
2138 PG_FUNCTION_INFO_V1(geos_intersects);
2139 Datum geos_intersects(PG_FUNCTION_ARGS)
2140 {
2141  GSERIALIZED *geom1;
2142  GSERIALIZED *geom2;
2143  int result;
2144  GBOX box1, box2;
2145  PrepGeomCache *prep_cache;
2146 
2147  geom1 = PG_GETARG_GSERIALIZED_P(0);
2148  geom2 = PG_GETARG_GSERIALIZED_P(1);
2149 
2150  error_if_srid_mismatch(gserialized_get_srid(geom1), gserialized_get_srid(geom2));
2151 
2152  /* A.Intersects(Empty) == FALSE */
2153  if ( gserialized_is_empty(geom1) || gserialized_is_empty(geom2) )
2154  PG_RETURN_BOOL(false);
2155 
2156  /*
2157  * short-circuit 1: if geom2 bounding box does not overlap
2158  * geom1 bounding box we can return FALSE.
2159  */
2160  if ( gserialized_get_gbox_p(geom1, &box1) &&
2161  gserialized_get_gbox_p(geom2, &box2) )
2162  {
2163  if ( gbox_overlaps_2d(&box1, &box2) == LW_FALSE )
2164  {
2165  PG_RETURN_BOOL(false);
2166  }
2167  }
2168 
2169  /*
2170  * short-circuit 2: if the geoms are a point and a polygon,
2171  * call the point_outside_polygon function.
2172  */
2173  if ((is_point(geom1) && is_poly(geom2)) || (is_poly(geom1) && is_point(geom2)))
2174  {
2175  GSERIALIZED* gpoly = is_poly(geom1) ? geom1 : geom2;
2176  GSERIALIZED* gpoint = is_point(geom1) ? geom1 : geom2;
2177  RTREE_POLY_CACHE* cache = GetRtreeCache(fcinfo, gpoly);
2178  int retval;
2179 
2180  POSTGIS_DEBUG(3, "Point in Polygon test requested...short-circuiting.");
2181  if (gserialized_get_type(gpoint) == POINTTYPE)
2182  {
2183  LWGEOM* point = lwgeom_from_gserialized(gpoint);
2184  int pip_result = pip_short_circuit(cache, lwgeom_as_lwpoint(point), gpoly);
2185  lwgeom_free(point);
2186 
2187  retval = (pip_result != -1); /* not outside */
2188  }
2189  else if (gserialized_get_type(gpoint) == MULTIPOINTTYPE)
2190  {
2191  LWMPOINT* mpoint = lwgeom_as_lwmpoint(lwgeom_from_gserialized(gpoint));
2192  uint32_t i;
2193 
2194  retval = LW_FALSE;
2195  for (i = 0; i < mpoint->ngeoms; i++)
2196  {
2197  int pip_result = pip_short_circuit(cache, mpoint->geoms[i], gpoly);
2198  if (pip_result != -1) /* not outside */
2199  {
2200  retval = LW_TRUE;
2201  break;
2202  }
2203  }
2204 
2205  lwmpoint_free(mpoint);
2206  }
2207  else
2208  {
2209  /* Never get here */
2210  elog(ERROR,"Type isn't point or multipoint!");
2211  PG_RETURN_NULL();
2212  }
2213 
2214  PG_FREE_IF_COPY(geom1, 0);
2215  PG_FREE_IF_COPY(geom2, 1);
2216  PG_RETURN_BOOL(retval);
2217  }
2218 
2219  initGEOS(lwpgnotice, lwgeom_geos_error);
2220  prep_cache = GetPrepGeomCache( fcinfo, geom1, geom2 );
2221 
2222  if ( prep_cache && prep_cache->prepared_geom )
2223  {
2224  if ( prep_cache->gcache.argnum == 1 )
2225  {
2226  GEOSGeometry *g = POSTGIS2GEOS(geom2);
2227  if (!g) HANDLE_GEOS_ERROR("Geometry could not be converted to GEOS");
2228  result = GEOSPreparedIntersects( prep_cache->prepared_geom, g);
2229  GEOSGeom_destroy(g);
2230  }
2231  else
2232  {
2233  GEOSGeometry *g = POSTGIS2GEOS(geom1);
2234  if (!g)
2235  HANDLE_GEOS_ERROR("Geometry could not be converted to GEOS");
2236  result = GEOSPreparedIntersects( prep_cache->prepared_geom, g);
2237  GEOSGeom_destroy(g);
2238  }
2239  }
2240  else
2241  {
2242  GEOSGeometry *g1;
2243  GEOSGeometry *g2;
2244  g1 = POSTGIS2GEOS(geom1);
2245  if (!g1) HANDLE_GEOS_ERROR("First argument geometry could not be converted to GEOS");
2246  g2 = POSTGIS2GEOS(geom2);
2247  if (!g2)
2248  {
2249  GEOSGeom_destroy(g1);
2250  HANDLE_GEOS_ERROR("Second argument geometry could not be converted to GEOS");
2251  }
2252  result = GEOSIntersects( g1, g2);
2253  GEOSGeom_destroy(g1);
2254  GEOSGeom_destroy(g2);
2255  }
2256 
2257  if (result == 2) HANDLE_GEOS_ERROR("GEOSIntersects");
2258 
2259  PG_FREE_IF_COPY(geom1, 0);
2260  PG_FREE_IF_COPY(geom2, 1);
2261 
2262  PG_RETURN_BOOL(result);
2263 }
2264 
2265 
2266 PG_FUNCTION_INFO_V1(touches);
2267 Datum touches(PG_FUNCTION_ARGS)
2268 {
2269  GSERIALIZED *geom1;
2270  GSERIALIZED *geom2;
2271  GEOSGeometry *g1, *g2;
2272  char result;
2273  GBOX box1, box2;
2274 
2275  geom1 = PG_GETARG_GSERIALIZED_P(0);
2276  geom2 = PG_GETARG_GSERIALIZED_P(1);
2277 
2278  errorIfGeometryCollection(geom1,geom2);
2279  error_if_srid_mismatch(gserialized_get_srid(geom1), gserialized_get_srid(geom2));
2280 
2281  /* A.Touches(Empty) == FALSE */
2282  if ( gserialized_is_empty(geom1) || gserialized_is_empty(geom2) )
2283  PG_RETURN_BOOL(false);
2284 
2285  /*
2286  * short-circuit 1: if geom2 bounding box does not overlap
2287  * geom1 bounding box we can return FALSE.
2288  */
2289  if ( gserialized_get_gbox_p(geom1, &box1) &&
2290  gserialized_get_gbox_p(geom2, &box2) )
2291  {
2292  if ( gbox_overlaps_2d(&box1, &box2) == LW_FALSE )
2293  {
2294  PG_RETURN_BOOL(false);
2295  }
2296  }
2297 
2298  initGEOS(lwpgnotice, lwgeom_geos_error);
2299 
2300  g1 = POSTGIS2GEOS(geom1 );
2301  if (!g1)
2302  HANDLE_GEOS_ERROR("First argument geometry could not be converted to GEOS");
2303 
2304  g2 = POSTGIS2GEOS(geom2 );
2305  if (!g2)
2306  {
2307  GEOSGeom_destroy(g1);
2308  HANDLE_GEOS_ERROR("Second argument geometry could not be converted to GEOS");
2309  }
2310 
2311  result = GEOSTouches(g1,g2);
2312 
2313  GEOSGeom_destroy(g1);
2314  GEOSGeom_destroy(g2);
2315 
2316  if (result == 2) HANDLE_GEOS_ERROR("GEOSTouches");
2317 
2318  PG_FREE_IF_COPY(geom1, 0);
2319  PG_FREE_IF_COPY(geom2, 1);
2320 
2321  PG_RETURN_BOOL(result);
2322 }
2323 
2324 
2325 PG_FUNCTION_INFO_V1(disjoint);
2326 Datum disjoint(PG_FUNCTION_ARGS)
2327 {
2328  GSERIALIZED *geom1;
2329  GSERIALIZED *geom2;
2330  GEOSGeometry *g1, *g2;
2331  char result;
2332  GBOX box1, box2;
2333 
2334  geom1 = PG_GETARG_GSERIALIZED_P(0);
2335  geom2 = PG_GETARG_GSERIALIZED_P(1);
2336 
2337  errorIfGeometryCollection(geom1,geom2);
2338  error_if_srid_mismatch(gserialized_get_srid(geom1), gserialized_get_srid(geom2));
2339 
2340  /* A.Disjoint(Empty) == TRUE */
2341  if ( gserialized_is_empty(geom1) || gserialized_is_empty(geom2) )
2342  PG_RETURN_BOOL(true);
2343 
2344  /*
2345  * short-circuit 1: if geom2 bounding box does not overlap
2346  * geom1 bounding box we can return TRUE.
2347  */
2348  if ( gserialized_get_gbox_p(geom1, &box1) &&
2349  gserialized_get_gbox_p(geom2, &box2) )
2350  {
2351  if ( gbox_overlaps_2d(&box1, &box2) == LW_FALSE )
2352  {
2353  PG_RETURN_BOOL(true);
2354  }
2355  }
2356 
2357  initGEOS(lwpgnotice, lwgeom_geos_error);
2358 
2359  g1 = POSTGIS2GEOS(geom1);
2360  if (!g1)
2361  HANDLE_GEOS_ERROR("First argument geometry could not be converted to GEOS");
2362 
2363  g2 = POSTGIS2GEOS(geom2);
2364  if (!g2)
2365  {
2366  GEOSGeom_destroy(g1);
2367  HANDLE_GEOS_ERROR("Second argument geometry could not be converted to GEOS");
2368  }
2369 
2370  result = GEOSDisjoint(g1,g2);
2371 
2372  GEOSGeom_destroy(g1);
2373  GEOSGeom_destroy(g2);
2374 
2375  if (result == 2) HANDLE_GEOS_ERROR("GEOSDisjoint");
2376 
2377  PG_FREE_IF_COPY(geom1, 0);
2378  PG_FREE_IF_COPY(geom2, 1);
2379 
2380  PG_RETURN_BOOL(result);
2381 }
2382 
2383 
2384 PG_FUNCTION_INFO_V1(relate_pattern);
2385 Datum relate_pattern(PG_FUNCTION_ARGS)
2386 {
2387  GSERIALIZED *geom1;
2388  GSERIALIZED *geom2;
2389  char *patt;
2390  char result;
2391  GEOSGeometry *g1, *g2;
2392  size_t i;
2393 
2394  geom1 = PG_GETARG_GSERIALIZED_P(0);
2395  geom2 = PG_GETARG_GSERIALIZED_P(1);
2396 
2397  /* TODO handle empty */
2398 
2399  errorIfGeometryCollection(geom1,geom2);
2400  error_if_srid_mismatch(gserialized_get_srid(geom1), gserialized_get_srid(geom2));
2401 
2402  initGEOS(lwpgnotice, lwgeom_geos_error);
2403 
2404  g1 = POSTGIS2GEOS(geom1);
2405  if (!g1)
2406  HANDLE_GEOS_ERROR("First argument geometry could not be converted to GEOS");
2407  g2 = POSTGIS2GEOS(geom2);
2408  if (!g2)
2409  {
2410  GEOSGeom_destroy(g1);
2411  HANDLE_GEOS_ERROR("Second argument geometry could not be converted to GEOS");
2412  }
2413 
2414  patt = DatumGetCString(DirectFunctionCall1(textout,
2415  PointerGetDatum(PG_GETARG_DATUM(2))));
2416 
2417  /*
2418  ** Need to make sure 't' and 'f' are upper-case before handing to GEOS
2419  */
2420  for ( i = 0; i < strlen(patt); i++ )
2421  {
2422  if ( patt[i] == 't' ) patt[i] = 'T';
2423  if ( patt[i] == 'f' ) patt[i] = 'F';
2424  }
2425 
2426  result = GEOSRelatePattern(g1,g2,patt);
2427  GEOSGeom_destroy(g1);
2428  GEOSGeom_destroy(g2);
2429  pfree(patt);
2430 
2431  if (result == 2) HANDLE_GEOS_ERROR("GEOSRelatePattern");
2432 
2433  PG_FREE_IF_COPY(geom1, 0);
2434  PG_FREE_IF_COPY(geom2, 1);
2435 
2436  PG_RETURN_BOOL(result);
2437 }
2438 
2439 
2440 
2441 PG_FUNCTION_INFO_V1(relate_full);
2442 Datum relate_full(PG_FUNCTION_ARGS)
2443 {
2444  GSERIALIZED *geom1;
2445  GSERIALIZED *geom2;
2446  GEOSGeometry *g1, *g2;
2447  char *relate_str;
2448  text *result;
2449  int bnr = GEOSRELATE_BNR_OGC;
2450 
2451  POSTGIS_DEBUG(2, "in relate_full()");
2452 
2453  /* TODO handle empty */
2454 
2455  geom1 = PG_GETARG_GSERIALIZED_P(0);
2456  geom2 = PG_GETARG_GSERIALIZED_P(1);
2457 
2458  if ( PG_NARGS() > 2 )
2459  {
2460  bnr = PG_GETARG_INT32(2);
2461  }
2462 
2463  errorIfGeometryCollection(geom1,geom2);
2464  error_if_srid_mismatch(gserialized_get_srid(geom1), gserialized_get_srid(geom2));
2465 
2466  initGEOS(lwpgnotice, lwgeom_geos_error);
2467 
2468  g1 = POSTGIS2GEOS(geom1 );
2469  if (!g1)
2470  HANDLE_GEOS_ERROR("First argument geometry could not be converted to GEOS");
2471  g2 = POSTGIS2GEOS(geom2 );
2472  if (!g2)
2473  {
2474  GEOSGeom_destroy(g1);
2475  HANDLE_GEOS_ERROR("Second argument geometry could not be converted to GEOS");
2476  }
2477 
2478  POSTGIS_DEBUG(3, "constructed geometries ");
2479 
2480  POSTGIS_DEBUGF(3, "%s", GEOSGeomToWKT(g1));
2481  POSTGIS_DEBUGF(3, "%s", GEOSGeomToWKT(g2));
2482 
2483  relate_str = GEOSRelateBoundaryNodeRule(g1, g2, bnr);
2484 
2485  GEOSGeom_destroy(g1);
2486  GEOSGeom_destroy(g2);
2487 
2488  if (!relate_str) HANDLE_GEOS_ERROR("GEOSRelate");
2489 
2490  result = cstring_to_text(relate_str);
2491  GEOSFree(relate_str);
2492 
2493  PG_FREE_IF_COPY(geom1, 0);
2494  PG_FREE_IF_COPY(geom2, 1);
2495 
2496  PG_RETURN_TEXT_P(result);
2497 }
2498 
2499 
2500 PG_FUNCTION_INFO_V1(ST_Equals);
2501 Datum ST_Equals(PG_FUNCTION_ARGS)
2502 {
2503  GSERIALIZED *geom1;
2504  GSERIALIZED *geom2;
2505  GEOSGeometry *g1, *g2;
2506  char result;
2507  GBOX box1, box2;
2508 
2509  geom1 = PG_GETARG_GSERIALIZED_P(0);
2510  geom2 = PG_GETARG_GSERIALIZED_P(1);
2511 
2512  errorIfGeometryCollection(geom1,geom2);
2513  error_if_srid_mismatch(gserialized_get_srid(geom1), gserialized_get_srid(geom2));
2514 
2515  /* Empty == Empty */
2516  if ( gserialized_is_empty(geom1) && gserialized_is_empty(geom2) )
2517  PG_RETURN_BOOL(true);
2518 
2519  /*
2520  * short-circuit: If geom1 and geom2 do not have the same bounding box
2521  * we can return FALSE.
2522  */
2523  if ( gserialized_get_gbox_p(geom1, &box1) &&
2524  gserialized_get_gbox_p(geom2, &box2) )
2525  {
2526  if ( gbox_same_2d_float(&box1, &box2) == LW_FALSE )
2527  {
2528  PG_RETURN_BOOL(false);
2529  }
2530  }
2531 
2532  /*
2533  * short-circuit: if geom1 and geom2 are binary-equivalent, we can return
2534  * TRUE. This is much faster than doing the comparison using GEOS.
2535  */
2536  if (VARSIZE(geom1) == VARSIZE(geom2) && !memcmp(geom1, geom2, VARSIZE(geom1))) {
2537  PG_RETURN_BOOL(true);
2538  }
2539 
2540  initGEOS(lwpgnotice, lwgeom_geos_error);
2541 
2542  g1 = POSTGIS2GEOS(geom1);
2543 
2544  if (!g1)
2545  HANDLE_GEOS_ERROR("First argument geometry could not be converted to GEOS");
2546 
2547  g2 = POSTGIS2GEOS(geom2);
2548 
2549  if (!g2)
2550  {
2551  GEOSGeom_destroy(g1);
2552  HANDLE_GEOS_ERROR("Second argument geometry could not be converted to GEOS");
2553  }
2554 
2555  result = GEOSEquals(g1,g2);
2556 
2557  GEOSGeom_destroy(g1);
2558  GEOSGeom_destroy(g2);
2559 
2560  if (result == 2) HANDLE_GEOS_ERROR("GEOSEquals");
2561 
2562  PG_FREE_IF_COPY(geom1, 0);
2563  PG_FREE_IF_COPY(geom2, 1);
2564 
2565  PG_RETURN_BOOL(result);
2566 }
2567 
2568 PG_FUNCTION_INFO_V1(issimple);
2569 Datum issimple(PG_FUNCTION_ARGS)
2570 {
2571  GSERIALIZED *geom;
2572  LWGEOM *lwgeom_in;
2573  int result;
2574 
2575  POSTGIS_DEBUG(2, "issimple called");
2576 
2577  geom = PG_GETARG_GSERIALIZED_P(0);
2578 
2579  if ( gserialized_is_empty(geom) )
2580  PG_RETURN_BOOL(true);
2581 
2582  lwgeom_in = lwgeom_from_gserialized(geom);
2583  result = lwgeom_is_simple(lwgeom_in);
2584  lwgeom_free(lwgeom_in) ;
2585  PG_FREE_IF_COPY(geom, 0);
2586 
2587  if (result == -1) {
2588  PG_RETURN_NULL(); /*never get here */
2589  }
2590 
2591  PG_RETURN_BOOL(result);
2592 }
2593 
2594 PG_FUNCTION_INFO_V1(isring);
2595 Datum isring(PG_FUNCTION_ARGS)
2596 {
2597  GSERIALIZED *geom;
2598  GEOSGeometry *g1;
2599  int result;
2600 
2601  geom = PG_GETARG_GSERIALIZED_P(0);
2602 
2603  /* Empty things can't close */
2604  if ( gserialized_is_empty(geom) )
2605  PG_RETURN_BOOL(false);
2606 
2607  initGEOS(lwpgnotice, lwgeom_geos_error);
2608 
2609  g1 = POSTGIS2GEOS(geom);
2610  if (!g1)
2611  HANDLE_GEOS_ERROR("First argument geometry could not be converted to GEOS");
2612 
2613  if ( GEOSGeomTypeId(g1) != GEOS_LINESTRING )
2614  {
2615  GEOSGeom_destroy(g1);
2616  elog(ERROR, "ST_IsRing() should only be called on a linear feature");
2617  }
2618 
2619  result = GEOSisRing(g1);
2620  GEOSGeom_destroy(g1);
2621 
2622  if (result == 2) HANDLE_GEOS_ERROR("GEOSisRing");
2623 
2624  PG_FREE_IF_COPY(geom, 0);
2625  PG_RETURN_BOOL(result);
2626 }
2627 
2628 GSERIALIZED *
2629 GEOS2POSTGIS(GEOSGeom geom, char want3d)
2630 {
2631  LWGEOM *lwgeom;
2632  GSERIALIZED *result;
2633 
2634  lwgeom = GEOS2LWGEOM(geom, want3d);
2635  if ( ! lwgeom )
2636  {
2637  lwpgerror("%s: GEOS2LWGEOM returned NULL", __func__);
2638  return NULL;
2639  }
2640 
2641  POSTGIS_DEBUGF(4, "%s: GEOS2LWGEOM returned a %s", __func__, lwgeom_summary(lwgeom, 0));
2642 
2643  if (lwgeom_needs_bbox(lwgeom)) lwgeom_add_bbox(lwgeom);
2644 
2645  result = geometry_serialize(lwgeom);
2646  lwgeom_free(lwgeom);
2647 
2648  return result;
2649 }
2650 
2651 /*-----=POSTGIS2GEOS= */
2652 
2653 
2654 GEOSGeometry *
2655 POSTGIS2GEOS(GSERIALIZED *pglwgeom)
2656 {
2657  GEOSGeometry *ret;
2658  LWGEOM *lwgeom = lwgeom_from_gserialized(pglwgeom);
2659  if ( ! lwgeom )
2660  {
2661  lwpgerror("POSTGIS2GEOS: unable to deserialize input");
2662  return NULL;
2663  }
2664  ret = LWGEOM2GEOS(lwgeom, 0);
2665  lwgeom_free(lwgeom);
2666 
2667  return ret;
2668 }
2669 
2670 uint32_t array_nelems_not_null(ArrayType* array) {
2671  ArrayIterator iterator;
2672  Datum value;
2673  bool isnull;
2674  uint32_t nelems_not_null = 0;
2675 
2676 #if POSTGIS_PGSQL_VERSION >= 95
2677  iterator = array_create_iterator(array, 0, NULL);
2678 #else
2679  iterator = array_create_iterator(array, 0);
2680 #endif
2681  while(array_iterate(iterator, &value, &isnull) )
2682  {
2683  if (!isnull)
2684  {
2685  nelems_not_null++;
2686  }
2687  }
2688  array_free_iterator(iterator);
2689 
2690  return nelems_not_null;
2691 }
2692 
2693 /* ARRAY2LWGEOM: Converts the non-null elements of a Postgres array into a LWGEOM* array */
2694 LWGEOM** ARRAY2LWGEOM(ArrayType* array, uint32_t nelems, int* is3d, int* srid)
2695 {
2696  ArrayIterator iterator;
2697  Datum value;
2698  bool isnull;
2699  bool gotsrid = false;
2700  uint32_t i = 0;
2701 
2702  LWGEOM** lw_geoms = palloc(nelems * sizeof(LWGEOM*));
2703 
2704 #if POSTGIS_PGSQL_VERSION >= 95
2705  iterator = array_create_iterator(array, 0, NULL);
2706 #else
2707  iterator = array_create_iterator(array, 0);
2708 #endif
2709 
2710  while(array_iterate(iterator, &value, &isnull))
2711  {
2712  GSERIALIZED *geom = (GSERIALIZED*) DatumGetPointer(value);
2713 
2714  if (isnull)
2715  {
2716  continue;
2717  }
2718 
2719  *is3d = *is3d || gserialized_has_z(geom);
2720 
2721  lw_geoms[i] = lwgeom_from_gserialized(geom);
2722  if (!lw_geoms[i]) /* error in creation */
2723  {
2724  lwpgerror("Geometry deserializing geometry");
2725  return NULL;
2726  }
2727  if (!gotsrid)
2728  {
2729  gotsrid = true;
2730  *srid = gserialized_get_srid(geom);
2731  }
2732  else if (*srid != gserialized_get_srid(geom))
2733  {
2734  error_if_srid_mismatch(*srid, gserialized_get_srid(geom));
2735  return NULL;
2736  }
2737 
2738  i++;
2739  }
2740 
2741  return lw_geoms;
2742 }
2743 
2744 /* ARRAY2GEOS: Converts the non-null elements of a Postgres array into a GEOSGeometry* array */
2745 GEOSGeometry** ARRAY2GEOS(ArrayType* array, uint32_t nelems, int* is3d, int* srid)
2746 {
2747  ArrayIterator iterator;
2748  Datum value;
2749  bool isnull;
2750  bool gotsrid = false;
2751  uint32_t i = 0;
2752 
2753  GEOSGeometry** geos_geoms = palloc(nelems * sizeof(GEOSGeometry*));
2754 
2755 #if POSTGIS_PGSQL_VERSION >= 95
2756  iterator = array_create_iterator(array, 0, NULL);
2757 #else
2758  iterator = array_create_iterator(array, 0);
2759 #endif
2760 
2761  while(array_iterate(iterator, &value, &isnull))
2762  {
2763  GSERIALIZED *geom = (GSERIALIZED*) DatumGetPointer(value);
2764 
2765  if (isnull)
2766  {
2767  continue;
2768  }
2769 
2770  *is3d = *is3d || gserialized_has_z(geom);
2771 
2772  geos_geoms[i] = POSTGIS2GEOS(geom);
2773  if (!geos_geoms[i])
2774  {
2775  uint32_t j;
2776  lwpgerror("Geometry could not be converted to GEOS");
2777 
2778  for (j = 0; j < i; j++) {
2779  GEOSGeom_destroy(geos_geoms[j]);
2780  }
2781  return NULL;
2782  }
2783 
2784  if (!gotsrid)
2785  {
2786  *srid = gserialized_get_srid(geom);
2787  gotsrid = true;
2788  }
2789  else if (*srid != gserialized_get_srid(geom))
2790  {
2791  uint32_t j;
2792  error_if_srid_mismatch(*srid, gserialized_get_srid(geom));
2793 
2794  for (j = 0; j <= i; j++) {
2795  GEOSGeom_destroy(geos_geoms[j]);
2796  }
2797  return NULL;
2798  }
2799 
2800  i++;
2801  }
2802 
2803  array_free_iterator(iterator);
2804  return geos_geoms;
2805 }
2806 
2807 PG_FUNCTION_INFO_V1(GEOSnoop);
2808 Datum GEOSnoop(PG_FUNCTION_ARGS)
2809 {
2810  GSERIALIZED *geom;
2811  GEOSGeometry *geosgeom;
2812  GSERIALIZED *lwgeom_result;
2813 
2814  initGEOS(lwpgnotice, lwgeom_geos_error);
2815 
2816  geom = PG_GETARG_GSERIALIZED_P(0);
2817  geosgeom = POSTGIS2GEOS(geom);
2818  if ( ! geosgeom ) PG_RETURN_NULL();
2819 
2820  lwgeom_result = GEOS2POSTGIS(geosgeom, gserialized_has_z(geom));
2821  GEOSGeom_destroy(geosgeom);
2822 
2823  PG_FREE_IF_COPY(geom, 0);
2824 
2825  PG_RETURN_POINTER(lwgeom_result);
2826 }
2827 
2828 PG_FUNCTION_INFO_V1(polygonize_garray);
2829 Datum polygonize_garray(PG_FUNCTION_ARGS)
2830 {
2831  ArrayType *array;
2832  int is3d = 0;
2833  uint32 nelems, i;
2834  GSERIALIZED *result;
2835  GEOSGeometry *geos_result;
2836  const GEOSGeometry **vgeoms;
2837  int srid=SRID_UNKNOWN;
2838 #if POSTGIS_DEBUG_LEVEL >= 3
2839  static int call=1;
2840 #endif
2841 
2842 #if POSTGIS_DEBUG_LEVEL >= 3
2843  call++;
2844 #endif
2845 
2846  if (PG_ARGISNULL(0))
2847  PG_RETURN_NULL();
2848 
2849  array = PG_GETARG_ARRAYTYPE_P(0);
2850  nelems = array_nelems_not_null(array);
2851 
2852  if (nelems == 0)
2853  PG_RETURN_NULL();
2854 
2855  POSTGIS_DEBUGF(3, "polygonize_garray: number of non-null elements: %d", nelems);
2856 
2857  /* Ok, we really need geos now ;) */
2858  initGEOS(lwpgnotice, lwgeom_geos_error);
2859 
2860  vgeoms = (const GEOSGeometry**) ARRAY2GEOS(array, nelems, &is3d, &srid);
2861 
2862  POSTGIS_DEBUG(3, "polygonize_garray: invoking GEOSpolygonize");
2863 
2864  geos_result = GEOSPolygonize(vgeoms, nelems);
2865 
2866  POSTGIS_DEBUG(3, "polygonize_garray: GEOSpolygonize returned");
2867 
2868  for (i=0; i<nelems; ++i) GEOSGeom_destroy((GEOSGeometry *)vgeoms[i]);
2869  pfree(vgeoms);
2870 
2871  if ( ! geos_result ) PG_RETURN_NULL();
2872 
2873  GEOSSetSRID(geos_result, srid);
2874  result = GEOS2POSTGIS(geos_result, is3d);
2875  GEOSGeom_destroy(geos_result);
2876  if (!result)
2877  {
2878  elog(ERROR, "%s returned an error", __func__);
2879  PG_RETURN_NULL(); /*never get here */
2880  }
2881 
2882  PG_RETURN_POINTER(result);
2883 }
2884 
2885 
2886 PG_FUNCTION_INFO_V1(clusterintersecting_garray);
2887 Datum clusterintersecting_garray(PG_FUNCTION_ARGS)
2888 {
2889  Datum* result_array_data;
2890  ArrayType *array, *result;
2891  int is3d = 0;
2892  uint32 nelems, nclusters, i;
2893  GEOSGeometry **geos_inputs, **geos_results;
2894  int srid=SRID_UNKNOWN;
2895 
2896  /* Parameters used to construct a result array */
2897  int16 elmlen;
2898  bool elmbyval;
2899  char elmalign;
2900 
2901  /* Null array, null geometry (should be empty?) */
2902  if (PG_ARGISNULL(0))
2903  PG_RETURN_NULL();
2904 
2905  array = PG_GETARG_ARRAYTYPE_P(0);
2906  nelems = array_nelems_not_null(array);
2907 
2908  POSTGIS_DEBUGF(3, "clusterintersecting_garray: number of non-null elements: %d", nelems);
2909 
2910  if ( nelems == 0 ) PG_RETURN_NULL();
2911 
2912  /* TODO short-circuit for one element? */
2913 
2914  /* Ok, we really need geos now ;) */
2915  initGEOS(lwpgnotice, lwgeom_geos_error);
2916 
2917  geos_inputs = ARRAY2GEOS(array, nelems, &is3d, &srid);
2918  if(!geos_inputs)
2919  {
2920  PG_RETURN_NULL();
2921  }
2922 
2923  if (cluster_intersecting(geos_inputs, nelems, &geos_results, &nclusters) != LW_SUCCESS)
2924  {
2925  elog(ERROR, "clusterintersecting: Error performing clustering");
2926  PG_RETURN_NULL();
2927  }
2928  pfree(geos_inputs); /* don't need to destroy items because GeometryCollections have taken ownership */
2929 
2930  if (!geos_results) PG_RETURN_NULL();
2931 
2932  result_array_data = palloc(nclusters * sizeof(Datum));
2933  for (i=0; i<nclusters; ++i)
2934  {
2935  result_array_data[i] = PointerGetDatum(GEOS2POSTGIS(geos_results[i], is3d));
2936  GEOSGeom_destroy(geos_results[i]);
2937  }
2938  pfree(geos_results);
2939 
2940  get_typlenbyvalalign(array->elemtype, &elmlen, &elmbyval, &elmalign);
2941  result = construct_array(result_array_data, nclusters, array->elemtype, elmlen, elmbyval, elmalign);
2942 
2943  if (!result)
2944  {
2945  elog(ERROR, "clusterintersecting: Error constructing return-array");
2946  PG_RETURN_NULL();
2947  }
2948 
2949  PG_RETURN_POINTER(result);
2950 }
2951 
2952 PG_FUNCTION_INFO_V1(cluster_within_distance_garray);
2953 Datum cluster_within_distance_garray(PG_FUNCTION_ARGS)
2954 {
2955  Datum* result_array_data;
2956  ArrayType *array, *result;
2957  int is3d = 0;
2958  uint32 nelems, nclusters, i;
2959  LWGEOM** lw_inputs;
2960  LWGEOM** lw_results;
2961  double tolerance;
2962  int srid=SRID_UNKNOWN;
2963 
2964  /* Parameters used to construct a result array */
2965  int16 elmlen;
2966  bool elmbyval;
2967  char elmalign;
2968 
2969  /* Null array, null geometry (should be empty?) */
2970  if (PG_ARGISNULL(0))
2971  PG_RETURN_NULL();
2972 
2973  array = PG_GETARG_ARRAYTYPE_P(0);
2974 
2975  tolerance = PG_GETARG_FLOAT8(1);
2976  if (tolerance < 0)
2977  {
2978  lwpgerror("Tolerance must be a positive number.");
2979  PG_RETURN_NULL();
2980  }
2981 
2982  nelems = array_nelems_not_null(array);
2983 
2984  POSTGIS_DEBUGF(3, "cluster_within_distance_garray: number of non-null elements: %d", nelems);
2985 
2986  if ( nelems == 0 ) PG_RETURN_NULL();
2987 
2988  /* TODO short-circuit for one element? */
2989 
2990  /* Ok, we really need geos now ;) */
2991  initGEOS(lwpgnotice, lwgeom_geos_error);
2992 
2993  lw_inputs = ARRAY2LWGEOM(array, nelems, &is3d, &srid);
2994  if (!lw_inputs)
2995  {
2996  PG_RETURN_NULL();
2997  }
2998 
2999  if (cluster_within_distance(lw_inputs, nelems, tolerance, &lw_results, &nclusters) != LW_SUCCESS)
3000  {
3001  elog(ERROR, "cluster_within: Error performing clustering");
3002  PG_RETURN_NULL();
3003  }
3004  pfree(lw_inputs); /* don't need to destroy items because GeometryCollections have taken ownership */
3005 
3006  if (!lw_results) PG_RETURN_NULL();
3007 
3008  result_array_data = palloc(nclusters * sizeof(Datum));
3009  for (i=0; i<nclusters; ++i)
3010  {
3011  result_array_data[i] = PointerGetDatum(gserialized_from_lwgeom(lw_results[i], NULL));
3012  lwgeom_free(lw_results[i]);
3013  }
3014  pfree(lw_results);
3015 
3016  get_typlenbyvalalign(array->elemtype, &elmlen, &elmbyval, &elmalign);
3017  result = construct_array(result_array_data, nclusters, array->elemtype, elmlen, elmbyval, elmalign);
3018 
3019  if (!result)
3020  {
3021  elog(ERROR, "clusterwithin: Error constructing return-array");
3022  PG_RETURN_NULL();
3023  }
3024 
3025  PG_RETURN_POINTER(result);
3026 }
3027 
3028 PG_FUNCTION_INFO_V1(linemerge);
3029 Datum linemerge(PG_FUNCTION_ARGS)
3030 {
3031  GSERIALIZED *geom1;
3032  GSERIALIZED *result;
3033  LWGEOM *lwgeom1, *lwresult ;
3034 
3035  geom1 = PG_GETARG_GSERIALIZED_P(0);
3036 
3037 
3038  lwgeom1 = lwgeom_from_gserialized(geom1) ;
3039 
3040  lwresult = lwgeom_linemerge(lwgeom1);
3041  result = geometry_serialize(lwresult) ;
3042 
3043  lwgeom_free(lwgeom1) ;
3044  lwgeom_free(lwresult) ;
3045 
3046  PG_FREE_IF_COPY(geom1, 0);
3047 
3048  PG_RETURN_POINTER(result);
3049 }
3050 
3051 /*
3052  * Take a geometry and return an areal geometry
3053  * (Polygon or MultiPolygon).
3054  * Actually a wrapper around GEOSpolygonize,
3055  * transforming the resulting collection into
3056  * a valid polygon Geometry.
3057  */
3058 PG_FUNCTION_INFO_V1(ST_BuildArea);
3059 Datum ST_BuildArea(PG_FUNCTION_ARGS)
3060 {
3061  GSERIALIZED *result;
3062  GSERIALIZED *geom;
3063  LWGEOM *lwgeom_in, *lwgeom_out;
3064 
3065  geom = PG_GETARG_GSERIALIZED_P(0);
3066  lwgeom_in = lwgeom_from_gserialized(geom);
3067 
3068  lwgeom_out = lwgeom_buildarea(lwgeom_in);
3069  lwgeom_free(lwgeom_in) ;
3070 
3071  if ( ! lwgeom_out )
3072  {
3073  PG_FREE_IF_COPY(geom, 0);
3074  PG_RETURN_NULL();
3075  }
3076 
3077  result = geometry_serialize(lwgeom_out) ;
3078  lwgeom_free(lwgeom_out) ;
3079 
3080  PG_FREE_IF_COPY(geom, 0);
3081  PG_RETURN_POINTER(result);
3082 }
3083 
3084 /*
3085  * Take the vertices of a geometry and builds
3086  * Delaunay triangles around them.
3087  */
3088 PG_FUNCTION_INFO_V1(ST_DelaunayTriangles);
3089 Datum ST_DelaunayTriangles(PG_FUNCTION_ARGS)
3090 {
3091  GSERIALIZED *result;
3092  GSERIALIZED *geom;
3093  LWGEOM *lwgeom_in, *lwgeom_out;
3094  double tolerance = 0.0;
3095  int flags = 0;
3096 
3097  geom = PG_GETARG_GSERIALIZED_P(0);
3098  tolerance = PG_GETARG_FLOAT8(1);
3099  flags = PG_GETARG_INT32(2);
3100 
3101  lwgeom_in = lwgeom_from_gserialized(geom);
3102  lwgeom_out = lwgeom_delaunay_triangulation(lwgeom_in, tolerance, flags);
3103  lwgeom_free(lwgeom_in) ;
3104 
3105  if ( ! lwgeom_out )
3106  {
3107  PG_FREE_IF_COPY(geom, 0);
3108  PG_RETURN_NULL();
3109  }
3110 
3111  result = geometry_serialize(lwgeom_out) ;
3112  lwgeom_free(lwgeom_out) ;
3113 
3114  PG_FREE_IF_COPY(geom, 0);
3115  PG_RETURN_POINTER(result);
3116 }
3117 
3118 /*
3119  * ST_Snap
3120  *
3121  * Snap a geometry to another with a given tolerance
3122  */
3123 Datum ST_Snap(PG_FUNCTION_ARGS);
3124 PG_FUNCTION_INFO_V1(ST_Snap);
3125 Datum ST_Snap(PG_FUNCTION_ARGS)
3126 {
3127  GSERIALIZED *geom1, *geom2, *result;
3128  LWGEOM *lwgeom1, *lwgeom2, *lwresult;
3129  double tolerance;
3130 
3131  geom1 = PG_GETARG_GSERIALIZED_P(0);
3132  geom2 = PG_GETARG_GSERIALIZED_P(1);
3133  tolerance = PG_GETARG_FLOAT8(2);
3134 
3135  lwgeom1 = lwgeom_from_gserialized(geom1);
3136  lwgeom2 = lwgeom_from_gserialized(geom2);
3137 
3138  lwresult = lwgeom_snap(lwgeom1, lwgeom2, tolerance);
3139  lwgeom_free(lwgeom1);
3140  lwgeom_free(lwgeom2);
3141  PG_FREE_IF_COPY(geom1, 0);
3142  PG_FREE_IF_COPY(geom2, 1);
3143 
3144  result = geometry_serialize(lwresult);
3145  lwgeom_free(lwresult);
3146 
3147  PG_RETURN_POINTER(result);
3148 }
3149 
3150 /*
3151  * ST_Split
3152  *
3153  * Split polygon by line, line by line, line by point.
3154  * Returns at most components as a collection.
3155  * First element of the collection is always the part which
3156  * remains after the cut, while the second element is the
3157  * part which has been cut out. We arbitrarily take the part
3158  * on the *right* of cut lines as the part which has been cut out.
3159  * For a line cut by a point the part which remains is the one
3160  * from start of the line to the cut point.
3161  *
3162  *
3163  * Author: Sandro Santilli <strk@kbt.io>
3164  *
3165  * Work done for Faunalia (http://www.faunalia.it) with fundings
3166  * from Regione Toscana - Sistema Informativo per il Governo
3167  * del Territorio e dell'Ambiente (RT-SIGTA).
3168  *
3169  * Thanks to the PostGIS community for sharing poly/line ideas [1]
3170  *
3171  * [1] http://trac.osgeo.org/postgis/wiki/UsersWikiSplitPolygonWithLineString
3172  *
3173  */
3174 Datum ST_Split(PG_FUNCTION_ARGS);
3175 PG_FUNCTION_INFO_V1(ST_Split);
3176 Datum ST_Split(PG_FUNCTION_ARGS)
3177 {
3178  GSERIALIZED *in, *blade_in, *out;
3179  LWGEOM *lwgeom_in, *lwblade_in, *lwgeom_out;
3180 
3181  in = PG_GETARG_GSERIALIZED_P(0);
3182  lwgeom_in = lwgeom_from_gserialized(in);
3183 
3184  blade_in = PG_GETARG_GSERIALIZED_P(1);
3185  lwblade_in = lwgeom_from_gserialized(blade_in);
3186 
3187  error_if_srid_mismatch(lwgeom_in->srid, lwblade_in->srid);
3188 
3189  lwgeom_out = lwgeom_split(lwgeom_in, lwblade_in);
3190  lwgeom_free(lwgeom_in);
3191  lwgeom_free(lwblade_in);
3192 
3193  if ( ! lwgeom_out )
3194  {
3195  PG_FREE_IF_COPY(in, 0); /* possibly referenced by lwgeom_out */
3196  PG_FREE_IF_COPY(blade_in, 1);
3197  PG_RETURN_NULL();
3198  }
3199 
3200  out = geometry_serialize(lwgeom_out);
3201  lwgeom_free(lwgeom_out);
3202  PG_FREE_IF_COPY(in, 0); /* possibly referenced by lwgeom_out */
3203  PG_FREE_IF_COPY(blade_in, 1);
3204 
3205  PG_RETURN_POINTER(out);
3206 }
3207 
3208 /**********************************************************************
3209  *
3210  * ST_SharedPaths
3211  *
3212  * Return the set of paths shared between two linear geometries,
3213  * and their direction (same or opposite).
3214  *
3215  * Developed by Sandro Santilli (strk@kbt.io) for Faunalia
3216  * (http://www.faunalia.it) with funding from Regione Toscana - Sistema
3217  * Informativo per la Gestione del Territorio e dell' Ambiente
3218  * [RT-SIGTA]". For the project: "Sviluppo strumenti software per il
3219  * trattamento di dati geografici basati su QuantumGIS e Postgis (CIG
3220  * 0494241492)"
3221  *
3222  **********************************************************************/
3223 Datum ST_SharedPaths(PG_FUNCTION_ARGS);
3224 PG_FUNCTION_INFO_V1(ST_SharedPaths);
3225 Datum ST_SharedPaths(PG_FUNCTION_ARGS)
3226 {
3227  GSERIALIZED *geom1, *geom2, *out;
3228  LWGEOM *g1, *g2, *lwgeom_out;
3229 
3230  geom1 = PG_GETARG_GSERIALIZED_P(0);
3231  geom2 = PG_GETARG_GSERIALIZED_P(1);
3232 
3233  g1 = lwgeom_from_gserialized(geom1);
3234  g2 = lwgeom_from_gserialized(geom2);
3235 
3236  lwgeom_out = lwgeom_sharedpaths(g1, g2);
3237  lwgeom_free(g1);
3238  lwgeom_free(g2);
3239 
3240  if ( ! lwgeom_out )
3241  {
3242  PG_FREE_IF_COPY(geom1, 0);
3243  PG_FREE_IF_COPY(geom2, 1);
3244  PG_RETURN_NULL();
3245  }
3246 
3247  out = geometry_serialize(lwgeom_out);
3248  lwgeom_free(lwgeom_out);
3249 
3250  PG_FREE_IF_COPY(geom1, 0);
3251  PG_FREE_IF_COPY(geom2, 1);
3252  PG_RETURN_POINTER(out);
3253 }
3254 
3255 
3256 /**********************************************************************
3257  *
3258  * ST_Node
3259  *
3260  * Fully node a set of lines using the least possible nodes while
3261  * preserving all of the input ones.
3262  *
3263  **********************************************************************/
3264 Datum ST_Node(PG_FUNCTION_ARGS);
3265 PG_FUNCTION_INFO_V1(ST_Node);
3266 Datum ST_Node(PG_FUNCTION_ARGS)
3267 {
3268  GSERIALIZED *geom1, *out;
3269  LWGEOM *g1, *lwgeom_out;
3270 
3271  geom1 = PG_GETARG_GSERIALIZED_P(0);
3272 
3273  g1 = lwgeom_from_gserialized(geom1);
3274 
3275  lwgeom_out = lwgeom_node(g1);
3276  lwgeom_free(g1);
3277 
3278  if ( ! lwgeom_out )
3279  {
3280  PG_FREE_IF_COPY(geom1, 0);
3281  PG_RETURN_NULL();
3282  }
3283 
3284  out = geometry_serialize(lwgeom_out);
3285  lwgeom_free(lwgeom_out);
3286 
3287  PG_FREE_IF_COPY(geom1, 0);
3288  PG_RETURN_POINTER(out);
3289 }
3290 
3291 /******************************************
3292  *
3293  * ST_Voronoi
3294  *
3295  * Returns a Voronoi diagram constructed
3296  * from the points of the input geometry.
3297  *
3298  ******************************************/
3299 Datum ST_Voronoi(PG_FUNCTION_ARGS);
3300 PG_FUNCTION_INFO_V1(ST_Voronoi);
3301 Datum ST_Voronoi(PG_FUNCTION_ARGS)
3302 {
3303  GSERIALIZED* input;
3304  GSERIALIZED* clip;
3305  GSERIALIZED* result;
3306  LWGEOM* lwgeom_input;
3307  LWGEOM* lwgeom_result;
3308  double tolerance;
3309  GBOX clip_envelope;
3310  int custom_clip_envelope;
3311  int return_polygons;
3312 
3313  /* Return NULL on NULL geometry */
3314  if (PG_ARGISNULL(0))
3315  PG_RETURN_NULL();
3316 
3317  /* Read our tolerance value */
3318  if (PG_ARGISNULL(2))
3319  {
3320  lwpgerror("Tolerance must be a positive number.");
3321  PG_RETURN_NULL();
3322  }
3323 
3324  tolerance = PG_GETARG_FLOAT8(2);
3325 
3326  if (tolerance < 0)
3327  {
3328  lwpgerror("Tolerance must be a positive number.");
3329  PG_RETURN_NULL();
3330  }
3331 
3332  /* Are we returning lines or polygons? */
3333  if (PG_ARGISNULL(3))
3334  {
3335  lwpgerror("return_polygons must be true or false.");
3336  PG_RETURN_NULL();
3337  }
3338  return_polygons = PG_GETARG_BOOL(3);
3339 
3340  /* Read our clipping envelope, if applicable. */
3341  custom_clip_envelope = !PG_ARGISNULL(1);
3342  if (custom_clip_envelope) {
3343  clip = PG_GETARG_GSERIALIZED_P(1);
3344  if (!gserialized_get_gbox_p(clip, &clip_envelope))
3345  {
3346  lwpgerror("Could not determine envelope of clipping geometry.");
3347  PG_FREE_IF_COPY(clip, 1);
3348  PG_RETURN_NULL();
3349  }
3350  PG_FREE_IF_COPY(clip, 1);
3351  }
3352 
3353  /* Read our input geometry */
3354  input = PG_GETARG_GSERIALIZED_P(0);
3355 
3356  lwgeom_input = lwgeom_from_gserialized(input);
3357 
3358  if(!lwgeom_input)
3359  {
3360  lwpgerror("Could not read input geometry.");
3361  PG_FREE_IF_COPY(input, 0);
3362  PG_RETURN_NULL();
3363  }
3364 
3365  lwgeom_result = lwgeom_voronoi_diagram(lwgeom_input, custom_clip_envelope ? &clip_envelope : NULL, tolerance, !return_polygons);
3366  lwgeom_free(lwgeom_input);
3367 
3368  if (!lwgeom_result)
3369  {
3370  lwpgerror("Error computing Voronoi diagram.");
3371  PG_FREE_IF_COPY(input, 0);
3372  PG_RETURN_NULL();
3373  }
3374 
3375  result = geometry_serialize(lwgeom_result);
3376  lwgeom_free(lwgeom_result);
3377 
3378  PG_FREE_IF_COPY(input, 0);
3379  PG_RETURN_POINTER(result);
3380 }
3381 
3382 /******************************************
3383  *
3384  * ST_MinimumClearance
3385  *
3386  * Returns the minimum clearance of a geometry.
3387  *
3388  ******************************************/
3389 Datum ST_MinimumClearance(PG_FUNCTION_ARGS);
3390 PG_FUNCTION_INFO_V1(ST_MinimumClearance);
3391 Datum ST_MinimumClearance(PG_FUNCTION_ARGS)
3392 {
3393 #if POSTGIS_GEOS_VERSION < 36
3394  lwpgerror("The GEOS version this PostGIS binary "
3395  "was compiled against (%d) doesn't support "
3396  "'ST_MinimumClearance' function (3.6.0+ required)",
3397  POSTGIS_GEOS_VERSION);
3398  PG_RETURN_NULL();
3399 #else /* POSTGIS_GEOS_VERSION >= 36 */
3400  GSERIALIZED* input;
3401  GEOSGeometry* input_geos;
3402  int error;
3403  double result;
3404 
3405  initGEOS(lwpgnotice, lwgeom_geos_error);
3406 
3407  input = PG_GETARG_GSERIALIZED_P(0);
3408  input_geos = POSTGIS2GEOS(input);
3409  if (!input_geos)
3410  HANDLE_GEOS_ERROR("Geometry could not be converted to GEOS");
3411 
3412  error = GEOSMinimumClearance(input_geos, &result);
3413  GEOSGeom_destroy(input_geos);
3414  if (error) HANDLE_GEOS_ERROR("Error computing minimum clearance");
3415 
3416  PG_FREE_IF_COPY(input, 0);
3417  PG_RETURN_FLOAT8(result);
3418 #endif
3419 }
3420 
3421 /******************************************
3422  *
3423  * ST_MinimumClearanceLine
3424  *
3425  * Returns the minimum clearance line of a geometry.
3426  *
3427  ******************************************/
3428 Datum ST_MinimumClearanceLine(PG_FUNCTION_ARGS);
3429 PG_FUNCTION_INFO_V1(ST_MinimumClearanceLine);
3430 Datum ST_MinimumClearanceLine(PG_FUNCTION_ARGS)
3431 {
3432 #if POSTGIS_GEOS_VERSION < 36
3433  lwpgerror("The GEOS version this PostGIS binary "
3434  "was compiled against (%d) doesn't support "
3435  "'ST_MinimumClearanceLine' function (3.6.0+ required)",
3436  POSTGIS_GEOS_VERSION);
3437  PG_RETURN_NULL();
3438 #else /* POSTGIS_GEOS_VERSION >= 36 */
3439  GSERIALIZED* input;
3440  GSERIALIZED* result;
3441  GEOSGeometry* input_geos;
3442  GEOSGeometry* result_geos;
3443  int srid;
3444 
3445  initGEOS(lwpgnotice, lwgeom_geos_error);
3446 
3447  input = PG_GETARG_GSERIALIZED_P(0);
3448  srid = gserialized_get_srid(input);
3449  input_geos = POSTGIS2GEOS(input);
3450  if (!input_geos)
3451  HANDLE_GEOS_ERROR("Geometry could not be converted to GEOS");
3452 
3453  result_geos = GEOSMinimumClearanceLine(input_geos);
3454  GEOSGeom_destroy(input_geos);
3455  if (!result_geos)
3456  HANDLE_GEOS_ERROR("Error computing minimum clearance");
3457 
3458  GEOSSetSRID(result_geos, srid);
3459  result = GEOS2POSTGIS(result_geos, LW_FALSE);
3460  GEOSGeom_destroy(result_geos);
3461 
3462  PG_FREE_IF_COPY(input, 0);
3463  PG_RETURN_POINTER(result);
3464 #endif
3465 }
3466 
3467 /******************************************
3468  *
3469  * ST_OrientedEnvelope
3470  *
3471  ******************************************/
3472 Datum ST_OrientedEnvelope(PG_FUNCTION_ARGS);
3473 PG_FUNCTION_INFO_V1(ST_OrientedEnvelope);
3474 Datum ST_OrientedEnvelope(PG_FUNCTION_ARGS)
3475 {
3476 #if POSTGIS_GEOS_VERSION < 36
3477  lwpgerror("The GEOS version this PostGIS binary "
3478  "was compiled against (%d) doesn't support "
3479  "'ST_OrientedEnvelope' function (3.6.0+ required)",
3480  POSTGIS_GEOS_VERSION);
3481  PG_RETURN_NULL();
3482 #else
3483  GSERIALIZED* input;
3484  GSERIALIZED* result;
3485  GEOSGeometry* input_geos;
3486  GEOSGeometry* result_geos;
3487  int srid;
3488 
3489  initGEOS(lwpgnotice, lwgeom_geos_error);
3490 
3491  input = PG_GETARG_GSERIALIZED_P(0);
3492  srid = gserialized_get_srid(input);
3493  input_geos = POSTGIS2GEOS(input);
3494  if (!input_geos)
3495  HANDLE_GEOS_ERROR("Geometry could not be converted to GEOS");
3496 
3497  result_geos = GEOSMinimumRotatedRectangle(input_geos);
3498  GEOSGeom_destroy(input_geos);
3499  if (!result_geos)
3500  HANDLE_GEOS_ERROR("Error computing oriented envelope");
3501 
3502  GEOSSetSRID(result_geos, srid);
3503  result = GEOS2POSTGIS(result_geos, LW_FALSE);
3504  GEOSGeom_destroy(result_geos);
3505 
3506  PG_FREE_IF_COPY(input, 0);
3507  PG_RETURN_POINTER(result);
3508 #endif
3509 }
gbox_overlaps_2d
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: g_box.c:330
gbox_copy
GBOX * gbox_copy(const GBOX *box)
Return a copy of the GBOX, based on dimensionality of flags.
Definition: g_box.c:433
gbox_same_2d_float
int gbox_same_2d_float(const GBOX *g1, const GBOX *g2)
Check if two given GBOX are the same in x and y, or would round to the same GBOX in x and if serializ...
Definition: g_box.c:194
gbox_contains_2d
int gbox_contains_2d(const GBOX *g1, const GBOX *g2)
Return LW_TRUE if the first GBOX contains the second on the 2d plane, LW_FALSE otherwise.
Definition: g_box.c:346
gserialized_get_srid
int32_t gserialized_get_srid(const GSERIALIZED *s)
Extract the SRID from the serialized form (it is packed into three bytes so this is a handy function)...
Definition: g_serialized.c:100
gserialized_has_z
int gserialized_has_z(const GSERIALIZED *gser)
Check if a GSERIALIZED has a Z ordinate.
Definition: g_serialized.c:45
lwgeom_from_gserialized
LWGEOM * lwgeom_from_gserialized(const GSERIALIZED *g)
Allocate a new LWGEOM from a GSERIALIZED.
Definition: g_serialized.c:1530
gserialized_is_empty
int gserialized_is_empty(const GSERIALIZED *g)
Check if a GSERIALIZED is empty without deserializing first.
Definition: g_serialized.c:179
gserialized_get_gbox_p
int gserialized_get_gbox_p(const GSERIALIZED *g, GBOX *box)
Read the bounding box off a serialization and calculate one if it is not already there.
Definition: g_serialized.c:640
gserialized_get_type
uint32_t gserialized_get_type(const GSERIALIZED *s)
Extract the geometry type from the serialized form (it hides in the anonymous data area,...
Definition: g_serialized.c:86
gserialized_from_lwgeom
GSERIALIZED * gserialized_from_lwgeom(LWGEOM *geom, size_t *size)
Allocate a new GSERIALIZED from an LWGEOM.
Definition: g_serialized.c:1176
lwgeom_geos_errmsg
char lwgeom_geos_errmsg[LWGEOM_GEOS_ERRMSG_MAXSIZE]
Definition: liblwgeom/lwgeom_geos.c:40
LWGEOM2GEOS
GEOSGeometry * LWGEOM2GEOS(const LWGEOM *lwgeom, uint8_t autofix)
Definition: liblwgeom/lwgeom_geos.c:390
GEOS2LWGEOM
LWGEOM * GEOS2LWGEOM(const GEOSGeometry *geom, uint8_t want3d)
Definition: liblwgeom/lwgeom_geos.c:166
lwgeom_geos_error
void lwgeom_geos_error(const char *fmt,...)
Definition: liblwgeom/lwgeom_geos.c:43
cluster_within_distance
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...
Definition: lwgeom_geos_cluster.c:526
cluster_intersecting
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...
Definition: lwgeom_geos_cluster.c:247
lwgeom_geos.h
lwgeom_centroid
LWGEOM * lwgeom_centroid(const LWGEOM *geom)
Definition: liblwgeom/lwgeom_geos.c:822
LW_FALSE
#define LW_FALSE
Definition: liblwgeom.h:77
WKT_SFSQL
#define WKT_SFSQL
Definition: liblwgeom.h:2076
COLLECTIONTYPE
#define COLLECTIONTYPE
Definition: liblwgeom.h:91
lwgeom_delaunay_triangulation
LWGEOM * lwgeom_delaunay_triangulation(const LWGEOM *geom, double tolerance, int32_t edgeOnly)
Take vertices of a geometry and build a delaunay triangulation on them.
Definition: liblwgeom/lwgeom_geos.c:1788
lwgeom_node
LWGEOM * lwgeom_node(const LWGEOM *lwgeom_in)
Definition: lwgeom_geos_node.c:130
lwmpoint_free
void lwmpoint_free(LWMPOINT *mpt)
Definition: lwmpoint.c:72
lwgeom_symdifference
LWGEOM * lwgeom_symdifference(const LWGEOM *geom1, const LWGEOM *geom2)
Definition: liblwgeom/lwgeom_geos.c:789
lwgeom_as_lwmpoint
LWMPOINT * lwgeom_as_lwmpoint(const LWGEOM *lwgeom)
Definition: lwgeom.c:233
lwgeom_free
void lwgeom_free(LWGEOM *geom)
Definition: lwgeom.c:1144
lwgeom_split
LWGEOM * lwgeom_split(const LWGEOM *lwgeom_in, const LWGEOM *blade_in)
Definition: lwgeom_geos_split.c:564
lwgeom_intersection
LWGEOM * lwgeom_intersection(const LWGEOM *geom1, const LWGEOM *geom2)
Definition: liblwgeom/lwgeom_geos.c:660
LINETYPE
#define LINETYPE
Definition: liblwgeom.h:86
lwgeom_offsetcurve
LWGEOM * lwgeom_offsetcurve(const LWGEOM *geom, double size, int quadsegs, int joinStyle, double mitreLimit)
Definition: liblwgeom/lwgeom_geos.c:1430
LW_SUCCESS
#define LW_SUCCESS
Definition: liblwgeom.h:80
lwpoly_as_lwgeom
LWGEOM * lwpoly_as_lwgeom(const LWPOLY *obj)
Definition: lwgeom.c:320
MULTIPOINTTYPE
#define MULTIPOINTTYPE
Definition: liblwgeom.h:88
lwmessage_truncate
char * lwmessage_truncate(char *str, int startpos, int endpos, int maxlength, int truncdirection)
Definition: lwutil.c:261
lwgeom_is_simple
int lwgeom_is_simple(const LWGEOM *lwgeom)
Definition: liblwgeom/lwgeom_geos.c:1246
lwgeom_snap
LWGEOM * lwgeom_snap(const LWGEOM *geom1, const LWGEOM *geom2, double tolerance)
Snap vertices and segments of a geometry to another using a given tolerance.
Definition: liblwgeom/lwgeom_geos.c:1296
lwgeom_get_type
uint32_t lwgeom_get_type(const LWGEOM *geom)
Return LWTYPE number.
Definition: lwgeom.c:923
lwgeom_difference
LWGEOM * lwgeom_difference(const LWGEOM *geom1, const LWGEOM *geom2)
Definition: liblwgeom/lwgeom_geos.c:756
lwgeom_geos_version
const char * lwgeom_geos_version(void)
Return GEOS version string (not to be freed)
Definition: liblwgeom/lwgeom_geos.c:597
lwgeom_as_lwmpoly
LWMPOLY * lwgeom_as_lwmpoly(const LWGEOM *lwgeom)
Definition: lwgeom.c:251
lwgeom_needs_bbox
int lwgeom_needs_bbox(const LWGEOM *geom)
Check whether or not a lwgeom is big enough to warrant a bounding box.
Definition: lwgeom.c:1197
POINTTYPE
#define POINTTYPE
LWTYPE numbers, used internally by PostGIS.
Definition: liblwgeom.h:85
lwgeom_summary
char * lwgeom_summary(const LWGEOM *lwgeom, int offset)
Definition: lwgeom_debug.c:158
lwgeom_pointonsurface
LWGEOM * lwgeom_pointonsurface(const LWGEOM *geom)
Definition: liblwgeom/lwgeom_geos.c:855
lwgeom_sharedpaths
LWGEOM * lwgeom_sharedpaths(const LWGEOM *geom1, const LWGEOM *geom2)
Definition: liblwgeom/lwgeom_geos.c:1323
lwgeom_as_lwpoint
LWPOINT * lwgeom_as_lwpoint(const LWGEOM *lwgeom)
Definition: lwgeom.c:161
TINTYPE
#define TINTYPE
Definition: liblwgeom.h:99
lwgeom_voronoi_diagram
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.
Definition: liblwgeom/lwgeom_geos.c:1872
MULTIPOLYGONTYPE
#define MULTIPOLYGONTYPE
Definition: liblwgeom.h:90
lwgeom_to_points
LWMPOINT * lwgeom_to_points(const LWGEOM *lwgeom, uint32_t npoints)
Definition: liblwgeom/lwgeom_geos.c:1705
lwgeom_unaryunion
LWGEOM * lwgeom_unaryunion(const LWGEOM *geom1)
Definition: liblwgeom/lwgeom_geos.c:725
POLYGONTYPE
#define POLYGONTYPE
Definition: liblwgeom.h:87
lwgeom_buildarea
LWGEOM * lwgeom_buildarea(const LWGEOM *geom)
Take a geometry and return an areal geometry (Polygon or MultiPolygon).
Definition: liblwgeom/lwgeom_geos.c:1203
WKB_EXTENDED
#define WKB_EXTENDED
Definition: liblwgeom.h:2068
lwgeom_to_hexwkb
char * lwgeom_to_hexwkb(const LWGEOM *geom, uint8_t variant, size_t *size_out)
Definition: lwout_wkb.c:848
lwgeom_is_empty
int lwgeom_is_empty(const LWGEOM *geom)
Return true or false depending on whether a geometry is an "empty" geometry (no vertices members)
Definition: lwgeom.c:1393
lwgeom_to_wkt
char * lwgeom_to_wkt(const LWGEOM *geom, uint8_t variant, int precision, size_t *size_out)
WKT emitter function.
Definition: lwout_wkt.c:676
lwgeom_get_bbox
const GBOX * lwgeom_get_bbox(const LWGEOM *lwgeom)
Get a non-empty geometry bounding box, computing and caching it if not already there.
Definition: lwgeom.c:734
TRIANGLETYPE
#define TRIANGLETYPE
Definition: liblwgeom.h:98
lwgeom_linemerge
LWGEOM * lwgeom_linemerge(const LWGEOM *geom1)
Definition: liblwgeom/lwgeom_geos.c:694
lwgeom_clip_by_rect
LWGEOM * lwgeom_clip_by_rect(const LWGEOM *geom1, double x0, double y0, double x1, double y1)
Definition: liblwgeom/lwgeom_geos.c:921
error_if_srid_mismatch
void error_if_srid_mismatch(int srid1, int srid2)
Definition: lwutil.c:338
LW_TRUE
#define LW_TRUE
Return types for functions with status returns.
Definition: liblwgeom.h:76
lwgeom_as_lwpoly
LWPOLY * lwgeom_as_lwpoly(const LWGEOM *lwgeom)
Definition: lwgeom.c:206
SRID_UNKNOWN
#define SRID_UNKNOWN
Unknown SRID value.
Definition: liblwgeom.h:188
lwgeom_construct_empty
LWGEOM * lwgeom_construct_empty(uint8_t type, int srid, char hasz, char hasm)
Definition: lwgeom.c:2092
lwgeom_add_bbox
void lwgeom_add_bbox(LWGEOM *lwgeom)
Compute a bbox if not already computed.
Definition: lwgeom.c:686
lwpoly_construct_empty
LWPOLY * lwpoly_construct_empty(int srid, char hasz, char hasm)
Definition: lwpoly.c:161
lwgeom_union
LWGEOM * lwgeom_union(const LWGEOM *geom1, const LWGEOM *geom2)
Definition: liblwgeom/lwgeom_geos.c:888
liblwgeom.h
This library is the generic geometry handling section of PostGIS.
point_in_multipolygon
int point_in_multipolygon(LWMPOLY *mpolygon, LWPOINT *point)
Definition: lwgeom_functions_analytic.c:1000
point_in_multipolygon_rtree
int point_in_multipolygon_rtree(RTREE_NODE **root, int polyCount, int *ringCounts, LWPOINT *point)
Definition: lwgeom_functions_analytic.c:891
point_in_polygon
int point_in_polygon(LWPOLY *polygon, LWPOINT *point)
Definition: lwgeom_functions_analytic.c:956
lwgeom_functions_analytic.h
GetPrepGeomCache
PrepGeomCache * GetPrepGeomCache(FunctionCallInfo fcinfo, GSERIALIZED *g1, GSERIALIZED *g2)
Given a couple potential geometries and a function call context, return a prepared structure for one ...
Definition: lwgeom_geos_prepared.c:462
lwgeom_geos_prepared.h
GetRtreeCache
RTREE_POLY_CACHE * GetRtreeCache(FunctionCallInfo fcinfo, GSERIALIZED *g1)
Checks for a cache hit against the provided geometry and returns a pre-built index structure (RTREE_P...
Definition: lwgeom_rtree.c:432
lwgeom_rtree.h
genraster.value
int value
Definition: genraster.py:61
genraster.count
int count
Definition: genraster.py:56
ovdump.type
type
Definition: ovdump.py:41
array_nelems_not_null
uint32_t array_nelems_not_null(ArrayType *array)
Definition: postgis/lwgeom_geos.c:2670
polygonize_garray
Datum polygonize_garray(PG_FUNCTION_ARGS)
Definition: postgis/lwgeom_geos.c:2829
pip_short_circuit
static int pip_short_circuit(RTREE_POLY_CACHE *poly_cache, LWPOINT *point, GSERIALIZED *gpoly)
Definition: postgis/lwgeom_geos.c:142
contains
Datum contains(PG_FUNCTION_ARGS)
Definition: postgis/lwgeom_geos.c:1624
coveredby
Datum coveredby(PG_FUNCTION_ARGS)
Definition: postgis/lwgeom_geos.c:1965
isvaliddetail
Datum isvaliddetail(PG_FUNCTION_ARGS)
Definition: postgis/lwgeom_geos.c:1472
geos_geomunion
Datum geos_geomunion(PG_FUNCTION_ARGS)
Definition: postgis/lwgeom_geos.c:554
HANDLE_GEOS_ERROR
#define HANDLE_GEOS_ERROR(label)
Definition: postgis/lwgeom_geos.c:56
isvalid
Datum isvalid(PG_FUNCTION_ARGS)
Definition: postgis/lwgeom_geos.c:1396
is_point
static char is_point(const GSERIALIZED *g)
Definition: postgis/lwgeom_geos.c:132
ST_Snap
Datum ST_Snap(PG_FUNCTION_ARGS)
Definition: postgis/lwgeom_geos.c:3125
ST_BuildArea
Datum ST_BuildArea(PG_FUNCTION_ARGS)
Definition: postgis/lwgeom_geos.c:3059
errorIfGeometryCollection
void errorIfGeometryCollection(GSERIALIZED *g1, GSERIALIZED *g2)
Throws an ereport ERROR if either geometry is a COLLECTIONTYPE.
Definition: postgis/lwgeom_geos.c:1357
ST_Split
Datum ST_Split(PG_FUNCTION_ARGS)
Definition: postgis/lwgeom_geos.c:3176
covers
Datum covers(PG_FUNCTION_ARGS)
Definition: postgis/lwgeom_geos.c:1830
ST_Node
Datum ST_Node(PG_FUNCTION_ARGS)
Definition: postgis/lwgeom_geos.c:3266
ST_GeneratePoints
Datum ST_GeneratePoints(PG_FUNCTION_ARGS)
Definition: postgis/lwgeom_geos.c:1028
hausdorffdistancedensify
Datum hausdorffdistancedensify(PG_FUNCTION_ARGS)
Definition: postgis/lwgeom_geos.c:224
ST_OffsetCurve
Datum ST_OffsetCurve(PG_FUNCTION_ARGS)
Definition: postgis/lwgeom_geos.c:1064
cluster_within_distance_garray
Datum cluster_within_distance_garray(PG_FUNCTION_ARGS)
Definition: postgis/lwgeom_geos.c:2953
centroid
Datum centroid(PG_FUNCTION_ARGS)
Definition: postgis/lwgeom_geos.c:1273
touches
Datum touches(PG_FUNCTION_ARGS)
Definition: postgis/lwgeom_geos.c:2267
disjoint
Datum disjoint(PG_FUNCTION_ARGS)
Definition: postgis/lwgeom_geos.c:2326
geos_intersects
Datum geos_intersects(PG_FUNCTION_ARGS)
Definition: postgis/lwgeom_geos.c:2139
ST_MinimumClearance
Datum ST_MinimumClearance(PG_FUNCTION_ARGS)
Definition: postgis/lwgeom_geos.c:3391
GEOSnoop
Datum GEOSnoop(PG_FUNCTION_ARGS)
Definition: postgis/lwgeom_geos.c:2808
ST_FrechetDistance
Datum ST_FrechetDistance(PG_FUNCTION_ARGS)
Definition: postgis/lwgeom_geos.c:274
buffer
Datum buffer(PG_FUNCTION_ARGS)
Definition: postgis/lwgeom_geos.c:796
convexhull
Datum convexhull(PG_FUNCTION_ARGS)
Definition: postgis/lwgeom_geos.c:686
POSTGIS2GEOS
GEOSGeometry * POSTGIS2GEOS(GSERIALIZED *pglwgeom)
Definition: postgis/lwgeom_geos.c:2655
ST_ClipByBox2d
Datum ST_ClipByBox2d(PG_FUNCTION_ARGS)
Definition: postgis/lwgeom_geos.c:1294
crosses
Datum crosses(PG_FUNCTION_ARGS)
Definition: postgis/lwgeom_geos.c:2081
relate_full
Datum relate_full(PG_FUNCTION_ARGS)
Definition: postgis/lwgeom_geos.c:2442
relate_pattern
Datum relate_pattern(PG_FUNCTION_ARGS)
Definition: postgis/lwgeom_geos.c:2385
issimple
Datum issimple(PG_FUNCTION_ARGS)
Definition: postgis/lwgeom_geos.c:2569
symdifference
Datum symdifference(PG_FUNCTION_ARGS)
Definition: postgis/lwgeom_geos.c:586
geos_difference
Datum geos_difference(PG_FUNCTION_ARGS)
Definition: postgis/lwgeom_geos.c:1222
ARRAY2LWGEOM
LWGEOM ** ARRAY2LWGEOM(ArrayType *array, uint32_t nelems, int *is3d, int *srid)
Definition: postgis/lwgeom_geos.c:2694
topologypreservesimplify
Datum topologypreservesimplify(PG_FUNCTION_ARGS)
Definition: postgis/lwgeom_geos.c:751
geos_intersection
Datum geos_intersection(PG_FUNCTION_ARGS)
Definition: postgis/lwgeom_geos.c:1190
PG_FUNCTION_INFO_V1
PG_FUNCTION_INFO_V1(postgis_geos_version)
GEOS2POSTGIS
GSERIALIZED * GEOS2POSTGIS(GEOSGeom geom, char want3d)
Definition: postgis/lwgeom_geos.c:2629
ST_DelaunayTriangles
Datum ST_DelaunayTriangles(PG_FUNCTION_ARGS)
Definition: postgis/lwgeom_geos.c:3089
is_poly
static char is_poly(const GSERIALIZED *g)
Definition: postgis/lwgeom_geos.c:125
boundary
Datum boundary(PG_FUNCTION_ARGS)
Definition: postgis/lwgeom_geos.c:614
postgis_geos_version
Datum postgis_geos_version(PG_FUNCTION_ARGS)
Definition: postgis/lwgeom_geos.c:117
ST_UnaryUnion
Datum ST_UnaryUnion(PG_FUNCTION_ARGS)
Definition: postgis/lwgeom_geos.c:525
overlaps
Datum overlaps(PG_FUNCTION_ARGS)
Definition: postgis/lwgeom_geos.c:1565
ST_Equals
Datum ST_Equals(PG_FUNCTION_ARGS)
Definition: postgis/lwgeom_geos.c:2501
ST_OrientedEnvelope
Datum ST_OrientedEnvelope(PG_FUNCTION_ARGS)
Definition: postgis/lwgeom_geos.c:3474
ST_SharedPaths
Datum ST_SharedPaths(PG_FUNCTION_ARGS)
Definition: postgis/lwgeom_geos.c:3225
ST_MinimumClearanceLine
Datum ST_MinimumClearanceLine(PG_FUNCTION_ARGS)
Definition: postgis/lwgeom_geos.c:3430
ST_Voronoi
Datum ST_Voronoi(PG_FUNCTION_ARGS)
Definition: postgis/lwgeom_geos.c:3301
clusterintersecting_garray
Datum clusterintersecting_garray(PG_FUNCTION_ARGS)
Definition: postgis/lwgeom_geos.c:2887
ARRAY2GEOS
GEOSGeometry ** ARRAY2GEOS(ArrayType *array, uint32_t nelems, int *is3d, int *srid)
Definition: postgis/lwgeom_geos.c:2745
pgis_union_geometry_array
Datum pgis_union_geometry_array(PG_FUNCTION_ARGS)
Definition: postgis/lwgeom_geos.c:343
hausdorffdistance
Datum hausdorffdistance(PG_FUNCTION_ARGS)
Definition: postgis/lwgeom_geos.c:175
isvalidreason
Datum isvalidreason(PG_FUNCTION_ARGS)
Definition: postgis/lwgeom_geos.c:1442
isring
Datum isring(PG_FUNCTION_ARGS)
Definition: postgis/lwgeom_geos.c:2595
linemerge
Datum linemerge(PG_FUNCTION_ARGS)
Definition: postgis/lwgeom_geos.c:3029
pointonsurface
Datum pointonsurface(PG_FUNCTION_ARGS)
Definition: postgis/lwgeom_geos.c:1253
containsproperly
Datum containsproperly(PG_FUNCTION_ARGS)
Definition: postgis/lwgeom_geos.c:1759
text_to_cstring
char * text_to_cstring(const text *textptr)
geometry_serialize
GSERIALIZED * geometry_serialize(LWGEOM *lwgeom)
int32
unsigned int int32
Definition: shpopen.c:273
POSTGIS_GEOS_VERSION
#define POSTGIS_GEOS_VERSION
Definition: sqldefines.h:10
GBOX::ymax
double ymax
Definition: liblwgeom.h:298
GBOX::xmax
double xmax
Definition: liblwgeom.h:296
GBOX::ymin
double ymin
Definition: liblwgeom.h:297
GBOX::xmin
double xmin
Definition: liblwgeom.h:295
GBOX::flags
uint8_t flags
Definition: liblwgeom.h:294
LWGEOM::type
uint8_t type
Definition: liblwgeom.h:399
LWGEOM::bbox
GBOX * bbox
Definition: liblwgeom.h:401
LWGEOM::flags
uint8_t flags
Definition: liblwgeom.h:400
LWGEOM::srid
int32_t srid
Definition: liblwgeom.h:402
LWMPOINT::ngeoms
uint32_t ngeoms
Definition: liblwgeom.h:471
LWMPOINT::geoms
LWPOINT ** geoms
Definition: liblwgeom.h:473
PrepGeomCache::prepared_geom
const GEOSPreparedGeometry * prepared_geom
Definition: lwgeom_geos_prepared.h:60
RTREE_POLY_CACHE::polyCount
int polyCount
Definition: lwgeom_rtree.h:62
RTREE_POLY_CACHE::ringCounts
int * ringCounts
Definition: lwgeom_rtree.h:61
RTREE_POLY_CACHE::ringIndices
RTREE_NODE ** ringIndices
Definition: lwgeom_rtree.h:60
RTREE_POLY_CACHE
The tree structure used for fast P-i-P tests by point_in_multipolygon_rtree()
Definition: lwgeom_rtree.h:59
uint32_t
unsigned int uint32_t
Definition: uthash.h:78