1

According to Google, these two points are about 137 meters apart, yet the result of ST_Distance reports just a fraction of a meter? The units of the result should be in meters, as I understand it, because the reference system (3857) is in meters.

# SELECT ST_Distance(
  ST_Transform('SRID=3857;POINT(-113.936911 43.3150997)'::geometry, 3857), 
  ST_Transform('SRID=3857;POINT(-113.9369126 43.3139828)'::geometry, 3857));
      st_distance
-----------------------
 0.0011169011460281386

Why does ST_Distance report this unexpected result?

2
  • 4
    Because you meant to use SRID=4326 as the EWKT prefix? Right now you're mapping in quadrant II on Null Island.
    – Vince
    Jan 5 at 1:14
  • 5
    BTW: You should not ever trust distances generated in 3857 (they're always wrong) -- SELECT ST_Distance( 'SRID=4326;POINT(-113.936911 43.3150997)'::geography, 'SRID=4326;POINT(-113.9369126 43.3139828)'::geography) as dist returns 124.0863924
    – Vince
    Jan 5 at 1:20

1 Answer 1

5

There's quite a few problems here.

First off, your ST_Transform from 3857 to 3857 doesn't actually do anything. The Extended Well-Known Text prefix incorporates the SRID of the coordinates that follow, so you asserted that the {-113,43} values were already in Web Mercator, and the values didn't change. The "distance" was therefore the Cartesian distance between the decimal degrees (DD) values.

If you had used the SRID for WGS 1984 (4326) in the EWKT, then the result would have looked like:

SELECT ST_Distance(
  ST_Transform('SRID=4326;POINT(-113.936911 43.3150997)'::geometry, 3857), 
  ST_Transform('SRID=4326;POINT(-113.9369126 43.3139828)'::geometry, 3857));
st_distance
----------------
170.881023396406

The next problem is that distance is pretty much always wrong in Web Mercator. The projection was designed to support lines of constant bearing, so distance from the Equator increases to infinity at the poles. The correct way to compute distance starting with decimal degrees is to use the geography type (without transform):

SELECT ST_Distance(
  'SRID=4326;POINT(-113.936911 43.3150997)'::geography, 
  'SRID=4326;POINT(-113.9369126 43.3139828)'::geography);
st_distance
----------------
124.0863924

Now, you'll note that this value is different than your reference value of "137 meters", and this is the last problem. Looking for confirmation of the PostGIS value, I tried using ArcPy (of ArcMap 10.5.1) with three different methods,

>>> import arcpy
>>> sr = arcpy.SpatialReference(4326)
>>> p1 = arcpy.PointGeometry(arcpy.Point(-113.936911,43.3150997),sr)
>>> p2 = arcpy.PointGeometry(arcpy.Point(-113.9369126,43.313982),sr)
>>> p1.angleAndDistanceTo(p2,'GEODESIC')
(-179.94010923675944, 124.17527144338943)
>>> p1.angleAndDistanceTo(p2,'GREAT_ELLIPTIC')
(-179.9401092387172, 124.17527144402266)
>>> p1.angleAndDistanceTo(p2,'LOXODROME')
(-179.94010978556503, 124.1752714433236)

I also tried ST_Distance with a False in the use_spheroid parameter:

SELECT ST_Distance(
  'SRID=4326;POINT(-113.936911 43.3150997)'::geography, 
  'SRID=4326;POINT(-113.9369126 43.3139828)'::geography,False);
st_distance
----------------
124.19385202

So it seems that 124 is more correct than 170 or 137, but the correctness of 124.086 is undetermined.

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.