PostGIS  2.4.9dev-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  double tolerance = 0.0;
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  /* Read our tolerance value. */
219  if ( PG_NARGS() > 2 && ! PG_ARGISNULL(2) )
220  tolerance = PG_GETARG_FLOAT8(2);
221 
222  /* Read our calculation type. */
223  if ( PG_NARGS() > 3 && ! PG_ARGISNULL(3) )
224  use_spheroid = PG_GETARG_BOOL(3);
225 
227 
228  /* Initialize spheroid */
229  spheroid_init_from_srid(fcinfo, gserialized_get_srid(g1), &s);
230 
231  /* Set to sphere if requested */
232  if ( ! use_spheroid )
233  s.a = s.b = s.radius;
234 
235  /* Return NULL on empty arguments. */
237  {
238  PG_FREE_IF_COPY(g1, 0);
239  PG_FREE_IF_COPY(g2, 1);
240  PG_RETURN_NULL();
241  }
242 
243  /* Do the brute force calculation if the cached calculation doesn't tick over */
244  if ( LW_FAILURE == geography_distance_cache(fcinfo, g1, g2, &s, &distance) )
245  {
246  LWGEOM* lwgeom1 = lwgeom_from_gserialized(g1);
247  LWGEOM* lwgeom2 = lwgeom_from_gserialized(g2);
248  distance = lwgeom_distance_spheroid(lwgeom1, lwgeom2, &s, tolerance);
249  lwgeom_free(lwgeom1);
250  lwgeom_free(lwgeom2);
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 type1, type2;
740  int result = LW_FALSE;
741 
742  /* Get our geometry objects loaded into memory. */
743  g1 = PG_GETARG_GSERIALIZED_P(0);
744  g2 = PG_GETARG_GSERIALIZED_P(1);
745 
746  type1 = gserialized_get_type(g1);
747  type2 = gserialized_get_type(g2);
748 
749  /* Construct our working geometries */
750  lwgeom1 = lwgeom_from_gserialized(g1);
751  lwgeom2 = lwgeom_from_gserialized(g2);
752 
753  error_if_srid_mismatch(lwgeom1->srid, lwgeom2->srid);
754 
755  /* EMPTY never intersects with another geometry */
756  if ( lwgeom_is_empty(lwgeom1) || lwgeom_is_empty(lwgeom2) )
757  {
758  lwgeom_free(lwgeom1);
759  lwgeom_free(lwgeom2);
760  PG_FREE_IF_COPY(g1, 0);
761  PG_FREE_IF_COPY(g2, 1);
762  PG_RETURN_BOOL(false);
763  }
764 
765  /* Calculate answer */
766  result = lwgeom_covers_lwgeom_sphere(lwgeom1, lwgeom2);
767 
768  /* Clean up */
769  lwgeom_free(lwgeom1);
770  lwgeom_free(lwgeom2);
771  PG_FREE_IF_COPY(g1, 0);
772  PG_FREE_IF_COPY(g2, 1);
773 
774  PG_RETURN_BOOL(result);
775 }
776 
777 /*
778 ** geography_bestsrid(GSERIALIZED *g, GSERIALIZED *g) returns int
779 ** Utility function. Returns negative SRID numbers that match to the
780 ** numbers handled in code by the transform(lwgeom, srid) function.
781 ** UTM, polar stereographic and mercator as fallback. To be used
782 ** in wrapping existing geometry functions in SQL to provide access
783 ** to them in the geography module.
784 */
786 Datum geography_bestsrid(PG_FUNCTION_ARGS)
787 {
788  GBOX gbox, gbox1, gbox2;
789  GSERIALIZED *g1 = NULL;
790  GSERIALIZED *g2 = NULL;
791  int empty1 = LW_FALSE;
792  int empty2 = LW_FALSE;
793  double xwidth, ywidth;
794  POINT2D center;
795 
796  Datum d1 = PG_GETARG_DATUM(0);
797  Datum d2 = PG_GETARG_DATUM(1);
798 
799  /* Get our geometry objects loaded into memory. */
800  g1 = (GSERIALIZED*)PG_DETOAST_DATUM(d1);
801  /* Synchronize our box types */
802  gbox1.flags = g1->flags;
803  /* Calculate if the geometry is empty. */
804  empty1 = gserialized_is_empty(g1);
805  /* Calculate a geocentric bounds for the objects */
806  if ( ! empty1 && gserialized_get_gbox_p(g1, &gbox1) == LW_FAILURE )
807  elog(ERROR, "Error in geography_bestsrid calling gserialized_get_gbox_p(g1, &gbox1)");
808 
809  POSTGIS_DEBUGF(4, "calculated gbox = %s", gbox_to_string(&gbox1));
810 
811  /* If we have a unique second argument, fill in all the necessary variables. */
812  if ( d1 != d2 )
813  {
814  g2 = (GSERIALIZED*)PG_DETOAST_DATUM(d2);
815  gbox2.flags = g2->flags;
816  empty2 = gserialized_is_empty(g2);
817  if ( ! empty2 && gserialized_get_gbox_p(g2, &gbox2) == LW_FAILURE )
818  elog(ERROR, "Error in geography_bestsrid calling gserialized_get_gbox_p(g2, &gbox2)");
819  }
820  /*
821  ** If no unique second argument, copying the box for the first
822  ** argument will give us the right answer for all subsequent tests.
823  */
824  else
825  {
826  gbox = gbox2 = gbox1;
827  }
828 
829  /* Both empty? We don't have an answer. */
830  if ( empty1 && empty2 )
831  PG_RETURN_NULL();
832 
833  /* One empty? We can use the other argument values as infill. Otherwise merge the boxen */
834  if ( empty1 )
835  gbox = gbox2;
836  else if ( empty2 )
837  gbox = gbox1;
838  else
839  gbox_union(&gbox1, &gbox2, &gbox);
840 
841  gbox_centroid(&gbox, &center);
842 
843  /* Width and height in degrees */
844  xwidth = 180.0 * gbox_angular_width(&gbox) / M_PI;
845  ywidth = 180.0 * gbox_angular_height(&gbox) / M_PI;
846 
847  POSTGIS_DEBUGF(2, "xwidth %g", xwidth);
848  POSTGIS_DEBUGF(2, "ywidth %g", ywidth);
849  POSTGIS_DEBUGF(2, "center POINT(%g %g)", center.x, center.y);
850 
851  /* Are these data arctic? Lambert Azimuthal Equal Area North. */
852  if ( center.y > 70.0 && ywidth < 45.0 )
853  {
854  PG_RETURN_INT32(SRID_NORTH_LAMBERT);
855  }
856 
857  /* Are these data antarctic? Lambert Azimuthal Equal Area South. */
858  if ( center.y < -70.0 && ywidth < 45.0 )
859  {
860  PG_RETURN_INT32(SRID_SOUTH_LAMBERT);
861  }
862 
863  /*
864  ** Can we fit these data into one UTM zone?
865  ** We will assume we can push things as
866  ** far as a half zone past a zone boundary.
867  ** Note we have no handling for the date line in here.
868  */
869  if ( xwidth < 6.0 )
870  {
871  int zone = floor((center.x + 180.0) / 6.0);
872 
873  if ( zone > 59 ) zone = 59;
874 
875  /* Are these data below the equator? UTM South. */
876  if ( center.y < 0.0 )
877  {
878  PG_RETURN_INT32( SRID_SOUTH_UTM_START + zone );
879  }
880  /* Are these data above the equator? UTM North. */
881  else
882  {
883  PG_RETURN_INT32( SRID_NORTH_UTM_START + zone );
884  }
885  }
886 
887  /*
888  ** Can we fit into a custom LAEA area? (30 degrees high, variable width)
889  ** We will allow overlap into adjoining areas, but use a slightly narrower test (25) to try
890  ** and minimize the worst case.
891  ** Again, we are hoping the dateline doesn't trip us up much
892  */
893  if ( ywidth < 25.0 )
894  {
895  int xzone = -1;
896  int yzone = 3 + floor(center.y / 30.0); /* (range of 0-5) */
897 
898  /* Equatorial band, 12 zones, 30 degrees wide */
899  if ( (yzone == 2 || yzone == 3) && xwidth < 30.0 )
900  {
901  xzone = 6 + floor(center.x / 30.0);
902  }
903  /* Temperate band, 8 zones, 45 degrees wide */
904  else if ( (yzone == 1 || yzone == 4) && xwidth < 45.0 )
905  {
906  xzone = 4 + floor(center.x / 45.0);
907  }
908  /* Arctic band, 4 zones, 90 degrees wide */
909  else if ( (yzone == 0 || yzone == 5) && xwidth < 90.0 )
910  {
911  xzone = 2 + floor(center.x / 90.0);
912  }
913 
914  /* Did we fit into an appropriate xzone? */
915  if ( xzone != -1 )
916  {
917  PG_RETURN_INT32(SRID_LAEA_START + 20 * yzone + xzone);
918  }
919  }
920 
921  /*
922  ** Running out of options... fall-back to Mercator
923  ** and hope for the best.
924  */
925  PG_RETURN_INT32(SRID_WORLD_MERCATOR);
926 
927 }
928 
929 /*
930 ** geography_project(GSERIALIZED *g, distance, azimuth)
931 ** returns point of projection given start point,
932 ** azimuth in radians (bearing) and distance in meters
933 */
935 Datum geography_project(PG_FUNCTION_ARGS)
936 {
937  LWGEOM *lwgeom = NULL;
938  LWPOINT *lwp_projected;
939  GSERIALIZED *g = NULL;
940  GSERIALIZED *g_out = NULL;
941  double azimuth = 0.0;
942  double distance;
943  SPHEROID s;
944  uint32_t type;
945 
946  /* Return NULL on NULL distance or geography */
947  if ( PG_NARGS() < 2 || PG_ARGISNULL(0) || PG_ARGISNULL(1) )
948  PG_RETURN_NULL();
949 
950  /* Get our geometry object loaded into memory. */
951  g = PG_GETARG_GSERIALIZED_P(0);
952 
953  /* Only return for points. */
954  type = gserialized_get_type(g);
955  if ( type != POINTTYPE )
956  {
957  elog(ERROR, "ST_Project(geography) is only valid for point inputs");
958  PG_RETURN_NULL();
959  }
960 
961  distance = PG_GETARG_FLOAT8(1); /* Distance in Meters */
962  lwgeom = lwgeom_from_gserialized(g);
963 
964  /* EMPTY things cannot be projected from */
965  if ( lwgeom_is_empty(lwgeom) )
966  {
967  lwgeom_free(lwgeom);
968  elog(ERROR, "ST_Project(geography) cannot project from an empty start point");
969  PG_RETURN_NULL();
970  }
971 
972  if ( PG_NARGS() > 2 && ! PG_ARGISNULL(2) )
973  azimuth = PG_GETARG_FLOAT8(2); /* Azimuth in Radians */
974 
975  /* Initialize spheroid */
976  spheroid_init_from_srid(fcinfo, gserialized_get_srid(g), &s);
977 
978  /* Handle the zero distance case */
979  if( FP_EQUALS(distance, 0.0) )
980  {
981  PG_RETURN_POINTER(g);
982  }
983 
984  /* Calculate the length */
985  lwp_projected = lwgeom_project_spheroid(lwgeom_as_lwpoint(lwgeom), &s, distance, azimuth);
986 
987  /* Something went wrong... */
988  if ( lwp_projected == NULL )
989  {
990  elog(ERROR, "lwgeom_project_spheroid returned null");
991  PG_RETURN_NULL();
992  }
993 
994  /* Clean up, but not all the way to the point arrays */
995  lwgeom_free(lwgeom);
996  g_out = geography_serialize(lwpoint_as_lwgeom(lwp_projected));
997  lwpoint_free(lwp_projected);
998 
999  PG_FREE_IF_COPY(g, 0);
1000  PG_RETURN_POINTER(g_out);
1001 }
1002 
1003 
1004 /*
1005 ** geography_azimuth(GSERIALIZED *g1, GSERIALIZED *g2)
1006 ** returns direction between points (north = 0)
1007 ** azimuth (bearing) and distance
1008 */
1010 Datum geography_azimuth(PG_FUNCTION_ARGS)
1011 {
1012  LWGEOM *lwgeom1 = NULL;
1013  LWGEOM *lwgeom2 = NULL;
1014  GSERIALIZED *g1 = NULL;
1015  GSERIALIZED *g2 = NULL;
1016  double azimuth;
1017  SPHEROID s;
1018  uint32_t type1, type2;
1019 
1020  /* Get our geometry object loaded into memory. */
1021  g1 = PG_GETARG_GSERIALIZED_P(0);
1022  g2 = PG_GETARG_GSERIALIZED_P(1);
1023 
1024  /* Only return for points. */
1025  type1 = gserialized_get_type(g1);
1026  type2 = gserialized_get_type(g2);
1027  if ( type1 != POINTTYPE || type2 != POINTTYPE )
1028  {
1029  elog(ERROR, "ST_Azimuth(geography, geography) is only valid for point inputs");
1030  PG_RETURN_NULL();
1031  }
1032 
1033  lwgeom1 = lwgeom_from_gserialized(g1);
1034  lwgeom2 = lwgeom_from_gserialized(g2);
1035 
1036  /* EMPTY things cannot be used */
1037  if ( lwgeom_is_empty(lwgeom1) || lwgeom_is_empty(lwgeom2) )
1038  {
1039  lwgeom_free(lwgeom1);
1040  lwgeom_free(lwgeom2);
1041  elog(ERROR, "ST_Azimuth(geography, geography) cannot work with empty points");
1042  PG_RETURN_NULL();
1043  }
1044 
1045  /* Initialize spheroid */
1046  spheroid_init_from_srid(fcinfo, gserialized_get_srid(g1), &s);
1047 
1048  /* Calculate the direction */
1049  azimuth = lwgeom_azumith_spheroid(lwgeom_as_lwpoint(lwgeom1), lwgeom_as_lwpoint(lwgeom2), &s);
1050 
1051  /* Clean up */
1052  lwgeom_free(lwgeom1);
1053  lwgeom_free(lwgeom2);
1054 
1055  PG_FREE_IF_COPY(g1, 0);
1056  PG_FREE_IF_COPY(g2, 1);
1057 
1058  /* Return NULL for unknown (same point) azimuth */
1059  if( isnan(azimuth) )
1060  {
1061  PG_RETURN_NULL();
1062  }
1063 
1064  PG_RETURN_FLOAT8(azimuth);
1065 }
1066 
1067 
1068 
1069 /*
1070 ** geography_segmentize(GSERIALIZED *g1, double max_seg_length)
1071 ** returns densified geometry with no segment longer than max
1072 */
1074 Datum geography_segmentize(PG_FUNCTION_ARGS)
1075 {
1076  LWGEOM *lwgeom1 = NULL;
1077  LWGEOM *lwgeom2 = NULL;
1078  GSERIALIZED *g1 = NULL;
1079  GSERIALIZED *g2 = NULL;
1080  double max_seg_length;
1081  uint32_t type1;
1082 
1083  /* Get our geometry object loaded into memory. */
1084  g1 = PG_GETARG_GSERIALIZED_P(0);
1085  type1 = gserialized_get_type(g1);
1086 
1087  /* Convert max_seg_length from metric units to radians */
1088  max_seg_length = PG_GETARG_FLOAT8(1) / WGS84_RADIUS;
1089 
1090  /* We can't densify points or points, reflect them back */
1091  if ( type1 == POINTTYPE || type1 == MULTIPOINTTYPE || gserialized_is_empty(g1) )
1092  PG_RETURN_POINTER(g1);
1093 
1094  /* Deserialize */
1095  lwgeom1 = lwgeom_from_gserialized(g1);
1096 
1097  /* Calculate the densified geometry */
1098  lwgeom2 = lwgeom_segmentize_sphere(lwgeom1, max_seg_length);
1099 
1100  /*
1101  ** Set the geodetic flag so subsequent
1102  ** functions do the right thing.
1103  */
1104  lwgeom_set_geodetic(lwgeom2, true);
1105 
1106  /* Recalculate the boxes after re-setting the geodetic bit */
1107  lwgeom_drop_bbox(lwgeom2);
1108 
1109  /* We are trusting geography_serialize will add a box if needed */
1110  g2 = geography_serialize(lwgeom2);
1111 
1112  /* Clean up */
1113  lwgeom_free(lwgeom1);
1114  lwgeom_free(lwgeom2);
1115  PG_FREE_IF_COPY(g1, 0);
1116 
1117  PG_RETURN_POINTER(g2);
1118 }
1119 
1120 
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:642
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:3012
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:660
GBOX * bbox
Definition: liblwgeom.h:398
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:2152
int geography_distance_cache(FunctionCallInfo fcinfo, const GSERIALIZED *g1, const GSERIALIZED *g2, const SPHEROID *s, double *distance)
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:404
#define WGS84_RADIUS
Definition: liblwgeom.h:131
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:2405
#define POLYGONTYPE
Definition: liblwgeom.h:87
LWPOINT * lwpoint_make2d(int srid, double x, double y)
Definition: lwpoint.c:163
Datum area(PG_FUNCTION_ARGS)
double b
Definition: liblwgeom.h:314
Datum geography_distance_knn(PG_FUNCTION_ARGS)
void lwpoint_free(LWPOINT *pt)
Definition: lwpoint.c:213
void lwgeom_free(LWGEOM *geom)
Definition: lwgeom.c:1099
#define INVMINDIST
Datum geography_dwithin_uncached(PG_FUNCTION_ARGS)
#define MULTIPOINTTYPE
Definition: liblwgeom.h:88
double radius
Definition: liblwgeom.h:318
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)
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:2095
void error_if_srid_mismatch(int srid1, int srid2)
Definition: lwutil.c:371
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:129
double lwgeom_length_spheroid(const LWGEOM *geom, const SPHEROID *s)
Calculate the geodetic length of a lwgeom on the unit sphere.
Definition: lwgeodetic.c:3281
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:635
int32_t srid
Definition: liblwgeom.h:399
int gserialized_is_empty(const GSERIALIZED *g)
Check if a GSERIALIZED is empty without deserializing first.
Definition: g_serialized.c:179
#define LW_FAILURE
Definition: liblwgeom.h:79
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:328
Datum geography_bestsrid(PG_FUNCTION_ARGS)
double zmax
Definition: liblwgeom.h:297
#define LW_FALSE
Definition: liblwgeom.h:77
#define LW_TRUE
Return types for functions with status returns.
Definition: liblwgeom.h:76
Datum geography_point_outside(PG_FUNCTION_ARGS)
#define FP_TOLERANCE
Floating point comparators.
double y
Definition: liblwgeom.h:328
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:146
double lwgeom_area_sphere(const LWGEOM *lwgeom, const SPHEROID *spheroid)
Calculate the geodetic area of a lwgeom on the sphere.
Definition: lwgeodetic.c:2027
void lwgeom_set_geodetic(LWGEOM *geom, int value)
Set the FLAGS geodetic bit on geometry an all sub-geometries and pointlists.
Definition: lwgeom.c:907
uint8_t flags
Definition: liblwgeom.h:291
#define MULTIPOLYGONTYPE
Definition: liblwgeom.h:90
int geography_dwithin_cache(FunctionCallInfo fcinfo, const GSERIALIZED *g1, const GSERIALIZED *g2, const SPHEROID *s, double tolerance, int *dwithin)
double a
Definition: liblwgeom.h:313
#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:1739
#define FALSE
Definition: dbfopen.c:168
Datum geography_expand(PG_FUNCTION_ARGS)
Datum geography_distance_tree(PG_FUNCTION_ARGS)
double zmin
Definition: liblwgeom.h:296
double gbox_angular_height(const GBOX *gbox)
GBOX utility functions to figure out coverage/location on the globe.
Definition: lwgeodetic.c:188
#define POINTTYPE
LWTYPE numbers, used internally by PostGIS.
Definition: liblwgeom.h:85
Datum geography_covers(PG_FUNCTION_ARGS)
Datum geography_azimuth(PG_FUNCTION_ARGS)
uint8_t type
Definition: liblwgeom.h:396
int gbox_pt_outside(const GBOX *gbox, POINT2D *pt_outside)
Calculate a spherical point that falls outside the geocentric gbox.
Definition: lwgeodetic.c:1548
type
Definition: ovdump.py:41
#define FP_EQUALS(A, B)
LWGEOM * lwpoint_as_lwgeom(const LWPOINT *obj)
Definition: lwgeom.c:303
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:2183
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:1346
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:100
#define COLLECTIONTYPE
Definition: liblwgeom.h:91
This library is the generic geometry handling section of PostGIS.
uint8_t flags
Definition: liblwgeom.h:383
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:215