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