Heat waves & public health — apparent temperature¶

Heat-related health risk depends on more than air temperature: humidity amplifies the physiological stress because evaporative cooling becomes less effective. The heat index (apparent temperature) combines air temperature and dewpoint into a single number that maps to discomfort categories.

Domain context. A common formula (Rothfusz, NOAA NWS) is

$$ HI = -42.379 + 2.04901523 T + 10.14333127 RH - 0.22475541 \, T \cdot RH - \dots $$

where $T$ is air temperature (°F) and $RH$ is relative humidity (%). Below $T \approx 80$°F the formula returns roughly $T$ itself; above it, humidity drives the index up sharply. A simpler proxy used here: the wet-bulb temperature, derivable from $T_{2m}$ and dewpoint $T_d$ via Davies-Jones / Stull approximations.

Step 1 — pull a summer's daily data¶

Box: 1° around Athens, Greece (37°–38°N, 23°–24°E) — well-known European heat-wave hotspot. Range June–August 2022 covers the infamous summer heat wave. Daily resolution captures the day-to-day variability that monthly aggregation would smooth out.

In [1]:

Copied!





from pathlib import Path
from earthlens import EarthLens, AggregationConfig

OUT = Path("data/era5-athens-summer")
OUT.mkdir(parents=True, exist_ok=True)

earthlens = EarthLens(
    data_source="ecmwf",
    temporal_resolution="daily",
    start="2022-06-01",
    end="2022-08-31",
    variables={
        "reanalysis-era5-single-levels": [
            "2m-temperature",
            "2m-dewpoint-temperature",
        ],
    },
    lat_lim=[37.0, 38.0],
    lon_lim=[23.0, 24.0],
    path=str(OUT),
)
earthlens.download(aggregate=AggregationConfig(freq="1D", op="mean"))
from pathlib import Path
from earthlens import EarthLens, AggregationConfig

OUT = Path("data/era5-athens-summer")
OUT.mkdir(parents=True, exist_ok=True)

earthlens = EarthLens(
    data_source="ecmwf",
    temporal_resolution="daily",
    start="2022-06-01",
    end="2022-08-31",
    variables={
        "reanalysis-era5-single-levels": [
            "2m-temperature",
            "2m-dewpoint-temperature",
        ],
    },
    lat_lim=[37.0, 38.0],
    lon_lim=[23.0, 24.0],
    path=str(OUT),
)
earthlens.download(aggregate=AggregationConfig(freq="1D", op="mean"))

2026-05-10 01:38:55.920 | INFO     | earthlens.ecmwf.backend:download:536 - Download ECMWF reanalysis-era5-single-levels/2m-temperature data for period 2022-06-01 00:00:00 till 2022-08-31 00:00:00

2026-05-10 01:38:57.079 | INFO     | earthlens.ecmwf.backend:_api:724 - Requesting reanalysis-era5-single-levels from CDS; this may take several minutes

2026-05-10 01:38:57,372 INFO [2025-12-11T00:00:00] Please note that a dedicated catalogue entry for this dataset, post-processed and stored in Analysis Ready Cloud Optimized (ARCO) format (Zarr), is available for optimised time-series retrievals (i.e. for retrieving data from selected variables for a single point over an extended period of time in an efficient way). You can discover it [here](https://cds.climate.copernicus.eu/datasets/reanalysis-era5-single-levels-timeseries?tab=overview)

2026-05-10 01:38:57,373 INFO Request ID is d5fd4922-d429-41b0-a71c-134c2d8ace0c

2026-05-10 01:38:57,431 INFO status has been updated to accepted

2026-05-10 01:39:18,624 INFO status has been updated to running

2026-05-10 01:39:30,143 INFO status has been updated to successful

2026-05-10 01:39:31 | INFO | pyramids.base.config | Logging is configured.

2026-05-10 01:39:32.480 | INFO     | earthlens.ecmwf.backend:download:536 - Download ECMWF reanalysis-era5-single-levels/2m-dewpoint-temperature data for period 2022-06-01 00:00:00 till 2022-08-31 00:00:00

2026-05-10 01:39:32.488 | INFO     | earthlens.ecmwf.backend:_api:724 - Requesting reanalysis-era5-single-levels from CDS; this may take several minutes

2026-05-10 01:39:32,983 INFO [2025-12-11T00:00:00] Please note that a dedicated catalogue entry for this dataset, post-processed and stored in Analysis Ready Cloud Optimized (ARCO) format (Zarr), is available for optimised time-series retrievals (i.e. for retrieving data from selected variables for a single point over an extended period of time in an efficient way). You can discover it [here](https://cds.climate.copernicus.eu/datasets/reanalysis-era5-single-levels-timeseries?tab=overview)

2026-05-10 01:39:32 | INFO | ecmwf.datastores.legacy_client | [2025-12-11T00:00:00] Please note that a dedicated catalogue entry for this dataset, post-processed and stored in Analysis Ready Cloud Optimized (ARCO) format (Zarr), is available for optimised time-series retrievals (i.e. for retrieving data from selected variables for a single point over an extended period of time in an efficient way). You can discover it [here](https://cds.climate.copernicus.eu/datasets/reanalysis-era5-single-levels-timeseries?tab=overview)

2026-05-10 01:39:32,984 INFO Request ID is 04ef6ffe-a609-4c02-a280-0f5088c4b55d

2026-05-10 01:39:32 | INFO | ecmwf.datastores.legacy_client | Request ID is 04ef6ffe-a609-4c02-a280-0f5088c4b55d

2026-05-10 01:39:33,088 INFO status has been updated to accepted

2026-05-10 01:39:33 | INFO | ecmwf.datastores.legacy_client | status has been updated to accepted

2026-05-10 01:40:06,335 INFO status has been updated to successful

2026-05-10 01:40:06 | INFO | ecmwf.datastores.legacy_client | status has been updated to successful

2026-05-10 01:40:06 | INFO | multiurl.base | Downloading https://object-store.os-api.cci2.ecmwf.int:443/cci2-prod-cache-1/2026-05-09/a9ed827790516684fd750f65857c3f9e.nc

2026-05-10 01:40:08.223 | INFO     | earthlens.ecmwf.backend:download:575 - ECMWF download summary: all 2 variables succeeded ([('reanalysis-era5-single-levels', '2m-temperature'), ('reanalysis-era5-single-levels', '2m-dewpoint-temperature')])

Step 2 — derive RH and wet-bulb from T and Td¶

Magnus-Tetens for saturation vapour pressure, then RH = e/es. Stull (2011) for wet-bulb temperature.

In [2]:

Copied!





import numpy as np
import pandas as pd
from pyramids.dataset import Dataset

agg = OUT / "aggregated"
t_air = np.stack([
    Dataset.read_file(str(p)).read_array()
    for p in sorted(agg.glob("2m_temperature_1D_*.tif"))
]) - 273.15
t_dew = np.stack([
    Dataset.read_file(str(p)).read_array()
    for p in sorted(agg.glob("2m_dewpoint_temperature_1D_*.tif"))
]) - 273.15

# Magnus formula for saturation vapour pressure (hPa).
def es(t_celsius):
    return 6.112 * np.exp(17.67 * t_celsius / (t_celsius + 243.5))

rh = 100.0 * es(t_dew) / es(t_air)

# Stull (2011) wet-bulb approximation, valid for rh > 5% and t -20..50 °C.
def wet_bulb(t, rh):
    return (
        t * np.arctan(0.151977 * (rh + 8.313659) ** 0.5)
        + np.arctan(t + rh)
        - np.arctan(rh - 1.676331)
        + 0.00391838 * (rh ** 1.5) * np.arctan(0.023101 * rh)
        - 4.686035
    )

t_wet_bulb = wet_bulb(t_air, rh)

site_t_air      = np.nanmean(t_air,      axis=(1, 2))
site_rh         = np.nanmean(rh,         axis=(1, 2))
site_t_wet_bulb = np.nanmean(t_wet_bulb, axis=(1, 2))

days = pd.date_range("2022-06-01", periods=len(site_t_air), freq="D")
summary = pd.DataFrame(
    {"T_2m [°C]": site_t_air.round(1), "RH [%]": site_rh.round(0),
     "T_wet-bulb [°C]": site_t_wet_bulb.round(1)},
    index=days,
)
summary.head()
import numpy as np
import pandas as pd
from pyramids.dataset import Dataset

agg = OUT / "aggregated"
t_air = np.stack([
    Dataset.read_file(str(p)).read_array()
    for p in sorted(agg.glob("2m_temperature_1D_*.tif"))
]) - 273.15
t_dew = np.stack([
    Dataset.read_file(str(p)).read_array()
    for p in sorted(agg.glob("2m_dewpoint_temperature_1D_*.tif"))
]) - 273.15

# Magnus formula for saturation vapour pressure (hPa).
def es(t_celsius):
    return 6.112 * np.exp(17.67 * t_celsius / (t_celsius + 243.5))

rh = 100.0 * es(t_dew) / es(t_air)

# Stull (2011) wet-bulb approximation, valid for rh > 5% and t -20..50 °C.
def wet_bulb(t, rh):
    return (
        t * np.arctan(0.151977 * (rh + 8.313659) ** 0.5)
        + np.arctan(t + rh)
        - np.arctan(rh - 1.676331)
        + 0.00391838 * (rh ** 1.5) * np.arctan(0.023101 * rh)
        - 4.686035
    )

t_wet_bulb = wet_bulb(t_air, rh)

site_t_air      = np.nanmean(t_air,      axis=(1, 2))
site_rh         = np.nanmean(rh,         axis=(1, 2))
site_t_wet_bulb = np.nanmean(t_wet_bulb, axis=(1, 2))

days = pd.date_range("2022-06-01", periods=len(site_t_air), freq="D")
summary = pd.DataFrame(
    {"T_2m [°C]": site_t_air.round(1), "RH [%]": site_rh.round(0),
     "T_wet-bulb [°C]": site_t_wet_bulb.round(1)},
    index=days,
)
summary.head()

Out[2]:

	T_2m [°C]	RH [%]	T_wet-bulb [°C]
2022-06-01	24.200001	66.0	19.600000
2022-06-02	25.600000	56.0	19.400000
2022-06-03	26.200001	53.0	19.400000
2022-06-04	26.000000	54.0	19.299999
2022-06-05	25.700001	57.0	19.600000

Step 3 — plot the heat-stress time series¶

Wet-bulb temperatures > 30 °C are physiologically stressful for healthy adults; > 35 °C is the theoretical limit of human survival (no evaporative cooling possible at all).

In [3]:

Copied!





import matplotlib.pyplot as plt

fig, ax = plt.subplots(figsize=(11, 5))
ax.plot(days, site_t_air, label="Air T_2m", color="tab:red", lw=1.5)
ax.plot(days, site_t_wet_bulb, label="Wet-bulb T", color="tab:purple", lw=2)
ax.axhline(30, color="orange", lw=0.8, ls="--", label="30 °C wet-bulb (stressful)")
ax.axhline(35, color="red", lw=0.8, ls="--", label="35 °C wet-bulb (survival limit)")
ax.fill_between(days, 30, np.maximum(site_t_wet_bulb, 30), where=site_t_wet_bulb >= 30,
                color="orange", alpha=0.25)
ax.set_ylabel("Temperature [°C]")
ax.set_title("Athens summer 2022 — daily-mean air and wet-bulb temperatures")
ax.legend()
ax.grid(alpha=0.3)
plt.tight_layout()
plt.show()
import matplotlib.pyplot as plt

fig, ax = plt.subplots(figsize=(11, 5))
ax.plot(days, site_t_air, label="Air T_2m", color="tab:red", lw=1.5)
ax.plot(days, site_t_wet_bulb, label="Wet-bulb T", color="tab:purple", lw=2)
ax.axhline(30, color="orange", lw=0.8, ls="--", label="30 °C wet-bulb (stressful)")
ax.axhline(35, color="red", lw=0.8, ls="--", label="35 °C wet-bulb (survival limit)")
ax.fill_between(days, 30, np.maximum(site_t_wet_bulb, 30), where=site_t_wet_bulb >= 30,
                color="orange", alpha=0.25)
ax.set_ylabel("Temperature [°C]")
ax.set_title("Athens summer 2022 — daily-mean air and wet-bulb temperatures")
ax.legend()
ax.grid(alpha=0.3)
plt.tight_layout()
plt.show()

No description has been provided for this image

Notes¶

Daily mean understates peak heat. Real heat-wave health exposure depends on daily maximum wet-bulb. To get it, retrieve hourly data (temporal_resolution="daily") and run aggregate_netcdf with op="max" instead of op="mean".
CDS-Beta has hourly products. Use the derived-era5-single- levels-daily-statistics family for pre-aggregated daily max / min statistics — but those NetCDFs are already daily aggregates, so pass op="mean" explicitly when re-aggregating to coarser steps to avoid double-summing.
Health categories. US NWS heat index categories: Caution > 27 °C, Extreme Caution > 32 °C, Danger > 41 °C, Extreme Danger > 54 °C (apparent temperature, not wet-bulb). The Stull wet-bulb here is an alternative — closer to the thermodynamic limit on evaporative cooling.
Stull approximation limits. Valid for $T \in [-20, 50]$°C and $RH > 5\%$. Outside that range use the iterative Davies-Jones (2008) formulation.