Lines 29
##### Keywords
gis (1) gps distance (1) PostgreSQL (7)
##### Permissions
Owner: Stou S.
Viewable by Everyone
Editable by All Siafoo Users

# PL/pgSQL Earth (GPS) distance formula 0

 In Brief A function for calculating the distance (in Km) between two GPS coordinates, along with some support functions like Heaver Sin/Arcsin and Earth radius calculation.... more
 Language SQL
# 's
` 1-- Heaver Sine 2CREATE OR REPLACE FUNCTION heaver_sine(FLOAT) RETURNS FLOAT AS \$\$ 3    SELECT power(sin(\$1/2.0), 2); 4\$\$ LANGUAGE SQL; 5 6-- Heaver Arcsine 7CREATE OR REPLACE FUNCTION arc_heaver_sin(FLOAT) RETURNS FLOAT AS \$\$ 8    SELECT 2.0*asin(sqrt(\$1)); 9\$\$ LANGUAGE SQL;1011-- Provides the Earth radius 12CREATE OR REPLACE FUNCTION earth_radius(lat FLOAT) RETURNS float AS \$\$13    DECLARE14        rad_equatorial FLOAT := 6378.137;15        rad_polar FLOAT := 6356.752;16    BEGIN17    RETURN sqrt(power(rad_equatorial, 2) * power(cos(lat), 2) + power(rad_polar, 2)* power(sin(lat), 2)/rad_equatorial*power(cos(lat), 2) + rad_polar*power(sin(lat), 2));18    END;19\$\$ LANGUAGE plpgsql;2021-- Give two sets of coordinates in degree form22CREATE OR REPLACE FUNCTION distance_in_km(lat1 float, lon1 float, lat2 float, long2 float) 23    RETURNS float AS \$\$24DECLARE25    lat1_r FLOAT := radians(lat1);26    lon1_r FLOAT := radians(lon1);27    lat2_r FLOAT := radians(lat2);28    lon2_r FLOAT := radians(lon2);29BEGIN3031    RETURN earth_radius((lat1_r + lat2_r) / 2.0) * arc_heaver_sin(heaver_sin(abs(lat1_r-lat2_r)) + cos(lat1_r)*cos(lat2_r)*heaver_sin(abs(lon1_r - lon2_r)));32    END;33\$\$ LANGUAGE plpgsql;`

A function for calculating the distance (in Km) between two GPS coordinates, along with some support functions like Heaver Sin/Arcsin and Earth radius calculation.

You should check out the PostgreSQL Earth Distance contrib module http://www.postgresql.org/docs/8.3/interactive/earthdistance.html.

This is based on A couple of GIS functions in Python.

over 12 years ago (06 May 2008 at 01:23 PM) by David Isaacson
Awesome.
over 9 years ago (08 Dec 2010 at 08:41 AM) by bricklen
Here is a slightly different version making use of Perl's Math::Trig module, and gives roughly the same results as the contrib Earthdistance & Cube version.

Note: You need to have plperlu installed as a language in your postgresql db.

-- this example using the contrib module "Earthdistance" yields the same results as the Math::Trig value below
select ROUND( ('(33.0, 97.1)'::point <@> '(24.0, 88.6)')::NUMERIC * 1.609344 ) as km;

create or replace function distance_in_km (numeric,numeric,numeric,numeric) returns numeric as
\$body\$
use strict;
# from http://perldoc.perl.org/Math/Trig.html
my \$lat1 = shift(@_);
my \$lon1 = shift(@_);
my \$lat2 = shift(@_);
my \$lon2 = shift(@_);

#elog WARNING, "\$lat1, \$lon1, \$lat2, \$lon2\n";

# Notice the 90 - latitude: phi zero is at the North Pole.