From b684748c6f36ccc35643496b080336a367241d3a Mon Sep 17 00:00:00 2001 From: Andreas Date: Wed, 16 Sep 2026 18:23:35 +0200 Subject: [PATCH] feat(pvforecast): add calibrated local Akkudoktor backend Port local PV modeling and outage calibration from feature commits f976335, 6dc58c3 and faed0fd by Andreas. Keep PVForecastAkkudoktor identity and remote default, adapt to async storage, and migrate legacy provider settings. --- docs/akkudoktoreos/prediction.md | 56 + src/akkudoktoreos/config/configmigrate.py | 18 + src/akkudoktoreos/prediction/pvforecast.py | 17 + .../prediction/pvforecastakkudoktor.py | 12 +- .../prediction/pvforecastakkudoktorlocal.py | 1248 +++++++++++++++++ tests/test_pvforecastakkudoktorlocal.py | 725 ++++++++++ 6 files changed, 2071 insertions(+), 5 deletions(-) create mode 100644 src/akkudoktoreos/prediction/pvforecastakkudoktorlocal.py create mode 100644 tests/test_pvforecastakkudoktorlocal.py diff --git a/docs/akkudoktoreos/prediction.md b/docs/akkudoktoreos/prediction.md index a18a3a7f..7aabbd67 100644 --- a/docs/akkudoktoreos/prediction.md +++ b/docs/akkudoktoreos/prediction.md @@ -839,6 +839,62 @@ Example: } ``` +#### Local Open-Meteo/pvlib backend + +`PVForecastAkkudoktor` keeps the remote API as its default backend. To calculate +PV output locally, keep the same public provider ID and set: + +```json +{ + "pvforecast": { + "provider": "PVForecastAkkudoktor", + "akkudoktor": { + "backend": "local", + "resolution_minutes": 15, + "calibration_enabled": false + } + } +} +``` + +Configure `general.latitude`, `general.longitude` and `pvforecast.planes` as +usual. The local PVWatts model uses plane peak power, tilt, azimuth, horizon, +tracking, mounting, losses and inverter power. It does not use the CEC module +and inverter names required by `PVForecastPVLib`; that provider remains a +separate option backed by the selected weather provider. + +The local backend requests irradiance and weather directly from Open-Meteo, +then computes DC and AC power with pvlib. Both prediction keys retain their +units of **W**. Fifteen-minute power samples represent 0.25 hours when converted +to Wh by an optimizer. Hourly consumers obtain the hourly mean power. Neither +backend selects or changes the optimization algorithm or the legacy `/optimize` +request schema. + +Open-Meteo radiation timestamps label the preceding interval. By default EOS +shifts them to interval starts and uses the interval midpoint for solar position. +Native 15-minute weather is available in Central Europe and North America; +elsewhere Open-Meteo interpolates hourly weather. See the +[Open-Meteo API documentation](https://open-meteo.com/en/docs). + +Set `calibration_enabled` to true and configure cumulative PV production meter +keys in `measurement.pv_production_emr_keys` to fit the model against measurements. +Calibration preserves daily energy while optionally correcting the intraday +azimuth shape, clamps correction factors and respects inverter capacity. +Hourly meters are fitted hourly; native quarter-hour meters retain their cadence. +Probable outage or curtailment days are excluded by default. Without sufficient +usable measurements the uncalibrated physical model is used. + +Transient weather errors are retried. If an update still fails, existing usable +AC forecasts are retained; a cold start or expired-only forecast reports the +error. Retention does not extend the available forecast horizon. Calibration and +forecast settings are listed in the generated PV configuration reference. + +Feature-branch configurations using `PVForecastAkkudoktorLocal` and +`provider_settings.PVForecastAkkudoktorLocal` migrate to the public provider ID +and `akkudoktor.backend = "local"`. Explicit current `akkudoktor` values take +precedence over legacy values. Invalid local settings fail validation instead +of being silently dropped during file migration. + ### PVForecastVrm Provider The `PVForecastVrm` provider retrieves pv power forecast data from the Victron Remote Management diff --git a/src/akkudoktoreos/config/configmigrate.py b/src/akkudoktoreos/config/configmigrate.py index 30291e70..f400dc37 100644 --- a/src/akkudoktoreos/config/configmigrate.py +++ b/src/akkudoktoreos/config/configmigrate.py @@ -255,8 +255,26 @@ def migrate_config_data(config_data: Dict[str, Any]) -> "SettingsEOSDefaults": skipped_paths = [] from akkudoktoreos.config.config import SettingsEOSDefaults + from akkudoktoreos.prediction.pvforecastakkudoktorlocal import ( + PVForecastAkkudoktorLocalCommonSettings, + normalize_akkudoktor_settings, + ) + + # Normalize this provider before the generic field-by-field transfer. Validate + # coupled bounds together so a transient intermediate default cannot lose them. + config_data = dict(config_data) + pv_settings = normalize_akkudoktor_settings(config_data.get("pvforecast")) + if isinstance(pv_settings, dict): + config_data["pvforecast"] = pv_settings new_config = SettingsEOSDefaults() + if isinstance(pv_settings, dict) and "akkudoktor" in pv_settings: + local_settings = PVForecastAkkudoktorLocalCommonSettings.model_validate( + pv_settings["akkudoktor"] + ) + new_config.set_nested_value("pvforecast/akkudoktor", local_settings) + migrated_source_paths.add("pvforecast/akkudoktor") + mapped_count += 1 # 1) Apply explicit migration map for old_path, mapping in MIGRATION_MAP.items(): diff --git a/src/akkudoktoreos/prediction/pvforecast.py b/src/akkudoktoreos/prediction/pvforecast.py index 7d0d6650..ed2dc3ff 100644 --- a/src/akkudoktoreos/prediction/pvforecast.py +++ b/src/akkudoktoreos/prediction/pvforecast.py @@ -7,6 +7,10 @@ from pydantic import Field, computed_field, field_validator, model_validator from akkudoktoreos.config.configabc import SettingsBaseModel from akkudoktoreos.core.coreabc import get_prediction from akkudoktoreos.prediction.pvforecastabc import PVForecastProvider +from akkudoktoreos.prediction.pvforecastakkudoktorlocal import ( + PVForecastAkkudoktorLocalCommonSettings, + normalize_akkudoktor_settings, +) from akkudoktoreos.prediction.pvforecastforecastsolar import ( PVForecastForecastSolarCommonSettings, ) @@ -202,6 +206,13 @@ class PVForecastCommonSettings(SettingsBaseModel): }, ) + akkudoktor: PVForecastAkkudoktorLocalCommonSettings = Field( + default_factory=PVForecastAkkudoktorLocalCommonSettings, + json_schema_extra={ + "description": "Akkudoktor forecast backend and local calibration settings", + }, + ) + pvforecastimport: PVForecastImportCommonSettings = Field( default_factory=PVForecastImportCommonSettings, json_schema_extra={"description": "PV forecast import provider settings"}, @@ -300,6 +311,12 @@ class PVForecastCommonSettings(SettingsBaseModel): return pvforecast_provider_ids() # Validators + @model_validator(mode="before") + @classmethod + def migrate_local_akkudoktor_settings(cls, data: Any) -> Any: + """Accept the feature branch's local backend configuration.""" + return normalize_akkudoktor_settings(data) + @field_validator("provider", mode="after") @classmethod def validate_provider(cls, value: Optional[str]) -> Optional[str]: diff --git a/src/akkudoktoreos/prediction/pvforecastakkudoktor.py b/src/akkudoktoreos/prediction/pvforecastakkudoktor.py index 97e0eee4..94b40b2c 100644 --- a/src/akkudoktoreos/prediction/pvforecastakkudoktor.py +++ b/src/akkudoktoreos/prediction/pvforecastakkudoktor.py @@ -87,10 +87,8 @@ from pydantic import Field, ValidationError, computed_field, field_validator from akkudoktoreos.core.cache import cache_in_file from akkudoktoreos.core.pydantic import PydanticBaseModel from akkudoktoreos.prediction.pvforecast import PVForecastPlaneSetting -from akkudoktoreos.prediction.pvforecastabc import ( - PVForecastDataRecord, - PVForecastProvider, -) +from akkudoktoreos.prediction.pvforecastabc import PVForecastDataRecord +from akkudoktoreos.prediction.pvforecastakkudoktorlocal import PVForecastAkkudoktorLocal from akkudoktoreos.utils.datetimeutil import compare_datetimes, to_datetime @@ -188,7 +186,7 @@ class PVForecastAkkudoktorDataRecord(PVForecastDataRecord): return self.pvforecast_ac_power -class PVForecastAkkudoktor(PVForecastProvider[PVForecastAkkudoktorDataRecord]): +class PVForecastAkkudoktor(PVForecastAkkudoktorLocal[PVForecastAkkudoktorDataRecord]): """Fetch and process PV forecast data from akkudoktor.net. PVForecastAkkudoktor is a singleton-based class that retrieves weather forecast data @@ -385,6 +383,10 @@ class PVForecastAkkudoktor(PVForecastProvider[PVForecastAkkudoktorDataRecord]): # Get Akkudoktor PV Forecast data for the given configuration. if force_update: logger.info("[PVForecastAkkudoktor] force update.") + if self.config.pvforecast.akkudoktor.backend == "local": + await self._update_local_data(force_update=force_update) + return + akkudoktor_data = self._request_forecast(force_update=force_update) # type: ignore # Timezone of the PV system diff --git a/src/akkudoktoreos/prediction/pvforecastakkudoktorlocal.py b/src/akkudoktoreos/prediction/pvforecastakkudoktorlocal.py new file mode 100644 index 00000000..4a6abc14 --- /dev/null +++ b/src/akkudoktoreos/prediction/pvforecastakkudoktorlocal.py @@ -0,0 +1,1248 @@ +"""Native PV power forecast computed inside EOS from Open-Meteo weather with pvlib. + +This optional backend of PVForecastAkkudoktor does not ask a third-party service +for PV power. It fetches raw irradiance and weather from Open-Meteo and runs the +whole modelling chain locally with pvlib: + + solar position -> horizon shading -> transposition to the module plane -> + incidence-angle modifier -> cell temperature -> PVWatts DC -> inverter AC + +Why this exists: + +* **Horizon.** Open-Meteo serves up to 16 forecast days at 15-minute resolution in a + single request. Services that wrap it (including api.akkudoktor.net) cut the + horizon much shorter, which limits the optimizer forecast tail. +* **Call budget.** One request per hour against a ~10k/day non-commercial budget, + with a cached weather response instead of a remote PV model. +* **Honest parameters.** ``albedo``, inverter efficiency and the module temperature + coefficient become real configuration instead of constants baked into a URL. + +Two conventions matter and are handled explicitly here: + +* Open-Meteo radiation values are the mean over the **preceding** interval, so the + representative sun position for a value stamped ``t`` is ``t - interval/2``. +* EOS records label an interval by its **start** (``key_to_array`` resamples with + left-labelled buckets), so a value stamped ``t`` by Open-Meteo is stored at + ``t - interval``. Set ``shift_to_interval_start`` to False to keep the raw stamps. + +Note also that ``direct_radiation`` in the Open-Meteo API is beam irradiance on the +*horizontal* plane; the DNI this chain needs is ``direct_normal_irradiance``. +""" + +import asyncio +import math +import time +from collections.abc import Callable, Iterable +from typing import Any, Literal, Optional, Self + +import numpy as np +import pandas as pd +import pendulum +import pvlib +import requests +from loguru import logger +from pydantic import Field, field_validator, model_validator + +from akkudoktoreos.config.configabc import SettingsBaseModel +from akkudoktoreos.core.cache import cache_in_file +from akkudoktoreos.prediction.pvforecastabc import PVForecastDataRecordT, PVForecastProvider +from akkudoktoreos.utils.datetimeutil import compare_datetimes, to_datetime, to_duration + +OPENMETEO_URL = "https://api.open-meteo.com/v1/forecast" + +# Open-Meteo variables the pvlib chain needs. `direct_normal_irradiance` is the DNI; +# `direct_radiation` (beam on the horizontal) would be wrong here. +OPENMETEO_VARIABLES = ( + "temperature_2m", + "relative_humidity_2m", + "wind_speed_10m", + "shortwave_radiation", + "diffuse_radiation", + "direct_normal_irradiance", +) + +# Open-Meteo hard limits for a single forecast request. +MAX_FORECAST_DAYS = 16 +MAX_PAST_DAYS = 92 + +TRANSPOSITION_MODELS = ("isotropic", "klucher", "haydavies", "reindl", "king", "perez") + +# pvlib SAPM cell temperature parameter set per EOS `mountingplace`. +MOUNTING_TEMPERATURE_MODEL = { + "free": "open_rack_glass_glass", + "building": "close_mount_glass_glass", +} + + +def normalize_akkudoktor_settings(data: Any) -> Any: + """Translate feature-branch PV settings without modifying the caller's data. + + Current settings win over legacy entries. Invalid local settings are validated + by the current settings model; they must never disappear during file migration. + """ + if not isinstance(data, dict): + return data + providers = data.get("provider_settings") + legacy = providers.get("PVForecastAkkudoktorLocal") if isinstance(providers, dict) else None + selected_legacy = data.get("provider") == "PVForecastAkkudoktorLocal" + if legacy is None and not selected_legacy: + return data + if legacy is not None and not isinstance(legacy, dict): + raise ValueError("PVForecastAkkudoktorLocal settings must be an object") + current = data.get("akkudoktor", {}) + if not isinstance(current, dict): + raise ValueError("pvforecast.akkudoktor settings must be an object") + normalized = dict(data) + settings = dict(legacy or {}) + if selected_legacy: + normalized["provider"] = "PVForecastAkkudoktor" + settings["backend"] = "local" + settings.update(current) + normalized["akkudoktor"] = settings + if isinstance(providers, dict): + remaining = dict(providers) + remaining.pop("PVForecastAkkudoktorLocal", None) + if remaining: + normalized["provider_settings"] = remaining + else: + normalized.pop("provider_settings", None) + return normalized + + +class PVForecastAkkudoktorLocalCommonSettings(SettingsBaseModel): + """Common settings for the local (pvlib) PV forecast provider.""" + + backend: Literal["remote", "local"] = Field( + default="remote", + json_schema_extra={ + "description": "Akkudoktor forecast backend: remote API or local Open-Meteo/pvlib model.", + "examples": ["remote", "local"], + }, + ) + + @model_validator(mode="after") + def validate_calibration_bounds(self) -> Self: + """Reject inverted correction bounds.""" + if self.calibration_min_factor > self.calibration_max_factor: + raise ValueError("calibration_min_factor must not exceed calibration_max_factor") + return self + + resolution_minutes: int = Field( + default=15, + json_schema_extra={ + "description": ( + "Forecast resolution in minutes. 15 requests Open-Meteo's `minutely_15` " + "block (natively resolved over Central Europe and North America, " + "interpolated from hourly elsewhere); 60 requests the `hourly` block." + ), + "examples": [15, 60], + }, + ) + forecast_days: Optional[int] = Field( + default=None, + ge=1, + le=MAX_FORECAST_DAYS, + json_schema_extra={ + "description": ( + "Forecast horizon in days (1-16). Leave empty to derive it from " + "`prediction.hours`, which is what keeps the optimizer's tail horizon fed." + ), + "examples": [None, 7], + }, + ) + past_days: Optional[int] = Field( + default=None, + ge=0, + le=MAX_PAST_DAYS, + json_schema_extra={ + "description": ( + "Days of past data to request (0-92). Leave empty to derive it from " + "`prediction.historic_hours`." + ), + "examples": [None, 3], + }, + ) + weather_models: list[str] = Field( + default=["best_match"], + min_length=1, + json_schema_extra={ + "description": ( + "Open-Meteo weather models to request. Listing more than one turns the " + "input into a poor-man's ensemble: the members are averaged per variable, " + "which is the cheapest reliable way to cut irradiance forecast error. " + "Costs no extra API calls." + ), + "examples": [ + ["best_match"], + ["icon_seamless", "ecmwf_ifs025", "gfs_seamless"], + ], + }, + ) + transposition_model: str = Field( + default="perez", + json_schema_extra={ + "description": ( + "pvlib sky-diffuse transposition model: isotropic, klucher, haydavies, " + "reindl, king or perez." + ), + "examples": ["perez", "haydavies"], + }, + ) + albedo: float = Field( + default=0.25, + ge=0.0, + le=1.0, + json_schema_extra={ + "description": "Ground albedo used for planes that do not set their own.", + "examples": [0.25, 0.2], + }, + ) + inverter_efficiency: float = Field( + default=0.96, + gt=0.0, + le=1.0, + json_schema_extra={ + "description": "Nominal inverter efficiency (PVWatts eta_inv_nom).", + "examples": [0.96, 0.94], + }, + ) + temperature_coefficient: float = Field( + default=-0.36, + json_schema_extra={ + "description": ( + "Module power temperature coefficient in %/degC (negative). Matches the " + "`cellCoEff` the akkudoktor.net forecast uses." + ), + "examples": [-0.36, -0.29], + }, + ) + apply_iam: bool = Field( + default=True, + json_schema_extra={ + "description": "Apply the ASHRAE incidence-angle modifier to the beam component.", + "examples": [True], + }, + ) + shift_to_interval_start: bool = Field( + default=True, + json_schema_extra={ + "description": ( + "Open-Meteo stamps an interval mean with the interval END. EOS labels an " + "interval by its START, so records are shifted back by one interval. " + "Disable only to compare like-for-like against a provider that does not." + ), + "examples": [True], + }, + ) + + calibration_enabled: bool = Field( + default=False, + json_schema_extra={ + "description": ( + "Correct systematic model error against measured PV production. Requires " + "`measurement.pv_production_emr_keys` to be configured and fed. Fits a " + "global scale factor plus per-solar-azimuth factors, which is what catches " + "near-field shading the horizon profile misses." + ), + "examples": [True], + }, + ) + calibration_days: int = Field( + default=30, + ge=1, + le=MAX_PAST_DAYS, + json_schema_extra={ + "description": "Length of the measurement window used to fit the correction.", + "examples": [30, 14], + }, + ) + calibration_reference_days: int = Field( + default=30, + ge=3, + le=MAX_PAST_DAYS, + json_schema_extra={ + "description": ( + "Lookback used to distinguish healthy production from outages or " + "curtailment. If the calibration window contains too few healthy days, " + "the most recent healthy days from this reference window are used." + ), + "examples": [30, 14], + }, + ) + calibration_outage_filter_enabled: bool = Field( + default=True, + json_schema_extra={ + "description": ( + "Exclude days whose measured production is far below the recent healthy " + "plant level. This prevents inverter, battery and curtailment events from " + "being learned as permanent PV model losses." + ), + "examples": [True], + }, + ) + calibration_outage_threshold: float = Field( + default=0.55, + gt=0.0, + lt=1.0, + json_schema_extra={ + "description": ( + "A day is treated as unavailable when its measured/modelled energy ratio " + "is below this fraction of the robust healthy reference ratio." + ), + "examples": [0.55, 0.5], + }, + ) + calibration_min_healthy_days: int = Field( + default=3, + ge=1, + le=31, + json_schema_extra={ + "description": ( + "Minimum number of healthy days used for a fit. Older healthy days from " + "the reference window are added when the recent window contains fewer." + ), + "examples": [3], + }, + ) + calibration_azimuth_bin_degrees: int = Field( + default=45, + ge=0, + le=180, + json_schema_extra={ + "description": ( + "Width of the solar-azimuth bins for the correction. 0 fits a single " + "global factor only." + ), + "examples": [45, 30, 15, 0], + }, + ) + calibration_prior_kwh: float = Field( + default=5.0, + ge=0.0, + json_schema_extra={ + "description": ( + "Shrinkage strength: a bin needs this much modelled energy before its own " + "factor outweighs the global one. Higher is more conservative." + ), + "examples": [5.0, 20.0], + }, + ) + calibration_min_factor: float = Field( + default=0.5, + gt=0.0, + json_schema_extra={ + "description": "Lower clamp on any fitted correction factor.", + "examples": [0.5], + }, + ) + calibration_max_factor: float = Field( + default=1.5, + gt=0.0, + json_schema_extra={ + "description": "Upper clamp on any fitted correction factor.", + "examples": [1.5], + }, + ) + + @field_validator("resolution_minutes") + @classmethod + def validate_resolution(cls, value: int) -> int: + if value not in (15, 60): + raise ValueError(f"resolution_minutes must be 15 or 60, got {value}") + return value + + @field_validator("transposition_model") + @classmethod + def validate_transposition_model(cls, value: str) -> str: + if value not in TRANSPOSITION_MODELS: + raise ValueError( + f"Invalid transposition_model '{value}', expected one of {TRANSPOSITION_MODELS}" + ) + return value + + +class PVForecastAkkudoktorLocal(PVForecastProvider[PVForecastDataRecordT]): + """Compute the PV forecast locally from Open-Meteo irradiance using pvlib.""" + + @classmethod + def provider_id(cls) -> str: + """Return the unique identifier for the PV-Forecast-Provider.""" + return "PVForecastAkkudoktor" + + @property + def _settings(self) -> PVForecastAkkudoktorLocalCommonSettings: + settings = self.config.pvforecast.akkudoktor + if settings is None: + settings = PVForecastAkkudoktorLocalCommonSettings() + return settings + + # ------------------------------------------------------------------ request + + def _horizon_days(self) -> tuple[int, int]: + """Resolve (forecast_days, past_days), deriving them from prediction config.""" + settings = self._settings + + forecast_days = settings.forecast_days + if forecast_days is None: + # +1 day so the last requested hour is still covered when the run starts + # late in the day. + hours = self.config.prediction.hours or 48 + forecast_days = min(MAX_FORECAST_DAYS, max(1, math.ceil(hours / 24) + 1)) + + past_days = settings.past_days + if past_days is None: + historic_hours = self.config.prediction.historic_hours or 0 + past_days = math.ceil(historic_hours / 24) + if settings.calibration_enabled: + # The fit compares modelled against measured power over the same past + # intervals, so the weather for that window has to come back with the request. + calibration_lookback = settings.calibration_days + if settings.calibration_outage_filter_enabled: + calibration_lookback = max( + calibration_lookback, settings.calibration_reference_days + ) + past_days = max(past_days, calibration_lookback) + + return forecast_days, min(MAX_PAST_DAYS, past_days) + + @cache_in_file(with_ttl="1 hour") + def _request_local_forecast(self) -> Any: + """Fetch raw irradiance and weather from Open-Meteo.""" + latitude = self.config.general.latitude + longitude = self.config.general.longitude + if latitude is None or longitude is None: + raise ValueError( + "PVForecastAkkudoktorLocal needs general.latitude and general.longitude" + ) + + settings = self._settings + block = "minutely_15" if settings.resolution_minutes == 15 else "hourly" + forecast_days, past_days = self._horizon_days() + + params = { + "latitude": latitude, + "longitude": longitude, + block: ",".join(OPENMETEO_VARIABLES), + # Ask for UTC so the returned stamps are unambiguous across DST changes. + "timezone": "UTC", + # pvlib's SAPM cell temperature model wants m/s; Open-Meteo defaults to km/h. + "wind_speed_unit": "ms", + "forecast_days": forecast_days, + "past_days": past_days, + "models": ",".join(settings.weather_models), + } + + response = None + for attempt in range(1, 4): + try: + response = requests.get(OPENMETEO_URL, params=params, timeout=(5, 30)) + logger.debug(f"Requesting Open-Meteo forecast: {response.url}") + response.raise_for_status() + break + except requests.RequestException as e: + response = None + status = getattr(e.response, "status_code", None) + # Open-Meteo answers 503 while it rotates its model runs and 429 + # when the free tier is briefly saturated. Both clear in seconds. + retryable = status is None or status in (429, 500, 502, 503, 504) + if not retryable or attempt == 3: + logger.error(f"Failed to fetch weather for local pvforecast: {e}") + raise RuntimeError("Failed to fetch weather from Open-Meteo API") from e + logger.warning( + "Open-Meteo request attempt {}/3 failed for local pvforecast: {}", attempt, e + ) + time.sleep(2 * attempt) + + if response is None: + raise RuntimeError("No response from Open-Meteo API") + data = response.json() + if block not in data: + raise ValueError( + f"Open-Meteo response is missing the '{block}' block: {list(data.keys())}" + ) + + return data + + # ------------------------------------------------------------------ weather + + def _weather_frame(self, data: Any) -> pd.DataFrame: + """Build a UTC-indexed weather frame from the Open-Meteo response.""" + block = "minutely_15" if self._settings.resolution_minutes == 15 else "hourly" + raw = data[block] + + index = pd.DatetimeIndex( + [pendulum.parse(str(t), tz="UTC") for t in raw["time"]], name="time" + ).tz_convert("UTC") + + frame = pd.DataFrame(index=index) + for variable in OPENMETEO_VARIABLES: + # With a single model Open-Meteo returns the bare variable name; with several + # it suffixes each member (`shortwave_radiation_icon_seamless`). Average the + # members that actually came back - not every model carries every variable. + members = [ + pd.to_numeric(pd.Series(raw[key], index=index), errors="coerce") + for key in raw + if key == variable or key.startswith(f"{variable}_") + ] + if not members: + raise ValueError(f"Open-Meteo response is missing '{variable}'") + frame[variable] = pd.concat(members, axis=1).mean(axis=1, skipna=True) + + # Night rows and occasional gaps arrive as null. Irradiance is genuinely zero + # then; temperature and wind are interpolated so the temperature model stays + # defined instead of poisoning the whole row with NaN. + for variable in ("shortwave_radiation", "diffuse_radiation", "direct_normal_irradiance"): + frame[variable] = frame[variable].fillna(0.0).clip(lower=0.0) + for variable in ("temperature_2m", "relative_humidity_2m", "wind_speed_10m"): + frame[variable] = frame[variable].interpolate(method="time").ffill().bfill() + frame["wind_speed_10m"] = frame["wind_speed_10m"].fillna(1.0).clip(lower=0.0) + frame["temperature_2m"] = frame["temperature_2m"].fillna(15.0) + + return frame + + # ------------------------------------------------------------------ geometry + + @staticmethod + def _horizon_elevation(userhorizon: list[float], solar_azimuth: np.ndarray) -> np.ndarray: + """Interpolate the horizon elevation at each solar azimuth. + + ``userhorizon`` follows the PVGIS convention: elevations in degrees at equally + spaced azimuths clockwise from north, the first entry being due north. The + profile wraps around, so the first entry is repeated at 360 degrees. + """ + horizon = np.asarray(userhorizon, dtype=float) + count = len(horizon) + azimuths = np.arange(count, dtype=float) * (360.0 / count) + return np.interp( + np.asarray(solar_azimuth, dtype=float) % 360.0, + np.append(azimuths, 360.0), + np.append(horizon, horizon[0]), + ) + + @staticmethod + def _tracked_orientation(plane: Any, solpos: pd.DataFrame) -> tuple[Any, Any]: + """Return (surface_tilt, surface_azimuth) honouring the plane's tracking type.""" + tilt = float(plane.surface_tilt if plane.surface_tilt is not None else 30.0) + azimuth = float(plane.surface_azimuth if plane.surface_azimuth is not None else 180.0) + tracking = plane.trackingtype + + if tracking in (None, 0): + return tilt, azimuth + + apparent_zenith = solpos["apparent_zenith"] + solar_azimuth = solpos["azimuth"] + + if tracking == 2: + # Two-axis: the plane always faces the sun. Below the horizon the angles are + # meaningless, so park the plane flat and let the zero irradiance do the rest. + return apparent_zenith.clip(lower=0.0, upper=90.0), solar_azimuth + if tracking == 3: + # Vertical axis: fixed tilt, azimuth follows the sun. + return tilt, solar_azimuth + if tracking in (1, 4, 5): + # Horizontal N-S (1), horizontal E-W (4), inclined N-S (5). + axis_tilt = tilt if tracking == 5 else 0.0 + axis_azimuth = 90.0 if tracking == 4 else 0.0 + tracker = pvlib.tracking.singleaxis( + apparent_zenith=apparent_zenith, + solar_azimuth=solar_azimuth, + axis_tilt=axis_tilt, + axis_azimuth=axis_azimuth, + max_angle=90, + backtrack=False, + ) + return ( + tracker["surface_tilt"].fillna(axis_tilt), + tracker["surface_azimuth"].fillna(axis_azimuth), + ) + + logger.warning( + f"Unsupported trackingtype {tracking} for local pvforecast, treating plane as fixed." + ) + return tilt, azimuth + + # ------------------------------------------------------------------ pv model + + def _plane_power( + self, + plane: Any, + weather: pd.DataFrame, + solpos: pd.DataFrame, + dni_extra: pd.Series, + airmass: pd.Series, + ) -> tuple[pd.Series, pd.Series]: + """Run the pvlib chain for one plane, returning (dc_power_w, ac_power_w).""" + settings = self._settings + + peakpower_kw = plane.peakpower + if peakpower_kw is None: + logger.warning("Plane without peakpower skipped by local pvforecast.") + zero = pd.Series(0.0, index=weather.index) + return zero, zero + pdc0 = float(peakpower_kw) * 1000.0 + + ghi = weather["shortwave_radiation"] + dhi = weather["diffuse_radiation"] + dni = weather["direct_normal_irradiance"] + + # Horizon shading kills the beam component; what is left of the global is the + # diffuse. Do this before transposition so the sky model sees consistent inputs. + if plane.userhorizon: + horizon = self._horizon_elevation(plane.userhorizon, solpos["azimuth"].to_numpy()) + shaded = solpos["apparent_elevation"].to_numpy() < horizon + dni = dni.where(~shaded, 0.0) + ghi = ghi.where(~shaded, dhi) + + surface_tilt, surface_azimuth = self._tracked_orientation(plane, solpos) + albedo = plane.albedo if plane.albedo is not None else settings.albedo + + poa = pvlib.irradiance.get_total_irradiance( + surface_tilt=surface_tilt, + surface_azimuth=surface_azimuth, + solar_zenith=solpos["apparent_zenith"], + solar_azimuth=solpos["azimuth"], + dni=dni, + ghi=ghi, + dhi=dhi, + dni_extra=dni_extra, + airmass=airmass, + albedo=float(albedo), + model=settings.transposition_model, + ) + poa_global = poa["poa_global"].fillna(0.0).clip(lower=0.0) + poa_direct = poa["poa_direct"].fillna(0.0).clip(lower=0.0) + poa_diffuse = poa["poa_diffuse"].fillna(0.0).clip(lower=0.0) + + if settings.apply_iam: + aoi = pvlib.irradiance.aoi( + surface_tilt, surface_azimuth, solpos["apparent_zenith"], solpos["azimuth"] + ) + iam = pvlib.iam.ashrae(aoi).fillna(0.0) + effective_irradiance = poa_direct * iam + poa_diffuse + else: + effective_irradiance = poa_global + + mounting = plane.mountingplace or "free" + temperature_params = pvlib.temperature.TEMPERATURE_MODEL_PARAMETERS["sapm"][ + MOUNTING_TEMPERATURE_MODEL.get(mounting, "open_rack_glass_glass") + ] + temp_cell = pvlib.temperature.sapm_cell( + poa_global=poa_global, + temp_air=weather["temperature_2m"], + wind_speed=weather["wind_speed_10m"], + **temperature_params, + ) + + dc_power = pvlib.pvsystem.pvwatts_dc( + effective_irradiance=effective_irradiance, + temp_cell=temp_cell, + pdc0=pdc0, + gamma_pdc=settings.temperature_coefficient / 100.0, + ) + # `loss` is the PVGIS-style lump of soiling, mismatch, wiring and ageing. + loss = plane.loss if plane.loss is not None else 0.0 + dc_power = (dc_power * (1.0 - float(loss) / 100.0)).fillna(0.0).clip(lower=0.0) + + paco = plane.inverter_paco + eta = settings.inverter_efficiency + if paco is None: + # No inverter rating configured: apply efficiency but do not clip. + ac_power = dc_power * eta + else: + ac_power = pd.Series( + pvlib.inverter.pvwatts( + pdc=dc_power.to_numpy(), + pdc0=float(paco) / eta, + eta_inv_nom=eta, + ), + index=dc_power.index, + ) + ac_power = ac_power.fillna(0.0).clip(lower=0.0) + + return dc_power, ac_power + + async def _forecast_frame(self, data: Any, calibrate: bool = True) -> pd.DataFrame: + """Run the full chain and return a frame with dc/ac power indexed as EOS records. + + Args: + data: The Open-Meteo response. + calibrate: Apply the measurement-fitted correction when it is enabled and + fittable. Pass False to obtain the raw model output, which is what the + fit itself compares against. + """ + settings = self._settings + weather = self._weather_frame(data) + if weather.empty: + return pd.DataFrame(columns=["dc_power", "ac_power"]) + + location = pvlib.location.Location( + latitude=float(self.config.general.latitude), + longitude=float(self.config.general.longitude), + tz="UTC", + altitude=data.get("elevation"), + ) + + # Open-Meteo stamps an interval mean with the interval end, so the sun position + # that produced it sits half an interval earlier. + interval = pd.Timedelta(minutes=settings.resolution_minutes) + solar_times = pd.DatetimeIndex(weather.index) - interval / 2 + + solpos_solar = location.get_solarposition(solar_times) + solpos = solpos_solar.set_axis(weather.index) + dni_extra = pd.Series( + np.asarray(pvlib.irradiance.get_extra_radiation(solar_times), dtype=float), + index=weather.index, + ) + airmass = pd.Series( + location.get_airmass(solar_times, solar_position=solpos_solar)[ + "airmass_relative" + ].to_numpy(), + index=weather.index, + ) + + total_dc = pd.Series(0.0, index=weather.index) + total_ac = pd.Series(0.0, index=weather.index) + for plane in self.config.pvforecast.planes or []: + dc_power, ac_power = self._plane_power(plane, weather, solpos, dni_extra, airmass) + total_dc = total_dc.add(dc_power, fill_value=0.0) + total_ac = total_ac.add(ac_power, fill_value=0.0) + + frame = pd.DataFrame( + { + "dc_power": total_dc, + "ac_power": total_ac, + "solar_elevation": solpos["apparent_elevation"].to_numpy(), + "solar_azimuth": solpos["azimuth"].to_numpy(), + } + ) + # Relabel from Open-Meteo's interval-end stamps to the interval-start stamps + # that EOS records use. + if settings.shift_to_interval_start: + frame.index = pd.DatetimeIndex(frame.index) - interval + + if calibrate: + # The fit needs the raw model to compare against, which is exactly `frame`. + calibration = await self._fit_calibration(frame) + if calibration is not None: + global_factor, factors = calibration + frame = self._apply_calibration( + frame, + factors, + self._installed_ac_capacity_w(), + global_factor=global_factor, + timezone=self.config.general.timezone, + ) + return frame + + # ------------------------------------------------------------ calibration + + def _installed_ac_capacity_w(self) -> float: + """Rough installed AC capacity, used only to threshold near-zero intervals.""" + total = 0.0 + for plane in self.config.pvforecast.planes or []: + if plane.inverter_paco is not None: + total += float(plane.inverter_paco) + elif plane.peakpower is not None: + total += float(plane.peakpower) * 1000.0 + return total + + async def _pv_measurement_window(self) -> Optional[tuple[Any, Any]]: + """First and last timestamp that actually carries PV production readings. + + The measurement store holds every meter, not just the PV ones. Load meters + routinely reach further than the PV meter in both directions, so deriving + the calibration window from the store as a whole would place it where no + PV reading exists - which silently skips calibration or downgrades it to + hourly fitting. + """ + measurement = self.measurement + first = await measurement.min_datetime() + last = await measurement.max_datetime() + if first is None or last is None: + return None + earliest: Any = None + latest: Any = None + for key in self.config.measurement.pv_production_emr_keys or []: + dates, _ = await measurement.key_to_lists( + key=key, + start_datetime=first, + end_datetime=last.add(minutes=1), + ) + if not dates: + continue + if earliest is None or compare_datetimes(dates[0], earliest).lt: + earliest = dates[0] + if latest is None or compare_datetimes(dates[-1], latest).gt: + latest = dates[-1] + if earliest is None or latest is None: + return None + return earliest, latest + + async def _calibration_interval_minutes(self, start: Any, end: Any) -> int: + """Use native forecast slots only when every PV meter resolves them. + + Interpolating an hourly cumulative meter onto quarter hours would create a + perfectly flat, but invented, intrahour profile. Fall back to hourly fitting + until all configured production meters actually provide native slot readings. + """ + native_minutes = self._settings.resolution_minutes + meter_resolutions: list[float] = [] + for key in self.config.measurement.pv_production_emr_keys or []: + dates, _ = await self.measurement.key_to_lists( + key=key, start_datetime=start, end_datetime=end + ) + if len(dates) < 3: + return 60 + deltas = np.asarray( + [ + (dates[index] - dates[index - 1]).total_seconds() / 60.0 + for index in range(1, len(dates)) + if dates[index] > dates[index - 1] + ], + dtype=float, + ) + if deltas.size == 0: + return 60 + meter_resolutions.append(float(np.median(deltas))) + + if not meter_resolutions: + return 60 + if max(meter_resolutions) <= native_minutes * 1.5: + return native_minutes + return 60 + + def _fit_azimuth_factors( + self, + modelled_kwh: np.ndarray, + measured_kwh: np.ndarray, + azimuth: np.ndarray, + global_factor: float, + ) -> np.ndarray: + """Fit an energy-weighted intraday shape while preserving global energy.""" + settings = self._settings + bin_degrees = settings.calibration_azimuth_bin_degrees + if bin_degrees <= 0: + return np.array([global_factor]) + + bin_count = max(1, int(round(360 / bin_degrees))) + bin_index = np.clip((azimuth % 360.0) / (360.0 / bin_count), 0, bin_count - 1).astype(int) + + # Fit the shape as a residual around the independently determined global + # energy factor. Day-level availability filtering has already removed outages; + # energy sums now give productive intervals the influence relevant to the EMS. + # The prior shrinks sparse bins back toward a neutral relative factor of one. + relative_shape = np.ones(bin_count, dtype=float) + prior = settings.calibration_prior_kwh + for bin_number in range(bin_count): + in_bin = bin_index == bin_number + weight = float(modelled_kwh[in_bin].sum()) + if weight <= 0.0: + continue + measured_sum = float(measured_kwh[in_bin].sum()) + raw_shape = measured_sum / (weight * global_factor) + relative_shape[bin_number] = (weight * raw_shape + prior) / (weight + prior) + + factors = np.clip( + global_factor * relative_shape, + settings.calibration_min_factor, + settings.calibration_max_factor, + ) + + # Smooth interpolation changes the exact weighted mean of the bin-centre + # values. Renormalize after interpolation so shape correction cannot silently + # change the global kWh calibration. Re-clipping is iterated to respect bounds. + target_energy = global_factor * float(modelled_kwh.sum()) + for _ in range(8): + scale = self._interpolate_azimuth_factors(azimuth, factors) + corrected_energy = float(np.dot(modelled_kwh, scale)) + if corrected_energy <= 0.0: + break + correction = target_energy / corrected_energy + updated = np.clip( + factors * correction, + settings.calibration_min_factor, + settings.calibration_max_factor, + ) + if np.allclose(updated, factors, rtol=1e-6, atol=1e-8): + factors = updated + break + factors = updated + return factors + + async def _pv_production_total_kwh( + self, start_datetime: Any, end_datetime: Any, interval: Any + ) -> np.ndarray: + """Sum configured cumulative PV meters using main's async energy conversion.""" + totals = None + for key in self.config.measurement.pv_production_emr_keys or []: + values = await self.measurement._energy_from_meter_readings( + key, start_datetime, end_datetime, interval + ) + totals = values if totals is None else totals + values + return np.asarray([] if totals is None else totals, dtype=float) + + async def _fit_calibration(self, frame: pd.DataFrame) -> Optional[tuple[float, np.ndarray]]: + """Fit correction factors from measured PV production against the model. + + The comparison runs on past intervals, where the Open-Meteo rows are analysed + rather than forecast weather. That is deliberate: it isolates the error of the + *PV model* (wrong kWp, soiling, degradation, shading the horizon profile misses) + from the error of the *weather forecast*, and only the former is systematic + enough to correct. + + Returns: + (global_factor, per_azimuth_bin_factors) or None when there is not enough + data to fit anything. + """ + settings = self._settings + if not settings.calibration_enabled: + return None + if not self.config.measurement.pv_production_emr_keys: + logger.info( + "PVForecastAkkudoktorLocal calibration is enabled but " + "measurement.pv_production_emr_keys is not configured - skipping." + ) + return None + + pv_window = await self._pv_measurement_window() + if pv_window is None: + logger.info("PVForecastAkkudoktorLocal calibration: no PV measurements yet - skipping.") + return None + pv_min_datetime, reference_end = pv_window + + lookback_days = settings.calibration_days + if settings.calibration_outage_filter_enabled: + lookback_days = max(lookback_days, settings.calibration_reference_days) + reference_start = reference_end.subtract(days=lookback_days) + interval_minutes = await self._calibration_interval_minutes(reference_start, reference_end) + interval_hours = interval_minutes / 60.0 + interval = to_duration(f"{interval_minutes} minutes") + end = reference_end.start_of("hour") + if interval_minutes < 60: + end = reference_end.start_of("minute").set( + minute=(reference_end.minute // interval_minutes) * interval_minutes + ) + start = end.subtract(days=lookback_days) + if compare_datetimes(start, pv_min_datetime).lt: + start = pv_min_datetime.start_of("minute").add(minutes=interval_minutes) + # The model side only exists for the weather window that was requested. + model_start = to_datetime(frame.index[0].to_pydatetime()) + if compare_datetimes(start, model_start).lt: + start = model_start.start_of("minute").add(minutes=interval_minutes) + if compare_datetimes(start, end).ge: + logger.info( + "PVForecastAkkudoktorLocal calibration: measurement window too short - skipping." + ) + return None + + measured_kwh = np.asarray( + await self._pv_production_total_kwh( + start_datetime=start, end_datetime=end, interval=interval + ), + dtype=float, + ) + if measured_kwh.size == 0 or not np.isfinite(measured_kwh).any(): + logger.info( + "PVForecastAkkudoktorLocal calibration: no usable PV measurements - skipping." + ) + return None + + # Model side on the same grid as the real meter. Mean power is converted to + # interval energy below; hourly meters stay hourly and native 15-minute meters + # retain the shape that matters to the EMS. + samples_frame = ( + frame[["ac_power", "solar_azimuth"]].resample(f"{interval_minutes}min").mean() + ) + grid = pd.date_range( + start=pd.Timestamp(start.in_timezone("UTC").isoformat()), + periods=len(measured_kwh), + freq=f"{interval_minutes}min", + ) + samples_frame = samples_frame.reindex(grid) + modelled_kwh = samples_frame["ac_power"].to_numpy(dtype=float) / 1000.0 * interval_hours + azimuth = samples_frame["solar_azimuth"].to_numpy(dtype=float) + + # Only fit where the model says something meaningful is being produced. Dawn and + # dusk intervals otherwise dominate the ratio with noise. + floor_kwh = max( + 0.02 * self._installed_ac_capacity_w() / 1000.0 * interval_hours, + 0.05 * interval_hours, + ) + usable = ( + np.isfinite(modelled_kwh) + & np.isfinite(measured_kwh) + & np.isfinite(azimuth) + & (modelled_kwh > floor_kwh) + & (measured_kwh >= 0.0) + ) + shape_usable = usable.copy() + + # Calibration represents the available PV potential. A battery or inverter + # outage can make an otherwise healthy plant cover only local demand; those + # intervals must not be learned as a permanent model loss. Detect this at day + # level, because individual cloudy hours are much too noisy for a reliable + # availability decision. + local_days = ( + pd.DatetimeIndex(samples_frame.index) + .tz_convert(self.config.general.timezone) + .normalize() + ) + fit_start = pd.Timestamp( + end.subtract(days=settings.calibration_days) + .in_timezone(self.config.general.timezone) + .isoformat() + ).normalize() + + if settings.calibration_outage_filter_enabled and usable.any(): + samples = pd.DataFrame( + { + "modelled_kwh": modelled_kwh[usable], + "measured_kwh": measured_kwh[usable], + "local_day": local_days[usable], + } + ) + daily = samples.groupby("local_day").agg( + modelled_kwh=("modelled_kwh", "sum"), + measured_kwh=("measured_kwh", "sum"), + usable_intervals=("modelled_kwh", "size"), + ) + daily["ratio"] = daily["measured_kwh"] / daily["modelled_kwh"] + + # Low-yield weather days do not carry enough evidence to call an outage. + # Half an equivalent full-load hour scales naturally with plant size. + minimum_day_kwh = max(0.5 * self._installed_ac_capacity_w() / 1000.0, 1.0) + reference_candidates = daily[ + (daily["modelled_kwh"] >= minimum_day_kwh) & np.isfinite(daily["ratio"]) + ] + + outage_days = pd.DatetimeIndex([]) + reference_ratio = float("nan") + if len(reference_candidates) >= settings.calibration_min_healthy_days: + # The upper quartile is a robust estimate of the available plant level: + # outages and curtailment only pull the ratio down, while a few weather + # outliers cannot dominate it as a maximum would. + reference_ratio = float(reference_candidates["ratio"].quantile(0.75)) + outage_limit = reference_ratio * settings.calibration_outage_threshold + outage_days = pd.DatetimeIndex( + reference_candidates.index[reference_candidates["ratio"] < outage_limit] + ) + + # A cloudy day inside a known low-production block can accidentally + # resemble a healthy ratio because both numerator and denominator are + # small. Bridge a single-day gap between two detected outage days so a + # continuous plant event is not partly admitted into the fit. + ordered_days = pd.DatetimeIndex(daily.index).sort_values() + bridged_days = [] + for index in range(1, len(ordered_days) - 1): + previous_day = ordered_days[index - 1] + day = ordered_days[index] + next_day = ordered_days[index + 1] + if ( + previous_day in outage_days + and next_day in outage_days + and (day.date() - previous_day.date()).days == 1 + and (next_day.date() - day.date()).days == 1 + ): + bridged_days.append(day) + if bridged_days: + outage_days = outage_days.union(pd.DatetimeIndex(bridged_days)).sort_values() + + healthy_days = pd.DatetimeIndex(daily.index).difference(outage_days).sort_values() + recent_healthy_days = healthy_days[healthy_days >= fit_start] + if len(recent_healthy_days) < settings.calibration_min_healthy_days: + fit_days = healthy_days[-settings.calibration_min_healthy_days :] + else: + fit_days = recent_healthy_days + + # The recent healthy window tracks the current energy level. The stable + # intraday signature uses the full healthy reference window, avoiding noisy + # shape factors learned from only a handful of days. + shape_usable &= np.asarray(local_days.isin(healthy_days), dtype=bool) + usable &= np.asarray(local_days.isin(fit_days), dtype=bool) + if len(outage_days) > 0: + day_list = ", ".join(day.strftime("%Y-%m-%d") for day in outage_days) + logger.info( + "PVForecastAkkudoktorLocal calibration: excluded probable outage/" + f"curtailment days [{day_list}] (healthy reference " + f"{reference_ratio:.3f})." + ) + else: + usable &= np.asarray(local_days >= fit_start, dtype=bool) + shape_usable = usable.copy() + + minimum_samples = math.ceil(12 / interval_hours) + if usable.sum() < minimum_samples: + logger.info( + "PVForecastAkkudoktorLocal calibration: only " + f"{int(usable.sum())} usable {interval_minutes}-minute intervals - skipping." + ) + return None + + shape_modelled_kwh = modelled_kwh[shape_usable] + shape_measured_kwh = measured_kwh[shape_usable] + shape_azimuth = azimuth[shape_usable] + modelled_kwh = modelled_kwh[usable] + measured_kwh = measured_kwh[usable] + azimuth = azimuth[usable] + + model_total = float(modelled_kwh.sum()) + if model_total <= 0.0: + return None + global_factor = float( + np.clip( + measured_kwh.sum() / model_total, + settings.calibration_min_factor, + settings.calibration_max_factor, + ) + ) + + factors = self._fit_azimuth_factors( + shape_modelled_kwh, + shape_measured_kwh, + shape_azimuth, + global_factor, + ) + if len(factors) == 1: + logger.info( + f"PVForecastAkkudoktorLocal calibration: global factor {global_factor:.3f} " + f"from {usable.sum()} {interval_minutes}-minute intervals." + ) + return global_factor, factors + + # Report how much of the bias the fit actually removes on its own training window. + before = float(np.abs(modelled_kwh - measured_kwh).mean()) + fitted_scale = self._interpolate_azimuth_factors(azimuth, factors) + after = float(np.abs(modelled_kwh * fitted_scale - measured_kwh).mean()) + logger.info( + f"PVForecastAkkudoktorLocal calibration over {usable.sum()} intervals: global factor " + f"{global_factor:.3f}, {len(factors)} azimuth bins at {interval_minutes}-minute " + "measurement resolution, " + f"MAE {before:.3f} -> {after:.3f} kWh per interval" + ) + return global_factor, factors + + @staticmethod + def _interpolate_azimuth_factors(azimuth: np.ndarray, factors: np.ndarray) -> np.ndarray: + """Interpolate fitted bin-centre factors without a discontinuity at north.""" + bin_count = len(factors) + if bin_count == 1: + return np.full(len(azimuth), factors[0], dtype=float) + + bin_width = 360.0 / bin_count + centers = (np.arange(bin_count, dtype=float) + 0.5) * bin_width + interpolation_azimuths = np.concatenate( + ([centers[-1] - 360.0], centers, [centers[0] + 360.0]) + ) + interpolation_factors = np.concatenate(([factors[-1]], factors, [factors[0]])) + return np.interp( + np.nan_to_num(azimuth % 360.0), + interpolation_azimuths, + interpolation_factors, + ) + + @staticmethod + def _apply_calibration( + frame: pd.DataFrame, + factors: np.ndarray, + ac_cap_w: float, + global_factor: Optional[float] = None, + timezone: Optional[str] = None, + ) -> pd.DataFrame: + """Scale power with a smooth, circular interpolation of azimuth factors.""" + azimuth = frame["solar_azimuth"].to_numpy(dtype=float) + scale = PVForecastAkkudoktorLocal._interpolate_azimuth_factors(azimuth, factors) + + # Preserve the independently fitted kWh correction for every forecast day. + # Azimuth factors may redistribute energy within a day, but cannot change its + # calibrated total (apart from the physical inverter cap applied below). + if global_factor is not None and len(frame) > 0: + ac_power = frame["ac_power"].to_numpy(dtype=float) + group_indices: Iterable[np.ndarray] + if isinstance(frame.index, pd.DatetimeIndex): + day_index = frame.index + if timezone is not None and day_index.tz is not None: + day_index = day_index.tz_convert(timezone) + groups = pd.Series(np.arange(len(frame)), index=frame.index).groupby( + day_index.normalize() + ) + group_indices = (group.to_numpy() for _, group in groups) + else: + group_indices = (np.arange(len(frame)),) + for indices in group_indices: + model_energy = float(ac_power[indices].sum()) + shaped_energy = float(np.dot(ac_power[indices], scale[indices])) + if model_energy > 0.0 and shaped_energy > 0.0: + scale[indices] *= global_factor * model_energy / shaped_energy + + frame = frame.copy() + frame["dc_power"] = frame["dc_power"] * scale + frame["ac_power"] = frame["ac_power"] * scale + if ac_cap_w > 0.0: + # A factor above 1 must not push the plant past its inverters. + frame["ac_power"] = frame["ac_power"].clip(upper=ac_cap_w) + return frame + + # ------------------------------------------------------------------ update + + async def _holds_usable_forecast(self) -> bool: + """Whether at least one finite AC forecast remains at or after the run start.""" + dates, values = await self.key_to_lists( + "pvforecast_ac_power", start_datetime=self.ems_start_datetime + ) + return any( + compare_datetimes(date, self.ems_start_datetime).ge + and value is not None + and math.isfinite(float(value)) + for date, value in zip(dates, values) + ) + + async def _update_local_data(self, force_update: Optional[bool] = False) -> None: + """Compute the PV forecast and store it as PVForecastDataRecord entries.""" + if not self.enabled(): + logger.info("PVForecastAkkudoktorLocal is disabled, skipping update.") + return + + if not self.config.pvforecast.planes: + error_msg = "Requested PV forecast, but no planes configured." + logger.error(f"Configuration error: {error_msg}") + raise ValueError(error_msg) + + try: + request: Callable[..., Any] = self._request_local_forecast + data = await asyncio.to_thread(request, force_update=force_update) + except Exception as exc: + if not await self._holds_usable_forecast(): + # No usable future AC prediction remains; the caller must see + # the source failure. + raise + # A momentary weather-API outage must not fail the whole prediction + # update and take every provider after this one down with it. The + # previous run retains some usable future data. This does not prove + # full horizon coverage: consumers must check their required window. + logger.warning( + "PVForecastAkkudoktorLocal update failed ({}); keeping the forecast from the " + "previous run until {}.", + exc, + await self.max_datetime(), + ) + return + frame = await self._forecast_frame(data) + if frame.empty: + logger.warning("Open-Meteo returned no weather rows for local pvforecast.") + return + + for timestamp, row in frame.iterrows(): + await self.update_value( + to_datetime(pd.Timestamp(str(timestamp)).to_pydatetime()), + { + "pvforecast_dc_power": round(float(row["dc_power"]), 1), + "pvforecast_ac_power": round(float(row["ac_power"]), 1), + }, + ) + + logger.debug( + f"Updated local pvforecast: {len(frame)} records at " + f"{self._settings.resolution_minutes} min over " + f"{len(self.config.pvforecast.planes)} plane(s)." + ) + self.update_datetime = to_datetime(in_timezone=self.config.general.timezone) diff --git a/tests/test_pvforecastakkudoktorlocal.py b/tests/test_pvforecastakkudoktorlocal.py new file mode 100644 index 00000000..58d11afe --- /dev/null +++ b/tests/test_pvforecastakkudoktorlocal.py @@ -0,0 +1,725 @@ +"""Tests for the native (pvlib) PV forecast provider.""" + +from unittest.mock import AsyncMock, Mock, patch + +import numpy as np +import pandas as pd +import pendulum +import pvlib +import pytest +import pytest_asyncio +import requests + +from akkudoktoreos.core.coreabc import get_ems, get_measurement, get_prediction +from akkudoktoreos.prediction.pvforecastakkudoktor import PVForecastAkkudoktor +from akkudoktoreos.prediction.pvforecastakkudoktorlocal import ( + PVForecastAkkudoktorLocal, + PVForecastAkkudoktorLocalCommonSettings, +) + +LATITUDE = 52.52 +LONGITUDE = 13.405 + +# A window that starts well before `START` so the calibration fit has past data. +WINDOW_START = pendulum.datetime(2025, 6, 1, 0, 0, tz="UTC") +WINDOW_END = pendulum.datetime(2025, 6, 20, 0, 0, tz="UTC") +START = pendulum.datetime(2025, 6, 15, 0, 0, tz="UTC") + + +def synthetic_openmeteo(resolution_minutes: int = 15, models: list[str] | None = None) -> dict: + """Build an Open-Meteo-shaped response from a pvlib clear-sky series. + + Open-Meteo stamps an interval mean with the interval END, and that is what the + provider expects, so the values are generated at those stamps directly. + """ + freq = f"{resolution_minutes}min" + index = pd.date_range( + start=WINDOW_START.format("YYYY-MM-DD HH:mm"), + end=WINDOW_END.format("YYYY-MM-DD HH:mm"), + freq=freq, + tz="UTC", + ) + location = pvlib.location.Location(LATITUDE, LONGITUDE, tz="UTC", altitude=37.0) + clearsky = location.get_clearsky(index, model="ineichen") + + block = "minutely_15" if resolution_minutes == 15 else "hourly" + values = { + "shortwave_radiation": clearsky["ghi"].round(1).tolist(), + "diffuse_radiation": clearsky["dhi"].round(1).tolist(), + "direct_normal_irradiance": clearsky["dni"].round(1).tolist(), + "temperature_2m": [20.0] * len(index), + "relative_humidity_2m": [50.0] * len(index), + "wind_speed_10m": [2.0] * len(index), + } + + data: dict = {"elevation": 37.0} + payload = {"time": [t.strftime("%Y-%m-%dT%H:%M") for t in index]} + if models: + # Multi-model requests come back with one suffixed series per member. + for name in models: + for key, series in values.items(): + payload[f"{key}_{name}"] = series + else: + payload.update(values) + data[block] = payload + return data + + +@pytest_asyncio.fixture +async def pvforecast_instance(config_eos): + config_eos.merge_settings_from_dict( + { + "general": {"latitude": LATITUDE, "longitude": LONGITUDE, "timezone": "UTC"}, + "prediction": {"hours": 96, "historic_hours": 48}, + "pvforecast": { + "provider": "PVForecastAkkudoktor", + "planes": [ + { + "surface_tilt": 30.0, + "surface_azimuth": 180.0, + "peakpower": 10.0, + "inverter_paco": 10000, + "loss": 14.0, + } + ], + "akkudoktor": {"backend": "local", "resolution_minutes": 15}, + }, + } + ) + provider = PVForecastAkkudoktor() + get_ems().set_start_datetime(START) + await provider.delete_by_datetime() + await get_measurement().delete_by_datetime() + yield provider + await provider.delete_by_datetime() + await get_measurement().delete_by_datetime() + + +def test_provider_id(pvforecast_instance): + assert PVForecastAkkudoktorLocal.provider_id() == "PVForecastAkkudoktor" + assert pvforecast_instance.enabled() is True + + +@pytest.mark.parametrize("value", [0, 5, 30, 61]) +def test_resolution_must_be_15_or_60(value): + with pytest.raises(ValueError, match="resolution_minutes must be 15 or 60"): + PVForecastAkkudoktorLocalCommonSettings(resolution_minutes=value) + + +def test_invalid_transposition_model(): + with pytest.raises(ValueError, match="Invalid transposition_model"): + PVForecastAkkudoktorLocalCommonSettings(transposition_model="nonsense") + + +@pytest.mark.asyncio +async def test_forecast_frame_is_quarter_hourly_and_plausible(pvforecast_instance): + frame = await pvforecast_instance._forecast_frame(synthetic_openmeteo()) + + assert not frame.empty + deltas = frame.index.to_series().diff().dropna().unique() + assert list(deltas) == [pd.Timedelta(minutes=15)] + + # A 10 kWp south-facing roof under clear June skies: below the inverter cap, + # but a substantial fraction of it. + peak = frame["ac_power"].max() + assert 5000.0 < peak <= 10000.0 + assert (frame["ac_power"] >= 0.0).all() + assert (frame["ac_power"] <= frame["dc_power"] + 1e-6).all() + + # Nights are dark. + midnight = frame.between_time("00:00", "01:00")["ac_power"] + assert midnight.max() == pytest.approx(0.0) + + +@pytest.mark.asyncio +async def test_records_are_shifted_to_interval_start(pvforecast_instance): + """Open-Meteo labels an interval by its end; EOS labels it by its start.""" + shifted = await pvforecast_instance._forecast_frame(synthetic_openmeteo()) + + pvforecast_instance.config.pvforecast.akkudoktor = PVForecastAkkudoktorLocalCommonSettings( + resolution_minutes=15, shift_to_interval_start=False + ) + raw = await pvforecast_instance._forecast_frame(synthetic_openmeteo()) + + assert raw.index[0] - shifted.index[0] == pd.Timedelta(minutes=15) + assert raw["ac_power"].to_numpy() == pytest.approx(shifted["ac_power"].to_numpy()) + + +@pytest.mark.asyncio +async def test_ensemble_members_are_averaged(pvforecast_instance): + """Several models in one request must be averaged, not dropped or duplicated.""" + single = await pvforecast_instance._forecast_frame(synthetic_openmeteo()) + ensemble = await pvforecast_instance._forecast_frame( + synthetic_openmeteo(models=["icon_seamless", "gfs_seamless", "ecmwf_ifs025"]) + ) + + # The synthetic members are identical, so the mean must reproduce the single run. + assert ensemble["ac_power"].to_numpy() == pytest.approx(single["ac_power"].to_numpy()) + + +def test_horizon_elevation_wraps_around_north(): + horizon = PVForecastAkkudoktorLocal._horizon_elevation( + [0.0, 10.0, 20.0, 30.0], np.array([0.0, 90.0, 180.0, 270.0, 359.999]) + ) + assert horizon[:4] == pytest.approx([0.0, 10.0, 20.0, 30.0]) + # Wrapping back to due north interpolates from 30 deg towards 0 deg. + assert horizon[4] == pytest.approx(0.0, abs=0.01) + + +@pytest.mark.asyncio +async def test_horizon_shading_reduces_yield(pvforecast_instance): + baseline = await pvforecast_instance._forecast_frame(synthetic_openmeteo()) + + # A 40 deg wall all around blocks the beam for most of the day. + pvforecast_instance.config.pvforecast.planes[0].userhorizon = [40.0] * 12 + shaded = await pvforecast_instance._forecast_frame(synthetic_openmeteo()) + + assert shaded["ac_power"].sum() < baseline["ac_power"].sum() * 0.9 + assert shaded["ac_power"].min() >= 0.0 + + # The low morning sun (below 40 deg elevation until ~07:00 UTC in June) is behind + # the wall, so its beam is gone entirely and only diffuse is left. + assert ( + shaded["ac_power"].between_time("05:00", "07:00").sum() + < baseline["ac_power"].between_time("05:00", "07:00").sum() * 0.5 + ) + + +@pytest.mark.asyncio +async def test_update_data_writes_records(pvforecast_instance): + with patch.object( + PVForecastAkkudoktorLocal, "_request_local_forecast", return_value=synthetic_openmeteo() + ): + await pvforecast_instance._update_data(force_update=True) + + dates, values = await pvforecast_instance.key_to_lists("pvforecast_ac_power") + assert len(dates) > 0 + assert all(value is not None for value in values) + _, dc_values = await pvforecast_instance.key_to_lists("pvforecast_dc_power") + assert all(value is not None for value in dc_values) + + +async def _feed_measurements( + instance: PVForecastAkkudoktorLocal, frame: pd.DataFrame, bias: float, key: str +) -> None: + """Write cumulative PV meter readings that are `bias` times the modelled power. + + Each test passes its own `key`. `Measurement` is database-backed, so clearing the + in-memory record list would not remove readings another test already stored. + """ + instance.config.measurement.pv_production_emr_keys = [key] + measurement = get_measurement() + + hourly = frame["ac_power"].resample("1h").mean() + hourly = hourly.loc[hourly.index < START] + + cumulative = 0.0 + for timestamp, power_w in hourly.items(): + timestamp = pd.Timestamp(str(timestamp)) + await measurement.update_value( + pendulum.instance(timestamp.to_pydatetime()), key, round(cumulative, 6) + ) + cumulative += float(power_w) * bias / 1000.0 + # Closing reading so the last interval has a difference to work with. + await measurement.update_value( + pendulum.instance(hourly.index[-1].to_pydatetime()).add(hours=1), + key, + round(cumulative, 6), + ) + + +async def _feed_measurements_with_recent_outage( + instance: PVForecastAkkudoktorLocal, + frame: pd.DataFrame, + healthy_bias: float, + outage_bias: float, + outage_days: int, + key: str, +) -> None: + """Write a healthy meter history followed by demand-limited PV production.""" + instance.config.measurement.pv_production_emr_keys = [key] + measurement = get_measurement() + + hourly = frame["ac_power"].resample("1h").mean() + hourly = hourly.loc[hourly.index < START] + outage_start = START.subtract(days=outage_days) + healthy_looking_gap = outage_start.add(days=2).date() + + cumulative = 0.0 + for timestamp, power_w in hourly.items(): + timestamp = pd.Timestamp(str(timestamp)) + await measurement.update_value( + pendulum.instance(timestamp.to_pydatetime()), key, round(cumulative, 6) + ) + in_outage = timestamp >= outage_start and timestamp.date() != healthy_looking_gap + bias = outage_bias if in_outage else healthy_bias + cumulative += float(power_w) * bias / 1000.0 + await measurement.update_value( + pendulum.instance(hourly.index[-1].to_pydatetime()).add(hours=1), + key, + round(cumulative, 6), + ) + + +async def _feed_native_quarter_hour_measurements( + instance: PVForecastAkkudoktorLocal, frame: pd.DataFrame, bias: float, key: str +) -> None: + """Write cumulative PV readings at the provider's native 15-minute cadence.""" + instance.config.measurement.pv_production_emr_keys = [key] + measurement = get_measurement() + slots = frame.loc[frame.index < START, "ac_power"] + + cumulative = 0.0 + for timestamp, power_w in slots.items(): + timestamp = pd.Timestamp(str(timestamp)) + await measurement.update_value( + pendulum.instance(timestamp.to_pydatetime()), key, round(cumulative, 6) + ) + cumulative += float(power_w) * 0.25 * bias / 1000.0 + await measurement.update_value( + pendulum.instance(slots.index[-1].to_pydatetime()).add(minutes=15), + key, + round(cumulative, 6), + ) + + +@pytest.mark.asyncio +async def test_calibration_is_off_by_default(pvforecast_instance): + frame = await pvforecast_instance._forecast_frame(synthetic_openmeteo(), calibrate=False) + assert await pvforecast_instance._fit_calibration(frame) is None + + +@pytest.mark.asyncio +async def test_calibration_skips_without_measurement_keys(pvforecast_instance): + pvforecast_instance.config.pvforecast.akkudoktor = PVForecastAkkudoktorLocalCommonSettings( + calibration_enabled=True + ) + frame = await pvforecast_instance._forecast_frame(synthetic_openmeteo(), calibrate=False) + assert await pvforecast_instance._fit_calibration(frame) is None + + +@pytest.mark.asyncio +async def test_calibration_recovers_a_systematic_bias(pvforecast_instance): + """A plant that consistently delivers 80% of the model must be corrected to 0.8.""" + pvforecast_instance.config.pvforecast.akkudoktor = PVForecastAkkudoktorLocalCommonSettings( + calibration_enabled=True, + calibration_days=14, + calibration_azimuth_bin_degrees=0, + ) + + frame = await pvforecast_instance._forecast_frame(synthetic_openmeteo(), calibrate=False) + await _feed_measurements(pvforecast_instance, frame, bias=0.8, key="pv_bias_emr") + + calibration = await pvforecast_instance._fit_calibration(frame) + assert calibration is not None + global_factor, factors = calibration + assert global_factor == pytest.approx(0.8, abs=0.03) + assert factors == pytest.approx([global_factor]) + + corrected = pvforecast_instance._apply_calibration(frame, factors, 10000.0) + assert corrected["ac_power"].sum() == pytest.approx(frame["ac_power"].sum() * global_factor) + + +@pytest.mark.asyncio +async def test_calibration_excludes_recent_demand_limited_outage(pvforecast_instance, caplog): + """A battery outage must not teach demand-limited PV as available generation.""" + pvforecast_instance.config.pvforecast.akkudoktor = PVForecastAkkudoktorLocalCommonSettings( + calibration_enabled=True, + calibration_days=5, + calibration_reference_days=14, + calibration_azimuth_bin_degrees=0, + calibration_min_factor=0.2, + ) + + frame = await pvforecast_instance._forecast_frame(synthetic_openmeteo(), calibrate=False) + await _feed_measurements_with_recent_outage( + pvforecast_instance, + frame, + healthy_bias=0.8, + outage_bias=0.2, + outage_days=5, + key="pv_demand_limited_emr", + ) + + with caplog.at_level("INFO"): + calibration = await pvforecast_instance._fit_calibration(frame) + + assert calibration is not None + global_factor, factors = calibration + assert global_factor == pytest.approx(0.8, abs=0.03) + assert factors == pytest.approx([global_factor]) + assert "excluded probable outage/curtailment days" in caplog.text + # A single statistically healthy-looking day inside the outage is bridged. + assert "2025-06-12" in caplog.text + + # Filtering only selects training data. Applying a global factor preserves the + # native quarter-hour shape instead of replacing it with hourly bucket values. + corrected = pvforecast_instance._apply_calibration(frame, factors, 10000.0) + producing = frame["ac_power"] > 0.0 + assert corrected.index.to_series().diff().dropna().unique().tolist() == [ + pd.Timedelta(minutes=15) + ] + assert ( + corrected.loc[producing, "ac_power"] / frame.loc[producing, "ac_power"] + ).to_numpy() == pytest.approx(np.full(producing.sum(), global_factor)) + + +@pytest.mark.asyncio +async def test_calibration_uses_real_quarter_hour_measurements(pvforecast_instance, caplog): + """Native meter slots permit shape calibration without inventing intrahour data.""" + pvforecast_instance.config.pvforecast.akkudoktor = PVForecastAkkudoktorLocalCommonSettings( + calibration_enabled=True, + calibration_days=14, + calibration_reference_days=14, + calibration_azimuth_bin_degrees=0, + ) + frame = await pvforecast_instance._forecast_frame(synthetic_openmeteo(), calibrate=False) + await _feed_native_quarter_hour_measurements( + pvforecast_instance, frame, bias=0.8, key="pv_quarter_hour_emr" + ) + + assert ( + await pvforecast_instance._calibration_interval_minutes(START.subtract(days=14), START) + == 15 + ) + with caplog.at_level("INFO"): + calibration = await pvforecast_instance._fit_calibration(frame) + + assert calibration is not None + global_factor, _ = calibration + assert global_factor == pytest.approx(0.8, abs=0.03) + assert "15-minute measurement resolution" in caplog.text or "15-minute intervals" in caplog.text + + +@pytest.mark.asyncio +async def test_calibration_factor_is_clamped(pvforecast_instance): + """A wildly wrong meter must not be allowed to swing the forecast.""" + pvforecast_instance.config.pvforecast.akkudoktor = PVForecastAkkudoktorLocalCommonSettings( + calibration_enabled=True, + calibration_days=14, + calibration_azimuth_bin_degrees=0, + calibration_min_factor=0.9, + calibration_max_factor=1.1, + ) + + frame = await pvforecast_instance._forecast_frame(synthetic_openmeteo(), calibrate=False) + await _feed_measurements(pvforecast_instance, frame, bias=0.2, key="pv_clamp_emr") + + calibration = await pvforecast_instance._fit_calibration(frame) + assert calibration is not None + global_factor, _ = calibration + assert global_factor == pytest.approx(0.9) + + +@pytest.mark.asyncio +async def test_calibration_respects_the_inverter_cap(pvforecast_instance): + frame = await pvforecast_instance._forecast_frame(synthetic_openmeteo(), calibrate=False) + corrected = PVForecastAkkudoktorLocal._apply_calibration(frame, np.array([1.5]), 10000.0) + assert corrected["ac_power"].max() <= 10000.0 + 1e-6 + + +def test_azimuth_calibration_is_interpolated_smoothly(): + """Azimuth correction must not introduce steps into the quarter-hour plan.""" + frame = pd.DataFrame( + { + "solar_azimuth": [44.9, 45.0, 45.1, 134.9, 135.0, 135.1], + "dc_power": [1000.0] * 6, + "ac_power": [1000.0] * 6, + } + ) + corrected = PVForecastAkkudoktorLocal._apply_calibration( + frame, np.array([0.5, 1.0, 1.5, 1.0]), ac_cap_w=2000.0 + ) + + assert corrected.loc[1, "ac_power"] == pytest.approx(500.0) + assert corrected.loc[4, "ac_power"] == pytest.approx(1000.0) + values = corrected["ac_power"].to_numpy(dtype=float) + assert abs(values[2] - values[0]) < 2.0 + assert abs(values[5] - values[3]) < 2.0 + + +def test_azimuth_shape_preserves_each_days_global_energy(): + """The EMS gets a changed shape without a changed daily energy budget.""" + index = pd.date_range("2025-06-01", periods=8, freq="12h", tz="UTC") + frame = pd.DataFrame( + { + "solar_azimuth": [45.0, 225.0, 45.0, 225.0, 45.0, 225.0, 45.0, 225.0], + "dc_power": [100.0, 300.0, 200.0, 200.0, 300.0, 100.0, 150.0, 250.0], + "ac_power": [100.0, 300.0, 200.0, 200.0, 300.0, 100.0, 150.0, 250.0], + }, + index=index, + ) + global_factor = 0.8 + corrected = PVForecastAkkudoktorLocal._apply_calibration( + frame, + np.array([0.5, 1.0, 1.5, 1.0]), + ac_cap_w=10_000.0, + global_factor=global_factor, + timezone="UTC", + ) + + raw_daily = frame["ac_power"].resample("1D").sum() + corrected_daily = corrected["ac_power"].resample("1D").sum() + assert corrected_daily.to_numpy() == pytest.approx(raw_daily.to_numpy() * global_factor) + assert np.std(corrected["ac_power"] / frame["ac_power"]) > 0.01 + + +def test_azimuth_shape_fit_preserves_global_energy(pvforecast_instance): + """Intraday correction must not undo the independently fitted daily kWh.""" + centers = np.arange(22.5, 360.0, 45.0) + azimuth = np.repeat(centers, 20) + modelled_kwh = np.ones(len(azimuth)) + expected_shape = np.repeat([0.8, 0.9, 1.0, 1.1, 1.2, 1.1, 1.0, 0.9], 20) + global_factor = 0.8 + measured_kwh = modelled_kwh * global_factor * expected_shape + + factors = pvforecast_instance._fit_azimuth_factors( + modelled_kwh, measured_kwh, azimuth, global_factor + ) + fitted_scale = pvforecast_instance._interpolate_azimuth_factors(azimuth, factors) + + assert np.std(factors) > 0.01 + assert np.dot(modelled_kwh, fitted_scale) == pytest.approx( + global_factor * modelled_kwh.sum(), rel=1e-6 + ) + + +@pytest.mark.asyncio +async def test_forecast_frame_applies_the_calibration(pvforecast_instance): + """The correction must reach every caller of the chain, not just `_update_data`.""" + pvforecast_instance.config.pvforecast.akkudoktor = PVForecastAkkudoktorLocalCommonSettings( + calibration_enabled=True, + calibration_days=14, + calibration_azimuth_bin_degrees=0, + ) + + raw = await pvforecast_instance._forecast_frame(synthetic_openmeteo(), calibrate=False) + await _feed_measurements(pvforecast_instance, raw, bias=0.8, key="pv_chain_emr") + + calibrated = await pvforecast_instance._forecast_frame(synthetic_openmeteo()) + ratio = calibrated["ac_power"].sum() / raw["ac_power"].sum() + assert ratio == pytest.approx(0.8, abs=0.03) + + +def test_transient_weather_outage_is_retried(pvforecast_instance): + """Open-Meteo answers 503 while rotating model runs; that clears in seconds.""" + error = requests.exceptions.HTTPError("503 Server Error") + error.response = Mock(status_code=503) + good = Mock() + good.url = "https://api.open-meteo.com/v1/forecast" + good.raise_for_status.return_value = None + good.json.return_value = synthetic_openmeteo() + + failing = Mock() + failing.url = good.url + failing.raise_for_status.side_effect = error + + with ( + patch( + "akkudoktoreos.prediction.pvforecastakkudoktorlocal.requests.get", + side_effect=[failing, good], + ) as request, + patch("akkudoktoreos.prediction.pvforecastakkudoktorlocal.time.sleep"), + ): + data = pvforecast_instance._request_local_forecast(force_update=True) + + assert request.call_count == 2 + assert "minutely_15" in data + + +@pytest.mark.asyncio +async def test_weather_outage_keeps_the_previous_forecast(pvforecast_instance): + """One dead weather API must not take the whole prediction update down. + + `PredictionContainer.update_data` re-raises whatever an enabled provider + raises, so every provider after this one would be skipped and the endpoint + would answer 400. This fixture retains future samples from the previous run. + """ + # The stored forecast has to reach past the run start for the fallback to be + # worth anything, so put the run inside the synthetic window. + get_ems().set_start_datetime(WINDOW_START.add(days=1)) + with patch.object( + pvforecast_instance, + "_request_local_forecast", + return_value=synthetic_openmeteo(), + ): + await pvforecast_instance._update_data(force_update=True) + stored = await pvforecast_instance.max_datetime() + assert stored is not None + records_before = len(pvforecast_instance) + + with patch.object( + pvforecast_instance, + "_request_local_forecast", + side_effect=RuntimeError("Failed to fetch weather from Open-Meteo API"), + ): + await pvforecast_instance._update_data(force_update=True) + + assert await pvforecast_instance.max_datetime() == stored + assert len(pvforecast_instance) == records_before + + +@pytest.mark.asyncio +async def test_weather_outage_without_any_forecast_still_fails(pvforecast_instance): + """A cold start has nothing to fall back to, so the caller must hear about it. + + The stored forecast lives in the backing store, not in `records`, so an empty + store is what a cold start actually looks like. + """ + await pvforecast_instance.delete_by_datetime() + with patch.object( + pvforecast_instance, + "_request_local_forecast", + side_effect=RuntimeError("Failed to fetch weather from Open-Meteo API"), + ): + with pytest.raises(RuntimeError, match="Failed to fetch weather"): + await pvforecast_instance._update_data(force_update=True) + + +def test_remote_backend_is_the_compatible_default(): + assert PVForecastAkkudoktorLocalCommonSettings().backend == "remote" + + +def test_calibration_rejects_inverted_bounds(): + with pytest.raises(ValueError, match="must not exceed"): + PVForecastAkkudoktorLocalCommonSettings( + calibration_min_factor=1.5, calibration_max_factor=0.5 + ) + + +def test_existing_provider_id_is_registered_once(pvforecast_instance): + providers = get_prediction().providers + selected = [p for p in providers if p.provider_id() == "PVForecastAkkudoktor"] + assert selected == [pvforecast_instance] + assert all(p.provider_id() != "PVForecastAkkudoktorLocal" for p in providers) + + +@pytest.mark.asyncio +async def test_local_power_contract_for_hourly_and_quarter_hour_consumers(pvforecast_instance): + """Both optimizer cadences receive W, ready for their W-to-Wh conversion.""" + start = START.add(hours=12) + frame = pd.DataFrame( + { + "ac_power": [1000.0, 2000.0, 3000.0, 4000.0], + "dc_power": [1100.0, 2200.0, 3300.0, 4400.0], + }, + index=pd.date_range(start=start, periods=4, freq="15min"), + ) + with ( + patch.object(pvforecast_instance, "_request_local_forecast", return_value={}), + patch.object(pvforecast_instance, "_forecast_frame", new=AsyncMock(return_value=frame)), + patch.object(pvforecast_instance, "_request_forecast") as remote, + ): + await pvforecast_instance._update_data(force_update=True) + remote.assert_not_called() + quarter_power = await get_prediction().key_to_array( + "pvforecast_ac_power", + start_datetime=start, + end_datetime=start.add(hours=1), + interval=pendulum.duration(minutes=15), + ) + hourly_power = await get_prediction().key_to_array( + "pvforecast_ac_power", + start_datetime=start, + end_datetime=start.add(hours=1), + interval=pendulum.duration(hours=1), + ) + assert quarter_power == pytest.approx([1000.0, 2000.0, 3000.0, 4000.0]) + assert hourly_power == pytest.approx([2500.0]) + assert float(sum(quarter_power * 0.25)) == pytest.approx(float(sum(hourly_power))) + + +@pytest.mark.asyncio +@pytest.mark.parametrize("record", ["expired_power", "future_measurement_only"]) +async def test_outage_does_not_accept_unusable_stored_data(pvforecast_instance, record): + if record == "expired_power": + await pvforecast_instance.update_value( + START.subtract(hours=1), "pvforecast_ac_power", 1000.0 + ) + else: + await pvforecast_instance.update_value( + START.add(hours=1), "pvforecastakkudoktor_ac_power_measured", 1000.0 + ) + with patch.object( + pvforecast_instance, "_request_local_forecast", side_effect=RuntimeError("outage") + ): + with pytest.raises(RuntimeError, match="outage"): + await pvforecast_instance._update_data(force_update=True) + + +@pytest.mark.parametrize("through_file_migration", [False, True]) +def test_feature_settings_migrate_without_loss_or_input_mutation( + config_eos, through_file_migration +): + import copy + from akkudoktoreos.config.configmigrate import migrate_config_data + from akkudoktoreos.prediction.pvforecast import PVForecastCommonSettings + + source = { + "provider": "PVForecastAkkudoktorLocal", + "provider_settings": { + "PVForecastAkkudoktorLocal": { + "calibration_enabled": True, + "calibration_days": 21, + "calibration_min_factor": 1.7, + "calibration_max_factor": 2.1, + } + }, + } + before = copy.deepcopy(source) + if through_file_migration: + migrated = migrate_config_data({"pvforecast": source}).pvforecast + else: + migrated = PVForecastCommonSettings.model_validate(source) + assert source == before + assert migrated.provider == "PVForecastAkkudoktor" + assert migrated.akkudoktor.backend == "local" + assert migrated.akkudoktor.calibration_enabled + assert migrated.akkudoktor.calibration_days == 21 + assert migrated.akkudoktor.calibration_min_factor == 1.7 + assert migrated.akkudoktor.calibration_max_factor == 2.1 + + +def test_current_pv_settings_take_precedence_during_migration(config_eos): + from akkudoktoreos.config.configmigrate import migrate_config_data + + migrated = migrate_config_data( + { + "pvforecast": { + "provider": "PVForecastAkkudoktorLocal", + "provider_settings": { + "PVForecastAkkudoktorLocal": { + "calibration_enabled": True, + "calibration_days": 21, + } + }, + "akkudoktor": {"backend": "remote", "calibration_days": 7}, + } + } + ).pvforecast + assert migrated.provider == "PVForecastAkkudoktor" + assert migrated.akkudoktor.backend == "remote" + assert migrated.akkudoktor.calibration_days == 7 + assert migrated.akkudoktor.calibration_enabled + + +@pytest.mark.parametrize( + "invalid", + [ + {"resolution_minutes": 5}, + {"calibration_min_factor": 2.0, "calibration_max_factor": 1.0}, + "invalid", + ], +) +def test_file_migration_rejects_invalid_local_settings(config_eos, invalid): + from akkudoktoreos.config.configmigrate import migrate_config_data + + with pytest.raises(ValueError): + migrate_config_data( + { + "pvforecast": { + "provider": "PVForecastAkkudoktorLocal", + "provider_settings": {"PVForecastAkkudoktorLocal": invalid}, + } + } + )