Skip to content

Commit b5ea302

Browse files
chross22claude
andcommitted
Request Copernicus variables by name, with a printable dictionary
Copernicus variable codes are terse and easy to misremember - thetao for temperature, mlotst for mixed layer depth, zos for sea surface height - and getting one wrong produces a failed download rather than an obvious mistake. Variables can now be requested by name instead. - copernicus_variables(): the catalog, 14 physical and biogeochemical variables mapped to their product, dataset, and code. - variable_dictionary(): the catalog as a data frame with a print method, so `variable_dictionary()` at the console shows what is available. Filterable by product; full descriptions stay in the returned object rather than wrapping unreadably in the printed view. - variable_dataset() and infer_dataset(): which dataset a variable lives in, and the product/dataset implied by a set of names. accessEnvDat() accepts either names or raw codes, so existing calls keep working unchanged, and result columns take the names the caller asked for - requesting "SST" yields a column called SST rather than thetao. product_id and dataset_id may now be omitted when every variable is in the catalog. Two failure modes are now caught rather than producing quietly wrong output: - Variables from different datasets cannot be fetched in one request. Mixing them is refused before any download is attempted, naming which variable belongs to which dataset, instead of failing obscurely at the API. - Column names are assigned positionally, so a raster with more layers than requested variables - which a depth range spanning several model levels produces - would silently mislabel every column. That now errors. Anything outside the dictionary is passed through to the API with a warning, since Copernicus serves far more than this catalog covers, but a typo looks identical to a real code and should not pass silently. README.Rmd gains real documentation in place of the scaffold placeholders. Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
1 parent 2a9a082 commit b5ea302

13 files changed

Lines changed: 773 additions & 34 deletions

NAMESPACE

Lines changed: 4 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -1,4 +1,8 @@
11
# Generated by roxygen2: do not edit by hand
22

3+
S3method(print,datamatch_dictionary)
34
export(accessEnvDat)
5+
export(copernicus_variables)
46
export(matchData)
7+
export(variable_dataset)
8+
export(variable_dictionary)

R/accessEnvDat.R

Lines changed: 32 additions & 5 deletions
Original file line numberDiff line numberDiff line change
@@ -22,10 +22,25 @@
2222
#' network/API, not CPU.
2323
#' @return envDat <sf object> sf object containing requested environmental data from Copernicus Marine Service
2424
#' @export
25-
accessEnvDat <- function(product_id, dataset_id, vars, years, months,
25+
accessEnvDat <- function(product_id = NULL, dataset_id = NULL, vars, years, months,
2626
bounding_box, depth = c(0,1),
2727
overwrite = FALSE, n_workers = 1) {
2828

29+
# `vars` may be catalog names ("SST") or raw Copernicus codes ("thetao").
30+
# Codes go to the API; names come back as the column names, so a caller who
31+
# asked for SST gets a column called SST rather than thetao.
32+
resolved <- resolve_variables(vars)
33+
var_codes <- resolved$codes
34+
var_names <- resolved$names
35+
36+
# With every variable in the catalog, the product and dataset are implied and
37+
# need not be repeated at the call site.
38+
if (is.null(product_id) || is.null(dataset_id)) {
39+
inferred <- infer_dataset(vars)
40+
product_id <- product_id %||% inferred$product_id
41+
dataset_id <- dataset_id %||% inferred$dataset_id
42+
}
43+
2944
# Build the full list of (year, month, day) combinations to fetch up front,
3045
# so they can be dispatched in parallel instead of three nested serial loops.
3146
is_daily <- substr(dataset_id, nchar(dataset_id) - 2, nchar(dataset_id) - 2) == "D"
@@ -75,20 +90,32 @@ accessEnvDat <- function(product_id, dataset_id, vars, years, months,
7590
})
7691
results <- tryCatch(
7792
parallel::parLapply(cl, work_items, fetch_one_day,
78-
product_id = product_id, dataset_id = dataset_id, vars = vars,
93+
product_id = product_id, dataset_id = dataset_id, vars = var_codes,
7994
bounding_box = bounding_box, depth = depth, overwrite = overwrite),
8095
error = function(e) stop("accessEnvDat: parallel fetch failed - ", conditionMessage(e), call. = FALSE)
8196
)
8297
} else {
8398
results <- lapply(work_items, fetch_one_day,
84-
product_id = product_id, dataset_id = dataset_id, vars = vars,
99+
product_id = product_id, dataset_id = dataset_id, vars = var_codes,
85100
bounding_box = bounding_box, depth = depth, overwrite = overwrite)
86101
}
87102

88103
covars <- dplyr::bind_rows(results)
89104

90-
# Define column names
91-
names(covars) <- c("x", "y", vars, "YEAR", "MONTH", "DAY")
105+
# Names are assigned positionally, so the raster must have exactly one layer
106+
# per requested variable. A depth range spanning several model levels returns
107+
# more, and silently mislabelling those columns would be worse than stopping.
108+
expected <- length(var_names) + 5 # x, y, vars..., YEAR, MONTH, DAY
109+
if (ncol(covars) != expected) {
110+
stop("Expected ", length(var_names), " variable column(s) but the download ",
111+
"returned ", ncol(covars) - 5, ". This usually means the depth range ",
112+
"spans several model levels; request a single level, or one variable ",
113+
"at a time.", call. = FALSE)
114+
}
115+
116+
# Columns take the names the caller asked for, so requesting "SST" yields a
117+
# column called SST rather than thetao.
118+
names(covars) <- c("x", "y", var_names, "YEAR", "MONTH", "DAY")
92119

93120
# Convert data to sf and return
94121
sf::st_as_sf(covars,

R/variables.R

Lines changed: 232 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,232 @@
1+
# Null-coalescing operator. Defined here rather than relying on base R's, which
2+
# only exists from R 4.4 onward.
3+
`%||%` <- function(x, y) if (is.null(x)) y else x
4+
5+
#' Catalog of Copernicus variables under familiar names
6+
#'
7+
#' Maps short names people actually use (`SST`, `CHL`, ...) onto the Copernicus
8+
#' Marine product, dataset, and variable code that supply them. Copernicus codes
9+
#' are terse and easy to misremember — `thetao` for temperature, `mlotst` for
10+
#' mixed layer depth, `zos` for sea surface height — and getting one wrong
11+
#' produces a failed download rather than an obvious mistake.
12+
#'
13+
#' All entries are monthly means from the global reanalyses: physical variables
14+
#' from `GLOBAL_MULTIYEAR_PHY_001_030` (GLORYS12V1) and biogeochemical ones from
15+
#' `GLOBAL_MULTIYEAR_BGC_001_029`.
16+
#'
17+
#' Copernicus revises dataset identifiers periodically. If a fetch fails with an
18+
#' unknown-dataset error, check the current identifier on the Copernicus Marine
19+
#' Data Store and pass `dataset_id` explicitly.
20+
#'
21+
#' @return a named list, one entry per variable, each with `variable`, `label`,
22+
#' `units`, `product_id`, `dataset_id`, and `description`
23+
#' @examples
24+
#' names(copernicus_variables())
25+
#' copernicus_variables()$SST$variable
26+
#' @seealso [variable_dictionary()] for a printable table
27+
#' @export
28+
copernicus_variables <- function() {
29+
phy_product <- "GLOBAL_MULTIYEAR_PHY_001_030"
30+
phy_dataset <- "cmems_mod_glo_phy_my_0.083deg_P1M-m"
31+
bgc_product <- "GLOBAL_MULTIYEAR_BGC_001_029"
32+
bgc_dataset <- "cmems_mod_glo_bgc_my_0.25deg_P1M-m"
33+
34+
physical <- function(variable, label, units, description) {
35+
list(variable = variable, label = label, units = units,
36+
product_id = phy_product, dataset_id = phy_dataset,
37+
description = description)
38+
}
39+
biogeochemical <- function(variable, label, units, description) {
40+
list(variable = variable, label = label, units = units,
41+
product_id = bgc_product, dataset_id = bgc_dataset,
42+
description = description)
43+
}
44+
45+
list(
46+
SST = physical("thetao", "Sea surface temperature", "degrees C",
47+
"Sea water potential temperature at the surface."),
48+
SSS = physical("so", "Sea surface salinity", "PSU",
49+
"Sea water salinity at the surface."),
50+
BOTT = physical("bottomT", "Bottom temperature", "degrees C",
51+
paste("Sea water potential temperature at the sea floor.",
52+
"Relevant to overwintering copepod stages.")),
53+
UO = physical("uo", "Eastward current velocity", "m/s",
54+
"Eastward component of sea water velocity."),
55+
VO = physical("vo", "Northward current velocity", "m/s",
56+
"Northward component of sea water velocity."),
57+
SSH = physical("zos", "Sea surface height", "m",
58+
paste("Sea surface height above geoid. A proxy for",
59+
"mesoscale circulation features.")),
60+
MLD = physical("mlotst", "Mixed layer depth", "m",
61+
paste("Ocean mixed layer thickness by a sigma-theta",
62+
"criterion. Controls how deeply plankton are mixed.")),
63+
SIC = physical("siconc", "Sea ice concentration", "fraction",
64+
"Fraction of the cell covered by sea ice."),
65+
CHL = biogeochemical("chl", "Chlorophyll-a concentration", "mg/m3",
66+
paste("Mass concentration of chlorophyll-a. A food",
67+
"availability proxy; usually worth",
68+
"log-transforming.")),
69+
NO3 = biogeochemical("no3", "Nitrate concentration", "mmol/m3",
70+
"Mole concentration of nitrate, a limiting nutrient."),
71+
PO4 = biogeochemical("po4", "Phosphate concentration", "mmol/m3",
72+
"Mole concentration of phosphate."),
73+
O2 = biogeochemical("o2", "Dissolved oxygen", "mmol/m3",
74+
"Mole concentration of dissolved molecular oxygen."),
75+
NPP = biogeochemical("nppv", "Net primary production", "mg/m3/day",
76+
paste("Net primary production of biomass expressed as",
77+
"carbon. A direct productivity measure rather",
78+
"than the standing stock chlorophyll reports.")),
79+
PH = biogeochemical("ph", "pH", "1", "Sea water pH reported on total scale.")
80+
)
81+
}
82+
83+
#' Printable dictionary of variable names
84+
#'
85+
#' The catalog as a data frame: what each short name means, its units, and the
86+
#' Copernicus code and dataset behind it. Print it to see what is available
87+
#' without leaving the console.
88+
#'
89+
#' @param product filter to `"physical"`, `"biogeochemical"`, or `"all"`
90+
#' @return a data frame of class `datamatch_dictionary` with columns `name`,
91+
#' `variable`, `label`, `units`, and `dataset`
92+
#' @examples
93+
#' variable_dictionary()
94+
#' variable_dictionary("biogeochemical")
95+
#'
96+
#' # As a plain data frame, for programmatic use
97+
#' as.data.frame(variable_dictionary())
98+
#' @export
99+
variable_dictionary <- function(product = c("all", "physical", "biogeochemical")) {
100+
product <- match.arg(product)
101+
catalog <- copernicus_variables()
102+
103+
dictionary <- do.call(rbind, lapply(names(catalog), function(name) {
104+
entry <- catalog[[name]]
105+
data.frame(
106+
name = name, variable = entry$variable, label = entry$label,
107+
units = entry$units, dataset = entry$dataset_id,
108+
description = entry$description, stringsAsFactors = FALSE
109+
)
110+
}))
111+
112+
if (product != "all") {
113+
is_physical <- grepl("_phy_", dictionary$dataset, fixed = TRUE)
114+
dictionary <- dictionary[if (product == "physical") is_physical else !is_physical, ]
115+
rownames(dictionary) <- NULL
116+
}
117+
118+
class(dictionary) <- c("datamatch_dictionary", "data.frame")
119+
dictionary
120+
}
121+
122+
#' @param x a `datamatch_dictionary`
123+
#' @param ... ignored
124+
#' @rdname variable_dictionary
125+
#' @export
126+
print.datamatch_dictionary <- function(x, ...) {
127+
cat("Copernicus variables available by name\n")
128+
cat(strrep("-", 62), "\n", sep = "")
129+
130+
# The description is the widest column by far and would wrap unreadably, so
131+
# the printed view drops it; it remains in the returned object.
132+
visible <- as.data.frame(x)[c("name", "variable", "label", "units")]
133+
print(visible, row.names = FALSE, right = FALSE)
134+
135+
cat("\nPass a name to accessEnvDat(vars = ...), or the Copernicus code.\n")
136+
cat("Full descriptions: as.data.frame(variable_dictionary())$description\n")
137+
invisible(x)
138+
}
139+
140+
#' Resolve variable names to Copernicus codes
141+
#'
142+
#' Accepts either a catalog name (`"SST"`) or a raw Copernicus code
143+
#' (`"thetao"`), so existing calls that pass codes keep working unchanged.
144+
#'
145+
#' @param vars variable names or codes
146+
#' @return a list with `codes` (Copernicus codes, in the given order) and
147+
#' `names` (what each should be called in the result)
148+
#' @keywords internal
149+
resolve_variables <- function(vars) {
150+
catalog <- copernicus_variables()
151+
known_codes <- vapply(catalog, function(entry) entry$variable, character(1))
152+
153+
codes <- character(length(vars))
154+
for (i in seq_along(vars)) {
155+
if (vars[i] %in% names(catalog)) {
156+
codes[i] <- catalog[[vars[i]]]$variable
157+
} else if (vars[i] %in% known_codes) {
158+
# Already a Copernicus code.
159+
codes[i] <- vars[i]
160+
} else {
161+
# Not in the catalog at all. Copernicus serves far more than this catalog
162+
# covers, so an unrecognized string is passed through as a code rather
163+
# than rejected - but say so, since a typo looks identical.
164+
warning("'", vars[i], "' is not in the variable dictionary; passing it to ",
165+
"Copernicus as a variable code. See variable_dictionary() for ",
166+
"known names.", call. = FALSE)
167+
codes[i] <- vars[i]
168+
}
169+
}
170+
171+
list(codes = unname(codes), names = vars)
172+
}
173+
174+
#' Infer the product and dataset from a set of variable names
175+
#'
176+
#' When every requested variable is in the catalog, the product and dataset are
177+
#' implied, so a call need not repeat identifiers that are already known.
178+
#'
179+
#' Variables from different datasets cannot be fetched in one request, so mixing
180+
#' them is an error here rather than a confusing failure at the API.
181+
#'
182+
#' @param vars variable names
183+
#' @return `list(product_id =, dataset_id =)`
184+
#' @keywords internal
185+
infer_dataset <- function(vars) {
186+
datasets <- variable_dataset(vars)
187+
188+
unknown <- vars[is.na(datasets)]
189+
if (length(unknown) > 0) {
190+
stop("Cannot infer the dataset for: ", paste(unknown, collapse = ", "),
191+
"\nEither use a name from variable_dictionary(), or pass product_id ",
192+
"and dataset_id explicitly.", call. = FALSE)
193+
}
194+
195+
distinct <- unique(unname(datasets))
196+
if (length(distinct) > 1) {
197+
catalog <- copernicus_variables()
198+
grouped <- vapply(distinct, function(d) {
199+
paste(vars[datasets == d], collapse = ", ")
200+
}, character(1))
201+
stop("These variables come from different Copernicus datasets and cannot ",
202+
"be fetched together:\n ",
203+
paste(paste0(grouped, " -> ", distinct), collapse = "\n "),
204+
"\nCall accessEnvDat() once per dataset.", call. = FALSE)
205+
}
206+
207+
catalog <- copernicus_variables()
208+
entry <- catalog[[vars[1]]]
209+
list(product_id = entry$product_id, dataset_id = entry$dataset_id)
210+
}
211+
212+
#' Look up the dataset a set of variables comes from
213+
#'
214+
#' Variables in one Copernicus dataset can be fetched together; variables from
215+
#' different datasets cannot. This reports which dataset each name belongs to so
216+
#' a caller can group them.
217+
#'
218+
#' @param vars variable names from the catalog
219+
#' @return a named character vector of dataset identifiers, `NA` for names not in
220+
#' the catalog
221+
#' @examples
222+
#' variable_dataset(c("SST", "SSS", "CHL"))
223+
#' @export
224+
variable_dataset <- function(vars) {
225+
catalog <- copernicus_variables()
226+
stats::setNames(
227+
vapply(vars, function(v) {
228+
if (v %in% names(catalog)) catalog[[v]]$dataset_id else NA_character_
229+
}, character(1)),
230+
vars
231+
)
232+
}

README.Rmd

Lines changed: 63 additions & 8 deletions
Original file line numberDiff line numberDiff line change
@@ -36,20 +36,75 @@ Before installation, you must install the Copernicus Marine toolbox from Coperni
3636

3737
To install this toolbox, visit the Copernicus website. In brief, ...
3838

39-
## Example
39+
## Variable names
4040

41-
This is a basic example which shows you how to solve a common problem:
41+
Copernicus variable codes are terse and easy to misremember — `thetao` for
42+
temperature, `mlotst` for mixed layer depth, `zos` for sea surface height — and
43+
getting one wrong produces a failed download rather than an obvious mistake. So
44+
variables can be requested by name instead, and the dictionary lists what is
45+
available:
4246

43-
```{r example}
47+
```{r dictionary}
4448
library(datamatch)
45-
## basic example code
49+
50+
variable_dictionary()
51+
```
52+
53+
Pass those names to `accessEnvDat()` and the result comes back with them as
54+
column names, rather than the Copernicus codes:
55+
56+
```{r named-fetch, eval = FALSE}
57+
env <- accessEnvDat(
58+
vars = c("SST", "SSS", "MLD"),
59+
years = 2003:2017,
60+
months = 1:12,
61+
bounding_box = list(xmin = -76, xmax = -65, ymin = 35, ymax = 45)
62+
)
63+
names(env)
64+
#> "SST" "SSS" "MLD" "YEAR" "MONTH" "DAY" "geometry"
65+
```
66+
67+
`product_id` and `dataset_id` can be omitted when every variable is in the
68+
dictionary, since the catalog already knows where they live. Variables from
69+
different datasets cannot be fetched in one request, and mixing them is refused
70+
before any download is attempted rather than failing obscurely at the API:
71+
72+
```{r mixed, eval = FALSE}
73+
accessEnvDat(vars = c("SST", "CHL"), ...)
74+
#> Error: These variables come from different Copernicus datasets and cannot be
75+
#> fetched together:
76+
#> SST -> cmems_mod_glo_phy_my_0.083deg_P1M-m
77+
#> CHL -> cmems_mod_glo_bgc_my_0.25deg_P1M-m
78+
#> Call accessEnvDat() once per dataset.
79+
```
80+
81+
Raw Copernicus codes still work exactly as before, so existing calls need no
82+
change, and anything outside the dictionary is passed through to the API with a
83+
warning — Copernicus serves far more than this catalog covers.
84+
85+
```{r lookup, eval = FALSE}
86+
variable_dictionary("biogeochemical") # filter by product
87+
variable_dataset(c("SST", "CHL")) # which dataset each comes from
88+
as.data.frame(variable_dictionary())$description # full descriptions
4689
```
4790

48-
What is special about using `README.Rmd` instead of just `README.md`? You can include R chunks like so:
91+
## Matching to observations
92+
93+
`matchData()` joins environmental data to species observations at the
94+
environmental data's own temporal resolution, inferred from its time steps. This
95+
matters for monthly products: a monthly mean carries one time step per month
96+
while observations fall on arbitrary days, so matching on exact dates would
97+
match nothing.
4998

50-
```{r cars}
51-
# ADD CODE HERE
99+
```{r match, eval = FALSE}
100+
matched <- matchData(speciesDat = observations, envDat = env)
52101
```
53102

54-
You'll still need to render `README.Rmd` regularly, to keep `README.md` up-to-date. `devtools::build_readme()` is handy for this.
103+
Observations falling in a period with no environmental data are returned with
104+
`NA` values and a warning naming the periods, rather than being dropped silently.
105+
106+
## Related packages
55107

108+
- [derivoce](https://github.com/chross22/derivoce) — derived covariates
109+
(gradients, FTLE/FSLE, front and isobath distances, lags, integrals) computed
110+
from what `accessEnvDat()` returns

0 commit comments

Comments
 (0)