st_read(): WFS data importeren in R met behulp van GeoJSON


Als een WFS ruimtelijke gegevens aanbiedt in het GeoJSON formaat, dan kan deze informatie eenvoudig in R ingelezen worden. Vervolgens kunnen deze data binnen R in allerlei ruimtelijke analyses gebruikt worden. WFS staat voor Web Feature Service, dit is een standaard van het Open Geospatial Consortium (OGC).

De webkaart hieronder - met daarop alle 20 nationale parken in Nederland - is gemaakt met behulp van R, met slechts enkele regels code. En in deze blog laten we zien hoe dat gedaan is.


Veel Open Data wordt aangeboden als WFS

Veel openbare geo-informatie die in Nederland wordt aangeboden op platforms als PDOK en het NGR, is beschikbaar als WFS. Voor deze demonstratie is oog ons gevallen op de dataset met nationale parken (via dit GetCapabilities request).

Het standaard output formaat van een WFS service is GML versie 3.2. Het nadeel van dit GML formaat is, dat het erg lastig verder te verwerken is. Gelukkig maken veel WFS severs - in Nederland - gebruik van de optie om ook andere output formaten, zoals GeoJOSN, aan te bieden. GeoJSON is een geografisch datauitwisselingsformaat dat gebaseerd is op JavaScript Object Notation (JSON). En GeoJSON kan eenvoudig ingelezen worden door de functie st_read() in R.

Klik op de volgende links om het verschil tussen de twee outputformaten te zien:

(Voor als het NGR ooit besluit om deze service uit de lucht te halen, staat hier een backup van de GeoJSON output, zoals die op het moment van schrijven van deze blog gegeven werd.)

R Code

We beginnen met het laden van de benodigde packages:

library(sf)
## Linking to GEOS 3.6.1, GDAL 2.2.3, proj.4 4.9.3
library(tmap)

Vervolgens definiëren we de URL - voor de leesbaarheid gebeurt dat hier in twee losse delen die vervolgens aan elkaar worden geplakt:

baseurl <- "https://geodata.nationaalgeoregister.nl/nationaleparken/wfs?"
# Verzoek om output in GeoJSON
wfs_request <- "request=GetFeature&service=WFS&version=2.0.0&typeName=nationaleparken&outputFormat=application%2Fjson"

nl_nationale_parken_wfs <- paste0(baseurl,wfs_request)

Deze URL geven we als input aan de functie st_read() om de Simple feature collection nl_nationale_parken aan te maken. (Een Simple feature collection is - kort gezegd - een R data frame met een geometrie kolom.)

nl_nationale_parken <- st_read(nl_nationale_parken_wfs)
## Reading layer `OGRGeoJSON' from data source `https://geodata.nationaalgeoregister.nl/nationaleparken/wfs?request=GetFeature&service=WFS&version=2.0.0&typeName=nationaleparken&outputFormat=application%2Fjson' using driver `GeoJSON'
## Simple feature collection with 20 features and 8 fields
## geometry type:  MULTIPOLYGON
## dimension:      XY
## bbox:           xmin: 36781.02 ymin: 350383.2 xmax: 246437.1 ymax: 615251.4
## epsg (SRID):    28992
## proj4string:    +proj=sterea +lat_0=52.15616055555555 +lon_0=5.38763888888889 +k=0.9999079 +x_0=155000 +y_0=463000 +ellps=bessel +towgs84=565.2369,50.0087,465.658,-0.406857,0.350733,-1.87035,4.0812 +units=m +no_defs

Je kunt nl_nationale_parken eventueel plotten of viewen (maar dat hebben we hier niet gedaan…)

plot(st_geometry(nl_nationale_parken))
View(nl_nationale_parken)

Het leukste is natuurlijk om de dataset in RStudio in een interactief kaartvenster te bekijken, en dat doen we met tmap:

tmap_mode("view")
qtm(nl_nationale_parken, borders = "red", fill = "darkgreen", basemaps = c("Esri.WorldTopoMap","Esri.NatGeoWorldMap","CartoDB.Positron", "OpenStreetMap"))

Deze laatste 2 regels code geven de interactieve kaart die hierboven op deze pagina getoond wordt.


GeoJSON is te verkiezen boven GML

Hierboven werd gesteld dat GML lastig verder te verwerken is. Om dat te demonstreren, proberen wij hieronder de GML output in te lezen met st_read():

baseurl <- "https://geodata.nationaalgeoregister.nl/nationaleparken/wfs?"
# Verzoek om output in GML
wfs_request2 <- "request=GetFeature&service=WFS&version=2.0.0&typeName=nationaleparken"

nl_nationale_parken_wfs2 <- paste0(baseurl,wfs_request2)

nl_nationale_parken2 <- st_read(nl_nationale_parken_wfs2)
## Reading layer `nationaleparken' from data source `https://geodata.nationaalgeoregister.nl/nationaleparken/wfs?request=GetFeature&service=WFS&version=2.0.0&typeName=nationaleparken' using driver `GML'
## Simple feature collection with 20 features and 8 fields
## geometry type:  MULTISURFACE
## dimension:      XY
## bbox:           xmin: 36781.02 ymin: 350383.2 xmax: 246437.1 ymax: 615251.4
## epsg (SRID):    NA
## proj4string:    NA

We zien dat GDAL de GML wel degelijk kan lezen en converteren. Er zij echter twee problemen met deze output:

Dus de boodschap van deze blog blijft: als een WFS server GeoJSON als outputformaat aanbiedt, dan verdient dat sterk de voorkeur boven GML. Daarom is het - zoals hierboven al aangegeven - raadzaam om bij de bevraging van een WFS server de parameter outputFormat=application/json mee te geven.


Egge-Jan Pollé