PostGIS 3.7.0dev-r@@SVN_REVISION@@
Loading...
Searching...
No Matches
lwgeom_geos_predicates.c
Go to the documentation of this file.
1/**********************************************************************
2 *
3 * PostGIS - Spatial Types for PostgreSQL
4 * http://postgis.net
5 *
6 * PostGIS is free software: you can redistribute it and/or modify
7 * it under the terms of the GNU General Public License as published by
8 * the Free Software Foundation, either version 2 of the License, or
9 * (at your option) any later version.
10 *
11 * PostGIS is distributed in the hope that it will be useful,
12 * but WITHOUT ANY WARRANTY; without even the implied warranty of
13 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14 * GNU General Public License for more details.
15 *
16 * You should have received a copy of the GNU General Public License
17 * along with PostGIS. If not, see <http://www.gnu.org/licenses/>.
18 *
19 **********************************************************************
20 *
21 * Copyright 2009-2014 Sandro Santilli <strk@kbt.io>
22 * Copyright 2008 Paul Ramsey <pramsey@cleverelephant.ca>
23 * Copyright 2001-2003 Refractions Research Inc.
24 *
25 **********************************************************************/
26
27
28#include "../postgis_config.h"
29
30/* PostgreSQL */
31#include "postgres.h"
32#include "funcapi.h"
33#include "utils/array.h"
34#include "utils/builtins.h"
35#include "utils/lsyscache.h"
36#include "utils/numeric.h"
37#include "access/htup_details.h"
38
39/* PostGIS */
40#include "lwgeom_geos.h"
41#include "liblwgeom.h"
42#include "liblwgeom_internal.h"
43#include "lwgeom_itree.h"
44#include "lwgeom_geos_prepared.h"
45#include "lwgeom_accum.h"
46
47
48/* Prototypes for SQL-bound functions */
49Datum relate_full(PG_FUNCTION_ARGS);
50Datum relate_pattern(PG_FUNCTION_ARGS);
51Datum disjoint(PG_FUNCTION_ARGS);
52Datum touches(PG_FUNCTION_ARGS);
53Datum ST_Intersects(PG_FUNCTION_ARGS);
54Datum crosses(PG_FUNCTION_ARGS);
55Datum contains(PG_FUNCTION_ARGS);
56Datum within(PG_FUNCTION_ARGS);
57Datum containsproperly(PG_FUNCTION_ARGS);
58Datum covers(PG_FUNCTION_ARGS);
59Datum overlaps(PG_FUNCTION_ARGS);
60Datum coveredby(PG_FUNCTION_ARGS);
61Datum ST_Equals(PG_FUNCTION_ARGS);
62
63
64/*
65 * Utility to quickly check for polygonal geometries
66 */
67static inline char
68is_poly(const GSERIALIZED* g)
69{
70 int type = gserialized_get_type(g);
71 return type == POLYGONTYPE || type == MULTIPOLYGONTYPE;
72}
73
74/*
75 * Utility to quickly check for point geometries
76 */
77static inline char
78is_point(const GSERIALIZED* g)
79{
80 int type = gserialized_get_type(g);
81 return type == POINTTYPE || type == MULTIPOINTTYPE;
82}
83
84
85
86
87
88
89
90PG_FUNCTION_INFO_V1(ST_Intersects);
91Datum ST_Intersects(PG_FUNCTION_ARGS)
92{
93 SHARED_GSERIALIZED *shared_geom1 = ToastCacheGetGeometry(fcinfo, 0);
94 SHARED_GSERIALIZED *shared_geom2 = ToastCacheGetGeometry(fcinfo, 1);
95 const GSERIALIZED *geom1 = shared_gserialized_get(shared_geom1);
96 const GSERIALIZED *geom2 = shared_gserialized_get(shared_geom2);
97 int result;
98 GBOX box1, box2;
99 PrepGeomCache *prep_cache;
100
101 gserialized_error_if_srid_mismatch(geom1, geom2, __func__);
102
103 /* A.Intersects(Empty) == FALSE */
104 if ( gserialized_is_empty(geom1) || gserialized_is_empty(geom2) )
105 PG_RETURN_BOOL(false);
106
107 /*
108 * Short-circuit 1: if geom2 bounding box does not overlap
109 * geom1 bounding box we can return FALSE.
110 */
111 if ( gserialized_get_gbox_p(geom1, &box1) &&
112 gserialized_get_gbox_p(geom2, &box2) )
113 {
114 if ( gbox_overlaps_2d(&box1, &box2) == LW_FALSE )
115 PG_RETURN_BOOL(false);
116 }
117
118 /*
119 * Short-circuit 2: if the geoms are a point and a polygon,
120 * call the itree_pip_intersects function.
121 */
122 if ((is_point(geom1) && is_poly(geom2)) ||
123 (is_point(geom2) && is_poly(geom1)))
124 {
125 SHARED_GSERIALIZED *shared_gpoly = is_poly(geom1) ? shared_geom1 : shared_geom2;
126 SHARED_GSERIALIZED *shared_gpoint = is_point(geom1) ? shared_geom1 : shared_geom2;
127 const GSERIALIZED *gpoint = shared_gserialized_get(shared_gpoint);
128 LWGEOM *lwpt = lwgeom_from_gserialized(gpoint);
129 IntervalTree *itree = GetIntervalTree(fcinfo, shared_gpoly);
130 bool result = itree_pip_intersects(itree, lwpt);
131 lwgeom_free(lwpt);
132 PG_RETURN_BOOL(result);
133 }
134
135 initGEOS(lwpgnotice, lwgeom_geos_error);
136
137 prep_cache = GetPrepGeomCache(fcinfo, shared_geom1, shared_geom2);
138 if ( prep_cache && prep_cache->prepared_geom )
139 {
140 GEOSGeometry *g = prep_cache->gcache.argnum == 1
141 ? POSTGIS2GEOS(geom2)
142 : POSTGIS2GEOS(geom1);
143 if (!g) HANDLE_GEOS_ERROR("Geometry could not be converted to GEOS");
144 result = GEOSPreparedIntersects(prep_cache->prepared_geom, g);
145 GEOSGeom_destroy(g);
146 }
147 else
148 {
149 GEOSGeometry *g1, *g2;
150 g1 = POSTGIS2GEOS(geom1);
151 if (!g1) HANDLE_GEOS_ERROR("First argument geometry could not be converted to GEOS");
152 g2 = POSTGIS2GEOS(geom2);
153 if (!g2)
154 {
155 GEOSGeom_destroy(g1);
156 HANDLE_GEOS_ERROR("Second argument geometry could not be converted to GEOS");
157 }
158 result = GEOSIntersects(g1, g2);
159 GEOSGeom_destroy(g1);
160 GEOSGeom_destroy(g2);
161 }
162
163 if (result == 2) HANDLE_GEOS_ERROR("GEOSIntersects");
164
165 PG_RETURN_BOOL(result);
166}
167
168
169PG_FUNCTION_INFO_V1(ST_Equals);
170Datum ST_Equals(PG_FUNCTION_ARGS)
171{
172 SHARED_GSERIALIZED *shared_geom1 = ToastCacheGetGeometry(fcinfo, 0);
173 SHARED_GSERIALIZED *shared_geom2 = ToastCacheGetGeometry(fcinfo, 1);
174 const GSERIALIZED *geom1 = shared_gserialized_get(shared_geom1);
175 const GSERIALIZED *geom2 = shared_gserialized_get(shared_geom2);
176 char result;
177 GBOX box1, box2;
178 GEOSGeometry *g1, *g2;
179
180 gserialized_error_if_srid_mismatch(geom1, geom2, __func__);
181
182 /* Empty == Empty */
183 if ( gserialized_is_empty(geom1) && gserialized_is_empty(geom2) )
184 PG_RETURN_BOOL(true);
185
186 /*
187 * Short-circuit: If geom1 and geom2 do not have the same bounding box
188 * we can return FALSE.
189 */
190 if ( gserialized_get_gbox_p(geom1, &box1) &&
191 gserialized_get_gbox_p(geom2, &box2) )
192 {
193 if ( gbox_same_2d_float(&box1, &box2) == LW_FALSE )
194 {
195 PG_RETURN_BOOL(false);
196 }
197 }
198
199 /*
200 * Short-circuit: if geom1 and geom2 are binary-equivalent, we can return
201 * TRUE. This is much faster than doing the comparison using GEOS.
202 */
203 if (VARSIZE(geom1) == VARSIZE(geom2) && !memcmp(geom1, geom2, VARSIZE(geom1))) {
204 PG_RETURN_BOOL(true);
205 }
206
207 initGEOS(lwpgnotice, lwgeom_geos_error);
208
209 g1 = POSTGIS2GEOS(geom1);
210 if (!g1) HANDLE_GEOS_ERROR("First argument geometry could not be converted to GEOS");
211 g2 = POSTGIS2GEOS(geom2);
212 if (!g2)
213 {
214 GEOSGeom_destroy(g1);
215 HANDLE_GEOS_ERROR("Second argument geometry could not be converted to GEOS");
216 }
217 result = GEOSEquals(g1, g2);
218 GEOSGeom_destroy(g1);
219 GEOSGeom_destroy(g2);
220
221 if (result == 2) HANDLE_GEOS_ERROR("GEOSEquals");
222
223 PG_RETURN_BOOL(result);
224}
225
226
227
228PG_FUNCTION_INFO_V1(touches);
229Datum touches(PG_FUNCTION_ARGS)
230{
231 SHARED_GSERIALIZED *shared_geom1 = ToastCacheGetGeometry(fcinfo, 0);
232 SHARED_GSERIALIZED *shared_geom2 = ToastCacheGetGeometry(fcinfo, 1);
233 const GSERIALIZED *geom1 = shared_gserialized_get(shared_geom1);
234 const GSERIALIZED *geom2 = shared_gserialized_get(shared_geom2);
235 char result;
236 GBOX box1, box2;
237 PrepGeomCache *prep_cache;
238
239 gserialized_error_if_srid_mismatch(geom1, geom2, __func__);
240
241 /* A.Touches(Empty) == FALSE */
242 if ( gserialized_is_empty(geom1) || gserialized_is_empty(geom2) )
243 PG_RETURN_BOOL(false);
244
245 /*
246 * Short-circuit 1: if geom2 bounding box does not overlap
247 * geom1 bounding box we can return FALSE.
248 */
249 if ( gserialized_get_gbox_p(geom1, &box1) &&
250 gserialized_get_gbox_p(geom2, &box2) )
251 {
252 if ( gbox_overlaps_2d(&box1, &box2) == LW_FALSE )
253 {
254 PG_RETURN_BOOL(false);
255 }
256 }
257
258 initGEOS(lwpgnotice, lwgeom_geos_error);
259
260 prep_cache = GetPrepGeomCache(fcinfo, shared_geom1, shared_geom2);
261 if ( prep_cache && prep_cache->prepared_geom )
262 {
263 GEOSGeometry *g = prep_cache->gcache.argnum == 1
264 ? POSTGIS2GEOS(geom2)
265 : POSTGIS2GEOS(geom1);
266 if (!g) HANDLE_GEOS_ERROR("Geometry could not be converted to GEOS");
267 result = GEOSPreparedTouches(prep_cache->prepared_geom, g);
268 GEOSGeom_destroy(g);
269 }
270 else
271 {
272 GEOSGeometry *g1, *g2;
273 g1 = POSTGIS2GEOS(geom1);
274 if (!g1) HANDLE_GEOS_ERROR("First argument geometry could not be converted to GEOS");
275 g2 = POSTGIS2GEOS(geom2);
276 if (!g2)
277 {
278 GEOSGeom_destroy(g1);
279 HANDLE_GEOS_ERROR("Second argument geometry could not be converted to GEOS");
280 }
281 result = GEOSTouches(g1, g2);
282 GEOSGeom_destroy(g1);
283 GEOSGeom_destroy(g2);
284 }
285
286 if (result == 2) HANDLE_GEOS_ERROR("GEOSTouches");
287
288 PG_RETURN_BOOL(result);
289}
290
291
292PG_FUNCTION_INFO_V1(disjoint);
293Datum disjoint(PG_FUNCTION_ARGS)
294{
295 SHARED_GSERIALIZED *shared_geom1 = ToastCacheGetGeometry(fcinfo, 0);
296 SHARED_GSERIALIZED *shared_geom2 = ToastCacheGetGeometry(fcinfo, 1);
297 const GSERIALIZED *geom1 = shared_gserialized_get(shared_geom1);
298 const GSERIALIZED *geom2 = shared_gserialized_get(shared_geom2);
299 char result;
300 GBOX box1, box2;
301 PrepGeomCache *prep_cache;
302
303 gserialized_error_if_srid_mismatch(geom1, geom2, __func__);
304
305 /* A.Disjoint(Empty) == TRUE */
306 if ( gserialized_is_empty(geom1) || gserialized_is_empty(geom2) )
307 PG_RETURN_BOOL(true);
308
309 /*
310 * Short-circuit 1: if geom2 bounding box does not overlap
311 * geom1 bounding box we can return TRUE.
312 */
313 if ( gserialized_get_gbox_p(geom1, &box1) &&
314 gserialized_get_gbox_p(geom2, &box2) )
315 {
316 if ( gbox_overlaps_2d(&box1, &box2) == LW_FALSE )
317 {
318 PG_RETURN_BOOL(true);
319 }
320 }
321
322 initGEOS(lwpgnotice, lwgeom_geos_error);
323
324 prep_cache = GetPrepGeomCache(fcinfo, shared_geom1, shared_geom2);
325 if ( prep_cache && prep_cache->prepared_geom )
326 {
327 GEOSGeometry *g = prep_cache->gcache.argnum == 1
328 ? POSTGIS2GEOS(geom2)
329 : POSTGIS2GEOS(geom1);
330 if (!g) HANDLE_GEOS_ERROR("Geometry could not be converted to GEOS");
331 result = GEOSPreparedDisjoint(prep_cache->prepared_geom, g);
332 GEOSGeom_destroy(g);
333 }
334 else
335 {
336 GEOSGeometry *g1, *g2;
337 g1 = POSTGIS2GEOS(geom1);
338 if (!g1) HANDLE_GEOS_ERROR("First argument geometry could not be converted to GEOS");
339 g2 = POSTGIS2GEOS(geom2);
340 if (!g2)
341 {
342 GEOSGeom_destroy(g1);
343 HANDLE_GEOS_ERROR("Second argument geometry could not be converted to GEOS");
344 }
345 result = GEOSDisjoint(g1, g2);
346 GEOSGeom_destroy(g1);
347 GEOSGeom_destroy(g2);
348 }
349
350 if (result == 2) HANDLE_GEOS_ERROR("GEOSDisjoint");
351
352 PG_RETURN_BOOL(result);
353}
354
355
356
360PG_FUNCTION_INFO_V1(overlaps);
361Datum overlaps(PG_FUNCTION_ARGS)
362{
363 SHARED_GSERIALIZED *shared_geom1 = ToastCacheGetGeometry(fcinfo, 0);
364 SHARED_GSERIALIZED *shared_geom2 = ToastCacheGetGeometry(fcinfo, 1);
365 const GSERIALIZED *geom1 = shared_gserialized_get(shared_geom1);
366 const GSERIALIZED *geom2 = shared_gserialized_get(shared_geom2);
367 char result;
368 GBOX box1, box2;
369 PrepGeomCache *prep_cache;
370
371 gserialized_error_if_srid_mismatch(geom1, geom2, __func__);
372
373 /* A.Overlaps(Empty) == FALSE */
374 if ( gserialized_is_empty(geom1) || gserialized_is_empty(geom2) )
375 PG_RETURN_BOOL(false);
376
377 /*
378 * Short-circuit 1: if geom2 bounding box does not overlap
379 * geom1 bounding box we can return FALSE.
380 */
381 if ( gserialized_get_gbox_p(geom1, &box1) &&
382 gserialized_get_gbox_p(geom2, &box2) )
383 {
384 if ( ! gbox_overlaps_2d(&box1, &box2) )
385 {
386 PG_RETURN_BOOL(false);
387 }
388 }
389
390 initGEOS(lwpgnotice, lwgeom_geos_error);
391
392 prep_cache = GetPrepGeomCache(fcinfo, shared_geom1, shared_geom2);
393 if ( prep_cache && prep_cache->prepared_geom )
394 {
395 GEOSGeometry *g = prep_cache->gcache.argnum == 1
396 ? POSTGIS2GEOS(geom2)
397 : POSTGIS2GEOS(geom1);
398 if (!g) HANDLE_GEOS_ERROR("Geometry could not be converted to GEOS");
399 result = GEOSPreparedOverlaps(prep_cache->prepared_geom, g);
400 GEOSGeom_destroy(g);
401 }
402 else
403 {
404 GEOSGeometry *g1, *g2;
405 g1 = POSTGIS2GEOS(geom1);
406 if (!g1) HANDLE_GEOS_ERROR("First argument geometry could not be converted to GEOS");
407 g2 = POSTGIS2GEOS(geom2);
408 if (!g2)
409 {
410 GEOSGeom_destroy(g1);
411 HANDLE_GEOS_ERROR("Second argument geometry could not be converted to GEOS");
412 }
413 result = GEOSOverlaps(g1, g2);
414 GEOSGeom_destroy(g1);
415 GEOSGeom_destroy(g2);
416 }
417
418 if (result == 2) HANDLE_GEOS_ERROR("GEOSOverlaps");
419
420 PG_RETURN_BOOL(result);
421}
422
423
427PG_FUNCTION_INFO_V1(crosses);
428Datum crosses(PG_FUNCTION_ARGS)
429{
430 SHARED_GSERIALIZED *shared_geom1 = ToastCacheGetGeometry(fcinfo, 0);
431 SHARED_GSERIALIZED *shared_geom2 = ToastCacheGetGeometry(fcinfo, 1);
432 const GSERIALIZED *geom1 = shared_gserialized_get(shared_geom1);
433 const GSERIALIZED *geom2 = shared_gserialized_get(shared_geom2);
434 int result;
435 GBOX box1, box2;
436 PrepGeomCache *prep_cache;
437
438 gserialized_error_if_srid_mismatch(geom1, geom2, __func__);
439
440 /* A.Crosses(Empty) == FALSE */
441 if ( gserialized_is_empty(geom1) || gserialized_is_empty(geom2) )
442 PG_RETURN_BOOL(false);
443
444 /*
445 * Short-circuit 1: if geom2 bounding box does not overlap
446 * geom1 bounding box we can return FALSE.
447 */
448 if ( gserialized_get_gbox_p(geom1, &box1) &&
449 gserialized_get_gbox_p(geom2, &box2) )
450 {
451 if ( gbox_overlaps_2d(&box1, &box2) == LW_FALSE )
452 {
453 PG_RETURN_BOOL(false);
454 }
455 }
456
457 initGEOS(lwpgnotice, lwgeom_geos_error);
458
459 prep_cache = GetPrepGeomCache(fcinfo, shared_geom1, shared_geom2);
460 if ( prep_cache && prep_cache->prepared_geom )
461 {
462 GEOSGeometry *g = prep_cache->gcache.argnum == 1
463 ? POSTGIS2GEOS(geom2)
464 : POSTGIS2GEOS(geom1);
465 if (!g) HANDLE_GEOS_ERROR("Geometry could not be converted to GEOS");
466 result = GEOSPreparedCrosses(prep_cache->prepared_geom, g);
467 GEOSGeom_destroy(g);
468 }
469 else
470 {
471 GEOSGeometry *g1, *g2;
472 g1 = POSTGIS2GEOS(geom1);
473 if (!g1) HANDLE_GEOS_ERROR("First argument geometry could not be converted to GEOS");
474 g2 = POSTGIS2GEOS(geom2);
475 if (!g2)
476 {
477 GEOSGeom_destroy(g1);
478 HANDLE_GEOS_ERROR("Second argument geometry could not be converted to GEOS");
479 }
480 result = GEOSCrosses(g1, g2);
481 GEOSGeom_destroy(g1);
482 GEOSGeom_destroy(g2);
483 }
484
485 if (result == 2) HANDLE_GEOS_ERROR("GEOSCrosses");
486
487 PG_RETURN_BOOL(result);
488}
489
490
491
495PG_FUNCTION_INFO_V1(contains);
496Datum contains(PG_FUNCTION_ARGS)
497{
498 SHARED_GSERIALIZED *shared_geom1 = ToastCacheGetGeometry(fcinfo, 0);
499 SHARED_GSERIALIZED *shared_geom2 = ToastCacheGetGeometry(fcinfo, 1);
500 const GSERIALIZED *geom1 = shared_gserialized_get(shared_geom1);
501 const GSERIALIZED *geom2 = shared_gserialized_get(shared_geom2);
502 int result;
503 GEOSGeometry *g1, *g2;
504 GBOX box1, box2;
505 PrepGeomCache *prep_cache;
506
507 gserialized_error_if_srid_mismatch(geom1, geom2, __func__);
508
509 /* A.Contains(Empty) == FALSE */
510 if (gserialized_is_empty(geom1) || gserialized_is_empty(geom2))
511 PG_RETURN_BOOL(false);
512
513 POSTGIS_DEBUGF(3, "Contains: type1: %d, type2: %d", gserialized_get_type(geom1), gserialized_get_type(geom2));
514
515 /*
516 ** Short-circuit 1: if geom2 bounding box is not completely inside
517 ** geom1 bounding box we can return FALSE.
518 */
519 if (gserialized_get_gbox_p(geom1, &box1) &&
520 gserialized_get_gbox_p(geom2, &box2))
521 {
522 if (!gbox_contains_2d(&box1, &box2))
523 PG_RETURN_BOOL(false);
524 }
525
526 /*
527 ** Short-circuit 2: if geom2 is a point and geom1 is a polygon
528 ** call the point-in-polygon function.
529 */
530 if (is_poly(geom1) && is_point(geom2))
531 {
532 const GSERIALIZED *gpoint = shared_gserialized_get(shared_geom2);
533 LWGEOM *lwpt = lwgeom_from_gserialized(gpoint);
534 IntervalTree *itree = GetIntervalTree(fcinfo, shared_geom1);
535 bool result = itree_pip_contains(itree, lwpt);
536 lwgeom_free(lwpt);
537 PG_RETURN_BOOL(result);
538 }
539
540 initGEOS(lwpgnotice, lwgeom_geos_error);
541
542 prep_cache = GetPrepGeomCache(fcinfo, shared_geom1, NULL);
543 if ( prep_cache && prep_cache->prepared_geom && prep_cache->gcache.argnum == 1 )
544 {
545 g1 = POSTGIS2GEOS(geom2);
546 if (!g1)
547 HANDLE_GEOS_ERROR("Geometry could not be converted to GEOS");
548
549 result = GEOSPreparedContains( prep_cache->prepared_geom, g1);
550 GEOSGeom_destroy(g1);
551 }
552 else
553 {
554 g1 = POSTGIS2GEOS(geom1);
555 if (!g1) HANDLE_GEOS_ERROR("First argument geometry could not be converted to GEOS");
556 g2 = POSTGIS2GEOS(geom2);
557 if (!g2)
558 {
559 HANDLE_GEOS_ERROR("Second argument geometry could not be converted to GEOS");
560 GEOSGeom_destroy(g1);
561 }
562 result = GEOSContains( g1, g2);
563 GEOSGeom_destroy(g1);
564 GEOSGeom_destroy(g2);
565 }
566
567 if (result == 2) HANDLE_GEOS_ERROR("GEOSContains");
568
569 PG_RETURN_BOOL(result > 0);
570}
571
572
576PG_FUNCTION_INFO_V1(within);
577Datum within(PG_FUNCTION_ARGS)
578{
579 SHARED_GSERIALIZED *shared_geom1 = ToastCacheGetGeometry(fcinfo, 0);
580 SHARED_GSERIALIZED *shared_geom2 = ToastCacheGetGeometry(fcinfo, 1);
581 const GSERIALIZED *geom1 = shared_gserialized_get(shared_geom1);
582 const GSERIALIZED *geom2 = shared_gserialized_get(shared_geom2);
583 int result;
584 GEOSGeometry *g1, *g2;
585 GBOX box1, box2;
586 PrepGeomCache *prep_cache;
587
588 gserialized_error_if_srid_mismatch(geom1, geom2, __func__);
589
590 /* A.Contains(Empty) == FALSE */
591 if (gserialized_is_empty(geom1) || gserialized_is_empty(geom2))
592 PG_RETURN_BOOL(false);
593
594 POSTGIS_DEBUGF(3, "Contains: type1: %d, type2: %d", gserialized_get_type(geom1), gserialized_get_type(geom2));
595
596 /*
597 ** Short-circuit 1: if geom1 bounding box is not completely inside
598 ** geom2 bounding box we can return FALSE.
599 */
600 if (gserialized_get_gbox_p(geom1, &box1) &&
601 gserialized_get_gbox_p(geom2, &box2))
602 {
603 if (!gbox_within_2d(&box1, &box2))
604 PG_RETURN_BOOL(false);
605 }
606
607 /*
608 ** Short-circuit 2: if geom2 is a polygon and geom1 is a point
609 ** call the point-in-polygon function.
610 */
611 if (is_poly(geom2) && is_point(geom1))
612 {
613 const GSERIALIZED *gpoint = shared_gserialized_get(shared_geom1);
614 LWGEOM *lwpt = lwgeom_from_gserialized(gpoint);
615 IntervalTree *itree = GetIntervalTree(fcinfo, shared_geom2);
616 bool result = itree_pip_contains(itree, lwpt);
617 lwgeom_free(lwpt);
618 PG_RETURN_BOOL(result);
619 }
620
621 initGEOS(lwpgnotice, lwgeom_geos_error);
622
623 prep_cache = GetPrepGeomCache(fcinfo, NULL, shared_geom2);
624 if ( prep_cache && prep_cache->prepared_geom && prep_cache->gcache.argnum == 2 )
625 {
626 g1 = POSTGIS2GEOS(geom1);
627 if (!g1)
628 HANDLE_GEOS_ERROR("Geometry could not be converted to GEOS");
629
630 result = GEOSPreparedContains(prep_cache->prepared_geom, g1);
631 GEOSGeom_destroy(g1);
632 }
633 else
634 {
635 g1 = POSTGIS2GEOS(geom1);
636 if (!g1) HANDLE_GEOS_ERROR("First argument geometry could not be converted to GEOS");
637 g2 = POSTGIS2GEOS(geom2);
638 if (!g2)
639 {
640 HANDLE_GEOS_ERROR("Second argument geometry could not be converted to GEOS");
641 GEOSGeom_destroy(g1);
642 }
643 result = GEOSWithin(g1, g2);
644 GEOSGeom_destroy(g1);
645 GEOSGeom_destroy(g2);
646 }
647
648 if ( result == 2 ) HANDLE_GEOS_ERROR("GEOSWithin");
649
650 PG_RETURN_BOOL(result > 0);
651}
652
656PG_FUNCTION_INFO_V1(containsproperly);
657Datum containsproperly(PG_FUNCTION_ARGS)
658{
659 SHARED_GSERIALIZED *shared_geom1 = ToastCacheGetGeometry(fcinfo, 0);
660 SHARED_GSERIALIZED *shared_geom2 = ToastCacheGetGeometry(fcinfo, 1);
661 const GSERIALIZED *geom1 = shared_gserialized_get(shared_geom1);
662 const GSERIALIZED *geom2 = shared_gserialized_get(shared_geom2);
663 char result;
664 GBOX box1, box2;
665 PrepGeomCache *prep_cache;
666
667 gserialized_error_if_srid_mismatch(geom1, geom2, __func__);
668
669 /* A.ContainsProperly(Empty) == FALSE */
670 if ( gserialized_is_empty(geom1) || gserialized_is_empty(geom2) )
671 PG_RETURN_BOOL(false);
672
673 /*
674 * Short-circuit: if geom2 bounding box is not completely inside
675 * geom1 bounding box we can return FALSE.
676 */
677 if ( gserialized_get_gbox_p(geom1, &box1) &&
678 gserialized_get_gbox_p(geom2, &box2) )
679 {
680 if ( ! gbox_contains_2d(&box1, &box2) )
681 PG_RETURN_BOOL(false);
682 }
683
684 initGEOS(lwpgnotice, lwgeom_geos_error);
685
686 prep_cache = GetPrepGeomCache(fcinfo, shared_geom1, 0);
687 if ( prep_cache && prep_cache->prepared_geom && prep_cache->gcache.argnum == 1 )
688 {
689 GEOSGeometry *g = POSTGIS2GEOS(geom2);
690 if (!g) HANDLE_GEOS_ERROR("First argument geometry could not be converted to GEOS");
691 result = GEOSPreparedContainsProperly( prep_cache->prepared_geom, g);
692 GEOSGeom_destroy(g);
693 }
694 else
695 {
696 GEOSGeometry *g1, *g2;
697 g1 = POSTGIS2GEOS(geom1);
698 if (!g1) HANDLE_GEOS_ERROR("First argument geometry could not be converted to GEOS");
699 g2 = POSTGIS2GEOS(geom2);
700 if (!g2)
701 {
702 GEOSGeom_destroy(g1);
703 HANDLE_GEOS_ERROR("Second argument geometry could not be converted to GEOS");
704 }
705 result = GEOSRelatePattern(g1, g2, "T**FF*FF*");
706 GEOSGeom_destroy(g1);
707 GEOSGeom_destroy(g2);
708 }
709
710 if (result == 2) HANDLE_GEOS_ERROR("GEOSContainProperly");
711
712 PG_RETURN_BOOL(result);
713}
714
715
721PG_FUNCTION_INFO_V1(covers);
722Datum covers(PG_FUNCTION_ARGS)
723{
724 SHARED_GSERIALIZED *shared_geom1 = ToastCacheGetGeometry(fcinfo, 0);
725 SHARED_GSERIALIZED *shared_geom2 = ToastCacheGetGeometry(fcinfo, 1);
726 const GSERIALIZED *geom1 = shared_gserialized_get(shared_geom1);
727 const GSERIALIZED *geom2 = shared_gserialized_get(shared_geom2);
728 int result;
729 GBOX box1, box2;
730 PrepGeomCache *prep_cache;
731
732 POSTGIS_DEBUGF(3, "Covers: type1: %d, type2: %d", gserialized_get_type(geom1), gserialized_get_type(geom2));
733
734 /* A.Covers(Empty) == FALSE */
735 if ( gserialized_is_empty(geom1) || gserialized_is_empty(geom2) )
736 PG_RETURN_BOOL(false);
737
738 gserialized_error_if_srid_mismatch(geom1, geom2, __func__);
739
740 /*
741 * Short-circuit 1: if geom2 bounding box is not completely inside
742 * geom1 bounding box we can return FALSE.
743 */
744 if ( gserialized_get_gbox_p(geom1, &box1) &&
745 gserialized_get_gbox_p(geom2, &box2) )
746 {
747 if ( ! gbox_contains_2d(&box1, &box2) )
748 {
749 PG_RETURN_BOOL(false);
750 }
751 }
752 /*
753 * Short-circuit 2: if geom2 is a point and geom1 is a polygon
754 * call the point-in-polygon function.
755 */
756 if (is_poly(geom1) && is_point(geom2))
757 {
758 const GSERIALIZED *gpoint = shared_gserialized_get(shared_geom2);
759 LWGEOM *lwpt = lwgeom_from_gserialized(gpoint);
760 IntervalTree *itree = GetIntervalTree(fcinfo, shared_geom1);
761 bool result = itree_pip_covers(itree, lwpt);
762 lwgeom_free(lwpt);
763 PG_RETURN_BOOL(result);
764 }
765
766 initGEOS(lwpgnotice, lwgeom_geos_error);
767
768 prep_cache = GetPrepGeomCache(fcinfo, shared_geom1, 0);
769 if ( prep_cache && prep_cache->prepared_geom && prep_cache->gcache.argnum == 1 )
770 {
771 GEOSGeometry *g1 = POSTGIS2GEOS(geom2);
772 if (!g1) HANDLE_GEOS_ERROR("First argument geometry could not be converted to GEOS");
773 result = GEOSPreparedCovers( prep_cache->prepared_geom, g1);
774 GEOSGeom_destroy(g1);
775 }
776 else
777 {
778 GEOSGeometry *g1, *g2;
779
780 g1 = POSTGIS2GEOS(geom1);
781 if (!g1) HANDLE_GEOS_ERROR("First argument geometry could not be converted to GEOS");
782 g2 = POSTGIS2GEOS(geom2);
783 if (!g2)
784 {
785 GEOSGeom_destroy(g1);
786 HANDLE_GEOS_ERROR("Second argument geometry could not be converted to GEOS");
787 }
788 result = GEOSRelatePattern( g1, g2, "******FF*" );
789 GEOSGeom_destroy(g1);
790 GEOSGeom_destroy(g2);
791 }
792
793 if (result == 2) HANDLE_GEOS_ERROR("GEOSCovers");
794
795 PG_RETURN_BOOL(result);
796}
797
798
805PG_FUNCTION_INFO_V1(coveredby);
806Datum coveredby(PG_FUNCTION_ARGS)
807{
808 SHARED_GSERIALIZED *shared_geom1 = ToastCacheGetGeometry(fcinfo, 0);
809 SHARED_GSERIALIZED *shared_geom2 = ToastCacheGetGeometry(fcinfo, 1);
810 const GSERIALIZED *geom1 = shared_gserialized_get(shared_geom1);
811 const GSERIALIZED *geom2 = shared_gserialized_get(shared_geom2);
812 int result;
813 GBOX box1, box2;
814 PrepGeomCache *prep_cache;
815
816 gserialized_error_if_srid_mismatch(geom1, geom2, __func__);
817
818 POSTGIS_DEBUGF(3, "CoveredBy: type1: %d, type2: %d", gserialized_get_type(geom1), gserialized_get_type(geom2));
819
820 /* A.CoveredBy(Empty) == FALSE */
821 if ( gserialized_is_empty(geom1) || gserialized_is_empty(geom2) )
822 PG_RETURN_BOOL(false);
823
824 /*
825 * Short-circuit 1: if geom1 bounding box is not completely inside
826 * geom2 bounding box we can return FALSE.
827 */
828 if ( gserialized_get_gbox_p(geom1, &box1) &&
829 gserialized_get_gbox_p(geom2, &box2) )
830 {
831 if ( ! gbox_within_2d(&box1, &box2) )
832 {
833 PG_RETURN_BOOL(false);
834 }
835 }
836 /*
837 * Short-circuit 2: if geom1 is a point and geom2 is a polygon
838 * call the point-in-polygon function.
839 */
840 if (is_point(geom1) && is_poly(geom2))
841 {
842 const GSERIALIZED *gpoint = shared_gserialized_get(shared_geom1);
843 IntervalTree *itree = GetIntervalTree(fcinfo, shared_geom2);
844 LWGEOM *lwpt = lwgeom_from_gserialized(gpoint);
845 bool result = itree_pip_covers(itree, lwpt);
846 lwgeom_free(lwpt);
847 PG_RETURN_BOOL(result);
848 }
849
850 initGEOS(lwpgnotice, lwgeom_geos_error);
851
852 prep_cache = GetPrepGeomCache(fcinfo, shared_geom1, 0);
853 if ( prep_cache && prep_cache->prepared_geom && prep_cache->gcache.argnum == 2)
854 {
855 GEOSGeometry *g1 = POSTGIS2GEOS(geom2);
856 if (!g1) HANDLE_GEOS_ERROR("First argument geometry could not be converted to GEOS");
857 result = GEOSPreparedCovers( prep_cache->prepared_geom, g1);
858 GEOSGeom_destroy(g1);
859 }
860 else
861 {
862 GEOSGeometry *g1, *g2;
863 g1 = POSTGIS2GEOS(geom1);
864 if (!g1) HANDLE_GEOS_ERROR("First argument geometry could not be converted to GEOS");
865
866 g2 = POSTGIS2GEOS(geom2);
867 if (!g2)
868 {
869 GEOSGeom_destroy(g1);
870 HANDLE_GEOS_ERROR("Second argument geometry could not be converted to GEOS");
871 }
872
873 result = GEOSCoveredBy(g1, g2);
874 GEOSGeom_destroy(g1);
875 GEOSGeom_destroy(g2);
876 }
877
878 if (result == 2) HANDLE_GEOS_ERROR("GEOSCoveredBy");
879
880 PG_RETURN_BOOL(result);
881}
882
883#if POSTGIS_GEOS_VERSION >= 31300
884/*
885 * Flip the DE9IM matrix around the diagonal to
886 * account for flipping the order of the operands
887 * to the GEOSPreparedRelatePattern function.
888 */
889static void
890imInvert(char *im)
891{
892 char t;
893 t = im[1]; im[1] = im[3]; im[3] = t;
894 t = im[2]; im[2] = im[6]; im[6] = t;
895 t = im[5]; im[5] = im[7]; im[7] = t;
896}
897#endif
898
899PG_FUNCTION_INFO_V1(relate_pattern);
900Datum relate_pattern(PG_FUNCTION_ARGS)
901{
902 SHARED_GSERIALIZED *shared_geom1 = ToastCacheGetGeometry(fcinfo, 0);
903 SHARED_GSERIALIZED *shared_geom2 = ToastCacheGetGeometry(fcinfo, 1);
904 const GSERIALIZED *geom1 = shared_gserialized_get(shared_geom1);
905 const GSERIALIZED *geom2 = shared_gserialized_get(shared_geom2);
906
907 /* Ensure DE9IM pattern is no more than 9 chars */
908 text *imPtr = DatumGetTextP(DirectFunctionCall2(text_left,
909 PG_GETARG_DATUM(2), Int32GetDatum(9)));
910 char *im = text_to_cstring(imPtr);
911 char result;
912 uint32_t i;
913#if POSTGIS_GEOS_VERSION >= 31300
914 PrepGeomCache *prep_cache;
915#endif
916
917 gserialized_error_if_srid_mismatch(geom1, geom2, __func__);
918
919 /*
920 ** Need to make sure 't' and 'f' are upper-case before handing to GEOS
921 */
922 for (i = 0; i < strlen(im); i++)
923 im[i] = toupper(im[i]);
924
925 initGEOS(lwpgnotice, lwgeom_geos_error);
926
927#if POSTGIS_GEOS_VERSION >= 31300
928 prep_cache = GetPrepGeomCache(fcinfo, shared_geom1, shared_geom2);
929 if ( prep_cache && prep_cache->prepared_geom )
930 {
931 GEOSGeometry *g = NULL;
932 if (prep_cache->gcache.argnum == 1)
933 {
934 g = POSTGIS2GEOS(geom2);
935 }
936 else
937 {
938 g = POSTGIS2GEOS(geom1);
939 imInvert(im);
940 }
941
942 if (!g) HANDLE_GEOS_ERROR("Geometry could not be converted to GEOS");
943 result = GEOSPreparedRelatePattern(prep_cache->prepared_geom, g, im);
944 GEOSGeom_destroy(g);
945 }
946 else
947#endif
948 {
949 GEOSGeometry *g1, *g2;
950 g1 = POSTGIS2GEOS(geom1);
951 if (!g1) HANDLE_GEOS_ERROR("First argument geometry could not be converted to GEOS");
952 g2 = POSTGIS2GEOS(geom2);
953 if (!g2)
954 {
955 HANDLE_GEOS_ERROR("Second argument geometry could not be converted to GEOS");
956 GEOSGeom_destroy(g1);
957 }
958 result = GEOSRelatePattern(g1, g2, im);
959 GEOSGeom_destroy(g1);
960 GEOSGeom_destroy(g2);
961 }
962
963 pfree(im);
964 if (result == 2) HANDLE_GEOS_ERROR("GEOSRelatePattern");
965
966 PG_RETURN_BOOL(result);
967}
968
969
970PG_FUNCTION_INFO_V1(relate_full);
971Datum relate_full(PG_FUNCTION_ARGS)
972{
973 GSERIALIZED *geom1 = PG_GETARG_GSERIALIZED_P(0);
974 GSERIALIZED *geom2 = PG_GETARG_GSERIALIZED_P(1);
975 GEOSGeometry *g1, *g2;
976 char *relate_str;
977 text *result;
978 int bnr = GEOSRELATE_BNR_OGC;
979
980 /* TODO handle empty */
981 gserialized_error_if_srid_mismatch(geom1, geom2, __func__);
982
983 if ( PG_NARGS() > 2 )
984 bnr = PG_GETARG_INT32(2);
985
986 initGEOS(lwpgnotice, lwgeom_geos_error);
987
988 g1 = POSTGIS2GEOS(geom1);
989 if (!g1)
990 HANDLE_GEOS_ERROR("First argument geometry could not be converted to GEOS");
991 g2 = POSTGIS2GEOS(geom2);
992 if (!g2)
993 {
994 GEOSGeom_destroy(g1);
995 HANDLE_GEOS_ERROR("Second argument geometry could not be converted to GEOS");
996 }
997
998 POSTGIS_DEBUG(3, "constructed geometries ");
999
1000 POSTGIS_DEBUGF(3, "%s", GEOSGeomToWKT(g1));
1001 POSTGIS_DEBUGF(3, "%s", GEOSGeomToWKT(g2));
1002
1003 relate_str = GEOSRelateBoundaryNodeRule(g1, g2, bnr);
1004 if (!relate_str) {
1005 GEOSGeom_destroy(g1);
1006 GEOSGeom_destroy(g2);
1007 HANDLE_GEOS_ERROR("GEOSRelate");
1008 }
1009 result = cstring_to_text(relate_str);
1010 GEOSFree(relate_str);
1011 GEOSGeom_destroy(g1);
1012 GEOSGeom_destroy(g2);
1013
1014 PG_FREE_IF_COPY(geom1, 0);
1015 PG_FREE_IF_COPY(geom2, 1);
1016
1017 PG_RETURN_TEXT_P(result);
1018}
1019
1020
result
char result[OUT_DOUBLE_BUFFER_SIZE]
Definition cu_print.c:267
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 gbox.c:323
gbox_within_2d
int gbox_within_2d(const GBOX *g1, const GBOX *g2)
Return LW_TRUE if the first GBOX is within the second on the 2d plane, LW_FALSE otherwise.
Definition gbox.c:351
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 gbox.c:187
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 gbox.c:339
gserialized_error_if_srid_mismatch
void gserialized_error_if_srid_mismatch(const GSERIALIZED *g1, const GSERIALIZED *g2, const char *funcname)
Definition gserialized.c:432
lwgeom_from_gserialized
LWGEOM * lwgeom_from_gserialized(const GSERIALIZED *g)
Allocate a new LWGEOM from a GSERIALIZED.
Definition gserialized.c:268
gserialized_get_gbox_p
int gserialized_get_gbox_p(const GSERIALIZED *g, GBOX *gbox)
Read the box from the GSERIALIZED or calculate it if necessary.
Definition gserialized.c:94
gserialized_is_empty
int gserialized_is_empty(const GSERIALIZED *g)
Check if a GSERIALIZED is empty without deserializing first.
Definition gserialized.c:181
gserialized_get_type
uint32_t gserialized_get_type(const GSERIALIZED *g)
Extract the geometry type from the serialized form (it hides in the anonymous data area,...
Definition gserialized.c:118
lwgeom_geos_error
void lwgeom_geos_error(const char *fmt,...)
Definition liblwgeom/lwgeom_geos.c:57
LW_FALSE
#define LW_FALSE
Definition liblwgeom.h:94
lwgeom_free
void lwgeom_free(LWGEOM *geom)
Definition lwgeom.c:1246
MULTIPOINTTYPE
#define MULTIPOINTTYPE
Definition liblwgeom.h:105
POINTTYPE
#define POINTTYPE
LWTYPE numbers, used internally by PostGIS.
Definition liblwgeom.h:102
MULTIPOLYGONTYPE
#define MULTIPOLYGONTYPE
Definition liblwgeom.h:107
POLYGONTYPE
#define POLYGONTYPE
Definition liblwgeom.h:104
liblwgeom.h
This library is the generic geometry handling section of PostGIS.
liblwgeom_internal.h
lwgeom_accum.h
contains
Datum contains(PG_FUNCTION_ARGS)
Definition lwgeom_geos_predicates.c:496
coveredby
Datum coveredby(PG_FUNCTION_ARGS)
Definition lwgeom_geos_predicates.c:806
is_point
static char is_point(const GSERIALIZED *g)
Definition lwgeom_geos_predicates.c:78
covers
Datum covers(PG_FUNCTION_ARGS)
Definition lwgeom_geos_predicates.c:722
PG_FUNCTION_INFO_V1
PG_FUNCTION_INFO_V1(ST_Intersects)
touches
Datum touches(PG_FUNCTION_ARGS)
Definition lwgeom_geos_predicates.c:229
disjoint
Datum disjoint(PG_FUNCTION_ARGS)
Definition lwgeom_geos_predicates.c:293
ST_Intersects
Datum ST_Intersects(PG_FUNCTION_ARGS)
Definition lwgeom_geos_predicates.c:91
crosses
Datum crosses(PG_FUNCTION_ARGS)
Definition lwgeom_geos_predicates.c:428
relate_full
Datum relate_full(PG_FUNCTION_ARGS)
Definition lwgeom_geos_predicates.c:971
relate_pattern
Datum relate_pattern(PG_FUNCTION_ARGS)
Definition lwgeom_geos_predicates.c:900
is_poly
static char is_poly(const GSERIALIZED *g)
Definition lwgeom_geos_predicates.c:68
overlaps
Datum overlaps(PG_FUNCTION_ARGS)
Definition lwgeom_geos_predicates.c:361
ST_Equals
Datum ST_Equals(PG_FUNCTION_ARGS)
Definition lwgeom_geos_predicates.c:170
within
Datum within(PG_FUNCTION_ARGS)
Definition lwgeom_geos_predicates.c:577
containsproperly
Datum containsproperly(PG_FUNCTION_ARGS)
Definition lwgeom_geos_predicates.c:657
GetPrepGeomCache
PrepGeomCache * GetPrepGeomCache(FunctionCallInfo fcinfo, SHARED_GSERIALIZED *g1, SHARED_GSERIALIZED *g2)
Given a couple potential geometries and a function call context, return a prepared structure for one ...
Definition lwgeom_geos_prepared.c:369
lwgeom_geos_prepared.h
itree_pip_intersects
bool itree_pip_intersects(const IntervalTree *itree, const LWGEOM *lwpoints)
Definition lwgeom_itree.c:217
GetIntervalTree
IntervalTree * GetIntervalTree(FunctionCallInfo fcinfo, SHARED_GSERIALIZED *g1)
Checks for a cache hit against the provided geometry and returns a pre-built index structure (RTREE_P...
Definition lwgeom_itree.c:114
itree_pip_covers
bool itree_pip_covers(const IntervalTree *itree, const LWGEOM *lwpoints)
Definition lwgeom_itree.c:185
itree_pip_contains
bool itree_pip_contains(const IntervalTree *itree, const LWGEOM *lwpoints)
Definition lwgeom_itree.c:140
lwgeom_itree.h
POSTGIS2GEOS
GEOSGeometry * POSTGIS2GEOS(const GSERIALIZED *pglwgeom)
Definition postgis/lwgeom_geos.c:1823
HANDLE_GEOS_ERROR
#define HANDLE_GEOS_ERROR(label)
Definition postgis/lwgeom_geos.h:56
lwgeom_geos.h
PrepGeomCache::prepared_geom
const GEOSPreparedGeometry * prepared_geom
Definition lwgeom_geos_prepared.h:54