calculate st_distance between the point in km (m)

Question

I try to calculate the distance between two geogprahical points, using st_distance. The points are stored in a data frame, but for simplicity, I take one pair of points to show the problem.

Here are the two points (an example of a pair of points):

a <-st_sfc(st_point(c(70.321117, 22.896517)))

b <-st_sfc(st_point(c(68.743738549303, 17.0129893190695)))

st_distance(a,b)

Which returns 6.09.

I understand this should be a distance in degrees, and to calculate the same in km I need to multiply by 111 km. But the distance seems to be too large (ca 666 km), compared to online calculators of great circle distance (I got 293 km).

Why is the difference so large? I do I calculate the distances in a wrong way?

Answer 1

sf package is tightly coupled with unitspackage, together they do really neat things, as giving you the units of the output of a function; in your case, the objects lack CRS, once set, you'll get:

a <- st_sfc(st_point(c(70.321117, 22.896517))) %>% 
  st_set_crs(4326)

b <- st_sfc(st_point(c(68.743738549303, 17.0129893190695))) %>% 
  st_set_crs(4326)

st_distance(a,b)

Units: [m]
         [,1]
[1,] 671910.7

From the function help: "Compute Euclidian or great circle distance between pairs of geometries; compute, the area or the length of a set of geometries."; on the parameter which :" character; for Cartesian coordinates only: one of Euclidean, Hausdorff or Frechet; for geodetic coordinates, great circle distances are computed; see details"

Answer 2

I think your coordinates might be reversed?

When I initially copied your coordinates, they put me out in the Indian ocean, and the distance was almost 700KM. But I thought that might be wrong, so I reversed them, and I think it works.

Here's the SQL I used (referencing a blog post below as well):

;with cte_geog as (
select 

st_makepoint(22.896517, 70.321117)::geography as geog_a
, st_makepoint(17.0129893190695, 68.743738549303)::geography as geog_b

)

select 

st_distance(d.geog_a, d.geog_b)/1000 as distance_km
, st_makeline(d.geog_a::geometry, d.geog_b::geometry)::geography as geog_line
, *
from cte_geog as d

The distance returned is 289.15KM

Here's the result in the DBeaver spatial viewer:

Referencing this blog post from Paul Ramsey (https://info.crunchydata.com/blog/postgis-and-the-geography-type), as long as your working with Geography the ST_Distance() function will do all the work of calculating the distances on the sphere.