PostGIS 3.7.0dev-r@@SVN_REVISION@@
Loading...
Searching...
No Matches
lwgeodetic.h
Go to the documentation of this file.
1/**********************************************************************
2 *
3 * PostGIS - Spatial Types for PostgreSQL
4 * http://postgis.net
5 *
6 * PostGIS is free software: you can redistribute it and/or modify
7 * it under the terms of the GNU General Public License as published by
8 * the Free Software Foundation, either version 2 of the License, or
9 * (at your option) any later version.
10 *
11 * PostGIS is distributed in the hope that it will be useful,
12 * but WITHOUT ANY WARRANTY; without even the implied warranty of
13 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14 * GNU General Public License for more details.
15 *
16 * You should have received a copy of the GNU General Public License
17 * along with PostGIS. If not, see <http://www.gnu.org/licenses/>.
18 *
19 **********************************************************************
20 *
21 * Copyright 2009 Paul Ramsey <pramsey@cleverelephant.ca>
22 *
23 **********************************************************************/
24
25
26#ifndef _LWGEODETIC_H
27#define _LWGEODETIC_H 1
28
29/* For NAN */
30#ifndef _GNU_SOURCE
31#define _GNU_SOURCE
32#endif
33
34#include <math.h>
35
36#ifndef NAN
37#define NAN 0.0/0.0
38#endif
39
40/* Override tolerance for geodetic */
41#ifdef FP_TOLERANCE
42#undef FP_TOLERANCE
43#define FP_TOLERANCE 5e-14
44#endif
45
46extern int gbox_geocentric_slow;
47
48#define POW2(x) ((x)*(x))
49
53typedef struct
54{
55 double lon;
56 double lat;
58
67
71typedef struct
72{
73 double measure;
74 uint32_t index;
76
80#define deg2rad(d) (M_PI * (d) / 180.0)
81#define rad2deg(r) (180.0 * (r) / M_PI)
82
83
87#define PIR_NO_INTERACT 0x00
88#define PIR_INTERSECTS 0x01
89#define PIR_COLINEAR 0x02
90#define PIR_A_TOUCH_RIGHT 0x04
91#define PIR_A_TOUCH_LEFT 0x08
92#define PIR_B_TOUCH_RIGHT 0x10
93#define PIR_B_TOUCH_LEFT 0x20
94
95
96/*
97* Geodetic calculations
98*/
99void geog2cart(const GEOGRAPHIC_POINT *g, POINT3D *p);
100void cart2geog(const POINT3D *p, GEOGRAPHIC_POINT *g);
102void x_to_z(POINT3D *p);
103void y_to_z(POINT3D *p);
108double z_to_latitude(double z, int top);
109int clairaut_cartesian(const POINT3D *start, const POINT3D *end, GEOGRAPHIC_POINT *g_top, GEOGRAPHIC_POINT *g_bottom);
110int clairaut_geographic(const GEOGRAPHIC_POINT *start, const GEOGRAPHIC_POINT *end, GEOGRAPHIC_POINT *g_top, GEOGRAPHIC_POINT *g_bottom);
111double sphere_distance(const GEOGRAPHIC_POINT *s, const GEOGRAPHIC_POINT *e);
112double sphere_distance_cartesian(const POINT3D *s, const POINT3D *e);
113int sphere_project(const GEOGRAPHIC_POINT *r, double distance, double azimuth, GEOGRAPHIC_POINT *n);
115int edge_calculate_gbox(const POINT3D *A1, const POINT3D *A2, GBOX *gbox);
117uint32_t edge_intersects(const POINT3D *A1, const POINT3D *A2, const POINT3D *B1, const POINT3D *B2);
118double edge_distance_to_point(const GEOGRAPHIC_EDGE *e, const GEOGRAPHIC_POINT *gp, GEOGRAPHIC_POINT *closest);
119double edge_distance_to_edge(const GEOGRAPHIC_EDGE *e1, const GEOGRAPHIC_EDGE *e2, GEOGRAPHIC_POINT *closest1, GEOGRAPHIC_POINT *closest2);
120void geographic_point_init(double lon, double lat, GEOGRAPHIC_POINT *g);
121int ptarray_contains_point_sphere(const POINTARRAY *pa, const POINT2D *pt_outside, const POINT2D *pt_to_test);
122int lwpoly_covers_point2d(const LWPOLY *poly, const POINT2D *pt_to_test);
123int lwpoly_covers_lwpoly(const LWPOLY *lwpoly1, const LWPOLY *lwpoly2);
124int lwpoly_covers_pointarray(const LWPOLY* lwpoly, const POINTARRAY* pta);
125int lwpoly_covers_lwline(const LWPOLY *poly, const LWLINE *line);
126int lwline_covers_lwline(const LWLINE* lwline1, const LWLINE* lwline2);
127int lwline_covers_lwpoint(const LWLINE* lwline, const LWPOINT* lwpoint);
128int lwpoly_intersects_line(const LWPOLY* lwpoly, const POINTARRAY* line);
129int lwpoly_pt_outside(const LWPOLY *poly, POINT2D *pt_outside);
130int ptarray_point_in_ring(const POINTARRAY *pa, const POINT2D *pt_outside, const POINT2D *pt_to_test);
131double latitude_degrees_normalize(double lat);
132double longitude_degrees_normalize(double lon);
133double ptarray_length_spheroid(const POINTARRAY *pa, const SPHEROID *s);
136void point_shift(GEOGRAPHIC_POINT *p, double shift);
137double longitude_radians_normalize(double lon);
138double latitude_radians_normalize(double lat);
139void vector_sum(const POINT3D *a, const POINT3D *b, POINT3D *n);
140void vector_scale(POINT3D *a, double s);
141double vector_angle(const POINT3D* v1, const POINT3D* v2);
142void vector_rotate(const POINT3D* v1, const POINT3D* v2, double angle, POINT3D* n);
143void normalize(POINT3D *p);
144void unit_normal(const POINT3D *P1, const POINT3D *P2, POINT3D *normal);
145double sphere_direction(const GEOGRAPHIC_POINT *s, const GEOGRAPHIC_POINT *e, double d);
146void ll2cart(const POINT2D *g, POINT3D *p);
147int gbox_geocentric_get_gbox_cartesian(const GBOX *gbox_geocentric, GBOX *gbox_planar);
148
149/*
150** Prototypes for spheroid functions.
151*/
152double spheroid_distance(const GEOGRAPHIC_POINT *a, const GEOGRAPHIC_POINT *b, const SPHEROID *spheroid);
153double spheroid_direction(const GEOGRAPHIC_POINT *r, const GEOGRAPHIC_POINT *s, const SPHEROID *spheroid);
154int spheroid_project(const GEOGRAPHIC_POINT *r, const SPHEROID *spheroid, double distance, double azimuth, GEOGRAPHIC_POINT *g);
155
156
157#endif /* _LWGEODETIC_H */
158
159
160
char * s
Definition cu_in_wkt.c:23
char * r
Definition cu_in_wkt.c:24
int clairaut_geographic(const GEOGRAPHIC_POINT *start, const GEOGRAPHIC_POINT *end, GEOGRAPHIC_POINT *g_top, GEOGRAPHIC_POINT *g_bottom)
Computes the pole of the great circle disk which is the intersection of the great circle with the lin...
int gbox_geocentric_get_gbox_cartesian(const GBOX *gbox_geocentric, GBOX *gbox_planar)
int lwpoly_covers_point2d(const LWPOLY *poly, const POINT2D *pt_to_test)
Given a polygon (lon/lat decimal degrees) and point (lon/lat decimal degrees) and a guaranteed outsid...
void normalize(POINT3D *p)
Normalize to a unit vector.
Definition lwgeodetic.c:615
int lwline_covers_lwpoint(const LWLINE *lwline, const LWPOINT *lwpoint)
return LW_TRUE if any of the line segments covers the point
int lwpoly_intersects_line(const LWPOLY *lwpoly, const POINTARRAY *line)
Checks if any edges of lwpoly intersect with the line formed by the pointarray return LW_TRUE if any ...
double longitude_radians_normalize(double lon)
Convert a longitude to the range of -PI,PI.
Definition lwgeodetic.c:50
int clairaut_cartesian(const POINT3D *start, const POINT3D *end, GEOGRAPHIC_POINT *g_top, GEOGRAPHIC_POINT *g_bottom)
Computes the pole of the great circle disk which is the intersection of the great circle with the lin...
void point_shift(GEOGRAPHIC_POINT *p, double shift)
Shift a point around by a number of radians.
Definition lwgeodetic.c:160
int ptarray_point_in_ring(const POINTARRAY *pa, const POINT2D *pt_outside, const POINT2D *pt_to_test)
double sphere_distance_cartesian(const POINT3D *s, const POINT3D *e)
Given two unit vectors, calculate their distance apart in radians.
Definition lwgeodetic.c:919
void vector_rotate(const POINT3D *v1, const POINT3D *v2, double angle, POINT3D *n)
Rotates v1 through an angle (in radians) within the plane defined by v1/v2, returns the rotated vecto...
Definition lwgeodetic.c:573
int lwpoly_covers_pointarray(const LWPOLY *lwpoly, const POINTARRAY *pta)
return LW_TRUE if all points are inside the polygon
void vector_scale(POINT3D *a, double s)
Scale a vector out by a factor.
Definition lwgeodetic.c:487
int ptarray_contains_point_sphere(const POINTARRAY *pa, const POINT2D *pt_outside, const POINT2D *pt_to_test)
This routine returns LW_TRUE if the stabline joining the pt_outside and pt_to_test crosses the ring a...
int edge_point_on_plane(const GEOGRAPHIC_EDGE *e, const GEOGRAPHIC_POINT *p)
Returns true if the point p is on the great circle plane.
Definition lwgeodetic.c:723
double latitude_radians_normalize(double lat)
Convert a latitude to the range of -PI/2,PI/2.
Definition lwgeodetic.c:78
double spheroid_distance(const GEOGRAPHIC_POINT *a, const GEOGRAPHIC_POINT *b, const SPHEROID *spheroid)
Computes the shortest distance along the surface of the spheroid between two points,...
Definition lwspheroid.c:79
void cart2geog(const POINT3D *p, GEOGRAPHIC_POINT *g)
Convert cartesian coordinates on unit sphere to spherical coordinates.
Definition lwgeodetic.c:414
void y_to_z(POINT3D *p)
Definition lwgeodetic.c:658
int lwpoly_covers_lwline(const LWPOLY *poly, const LWLINE *line)
int lwpoly_pt_outside(const LWPOLY *poly, POINT2D *pt_outside)
int sphere_project(const GEOGRAPHIC_POINT *r, double distance, double azimuth, GEOGRAPHIC_POINT *n)
Given a starting location r, a distance and an azimuth to the new point, compute the location of the ...
void ll2cart(const POINT2D *g, POINT3D *p)
Convert lon/lat coordinates to cartesian coordinates on unit sphere.
Definition lwgeodetic.c:423
int gbox_geocentric_slow
For testing geodetic bounding box, we have a magic global variable.
Definition lwgeodetic.c:36
int edge_point_in_cone(const GEOGRAPHIC_EDGE *e, const GEOGRAPHIC_POINT *p)
Returns true if the point p is inside the cone defined by the two ends of the edge e.
Definition lwgeodetic.c:736
double longitude_degrees_normalize(double lon)
Convert a longitude to the range of -180,180.
Definition lwgeodetic.c:106
double z_to_latitude(double z, int top)
Used in great circle to compute the pole of the great circle.
void robust_cross_product(const GEOGRAPHIC_POINT *p, const GEOGRAPHIC_POINT *q, POINT3D *a)
Computes the cross product of two vectors using their lat, lng representations.
Definition lwgeodetic.c:634
int lwpoly_covers_lwpoly(const LWPOLY *lwpoly1, const LWPOLY *lwpoly2)
Given a polygon1 check if all points of polygon2 are inside polygon1 and no intersections of the poly...
void geographic_point_init(double lon, double lat, GEOGRAPHIC_POINT *g)
Initialize a geographic point.
Definition lwgeodetic.c:180
double ptarray_length_spheroid(const POINTARRAY *pa, const SPHEROID *s)
void unit_normal(const POINT3D *P1, const POINT3D *P2, POINT3D *normal)
Calculates the unit normal to two vectors, trying to avoid problems with over-narrow or over-wide cas...
Definition lwgeodetic.c:541
int spheroid_project(const GEOGRAPHIC_POINT *r, const SPHEROID *spheroid, double distance, double azimuth, GEOGRAPHIC_POINT *g)
Given a location, an azimuth and a distance, computes the location of the projected point.
Definition lwspheroid.c:128
double sphere_distance(const GEOGRAPHIC_POINT *s, const GEOGRAPHIC_POINT *e)
Given two points on a unit sphere, calculate their distance apart in radians.
Definition lwgeodetic.c:896
double sphere_direction(const GEOGRAPHIC_POINT *s, const GEOGRAPHIC_POINT *e, double d)
Given two points on a unit sphere, calculate the direction from s to e.
Definition lwgeodetic.c:927
int edge_calculate_gbox_slow(const GEOGRAPHIC_EDGE *e, GBOX *gbox)
int edge_intersection(const GEOGRAPHIC_EDGE *e1, const GEOGRAPHIC_EDGE *e2, GEOGRAPHIC_POINT *g)
Returns true if an intersection can be calculated, and places it in *g.
void vector_sum(const POINT3D *a, const POINT3D *b, POINT3D *n)
Calculate the sum of two vectors.
Definition lwgeodetic.c:465
double spheroid_direction(const GEOGRAPHIC_POINT *r, const GEOGRAPHIC_POINT *s, const SPHEROID *spheroid)
Computes the forward azimuth of the geodesic joining two points on the spheroid, using the inverse ge...
Definition lwspheroid.c:105
uint32_t edge_intersects(const POINT3D *A1, const POINT3D *A2, const POINT3D *B1, const POINT3D *B2)
Returns non-zero if edges A and B interact.
int edge_calculate_gbox(const POINT3D *A1, const POINT3D *A2, GBOX *gbox)
The magic function, given an edge in spherical coordinates, calculate a 3D bounding box that fully co...
double edge_distance_to_point(const GEOGRAPHIC_EDGE *e, const GEOGRAPHIC_POINT *gp, GEOGRAPHIC_POINT *closest)
void geog2cart(const GEOGRAPHIC_POINT *g, POINT3D *p)
Convert spherical coordinates to cartesian coordinates on unit sphere.
Definition lwgeodetic.c:404
int lwline_covers_lwline(const LWLINE *lwline1, const LWLINE *lwline2)
Check if first and last point of line2 are covered by line1 and then each point in between has to be ...
int edge_contains_point(const GEOGRAPHIC_EDGE *e, const GEOGRAPHIC_POINT *p)
Returns true if the point p is on the minor edge defined by the end points of e.
Definition lwgeodetic.c:986
void x_to_z(POINT3D *p)
Definition lwgeodetic.c:651
double edge_distance_to_edge(const GEOGRAPHIC_EDGE *e1, const GEOGRAPHIC_EDGE *e2, GEOGRAPHIC_POINT *closest1, GEOGRAPHIC_POINT *closest2)
Calculate the distance between two edges.
int geographic_point_equals(const GEOGRAPHIC_POINT *g1, const GEOGRAPHIC_POINT *g2)
Definition lwgeodetic.c:170
int crosses_dateline(const GEOGRAPHIC_POINT *s, const GEOGRAPHIC_POINT *e)
Definition lwgeodetic.c:666
double latitude_degrees_normalize(double lat)
Convert a latitude to the range of -90,90.
Definition lwgeodetic.c:133
double vector_angle(const POINT3D *v1, const POINT3D *v2)
Angle between two unit vectors.
Definition lwgeodetic.c:505
int edge_contains_coplanar_point(const GEOGRAPHIC_EDGE *e, const GEOGRAPHIC_POINT *p)
True if the longitude of p is within the range of the longitude of the ends of e.
Definition lwgeodetic.c:783
static double distance(double x1, double y1, double x2, double y2)
Definition lwtree.c:1032
uint32_t index
Definition lwgeodetic.h:74
Holder for sorting points in distance algorithm.
Definition lwgeodetic.h:72
GEOGRAPHIC_POINT start
Definition lwgeodetic.h:64
GEOGRAPHIC_POINT end
Definition lwgeodetic.h:65
Two-point great circle segment from a to b.
Definition lwgeodetic.h:63
Point in spherical coordinates on the world.
Definition lwgeodetic.h:54