extract_era5.R

# Load necessary library
library(ecmwfr)
library(ncdf4)

# Define query parameters
latitude <- 48.666848
longitude <- 9.834224
years <- 2020:2022
months <- 1:12
days <- 1:31
hours <- sprintf("%02d:00", 0:23)  # Every hour of the day

# Define mandatory and optional variables
variables <- c("cloud_base_height", "low_cloud_cover", "total_cloud_cover")
optional_variables <- c("total_precipitation", "convective_precipitation", 
                        "convective_rain_rate", "precipitation_type")

# Set up the request
request <- list(
  dataset_short_name = "reanalysis-era5-single-levels",
  product_type = "reanalysis",
  variable = variables,
  year = years,
  month = months,
  day = days,
  time = hours,
  # Define the bounding box area (North, West, South, East)
  area = c(latitude + .001, longitude - .001, latitude - .001, longitude + .001),
  format = "netcdf",
  target = "test.nc"
)

# Execute the request (you'll need to have your API key set up)
wf_request(request)

# Url displayed by the previous line code
wf_transfer(
  url = 'https://cds.climate.copernicus.eu/api/retrieve/v1/jobs/a2ba87f5-adfc-4e79-b807-7969f77f5c49',
  path = "~/Downloads/",
  filename = "test.nc"
)

# Save data to file
nc <- nc_open("~/Downloads/test.nc")

# Extract time while converting to the correct format
time <- as.POSIXct(ncvar_get(nc, "valid_time"), origin = "1970-01-01", tz = "UTC")

# Extract variable
lcc <- ncvar_get(nc, "lcc")
tcc <- ncvar_get(nc, "tcc")
cbh <- ncvar_get(nc, "cbh")

# Plot variable for instance
plot(time, lcc)

# Rest of your analysis