Get DDL data
library(tidyverse) library(lubridate) # devtools::install_github("wstolte/rwsapi") library(rwsapi) library(rlist) library(httr) # A description (in dutch) of the RWS api is found here: # https://rahuls01.github.io/?json#introduction # I wrote a package to make it a bit easier (rwsapi) # Still not easy. Any suggestions or contributions are very welcome. # first, download the metadata catalogue from the RWS data service # read the metadata set # newMetadata = T # if(newMetadata){ metadata <- rwsapi::rws_metadata() parmeta <- metadata$content$AquoMetadataLijst %>% rlist::list.flatten()%>% as_tibble() %>% rename_all(tolower) locmeta <- metadata$content$LocatieLijst %>% rename_all(tolower) mapping <- metadata$content$AquoMetadataLocatieLijst %>% rename_all(tolower) ddlMetadata <- parmeta %>% full_join(mapping) %>% full_join(locmeta) # write_delim(ddlMetadata, "ddlMetadata.csv", delim = ";") # } else{ # ddlMetadata <- read_delim("ddlMetadata.csv", delim = ";") # } # find parameter # something with pH = "zuurgraad" in Dutch # This is in the used AQUO standard a quantity (grootheid) # nitrate, for example is a parameter, with grootheid == "CONCTTE" (concentration) # Salinity and temperature are also quantities grootheid = "zuurgraad" # find exact parameter description myPar <- ddlMetadata %>% filter(grepl(grootheid, tolower(ddlMetadata$parameter_wat_omschrijving))) %>% distinct(parameter_wat_omschrijving) %>% unlist() %>% unname() myPar # looks good # testing location = "NOORDWK20" # 20 km off the coast at Noordwijk firstYear = 2000 lastYear = 2020 getDDLdata <- function(parameter_wat_omschrijving, location, firstYear, lastYear, metadata){ myMetadata <- metadata %>% filter(parameter_wat_omschrijving == parameter, code == location) %>% rename(locatie.code = code) beginDatumTijd <- paste0(firstYear, "-01-01T00:00:00.000+01:00") eindDatumTijd <- paste0(lastYear, "-12-31T23:59:59.000+01:00") getlist <- rwsapi::rws_makeDDLapiList(mijnCatalogus = myMetadata, beginDatumTijd, eindDatumTijd, "OW") df <- rwsapi::rws_observations2(getlist[[1]])$content return(df) } # test df <- getDDLdata(myPar, location, firstYear , lastYear, ddlMetadata) ## all pH data for this location # vind locations where pH is measured (these are also fresh water locations!!!) # locs <- ddlMetadata %>% filter(parameter_wat_omschrijving == myPar) %>% distinct(code, x,y, coordinatenstelsel) # from here, you would want to do a spatial query and # get data while looping over getDDLdata require(sf) locations <- locs %>% st_as_sf(coords = c("x", "y"), crs = 25831) %>% st_transform(4326) # transform to wgs84 # Marine regions from marineregions.org (VLIZ) layerurl <- "http://geo.vliz.be/geoserver/MarineRegions/ows?service=WFS&version=1.0.0&request=GetFeature&typeName=MarineRegions:eez_iho_union_v2&outputFormat=application/json" marineRegions <- sf::st_read(layerurl) # play a bit with the buffer (decimial degrees). # if zero, many freshwater stations are also selected. marineLocs <- locations %>% sf::st_intersection(st_buffer(marineRegions, -0.05)) marineRegions %>% st_intersection(marineLocs) %>% ggplot() + geom_sf() myLocs <- marineLocs %>% st_drop_geometry() %>% distinct(code) %>% unlist() %>% unname() # then loop function getDDLdata over locs for(ii in 1:length(myLocs)){ df <- getDDLdata(myPar, myLocs[ii], firstYear , lastYear, ddlMetadata) write_delim(df, paste(myPar,myLocs[ii], ".csv")) }