Measure distance for spatialLine in R returns bad value

I have this spatialLine in R:

library(sp)
library(maptools)
library(SDraw)
library(leaflet)

x <- c(153.9, 152.67072,147.9879,144.25364,135.45557,116.03998,76.92708,45.10622,37.02522,36.76888,36.41,37.19583,41.55226,50.99487,53.18)

y <- c(73.52,75.80461,78.28467,80.82467,83.32397,85.49637,86.39683,85.16955,82.60437,78.01859,73.35001,72.87417,71.23226,71.55405,68.73)

d <- SpatialLines(list(Lines(Line(cbind(x,y)), ID='a')))
crs(d) <- "+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs"

leaflet(data = d) %>%
addPolylines() # just plot the line


I want to measure its length. I tried:

> lineLength(d, byid = T)
[1] 148.407


However, this is not good (the actual size of this would be aroung ~4000 km as you can see). How can I actually measure the path length?

• The measurement is correct, if you want the output in Cartesian degrees (which are meaningless). You probably want to sum on geodetic distance. – Vince Jan 12 at 23:37
• try sp::LineLength(d, longlat = TRUE) – mdsumner Jan 14 at 5:24

Edit: use geosphere...

There's a few packages that can compute true distance given lat-long coordinates, for example fields has the rdist.earth and rdist.earth.vec functions.

rdist.earth.vec takes two lat-long matrix arguments with the same number of rows and computes the distance from the points in matching rows. So for your simple single line feature we can get all the coordinates with:

> xy = coordinates(d)[[1]][[1]]


By dropping off the last row and the first row we can create two matrices where the rows are points 1 and 2, and then 2 and 3, and then 3 and 4 and so on, and then sum the total distance between those pairs:

> sum(fields::rdist.earth.vec(xy[-15,],xy[-1,], miles=FALSE))
[1] 4367.507


I think the package authors are from the USA, so don't forget to use miles=FALSE to get kilometres!

Or even easier, use geosphere package:

> library(geosphere)
> lengthLine(d)
[1] 4380327


This is in metres, and is slightly different to the rdist.earth measurement. It uses an ellipsoid model for the earth rather than a spherical model that rdist.earth uses. It also works directly on Spatial data objects.

Another solution is to convert lat-long to a cartesian coordinate system and then use pythagoras theorem, but that is tricky for data over a large area and is likely to be less accurate.

• thank you very much! I have multiple features in the spatialLine in fact but that should be easy to do with a loop. thanks for the help – Andrei Niță Jan 13 at 12:18
• I think geosphere is more accurate and convenient - see edit. – Spacedman Jan 13 at 13:13