[postgis-users] deriving distance between lat/long points.

[Brent Wood]

Hi y'all...

I am often being asked to work out the distance between two points (lat/long
WGS84), for genetic cline or mark recapture type studies. 

I figured to make my life easier by writing a shell script to invoke PostGIS to
do this. But it seems not, at least so far  :-)


The first problem is that I'm stuck (for now) with PostGIS v0.8, but this
should be able to do what I want, and proj version 4.4.7

The second is probably my failure to understand projections properly.

I looked for global equidistant projections, which seemed a sensible option and
found a few:
World_Azimuthal_Equidistant (54032)
World_Plate_Carree          (54001)
World_Equidistant_Conic     (54027)


So, my questions:

1. Is it feasible to have a simple script like this to return the distance?
2. Which is the "best" projection to use?
3. Which is the "best" Postgis function to use?
4. Any other advice/suggestions?

Thanks,

  Brent Wood





The results for 3 iterations, one for each projection, were:
(script called arc_dist with 2 coord pairs on the command line)

I have not been able to get distance_spheroid or length_spheroid to work, &
distance_sphere doesn't exist in postgis.sql. I have no say on the versions to
be installed, and can but hope for an upgrade someday.


+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
[woodb at otter ~]$ arc_dist -45 160 -45 180

Parameters:
Line:        GeometryFromText('LINESTRING(160 -45,180 -45)', 4326)
Start point: GeometryFromText('POINT(160 -45)', 4326)
End point:   GeometryFromText('POINT(180 -45)', 4326)
Spheroid = SPHEROID["WGS_1984",6378137,298.257223563]
SRID:      54001

distance(point,point)
ERROR:  tranform: couldnt parse proj4 output string
CONTEXT:  PL/pgSQL function "transform" line 2 at return
distance_spheroid(point, point, spheroid)
NaN
[woodb at otter ~]$ arc_dist -45 160 -45 180

Parameters:
Line:        GeometryFromText('LINESTRING(160 -45,180 -45)', 4326)
Start point: GeometryFromText('POINT(160 -45)', 4326)
End point:   GeometryFromText('POINT(180 -45)', 4326)
Spheroid = SPHEROID["WGS_1984",6378137,298.257223563]
SRID:      54027

distance(point,point)
4616.85
distance_spheroid(point, point, spheroid)
NaN
[woodb at otter ~]$ arc_dist -45 160 -45 180

Parameters:
Line:        GeometryFromText('LINESTRING(160 -45,180 -45)', 4326)
Start point: GeometryFromText('POINT(160 -45)', 4326)
End point:   GeometryFromText('POINT(180 -45)', 4326)
Spheroid = SPHEROID["WGS_1984",6378137,298.257223563]
SRID:      54032

distance(point,point)
1747.37
distance_spheroid(point, point, spheroid)
NaN
[woodb at otter ~]$

+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++

The postgis lines from the script are:

echo "distance(point,point)"
psql $DB -U woodb -q -A -t -c "select (distance(transform($START,$SRID ), 
                        transform($END,$SRID ))/1000)::decimal(9,2) as Km;"    
                              

echo "distance_spheroid(point, point, spheroid)"
psql $DB -U woodb -q -A -t -c "select (distance_spheroid($START , $END,
'$SPHEROID' )/1000)::decimal(9,2) as Km;"