PostGIS  2.3.7dev-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 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) / WGS84_RADIUS;
476 
477  /* Try the expansion */
478  g_out = gserialized_expand(g, distance);
479 
480  /* If the expansion fails, the return our input */
481  if ( g_out == NULL )
482  {
483  PG_RETURN_POINTER(g);
484  }
485 
486  if ( g_out != g )
487  {
488  pfree(g);
489  }
490 
491  PG_RETURN_POINTER(g_out);
492 }
493 
494 /*
495 ** geography_area(GSERIALIZED *g)
496 ** returns double area in meters square
497 */
499 Datum geography_area(PG_FUNCTION_ARGS)
500 {
501  LWGEOM *lwgeom = NULL;
502  GSERIALIZED *g = NULL;
503  GBOX gbox;
504  double area;
505  bool use_spheroid = LW_TRUE;
506  SPHEROID s;
507 
508  /* Get our geometry object loaded into memory. */
509  g = PG_GETARG_GSERIALIZED_P(0);
510 
511  /* Read our calculation type */
512  use_spheroid = PG_GETARG_BOOL(1);
513 
514  /* Initialize spheroid */
515  spheroid_init_from_srid(fcinfo, gserialized_get_srid(g), &s);
516 
517  lwgeom = lwgeom_from_gserialized(g);
518 
519  /* EMPTY things have no area */
520  if ( lwgeom_is_empty(lwgeom) )
521  {
522  lwgeom_free(lwgeom);
523  PG_RETURN_FLOAT8(0.0);
524  }
525 
526  if ( lwgeom->bbox )
527  gbox = *(lwgeom->bbox);
528  else
529  lwgeom_calculate_gbox_geodetic(lwgeom, &gbox);
530 
531 #if ! PROJ_GEODESIC
532  /* Test for cases that are currently not handled by spheroid code */
533  if ( use_spheroid )
534  {
535  /* We can't circle the poles right now */
536  if ( FP_GTEQ(gbox.zmax,1.0) || FP_LTEQ(gbox.zmin,-1.0) )
537  use_spheroid = LW_FALSE;
538  /* We can't cross the equator right now */
539  if ( gbox.zmax > 0.0 && gbox.zmin < 0.0 )
540  use_spheroid = LW_FALSE;
541  }
542 #endif /* if ! PROJ_GEODESIC */
543 
544  /* User requests spherical calculation, turn our spheroid into a sphere */
545  if ( ! use_spheroid )
546  s.a = s.b = s.radius;
547 
548  /* Calculate the area */
549  if ( use_spheroid )
550  area = lwgeom_area_spheroid(lwgeom, &s);
551  else
552  area = lwgeom_area_sphere(lwgeom, &s);
553 
554  /* Clean up */
555  lwgeom_free(lwgeom);
556  PG_FREE_IF_COPY(g, 0);
557 
558  /* Something went wrong... */
559  if ( area < 0.0 )
560  {
561  elog(ERROR, "lwgeom_area_spher(oid) returned area < 0.0");
562  PG_RETURN_NULL();
563  }
564 
565  PG_RETURN_FLOAT8(area);
566 }
567 
568 /*
569 ** geography_perimeter(GSERIALIZED *g)
570 ** returns double perimeter in meters for area features
571 */
573 Datum geography_perimeter(PG_FUNCTION_ARGS)
574 {
575  LWGEOM *lwgeom = NULL;
576  GSERIALIZED *g = NULL;
577  double length;
578  bool use_spheroid = LW_TRUE;
579  SPHEROID s;
580  int type;
581 
582  /* Get our geometry object loaded into memory. */
583  g = PG_GETARG_GSERIALIZED_P(0);
584 
585  /* Only return for area features. */
586  type = gserialized_get_type(g);
587  if ( ! (type == POLYGONTYPE || type == MULTIPOLYGONTYPE || type == COLLECTIONTYPE) )
588  {
589  PG_RETURN_FLOAT8(0.0);
590  }
591 
592  lwgeom = lwgeom_from_gserialized(g);
593 
594  /* EMPTY things have no perimeter */
595  if ( lwgeom_is_empty(lwgeom) )
596  {
597  lwgeom_free(lwgeom);
598  PG_RETURN_FLOAT8(0.0);
599  }
600 
601  /* Read our calculation type */
602  use_spheroid = PG_GETARG_BOOL(1);
603 
604  /* Initialize spheroid */
605  spheroid_init_from_srid(fcinfo, gserialized_get_srid(g), &s);
606 
607  /* User requests spherical calculation, turn our spheroid into a sphere */
608  if ( ! use_spheroid )
609  s.a = s.b = s.radius;
610 
611  /* Calculate the length */
612  length = lwgeom_length_spheroid(lwgeom, &s);
613 
614  /* Something went wrong... */
615  if ( length < 0.0 )
616  {
617  elog(ERROR, "lwgeom_length_spheroid returned length < 0.0");
618  PG_RETURN_NULL();
619  }
620 
621  /* Clean up, but not all the way to the point arrays */
622  lwgeom_free(lwgeom);
623 
624  PG_FREE_IF_COPY(g, 0);
625  PG_RETURN_FLOAT8(length);
626 }
627 
628 /*
629 ** geography_length(GSERIALIZED *g)
630 ** returns double length in meters
631 */
633 Datum geography_length(PG_FUNCTION_ARGS)
634 {
635  LWGEOM *lwgeom = NULL;
636  GSERIALIZED *g = NULL;
637  double length;
638  bool use_spheroid = LW_TRUE;
639  SPHEROID s;
640 
641  /* Get our geometry object loaded into memory. */
642  g = PG_GETARG_GSERIALIZED_P(0);
643  lwgeom = lwgeom_from_gserialized(g);
644 
645  /* EMPTY things have no length */
646  if ( lwgeom_is_empty(lwgeom) || lwgeom->type == POLYGONTYPE || lwgeom->type == MULTIPOLYGONTYPE )
647  {
648  lwgeom_free(lwgeom);
649  PG_RETURN_FLOAT8(0.0);
650  }
651 
652  /* Read our calculation type */
653  use_spheroid = PG_GETARG_BOOL(1);
654 
655  /* Initialize spheroid */
656  spheroid_init_from_srid(fcinfo, gserialized_get_srid(g), &s);
657 
658  /* User requests spherical calculation, turn our spheroid into a sphere */
659  if ( ! use_spheroid )
660  s.a = s.b = s.radius;
661 
662  /* Calculate the length */
663  length = lwgeom_length_spheroid(lwgeom, &s);
664 
665  /* Something went wrong... */
666  if ( length < 0.0 )
667  {
668  elog(ERROR, "lwgeom_length_spheroid returned length < 0.0");
669  PG_RETURN_NULL();
670  }
671 
672  /* Clean up */
673  lwgeom_free(lwgeom);
674 
675  PG_FREE_IF_COPY(g, 0);
676  PG_RETURN_FLOAT8(length);
677 }
678 
679 /*
680 ** geography_point_outside(GSERIALIZED *g)
681 ** returns point outside the object
682 */
684 Datum geography_point_outside(PG_FUNCTION_ARGS)
685 {
686  GBOX gbox;
687  GSERIALIZED *g = NULL;
688  GSERIALIZED *g_out = NULL;
689  size_t g_out_size;
690  LWGEOM *lwpoint = NULL;
691  POINT2D pt;
692 
693  /* Get our geometry object loaded into memory. */
694  g = PG_GETARG_GSERIALIZED_P(0);
695 
696  /* We need the bounding box to get an outside point for area algorithm */
697  if ( gserialized_get_gbox_p(g, &gbox) == LW_FAILURE )
698  {
699  POSTGIS_DEBUG(4,"gserialized_get_gbox_p returned LW_FAILURE");
700  elog(ERROR, "Error in gserialized_get_gbox_p calculation.");
701  PG_RETURN_NULL();
702  }
703  POSTGIS_DEBUGF(4, "got gbox %s", gbox_to_string(&gbox));
704 
705  /* Get an exterior point, based on this gbox */
706  gbox_pt_outside(&gbox, &pt);
707 
708  lwpoint = (LWGEOM*) lwpoint_make2d(4326, pt.x, pt.y);
709  /* TODO: Investigate where this is used, this was probably not
710  * returning a geography object before. How did this miss checking
711  */
712  lwgeom_set_geodetic(lwpoint, true);
713  g_out = gserialized_from_lwgeom(lwpoint, &g_out_size);
714  SET_VARSIZE(g_out, g_out_size);
715 
716  PG_FREE_IF_COPY(g, 0);
717  PG_RETURN_POINTER(g_out);
718 
719 }
720 
721 /*
722 ** geography_covers(GSERIALIZED *g, GSERIALIZED *g) returns boolean
723 ** Only works for (multi)points and (multi)polygons currently.
724 ** Attempts a simple point-in-polygon test on the polygon and point.
725 ** Current algorithm does not distinguish between points on edge
726 ** and points within.
727 */
729 Datum geography_covers(PG_FUNCTION_ARGS)
730 {
731  LWGEOM *lwgeom1 = NULL;
732  LWGEOM *lwgeom2 = NULL;
733  GSERIALIZED *g1 = NULL;
734  GSERIALIZED *g2 = NULL;
735  int type1, type2;
736  int result = LW_FALSE;
737 
738  /* Get our geometry objects loaded into memory. */
739  g1 = PG_GETARG_GSERIALIZED_P(0);
740  g2 = PG_GETARG_GSERIALIZED_P(1);
741 
742  type1 = gserialized_get_type(g1);
743  type2 = gserialized_get_type(g2);
744 
745  /* Right now we only handle points and polygons */
746  if ( ! ( (type1 == POLYGONTYPE || type1 == MULTIPOLYGONTYPE || type1 == COLLECTIONTYPE) &&
747  (type2 == POINTTYPE || type2 == MULTIPOINTTYPE || type2 == COLLECTIONTYPE) ) )
748  {
749  elog(ERROR, "geography_covers: only POLYGON and POINT types are currently supported");
750  PG_RETURN_NULL();
751  }
752 
753  /* Construct our working geometries */
754  lwgeom1 = lwgeom_from_gserialized(g1);
755  lwgeom2 = lwgeom_from_gserialized(g2);
756 
757  error_if_srid_mismatch(lwgeom1->srid, lwgeom2->srid);
758 
759  /* EMPTY never intersects with another geometry */
760  if ( lwgeom_is_empty(lwgeom1) || lwgeom_is_empty(lwgeom2) )
761  {
762  lwgeom_free(lwgeom1);
763  lwgeom_free(lwgeom2);
764  PG_FREE_IF_COPY(g1, 0);
765  PG_FREE_IF_COPY(g2, 1);
766  PG_RETURN_BOOL(false);
767  }
768 
769  /* Calculate answer */
770  result = lwgeom_covers_lwgeom_sphere(lwgeom1, lwgeom2);
771 
772  /* Clean up */
773  lwgeom_free(lwgeom1);
774  lwgeom_free(lwgeom2);
775  PG_FREE_IF_COPY(g1, 0);
776  PG_FREE_IF_COPY(g2, 1);
777 
778  PG_RETURN_BOOL(result);
779 }
780 
781 /*
782 ** geography_bestsrid(GSERIALIZED *g, GSERIALIZED *g) returns int
783 ** Utility function. Returns negative SRID numbers that match to the
784 ** numbers handled in code by the transform(lwgeom, srid) function.
785 ** UTM, polar stereographic and mercator as fallback. To be used
786 ** in wrapping existing geometry functions in SQL to provide access
787 ** to them in the geography module.
788 */
790 Datum geography_bestsrid(PG_FUNCTION_ARGS)
791 {
792  GBOX gbox, gbox1, gbox2;
793  GSERIALIZED *g1 = NULL;
794  GSERIALIZED *g2 = NULL;
795  int empty1 = LW_FALSE;
796  int empty2 = LW_FALSE;
797  double xwidth, ywidth;
798  POINT2D center;
799 
800  Datum d1 = PG_GETARG_DATUM(0);
801  Datum d2 = PG_GETARG_DATUM(1);
802 
803  /* Get our geometry objects loaded into memory. */
804  g1 = (GSERIALIZED*)PG_DETOAST_DATUM(d1);
805  /* Synchronize our box types */
806  gbox1.flags = g1->flags;
807  /* Calculate if the geometry is empty. */
808  empty1 = gserialized_is_empty(g1);
809  /* Calculate a geocentric bounds for the objects */
810  if ( ! empty1 && gserialized_get_gbox_p(g1, &gbox1) == LW_FAILURE )
811  elog(ERROR, "Error in geography_bestsrid calling gserialized_get_gbox_p(g1, &gbox1)");
812 
813  POSTGIS_DEBUGF(4, "calculated gbox = %s", gbox_to_string(&gbox1));
814 
815  /* If we have a unique second argument, fill in all the necessary variables. */
816  if ( d1 != d2 )
817  {
818  g2 = (GSERIALIZED*)PG_DETOAST_DATUM(d2);
819  gbox2.flags = g2->flags;
820  empty2 = gserialized_is_empty(g2);
821  if ( ! empty2 && gserialized_get_gbox_p(g2, &gbox2) == LW_FAILURE )
822  elog(ERROR, "Error in geography_bestsrid calling gserialized_get_gbox_p(g2, &gbox2)");
823  }
824  /*
825  ** If no unique second argument, copying the box for the first
826  ** argument will give us the right answer for all subsequent tests.
827  */
828  else
829  {
830  gbox = gbox2 = gbox1;
831  }
832 
833  /* Both empty? We don't have an answer. */
834  if ( empty1 && empty2 )
835  PG_RETURN_NULL();
836 
837  /* One empty? We can use the other argument values as infill. Otherwise merge the boxen */
838  if ( empty1 )
839  gbox = gbox2;
840  else if ( empty2 )
841  gbox = gbox1;
842  else
843  gbox_union(&gbox1, &gbox2, &gbox);
844 
845  gbox_centroid(&gbox, &center);
846 
847  /* Width and height in degrees */
848  xwidth = 180.0 * gbox_angular_width(&gbox) / M_PI;
849  ywidth = 180.0 * gbox_angular_height(&gbox) / M_PI;
850 
851  POSTGIS_DEBUGF(2, "xwidth %g", xwidth);
852  POSTGIS_DEBUGF(2, "ywidth %g", ywidth);
853  POSTGIS_DEBUGF(2, "center POINT(%g %g)", center.x, center.y);
854 
855  /* Are these data arctic? Lambert Azimuthal Equal Area North. */
856  if ( center.y > 70.0 && ywidth < 45.0 )
857  {
858  PG_RETURN_INT32(SRID_NORTH_LAMBERT);
859  }
860 
861  /* Are these data antarctic? Lambert Azimuthal Equal Area South. */
862  if ( center.y < -70.0 && ywidth < 45.0 )
863  {
864  PG_RETURN_INT32(SRID_SOUTH_LAMBERT);
865  }
866 
867  /*
868  ** Can we fit these data into one UTM zone?
869  ** We will assume we can push things as
870  ** far as a half zone past a zone boundary.
871  ** Note we have no handling for the date line in here.
872  */
873  if ( xwidth < 6.0 )
874  {
875  int zone = floor((center.x + 180.0) / 6.0);
876 
877  if ( zone > 59 ) zone = 59;
878 
879  /* Are these data below the equator? UTM South. */
880  if ( center.y < 0.0 )
881  {
882  PG_RETURN_INT32( SRID_SOUTH_UTM_START + zone );
883  }
884  /* Are these data above the equator? UTM North. */
885  else
886  {
887  PG_RETURN_INT32( SRID_NORTH_UTM_START + zone );
888  }
889  }
890 
891  /*
892  ** Can we fit into a custom LAEA area? (30 degrees high, variable width)
893  ** We will allow overlap into adjoining areas, but use a slightly narrower test (25) to try
894  ** and minimize the worst case.
895  ** Again, we are hoping the dateline doesn't trip us up much
896  */
897  if ( ywidth < 25.0 )
898  {
899  int xzone = -1;
900  int yzone = 3 + floor(center.y / 30.0); /* (range of 0-5) */
901 
902  /* Equatorial band, 12 zones, 30 degrees wide */
903  if ( (yzone == 2 || yzone == 3) && xwidth < 30.0 )
904  {
905  xzone = 6 + floor(center.x / 30.0);
906  }
907  /* Temperate band, 8 zones, 45 degrees wide */
908  else if ( (yzone == 1 || yzone == 4) && xwidth < 45.0 )
909  {
910  xzone = 4 + floor(center.x / 45.0);
911  }
912  /* Arctic band, 4 zones, 90 degrees wide */
913  else if ( (yzone == 0 || yzone == 5) && xwidth < 90.0 )
914  {
915  xzone = 2 + floor(center.x / 90.0);
916  }
917 
918  /* Did we fit into an appropriate xzone? */
919  if ( xzone != -1 )
920  {
921  PG_RETURN_INT32(SRID_LAEA_START + 20 * yzone + xzone);
922  }
923  }
924 
925  /*
926  ** Running out of options... fall-back to Mercator
927  ** and hope for the best.
928  */
929  PG_RETURN_INT32(SRID_WORLD_MERCATOR);
930 
931 }
932 
933 /*
934 ** geography_project(GSERIALIZED *g, distance, azimuth)
935 ** returns point of projection given start point,
936 ** azimuth in radians (bearing) and distance in meters
937 */
939 Datum geography_project(PG_FUNCTION_ARGS)
940 {
941  LWGEOM *lwgeom = NULL;
942  LWPOINT *lwp_projected;
943  GSERIALIZED *g = NULL;
944  GSERIALIZED *g_out = NULL;
945  double azimuth = 0.0;
946  double distance;
947  SPHEROID s;
948  uint32_t type;
949 
950  /* Return NULL on NULL distance or geography */
951  if ( PG_NARGS() < 2 || PG_ARGISNULL(0) || PG_ARGISNULL(1) )
952  PG_RETURN_NULL();
953 
954  /* Get our geometry object loaded into memory. */
955  g = PG_GETARG_GSERIALIZED_P(0);
956 
957  /* Only return for points. */
958  type = gserialized_get_type(g);
959  if ( type != POINTTYPE )
960  {
961  elog(ERROR, "ST_Project(geography) is only valid for point inputs");
962  PG_RETURN_NULL();
963  }
964 
965  distance = PG_GETARG_FLOAT8(1); /* Distance in Meters */
966  lwgeom = lwgeom_from_gserialized(g);
967 
968  /* EMPTY things cannot be projected from */
969  if ( lwgeom_is_empty(lwgeom) )
970  {
971  lwgeom_free(lwgeom);
972  elog(ERROR, "ST_Project(geography) cannot project from an empty start point");
973  PG_RETURN_NULL();
974  }
975 
976  if ( PG_NARGS() > 2 && ! PG_ARGISNULL(2) )
977  azimuth = PG_GETARG_FLOAT8(2); /* Azimuth in Radians */
978 
979  /* Initialize spheroid */
980  spheroid_init_from_srid(fcinfo, gserialized_get_srid(g), &s);
981 
982  /* Handle the zero distance case */
983  if( FP_EQUALS(distance, 0.0) )
984  {
985  PG_RETURN_POINTER(g);
986  }
987 
988  /* Calculate the length */
989  lwp_projected = lwgeom_project_spheroid(lwgeom_as_lwpoint(lwgeom), &s, distance, azimuth);
990 
991  /* Something went wrong... */
992  if ( lwp_projected == NULL )
993  {
994  elog(ERROR, "lwgeom_project_spheroid returned null");
995  PG_RETURN_NULL();
996  }
997 
998  /* Clean up, but not all the way to the point arrays */
999  lwgeom_free(lwgeom);
1000  g_out = geography_serialize(lwpoint_as_lwgeom(lwp_projected));
1001  lwpoint_free(lwp_projected);
1002 
1003  PG_FREE_IF_COPY(g, 0);
1004  PG_RETURN_POINTER(g_out);
1005 }
1006 
1007 
1008 /*
1009 ** geography_azimuth(GSERIALIZED *g1, GSERIALIZED *g2)
1010 ** returns direction between points (north = 0)
1011 ** azimuth (bearing) and distance
1012 */
1014 Datum geography_azimuth(PG_FUNCTION_ARGS)
1015 {
1016  LWGEOM *lwgeom1 = NULL;
1017  LWGEOM *lwgeom2 = NULL;
1018  GSERIALIZED *g1 = NULL;
1019  GSERIALIZED *g2 = NULL;
1020  double azimuth;
1021  SPHEROID s;
1022  uint32_t type1, type2;
1023 
1024  /* Get our geometry object loaded into memory. */
1025  g1 = PG_GETARG_GSERIALIZED_P(0);
1026  g2 = PG_GETARG_GSERIALIZED_P(1);
1027 
1028  /* Only return for points. */
1029  type1 = gserialized_get_type(g1);
1030  type2 = gserialized_get_type(g2);
1031  if ( type1 != POINTTYPE || type2 != POINTTYPE )
1032  {
1033  elog(ERROR, "ST_Azimuth(geography, geography) is only valid for point inputs");
1034  PG_RETURN_NULL();
1035  }
1036 
1037  lwgeom1 = lwgeom_from_gserialized(g1);
1038  lwgeom2 = lwgeom_from_gserialized(g2);
1039 
1040  /* EMPTY things cannot be used */
1041  if ( lwgeom_is_empty(lwgeom1) || lwgeom_is_empty(lwgeom2) )
1042  {
1043  lwgeom_free(lwgeom1);
1044  lwgeom_free(lwgeom2);
1045  elog(ERROR, "ST_Azimuth(geography, geography) cannot work with empty points");
1046  PG_RETURN_NULL();
1047  }
1048 
1049  /* Initialize spheroid */
1050  spheroid_init_from_srid(fcinfo, gserialized_get_srid(g1), &s);
1051 
1052  /* Calculate the direction */
1053  azimuth = lwgeom_azumith_spheroid(lwgeom_as_lwpoint(lwgeom1), lwgeom_as_lwpoint(lwgeom2), &s);
1054 
1055  /* Clean up */
1056  lwgeom_free(lwgeom1);
1057  lwgeom_free(lwgeom2);
1058 
1059  PG_FREE_IF_COPY(g1, 0);
1060  PG_FREE_IF_COPY(g2, 1);
1061 
1062  /* Return NULL for unknown (same point) azimuth */
1063  if( isnan(azimuth) )
1064  {
1065  PG_RETURN_NULL();
1066  }
1067 
1068  PG_RETURN_FLOAT8(azimuth);
1069 }
1070 
1071 
1072 
1073 /*
1074 ** geography_segmentize(GSERIALIZED *g1, double max_seg_length)
1075 ** returns densified geometry with no segment longer than max
1076 */
1078 Datum geography_segmentize(PG_FUNCTION_ARGS)
1079 {
1080  LWGEOM *lwgeom1 = NULL;
1081  LWGEOM *lwgeom2 = NULL;
1082  GSERIALIZED *g1 = NULL;
1083  GSERIALIZED *g2 = NULL;
1084  double max_seg_length;
1085  uint32_t type1;
1086 
1087  /* Get our geometry object loaded into memory. */
1088  g1 = PG_GETARG_GSERIALIZED_P(0);
1089  type1 = gserialized_get_type(g1);
1090 
1091  /* Convert max_seg_length from metric units to radians */
1092  max_seg_length = PG_GETARG_FLOAT8(1) / WGS84_RADIUS;
1093 
1094  /* We can't densify points or points, reflect them back */
1095  if ( type1 == POINTTYPE || type1 == MULTIPOINTTYPE || gserialized_is_empty(g1) )
1096  PG_RETURN_POINTER(g1);
1097 
1098  /* Deserialize */
1099  lwgeom1 = lwgeom_from_gserialized(g1);
1100 
1101  /* Calculate the densified geometry */
1102  lwgeom2 = lwgeom_segmentize_sphere(lwgeom1, max_seg_length);
1103 
1104  /*
1105  ** Set the geodetic flag so subsequent
1106  ** functions do the right thing.
1107  */
1108  lwgeom_set_geodetic(lwgeom2, true);
1109 
1110  /* Recalculate the boxes after re-setting the geodetic bit */
1111  lwgeom_drop_bbox(lwgeom2);
1112 
1113  /* We are trusting geography_serialize will add a box if needed */
1114  g2 = geography_serialize(lwgeom2);
1115 
1116  /* Clean up */
1117  lwgeom_free(lwgeom1);
1118  lwgeom_free(lwgeom2);
1119  PG_FREE_IF_COPY(g1, 0);
1120 
1121  PG_RETURN_POINTER(g2);
1122 }
1123 
1124 
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:398
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:69
int lwgeom_calculate_gbox_geodetic(const LWGEOM *geom, GBOX *gbox)
Calculate the geodetic bounding box for an LWGEOM.
Definition: lwgeodetic.c:2626
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:624
GBOX * bbox
Definition: liblwgeom.h:397
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:2061
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: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:2314
#define POLYGONTYPE
Definition: liblwgeom.h:86
LWPOINT * lwpoint_make2d(int srid, double x, double y)
Definition: lwpoint.c:145
Datum area(PG_FUNCTION_ARGS)
double b
Definition: liblwgeom.h:313
Datum geography_distance_knn(PG_FUNCTION_ARGS)
void lwpoint_free(LWPOINT *pt)
Definition: lwpoint.c:195
void lwgeom_free(LWGEOM *geom)
Definition: lwgeom.c:1063
#define INVMINDIST
Datum geography_dwithin_uncached(PG_FUNCTION_ARGS)
#define MULTIPOINTTYPE
Definition: liblwgeom.h:87
double radius
Definition: liblwgeom.h:317
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:2006
void error_if_srid_mismatch(int srid1, int srid2)
Definition: lwutil.c:369
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:93
double lwgeom_length_spheroid(const LWGEOM *geom, const SPHEROID *s)
Calculate the geodetic length of a lwgeom on the unit sphere.
Definition: lwgeodetic.c:2895
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:599
int32_t srid
Definition: liblwgeom.h:398
int gserialized_is_empty(const GSERIALIZED *g)
Check if a GSERIALIZED is empty without deserializing first.
Definition: g_serialized.c:153
#define LW_FAILURE
Definition: liblwgeom.h:78
GSERIALIZED * gserialized_expand(GSERIALIZED *g, double distance)
Return a GSERIALIZED with an expanded bounding box.
double x
Definition: liblwgeom.h:327
Datum geography_bestsrid(PG_FUNCTION_ARGS)
double zmax
Definition: liblwgeom.h:296
#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:327
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:1938
void lwgeom_set_geodetic(LWGEOM *geom, int value)
Set the FLAGS geodetic bit on geometry an all sub-geometries and pointlists.
Definition: lwgeom.c:871
uint8_t flags
Definition: liblwgeom.h:290
#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:1461
double a
Definition: liblwgeom.h:312
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:1650
#define FALSE
Definition: dbfopen.c:168
Datum geography_expand(PG_FUNCTION_ARGS)
Datum geography_distance_tree(PG_FUNCTION_ARGS)
double zmin
Definition: liblwgeom.h:295
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:395
#define FP_EQUALS(A, B)
LWGEOM * lwpoint_as_lwgeom(const LWPOINT *obj)
Definition: lwgeom.c:267
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:2092
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:1310
GSERIALIZED * gserialized_from_lwgeom(LWGEOM *geom, size_t *size)
Allocate a new GSERIALIZED from an LWGEOM.
Definition: g_serialized.c:933
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:83
#define COLLECTIONTYPE
Definition: liblwgeom.h:90
This library is the generic geometry handling section of PostGIS.
uint8_t flags
Definition: liblwgeom.h:382
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