0

I have been using osmdata to pull in OpenStreetMap data. The queries that I have been running are too large and some data sets are missing.

bbuk <- getbb("Great Britain")
x <- bbuk %>%
  opq(timeout = 25*1000)%>%
  add_osm_feature("plant:source","wind") %>%
  osmdata_sf()
y <- bbuk %>%
  opq(timeout = 25*100)%>%
  add_osm_feature("power","line") %>%
  osmdata_sf()

osmextract seems a better alternative but I cant figure out how to make the query specific for oe_get()

library(osmextract)
place_name = "Great Britain"
et = c("amenity")
q_points = "SELECT * FROM multipolygon WHERE power IN ('wind')"
oe_school_points = oe_get(place_name, query = q_points, extra_tags = et)

Is there a source to structure the query in oe_get() to pull in the specific data from the pdf that is needed>

Thanks,

user1966593
  • 69
  • 10

1 Answers1

1

The two queries described in the first chunk of code can be replicated using osmextract as follows:

# packages
remotes::install_github("ropensci/osmextract", quiet = TRUE, upgrade = "never")
library(osmextract)
#> Data (c) OpenStreetMap contributors, ODbL 1.0. https://www.openstreetmap.org/copyright.
#> Check the package website, https://docs.ropensci.org/osmextract/, for more details.
library(sf)
#> Linking to GEOS 3.10.2, GDAL 3.4.1, PROJ 7.2.1; sf_use_s2() is TRUE
options(timeout = 600L)

# Download and convert data. NB The following takes several minutes
oe_get("England", layer = "lines", extra_tags = "power", download_only = TRUE)
oe_get("England", layer = "multipolygons", extra_tags = c("power", "plant:source"), download_only = TRUE)

boundary <- oe_get(
  place = "England", 
  query = "SELECT * FROM multipolygons WHERE type = 'boundary' AND name = 'England'", 
  quiet = TRUE
)

# Analogous to add_osm_feature("plant:source","wind")
es1 <- oe_get(
  place = "England", 
  query = "SELECT * FROM multipolygons WHERE power = 'plant' AND plant_source = 'wind'", 
  quiet = TRUE
)

# Analogous to add_osm_feature("power","line")
es2 <- oe_get(
  place = "England", 
  query = "SELECT * FROM lines WHERE power = 'line'", 
  quiet = TRUE
)

Plot

par(mar = rep(0, 4))
plot(st_geometry(boundary), reset = FALSE)
plot(st_geometry(es1), add = TRUE, col = "blue")
plot(st_geometry(es2), add = TRUE, col = "grey")

Created on 2022-08-17 by the reprex package (v2.0.1)

agila
  • 3,289
  • 2
  • 9
  • 20