'Calculate length of LineString in meters using Python
I know this question has been asked before, but I just cannot figure it out. I'm just trying to find the length in meters from a LineString that I've made from two points ('EPSG:4326'). In R I get around 311 meters which is about the same as google maps straight line distance give using this code:
library(purrr)
library(sf)
library(dplyr)
xy <- list(c(59.94661244955431, 10.72052329576658), c(59.948009419166226, 10.725362000881619))
lat <- xy %>%
map(1) %>%
unlist()
lon <- xy %>%
map(2) %>%
unlist()
st_as_sf(data.frame(lon, lat), coords=c("lon", "lat"), crs="EPSG:4326") %>%
summarise(do_union = TRUE) %>%
st_cast("LINESTRING") %>%
st_length()
311.0174 [m]
When I do it in Python I instead get 621 even if I reproject to 'EPSG:3857'. What am I doing wrong here? Python code below:
from pyproj import Transformer
from shapely.geometry import LineString
import shapely.ops as sp_ops
latlon = [[59.94661244955431, 10.72052329576658], [59.948009419166226, 10.725362000881619]]
xy = [(lon, lat) for lat, lon in latlon]
line = LineString(xy)
transformer = Transformer.from_crs('EPSG:4326', 'EPSG:3857', always_xy=True)
line_transformed = sp_ops.transform(transformer.transform, line)
line_transformed.length
Out[147]: 621.7406578605734
Solution 1:[1]
EPSG:3857 is useless for any kind of computation. It is nominally in meters but, because of the brutal distortion of the Mercator projection, those meters don't really mean anything. Taking a guess from your name and your sample coordinates your are in Norway, meaning that some UTM projection is probably a good bet.
Solution 2:[2]
Update: As mentioned I need to use the right projection which in this case is EPSG:32632 (UTM 32N):
from pyproj import Transformer
from shapely.geometry import LineString
import shapely.ops as sp_ops
latlon = [[59.94661244955431, 10.72052329576658], [59.948009419166226, 10.725362000881619]]
xy = [(lon, lat) for lat, lon in latlon]
line = LineString(xy)
transformer = Transformer.from_crs('EPSG:4326', 'EPSG:32632', always_xy=True)
line_transformed = sp_ops.transform(transformer.transform, line)
line_transformed.length
Out[2]: 311.9284270277389
Sources
This article follows the attribution requirements of Stack Overflow and is licensed under CC BY-SA 3.0.
Source: Stack Overflow
| Solution | Source |
|---|---|
| Solution 1 | Ture PÄlsson |
| Solution 2 | Kjetil |
