sf error that column names longitude and latitude reversed?

[yaakovfeldman]

I was unsure whether to post this to SO instead, or to sf, but it is possibly an issue involving one of these libraries. Apologies if this is the wrong place.

opq("grantown on spey") %>%
  add_osm_feature(key = "railway", value = "construction") %>%
  osmdata_sf %>%
  .$osm_multilines %>%
  st_length

##3558.229 m
##Warning messages:
##1: In .pointsToMatrix(p1) :
##  Suspect column names (longitude and latitude reversed?)
##2: In .pointsToMatrix(p2) :
##  Suspect column names (longitude and latitude reversed?)
##3: In .pointsToMatrix(p1) :
##  Suspect column names (longitude and latitude reversed?)
##4: In .pointsToMatrix(p2) :
##  Suspect column names (longitude and latitude reversed?)
##5: In .pointsToMatrix(p1) :
##  Suspect column names (longitude and latitude reversed?)
##6: In .pointsToMatrix(p2) :
##  Suspect column names (longitude and latitude reversed?)
##7: In .pointsToMatrix(p1) :
##  Suspect column names (longitude and latitude reversed?)
##8: In .pointsToMatrix(p2) :
##  Suspect column names (longitude and latitude reversed?)

I do not understand the warning messages as the data comes straight from osmdata. This may be an error in sf or another package, and not in osmdata as the geometry looks fine to me, but has anyone come across this? Looking at https://github.com/cran/geosphere/blob/d7a52f11998ffc06ad1a9037809aba92bd5b2c6e/R/pointsToMatrix.R#L40 which I assume is where the warning comes from, it should only be triggered when I have columns with those names - which osmdata/sf don't.

[mpadge]

It seems like a geosphere issue, as you indicate, but you could nevertheless confirm whether it is an osmdata issue by using the GDAL driver to import the data via sf. Try:

library(sf)
opq(...) %>% add_osm_feature(...) %>%
    osmdata_xml (filename = "file.osm")
st_read("file.osm", layer="multilinestring") %>% st_length()

If that gives the same warning, then it's an sf or geosphere issue. If not, please let us know here and we'll dig deeper. (st_read won't give the same result as osmdata, but it will at least give you the lat/lon columns, so should reproduce.)

[yaakovfeldman]

Thank you for the suggestion. I tried

> opq("grantown on spey") %>%
+     add_osm_feature(key = "railway", value = "construction") %>%
+     osmdata_xml (filename = "file.osm")
> st_read("file.osm", layer="multilinestrings") %>% st_length()
Reading layer `multilinestrings' from data source `\\nao-lon-home\HomesA-H$\FELD067\file.osm' using driver `OSM'
Simple feature collection with 1 feature and 4 fields
geometry type:  MULTILINESTRING
dimension:      XY
bbox:           xmin: -3.656311 ymin: 57.29804 xmax: -3.62154 ymax: 57.32346
epsg (SRID):    4326
proj4string:    +proj=longlat +datum=WGS84 +no_defs
3558.007 m

So a very slightly different number (by a few cm) and no error messages! Could this point to an issue in osmdata_sf?

[mpadge]

Your first bit of code above works fine for me with latest versions of both sf and osmdata. The difference in distances is also okay - osmdata reads far more of the actual data than does GDAL (the reader for sf), so the osmdata values should generally be more accurate.