PostGIS 3.0.6dev-r@@SVN_REVISION@@
Loading...
Searching...
No Matches
geography_measurement.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 (C) 2009 Paul Ramsey <pramsey@cleverelephant.ca>
22 *
23 **********************************************************************/
24
25
26#include "postgres.h"
27#include "catalog/pg_type.h" /* for CSTRINGOID */
28
29#include "../postgis_config.h"
30
31#include <math.h>
32#include <float.h>
33#include <string.h>
34#include <stdio.h>
35
36#include "liblwgeom.h" /* For standard geometry types. */
37#include "liblwgeom_internal.h" /* For FP comparators. */
38#include "lwgeom_pg.h" /* For debugging macros. */
39#include "geography.h" /* For utility functions. */
40#include "geography_measurement_trees.h" /* For circ_tree caching */
41#include "lwgeom_transform.h" /* For SRID functions */
42
43#ifdef PROJ_GEODESIC
44/* round to 10 nm precision */
45#define INVMINDIST 1.0e8
46#else
47/* round to 100 nm precision */
48#define INVMINDIST 1.0e9
49#endif
50
51Datum geography_distance(PG_FUNCTION_ARGS);
52Datum geography_distance_uncached(PG_FUNCTION_ARGS);
53Datum geography_distance_knn(PG_FUNCTION_ARGS);
54Datum geography_distance_tree(PG_FUNCTION_ARGS);
55Datum geography_dwithin(PG_FUNCTION_ARGS);
56Datum geography_dwithin_uncached(PG_FUNCTION_ARGS);
57Datum geography_area(PG_FUNCTION_ARGS);
58Datum geography_length(PG_FUNCTION_ARGS);
59Datum geography_expand(PG_FUNCTION_ARGS);
60Datum geography_point_outside(PG_FUNCTION_ARGS);
61Datum geography_covers(PG_FUNCTION_ARGS);
62Datum geography_coveredby(PG_FUNCTION_ARGS);
63Datum geography_bestsrid(PG_FUNCTION_ARGS);
64Datum geography_perimeter(PG_FUNCTION_ARGS);
65Datum geography_project(PG_FUNCTION_ARGS);
66Datum geography_azimuth(PG_FUNCTION_ARGS);
67Datum geography_segmentize(PG_FUNCTION_ARGS);
68
69
71Datum geography_distance_knn(PG_FUNCTION_ARGS)
72{
73 LWGEOM *lwgeom1 = NULL;
74 LWGEOM *lwgeom2 = NULL;
75 GSERIALIZED *g1 = NULL;
76 GSERIALIZED *g2 = NULL;
77 double distance;
78 double tolerance = FP_TOLERANCE;
79 bool use_spheroid = false; /* must use sphere, can't get index to harmonize with spheroid */
80 SPHEROID s;
81
82 /* Get our geometry objects loaded into memory. */
83 g1 = PG_GETARG_GSERIALIZED_P(0);
84 g2 = PG_GETARG_GSERIALIZED_P(1);
85
86 gserialized_error_if_srid_mismatch(g1, g2, __func__);
87
88 /* Initialize spheroid */
89 spheroid_init_from_srid(fcinfo, gserialized_get_srid(g1), &s);
90
91 /* Set to sphere if requested */
92 if ( ! use_spheroid )
93 s.a = s.b = s.radius;
94
95 lwgeom1 = lwgeom_from_gserialized(g1);
96 lwgeom2 = lwgeom_from_gserialized(g2);
97
98 /* Return NULL on empty arguments. */
99 if ( lwgeom_is_empty(lwgeom1) || lwgeom_is_empty(lwgeom2) )
100 {
101 PG_FREE_IF_COPY(g1, 0);
102 PG_FREE_IF_COPY(g2, 1);
103 PG_RETURN_NULL();
104 }
105
106 /* Make sure we have boxes attached */
107 lwgeom_add_bbox_deep(lwgeom1, NULL);
108 lwgeom_add_bbox_deep(lwgeom2, NULL);
109
110 distance = lwgeom_distance_spheroid(lwgeom1, lwgeom2, &s, tolerance);
111
112 POSTGIS_DEBUGF(2, "[GIST] '%s' got distance %g", __func__, distance);
113
114 /* Clean up */
115 lwgeom_free(lwgeom1);
116 lwgeom_free(lwgeom2);
117 PG_FREE_IF_COPY(g1, 0);
118 PG_FREE_IF_COPY(g2, 1);
119
120 /* Something went wrong, negative return... should already be eloged, return NULL */
121 if ( distance < 0.0 )
122 {
123 PG_RETURN_NULL();
124 }
125
126 PG_RETURN_FLOAT8(distance);
127}
128
129/*
130** geography_distance_uncached(GSERIALIZED *g1, GSERIALIZED *g2, double tolerance, boolean use_spheroid)
131** returns double distance in meters
132*/
134Datum geography_distance_uncached(PG_FUNCTION_ARGS)
135{
136 LWGEOM *lwgeom1 = NULL;
137 LWGEOM *lwgeom2 = NULL;
138 GSERIALIZED *g1 = NULL;
139 GSERIALIZED *g2 = NULL;
140 double distance;
141 double tolerance = FP_TOLERANCE;
142 bool use_spheroid = true;
143 SPHEROID s;
144
145 /* Get our geometry objects loaded into memory. */
146 g1 = PG_GETARG_GSERIALIZED_P(0);
147 g2 = PG_GETARG_GSERIALIZED_P(1);
148
149 /* Read our tolerance value. */
150 if ( PG_NARGS() > 2 && ! PG_ARGISNULL(2) )
151 tolerance = PG_GETARG_FLOAT8(2);
152
153 /* Read our calculation type. */
154 if ( PG_NARGS() > 3 && ! PG_ARGISNULL(3) )
155 use_spheroid = PG_GETARG_BOOL(3);
156
157 gserialized_error_if_srid_mismatch(g1, g2, __func__);
158
159 /* Initialize spheroid */
160 spheroid_init_from_srid(fcinfo, gserialized_get_srid(g1), &s);
161
162 /* Set to sphere if requested */
163 if ( ! use_spheroid )
164 s.a = s.b = s.radius;
165
166 lwgeom1 = lwgeom_from_gserialized(g1);
167 lwgeom2 = lwgeom_from_gserialized(g2);
168
169 /* Return NULL on empty arguments. */
170 if ( lwgeom_is_empty(lwgeom1) || lwgeom_is_empty(lwgeom2) )
171 {
172 PG_FREE_IF_COPY(g1, 0);
173 PG_FREE_IF_COPY(g2, 1);
174 PG_RETURN_NULL();
175 }
176
177 /* Make sure we have boxes attached */
178 lwgeom_add_bbox_deep(lwgeom1, NULL);
179 lwgeom_add_bbox_deep(lwgeom2, NULL);
180
181 distance = lwgeom_distance_spheroid(lwgeom1, lwgeom2, &s, tolerance);
182
183 POSTGIS_DEBUGF(2, "[GIST] '%s' got distance %g", __func__, distance);
184
185 /* Clean up */
186 lwgeom_free(lwgeom1);
187 lwgeom_free(lwgeom2);
188 PG_FREE_IF_COPY(g1, 0);
189 PG_FREE_IF_COPY(g2, 1);
190
191 /* Something went wrong, negative return... should already be eloged, return NULL */
192 if ( distance < 0.0 )
193 {
194 PG_RETURN_NULL();
195 }
196
197 PG_RETURN_FLOAT8(distance);
198}
199
200
201/*
202** geography_distance(GSERIALIZED *g1, GSERIALIZED *g2, double tolerance, boolean use_spheroid)
203** returns double distance in meters
204*/
206Datum geography_distance(PG_FUNCTION_ARGS)
207{
208 GSERIALIZED* g1 = NULL;
209 GSERIALIZED* g2 = NULL;
210 double distance;
211 bool use_spheroid = true;
212 SPHEROID s;
213
214 /* Get our geometry objects loaded into memory. */
215 g1 = PG_GETARG_GSERIALIZED_P(0);
216 g2 = PG_GETARG_GSERIALIZED_P(1);
217
218 if (PG_NARGS() > 2)
219 use_spheroid = PG_GETARG_BOOL(2);
220
221 gserialized_error_if_srid_mismatch(g1, g2, __func__);
222
223 /* Initialize spheroid */
224 spheroid_init_from_srid(fcinfo, gserialized_get_srid(g1), &s);
225
226 /* Set to sphere if requested */
227 if ( ! use_spheroid )
228 s.a = s.b = s.radius;
229
230 /* Return NULL on empty arguments. */
232 {
233 PG_FREE_IF_COPY(g1, 0);
234 PG_FREE_IF_COPY(g2, 1);
235 PG_RETURN_NULL();
236 }
237
238 /* Do the brute force calculation if the cached calculation doesn't tick over */
239 if ( LW_FAILURE == geography_distance_cache(fcinfo, g1, g2, &s, &distance) )
240 {
241 /* default to using tree-based distance calculation at all times */
242 /* in standard distance call. */
244 /*
245 LWGEOM* lwgeom1 = lwgeom_from_gserialized(g1);
246 LWGEOM* lwgeom2 = lwgeom_from_gserialized(g2);
247 distance = lwgeom_distance_spheroid(lwgeom1, lwgeom2, &s, tolerance);
248 lwgeom_free(lwgeom1);
249 lwgeom_free(lwgeom2);
250 */
251 }
252
253 /* Clean up */
254 PG_FREE_IF_COPY(g1, 0);
255 PG_FREE_IF_COPY(g2, 1);
256
257 /* Knock off any funny business at the nanometer level, ticket #2168 */
259
260 /* Something went wrong, negative return... should already be eloged, return NULL */
261 if ( distance < 0.0 )
262 {
263 elog(ERROR, "distance returned negative!");
264 PG_RETURN_NULL();
265 }
266
267 PG_RETURN_FLOAT8(distance);
268}
269
270/*
271** geography_dwithin(GSERIALIZED *g1, GSERIALIZED *g2, double tolerance, boolean use_spheroid)
272** returns double distance in meters
273*/
275Datum geography_dwithin(PG_FUNCTION_ARGS)
276{
277 GSERIALIZED *g1 = PG_GETARG_GSERIALIZED_P(0);
278 GSERIALIZED *g2 = PG_GETARG_GSERIALIZED_P(1);
279 SPHEROID s;
280 double tolerance = 0.0;
281 bool use_spheroid = true;
282 double distance;
283 int dwithin = LW_FALSE;
284
285 gserialized_error_if_srid_mismatch(g1, g2, __func__);
286
287 /* Read our tolerance value. */
288 if ( PG_NARGS() > 2 && ! PG_ARGISNULL(2) )
289 tolerance = PG_GETARG_FLOAT8(2);
290
291 /* Read our calculation type. */
292 if ( PG_NARGS() > 3 && ! PG_ARGISNULL(3) )
293 use_spheroid = PG_GETARG_BOOL(3);
294
295 /* Initialize spheroid */
296 spheroid_init_from_srid(fcinfo, gserialized_get_srid(g1), &s);
297
298 /* Set to sphere if requested */
299 if ( ! use_spheroid )
300 s.a = s.b = s.radius;
301
302 /* Return FALSE on empty arguments. */
304 PG_RETURN_BOOL(false);
305
306 /* Do the brute force calculation if the cached calculation doesn't tick over */
307 if ( LW_FAILURE == geography_dwithin_cache(fcinfo, g1, g2, &s, tolerance, &dwithin) )
308 {
309 LWGEOM* lwgeom1 = lwgeom_from_gserialized(g1);
310 LWGEOM* lwgeom2 = lwgeom_from_gserialized(g2);
311 distance = lwgeom_distance_spheroid(lwgeom1, lwgeom2, &s, tolerance);
312 /* Something went wrong... */
313 if ( distance < 0.0 )
314 elog(ERROR, "lwgeom_distance_spheroid returned negative!");
315 dwithin = (distance <= tolerance);
316 lwgeom_free(lwgeom1);
317 lwgeom_free(lwgeom2);
318 }
319
320 PG_FREE_IF_COPY(g1, 0);
321 PG_FREE_IF_COPY(g2, 1);
322 PG_RETURN_BOOL(dwithin);
323}
324
326Datum geography_intersects(PG_FUNCTION_ARGS)
327{
328 PG_RETURN_BOOL(CallerFInfoFunctionCall2(
329 geography_dwithin, fcinfo->flinfo, InvalidOid,
330 PG_GETARG_DATUM(0), PG_GETARG_DATUM(1)));
331}
332
333/*
334** geography_distance_tree(GSERIALIZED *g1, GSERIALIZED *g2, double tolerance, boolean use_spheroid)
335** returns double distance in meters
336*/
338Datum geography_distance_tree(PG_FUNCTION_ARGS)
339{
340 GSERIALIZED *g1 = NULL;
341 GSERIALIZED *g2 = NULL;
342 double tolerance = 0.0;
343 double distance;
344 bool use_spheroid = true;
345 SPHEROID s;
346
347 /* Get our geometry objects loaded into memory. */
348 g1 = PG_GETARG_GSERIALIZED_P(0);
349 g2 = PG_GETARG_GSERIALIZED_P(1);
350
351 gserialized_error_if_srid_mismatch(g1, g2, __func__);
352
353 /* Return zero on empty arguments. */
355 {
356 PG_FREE_IF_COPY(g1, 0);
357 PG_FREE_IF_COPY(g2, 1);
358 PG_RETURN_FLOAT8(0.0);
359 }
360
361 /* Read our tolerance value. */
362 if ( PG_NARGS() > 2 && ! PG_ARGISNULL(2) )
363 tolerance = PG_GETARG_FLOAT8(2);
364
365 /* Read our calculation type. */
366 if ( PG_NARGS() > 3 && ! PG_ARGISNULL(3) )
367 use_spheroid = PG_GETARG_BOOL(3);
368
369 /* Initialize spheroid */
370 spheroid_init_from_srid(fcinfo, gserialized_get_srid(g1), &s);
371
372 /* Set to sphere if requested */
373 if ( ! use_spheroid )
374 s.a = s.b = s.radius;
375
376 if ( geography_tree_distance(g1, g2, &s, tolerance, &distance) == LW_FAILURE )
377 {
378 elog(ERROR, "geography_distance_tree failed!");
379 PG_RETURN_NULL();
380 }
381
382 /* Knock off any funny business at the nanometer level, ticket #2168 */
384
385 PG_RETURN_FLOAT8(distance);
386}
387
388
389
390/*
391** geography_dwithin_uncached(GSERIALIZED *g1, GSERIALIZED *g2, double tolerance, boolean use_spheroid)
392** returns double distance in meters
393*/
395Datum geography_dwithin_uncached(PG_FUNCTION_ARGS)
396{
397 LWGEOM *lwgeom1 = NULL;
398 LWGEOM *lwgeom2 = NULL;
399 GSERIALIZED *g1 = NULL;
400 GSERIALIZED *g2 = NULL;
401 double tolerance = 0.0;
402 double distance;
403 bool use_spheroid = true;
404 SPHEROID s;
405
406 /* Get our geometry objects loaded into memory. */
407 g1 = PG_GETARG_GSERIALIZED_P(0);
408 g2 = PG_GETARG_GSERIALIZED_P(1);
409 gserialized_error_if_srid_mismatch(g1, g2, __func__);
410
411 /* Read our tolerance value. */
412 if ( PG_NARGS() > 2 && ! PG_ARGISNULL(2) )
413 tolerance = PG_GETARG_FLOAT8(2);
414
415 /* Read our calculation type. */
416 if ( PG_NARGS() > 3 && ! PG_ARGISNULL(3) )
417 use_spheroid = PG_GETARG_BOOL(3);
418
419 /* Initialize spheroid */
420 spheroid_init_from_srid(fcinfo, gserialized_get_srid(g1), &s);
421
422 /* Set to sphere if requested */
423 if ( ! use_spheroid )
424 s.a = s.b = s.radius;
425
426 lwgeom1 = lwgeom_from_gserialized(g1);
427 lwgeom2 = lwgeom_from_gserialized(g2);
428
429 /* Return FALSE on empty arguments. */
430 if ( lwgeom_is_empty(lwgeom1) || lwgeom_is_empty(lwgeom2) )
431 {
432 PG_RETURN_BOOL(false);
433 }
434
435 distance = lwgeom_distance_spheroid(lwgeom1, lwgeom2, &s, tolerance);
436
437 /* Clean up */
438 lwgeom_free(lwgeom1);
439 lwgeom_free(lwgeom2);
440 PG_FREE_IF_COPY(g1, 0);
441 PG_FREE_IF_COPY(g2, 1);
442
443 /* Something went wrong... should already be eloged, return FALSE */
444 if ( distance < 0.0 )
445 {
446 elog(ERROR, "lwgeom_distance_spheroid returned negative!");
447 PG_RETURN_BOOL(false);
448 }
449
450 PG_RETURN_BOOL(distance <= tolerance);
451}
452
453
454/*
455** geography_expand(GSERIALIZED *g) returns *GSERIALIZED
456**
457** warning, this tricky little function does not expand the
458** geometry at all, just re-writes bounding box value to be
459** a bit bigger. only useful when passing the result along to
460** an index operator (&&)
461*/
463Datum geography_expand(PG_FUNCTION_ARGS)
464{
465 GSERIALIZED *g = NULL;
466 GSERIALIZED *g_out = NULL;
467 double unit_distance, distance;
468
469 /* Get a wholly-owned pointer to the geography */
470 g = PG_GETARG_GSERIALIZED_P_COPY(0);
471
472 /* Read our distance value and normalize to unit-sphere. */
473 distance = PG_GETARG_FLOAT8(1);
474 /* Magic 1% expansion is to bridge difference between potential */
475 /* spheroidal input distance and fact that expanded box filter is */
476 /* calculated on sphere */
477 unit_distance = 1.01 * distance / WGS84_RADIUS;
478
479 /* Try the expansion */
480 g_out = gserialized_expand(g, unit_distance);
481
482 /* If the expansion fails, the return our input */
483 if ( g_out == NULL )
484 {
485 PG_RETURN_POINTER(g);
486 }
487
488 if ( g_out != g )
489 {
490 pfree(g);
491 }
492
493 PG_RETURN_POINTER(g_out);
494}
495
496/*
497** geography_area(GSERIALIZED *g)
498** returns double area in meters square
499*/
501Datum geography_area(PG_FUNCTION_ARGS)
502{
503 LWGEOM *lwgeom = NULL;
504 GSERIALIZED *g = NULL;
505 GBOX gbox;
506 double area;
507 bool use_spheroid = LW_TRUE;
508 SPHEROID s;
509
510 /* Get our geometry object loaded into memory. */
511 g = PG_GETARG_GSERIALIZED_P(0);
512
513 /* Read our calculation type */
514 use_spheroid = PG_GETARG_BOOL(1);
515
516 /* Initialize spheroid */
517 spheroid_init_from_srid(fcinfo, gserialized_get_srid(g), &s);
518
519 lwgeom = lwgeom_from_gserialized(g);
520
521 /* EMPTY things have no area */
522 if ( lwgeom_is_empty(lwgeom) )
523 {
524 lwgeom_free(lwgeom);
525 PG_RETURN_FLOAT8(0.0);
526 }
527
528 if ( lwgeom->bbox )
529 gbox = *(lwgeom->bbox);
530 else
531 lwgeom_calculate_gbox_geodetic(lwgeom, &gbox);
532
533#ifndef PROJ_GEODESIC
534 /* Test for cases that are currently not handled by spheroid code */
535 if ( use_spheroid )
536 {
537 /* We can't circle the poles right now */
538 if ( FP_GTEQ(gbox.zmax,1.0) || FP_LTEQ(gbox.zmin,-1.0) )
539 use_spheroid = LW_FALSE;
540 /* We can't cross the equator right now */
541 if ( gbox.zmax > 0.0 && gbox.zmin < 0.0 )
542 use_spheroid = LW_FALSE;
543 }
544#endif /* ifndef PROJ_GEODESIC */
545
546 /* User requests spherical calculation, turn our spheroid into a sphere */
547 if ( ! use_spheroid )
548 s.a = s.b = s.radius;
549
550 /* Calculate the area */
551 if ( use_spheroid )
552 area = lwgeom_area_spheroid(lwgeom, &s);
553 else
554 area = lwgeom_area_sphere(lwgeom, &s);
555
556 /* Clean up */
557 lwgeom_free(lwgeom);
558 PG_FREE_IF_COPY(g, 0);
559
560 /* Something went wrong... */
561 if ( area < 0.0 )
562 {
563 elog(ERROR, "lwgeom_area_spher(oid) returned area < 0.0");
564 PG_RETURN_NULL();
565 }
566
567 PG_RETURN_FLOAT8(area);
568}
569
570/*
571** geography_perimeter(GSERIALIZED *g)
572** returns double perimeter in meters for area features
573*/
575Datum geography_perimeter(PG_FUNCTION_ARGS)
576{
577 LWGEOM *lwgeom = NULL;
578 GSERIALIZED *g = NULL;
579 double length;
580 bool use_spheroid = LW_TRUE;
581 SPHEROID s;
582 int type;
583
584 /* Get our geometry object loaded into memory. */
585 g = PG_GETARG_GSERIALIZED_P(0);
586
587 /* Only return for area features. */
588 type = gserialized_get_type(g);
589 if ( ! (type == POLYGONTYPE || type == MULTIPOLYGONTYPE || type == COLLECTIONTYPE) )
590 {
591 PG_RETURN_FLOAT8(0.0);
592 }
593
594 lwgeom = lwgeom_from_gserialized(g);
595
596 /* EMPTY things have no perimeter */
597 if ( lwgeom_is_empty(lwgeom) )
598 {
599 lwgeom_free(lwgeom);
600 PG_RETURN_FLOAT8(0.0);
601 }
602
603 /* Read our calculation type */
604 use_spheroid = PG_GETARG_BOOL(1);
605
606 /* Initialize spheroid */
607 spheroid_init_from_srid(fcinfo, gserialized_get_srid(g), &s);
608
609 /* User requests spherical calculation, turn our spheroid into a sphere */
610 if ( ! use_spheroid )
611 s.a = s.b = s.radius;
612
613 /* Calculate the length */
614 length = lwgeom_length_spheroid(lwgeom, &s);
615
616 /* Something went wrong... */
617 if ( length < 0.0 )
618 {
619 elog(ERROR, "lwgeom_length_spheroid returned length < 0.0");
620 PG_RETURN_NULL();
621 }
622
623 /* Clean up, but not all the way to the point arrays */
624 lwgeom_free(lwgeom);
625
626 PG_FREE_IF_COPY(g, 0);
627 PG_RETURN_FLOAT8(length);
628}
629
630/*
631** geography_length(GSERIALIZED *g)
632** returns double length in meters
633*/
635Datum geography_length(PG_FUNCTION_ARGS)
636{
637 LWGEOM *lwgeom = NULL;
638 GSERIALIZED *g = NULL;
639 double length;
640 bool use_spheroid = LW_TRUE;
641 SPHEROID s;
642
643 /* Get our geometry object loaded into memory. */
644 g = PG_GETARG_GSERIALIZED_P(0);
645 lwgeom = lwgeom_from_gserialized(g);
646
647 /* EMPTY things have no length */
648 if ( lwgeom_is_empty(lwgeom) || lwgeom->type == POLYGONTYPE || lwgeom->type == MULTIPOLYGONTYPE )
649 {
650 lwgeom_free(lwgeom);
651 PG_RETURN_FLOAT8(0.0);
652 }
653
654 /* Read our calculation type */
655 use_spheroid = PG_GETARG_BOOL(1);
656
657 /* Initialize spheroid */
658 spheroid_init_from_srid(fcinfo, gserialized_get_srid(g), &s);
659
660 /* User requests spherical calculation, turn our spheroid into a sphere */
661 if ( ! use_spheroid )
662 s.a = s.b = s.radius;
663
664 /* Calculate the length */
665 length = lwgeom_length_spheroid(lwgeom, &s);
666
667 /* Something went wrong... */
668 if ( length < 0.0 )
669 {
670 elog(ERROR, "lwgeom_length_spheroid returned length < 0.0");
671 PG_RETURN_NULL();
672 }
673
674 /* Clean up */
675 lwgeom_free(lwgeom);
676
677 PG_FREE_IF_COPY(g, 0);
678 PG_RETURN_FLOAT8(length);
679}
680
681/*
682** geography_point_outside(GSERIALIZED *g)
683** returns point outside the object
684*/
686Datum geography_point_outside(PG_FUNCTION_ARGS)
687{
688 GBOX gbox;
689 GSERIALIZED *g = NULL;
690 GSERIALIZED *g_out = NULL;
691 LWGEOM *lwpoint = NULL;
692 POINT2D pt;
693
694 /* Get our geometry object loaded into memory. */
695 g = PG_GETARG_GSERIALIZED_P(0);
696
697 /* We need the bounding box to get an outside point for area algorithm */
698 if ( gserialized_get_gbox_p(g, &gbox) == LW_FAILURE )
699 {
700 POSTGIS_DEBUG(4,"gserialized_get_gbox_p returned LW_FAILURE");
701 elog(ERROR, "Error in gserialized_get_gbox_p calculation.");
702 PG_RETURN_NULL();
703 }
704 POSTGIS_DEBUGF(4, "got gbox %s", gbox_to_string(&gbox));
705
706 /* Get an exterior point, based on this gbox */
707 gbox_pt_outside(&gbox, &pt);
708
709 lwpoint = (LWGEOM*) lwpoint_make2d(4326, pt.x, pt.y);
710
711 g_out = geography_serialize(lwpoint);
712
713 PG_FREE_IF_COPY(g, 0);
714 PG_RETURN_POINTER(g_out);
715
716}
717
718/*
719** geography_covers(GSERIALIZED *g, GSERIALIZED *g) returns boolean
720** Only works for (multi)points and (multi)polygons currently.
721** Attempts a simple point-in-polygon test on the polygon and point.
722** Current algorithm does not distinguish between points on edge
723** and points within.
724*/
726Datum geography_covers(PG_FUNCTION_ARGS)
727{
728 LWGEOM *lwgeom1 = NULL;
729 LWGEOM *lwgeom2 = NULL;
730 GSERIALIZED *g1 = NULL;
731 GSERIALIZED *g2 = NULL;
732 int result = LW_FALSE;
733
734 /* Get our geometry objects loaded into memory. */
735 g1 = PG_GETARG_GSERIALIZED_P(0);
736 g2 = PG_GETARG_GSERIALIZED_P(1);
737 gserialized_error_if_srid_mismatch(g1, g2, __func__);
738
739 /* Construct our working geometries */
740 lwgeom1 = lwgeom_from_gserialized(g1);
741 lwgeom2 = lwgeom_from_gserialized(g2);
742
743 /* EMPTY never intersects with another geometry */
744 if ( lwgeom_is_empty(lwgeom1) || lwgeom_is_empty(lwgeom2) )
745 {
746 lwgeom_free(lwgeom1);
747 lwgeom_free(lwgeom2);
748 PG_FREE_IF_COPY(g1, 0);
749 PG_FREE_IF_COPY(g2, 1);
750 PG_RETURN_BOOL(false);
751 }
752
753 /* Calculate answer */
754 result = lwgeom_covers_lwgeom_sphere(lwgeom1, lwgeom2);
755
756 /* Clean up */
757 lwgeom_free(lwgeom1);
758 lwgeom_free(lwgeom2);
759 PG_FREE_IF_COPY(g1, 0);
760 PG_FREE_IF_COPY(g2, 1);
761
762 PG_RETURN_BOOL(result);
763}
764
766Datum geography_coveredby(PG_FUNCTION_ARGS)
767{
768 LWGEOM *lwgeom1 = NULL;
769 LWGEOM *lwgeom2 = NULL;
770 GSERIALIZED *g1 = NULL;
771 GSERIALIZED *g2 = NULL;
772 int result = LW_FALSE;
773
774 /* Get our geometry objects loaded into memory. */
775 /* Pick them up in reverse order to covers */
776 g1 = PG_GETARG_GSERIALIZED_P(1);
777 g2 = PG_GETARG_GSERIALIZED_P(0);
778 gserialized_error_if_srid_mismatch(g1, g2, __func__);
779
780 /* Construct our working geometries */
781 lwgeom1 = lwgeom_from_gserialized(g1);
782 lwgeom2 = lwgeom_from_gserialized(g2);
783
784 /* EMPTY never intersects with another geometry */
785 if ( lwgeom_is_empty(lwgeom1) || lwgeom_is_empty(lwgeom2) )
786 {
787 lwgeom_free(lwgeom1);
788 lwgeom_free(lwgeom2);
789 PG_FREE_IF_COPY(g1, 1);
790 PG_FREE_IF_COPY(g2, 0);
791 PG_RETURN_BOOL(false);
792 }
793
794 /* Calculate answer */
795 result = lwgeom_covers_lwgeom_sphere(lwgeom1, lwgeom2);
796
797 /* Clean up */
798 lwgeom_free(lwgeom1);
799 lwgeom_free(lwgeom2);
800 PG_FREE_IF_COPY(g1, 1);
801 PG_FREE_IF_COPY(g2, 0);
802
803 PG_RETURN_BOOL(result);
804}
805
806/*
807** geography_bestsrid(GSERIALIZED *g, GSERIALIZED *g) returns int
808** Utility function. Returns negative SRID numbers that match to the
809** numbers handled in code by the transform(lwgeom, srid) function.
810** UTM, polar stereographic and mercator as fallback. To be used
811** in wrapping existing geometry functions in SQL to provide access
812** to them in the geography module.
813*/
815Datum geography_bestsrid(PG_FUNCTION_ARGS)
816{
817 GBOX gbox, gbox1, gbox2;
818 GSERIALIZED *g1 = NULL;
819 GSERIALIZED *g2 = NULL;
820 int empty1 = LW_FALSE;
821 int empty2 = LW_FALSE;
822 double xwidth, ywidth;
823 POINT2D center;
824 LWGEOM *lwgeom;
825
826 /* Get our geometry objects loaded into memory. */
827 g1 = PG_GETARG_GSERIALIZED_P(0);
828 /* Synchronize our box types */
829 gbox1.flags = gserialized_get_lwflags(g1);
830 /* Calculate if the geometry is empty. */
831 empty1 = gserialized_is_empty(g1);
832
833 /* Convert g1 to LWGEOM type */
834 lwgeom = lwgeom_from_gserialized(g1);
835
836 /* Calculate a geocentric bounds for the objects */
837 if ( ! empty1 && gserialized_get_gbox_p(g1, &gbox1) == LW_FAILURE )
838 elog(ERROR, "Error in geography_bestsrid calling gserialized_get_gbox_p(g1, &gbox1)");
839
840 POSTGIS_DEBUGF(4, "calculated gbox = %s", gbox_to_string(&gbox1));
841
842 if ( !lwgeom_isfinite(lwgeom) ) {
843 elog(ERROR, "Error in geography_bestsrid calling with infinite coordinate geographies");
844 }
845 lwgeom_free(lwgeom);
846
847 /* If we have a unique second argument, fill in all the necessary variables. */
848 if (PG_NARGS() > 1)
849 {
850 g2 = PG_GETARG_GSERIALIZED_P(1);
851 gbox2.flags = gserialized_get_lwflags(g2);
852 empty2 = gserialized_is_empty(g2);
853 if ( ! empty2 && gserialized_get_gbox_p(g2, &gbox2) == LW_FAILURE )
854 elog(ERROR, "Error in geography_bestsrid calling gserialized_get_gbox_p(g2, &gbox2)");
855
856 /* Convert g2 to LWGEOM type */
857 lwgeom = lwgeom_from_gserialized(g2);
858
859 if ( !lwgeom_isfinite(lwgeom) ) {
860 elog(ERROR, "Error in geography_bestsrid calling with second arg infinite coordinate geographies");
861 }
862 lwgeom_free(lwgeom);
863 }
864 /*
865 ** If no unique second argument, copying the box for the first
866 ** argument will give us the right answer for all subsequent tests.
867 */
868 else
869 {
870 gbox = gbox2 = gbox1;
871 }
872
873 /* Both empty? We don't have an answer. */
874 if ( empty1 && empty2 )
875 PG_RETURN_NULL();
876
877 /* One empty? We can use the other argument values as infill. Otherwise merge the boxen */
878 if ( empty1 )
879 gbox = gbox2;
880 else if ( empty2 )
881 gbox = gbox1;
882 else
883 gbox_union(&gbox1, &gbox2, &gbox);
884
885 gbox_centroid(&gbox, &center);
886
887 /* Width and height in degrees */
888 xwidth = 180.0 * gbox_angular_width(&gbox) / M_PI;
889 ywidth = 180.0 * gbox_angular_height(&gbox) / M_PI;
890
891 POSTGIS_DEBUGF(2, "xwidth %g", xwidth);
892 POSTGIS_DEBUGF(2, "ywidth %g", ywidth);
893 POSTGIS_DEBUGF(2, "center POINT(%g %g)", center.x, center.y);
894
895 /* Are these data arctic? Lambert Azimuthal Equal Area North. */
896 if ( center.y > 70.0 && ywidth < 45.0 )
897 {
898 PG_RETURN_INT32(SRID_NORTH_LAMBERT);
899 }
900
901 /* Are these data antarctic? Lambert Azimuthal Equal Area South. */
902 if ( center.y < -70.0 && ywidth < 45.0 )
903 {
904 PG_RETURN_INT32(SRID_SOUTH_LAMBERT);
905 }
906
907 /*
908 ** Can we fit these data into one UTM zone?
909 ** We will assume we can push things as
910 ** far as a half zone past a zone boundary.
911 ** Note we have no handling for the date line in here.
912 */
913 if ( xwidth < 6.0 )
914 {
915 int zone = floor((center.x + 180.0) / 6.0);
916
917 if ( zone > 59 ) zone = 59;
918
919 /* Are these data below the equator? UTM South. */
920 if ( center.y < 0.0 )
921 {
922 PG_RETURN_INT32( SRID_SOUTH_UTM_START + zone );
923 }
924 /* Are these data above the equator? UTM North. */
925 else
926 {
927 PG_RETURN_INT32( SRID_NORTH_UTM_START + zone );
928 }
929 }
930
931 /*
932 ** Can we fit into a custom LAEA area? (30 degrees high, variable width)
933 ** We will allow overlap into adjoining areas, but use a slightly narrower test (25) to try
934 ** and minimize the worst case.
935 ** Again, we are hoping the dateline doesn't trip us up much
936 */
937 if ( ywidth < 25.0 )
938 {
939 int xzone = -1;
940 int yzone = 3 + floor(center.y / 30.0); /* (range of 0-5) */
941
942 /* Equatorial band, 12 zones, 30 degrees wide */
943 if ( (yzone == 2 || yzone == 3) && xwidth < 30.0 )
944 {
945 xzone = 6 + floor(center.x / 30.0);
946 }
947 /* Temperate band, 8 zones, 45 degrees wide */
948 else if ( (yzone == 1 || yzone == 4) && xwidth < 45.0 )
949 {
950 xzone = 4 + floor(center.x / 45.0);
951 }
952 /* Arctic band, 4 zones, 90 degrees wide */
953 else if ( (yzone == 0 || yzone == 5) && xwidth < 90.0 )
954 {
955 xzone = 2 + floor(center.x / 90.0);
956 }
957
958 /* Did we fit into an appropriate xzone? */
959 if ( xzone != -1 )
960 {
961 PG_RETURN_INT32(SRID_LAEA_START + 20 * yzone + xzone);
962 }
963 }
964
965 /*
966 ** Running out of options... fall-back to Mercator
967 ** and hope for the best.
968 */
969 PG_RETURN_INT32(SRID_WORLD_MERCATOR);
970
971}
972
973/*
974** geography_project(GSERIALIZED *g, distance, azimuth)
975** returns point of projection given start point,
976** azimuth in radians (bearing) and distance in meters
977*/
979Datum geography_project(PG_FUNCTION_ARGS)
980{
981 LWGEOM *lwgeom = NULL;
982 LWPOINT *lwp_projected;
983 GSERIALIZED *g = NULL;
984 GSERIALIZED *g_out = NULL;
985 double azimuth = 0.0;
986 double distance;
987 SPHEROID s;
988 uint32_t type;
989
990 /* Return NULL on NULL distance or geography */
991 if ( PG_NARGS() < 2 || PG_ARGISNULL(0) || PG_ARGISNULL(1) )
992 PG_RETURN_NULL();
993
994 /* Get our geometry object loaded into memory. */
995 g = PG_GETARG_GSERIALIZED_P(0);
996
997 /* Only return for points. */
998 type = gserialized_get_type(g);
999 if ( type != POINTTYPE )
1000 {
1001 elog(ERROR, "ST_Project(geography) is only valid for point inputs");
1002 PG_RETURN_NULL();
1003 }
1004
1005 distance = PG_GETARG_FLOAT8(1); /* Distance in Meters */
1006 lwgeom = lwgeom_from_gserialized(g);
1007
1008 /* EMPTY things cannot be projected from */
1009 if ( lwgeom_is_empty(lwgeom) )
1010 {
1011 lwgeom_free(lwgeom);
1012 elog(ERROR, "ST_Project(geography) cannot project from an empty start point");
1013 PG_RETURN_NULL();
1014 }
1015
1016 if ( PG_NARGS() > 2 && ! PG_ARGISNULL(2) )
1017 azimuth = PG_GETARG_FLOAT8(2); /* Azimuth in Radians */
1018
1019 /* Initialize spheroid */
1020 spheroid_init_from_srid(fcinfo, gserialized_get_srid(g), &s);
1021
1022 /* Handle the zero distance case */
1023 if( FP_EQUALS(distance, 0.0) )
1024 {
1025 PG_RETURN_POINTER(g);
1026 }
1027
1028 /* Calculate the length */
1029 lwp_projected = lwgeom_project_spheroid(lwgeom_as_lwpoint(lwgeom), &s, distance, azimuth);
1030
1031 /* Something went wrong... */
1032 if ( lwp_projected == NULL )
1033 {
1034 elog(ERROR, "lwgeom_project_spheroid returned null");
1035 PG_RETURN_NULL();
1036 }
1037
1038 /* Clean up, but not all the way to the point arrays */
1039 lwgeom_free(lwgeom);
1040 g_out = geography_serialize(lwpoint_as_lwgeom(lwp_projected));
1041 lwpoint_free(lwp_projected);
1042
1043 PG_FREE_IF_COPY(g, 0);
1044 PG_RETURN_POINTER(g_out);
1045}
1046
1047
1048/*
1049** geography_azimuth(GSERIALIZED *g1, GSERIALIZED *g2)
1050** returns direction between points (north = 0)
1051** azimuth (bearing) and distance
1052*/
1054Datum geography_azimuth(PG_FUNCTION_ARGS)
1055{
1056 LWGEOM *lwgeom1 = NULL;
1057 LWGEOM *lwgeom2 = NULL;
1058 GSERIALIZED *g1 = NULL;
1059 GSERIALIZED *g2 = NULL;
1060 double azimuth;
1061 SPHEROID s;
1062 uint32_t type1, type2;
1063
1064 /* Get our geometry object loaded into memory. */
1065 g1 = PG_GETARG_GSERIALIZED_P(0);
1066 g2 = PG_GETARG_GSERIALIZED_P(1);
1067
1068 /* Only return for points. */
1069 type1 = gserialized_get_type(g1);
1070 type2 = gserialized_get_type(g2);
1071 if ( type1 != POINTTYPE || type2 != POINTTYPE )
1072 {
1073 elog(ERROR, "ST_Azimuth(geography, geography) is only valid for point inputs");
1074 PG_RETURN_NULL();
1075 }
1076
1077 lwgeom1 = lwgeom_from_gserialized(g1);
1078 lwgeom2 = lwgeom_from_gserialized(g2);
1079
1080 /* EMPTY things cannot be used */
1081 if ( lwgeom_is_empty(lwgeom1) || lwgeom_is_empty(lwgeom2) )
1082 {
1083 lwgeom_free(lwgeom1);
1084 lwgeom_free(lwgeom2);
1085 elog(ERROR, "ST_Azimuth(geography, geography) cannot work with empty points");
1086 PG_RETURN_NULL();
1087 }
1088
1089 /* Initialize spheroid */
1090 spheroid_init_from_srid(fcinfo, gserialized_get_srid(g1), &s);
1091
1092 /* Calculate the direction */
1093 azimuth = lwgeom_azumith_spheroid(lwgeom_as_lwpoint(lwgeom1), lwgeom_as_lwpoint(lwgeom2), &s);
1094
1095 /* Clean up */
1096 lwgeom_free(lwgeom1);
1097 lwgeom_free(lwgeom2);
1098
1099 PG_FREE_IF_COPY(g1, 0);
1100 PG_FREE_IF_COPY(g2, 1);
1101
1102 /* Return NULL for unknown (same point) azimuth */
1103 if( isnan(azimuth) )
1104 {
1105 PG_RETURN_NULL();
1106 }
1107
1108 PG_RETURN_FLOAT8(azimuth);
1109}
1110
1111
1112
1113/*
1114** geography_segmentize(GSERIALIZED *g1, double max_seg_length)
1115** returns densified geometry with no segment longer than max
1116*/
1118Datum geography_segmentize(PG_FUNCTION_ARGS)
1119{
1120 LWGEOM *lwgeom1 = NULL;
1121 LWGEOM *lwgeom2 = NULL;
1122 GSERIALIZED *g1 = NULL;
1123 GSERIALIZED *g2 = NULL;
1124 double max_seg_length;
1125 uint32_t type1;
1126
1127 /* Get our geometry object loaded into memory. */
1128 g1 = PG_GETARG_GSERIALIZED_P(0);
1129 type1 = gserialized_get_type(g1);
1130
1131 /* Convert max_seg_length from metric units to radians */
1132 max_seg_length = PG_GETARG_FLOAT8(1) / WGS84_RADIUS;
1133
1134 /* We can't densify points or points, reflect them back */
1135 if ( type1 == POINTTYPE || type1 == MULTIPOINTTYPE || gserialized_is_empty(g1) )
1136 PG_RETURN_POINTER(g1);
1137
1138 /* Deserialize */
1139 lwgeom1 = lwgeom_from_gserialized(g1);
1140
1141 /* Calculate the densified geometry */
1142 lwgeom2 = lwgeom_segmentize_sphere(lwgeom1, max_seg_length);
1143
1144 /*
1145 ** Set the geodetic flag so subsequent
1146 ** functions do the right thing.
1147 */
1148 lwgeom_set_geodetic(lwgeom2, true);
1149
1150 /* Recalculate the boxes after re-setting the geodetic bit */
1151 lwgeom_drop_bbox(lwgeom2);
1152
1153 /* We are trusting geography_serialize will add a box if needed */
1154 g2 = geography_serialize(lwgeom2);
1155
1156 /* Clean up */
1157 lwgeom_free(lwgeom1);
1158 lwgeom_free(lwgeom2);
1159 PG_FREE_IF_COPY(g1, 0);
1160
1161 PG_RETURN_POINTER(g2);
1162}
1163
1164
char * s
Definition cu_in_wkt.c:23
int gbox_union(const GBOX *g1, const GBOX *g2, GBOX *gout)
Update the output GBOX to be large enough to include both inputs.
Definition gbox.c:135
char * gbox_to_string(const GBOX *gbox)
Allocate a string representation of the GBOX, based on dimensionality of flags.
Definition gbox.c:392
GSERIALIZED * gserialized_expand(GSERIALIZED *g, double distance)
Return a GSERIALIZED with an expanded bounding box.
Datum geography_length(PG_FUNCTION_ARGS)
Datum geography_area(PG_FUNCTION_ARGS)
Datum geography_expand(PG_FUNCTION_ARGS)
PG_FUNCTION_INFO_V1(geography_distance_knn)
Datum geography_azimuth(PG_FUNCTION_ARGS)
Datum geography_distance_tree(PG_FUNCTION_ARGS)
Datum geography_point_outside(PG_FUNCTION_ARGS)
Datum geography_distance_uncached(PG_FUNCTION_ARGS)
Datum geography_distance(PG_FUNCTION_ARGS)
Datum geography_segmentize(PG_FUNCTION_ARGS)
Datum geography_covers(PG_FUNCTION_ARGS)
Datum geography_intersects(PG_FUNCTION_ARGS)
Datum geography_dwithin(PG_FUNCTION_ARGS)
Datum geography_coveredby(PG_FUNCTION_ARGS)
Datum geography_perimeter(PG_FUNCTION_ARGS)
Datum geography_bestsrid(PG_FUNCTION_ARGS)
Datum geography_dwithin_uncached(PG_FUNCTION_ARGS)
Datum geography_project(PG_FUNCTION_ARGS)
#define INVMINDIST
Datum geography_distance_knn(PG_FUNCTION_ARGS)
int geography_dwithin_cache(FunctionCallInfo fcinfo, const GSERIALIZED *g1, const GSERIALIZED *g2, const SPHEROID *s, double tolerance, int *dwithin)
int geography_tree_distance(const GSERIALIZED *g1, const GSERIALIZED *g2, const SPHEROID *s, double tolerance, double *distance)
int geography_distance_cache(FunctionCallInfo fcinfo, const GSERIALIZED *g1, const GSERIALIZED *g2, const SPHEROID *s, double *distance)
void gserialized_error_if_srid_mismatch(const GSERIALIZED *g1, const GSERIALIZED *g2, const char *funcname)
int32_t gserialized_get_srid(const GSERIALIZED *g)
Extract the SRID from the serialized form (it is packed into three bytes so this is a handy function)...
LWGEOM * lwgeom_from_gserialized(const GSERIALIZED *g)
Allocate a new LWGEOM from a GSERIALIZED.
int gserialized_get_gbox_p(const GSERIALIZED *g, GBOX *gbox)
Read the box from the GSERIALIZED or calculate it if necessary.
Definition gserialized.c:65
int gserialized_is_empty(const GSERIALIZED *g)
Check if a GSERIALIZED is empty without deserializing first.
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:89
lwflags_t gserialized_get_lwflags(const GSERIALIZED *g)
Read the flags from a GSERIALIZED and return a standard lwflag integer.
Definition gserialized.c:18
void lwgeom_set_geodetic(LWGEOM *geom, int value)
Set the FLAGS geodetic bit on geometry an all sub-geometries and pointlists.
Definition lwgeom.c:946
LWGEOM * lwpoint_as_lwgeom(const LWPOINT *obj)
Definition lwgeom.c:326
#define LW_FALSE
Definition liblwgeom.h:108
#define COLLECTIONTYPE
Definition liblwgeom.h:122
void lwpoint_free(LWPOINT *pt)
Definition lwpoint.c:213
double lwgeom_distance_spheroid(const LWGEOM *lwgeom1, const LWGEOM *lwgeom2, const SPHEROID *spheroid, double tolerance)
Calculate the geodetic distance from lwgeom1 to lwgeom2 on the spheroid.
#define LW_FAILURE
Definition liblwgeom.h:110
void lwgeom_free(LWGEOM *geom)
Definition lwgeom.c:1138
int gbox_pt_outside(const GBOX *gbox, POINT2D *pt_outside)
Calculate a spherical point that falls outside the geocentric gbox.
int lwgeom_calculate_gbox_geodetic(const LWGEOM *geom, GBOX *gbox)
Calculate the geodetic bounding box for an LWGEOM.
#define WGS84_RADIUS
Definition liblwgeom.h:162
double lwgeom_area_spheroid(const LWGEOM *lwgeom, const SPHEROID *spheroid)
Calculate the geodetic area of a lwgeom on the spheroid.
Definition lwspheroid.c:647
#define MULTIPOINTTYPE
Definition liblwgeom.h:119
int lwgeom_covers_lwgeom_sphere(const LWGEOM *lwgeom1, const LWGEOM *lwgeom2)
Calculate covers predicate for two lwgeoms on the sphere.
#define POINTTYPE
LWTYPE numbers, used internally by PostGIS.
Definition liblwgeom.h:116
void lwgeom_drop_bbox(LWGEOM *lwgeom)
Call this function to drop BBOX and SRID from LWGEOM.
Definition lwgeom.c:664
int lwgeom_isfinite(const LWGEOM *lwgeom)
Check if a LWGEOM has any non-finite (NaN or Inf) coordinates.
Definition lwgeom.c:2529
#define MULTIPOLYGONTYPE
Definition liblwgeom.h:121
LWPOINT * lwpoint_make2d(int32_t srid, double x, double y)
Definition lwpoint.c:163
#define POLYGONTYPE
Definition liblwgeom.h:118
void lwgeom_add_bbox_deep(LWGEOM *lwgeom, GBOX *gbox)
Compute a box for geom and all sub-geometries, if not already computed.
Definition lwgeom.c:696
double lwgeom_length_spheroid(const LWGEOM *geom, const SPHEROID *s)
Calculate the geodetic length of a lwgeom on the unit sphere.
LWPOINT * lwgeom_project_spheroid(const LWPOINT *r, const SPHEROID *spheroid, double distance, double azimuth)
Calculate the location of a point on a spheroid, give a start point, bearing and distance.
#define LW_TRUE
Return types for functions with status returns.
Definition liblwgeom.h:107
double lwgeom_azumith_spheroid(const LWPOINT *r, const LWPOINT *s, const SPHEROID *spheroid)
Calculate the bearing between two points on a spheroid.
LWGEOM * lwgeom_segmentize_sphere(const LWGEOM *lwg_in, double max_seg_length)
Derive a new geometry with vertices added to ensure no vertex is more than max_seg_length (in radians...
double lwgeom_area_sphere(const LWGEOM *lwgeom, const SPHEROID *spheroid)
Calculate the geodetic area of a lwgeom on the sphere.
This library is the generic geometry handling section of PostGIS.
#define FP_GTEQ(A, B)
double gbox_angular_height(const GBOX *gbox)
GBOX utility functions to figure out coverage/location on the globe.
Definition lwgeodetic.c:188
double gbox_angular_width(const GBOX *gbox)
Returns the angular width (longitudinal span) of the box in radians.
Definition lwgeodetic.c:215
#define FP_TOLERANCE
Floating point comparators.
#define FP_EQUALS(A, B)
int gbox_centroid(const GBOX *gbox, POINT2D *out)
Computes the average(ish) center of the box and returns success.
Definition lwgeodetic.c:267
#define FP_LTEQ(A, B)
static LWPOINT * lwgeom_as_lwpoint(const LWGEOM *lwgeom)
Definition lwinline.h:121
static int lwgeom_is_empty(const LWGEOM *geom)
Return true or false depending on whether a geometry is an "empty" geometry (no vertices members)
Definition lwinline.h:193
static double distance(double x1, double y1, double x2, double y2)
Definition lwtree.c:1032
double zmax
Definition liblwgeom.h:345
double zmin
Definition liblwgeom.h:344
lwflags_t flags
Definition liblwgeom.h:339
uint8_t type
Definition liblwgeom.h:448
GBOX * bbox
Definition liblwgeom.h:444
double y
Definition liblwgeom.h:376
double x
Definition liblwgeom.h:376