Issue #1292 mete grids from imod5

docs/api/changelog.rst

@@ -9,6 +9,12 @@ The format is based on `Keep a Changelog`_, and this project adheres to
 [Unreleased]
 ------------
 
+Added
+~~~~~
+
+- :class:`imod.msw.MeteoGridCopy` to copy existing `mete_grid.inp` files, so
+  ASCII grids in large existing meteo databases do not have to be read.
+
 Changed
 ~~~~~~~
 

imod/msw/__init__.py

@@ -9,7 +9,7 @@
     InitialConditionsSavedState,
 )
 from imod.msw.landuse import LanduseOptions
-from imod.msw.meteo_grid import MeteoGrid
+from imod.msw.meteo_grid import MeteoGrid, MeteoGridCopy
 from imod.msw.meteo_mapping import EvapotranspirationMapping, PrecipitationMapping
 from imod.msw.model import MetaSwapModel
 from imod.msw.output_control import TimeOutputControl, VariableOutputControl

imod/msw/meteo_grid.py

@@ -1,6 +1,7 @@
 import csv
 from pathlib import Path
-from typing import Optional, Union
+from shutil import copyfile
+from typing import Optional, Union, cast
 
 import numpy as np
 import pandas as pd
@@ -11,6 +12,9 @@
 from imod.msw.pkgbase import MetaSwapPackage
 from imod.msw.regrid.regrid_schemes import MeteoGridRegridMethod
 from imod.msw.timeutil import to_metaswap_timeformat
+from imod.msw.utilities.common import find_in_file_list
+from imod.typing import Imod5DataDict
+from imod.util.regrid_method_type import EmptyRegridMethod, RegridMethodType
 
 
 class MeteoGrid(MetaSwapPackage, IRegridPackage):
@@ -188,3 +192,38 @@ def _pkgcheck(self):
                     f"Received excess dims {excess_dims} in {self.__class__} for "
                     f"{varname}, please provide data with {allowed_dims}"
                 )
+
+
+class MeteoGridCopy(MetaSwapPackage, IRegridPackage):
+    """
+    Class to copy existing ``mete_grid.inp``, which contains the meteorological
+    grid data. Next to a MeteoGridCopy instance, instances of
+    PrecipitationMapping and EvapotranspirationMapping are required as well to
+    specify meteorological information to MetaSWAP.
+
+    Parameters
+    ----------
+    path: Path to mete_grid.inp file
+    """
+
+    _file_name = "mete_grid.inp"
+    _meteo_dirname = "meteo_grids"
+
+    _regrid_method: RegridMethodType = EmptyRegridMethod()
+
+    def __init__(self, path: Path | str):
+        super().__init__()
+        self.dataset["path"] = path
+
+    def write(self, directory: Path | str, *args):
+        directory = Path(directory)
+        path_metegrid = Path(str(self.dataset["path"].values[()]))
+        new_path = directory / self._file_name
+        copyfile(path_metegrid, new_path)
+
+    @classmethod
+    def from_imod5_data(cls, imod5_data: Imod5DataDict) -> "MeteoGridCopy":
+        paths = cast(list[str], imod5_data["extra"]["paths"])
+        filepath = find_in_file_list(cls._file_name, paths)
+
+        return cls(filepath)

imod/msw/meteo_mapping.py

@@ -1,18 +1,68 @@
 from copy import deepcopy
-from typing import Any, Optional, TextIO
+from pathlib import Path
+from textwrap import dedent
+from typing import Any, Optional, TextIO, cast
 
 import numpy as np
 import pandas as pd
 import xarray as xr
 
+import imod
 from imod.mf6.utilities.regrid import RegridderWeightsCache
 from imod.msw.fixed_format import VariableMetaData
 from imod.msw.pkgbase import MetaSwapPackage
+from imod.msw.utilities.common import find_in_file_list
 from imod.prepare import common
-from imod.typing import GridDataArray, IntArray
+from imod.typing import GridDataArray, Imod5DataDict, IntArray
 from imod.util.regrid_method_type import RegridMethodType
 
 
+def _is_parsable_and_existing_path(potential_path: str, mete_grid_path: Path) -> bool:
+    """
+    mete_grid.inp can contain values like "0.", which are converted to float by
+    MetaSWAP. String is converted to path and checked if existing path.
+    """
+    try:
+        float(potential_path)
+        return False
+    except ValueError:
+        # Resolve paths relative to mete_grid.inp path.
+        path = mete_grid_path / ".." / Path(potential_path)
+        return path.is_file()
+
+
+def open_first_meteo_grid(mete_grid_path: str | Path, column_nr: int) -> xr.DataArray:
+    """
+    Find and open first meteo grid path in mete_grid.inp. This grid is enough to
+    generate meteomappings. There can be floats before in the column which
+    should be skipped.
+    """
+    if column_nr not in [2, 3]:
+        raise ValueError("Column nr should be 2 or 3")
+
+    mete_grid_path = Path(mete_grid_path)
+    with open(mete_grid_path, "r") as f:
+        lines = f.readlines()
+
+    potential_paths = [line.split(",")[column_nr].replace('"', "") for line in lines]
+    for potential_path in potential_paths:
+        if _is_parsable_and_existing_path(potential_path, mete_grid_path):
+            resolved_path = mete_grid_path / ".." / Path(potential_path)
+            return imod.rasterio.open(resolved_path)
+
+    error_message = dedent(f"""    
+    Did not find parsable path to existing .ASC file in column {column_nr}. Got
+    values (printing first 10): {potential_paths[:10]}.""")
+
+    raise ValueError(error_message)
+
+
+def open_first_meteo_grid_from_imod5_data(imod5_data: Imod5DataDict, column_nr: int):
+    paths = cast(list[str], imod5_data["extra"]["paths"])
+    metegrid_path = find_in_file_list("mete_grid.inp", paths)
+    return open_first_meteo_grid(metegrid_path, column_nr=column_nr)
+
+
 class MeteoMapping(MetaSwapPackage):
     """
     This class provides common methods for creating mappings between
@@ -125,6 +175,18 @@ def __init__(
         super().__init__()
         self.meteo = precipitation
 
+    @classmethod
+    def from_imod5_data(cls, imod5_data: Imod5DataDict) -> "PrecipitationMapping":
+        """
+        Construct precipitation mapping from imod5 data. Opens first ascii grid
+        in mete_grid.inp, which is used to construct mappings to svats. The
+        grids should not change in dimension over time. No checks are done
+        whether cells switch from inactive to active or vice versa.
+        """
+        column_nr = 2
+        meteo_grid = open_first_meteo_grid_from_imod5_data(imod5_data, column_nr)
+        return cls(meteo_grid)
+
 
 class EvapotranspirationMapping(MeteoMapping):
     """
@@ -156,3 +218,15 @@ def __init__(
     ):
         super().__init__()
         self.meteo = evapotranspiration
+
+    @classmethod
+    def from_imod5_data(cls, imod5_data: Imod5DataDict) -> "EvapotranspirationMapping":
+        """
+        Construct evapotranspiration mapping from imod5 data. Opens first ascii
+        grid in mete_grid.inp, which is used to construct mappings to svats. The
+        grids should not change in dimension over time. No checks are done
+        whether cells switch from inactive to active or vice versa.
+        """
+        column_nr = 3
+        meteo_grid = open_first_meteo_grid_from_imod5_data(imod5_data, column_nr)
+        return cls(meteo_grid)

imod/msw/utilities/common.py

@@ -1,6 +1,15 @@
+from pathlib import Path
+
 from imod.typing import GridDataArray
 from imod.typing.grid import concat
 
 
 def concat_imod5(arg1: GridDataArray, arg2: GridDataArray) -> GridDataArray:
     return concat([arg1, arg2], dim="subunit").assign_coords(subunit=[0, 1])
+
+
+def find_in_file_list(filename: str, paths: list[str]) -> str:
+    for file in paths:
+        if filename == Path(file[0]).name.lower():
+            return file[0]
+    raise ValueError(f"could not find {filename} in list of paths: {paths}")