Add automated spot-market price threshold optimizer and DWD weather forecast modules

battery-control/README

@@ -1,4 +1,4 @@
-To run these script you need paho-mqtt, serial und urllib3 for python
+To run these scripts you need paho-mqtt, serial, urllib3 and requests for python
 
 Carefully check all scripts for sanity, e.g. voltage and power limits.
 
@@ -15,3 +15,67 @@ For backward compatibility, a single expression can still be configured via:
 - result_topic: topic where the computed result is published
 
 Broker connection uses "broker.address" and optionally "broker.port" and "broker.keepalive" from config.json.
+
+dwd_forecast.py fetches DWD MOSMIX_L weather forecasts once per hour and publishes solar and
+wind estimates via MQTT so that price_optimizer.py can decide whether PV production will cover
+the battery's charging needs the following day.
+Configure it in config.json under "dwd_forecast":
+- station_id          DWD MOSMIX station ID (e.g. "P755" for Augsburg)
+                      Find stations at https://wettwarn.de/mosmix/mosmix.html
+- panel_area_m2       total PV panel area in m²
+- panel_efficiency    module efficiency as a fraction (default 0.18)
+
+Published topics:
+  /forecast/pv_energy_today     – estimated PV yield for today (kWh)
+  /forecast/pv_energy_tomorrow  – estimated PV yield for tomorrow (kWh)
+  /forecast/avg_wind_ms         – mean wind speed over the next 24 h (m/s)
+  /forecast/wind_factor         – wind speed normalised to [0, 1]
+  /forecast/season              – "summer" (Apr–Sep) or "winter" (Oct–Mar)
+
+price_optimizer.py computes optimal charge and discharge price thresholds and publishes them so
+that netzero.py acts on them automatically.  It supports two storage configurations:
+  - EV (55 kWh default) when pyPlc/fsm_state == "CurrentDemand"
+  - Stationary battery (6 kWh default) at all other times
+
+Algorithm:
+  - Discharge threshold: determined by how cheaply the battery can be recharged after
+    discharging, accounting for domestic consumption (daily_consumption_kwh, EV excluded):
+    1. Estimate current stored energy = capacity - energy_still_needed.
+    2. Domestic load will drain the battery regardless; subtract it from stored energy.
+    3. Net kWh that must be bought back = max(0, stored - domestic_drain).
+    4. Find the cheapest N slots to cover net_recharge_kwh; the current (active) slot is
+       included because it is still available for recharging.  Use the Nth slot price as
+       the discharge threshold.
+    When tomorrow's prices are very cheap (e.g. −1 €/MWh) the threshold drops
+    accordingly, so energy is discharged freely today.
+  - Charge threshold (summer, sufficient PV forecast): set to -9999 to suppress grid charging
+    entirely – the battery will be filled by PV the following day.
+  - Charge threshold (winter or low PV): find the cheapest N hours needed to fill the storage
+    deficit and set the threshold to that boundary price, reduced by a wind discount
+    (wind_discount_eur_mwh * wind_factor) to capture wind-driven cheap hours.
+    The threshold is further capped at discharge_thresh × roundtrip_efficiency: if the
+    cheapest available charge price exceeds this limit the stored energy could not be
+    recouped after losses, so grid charging is suppressed automatically.
+
+MQTT resilience: the script uses paho loop_start() so the network loop runs in a background
+thread independent of the main sleep.  reconnect_delay_set() triggers automatic reconnection
+after any broker disconnect; on_connect re-subscribes and recalculates thresholds immediately.
+
+Configure it in config.json under "price_optimizer":
+- ev_capacity_kwh            EV battery capacity in kWh (default 55)
+- stat_capacity_kwh          stationary battery capacity in kWh (default 6)
+- stat_max_charge_power_w    stationary battery max charge power in W (default 1800)
+- wind_discount_eur_mwh      threshold discount per unit of wind_factor (default 20 EUR/MWh)
+- daily_consumption_kwh      expected household consumption per day excluding EV charging
+                             (default 10 kWh); used to net out domestic load from the
+                             recharge cost calculation for the discharge threshold
+- roundtrip_efficiency       storage round-trip efficiency as a fraction (default 0.8);
+                             grid charging is not triggered when its cost exceeds
+                             discharge_thresh × efficiency
+
+Published topics (retained):
+  /grid/chargethresh    – charge stationary battery when spot price is below this (EUR/MWh)
+  /grid/evchargethresh  – charge EV when spot price is below this (EUR/MWh)
+  /grid/dischargethresh – discharge storage only when spot price is above this (EUR/MWh)
+
+Recommended service startup order: spotmarket → dwd_forecast → price_optimizer → netzero

battery-control/config.json

@@ -54,6 +54,19 @@
 		"apitype": "tibber",
 		"apikey": "your tibber key"
 	},
+	"dwd_forecast": {
+		"station_id": "P755",
+		"panel_area_m2": 20.0,
+		"panel_efficiency": 0.18
+	},
+	"price_optimizer": {
+		"ev_capacity_kwh": 55,
+		"stat_capacity_kwh": 6,
+		"stat_max_charge_power_w": 1800,
+		"wind_discount_eur_mwh": 20,
+		"daily_consumption_kwh": 10,
+		"roundtrip_efficiency": 0.8
+	},
 	"mqtt_math": {
 		"aliases": [ {
 				"topic": "deye-dummycloud/3851343667/grid/active_power_w",

battery-control/dwd_forecast.py

@@ -0,0 +1,151 @@
+#!/usr/bin/python3
+#
+# Fetches DWD MOSMIX_L weather forecasts for a configurable station and
+# publishes estimates to MQTT so that the price_optimizer module can decide
+# whether PV production will cover the battery's charging needs.
+#
+# Published topics
+#   /forecast/pv_energy_today     – estimated PV yield for today (kWh)
+#   /forecast/pv_energy_tomorrow  – estimated PV yield for tomorrow (kWh)
+#   /forecast/avg_wind_ms         – mean wind speed over the next 24 h (m/s)
+#   /forecast/wind_factor         – wind speed normalised to [0, 1] (10 m/s = 1)
+#   /forecast/season              – "summer" (Apr–Sep) or "winter" (Oct–Mar)
+#
+# Config key: "dwd_forecast" in config.json
+#   station_id        DWD station identifier, e.g. "P755"
+#   panel_area_m2     total PV panel area in m²
+#   panel_efficiency  PV panel efficiency (default 0.18)
+
+import io
+import json
+import time
+import zipfile
+import xml.etree.ElementTree as ET
+from datetime import date, datetime, timedelta, timezone
+
+import paho.mqtt.client as mqtt
+import requests
+
+DWD_NS = 'https://opendata.dwd.de/weather/lib/pointforecast_dwd_extension_V1_0.xsd'
+DWD_BASE = (
+    'https://opendata.dwd.de/weather/local_forecasts/mos/MOSMIX_L'
+    '/single_stations/{sid}/kml/MOSMIX_L_LATEST_{sid}.kmz'
+)
+
+
+def fetch_dwd_forecast(station_id):
+    """Download and parse the MOSMIX_L KMZ file for *station_id*."""
+    url = DWD_BASE.format(sid=station_id)
+    r = requests.get(url, timeout=60)
+    r.raise_for_status()
+    with zipfile.ZipFile(io.BytesIO(r.content)) as z:
+        kml_name = next(n for n in z.namelist() if n.endswith('.kml'))
+        return ET.fromstring(z.read(kml_name))
+
+
+def parse_timestamps(root):
+    """Return a list of UTC-aware datetimes for each forecast step."""
+    steps_elem = root.find(f'.//{{{DWD_NS}}}ForecastTimeSteps')
+    if steps_elem is None:
+        return []
+    return [
+        datetime.fromisoformat(ts.text.replace('Z', '+00:00'))
+        for ts in steps_elem.findall(f'{{{DWD_NS}}}TimeStep')
+    ]
+
+
+def parse_parameter(root, name):
+    """Return a list of float values for DWD parameter *name*, or None where missing."""
+    attr_key = f'{{{DWD_NS}}}elementName'
+    for elem in root.iter(f'{{{DWD_NS}}}Forecast'):
+        if elem.get(attr_key) == name:
+            val_elem = elem.find(f'{{{DWD_NS}}}value')
+            if val_elem is not None and val_elem.text:
+                result = []
+                for token in val_elem.text.split():
+                    try:
+                        v = float(token)
+                        result.append(None if v != v else v)  # guard against NaN
+                    except ValueError:
+                        result.append(None)
+                return result
+    return []
+
+
+def estimate_pv_energy(timestamps, radiation, target_date, panel_area_m2, efficiency):
+    """Estimate total PV yield for *target_date* in kWh.
+
+    DWD Rad1h is global horizontal irradiance integrated over the hour in kJ/m².
+    Energy [kWh] = Rad1h [kJ/m²] * area [m²] * η / 3600
+    """
+    total = 0.0
+    for ts, rad in zip(timestamps, radiation):
+        if ts.date() == target_date and rad is not None and rad > 0:
+            total += rad * panel_area_m2 * efficiency / 3600.0
+    return total
+
+
+def estimate_avg_wind_speed(timestamps, wind, hours=24):
+    """Average wind speed (m/s) over the next *hours* hours."""
+    now = datetime.now(timezone.utc)
+    cutoff = now + timedelta(hours=hours)
+    values = [w for ts, w in zip(timestamps, wind)
+              if now <= ts <= cutoff and w is not None]
+    return sum(values) / len(values) if values else 0.0
+
+
+def is_summer():
+    return datetime.now().month in range(4, 10)
+
+
+with open('config.json') as f:
+    config = json.load(f)
+
+dwd_cfg = config.get('dwd_forecast', {})
+station_id = dwd_cfg.get('station_id', 'P755')
+panel_area = float(dwd_cfg.get('panel_area_m2', 20.0))
+efficiency = float(dwd_cfg.get('panel_efficiency', 0.18))
+
+client = mqtt.Client(client_id='dwd_forecast')
+client.connect(
+    config['broker']['address'],
+    config['broker'].get('port', 1883),
+    config['broker'].get('keepalive', 60),
+)
+
+last_fetch = 0
+
+while True:
+    if time.time() - last_fetch >= 3600:
+        try:
+            root = fetch_dwd_forecast(station_id)
+            timestamps = parse_timestamps(root)
+            radiation = parse_parameter(root, 'Rad1h')
+            wind = parse_parameter(root, 'FF')
+
+            today = datetime.now(timezone.utc).date()
+            tomorrow = today + timedelta(days=1)
+
+            pv_today = estimate_pv_energy(timestamps, radiation, today, panel_area, efficiency)
+            pv_tomorrow = estimate_pv_energy(timestamps, radiation, tomorrow, panel_area, efficiency)
+            avg_wind = estimate_avg_wind_speed(timestamps, wind)
+            wind_factor = min(1.0, avg_wind / 10.0)
+            season = 'summer' if is_summer() else 'winter'
+
+            client.publish('/forecast/pv_energy_today', round(pv_today, 2), retain=True)
+            client.publish('/forecast/pv_energy_tomorrow', round(pv_tomorrow, 2), retain=True)
+            client.publish('/forecast/avg_wind_ms', round(avg_wind, 1), retain=True)
+            client.publish('/forecast/wind_factor', round(wind_factor, 3), retain=True)
+            client.publish('/forecast/season', season, retain=True)
+
+            print(
+                f'DWD forecast updated: PV today={pv_today:.1f} kWh, '
+                f'tomorrow={pv_tomorrow:.1f} kWh, '
+                f'wind={avg_wind:.1f} m/s ({wind_factor:.2f}), season={season}'
+            )
+            last_fetch = time.time()
+        except Exception as e:
+            print(f'DWD forecast fetch failed: {e}')
+
+    client.loop(timeout=0.1)
+    time.sleep(60)