Add alternative source for ERA5 data to eliminate queueing

.github/workflows/test.yaml

@@ -55,6 +55,7 @@ jobs:
       MPLBACKEND: Agg  # https://github.com/orgs/community/discussions/26434
       CDSAPI_URL: ${{ vars.CDSAPI_URL }}
       CDSAPI_KEY: ${{ secrets.CDSAPI_TOKEN }}
+      EARTHDATAHUB_API_KEY: ${{ secrets.EARTHDATAHUB_API_KEY}}
 
     steps:
     - name: Setup cache and secrets (Linux & MacOS)

RELEASE_NOTES.rst

@@ -24,6 +24,12 @@ Upcoming Release
   buffer that computes a geometrically accurate Euclidean buffer, correct in
   diagonal directions.
 
+* Add an Earth Data Hub (EDH) source for ERA5 as a new dataset module ``era5-edh``,
+  which serves data without the CDS processing queue. Configure a DestinE API key
+  (see :doc:`installation`) and pass ``module="era5-edh"`` when creating a cutout.
+  Only the native 0.25°×0.25° grid is supported; use ``module="era5"`` (CDS) for
+  other resolutions.
+
 **Bug fixes**
 
 * Fix ``Cutout.line_rating`` passing line azimuth in radians while

atlite/datasets/__init__.py

@@ -4,6 +4,11 @@
 
 """atlite datasets."""
 
-from atlite.datasets import era5, gebco, sarah
+from atlite.datasets import era5, era5_edh, gebco, sarah
 
-modules = {"era5": era5, "sarah": sarah, "gebco": gebco}
+modules = {
+    "era5": era5,
+    "era5-edh": era5_edh,
+    "sarah": sarah,
+    "gebco": gebco,
+}

atlite/datasets/era5.py

@@ -23,7 +23,6 @@
 import pandas as pd
 import xarray as xr
 from dask import compute, delayed
-from dask.array import arctan2, sqrt
 from numpy import atleast_1d
 
 from atlite.gis import maybe_swap_spatial_dims
@@ -129,6 +128,40 @@ def _rename_and_clean_coords(ds: xr.Dataset, add_lon_lat: bool = True) -> xr.Dat
     return ds.drop_vars(["expver", "number"], errors="ignore")
 
 
+def _process_wind(ds: xr.Dataset, single_precision: bool = False) -> xr.Dataset:
+    """
+    Derive wind speed, shear exponent, azimuth and roughness from raw components.
+
+    Shared by the CDS (:mod:`atlite.datasets.era5`) and EDH
+    (:mod:`atlite.datasets.era5_edh`) backends. Operates on a dataset carrying
+    the raw ``u10``/``v10``/``u100``/``v100``/``fsr`` variables.
+    ``single_precision`` casts the float64-promoted shear and azimuth back to
+    float32 (used by EDH to halve on-disk size).
+
+    Returns
+    -------
+    xr.Dataset
+        Dataset with variables: wnd100m, wnd_shear_exp, wnd_azimuth, roughness.
+    """
+    for h in (10, 100):
+        units = ds[f"u{h}"].attrs.get("units", "m s**-1")
+        ds[f"wnd{h}m"] = np.sqrt(ds[f"u{h}"] ** 2 + ds[f"v{h}"] ** 2).assign_attrs(
+            units=units, long_name=f"{h} metre wind speed"
+        )
+    shear = (np.log(ds["wnd10m"] / ds["wnd100m"]) / np.log(10 / 100)).assign_attrs(
+        units="", long_name="wind shear exponent"
+    )
+    ds["wnd_shear_exp"] = shear.astype(np.float32) if single_precision else shear
+
+    # span the whole circle: 0 is north, π/2 is east, -π is south, 3π/2 is west
+    azimuth = np.arctan2(ds["u100"], ds["v100"])
+    azimuth = azimuth.where(azimuth >= 0, azimuth + 2 * np.pi)
+    ds["wnd_azimuth"] = azimuth.astype(np.float32) if single_precision else azimuth
+
+    ds = ds.drop_vars(["u100", "v100", "u10", "v10", "wnd10m"])
+    return ds.rename({"fsr": "roughness"})
+
+
 def get_data_wind(retrieval_params: dict[str, Any]) -> xr.Dataset:
     """
     Retrieve and compute wind speed variables from ERA5.
@@ -157,21 +190,7 @@ def get_data_wind(retrieval_params: dict[str, Any]) -> xr.Dataset:
         **retrieval_params,
     )
     ds = _rename_and_clean_coords(ds)
-
-    for h in [10, 100]:
-        ds[f"wnd{h}m"] = sqrt(ds[f"u{h}"] ** 2 + ds[f"v{h}"] ** 2).assign_attrs(
-            units=ds[f"u{h}"].attrs["units"], long_name=f"{h} metre wind speed"
-        )
-    ds["wnd_shear_exp"] = (
-        np.log(ds["wnd10m"] / ds["wnd100m"]) / np.log(10 / 100)
-    ).assign_attrs(units="", long_name="wind shear exponent")
-
-    # span the whole circle: 0 is north, π/2 is east, -π is south, 3π/2 is west
-    azimuth = arctan2(ds["u100"], ds["v100"])
-    ds["wnd_azimuth"] = azimuth.where(azimuth >= 0, azimuth + 2 * np.pi)
-
-    ds = ds.drop_vars(["u100", "v100", "u10", "v10", "wnd10m"])
-    return ds.rename({"fsr": "roughness"})
+    return _process_wind(ds)
 
 
 def sanitize_wind(ds: xr.Dataset) -> xr.Dataset:
@@ -192,36 +211,22 @@ def sanitize_wind(ds: xr.Dataset) -> xr.Dataset:
     return ds
 
 
-def get_data_influx(retrieval_params: dict[str, Any]) -> xr.Dataset:
+def _process_influx(ds: xr.Dataset, single_precision: bool = False) -> xr.Dataset:
     """
-    Retrieve and compute solar radiation variables from ERA5.
+    Derive influx variables and solar position from raw radiation fields.
 
-    Downloads radiation components, converts from J/m² to W/m², computes
-    albedo, diffuse radiation, and solar position (altitude/azimuth).
-
-    Parameters
-    ----------
-    retrieval_params : dict[str, Any]
-        CDS API retrieval parameters including area, time, and format.
+    Shared by the CDS (:mod:`atlite.datasets.era5`) and EDH
+    (:mod:`atlite.datasets.era5_edh`) backends. Operates on a dataset carrying
+    the raw ``ssrd``/``ssr``/``fdir``/``tisr`` variables. ``single_precision``
+    casts the solar-position fields back to float32 (used by EDH to halve
+    on-disk size).
 
     Returns
     -------
     xr.Dataset
         Dataset with variables: influx_toa, influx_direct, influx_diffuse,
         albedo, solar_altitude, solar_azimuth.
     """
-    ds = retrieve_data(
-        variable=[
-            "surface_net_solar_radiation",
-            "surface_solar_radiation_downwards",
-            "toa_incident_solar_radiation",
-            "total_sky_direct_solar_radiation_at_surface",
-        ],
-        **retrieval_params,
-    )
-
-    ds = _rename_and_clean_coords(ds)
-
     ds = ds.rename({"fdir": "influx_direct", "tisr": "influx_toa"})
     ds["albedo"] = (
         ((ds["ssrd"] - ds["ssr"]) / ds["ssrd"].where(ds["ssrd"] != 0))
@@ -245,13 +250,46 @@ def get_data_influx(retrieval_params: dict[str, Any]) -> xr.Dataset:
     # Suppress DeprecationWarning from new SolarPosition calculation (#199)
     with warnings.catch_warnings():
         warnings.simplefilter("ignore", DeprecationWarning)
-        time_shift = pd.to_timedelta("-30 minutes")
-        sp = SolarPosition(ds, time_shift=time_shift)
+        sp = SolarPosition(ds, time_shift=pd.to_timedelta("-30 minutes"))
     sp = sp.rename({v: f"solar_{v}" for v in sp.data_vars})
+    if single_precision:
+        sp = sp.astype(np.float32)
 
     return xr.merge([ds, sp])
 
 
+def get_data_influx(retrieval_params: dict[str, Any]) -> xr.Dataset:
+    """
+    Retrieve and compute solar radiation variables from ERA5.
+
+    Downloads radiation components, converts from J/m² to W/m², computes
+    albedo, diffuse radiation, and solar position (altitude/azimuth).
+
+    Parameters
+    ----------
+    retrieval_params : dict[str, Any]
+        CDS API retrieval parameters including area, time, and format.
+
+    Returns
+    -------
+    xr.Dataset
+        Dataset with variables: influx_toa, influx_direct, influx_diffuse,
+        albedo, solar_altitude, solar_azimuth.
+    """
+    ds = retrieve_data(
+        variable=[
+            "surface_net_solar_radiation",
+            "surface_solar_radiation_downwards",
+            "toa_incident_solar_radiation",
+            "total_sky_direct_solar_radiation_at_surface",
+        ],
+        **retrieval_params,
+    )
+
+    ds = _rename_and_clean_coords(ds)
+    return _process_influx(ds)
+
+
 def sanitize_influx(ds: xr.Dataset) -> xr.Dataset:
     """
     Clip negative radiation values to zero.