3

I have the following linestring which approximately marks an airport runway:

LINESTRING(-381578.819887 6574044.360475, -378513.534443 6574799.475787)  # EPSG 3857
LINESTRING(-3.427781 50.732117, -3.400245 50.73641)                       # EPSG 4326

If I use the GLength function, I get a length of

SELECT GLength(GeomFromText('LINESTRING(-381578.819887 6574044.360475, -378513.534443 6574799.475787)'))
3156.925m   # EPSG 3857
SELECT GLength(GeomFromText('LINESTRING(-3.427781 50.732117, -3.400245 50.73641)'))
0.027869°   # EPSG 4326 ~ 0.027869 x 111319.9 = 3102m

The actual length of this runway is 2076m. QGIS will display this correct distance if I check the ellipsoidal box using the measure tool, so I assume this performs a projection/trigonometry conversion. Yet, I can't seem to find the correct calculation to replicate this conversion and go from the (planar?) distance returned by the GLength function to the real distance.

Can anyone direct me to a tutorial for this?

This question is a continuation of my pervious question How can I get the ellipsoidal values for length and area using spatialite? which got side tracked into correcting a different issue.

3

ST_Length() (which is an alias for GLength(), but I prefer the consistent naming) can take a second parameter which is a flag to tell it to use the WGS84 ellipsoid rather than the sphere, at the cost of computation speed. Naturally, only WGS84 coordinates are allowed with this second function, so if your coordinates are in EPSG:3857, you'll need to transform them first.

So running the query:

SELECT ST_Length(ST_LineFromText('LINESTRING(-3.427781 50.732117, -3.400245 50.73641)', 4326), 1) AS WGS84,
ST_Length(ST_LineFromText('LINESTRING(-381578.819887 6574044.360475, -378513.534443 6574799.475787)', 3857)) AS WM_sphere,
ST_Length(ST_Transform(ST_LineFromText('LINESTRING(-381578.819887 6574044.360475, -378513.534443 6574799.475787)', 3857), 4326), 1) AS WM_ellipsoid;

Gives me the result:

WGS84             WM_sphere         WM_ellipsoid
----------------  ----------------  ----------------
2001.78486098528  3156.92476749369  2001.78741760296

Which is much closer, but not spot on. Using the measure tool in Google Maps for EGTE RWY 08/26 gives me:

Total distance: 2023.68 m

I'm not sure how they've arrived at the official figure of 2076 m, unless it's a particularly lumpy runway. But using the ellipsoid version of ST_Length() should get you a lot closer.

  • Thanks for your expert help MerseyViking. I get the response "wrong number of arguments to function ST_Length()" when using the above. Should this work with spatialite? Do I need a specific version of spatialite? – Dan Jan 21 '14 at 13:27
  • 1
    Hmmm... odd. Yes, I ran it on SpatiaLite 4.1.0 (and just tried it again in case I had a copy/paste error) command line as well as SpatiaLite GUI. The rest of my version info is: PROJ.4 version ......: Rel. 4.8.0, 6 March 2012 GEOS version ........: 3.3.8-CAPI-1.7.8 and SQLite version ......: 3.7.13 – MerseyViking Jan 21 '14 at 13:39
  • 2
    Seems it is a version issue. More recent versions of the spatialite DLL (e.g. 4.0.0) run the query without issue. – Dan Jan 21 '14 at 16:54

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.