Skip to content

fr_weather_reader module

simpest.models.fr_weather_reader

Weather readers that build an hourly forcing series for the model.

The disease sub-model operates on an hourly time step, so these readers ingest either daily or hourly weather CSV files and return a complete hourly series. When the input is daily, a physically based diurnal synthesis fills in the 24 hourly records for each day.

Key behaviour:

  • Header parsing. Column names are matched case-insensitively against a set of common aliases, so files from different sources are accepted without renaming.
  • Date columns. A single Date/Datetime column or separate Year/Month/Day (plus Hour for hourly files) are both accepted.
  • Radiation. Measured radiation is used when present; otherwise daily global solar radiation is estimated with the Hargreaves–Samani relationship from the daily temperature range and distributed over the day by clear-sky (extraterrestrial) fractions.
  • Humidity. Relative humidity is reconstructed from daily RHx/RHn extrema with a cosine diurnal curve, or estimated from a dew-point relation when extrema are absent.

read_daily(file, start_year, end_year, latitude=0.0)

Read a daily weather CSV and return synthesized hourly records.

Supports year/month/day columns OR a single date/datetime column. Radiation is estimated when absent (Hargreaves-Samani) if latitude is known.

Parameters:

Name Type Description Default
file str | Path

Path to the daily weather CSV.

required
start_year int

First year to include (inclusive).

required
end_year int

Last year to include (inclusive).

required
latitude float

Fallback latitude in decimal degrees when not in CSV.

0.0

Returns:

Type Description
Dict[datetime, InputsHourly]

Dict keyed by datetime(year, month, day, hour)InputsHourly.

Source code in simpest/models/fr_weather_reader.py
def read_daily(
    file: str | Path,
    start_year: int,
    end_year: int,
    latitude: float = 0.0,
) -> Dict[datetime, InputsHourly]:
    """Read a daily weather CSV and return synthesized hourly records.

    Supports year/month/day columns OR a single date/datetime column.
    Radiation is estimated when absent (Hargreaves-Samani) if latitude is known.

    Args:
        file:       Path to the daily weather CSV.
        start_year: First year to include (inclusive).
        end_year:   Last year to include (inclusive).
        latitude:   Fallback latitude in decimal degrees when not in CSV.

    Returns:
        Dict keyed by ``datetime(year, month, day, hour)`` → ``InputsHourly``.
    """
    result: Dict[datetime, InputsHourly] = {}
    path = Path(file)

    with path.open(newline="", encoding="utf-8-sig") as fh:
        sample = fh.read(2048)
        fh.seek(0)
        delimiter = "\t" if "\t" in sample.split("\n")[0] else ","

        reader = csv.reader(fh, delimiter=delimiter)
        headers = [_clean(h) for h in next(reader)]
        col = {name: i for i, name in enumerate(headers)}

        # Date columns
        date_idx = _idx(col, ["date", "datetime", "timestamp"], optional=True)
        year_idx = _idx(col, ["year"], optional=True)
        month_idx = _idx(col, ["month"], optional=True)
        day_idx = _idx(col, ["day"], optional=True)

        has_ymd = (year_idx >= 0 and month_idx >= 0 and day_idx >= 0)
        has_date = date_idx >= 0


        if not has_ymd and not has_date:
            raise ValueError(
                "Weather file must contain (year, month, day) OR a single date/datetime column."
            )

        # Meteorological columns
        rad_idx = _idx(col, ["rad", "radiation", "solar", "solarrad", "srad"], optional=True)
        lat_idx = _idx(col, ["lat", "latitude", "sitelat", "phi"], optional=True)
        tmax_idx = _idx(col, ["tmax", "t2mmax", "maxtemp", "tx"])
        tmin_idx = _idx(col, ["tmin", "t2mmin", "mintemp", "tn"])
        prec_idx = _idx(col, ["prec", "precip", "rain", "rainfall", "precipitation", "p"], optional=True)
        rhx_idx = _idx(col, ["rhmax", "humiditymax", "relativehumiditymax", "hummax", "rhx"], optional=True)
        rhn_idx = _idx(col, ["rhmin", "humiditymin", "relativehumiditymin", "hummin", "rhn"], optional=True)

        if rad_idx < 0 and lat_idx < 0 and latitude == 0.0:
            raise ValueError(
                "Weather file must contain a radiation column OR a latitude column "
                "(or pass latitude= to read_daily)."
            )

        for row in reader:
            if not row or all(c.strip() == "" for c in row):
                continue

            # --- Date ---
            try:
                if has_ymd:
                    yr = int(row[year_idx])
                    mo = int(row[month_idx])
                    dy = int(row[day_idx])
                    day_date = date(yr, mo, dy)
                else:
                    raw_date = row[date_idx].strip().strip('"')
                    day_date = datetime.fromisoformat(raw_date).date()
            except (ValueError, IndexError):
                continue

            if not (start_year <= day_date.year <= end_year):
                continue

            # --- Temperatures ---
            try:
                tmax = _pf(row, tmax_idx)
                tmin = _pf(row, tmin_idx)
            except (ValueError, IndexError):
                continue

            # --- Radiation ---
            rad = _pf(row, rad_idx) if rad_idx >= 0 else float("nan")
            rad_missing = not math.isfinite(rad) or rad <= 0.0

            # --- Latitude (for radiation estimation) ---
            lat = _pf(row, lat_idx) if lat_idx >= 0 else latitude

            if not math.isfinite(lat):
                lat = latitude

            if rad_missing:
                if lat == 0.0:
                    continue  # cannot estimate without latitude
                rd = _day_length(day_date, lat, tmax, tmin)
                rad = rd["gsr"]
                if not math.isfinite(rad) or rad <= 0.0:
                    continue

            # --- Build InputsDaily ---
            prec = _pf(row, prec_idx) if prec_idx >= 0 else 0.0
            if not math.isfinite(prec):
                prec = 0.0

            id_ = InputsDaily(
                date=datetime(day_date.year, day_date.month, day_date.day),
                tmax=tmax,
                tmin=tmin,
                rad=rad,
                precipitation=prec,
                latitude=lat,
            )
            id_.dew_point = _dew_point(id_.tmax, id_.tmin)

            if rhx_idx >= 0:
                v = _pf(row, rhx_idx)
                if math.isfinite(v):
                    id_.rhx = v
            if rhn_idx >= 0:
                v = _pf(row, rhn_idx)
                if math.isfinite(v):
                    id_.rhn = v

            # Calculate leaf wetness from RH and precipitation
            # High humidity (>80%) or significant precipitation (>0.5mm) indicate wet leaves
            avg_rh = (id_.rhx + id_.rhn) / 2.0 if id_.rhx > 0 else 0.0
            if prec >= 0.5 or avg_rh >= 80.0:
                # Estimate wetness hours: scale from 60% RH (minimal) to 100% RH (24 hours)
                if avg_rh >= 80.0:
                    id_.leaf_wetness = min(24.0, (avg_rh - 80.0) * 1.2)  # 80% -> 0h, 100% -> 24h
                else:
                    id_.leaf_wetness = 0.0
                # Add hours if precipitation present
                if prec >= 0.5:
                    id_.leaf_wetness = min(24.0, id_.leaf_wetness + prec * 2.0)  # Add ~2h per mm rainfall

            # --- Synthesize 24 hourly records ---
            hourly = _estimate_hourly(id_, day_date, lat)
            result.update(hourly)

    return result

read_hourly(file, start_year, end_year, site='', latitude=0.0)

Read an hourly weather CSV and return the hourly records.

Missing hourly radiation is filled per-day using the Hargreaves estimate distributed via clear-sky (ETR) fractions.

Parameters:

Name Type Description Default
file str | Path

Path to the hourly weather CSV.

required
start_year int

First year to include (inclusive).

required
end_year int

Last year to include (inclusive).

required
site str

Optional site filter (column 'site' or 'station').

''
latitude float

Fallback latitude when not in the file.

0.0

Returns:

Type Description
Dict[datetime, InputsHourly]

Dict keyed by datetime(year, month, day, hour)InputsHourly.

Source code in simpest/models/fr_weather_reader.py
def read_hourly(
    file: str | Path,
    start_year: int,
    end_year: int,
    site: str = "",
    latitude: float = 0.0,
) -> Dict[datetime, InputsHourly]:
    """Read an hourly weather CSV and return the hourly records.

    Missing hourly radiation is filled per-day using the Hargreaves estimate
    distributed via clear-sky (ETR) fractions.

    Args:
        file:       Path to the hourly weather CSV.
        start_year: First year to include (inclusive).
        end_year:   Last year to include (inclusive).
        site:       Optional site filter (column 'site' or 'station').
        latitude:   Fallback latitude when not in the file.

    Returns:
        Dict keyed by ``datetime(year, month, day, hour)`` → ``InputsHourly``.
    """
    result: Dict[datetime, InputsHourly] = {}
    daily_buckets: Dict[date, List[InputsHourly]] = {}

    path = Path(file)

    with path.open(newline="", encoding="utf-8-sig") as fh:
        sample = fh.read(2048)
        fh.seek(0)
        delimiter = "\t" if "\t" in sample.split("\n")[0] else ","
        reader = csv.reader(fh, delimiter=delimiter)
        headers = [_clean(h) for h in next(reader)]
        col = {name: i for i, name in enumerate(headers)}

        year_idx = _idx(col, ["year"], optional=True)
        month_idx = _idx(col, ["month", "mo"], optional=True)
        day_idx = _idx(col, ["day", "dd", "dy"], optional=True)
        hour_idx = _idx(col, ["hour", "hr", "h"], optional=True)
        date_idx = _idx(col, ["date", "datetime", "timestamp"], optional=True)

        has_ymdh = (year_idx >= 0 and month_idx >= 0 and day_idx >= 0 and hour_idx >= 0)
        has_dateh = (date_idx >= 0 and hour_idx >= 0)

        if not has_ymdh and not has_dateh:
            raise ValueError(
                "Hourly file must contain (year,month,day,hour) OR (date,hour)."
            )

        tmax_idx = _idx(col, ["tmax", "t2mmax", "maxtemp"], optional=True)
        tmin_idx = _idx(col, ["tmin", "t2mmin", "mintemp"], optional=True)
        temp_idx = _idx(col, ["temp", "temperature", "t2m"], optional=True)
        prec_idx = _idx(col, ["prec", "precip", "precipitation", "prectotcorr", "rain", "rainfall"], optional=True)
        rh_idx = _idx(col, ["rh", "humidity", "relhumidity", "relativehumidity"], optional=True)
        rad_idx = _idx(col, ["rad", "radiation", "solar", "solarrad"], optional=True)
        lat_idx = _idx(col, ["latitude", "lat"], optional=True)

        if rad_idx < 0 and lat_idx < 0 and latitude == 0.0:
            raise ValueError(
                "Hourly weather file must contain a radiation column OR a latitude column."
            )

        for row in reader:
            if not row or all(c.strip() == "" for c in row):
                continue

            # --- Timestamp ---
            try:
                if has_ymdh:
                    yr = int(row[year_idx])
                    mo = int(row[month_idx])
                    dy = int(row[day_idx])
                    hr = int(row[hour_idx])
                    ts = datetime(yr, mo, dy, hr)
                else:
                    raw = row[date_idx].strip().strip('"')
                    ts = datetime.fromisoformat(raw)
                    ts = ts.replace(minute=0, second=0, microsecond=0)
                    hr = int(row[hour_idx])
                    ts = ts.replace(hour=hr)
            except (ValueError, IndexError):
                continue

            if not (start_year <= ts.year <= end_year):
                continue

            gw = InputsHourly(date=ts)

            # Temperature
            if temp_idx >= 0:
                v = _pf(row, temp_idx)
                if math.isfinite(v):
                    gw.air_temperature = v
            elif tmax_idx >= 0 and tmin_idx >= 0:
                tx = _pf(row, tmax_idx)
                tn = _pf(row, tmin_idx)
                if math.isfinite(tx) and math.isfinite(tn):
                    gw.air_temperature = 0.5 * (tx + tn)

            # Precipitation
            if prec_idx >= 0:
                v = _pf(row, prec_idx)
                if math.isfinite(v) and v >= 0:
                    gw.precipitation = v

            # Relative humidity
            if rh_idx >= 0:
                v = _pf(row, rh_idx)
                if math.isfinite(v):
                    gw.relative_humidity = max(0.0, min(100.0, v))
            else:
                t = gw.air_temperature
                dp = _dew_point(t, t)
                es = 0.61121 * math.exp((17.502 * t) / (240.97 + t))
                ea = 0.61121 * math.exp((17.502 * dp) / (240.97 + dp))
                gw.relative_humidity = max(0.0, min(100.0, ea / es * 100.0)) if es > 0 else 0.0

            # Radiation
            if rad_idx >= 0:
                v = _pf(row, rad_idx)
                if math.isfinite(v) and v >= 0:
                    gw.rad = v

            # Latitude
            if lat_idx >= 0:
                v = _pf(row, lat_idx)
                if math.isfinite(v):
                    gw.latitude = v
            elif latitude != 0.0:
                gw.latitude = latitude

            # Leaf wetness
            gw.leaf_wetness = 1.0 if (gw.relative_humidity > 90.0 or gw.precipitation >= 0.2) else 0.0

            day_key = ts.date()
            daily_buckets.setdefault(day_key, []).append(gw)

    # --- Post-process: fill missing hourly radiation per day ---
    for day_key, records in daily_buckets.items():
        if not records:
            continue

        lat = records[0].latitude if records[0].latitude else latitude
        temps = [r.air_temperature for r in records]
        tmax_d = max(temps)
        tmin_d = min(temps)

        all_rad_zero = all(r.rad <= 0.0 for r in records)

        rd = None
        if all_rad_zero and lat:
            rd = _day_length(day_key, lat, tmax_d, tmin_d)

        for rec in records:
            hr = rec.date.hour
            if rec.rad <= 0.0 and rd is not None:
                rec.rad = rd["gsr_hourly"][hr]
            rec.leaf_wetness = 1.0 if (rec.relative_humidity > 90.0 or rec.precipitation >= 0.2) else 0.0
            result[rec.date] = rec

    return result