Client for OGC API - Environmental Data Retrieval (EDR) • edr4r

An R client for OGC API - Environmental Data Retrieval (EDR) services that expose JSON discovery metadata and CoverageJSON, GeoJSON, or CSV query responses. The spec is general, but in practice this package gets the most use against in-situ monitoring networks — stream gauges, weather stations, snow telemetry, reservoir telemetry — that expose their stations and time series through EDR.

Two known-good places to point it:

USGS waterdata OGC API — stream gauges and water-quality stations from the U.S. Geological Survey.
Western Water Datahub — a pygeoapi deployment that wraps RISE, SNOTEL, USACE, AWDB and other monitoring sources behind a single EDR endpoint.

For cross-server experiments, the Met Office Labs EDR demonstrator is another useful endpoint. It is a technical demonstrator, not an operational service: availability, collections, and response details can change without notice, so do not build production workflows around it. The cross-endpoint Lake Mead vignette shows its 2015 population grid alongside USGS river discharge and Western Water Datahub reservoir storage without mixing their provenance or units.

The goal is to take the tedious parts of EDR off your hands — URL construction, comma-separated parameter lists, WKT coordinate encoding, retries, content negotiation — and hand back something you can actually do data analysis with:

CoverageJSON → a long tibble (one row per coverage × parameter × domain position), via covjson_to_tibble(). Nonstandard dimensions are retained as .axis_* columns and normalized domain type, axis summaries, and effective referencing remain available in the edr_covjson_metadata attribute.
GeoJSON → an sf object, via geojson_to_sf().
CSV → a tibble, parsed directly by the query helper.

Installation

Install the current stable release from CRAN:

install.packages("edr4r")

The 0.3.0 source release is also available from GitHub:

# install.packages("pak")
pak::pak("ksonda/edr4r@v0.3.0")

# Follow the mutable development branch instead:
pak::pak("ksonda/edr4r")

# or
# install.packages("remotes")
remotes::install_github("ksonda/edr4r@v0.3.0")

The v0.2.0-rc.1 and v0.3.0-rc.1 tags were GitHub-only previews. Use the final v0.3.0 tag for a reproducible installation or the default branch to follow ongoing development.

For local development:

git clone https://github.com/ksonda/edr4r.git
cd edr4r
R -e 'devtools::install()'

Requires R >= 4.1. The sf package is optional but recommended (used to turn location lists and GeoJSON into spatial objects).

Quick start

Start by pointing a client at a server. The base URL is the only thing it really needs:

library(edr4r)

client <- edr_client("https://api.waterdata.usgs.gov/ogcapi/beta")
wwdh <- edr_client("https://api.wwdh.internetofwater.app")
# or use "http://localhost:5005" for a local pygeoapi deployment

collections <- edr_collections(client)
collections[, c("id", "title", "data_queries", "output_formats")]

Collection IDs are service-specific. The first thing to do against a new service is run edr_collections() and read the data_queries column to see which EDR endpoints each collection supports.

For a new or unfamiliar implementation, inspect its advertised support before issuing data queries:

edr_capabilities(client, "daily-edr")
edr_supports(client, "daily-edr", query = "locations")
edr_diagnose(client, "daily-edr")

edr_supports() reports what metadata advertises; FALSE is not proof that a partially conformant server cannot handle the request.

Discovery metadata is cached per client for a short, configurable period. Use refresh = TRUE when current server state matters, or edr_cache_clear(client) to clear it explicitly.

To try the non-operational Met Office demonstrator with a deliberately small request, query one terrain point rather than a forecast collection:

met <- edr_client(
  "https://labs.metoffice.gov.uk/edr",
  timeout = 10,
  max_tries = 1
)

terrain <- edr_position(
  met,
  "terrain_tiles",
  coords = c(-0.1276, 51.5072),
  parameter_name = "Height"
)
covjson_to_tibble(terrain)

This example is also exercised by a scheduled, non-blocking live smoke check; it is never run as part of CRAN checks or the regular test suite.

Collections representing model runs may advertise instances. The same query verbs work below an instance when instance_id is named explicitly:

runs <- edr_instances(met, "moglobal-station-level")
run_id <- runs$id[[1]]

run_capabilities <- edr_capabilities(
  met, "moglobal-station-level", instance_id = run_id
)
edr_supports(
  run_capabilities, query = "locations"
)
run_locations <- edr_locations(
  met, "moglobal-station-level",
  instance_id = run_id
)

Find stations

edr_locations() returns one server response by default. If the response is GeoJSON and sf is installed, it is promoted to an sf object automatically. For a complete result from a server that advertises rel = "next", opt into bounded pagination:

piedmont_bbox <- c(-78.60, 36.04, -78.28, 36.22)

locs <- edr_locations(
  client, "daily-edr",
  bbox = piedmont_bbox,
  limit = 100,              # server page size
  paginate = TRUE,
  max_pages = 10,
  max_features = 100
)
locs                            # sf POINTs with station attributes
plot(sf::st_geometry(locs))

Pagination follows the server’s next URL as an opaque cursor or offset. It stops with a typed error if a page/feature cap is reached while another page still exists, so a bounded result is never presented as complete.

Pull a time series for one station

Once you know a station ID, ask for its values. The server returns CoverageJSON; covjson_to_tibble() flattens it into one row per (coverage × parameter × domain position):

station_id <- locs$id[[1]]      # USGS ids include the "USGS-" prefix

resp <- edr_location(
  client, "daily-edr",
  location_id    = station_id,
  parameter_name = "00060",    # daily discharge
  limit          = 200
)

df <- covjson_to_tibble(resp)
df[, c("datetime", "value", "unit")]

The USGS beta location endpoint currently ignores datetime and returns its latest limit records, so filter the resulting tibble client-side when you need a shorter period. Other EDR implementations, including WWDH, honor bounded intervals.

Pull several station time series safely

When you already have station IDs, edr_location_batch() makes one bounded, sequential edr_location() request per ID and keeps request provenance and failures visible:

selected_ids <- head(locs$id, 10)

pull <- edr_location_batch(
  client, "daily-edr",
  location_id    = selected_ids,
  parameter_name = "00060",
  limit           = 200,
  max_requests   = length(selected_ids),
  on_error       = "collect"
)

pull$data                  # .request_id and .location_id identify every row
pull$errors                # typed, empty tibble when every request succeeded
pull$requests              # success / empty / error status for every ID

The helper deliberately does not discover stations or parallelize requests; the full request count is known and validated before network activity.

For a server that honors bounded datetime intervals, split a long pull into calendar windows with chunk. The request cap applies to the complete station-by-window plan, so ten stations over twelve monthly windows requires a cap of at least 120:

monthly_pull <- edr_location_batch(
  wwdh, "rise-edr",
  location_id    = "3514",
  datetime       = "2023-01-01/2024-01-01",
  parameter_name = "3",
  chunk           = "1 month",
  checkpoint      = "lake-mead-2023-checkpoint",
  resume          = TRUE,
  include_parameters = TRUE,
  max_requests    = 12,
  on_error        = "collect"
)

monthly_pull$parameters[
  monthly_pull$parameters$id == "3",
  c("id", "name", "description", "unit_symbol", "unit_definition")
]

Adjacent windows share their boundary because EDR intervals are closed and server inclusivity varies. By default, exact observations repeated across windows for the same station are retained only from the earliest request; deduplicate = FALSE preserves the raw responses. requests$n_rows records the pre-deduplication row count for each request.

resume = TRUE initializes a missing checkpoint on the first run. Later runs with the same effective request plan load successful and empty windows from disk and retry only unresolved work. The expanded plan still has to satisfy max_requests before the checkpoint is opened. Checkpoints contain parsed observations, but do not store client headers, query URLs, or error conditions; protect the directory like any other local data extract.

include_parameters = TRUE makes one explicit, cacheable discovery request and attaches the full collection parameter catalog once at monthly_pull$parameters; definitions and units are not repeated on every observation row. This discovery request is not counted by max_requests, is not persisted in the checkpoint, and is not collected by on_error if it fails.

The USGS beta location endpoint currently ignores datetime and returns its latest records. Chunking is therefore useful only for endpoints that honor the requested interval; it is not a workaround for retrieving USGS history.

Spatial filters — bbox and polygon

To grab everything inside a rectangle, use edr_cube():

cube <- edr_cube(
  wwdh, "rise-edr",
  bbox           = c(-115.5, 35.5, -114.5, 36.5),
  datetime       = "2023-01-01/2023-03-31",
  parameter_name = "3"
)
covjson_to_tibble(cube)

For an arbitrary polygon, edr_area() takes WKT, an sf polygon, or a matrix of (lon, lat) rows (it’ll close the ring for you):

ring <- matrix(
  c(-115.5, 35.5, -114.5, 35.5, -114.5, 36.5, -115.5, 36.5),
  ncol = 2, byrow = TRUE
)
area <- edr_area(
  wwdh, "rise-edr",
  coords = ring,
  datetime = "2023-01-01/2023-03-31",
  parameter_name = "3"
)
covjson_to_tibble(area)

Plot a time series

edr_plot() is a small ggplot2 wrapper over the tidy tibble:

edr_plot(resp)            # accepts an edr_response directly

Facets by parameter (so different units don’t share a y-axis) and colours by station. Add layers or themes like any other ggplot.

It also auto-detects common non-station shapes:

edr_plot(cube)            # x/y grid -> tile map
edr_plot(profile)         # varying z -> vertical profile

# or force the layout
edr_plot(profile, view = "profile")
edr_plot(cube, view = "grid")

Map stations with per-station popups

edr_map() puts the stations on a leaflet basemap. Pass data = as a named list keyed by station id (the shape [edr_explore()] produces) and each marker gets a popup with an inline plot and a “Download CSV” link for that station’s data — embedded as a data: URI so the saved HTML is selfcontained:

station <- locs[locs$id == station_id, ]
data_list <- stats::setNames(list(df), station_id)
m <- edr_map(station, data = data_list, popup = "plot+csv")
edr_save_html(m, "stations.html")

For a quick exploratory pass over a whole collection, edr_explore() does the fetch + plot + map in one call:

edr_explore(
  client, "daily-edr",
  bbox           = piedmont_bbox,
  parameter_name = "00060",
  limit          = 25,
  file           = "snapshot.html"
)

Gridded coverages and vertical profiles can be mapped too. edr_map() detects tidy CoverageJSON grids/profiles and puts slice selectors inside the leaflet widget when there are multiple parameters, datetimes, vertical levels, or nonstandard dimensions such as ensemble realisations:

grid <- covjson_to_tibble(cube)
edr_map(grid)

profile <- covjson_to_tibble(profile_resp)
edr_map(profile)

Coverage maps require WGS 84/CRS84 longitude/latitude coordinates. A declared projected or other known geographic CRS is rejected rather than passed to Leaflet as WGS 84 degrees. Missing or custom horizontal references are accepted only after coordinate and inferred cell-bound checks, with a warning; request a geographic response such as crs = "CRS84" when the endpoint supports it. Station sf geometries are transformed to WGS 84 for display while spatial matching and max_match_distance remain in the source CRS; missing sf CRS metadata warns before plausible degree-range coordinates are used as-is.

When a coverage has a custom dimension, its selector can be initialized by the original coordinate name:

edr_map(ensemble_grid, initial = list(realisations = "control"))

Coverage and station layers can share one widget. Start with the coverage map, then add independently styled station groups with their normal chart/CSV popups:

m <- edr_map(grid, grid_transform = "sqrt")
m <- edr_add_stations(
  m, stations,
  data = data_list,
  popup = "plot+csv",
  group = "USGS"
)

edr_explore() uses the same behavior for bulk coverage queries. Use output = "plot" when you want a ggplot instead of the interactive map:

edr_explore(client, "gridded-collection",
            bbox = c(-120, 39, -118, 41),
            method = "cube")

Here "gridded-collection" is intentionally schematic: substitute a collection whose data_queries metadata advertises cube. The concrete WWDH rise-edr cube above and the Met Office examples in the cross-endpoint vignette are runnable counterparts.

Weird IDs, CSV, and an escape hatch

Some monitoring networks use compound station IDs — colon-separated triplets are a common pattern. The client URL-encodes reserved characters for you. Because the identifier field varies by deployment, this call is schematic and must use an ID advertised by that collection:

edr_location(your_client, "station-network", "1185:CO:SNTL",
             datetime = "2024-01-01/2024-01-03")

If a collection advertises CSV, ask for it instead of CoverageJSON by setting the parser format. Collection IDs and format support are server-specific, so this low-level form is deliberately schematic:

edr_location(client, "collection-id", "location-id", format = "csv")

And if you need to hit an endpoint or encoding the package doesn’t wrap, edr_request() is the raw escape hatch:

service_description <- edr_request(
  client, "api", format = "raw", parse = FALSE
)

API at a glance

Function	EDR endpoint
`edr_client()`	construct a client
`edr_landing()` / `edr_conformance()`	`/`, `/conformance`
`edr_collections()` / `edr_collection()`	`/collections`
`edr_capabilities()` / `edr_supports()` / `edr_diagnose()`	inspect advertised support
`edr_cache_clear()`	clear cached discovery metadata
`edr_queryables()`	`/collections/{id}/queryables`
`edr_instances()` / `edr_instance()`	`/collections/{id}/instances[/{instance}]`
`edr_locations()` / `edr_location()`	`/collections/{id}/locations[/{loc}]`
`edr_location_batch()`	bounded sequential requests to explicit location IDs
`edr_items()` / `edr_item()`	`/collections/{id}/items[/{item}]`
`edr_position()`	`/collections/{id}/position`
`edr_area()`	`/collections/{id}/area`
`edr_cube()`	`/collections/{id}/cube`
`edr_radius()`	`/collections/{id}/radius`
`edr_trajectory()`	`/collections/{id}/trajectory`
`edr_corridor()`	`/collections/{id}/corridor`
`edr_request()`	low-level escape hatch
`covjson_to_tibble()` / `geojson_to_sf()`	response parsers

Every collection query helper also accepts named instance_id =; when set, the path becomes /collections/{id}/instances/{instance_id}/{query}.

What a server actually supports varies. Every query verb above is in the EDR spec and supported by the client, but most servers implement only a subset. On in-situ monitoring deployments, locations, position, cube, and area are common; radius, trajectory, and corridor less so. Hitting a verb the server doesn’t implement gives you an HTTP error. Check the data_queries column from edr_collections() before you assume a query will work.

See vignette("compatibility") for the precise supported subset, return formats, known limitations, and the distinction between verified and merely advertised endpoint behavior.

Common parameters

Every query verb accepts the standard EDR filters:

datetime — an ISO-8601 instant or interval. Accepts "2020-01-01/2020-12-31", an open interval "2020-01-01/..", or a length-2 character vector c("2020-01-01", "2020-12-31").
parameter_name — a character vector of parameter names; sent as a comma-separated parameter-name= query. Use edr_parameters() to discover valid names.
bbox — numeric length-4 (minx, miny, maxx, maxy) or length-6 (with z).
coords — for position/area/radius/trajectory/corridor: a WKT string, a numeric vector / 2-column matrix of lon-lat, or an sf/sfc geometry.
z, crs, limit — passed through when supplied.
f — an exact server-advertised EDR output token. Keep format as the client-side parser selector; for example, pair format = "covjson" with f = "CoverageJSON" when a strict endpoint requires that token.
instance_id — named, optional model-run/version identifier; inserts the standard /instances/{id} path segment before the query type.
... — any extra query parameter is forwarded verbatim.

License

MIT