When matching points to vertices, all points are matched to the same one

[marcboschmatas]

I am trying to compute the walking distance to all swimming pools in Barcelona, but when trying to match them to vertices in the street network, they all get matched to the same one. This is the code I'm using.

library(sf)
library(tidyverse)
library(osmdata)
library(dodgr)
library(tmap)
library(tmaptools)
library(httr)
library(jsonlite)
# swimming pools
# download from BCN open data API
sp <- httr::GET("https://opendata-ajuntament.barcelona.cat/data/api/action/datastore_search?resource_id=0508ed3c-362b-4f1b-8980-bbdc06358155&limit=10000") |> 
  content("text") |> 
  fromJSON()

sp <- sp$result

sp <- sp$records

# Transform into SF object
sp_sf <- st_as_sf(sp, coords = c("geo_epgs_4326_x", "geo_epgs_4326_y"), crs = "EPSG:4326")

q <- opq ("Barcelona") %>%
  add_osm_features(features = c ("\"highway\"=\"primary\"",
                                 "\"highway\"=\"secondary\"",
                                 "\"highway\"=\"tertiary\"",
                                 "\"highway\"=\"living_street\"",
                                 "\"highway\"=\"pedestrian\"",
                                 "\"highway\"=\"steps\"",
                                 "\"highway\"=\"footway\""))

streets <- osmdata::osmdata_sf(q)

streets <- streets$osm_lines

streets <- st_transform(streets, "EPSG:4326") # it doesn't work with projected coordinates

# transform sf to graph - create vertices
streets_graph <- dodgr::weight_streetnet(streets, wt_profile = "foot")
v <- dodgr_vertices(streets_graph)

# create "from" (every node in the network)

from_x <- v$x
from_y <- v$y

# create "to" (Every swimming pool)

to_x <- as.double(sp$geo_epgs_4326_x)
to_y <- as.double(sp$geo_epgs_4326_y)

# match points to the nearest vertex in network
tos <- match_points_to_verts(v, sp_sf, connected = TRUE)
froms <- match_points_to_verts(v, cbind(from_x, from_y), connected = TRUE)

froms <- v[froms, c("x", "y")]
tos <- v[tos, c("x", "y")]

then tos returns the following.

head(tos)
                x        y
179898   2.187494 41.36358
179898.1 2.187494 41.36358
179898.2 2.187494 41.36358
179898.3 2.187494 41.36358
179898.4 2.187494 41.36358
179898.5 2.187494 41.36358

Thanks!

[mpadge]

Thanks @marcboschmatas. That one was pretty easy to debug: your sp_sf data has the coordinates the wrong way around:

library(sf)
library(jsonlite)
library(httr)
sp <- httr::GET("https://opendata-ajuntament.barcelona.cat/data/api/action/datastore_search?resource_id=0508ed3c-362b-4f1b-8980-bbdc06358155&limit=10000") |> 
  content("text") |> 
  fromJSON()

sp <- sp$result$records
sp_sf <- st_as_sf(sp, coords = c("geo_epgs_4326_x", "geo_epgs_4326_y"), crs = "EPSG:4326")

print (sp_sf$geometry)
#> Geometry set for 132 features 
#> Geometry type: POINT
#> Dimension:     XY
#> Bounding box:  xmin: 41.36357 ymin: 2.104831 xmax: 41.45992 ymax: 2.21323
#> Geodetic CRS:  WGS 84
#> First 5 geometries:
#> POINT (41.45147 2.19313)
#> POINT (41.41037 2.171895)
#> POINT (41.39369 2.118332)
#> POINT (41.39736 2.142036)
#> POINT (41.41538 2.158172)

^{Created on 2022-09-20 with reprex v2.0.2}

Barcelona is at around long = 2.2, lat = 41.4, but that has long = 41.4, lat = 2.2. I recently encountered exactly this issue myself, and spent a while trying to debug before realising that was the cause. Switching those columns makes everything work fine, so i'm going to close this now, but feel free to ask any further questions.

Note also that the point matching algorithms were expanded in the previous release of dodgr, so the match_pts_to_graph() function now matches points to the nearest perpendicular point of intersection of graph edges, rather than just the nearest vertex. In case that helps with what you're trying to do here.