From 17f268c2128ea2da39a689308c248065e7abd459 Mon Sep 17 00:00:00 2001
From: nurkholi <afid.nur-kholis@ufz.de>
Date: Tue, 14 Jan 2025 14:58:38 +0100
Subject: [PATCH] modify component and constants

---
 src/finam_mhm/component.py | 127 ++++++++++++++++++++++++++++++-------
 src/finam_mhm/constants.py |   3 +
 2 files changed, 107 insertions(+), 23 deletions(-)

diff --git a/src/finam_mhm/component.py b/src/finam_mhm/component.py
index 631f2e6..fa5dda7 100644
--- a/src/finam_mhm/component.py
+++ b/src/finam_mhm/component.py
@@ -20,14 +20,19 @@ from .constants import (
     OUTPUT_META,
 )
 
-
 def _horizon_name(name, horizon):
     return name + "_L" + str(horizon).zfill(2)
 
 
 def _get_grid_name(var):
     grid_name = var.split("_")[0]
-    return "L1" if grid_name == "METEO" else grid_name
+    if grid_name == "METEO":
+        name = "L1"
+    elif grid_name == "GW":
+        name = "L0"
+    else:
+        name = grid_name
+    return name
 
 
 def _get_var_name(var):
@@ -39,6 +44,10 @@ def _get_meteo_inputs(inputs):
         _get_var_name(var).lower(): var for var in inputs if var.startswith("METEO")
     }
 
+def _get_GWvar_inputs(inputs):
+    return {
+        _get_var_name(var).lower(): var for var in inputs if var.startswith("GW")
+    }
 
 class MHM(fm.TimeComponent):
     """
@@ -80,8 +89,9 @@ class MHM(fm.TimeComponent):
         cwd=".",
         input_names=None,
         meteo_timestep=None,
+        gw_timestep=None,
         ignore_input_grid=False,
-    ):
+     ):
         super().__init__()
         self.gridspec = {}
         self.no_data = None
@@ -106,8 +116,15 @@ class MHM(fm.TimeComponent):
         self.cwd = cwd  # needed for @fm.tools.execute_in_cwd
         # mHM always has hourly stepping
         self.step = timedelta(hours=1)
+        
+        #>>> Meteo
         self.meteo_timestep = meteo_timestep
         self.meteo_inputs = _get_meteo_inputs(self.INPUT_NAMES)
+        
+        #>>> GW
+        self.gw_timestep = gw_timestep
+        self.gw_inputs = _get_GWvar_inputs(self.INPUT_NAMES)
+
         self.ignore_input_grid = ignore_input_grid
 
         if self.meteo_inputs and self.meteo_timestep not in HOURS_TO_TIMESTEP:
@@ -131,13 +148,15 @@ class MHM(fm.TimeComponent):
     def _initialize(self):
         # only show errors
         mhm.model.set_verbosity(level=1)
-        # configure coupling
+
+        #>> configure coupling
         if self.meteo_inputs:
             kwargs = {f"meteo_expect_{var}": True for var in self.meteo_inputs}
             kwargs["couple_case"] = 1
             kwargs["meteo_timestep"] = self.meteo_timestep
             mhm.model.config_coupling(**kwargs)
-        # init
+
+        #>> init
         mhm.model.init(
             namelist_mhm=self.namelist_mhm,
             namelist_mhm_param=self.namelist_mhm_param,
@@ -145,13 +164,18 @@ class MHM(fm.TimeComponent):
             namelist_mrm_output=self.namelist_mrm_output,
             cwd=".",
         )
-        # disable file output of mHM
+
+        #>> disable file output of mHM
         mhm.model.disable_output()
+
         mhm.run.prepare()
-        # only one domain possible
+
+        #>> only one domain possible
         mhm.run.prepare_domain()
+
         self.number_of_horizons = mhm.get.number_of_horizons()
-        # prepare outputs
+
+        #>> prepare outputs
         self.OUTPUT_NAMES = list(OUTPUT_META)
         if self.mrm_active:
             self.OUTPUT_NAMES += list(MRM_OUTPUT_META)
@@ -166,25 +190,28 @@ class MHM(fm.TimeComponent):
             for var in OUTPUT_CALC_HORIZONS_META
             for horizon in self.horizons
         ]
-        # get start time
+
+        #>> get start time
         year, month, day, hour = mhm.run.current_time()
         self.time = datetime(year=year, month=month, day=max(day, 0), hour=max(hour, 0))
         # first time step compensate by negative values in mHM
         if day < 0 or hour < 0:
             self.time += timedelta(days=min(day, 0), hours=min(hour, 0))
 
-        # store Grid specifications
+        #>> store Grid specifications
         # get grid info l0 (swap rows/cols to get "ij" indexing)
         nrows, ncols, __, xll, yll, cell_size, no_data = mhm.get.l0_domain_info()
         self.no_data = no_data
         self.gridspec["L0"] = fm.EsriGrid(
             ncols=ncols, nrows=nrows, cellsize=cell_size, xllcorner=xll, yllcorner=yll
         )
+
         # get grid info l1 (swap rows/cols to get "ij" indexing)
         nrows, ncols, __, xll, yll, cell_size, no_data = mhm.get.l1_domain_info()
         self.gridspec["L1"] = fm.EsriGrid(
             ncols=ncols, nrows=nrows, cellsize=cell_size, xllcorner=xll, yllcorner=yll
         )
+
         if self.mrm_active:
             # get grid info l11 (swap rows/cols to get "ij" indexing)
             nrows, ncols, __, xll, yll, cell_size, no_data = mhm.get.l11_domain_info()
@@ -195,6 +222,7 @@ class MHM(fm.TimeComponent):
                 xllcorner=xll,
                 yllcorner=yll,
             )
+
         # get grid info l2 (swap rows/cols to get "ij" indexing)
         nrows, ncols, __, xll, yll, cell_size, no_data = mhm.get.l2_domain_info()
         self.gridspec["L2"] = fm.EsriGrid(
@@ -263,18 +291,45 @@ class MHM(fm.TimeComponent):
                     _FillValue=self.no_data,
                     **n_meta,
                 )
+
+        #..............................
+        # MODIFICATION
+        #..............................
         for var in self.INPUT_NAMES:
-            grid_name = _get_grid_name(var)
-            self.inputs.add(
-                name=var,
-                time=self.time,
-                grid=None if self.ignore_input_grid else self.gridspec[grid_name],
-                missing_value=self.no_data,
-                _FillValue=self.no_data,
-                units=INPUT_UNITS[var].format(
-                    ts=HOURS_TO_TIMESTEP[self.meteo_timestep]
-                ),
-            )
+            if var == "GW_DEPTH":
+                grid_name = _get_grid_name(var)
+                self.inputs.add(
+                    name=var,
+                    time=self.time,
+                    grid=None if self.ignore_input_grid else self.gridspec[grid_name],
+                    missing_value=self.no_data,
+                    _FillValue=self.no_data,
+                    units=INPUT_UNITS[var],
+                )
+
+            elif var == "GW_BASEFLOW":
+                grid_name = _get_grid_name(var)
+                self.inputs.add(
+                    name=var,
+                    time=self.time,
+                    grid=None if self.ignore_input_grid else self.gridspec[grid_name],
+                    missing_value=self.no_data,
+                    _FillValue=self.no_data,
+                    units=INPUT_UNITS[var],
+                )
+
+            elif var == "METEO_PRE":
+                grid_name = _get_grid_name(var)
+                self.inputs.add(
+                    name=var,
+                    time=self.time,
+                    grid=None if self.ignore_input_grid else self.gridspec[grid_name],
+                    missing_value=self.no_data,
+                    _FillValue=self.no_data,
+                    units=INPUT_UNITS[var].format(
+                        ts=HOURS_TO_TIMESTEP[self.meteo_timestep]),
+                )
+
         self.create_connector()
 
     def _connect(self, start_time):
@@ -303,9 +358,9 @@ class MHM(fm.TimeComponent):
         # Don't run further than mHM can
         if mhm.run.finished():
             return
-        # set meteo data
+        
+        #>> set meteo data
         if self.meteo_inputs:
-            # every hour or every 24 hours
             if self.time.hour % self.meteo_timestep == 0:
                 kwargs = {
                     var: self.inputs[name].pull_data(self.next_time)[0].magnitude
@@ -313,8 +368,34 @@ class MHM(fm.TimeComponent):
                 }
                 kwargs["time"] = self.time
                 mhm.set_meteo(**kwargs)
+
+        #>> Set GW Coupling 
+        if self.gw_inputs:
+            if self.gw_timestep == 24 : # hours
+                if self.time.hour % self.gw_timestep == 0:
+                    print(f"time={self.time}")
+                    kwargs = {
+                        var: self.inputs[name].pull_data(self.next_time)[0].magnitude
+                        for var, name in self.gw_inputs.items()
+                    }
+                    kwargs["time"] = self.time
+                    mhm.set_gw_coupling(**kwargs)
+
+            if self.gw_timestep == 30 : # days
+                if self.time.day == 1 and self.time.hour == 0: 
+                    print(f"time={self.time}")
+                    kwargs = {
+                        var: self.inputs[name].pull_data(self.next_time)[0].magnitude
+                        for var, name in self.gw_inputs.items()
+                    }
+                    kwargs["time"] = self.time
+                    mhm.set_gw_coupling(**kwargs)
+
+
         # run mhm
         mhm.run.do_time_step()
+
+        
         # update time
         year, month, day, hour = mhm.run.current_time()
         self.time = datetime(year=year, month=month, day=day, hour=hour)
diff --git a/src/finam_mhm/constants.py b/src/finam_mhm/constants.py
index bcb9c68..fb63be9 100644
--- a/src/finam_mhm/constants.py
+++ b/src/finam_mhm/constants.py
@@ -163,6 +163,9 @@ INPUT_UNITS = {
     "METEO_SSRD": "W m-2",
     "METEO_STRD": "W m-2",
     "METEO_TANN": "degC",
+    "GW_DEPTH": "m",
+    "GW_BASEFLOW": "mm / d",
+
 }
 """units of the available inputs in mHM."""
 
-- 
GitLab