1/* contrib/earthdistance/earthdistance--1.1.sql */
2
3-- complain if script is sourced in psql, rather than via CREATE EXTENSION
4\echo Use "CREATE EXTENSION earthdistance" to load this file. \quit
5
6-- earth() returns the radius of the earth in meters. This is the only
7-- place you need to change things for the cube base distance functions
8-- in order to use different units (or a better value for the Earth's radius).
9
10CREATE FUNCTION earth() RETURNS float8
11LANGUAGE SQL IMMUTABLE PARALLEL SAFE
12AS 'SELECT ''6378168''::float8';
13
14-- Astronomers may want to change the earth function so that distances will be
15-- returned in degrees. To do this comment out the above definition and
16-- uncomment the one below. Note that doing this will break the regression
17-- tests.
18--
19-- CREATE FUNCTION earth() RETURNS float8
20-- LANGUAGE SQL IMMUTABLE
21-- AS 'SELECT 180/pi()';
22
23-- Define domain for locations on the surface of the earth using a cube
24-- datatype with constraints. cube provides 3D indexing.
25-- The cube is restricted to be a point, no more than 3 dimensions
26-- (for less than 3 dimensions 0 is assumed for the missing coordinates)
27-- and that the point must be very near the surface of the sphere
28-- centered about the origin with the radius of the earth.
29
30CREATE DOMAIN earth AS cube
31  CONSTRAINT not_point check(cube_is_point(value))
32  CONSTRAINT not_3d check(cube_dim(value) <= 3)
33  CONSTRAINT on_surface check(abs(cube_distance(value, '(0)'::cube) /
34  earth() - '1'::float8) < '10e-7'::float8);
35
36CREATE FUNCTION sec_to_gc(float8)
37RETURNS float8
38LANGUAGE SQL
39IMMUTABLE STRICT
40PARALLEL SAFE
41AS 'SELECT CASE WHEN $1 < 0 THEN 0::float8 WHEN $1/(2*earth()) > 1 THEN pi()*earth() ELSE 2*earth()*asin($1/(2*earth())) END';
42
43CREATE FUNCTION gc_to_sec(float8)
44RETURNS float8
45LANGUAGE SQL
46IMMUTABLE STRICT
47PARALLEL SAFE
48AS 'SELECT CASE WHEN $1 < 0 THEN 0::float8 WHEN $1/earth() > pi() THEN 2*earth() ELSE 2*earth()*sin($1/(2*earth())) END';
49
50CREATE FUNCTION ll_to_earth(float8, float8)
51RETURNS earth
52LANGUAGE SQL
53IMMUTABLE STRICT
54PARALLEL SAFE
55AS 'SELECT cube(cube(cube(earth()*cos(radians($1))*cos(radians($2))),earth()*cos(radians($1))*sin(radians($2))),earth()*sin(radians($1)))::earth';
56
57CREATE FUNCTION latitude(earth)
58RETURNS float8
59LANGUAGE SQL
60IMMUTABLE STRICT
61PARALLEL SAFE
62AS 'SELECT CASE WHEN cube_ll_coord($1, 3)/earth() < -1 THEN -90::float8 WHEN cube_ll_coord($1, 3)/earth() > 1 THEN 90::float8 ELSE degrees(asin(cube_ll_coord($1, 3)/earth())) END';
63
64CREATE FUNCTION longitude(earth)
65RETURNS float8
66LANGUAGE SQL
67IMMUTABLE STRICT
68PARALLEL SAFE
69AS 'SELECT degrees(atan2(cube_ll_coord($1, 2), cube_ll_coord($1, 1)))';
70
71CREATE FUNCTION earth_distance(earth, earth)
72RETURNS float8
73LANGUAGE SQL
74IMMUTABLE STRICT
75PARALLEL SAFE
76AS 'SELECT sec_to_gc(cube_distance($1, $2))';
77
78CREATE FUNCTION earth_box(earth, float8)
79RETURNS cube
80LANGUAGE SQL
81IMMUTABLE STRICT
82PARALLEL SAFE
83AS 'SELECT cube_enlarge($1, gc_to_sec($2), 3)';
84
85--------------- geo_distance
86
87CREATE FUNCTION geo_distance (point, point)
88RETURNS float8
89LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE AS 'MODULE_PATHNAME';
90
91--------------- geo_distance as operator <@>
92
93CREATE OPERATOR <@> (
94  LEFTARG = point,
95  RIGHTARG = point,
96  PROCEDURE = geo_distance,
97  COMMUTATOR = <@>
98);
99