Era5Prepare

Functions:

• era5_tiler –

  Splits processed ERA5 NetCDF files into the tiled directory structure expected

• preprocess_era5_data –

  Ssrd and fdir are hourly backward accumulated quantities in ERA5 with the unit: J m⁻².

era5_tiler

era5_tiler(
    processed_dir: str,
    tile_output_dir: str,
    zoom_level: int = 4,
    raw_nc: Optional[str] = None,
    raw_variables: Optional[List[str]] = None,
) -> str

Splits processed ERA5 NetCDF files into the tiled directory structure expected by Era5Source and execute_workflow_iteratively().

Output follows the shared-data naming convention: ///// reanalysis-era5-single-levels.z.x.y.y..nc where is taken from _ERA5_NC_TO_TILE_LABEL (e.g. "100m_wind_speed.processed"). One output file is written per variable per tile per year.

Parameters:

• processed_dir

  (str) –

  Directory containing processed NC files (output of preprocess_era5_data).

• tile_output_dir

  (str) –

  Root directory for the tiled output.

• zoom_level

  (int, default: 4 ) –

  Web Mercator zoom level (default: 4 → 16×16 global grid).

• raw_nc

  (str, default: None ) –

  Path to the raw ERA5 download file, used to tile raw_variables.

• raw_variables

  (list of str, default: None ) –

  NC short names to tile from raw_nc (e.g. ['t2m', 'sp', 'blh']). Ignored if raw_nc is None.

Source code in reskit/weather/Era5Source/Era5Prepare.py

def era5_tiler(
    processed_dir: str,
    tile_output_dir: str,
    zoom_level: int = 4,
    raw_nc: Optional[str] = None,
    raw_variables: Optional[List[str]] = None,
) -> str:
    """
    Splits processed ERA5 NetCDF files into the tiled directory structure expected
    by Era5Source and execute_workflow_iteratively().

    Output follows the shared-data naming convention:
        <tile_output_dir>/<zoom>/<xi>/<yi>/<year>/
            reanalysis-era5-single-levels.z<z>.x<xi>.y<yi>.y<year>.<label>.nc
    where <label> is taken from _ERA5_NC_TO_TILE_LABEL (e.g. "100m_wind_speed.processed").
    One output file is written per variable per tile per year.

    Parameters
    ----------
    processed_dir : str
        Directory containing processed NC files (output of preprocess_era5_data).
    tile_output_dir : str
        Root directory for the tiled output.
    zoom_level : int
        Web Mercator zoom level (default: 4 → 16×16 global grid).
    raw_nc : str, optional
        Path to the raw ERA5 download file, used to tile raw_variables.
    raw_variables : list of str, optional
        NC short names to tile from raw_nc (e.g. ['t2m', 'sp', 'blh']).
        Ignored if raw_nc is None.
    """
    source_group = "reanalysis-era5-single-levels"

    # (source_file, [nc_var_names_to_tile]) pairs
    file_var_pairs: list[tuple[str, list[str]]] = []

    for f in sorted(os.listdir(processed_dir)):
        if not f.endswith(".nc") or "_processed_" not in f:
            continue
        path = os.path.join(processed_dir, f)
        vars_in_file = _nc_data_var_names(path)
        vars_to_tile = [v for v in vars_in_file if v in _ERA5_NC_TO_TILE_LABEL]
        if vars_to_tile:
            file_var_pairs.append((path, vars_to_tile))

    if raw_nc and raw_variables:
        vars_in_raw = set(_nc_data_var_names(raw_nc))
        vars_to_tile = [v for v in raw_variables if v in vars_in_raw and v in _ERA5_NC_TO_TILE_LABEL]
        if vars_to_tile:
            file_var_pairs.append((raw_nc, vars_to_tile))

    for source_file, variables in file_var_pairs:
        with nc4.Dataset(source_file) as ds:
            lats = ds.variables["latitude"][:]
            lons = ds.variables["longitude"][:]
        lon_min, lon_max = float(lons.min()), float(lons.max())
        lat_min, lat_max = float(lats.min()), float(lats.max())
        years = _nc_years(source_file)

        # NW corner → SE corner (tile Y increases southward)
        xi_values = _iter_tile_x_indices(zoom_level=zoom_level, lon_west=lon_min, lon_east=lon_max)
        _, yi_nw = get_tile_XY(zoom_level, lon=_tile_lookup_lon(lon_min), lat=lat_max)
        _, yi_se = get_tile_XY(zoom_level, lon=_tile_lookup_lon(lon_max), lat=lat_min)

        for xi in xi_values:
            for yi in range(yi_nw, yi_se + 1):
                extent = gk.Extent.fromTile(xi, yi, zoom_level).castTo(gk.srs.EPSG4326).pad(2)
                lon_west, lon_east, lat_south, lat_north = extent.xXyY
                lon_boxes = _get_source_lon_boxes(
                    lon_west=lon_west,
                    lon_east=lon_east,
                    source_lon_min=lon_min,
                    source_lon_max=lon_max,
                )

                for year in years:
                    target_dir = os.path.join(tile_output_dir, str(zoom_level), str(xi), str(yi), str(year))
                    os.makedirs(target_dir, exist_ok=True)

                    for var in variables:
                        label = _ERA5_NC_TO_TILE_LABEL[var]
                        target_file = os.path.join(
                            target_dir,
                            f"{source_group}.z{zoom_level}.x{xi}.y{yi}.y{year}.{label}.nc",
                        )
                        if os.path.exists(target_file):
                            print(f"Skipping tile (exists): {target_file}")
                            continue

                        _tile_variable_to_file(
                            source_file=source_file,
                            var=var,
                            year=year,
                            lat_south=lat_south,
                            lat_north=lat_north,
                            lon_boxes=lon_boxes,
                            target_file=target_file,
                            source_lon_min=lon_min,
                            source_lon_max=lon_max,
                        )

    return tile_output_dir

preprocess_era5_data

preprocess_era5_data(
    focus_nc: str, processed_dir: Optional[str] = None
)

Ssrd and fdir are hourly backward accumulated quantities in ERA5 with the unit: J m⁻². Each value at time t represents the accumulated energy over the previous hour.

When you divide it by 3600 seconds, you convert the accumulated energy (J m⁻²) into an average power flux (W m⁻²) over that hour.

However, in solar observation and other models, this value represents the instantaneous mean over the next hour. This matches: PV modeling conventions Many energy system models atlite / PyPSA conventions

So, we need to do two things: 1. Convert the accumulated quantity to an average power flux by dividing by 3600 2. Shift the time axis forward by one hour to represent the average over the next hour.

Parameters:

• focus_nc

  (str) –

  Path to the raw ERA5 NetCDF file.

• processed_dir

  (str, default: None ) –

  Directory to write processed output files. Defaults to the same directory as focus_nc.

Source code in reskit/weather/Era5Source/Era5Prepare.py

def preprocess_era5_data(focus_nc: str, processed_dir: Optional[str] = None):
    """
    Ssrd and fdir are hourly backward accumulated quantities in ERA5 with the unit: J m⁻².
    Each value at time t represents the accumulated energy over the previous hour.

    When you divide it by 3600 seconds, you convert the accumulated energy (J m⁻²)
    into an average power flux (W m⁻²) over that hour.

    However, in solar observation and other models,
    this value represents the instantaneous mean over the next hour.
    This matches:
        PV modeling conventions
        Many energy system models
        atlite / PyPSA conventions

    So, we need to do two things:
    1. Convert the accumulated quantity to an average power flux by dividing by 3600
    2. Shift the time axis forward by one hour to represent the average over the next hour.

    Parameters
    ----------
    focus_nc : str
        Path to the raw ERA5 NetCDF file.
    processed_dir : str, optional
        Directory to write processed output files. Defaults to the same directory
        as focus_nc.
    """
    out_dir = processed_dir or os.path.dirname(focus_nc)
    if processed_dir:
        os.makedirs(processed_dir, exist_ok=True)
    f_name = os.path.basename(focus_nc)

    with _open_era5_dataset(focus_nc) as ds:
        varset = set(ds.data_vars)

        # process for solar radiation variables (time adjusted)
        # Build the variable list dynamically from whichever of ssrd/fdir is present so that
        # workflows requesting only one of them (e.g. CSP needs fdir but not ssrd) still work.
        solar_t_out = os.path.join(out_dir, f"{f_name.split('.')[0]}_processed_solar_t_adjusted.nc")
        solar_vars = [v for v in ("ssrd", "fdir") if v in varset]
        if solar_vars:
            out_names = [f"{v}_t_adj" for v in solar_vars]
            if not _nc_file_has_vars(solar_t_out, out_names):
                # 1. accumulated J m**-2 -> mean power flux W m**-2 (divide by 3600s)
                # 2. shift time axis +1h so each value is the mean over the *next* hour
                time_encoding = dict(ds["time"].encoding)
                out = ds[solar_vars] / 3600.0
                out = out.assign_coords(time=out["time"] + pd.Timedelta(hours=1))
                out["time"].encoding = time_encoding
                out = out.rename({v: f"{v}_t_adj" for v in solar_vars})
                for v in solar_vars:
                    attrs = dict(ds[v].attrs)
                    attrs["units"] = "W m**-2"
                    out[f"{v}_t_adj"].attrs = attrs
                _write_netcdf(out.load(), solar_t_out)
            else:
                print(f"Skipping process time-adjusted solar (exists): {solar_t_out}")

        # process for wind speed variables
        ws100_out = os.path.join(out_dir, f"{f_name.split('.')[0]}_processed_ws100.nc")
        if {"u100", "v100"} <= varset:
            if not _nc_file_has_vars(ws100_out, ["ws100"]):
                ws100 = np.sqrt(ds["u100"] ** 2 + ds["v100"] ** 2)
                ws100.name = "ws100"
                ws100.attrs = {"long_name": "100 metre wind speed", "units": "m s**-1"}
                _write_netcdf(ws100.to_dataset().load(), ws100_out)
            else:
                print(f"Skipping process ws100 (exists): {ws100_out}")

        ws10_out = os.path.join(out_dir, f"{f_name.split('.')[0]}_processed_ws10.nc")
        if {"u10", "v10"} <= varset:
            if not _nc_file_has_vars(ws10_out, ["ws10"]):
                ws10 = np.sqrt(ds["u10"] ** 2 + ds["v10"] ** 2)
                ws10.name = "ws10"
                ws10.attrs = {"long_name": "10 metre wind speed", "units": "m s**-1"}
                _write_netcdf(ws10.to_dataset().load(), ws10_out)
            else:
                print(f"Skipping process ws10 (exists): {ws10_out}")