From ddd87e16ef5c312c3d3e9c8ed21d4d9213326fb3 Mon Sep 17 00:00:00 2001 From: JoerivanEngelen Date: Mon, 2 Dec 2024 15:46:26 +0100 Subject: [PATCH 01/50] Add parser function for parasim.inp --- imod/msw/utilities/parse.py | 32 +++++++++++++++++++ .../test_msw/test_utilities/test_parse.py | 28 ++++++++++++++++ 2 files changed, 60 insertions(+) create mode 100644 imod/msw/utilities/parse.py create mode 100644 imod/tests/test_msw/test_utilities/test_parse.py diff --git a/imod/msw/utilities/parse.py b/imod/msw/utilities/parse.py new file mode 100644 index 000000000..a6cdba323 --- /dev/null +++ b/imod/msw/utilities/parse.py @@ -0,0 +1,32 @@ +from pathlib import Path + + +def _try_parsing_string_to_value(s: str): + """ + Convert string to value: + + "1" -> 1 + "1.0" -> 1.0 + "a" -> "a" + """ + try: + value = int(s) + except ValueError: + try: + value = float(s) + except ValueError: + value = s + return value + + +def read_para_sim(file: Path | str) -> dict[str, str]: + with open(file, "r") as f: + lines = f.readlines() + out = {} + for line in lines: + ll = line[0 : line.find("!")].split("=") + if len(ll) > 1: + key = ll[0].strip() + value = _try_parsing_string_to_value(ll[1].strip()) + out[key] = value + return out diff --git a/imod/tests/test_msw/test_utilities/test_parse.py b/imod/tests/test_msw/test_utilities/test_parse.py new file mode 100644 index 000000000..7d945db6e --- /dev/null +++ b/imod/tests/test_msw/test_utilities/test_parse.py @@ -0,0 +1,28 @@ +from pytest_cases import parametrize_with_cases + +from imod.msw.utilities.parse import _try_parsing_string_to_value + + +class ParseCases: + def case_int(self): + return "1", int + + def case_float(self): + return "1.0", float + + def case_string(self): + return "a", str + + def case_aterisk(self): + return "*", str + + def case_exclamation(self): + return "!", str + + +@parametrize_with_cases(["s", "expected_type"], cases=ParseCases) +def test_try_parsing_string_to_value(s, expected_type): + # Act + parsed = _try_parsing_string_to_value(s) + # Assert + assert type(parsed) is expected_type From 8734a5da02bc820d5516e9db93a867ccba962cfc Mon Sep 17 00:00:00 2001 From: JoerivanEngelen Date: Mon, 2 Dec 2024 16:40:06 +0100 Subject: [PATCH 02/50] Start adding a MetaSwapModel.from_imod5_data method --- imod/msw/model.py | 54 +++++++++++++++++++++++++++++++++++++++++++---- 1 file changed, 50 insertions(+), 4 deletions(-) diff --git a/imod/msw/model.py b/imod/msw/model.py index f349021bb..2939b86e5 100644 --- a/imod/msw/model.py +++ b/imod/msw/model.py @@ -1,13 +1,16 @@ import collections from copy import copy from pathlib import Path -from typing import Optional, Tuple, Union +from typing import Any, Optional, Tuple, Union import jinja2 import numpy as np from imod.mf6.dis import StructuredDiscretization from imod.mf6.mf6_wel_adapter import Mf6Wel +from imod.mf6.regrid.regrid_schemes import ( + RegridMethodType, +) from imod.mf6.utilities.regrid import RegridderWeightsCache from imod.msw.coupler_mapping import CouplerMapping from imod.msw.grid_data import GridData @@ -20,11 +23,15 @@ 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.output_control import TimeOutputControl +from imod.msw.ponding import Ponding from imod.msw.timeutil import to_metaswap_timeformat +from imod.msw.utilities.common import find_in_file_list +from imod.msw.utilities.parse import read_para_sim from imod.msw.vegetation import AnnualCropFactors +from imod.typing import Imod5DataDict from imod.util.regrid_method_type import RegridderType REQUIRED_PACKAGES = ( @@ -88,6 +95,7 @@ class MetaSwapModel(Model): ---------- unsaturated_database: Path-like or str Path to the MetaSWAP soil physical database folder. + settings: dict """ _pkg_id = "model" @@ -99,10 +107,16 @@ class MetaSwapModel(Model): "{%endfor%}" ) - def __init__(self, unsaturated_database): + def __init__( + self, unsaturated_database: Path | str, settings: dict[str, Any] = None + ): super().__init__() - self.simulation_settings = copy(DEFAULT_SETTINGS) + if settings is None: + self.simulation_settings = copy(DEFAULT_SETTINGS) + else: + self.simulation_settings = settings + self.simulation_settings["unsa_svat_path"] = ( self._render_unsaturated_database_path(unsaturated_database) ) @@ -282,3 +296,35 @@ def regrid_like( regridded_model[mod2svat_name] = CouplerMapping() return regridded_model + + @classmethod + def from_imod5_data( + cls, + imod5_data: Imod5DataDict, + target_dis: StructuredDiscretization, + regridder_types: Optional[RegridMethodType] = None, + regrid_cache: RegridderWeightsCache = RegridderWeightsCache(), + ): + extra_paths = imod5_data["extra"]["paths"] + path_to_parasim = find_in_file_list("para_sim.inp", extra_paths) + parasim_settings = read_para_sim(path_to_parasim) + + model = cls(parasim_settings["unsat_svat_path"], parasim_settings) + + model["grid"] = GridData.from_imod5_data( + imod5_data, target_dis, regridder_types, regrid_cache + ) + model["infiltration"] = Infiltration.from_imod5_data(imod5_data) + model["ponding"] = Ponding.from_imod5_data(imod5_data) + # model["sprinkling"] = Sprinkling.from_imod5_data(imod5_data) + model["meteo_grid"] = MeteoGridCopy.from_imod5_data(imod5_data) + model["prec_mapping"] = PrecipitationMapping.from_imod5_data(imod5_data) + model["evt_mapping"] = EvapotranspirationMapping.from_imod5_data(imod5_data) + model["idf_mapping"] = IdfMapping(model["grid"]["area"], np.nan) + + model["coupling"] = CouplerMapping() + + # TODO: + # Add CopyFiles package here + + return model From 435641c0d361ecc0ca047ad6bf1d25774921d453 Mon Sep 17 00:00:00 2001 From: JoerivanEngelen Date: Tue, 3 Dec 2024 17:28:34 +0100 Subject: [PATCH 03/50] Enable importing Sprinkling from_imod5_data --- imod/msw/model.py | 7 +++---- 1 file changed, 3 insertions(+), 4 deletions(-) diff --git a/imod/msw/model.py b/imod/msw/model.py index 2939b86e5..54dbf2693 100644 --- a/imod/msw/model.py +++ b/imod/msw/model.py @@ -27,6 +27,7 @@ from imod.msw.meteo_mapping import EvapotranspirationMapping, PrecipitationMapping from imod.msw.output_control import TimeOutputControl from imod.msw.ponding import Ponding +from imod.msw.sprinkling import Sprinkling from imod.msw.timeutil import to_metaswap_timeformat from imod.msw.utilities.common import find_in_file_list from imod.msw.utilities.parse import read_para_sim @@ -316,15 +317,13 @@ def from_imod5_data( ) model["infiltration"] = Infiltration.from_imod5_data(imod5_data) model["ponding"] = Ponding.from_imod5_data(imod5_data) - # model["sprinkling"] = Sprinkling.from_imod5_data(imod5_data) + model["sprinkling"] = Sprinkling.from_imod5_data(imod5_data) model["meteo_grid"] = MeteoGridCopy.from_imod5_data(imod5_data) model["prec_mapping"] = PrecipitationMapping.from_imod5_data(imod5_data) model["evt_mapping"] = EvapotranspirationMapping.from_imod5_data(imod5_data) model["idf_mapping"] = IdfMapping(model["grid"]["area"], np.nan) model["coupling"] = CouplerMapping() - - # TODO: - # Add CopyFiles package here + #model["extra_files"] = CopyFiles.from_imod5_data(imod5_data) return model From c0ad931e2cc94d82d0435212ab17478c75d65021 Mon Sep 17 00:00:00 2001 From: JoerivanEngelen Date: Wed, 4 Dec 2024 09:55:27 +0100 Subject: [PATCH 04/50] Fix input for IdfMapping --- imod/msw/model.py | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/imod/msw/model.py b/imod/msw/model.py index 54dbf2693..0da5e4dd9 100644 --- a/imod/msw/model.py +++ b/imod/msw/model.py @@ -321,7 +321,8 @@ def from_imod5_data( model["meteo_grid"] = MeteoGridCopy.from_imod5_data(imod5_data) model["prec_mapping"] = PrecipitationMapping.from_imod5_data(imod5_data) model["evt_mapping"] = EvapotranspirationMapping.from_imod5_data(imod5_data) - model["idf_mapping"] = IdfMapping(model["grid"]["area"], np.nan) + area = model["grid"]["area"].isel(subunit=0, drop=True) + model["idf_mapping"] = IdfMapping(area, -9999.0) model["coupling"] = CouplerMapping() #model["extra_files"] = CopyFiles.from_imod5_data(imod5_data) From 213bfc409b5ba941e4b9dd40c31c8165f5c4b12b Mon Sep 17 00:00:00 2001 From: JoerivanEngelen Date: Wed, 4 Dec 2024 13:42:37 +0100 Subject: [PATCH 05/50] Fix some mypy issues --- imod/msw/model.py | 17 ++++++++++------- imod/msw/utilities/parse.py | 13 ++++++++----- .../tests/test_msw/test_utilities/test_parse.py | 4 ++-- 3 files changed, 20 insertions(+), 14 deletions(-) diff --git a/imod/msw/model.py b/imod/msw/model.py index 0da5e4dd9..486c3c7b8 100644 --- a/imod/msw/model.py +++ b/imod/msw/model.py @@ -1,7 +1,7 @@ import collections from copy import copy from pathlib import Path -from typing import Any, Optional, Tuple, Union +from typing import Any, Optional, Tuple, Union, cast import jinja2 import numpy as np @@ -12,6 +12,7 @@ RegridMethodType, ) from imod.mf6.utilities.regrid import RegridderWeightsCache +from imod.msw.copy_files import FileCopier from imod.msw.coupler_mapping import CouplerMapping from imod.msw.grid_data import GridData from imod.msw.idf_mapping import IdfMapping @@ -109,7 +110,9 @@ class MetaSwapModel(Model): ) def __init__( - self, unsaturated_database: Path | str, settings: dict[str, Any] = None + self, + unsaturated_database: Path | str, + settings: Optional[dict[str, Any]] = None, ): super().__init__() @@ -273,7 +276,7 @@ def regrid_like( regrid_context: Optional[RegridderWeightsCache] = None, regridder_types: Optional[dict[str, Tuple[RegridderType, str]]] = None, ) -> "MetaSwapModel": - unsat_database = self.simulation_settings["unsa_svat_path"] + unsat_database = cast(str, self.simulation_settings["unsa_svat_path"]) regridded_model = MetaSwapModel(unsat_database) target_grid = mf6_regridded_dis["idomain"] @@ -306,11 +309,12 @@ def from_imod5_data( regridder_types: Optional[RegridMethodType] = None, regrid_cache: RegridderWeightsCache = RegridderWeightsCache(), ): - extra_paths = imod5_data["extra"]["paths"] + extra_paths = cast(list[str], imod5_data["extra"]["paths"]) path_to_parasim = find_in_file_list("para_sim.inp", extra_paths) parasim_settings = read_para_sim(path_to_parasim) + unsat_svat_path = cast(str, parasim_settings["unsat_svat_path"]) - model = cls(parasim_settings["unsat_svat_path"], parasim_settings) + model = cls(unsat_svat_path, parasim_settings) model["grid"] = GridData.from_imod5_data( imod5_data, target_dis, regridder_types, regrid_cache @@ -323,8 +327,7 @@ def from_imod5_data( model["evt_mapping"] = EvapotranspirationMapping.from_imod5_data(imod5_data) area = model["grid"]["area"].isel(subunit=0, drop=True) model["idf_mapping"] = IdfMapping(area, -9999.0) - model["coupling"] = CouplerMapping() - #model["extra_files"] = CopyFiles.from_imod5_data(imod5_data) + model["extra_files"] = FileCopier.from_imod5_data(imod5_data) return model diff --git a/imod/msw/utilities/parse.py b/imod/msw/utilities/parse.py index a6cdba323..6177775e6 100644 --- a/imod/msw/utilities/parse.py +++ b/imod/msw/utilities/parse.py @@ -1,16 +1,19 @@ from pathlib import Path +from typing import TypeAlias +ScalarType: TypeAlias = str | float | int -def _try_parsing_string_to_value(s: str): + +def _try_parsing_string_to_number(s: str) -> ScalarType: """ - Convert string to value: + Convert string to number: "1" -> 1 "1.0" -> 1.0 "a" -> "a" """ try: - value = int(s) + value: ScalarType = int(s) except ValueError: try: value = float(s) @@ -19,7 +22,7 @@ def _try_parsing_string_to_value(s: str): return value -def read_para_sim(file: Path | str) -> dict[str, str]: +def read_para_sim(file: Path | str) -> dict[str, ScalarType]: with open(file, "r") as f: lines = f.readlines() out = {} @@ -27,6 +30,6 @@ def read_para_sim(file: Path | str) -> dict[str, str]: ll = line[0 : line.find("!")].split("=") if len(ll) > 1: key = ll[0].strip() - value = _try_parsing_string_to_value(ll[1].strip()) + value = _try_parsing_string_to_number(ll[1].strip()) out[key] = value return out diff --git a/imod/tests/test_msw/test_utilities/test_parse.py b/imod/tests/test_msw/test_utilities/test_parse.py index 7d945db6e..a4293252a 100644 --- a/imod/tests/test_msw/test_utilities/test_parse.py +++ b/imod/tests/test_msw/test_utilities/test_parse.py @@ -1,6 +1,6 @@ from pytest_cases import parametrize_with_cases -from imod.msw.utilities.parse import _try_parsing_string_to_value +from imod.msw.utilities.parse import _try_parsing_string_to_number class ParseCases: @@ -23,6 +23,6 @@ def case_exclamation(self): @parametrize_with_cases(["s", "expected_type"], cases=ParseCases) def test_try_parsing_string_to_value(s, expected_type): # Act - parsed = _try_parsing_string_to_value(s) + parsed = _try_parsing_string_to_number(s) # Assert assert type(parsed) is expected_type From 82b593ffcb26eab2550553c278fbf1e2b41b4891 Mon Sep 17 00:00:00 2001 From: JoerivanEngelen Date: Thu, 5 Dec 2024 11:32:59 +0100 Subject: [PATCH 06/50] Add docstring, add timeoutputcontrol (so that startime can be set), and prepare for adding ScalingFactors --- imod/msw/model.py | 31 ++++++++++++++++++++++++++++--- 1 file changed, 28 insertions(+), 3 deletions(-) diff --git a/imod/msw/model.py b/imod/msw/model.py index 486c3c7b8..1a9c038bb 100644 --- a/imod/msw/model.py +++ b/imod/msw/model.py @@ -1,6 +1,7 @@ import collections from copy import copy from pathlib import Path +from datetime import datetime from typing import Any, Optional, Tuple, Union, cast import jinja2 @@ -28,6 +29,7 @@ from imod.msw.meteo_mapping import EvapotranspirationMapping, PrecipitationMapping from imod.msw.output_control import TimeOutputControl from imod.msw.ponding import Ponding +from imod.msw.scaling_factors import ScalingFactors from imod.msw.sprinkling import Sprinkling from imod.msw.timeutil import to_metaswap_timeformat from imod.msw.utilities.common import find_in_file_list @@ -306,9 +308,29 @@ def from_imod5_data( cls, imod5_data: Imod5DataDict, target_dis: StructuredDiscretization, - regridder_types: Optional[RegridMethodType] = None, - regrid_cache: RegridderWeightsCache = RegridderWeightsCache(), + times: list[datetime], ): + """ + Construct a MetaSWAP model from iMOD5 data in the CAP package, loaded + with the :func:`imod.formats.prj.open_projectfile_data` function. + + Parameters + ---------- + imod5_data: dict + Dictionary with iMOD5 data. This can be constructed from the + :func:`imod.formats.prj.open_projectfile_data` method. + target_dis: imod.mf6.StructuredDiscretization + Target discretization, cells where MODLOW6 is inactive will be + inactive in MetaSWAP as well. + times: list[datetime] + List of datetimes, will be used to set the output control times. + Is also used to infer the starttime of the simulation. + + Returns + ------- + MetaSwapModel + MetaSWAP model imported from imod5 data. + """ extra_paths = cast(list[str], imod5_data["extra"]["paths"]) path_to_parasim = find_in_file_list("para_sim.inp", extra_paths) parasim_settings = read_para_sim(path_to_parasim) @@ -317,7 +339,7 @@ def from_imod5_data( model = cls(unsat_svat_path, parasim_settings) model["grid"] = GridData.from_imod5_data( - imod5_data, target_dis, regridder_types, regrid_cache + imod5_data, target_dis ) model["infiltration"] = Infiltration.from_imod5_data(imod5_data) model["ponding"] = Ponding.from_imod5_data(imod5_data) @@ -325,9 +347,12 @@ def from_imod5_data( model["meteo_grid"] = MeteoGridCopy.from_imod5_data(imod5_data) model["prec_mapping"] = PrecipitationMapping.from_imod5_data(imod5_data) model["evt_mapping"] = EvapotranspirationMapping.from_imod5_data(imod5_data) + #model["scaling_factor"] = ScalingFactors.from_imod5_data(imod5_data) area = model["grid"]["area"].isel(subunit=0, drop=True) model["idf_mapping"] = IdfMapping(area, -9999.0) model["coupling"] = CouplerMapping() model["extra_files"] = FileCopier.from_imod5_data(imod5_data) + model["time_oc"] = TimeOutputControl(times) + return model From 7dfb56cb59f0111f13bfb928f6fdfe62fbee6a1a Mon Sep 17 00:00:00 2001 From: JoerivanEngelen Date: Thu, 5 Dec 2024 11:33:13 +0100 Subject: [PATCH 07/50] format --- imod/msw/model.py | 18 ++++++------------ 1 file changed, 6 insertions(+), 12 deletions(-) diff --git a/imod/msw/model.py b/imod/msw/model.py index 1a9c038bb..a6708065b 100644 --- a/imod/msw/model.py +++ b/imod/msw/model.py @@ -1,7 +1,7 @@ import collections from copy import copy -from pathlib import Path from datetime import datetime +from pathlib import Path from typing import Any, Optional, Tuple, Union, cast import jinja2 @@ -9,9 +9,6 @@ from imod.mf6.dis import StructuredDiscretization from imod.mf6.mf6_wel_adapter import Mf6Wel -from imod.mf6.regrid.regrid_schemes import ( - RegridMethodType, -) from imod.mf6.utilities.regrid import RegridderWeightsCache from imod.msw.copy_files import FileCopier from imod.msw.coupler_mapping import CouplerMapping @@ -29,7 +26,6 @@ from imod.msw.meteo_mapping import EvapotranspirationMapping, PrecipitationMapping from imod.msw.output_control import TimeOutputControl from imod.msw.ponding import Ponding -from imod.msw.scaling_factors import ScalingFactors from imod.msw.sprinkling import Sprinkling from imod.msw.timeutil import to_metaswap_timeformat from imod.msw.utilities.common import find_in_file_list @@ -313,7 +309,7 @@ def from_imod5_data( """ Construct a MetaSWAP model from iMOD5 data in the CAP package, loaded with the :func:`imod.formats.prj.open_projectfile_data` function. - + Parameters ---------- imod5_data: dict @@ -324,8 +320,8 @@ def from_imod5_data( inactive in MetaSWAP as well. times: list[datetime] List of datetimes, will be used to set the output control times. - Is also used to infer the starttime of the simulation. - + Is also used to infer the starttime of the simulation. + Returns ------- MetaSwapModel @@ -338,16 +334,14 @@ def from_imod5_data( model = cls(unsat_svat_path, parasim_settings) - model["grid"] = GridData.from_imod5_data( - imod5_data, target_dis - ) + model["grid"] = GridData.from_imod5_data(imod5_data, target_dis) model["infiltration"] = Infiltration.from_imod5_data(imod5_data) model["ponding"] = Ponding.from_imod5_data(imod5_data) model["sprinkling"] = Sprinkling.from_imod5_data(imod5_data) model["meteo_grid"] = MeteoGridCopy.from_imod5_data(imod5_data) model["prec_mapping"] = PrecipitationMapping.from_imod5_data(imod5_data) model["evt_mapping"] = EvapotranspirationMapping.from_imod5_data(imod5_data) - #model["scaling_factor"] = ScalingFactors.from_imod5_data(imod5_data) + # model["scaling_factor"] = ScalingFactors.from_imod5_data(imod5_data) area = model["grid"]["area"].isel(subunit=0, drop=True) model["idf_mapping"] = IdfMapping(area, -9999.0) model["coupling"] = CouplerMapping() From 6f15441cb34baf24a3a5a81bf9b4532bf2a73f15 Mon Sep 17 00:00:00 2001 From: JoerivanEngelen Date: Thu, 5 Dec 2024 13:14:33 +0100 Subject: [PATCH 08/50] Fix keyerror --- imod/msw/model.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/imod/msw/model.py b/imod/msw/model.py index a6708065b..b8c658b99 100644 --- a/imod/msw/model.py +++ b/imod/msw/model.py @@ -330,9 +330,9 @@ def from_imod5_data( extra_paths = cast(list[str], imod5_data["extra"]["paths"]) path_to_parasim = find_in_file_list("para_sim.inp", extra_paths) parasim_settings = read_para_sim(path_to_parasim) - unsat_svat_path = cast(str, parasim_settings["unsat_svat_path"]) + unsa_svat_path = cast(str, parasim_settings["unsa_svat_path"]) - model = cls(unsat_svat_path, parasim_settings) + model = cls(unsa_svat_path, parasim_settings) model["grid"] = GridData.from_imod5_data(imod5_data, target_dis) model["infiltration"] = Infiltration.from_imod5_data(imod5_data) From 19e48f3259d7b00880f4a131597d4c5e211dc1b8 Mon Sep 17 00:00:00 2001 From: JoerivanEngelen Date: Thu, 5 Dec 2024 14:32:49 +0100 Subject: [PATCH 09/50] Take boundary dataarray, as this is more likely to be the grid. --- imod/msw/utilities/imod5_converter.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/imod/msw/utilities/imod5_converter.py b/imod/msw/utilities/imod5_converter.py index 70f88c982..040f55c4a 100644 --- a/imod/msw/utilities/imod5_converter.py +++ b/imod/msw/utilities/imod5_converter.py @@ -10,7 +10,7 @@ def get_cell_area_from_imod5_data( imod5_cap: GridDataDict, ) -> GridDataArray: # area's per type of svats - mf6_area = get_cell_area(imod5_cap["wetted_area"]) + mf6_area = get_cell_area(imod5_cap["boundary"]) wetted_area = imod5_cap["wetted_area"] urban_area = imod5_cap["urban_area"] rural_area = mf6_area - (wetted_area + urban_area) From 78f9c61b373e0a6decdc53fc4d7af0ddc659c63d Mon Sep 17 00:00:00 2001 From: JoerivanEngelen Date: Thu, 5 Dec 2024 14:34:07 +0100 Subject: [PATCH 10/50] Drop layer coord and broadcast scalars to grid --- imod/msw/model.py | 46 ++++++++++++++++++++++++++++++++++------------ 1 file changed, 34 insertions(+), 12 deletions(-) diff --git a/imod/msw/model.py b/imod/msw/model.py index b8c658b99..5c741c2d1 100644 --- a/imod/msw/model.py +++ b/imod/msw/model.py @@ -26,12 +26,13 @@ from imod.msw.meteo_mapping import EvapotranspirationMapping, PrecipitationMapping from imod.msw.output_control import TimeOutputControl from imod.msw.ponding import Ponding +from imod.msw.scaling_factors import ScalingFactors from imod.msw.sprinkling import Sprinkling from imod.msw.timeutil import to_metaswap_timeformat from imod.msw.utilities.common import find_in_file_list from imod.msw.utilities.parse import read_para_sim from imod.msw.vegetation import AnnualCropFactors -from imod.typing import Imod5DataDict +from imod.typing import Imod5DataDict, SelSettingsType from imod.util.regrid_method_type import RegridderType REQUIRED_PACKAGES = ( @@ -76,6 +77,12 @@ "clocktime": 0, } +_DROP_LAYER_KWARGS: SelSettingsType = { + "layer": 0, + "drop": True, + "missing_dims": "ignore", +} + class Model(collections.UserDict): def __setitem__(self, key, value): @@ -331,21 +338,36 @@ def from_imod5_data( path_to_parasim = find_in_file_list("para_sim.inp", extra_paths) parasim_settings = read_para_sim(path_to_parasim) unsa_svat_path = cast(str, parasim_settings["unsa_svat_path"]) - + # Drop layer coord + imod5_cap_no_layer: Imod5DataDict = { + "cap": { + key: da.isel(**_DROP_LAYER_KWARGS).compute() + for key, da in imod5_data["cap"].items() + } + } model = cls(unsa_svat_path, parasim_settings) - - model["grid"] = GridData.from_imod5_data(imod5_data, target_dis) - model["infiltration"] = Infiltration.from_imod5_data(imod5_data) - model["ponding"] = Ponding.from_imod5_data(imod5_data) - model["sprinkling"] = Sprinkling.from_imod5_data(imod5_data) - model["meteo_grid"] = MeteoGridCopy.from_imod5_data(imod5_data) - model["prec_mapping"] = PrecipitationMapping.from_imod5_data(imod5_data) - model["evt_mapping"] = EvapotranspirationMapping.from_imod5_data(imod5_data) - # model["scaling_factor"] = ScalingFactors.from_imod5_data(imod5_data) + model["grid"] = GridData.from_imod5_data(imod5_cap_no_layer, target_dis) + active = model["grid"].dataset["active"] + # Mask grid cells and broadcast scalars to grid + imod5_cap_masked: Imod5DataDict = { + "cap": { + key: da.where(active) for key, da in imod5_cap_no_layer["cap"].items() + }, + "extra": {"paths": extra_paths}, + } + model["infiltration"] = Infiltration.from_imod5_data(imod5_cap_masked) + model["ponding"] = Ponding.from_imod5_data(imod5_cap_masked) + model["sprinkling"] = Sprinkling.from_imod5_data(imod5_cap_masked) + model["meteo_grid"] = MeteoGridCopy.from_imod5_data(imod5_cap_masked) + model["prec_mapping"] = PrecipitationMapping.from_imod5_data(imod5_cap_masked) + model["evt_mapping"] = EvapotranspirationMapping.from_imod5_data( + imod5_cap_masked + ) + model["scaling_factor"] = ScalingFactors.from_imod5_data(imod5_cap_masked) area = model["grid"]["area"].isel(subunit=0, drop=True) model["idf_mapping"] = IdfMapping(area, -9999.0) model["coupling"] = CouplerMapping() - model["extra_files"] = FileCopier.from_imod5_data(imod5_data) + model["extra_files"] = FileCopier.from_imod5_data(imod5_cap_masked) model["time_oc"] = TimeOutputControl(times) From e5acfb3314c057fd41ae7652fd104a1167d31734 Mon Sep 17 00:00:00 2001 From: JoerivanEngelen Date: Thu, 5 Dec 2024 16:30:27 +0100 Subject: [PATCH 11/50] Add masking utilities --- imod/msw/utilities/mask.py | 41 ++++++++++++++++++++++++++++++++++++++ 1 file changed, 41 insertions(+) create mode 100644 imod/msw/utilities/mask.py diff --git a/imod/msw/utilities/mask.py b/imod/msw/utilities/mask.py new file mode 100644 index 000000000..246127fe5 --- /dev/null +++ b/imod/msw/utilities/mask.py @@ -0,0 +1,41 @@ +import numbers +from dataclasses import dataclass + +from imod.msw.pkgbase import MetaSwapPackage +from imod.typing import GridDataArray, GridDataDict + + +@dataclass +class MetaSwapActive: + all: GridDataArray + per_subunit: GridDataArray + + +def mask_and_broadcast_grid_data( + package: MetaSwapPackage, + grid_data: GridDataDict, + msw_active: MetaSwapActive + ): + """ + Mask and broadcast grid data, carefully mask per subunit if variable needs + to contain subunits. + """ + + return { + key: ( + _mask_spatial_var(grid, msw_active.per_subunit) + if key in package._with_subunit + else _mask_spatial_var(grid, msw_active.all) + ) + for key, grid in grid_data.items() + } + +def _mask_spatial_var(self, da: GridDataArray, active: GridDataArray) -> GridDataArray: + if issubclass(da.dtype.type, numbers.Integral): + return da.where(active, other=0) + elif issubclass(da.dtype.type, numbers.Real): + return da.where(active) + else: + raise TypeError( + f"Expected dtype float or integer. Received instead: {da.dtype}" + ) \ No newline at end of file From a6878f06c6e2c506d9b79ae33feb30c0dc1352ef Mon Sep 17 00:00:00 2001 From: JoerivanEngelen Date: Thu, 5 Dec 2024 16:30:50 +0100 Subject: [PATCH 12/50] Call masking utilities --- imod/msw/grid_data.py | 21 ++++++--------------- 1 file changed, 6 insertions(+), 15 deletions(-) diff --git a/imod/msw/grid_data.py b/imod/msw/grid_data.py index ad7299206..d2090c7cc 100644 --- a/imod/msw/grid_data.py +++ b/imod/msw/grid_data.py @@ -18,6 +18,7 @@ get_rootzone_depth_from_imod5_data, is_msw_active_cell, ) +from imod.msw.utilities.mask import mask_and_broadcast_grid_data, MetaSwapActive from imod.typing import GridDataDict from imod.util.spatial import get_cell_area, spatial_reference @@ -131,9 +132,7 @@ def from_imod5_data( cls, imod5_data: dict[str, GridDataDict], target_dis: StructuredDiscretization, - regridder_types: Optional[RegridMethodType] = None, - regrid_cache: RegridderWeightsCache = RegridderWeightsCache(), - ) -> "GridData": + ) -> tuple["GridData", MetaSwapActive]: # Get iMOD5 capillary zone data imod5_cap = imod5_data["cap"] @@ -144,15 +143,7 @@ def from_imod5_data( data["surface_elevation"] = imod5_cap["surface_elevation"] data["soil_physical_unit"] = imod5_cap["soil_physical_unit"] - active, subunit_active = is_msw_active_cell(target_dis, imod5_cap, data["area"]) - - data_active = { - key: ( - griddata.where(subunit_active) - if key in cls._with_subunit - else griddata.where(active) - ) - for key, griddata in data.items() - } - data_active["active"] = active - return cls(**data_active) + msw_active = is_msw_active_cell(target_dis, imod5_cap, data["area"]) + data_active = mask_and_broadcast_grid_data(cls, data, msw_active) + data_active["active"] = msw_active.all + return cls(**data_active), msw_active From 486fa8469f4e3876e5ba5f34034ecf9d56449319 Mon Sep 17 00:00:00 2001 From: JoerivanEngelen Date: Thu, 5 Dec 2024 16:31:38 +0100 Subject: [PATCH 13/50] Use MetaSwapActive --- imod/msw/utilities/imod5_converter.py | 5 +++-- 1 file changed, 3 insertions(+), 2 deletions(-) diff --git a/imod/msw/utilities/imod5_converter.py b/imod/msw/utilities/imod5_converter.py index 040f55c4a..14c36451d 100644 --- a/imod/msw/utilities/imod5_converter.py +++ b/imod/msw/utilities/imod5_converter.py @@ -1,6 +1,7 @@ from imod.logging import LogLevel, logger from imod.mf6 import StructuredDiscretization from imod.msw.utilities.common import concat_imod5 +from imod.msw.utilities.mask import MetaSwapActive from imod.typing import GridDataArray, GridDataDict from imod.typing.grid import ones_like from imod.util.spatial import get_cell_area @@ -63,7 +64,7 @@ def is_msw_active_cell( target_dis: StructuredDiscretization, imod5_cap: GridDataDict, msw_area: GridDataArray, -) -> tuple[GridDataArray, GridDataArray]: +) -> MetaSwapActive: """ Return grid of cells that are active in the coupled computation, based on following criteria: @@ -84,4 +85,4 @@ def is_msw_active_cell( (imod5_cap["boundary"] == 1) & (msw_area > 0) & (mf6_top_active >= 1) ) active = subunit_active.any(dim="subunit") - return active, subunit_active + return MetaSwapActive(active, subunit_active) From f98111ec0a7e96cb53297ea644f09b456a013d6c Mon Sep 17 00:00:00 2001 From: JoerivanEngelen Date: Thu, 5 Dec 2024 16:33:30 +0100 Subject: [PATCH 14/50] Add validate option for write method and add TimeOutputControl, --- imod/msw/model.py | 15 ++++++++++----- 1 file changed, 10 insertions(+), 5 deletions(-) diff --git a/imod/msw/model.py b/imod/msw/model.py index 5c741c2d1..2088cb0f0 100644 --- a/imod/msw/model.py +++ b/imod/msw/model.py @@ -6,6 +6,7 @@ import jinja2 import numpy as np +import xarray as xr from imod.mf6.dis import StructuredDiscretization from imod.mf6.mf6_wel_adapter import Mf6Wel @@ -233,6 +234,7 @@ def write( directory: Union[str, Path], mf6_dis: StructuredDiscretization, mf6_wel: Mf6Wel, + validate: bool = True, ): """ Write packages and simulation settings (para_sim.inp). @@ -244,9 +246,10 @@ def write( """ # Model checks - self._check_required_packages() - self._check_vegetation_indices_in_annual_crop_factors() - self._check_landuse_indices_in_lookup_options() + if validate: + self._check_required_packages() + self._check_vegetation_indices_in_annual_crop_factors() + self._check_landuse_indices_in_lookup_options() # Force to Path directory = Path(directory) @@ -273,6 +276,7 @@ def write( # write package contents for pkgname in self: + print(pkgname) self[pkgname].write(directory, index, svat, mf6_dis, mf6_wel) def regrid_like( @@ -346,7 +350,7 @@ def from_imod5_data( } } model = cls(unsa_svat_path, parasim_settings) - model["grid"] = GridData.from_imod5_data(imod5_cap_no_layer, target_dis) + model["grid"], msw_active = GridData.from_imod5_data(imod5_cap_no_layer, target_dis) active = model["grid"].dataset["active"] # Mask grid cells and broadcast scalars to grid imod5_cap_masked: Imod5DataDict = { @@ -369,6 +373,7 @@ def from_imod5_data( model["coupling"] = CouplerMapping() model["extra_files"] = FileCopier.from_imod5_data(imod5_cap_masked) - model["time_oc"] = TimeOutputControl(times) + times_da = xr.DataArray(times, coords={"time": times}, dims=("time",)) + model["time_oc"] = TimeOutputControl(times_da) return model From d5db2710b23b9c2114ab86cf17952387347d4f2e Mon Sep 17 00:00:00 2001 From: JoerivanEngelen Date: Thu, 5 Dec 2024 16:34:34 +0100 Subject: [PATCH 15/50] Type annotate return type --- imod/msw/utilities/mask.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/imod/msw/utilities/mask.py b/imod/msw/utilities/mask.py index 246127fe5..3a7b4d5c2 100644 --- a/imod/msw/utilities/mask.py +++ b/imod/msw/utilities/mask.py @@ -15,7 +15,7 @@ def mask_and_broadcast_grid_data( package: MetaSwapPackage, grid_data: GridDataDict, msw_active: MetaSwapActive - ): + ) -> GridDataDict: """ Mask and broadcast grid data, carefully mask per subunit if variable needs to contain subunits. From 299d4616d04677ab5f043c7357d7a8281f120611 Mon Sep 17 00:00:00 2001 From: JoerivanEngelen Date: Thu, 5 Dec 2024 18:38:30 +0100 Subject: [PATCH 16/50] mask package data and force soil physical unit to int --- imod/msw/grid_data.py | 12 +++--------- 1 file changed, 3 insertions(+), 9 deletions(-) diff --git a/imod/msw/grid_data.py b/imod/msw/grid_data.py index d2090c7cc..d5dd3eda7 100644 --- a/imod/msw/grid_data.py +++ b/imod/msw/grid_data.py @@ -1,14 +1,8 @@ -from typing import Optional - import numpy as np import xarray as xr from imod.mf6 import StructuredDiscretization from imod.mf6.interfaces.iregridpackage import IRegridPackage -from imod.mf6.regrid.regrid_schemes import ( - RegridMethodType, -) -from imod.mf6.utilities.regrid import RegridderWeightsCache from imod.msw.fixed_format import VariableMetaData from imod.msw.pkgbase import MetaSwapPackage from imod.msw.regrid.regrid_schemes import GridDataRegridMethod @@ -18,7 +12,7 @@ get_rootzone_depth_from_imod5_data, is_msw_active_cell, ) -from imod.msw.utilities.mask import mask_and_broadcast_grid_data, MetaSwapActive +from imod.msw.utilities.mask import mask_package_data, MetaSwapActive from imod.typing import GridDataDict from imod.util.spatial import get_cell_area, spatial_reference @@ -141,9 +135,9 @@ def from_imod5_data( data["landuse"] = get_landuse_from_imod5_data(imod5_cap) data["rootzone_depth"] = get_rootzone_depth_from_imod5_data(imod5_cap) data["surface_elevation"] = imod5_cap["surface_elevation"] - data["soil_physical_unit"] = imod5_cap["soil_physical_unit"] + data["soil_physical_unit"] = imod5_cap["soil_physical_unit"].astype(int) msw_active = is_msw_active_cell(target_dis, imod5_cap, data["area"]) - data_active = mask_and_broadcast_grid_data(cls, data, msw_active) + data_active = mask_package_data(cls, data, msw_active) data_active["active"] = msw_active.all return cls(**data_active), msw_active From 5a347c6f05d518201fe4ebc550ae3ef0408e9cb0 Mon Sep 17 00:00:00 2001 From: JoerivanEngelen Date: Thu, 5 Dec 2024 18:38:40 +0100 Subject: [PATCH 17/50] Deactivate bottom resistance --- imod/msw/infiltration.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/imod/msw/infiltration.py b/imod/msw/infiltration.py index 830365681..ed4f27b23 100644 --- a/imod/msw/infiltration.py +++ b/imod/msw/infiltration.py @@ -120,7 +120,7 @@ def from_imod5_data(cls, imod5_data: dict[str, GridDataDict]) -> "Infiltration": data = deactivate_small_resistances_in_data(data) like = data["downward_resistance"].isel(subunit=0, drop=True) - data["bottom_resistance"] = ones_like(like) + data["bottom_resistance"] = ones_like(like) * -9999.0 data["extra_storage_coefficient"] = ones_like(like) return cls(**data) From 11a670784bae0d18bdf04f4446e3cf71eb838040 Mon Sep 17 00:00:00 2001 From: JoerivanEngelen Date: Thu, 5 Dec 2024 18:39:11 +0100 Subject: [PATCH 18/50] Mask cap data with new function --- imod/msw/model.py | 29 +++++++++++------------------ 1 file changed, 11 insertions(+), 18 deletions(-) diff --git a/imod/msw/model.py b/imod/msw/model.py index 2088cb0f0..ed132cb48 100644 --- a/imod/msw/model.py +++ b/imod/msw/model.py @@ -31,6 +31,7 @@ from imod.msw.sprinkling import Sprinkling from imod.msw.timeutil import to_metaswap_timeformat from imod.msw.utilities.common import find_in_file_list +from imod.msw.utilities.mask import mask_and_broadcast_grid_data from imod.msw.utilities.parse import read_para_sim from imod.msw.vegetation import AnnualCropFactors from imod.typing import Imod5DataDict, SelSettingsType @@ -351,27 +352,19 @@ def from_imod5_data( } model = cls(unsa_svat_path, parasim_settings) model["grid"], msw_active = GridData.from_imod5_data(imod5_cap_no_layer, target_dis) - active = model["grid"].dataset["active"] - # Mask grid cells and broadcast scalars to grid - imod5_cap_masked: Imod5DataDict = { - "cap": { - key: da.where(active) for key, da in imod5_cap_no_layer["cap"].items() - }, - "extra": {"paths": extra_paths}, - } - model["infiltration"] = Infiltration.from_imod5_data(imod5_cap_masked) - model["ponding"] = Ponding.from_imod5_data(imod5_cap_masked) - model["sprinkling"] = Sprinkling.from_imod5_data(imod5_cap_masked) - model["meteo_grid"] = MeteoGridCopy.from_imod5_data(imod5_cap_masked) - model["prec_mapping"] = PrecipitationMapping.from_imod5_data(imod5_cap_masked) - model["evt_mapping"] = EvapotranspirationMapping.from_imod5_data( - imod5_cap_masked - ) - model["scaling_factor"] = ScalingFactors.from_imod5_data(imod5_cap_masked) + cap_data_masked = mask_and_broadcast_grid_data(imod5_cap_no_layer["cap"], msw_active) + imod5_masked = {"cap": cap_data_masked, "extra": {"paths": extra_paths}} + model["infiltration"] = Infiltration.from_imod5_data(imod5_masked) + model["ponding"] = Ponding.from_imod5_data(imod5_masked) + model["sprinkling"] = Sprinkling.from_imod5_data(imod5_masked) + model["meteo_grid"] = MeteoGridCopy.from_imod5_data(imod5_masked) + model["prec_mapping"] = PrecipitationMapping.from_imod5_data(imod5_masked) + model["evt_mapping"] = EvapotranspirationMapping.from_imod5_data(imod5_masked) + # model["scaling_factor"] = ScalingFactors.from_imod5_data(imod5_masked) area = model["grid"]["area"].isel(subunit=0, drop=True) model["idf_mapping"] = IdfMapping(area, -9999.0) model["coupling"] = CouplerMapping() - model["extra_files"] = FileCopier.from_imod5_data(imod5_cap_masked) + model["extra_files"] = FileCopier.from_imod5_data(imod5_masked) times_da = xr.DataArray(times, coords={"time": times}, dims=("time",)) model["time_oc"] = TimeOutputControl(times_da) From 5f6d449bf6085972fb65b452eeb353b7c31589c0 Mon Sep 17 00:00:00 2001 From: JoerivanEngelen Date: Thu, 5 Dec 2024 18:39:23 +0100 Subject: [PATCH 19/50] force landuse to int dtype --- imod/msw/utilities/imod5_converter.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/imod/msw/utilities/imod5_converter.py b/imod/msw/utilities/imod5_converter.py index 14c36451d..fc5066a78 100644 --- a/imod/msw/utilities/imod5_converter.py +++ b/imod/msw/utilities/imod5_converter.py @@ -44,7 +44,7 @@ def get_landuse_from_imod5_data( rural_landuse = imod5_cap["landuse"] # Urban landuse = 18 urban_landuse = ones_like(rural_landuse) * 18 - return concat_imod5(rural_landuse, urban_landuse) + return concat_imod5(rural_landuse, urban_landuse).astype(int) def get_rootzone_depth_from_imod5_data( From e5df224390b75b2f496a59bd48ea6e5467c7aeca Mon Sep 17 00:00:00 2001 From: JoerivanEngelen Date: Thu, 5 Dec 2024 18:40:12 +0100 Subject: [PATCH 20/50] Rename to mask_package_data and implement simplified mask_and_broadcast_grid_data --- imod/msw/utilities/mask.py | 14 ++++++++++++-- 1 file changed, 12 insertions(+), 2 deletions(-) diff --git a/imod/msw/utilities/mask.py b/imod/msw/utilities/mask.py index 3a7b4d5c2..f229f037e 100644 --- a/imod/msw/utilities/mask.py +++ b/imod/msw/utilities/mask.py @@ -10,8 +10,18 @@ class MetaSwapActive: all: GridDataArray per_subunit: GridDataArray - def mask_and_broadcast_grid_data( + grid_data: GridDataDict, + msw_active: MetaSwapActive + ) -> GridDataDict: + """ + Mask and broadcast grid data, + """ + return { + key: _mask_spatial_var(grid, msw_active.all) for key, grid in grid_data.items() + } + +def mask_package_data( package: MetaSwapPackage, grid_data: GridDataDict, msw_active: MetaSwapActive @@ -30,7 +40,7 @@ def mask_and_broadcast_grid_data( for key, grid in grid_data.items() } -def _mask_spatial_var(self, da: GridDataArray, active: GridDataArray) -> GridDataArray: +def _mask_spatial_var(da: GridDataArray, active: GridDataArray) -> GridDataArray: if issubclass(da.dtype.type, numbers.Integral): return da.where(active, other=0) elif issubclass(da.dtype.type, numbers.Real): From 22b12a07a30e74a841197808b0ccea6c1689e687 Mon Sep 17 00:00:00 2001 From: JoerivanEngelen Date: Thu, 5 Dec 2024 18:40:27 +0100 Subject: [PATCH 21/50] format --- imod/msw/grid_data.py | 2 +- imod/msw/model.py | 9 ++++++--- imod/msw/utilities/mask.py | 18 +++++++++--------- 3 files changed, 16 insertions(+), 13 deletions(-) diff --git a/imod/msw/grid_data.py b/imod/msw/grid_data.py index d5dd3eda7..4c202697c 100644 --- a/imod/msw/grid_data.py +++ b/imod/msw/grid_data.py @@ -12,7 +12,7 @@ get_rootzone_depth_from_imod5_data, is_msw_active_cell, ) -from imod.msw.utilities.mask import mask_package_data, MetaSwapActive +from imod.msw.utilities.mask import MetaSwapActive, mask_package_data from imod.typing import GridDataDict from imod.util.spatial import get_cell_area, spatial_reference diff --git a/imod/msw/model.py b/imod/msw/model.py index ed132cb48..adae7c972 100644 --- a/imod/msw/model.py +++ b/imod/msw/model.py @@ -27,7 +27,6 @@ from imod.msw.meteo_mapping import EvapotranspirationMapping, PrecipitationMapping from imod.msw.output_control import TimeOutputControl from imod.msw.ponding import Ponding -from imod.msw.scaling_factors import ScalingFactors from imod.msw.sprinkling import Sprinkling from imod.msw.timeutil import to_metaswap_timeformat from imod.msw.utilities.common import find_in_file_list @@ -351,8 +350,12 @@ def from_imod5_data( } } model = cls(unsa_svat_path, parasim_settings) - model["grid"], msw_active = GridData.from_imod5_data(imod5_cap_no_layer, target_dis) - cap_data_masked = mask_and_broadcast_grid_data(imod5_cap_no_layer["cap"], msw_active) + model["grid"], msw_active = GridData.from_imod5_data( + imod5_cap_no_layer, target_dis + ) + cap_data_masked = mask_and_broadcast_grid_data( + imod5_cap_no_layer["cap"], msw_active + ) imod5_masked = {"cap": cap_data_masked, "extra": {"paths": extra_paths}} model["infiltration"] = Infiltration.from_imod5_data(imod5_masked) model["ponding"] = Ponding.from_imod5_data(imod5_masked) diff --git a/imod/msw/utilities/mask.py b/imod/msw/utilities/mask.py index f229f037e..62e2ef87e 100644 --- a/imod/msw/utilities/mask.py +++ b/imod/msw/utilities/mask.py @@ -10,22 +10,21 @@ class MetaSwapActive: all: GridDataArray per_subunit: GridDataArray + def mask_and_broadcast_grid_data( - grid_data: GridDataDict, - msw_active: MetaSwapActive - ) -> GridDataDict: + grid_data: GridDataDict, msw_active: MetaSwapActive +) -> GridDataDict: """ - Mask and broadcast grid data, + Mask and broadcast grid data, """ return { key: _mask_spatial_var(grid, msw_active.all) for key, grid in grid_data.items() } + def mask_package_data( - package: MetaSwapPackage, - grid_data: GridDataDict, - msw_active: MetaSwapActive - ) -> GridDataDict: + package: MetaSwapPackage, grid_data: GridDataDict, msw_active: MetaSwapActive +) -> GridDataDict: """ Mask and broadcast grid data, carefully mask per subunit if variable needs to contain subunits. @@ -40,6 +39,7 @@ def mask_package_data( for key, grid in grid_data.items() } + def _mask_spatial_var(da: GridDataArray, active: GridDataArray) -> GridDataArray: if issubclass(da.dtype.type, numbers.Integral): return da.where(active, other=0) @@ -48,4 +48,4 @@ def _mask_spatial_var(da: GridDataArray, active: GridDataArray) -> GridDataArray else: raise TypeError( f"Expected dtype float or integer. Received instead: {da.dtype}" - ) \ No newline at end of file + ) From e7920959a79d37a316c65fd2bbffe81e6c960c35 Mon Sep 17 00:00:00 2001 From: JoerivanEngelen Date: Fri, 6 Dec 2024 12:07:01 +0100 Subject: [PATCH 22/50] Move setup_meteo_grids to separate fixture in module --- imod/tests/fixtures/msw_meteo_fixture.py | 40 ++++++++++++ imod/tests/test_msw/test_meteo_grid.py | 80 +++++------------------- 2 files changed, 54 insertions(+), 66 deletions(-) create mode 100644 imod/tests/fixtures/msw_meteo_fixture.py diff --git a/imod/tests/fixtures/msw_meteo_fixture.py b/imod/tests/fixtures/msw_meteo_fixture.py new file mode 100644 index 000000000..0331ca889 --- /dev/null +++ b/imod/tests/fixtures/msw_meteo_fixture.py @@ -0,0 +1,40 @@ +import pandas as pd +import xarray as xr +import numpy as np +from numpy import nan +import pytest + +@pytest.fixture(scope="function") +def meteo_grids(): + x = [1.0, 2.0, 3.0] + y = [1.0, 2.0, 3.0] + time = pd.date_range(start="2000-01-01", end="2000-01-02", freq="D") + + dx = 1.0 + dy = -1.0 + # fmt: off + precipitation = xr.DataArray( + np.array( + [ + [[1.0, 1.0, 1.0], + [nan, nan, nan], + [1.0, 1.0, 1.0]], + + [[2.0, 2.0, 1.0], + [nan, nan, nan], + [1.0, 2.0, 1.0]], + ] + ), + dims=("time", "y", "x"), + coords={"time": time, "y": y, "x": x, "dx": dx, "dy": dy} + ) + + evapotranspiration = xr.DataArray( + np.array( + [1.0, 3.0] + ), + dims=("time",), + coords={"time": time} + ) + # fmt: on + return precipitation, evapotranspiration \ No newline at end of file diff --git a/imod/tests/test_msw/test_meteo_grid.py b/imod/tests/test_msw/test_meteo_grid.py index b5108842e..7d614f67a 100644 --- a/imod/tests/test_msw/test_meteo_grid.py +++ b/imod/tests/test_msw/test_meteo_grid.py @@ -15,46 +15,15 @@ from imod.msw import MeteoGrid, MeteoGridCopy -def setup_meteo_grid(): - x = [1.0, 2.0, 3.0] - y = [1.0, 2.0, 3.0] - time = pd.date_range(start="2000-01-01", end="2000-01-02", freq="D") - - dx = 1.0 - dy = -1.0 - # fmt: off - precipitation = xr.DataArray( - np.array( - [ - [[1.0, 1.0, 1.0], - [nan, nan, nan], - [1.0, 1.0, 1.0]], - - [[2.0, 2.0, 1.0], - [nan, nan, nan], - [1.0, 2.0, 1.0]], - ] - ), - dims=("time", "y", "x"), - coords={"time": time, "y": y, "x": x, "dx": dx, "dy": dy} - ) - - evapotranspiration = xr.DataArray( - np.array( - [1.0, 3.0] - ), - dims=("time",), - coords={"time": time} - ) - # fmt: on - - meteo_grid = MeteoGrid(precipitation, evapotranspiration) +def test_meteo_grid_init(meteo_grids): + meteo_grids = MeteoGrid(*meteo_grids) - return meteo_grid + assert meteo_grids.dataset["precipitation"].dims == ("time", "y", "x") + assert meteo_grids.dataset["evapotranspiration"].dims == ("time",) -def test_meteo_grid(): - meteo_grid = setup_meteo_grid() +def test_meteo_grid_write(meteo_grids): + meteo_grid = MeteoGrid(*meteo_grids) with tempfile.TemporaryDirectory() as output_dir: output_dir = Path(output_dir) @@ -83,29 +52,8 @@ def test_meteo_grid(): assert gridnames == expected_filenames -def test_meteo_no_time_grid(): - x = [1.0, 2.0, 3.0] - y = [1.0, 2.0, 3.0] - time = pd.date_range(start="2000-01-01", end="2000-01-02", freq="D") - - dx = 1.0 - dy = -1.0 - # fmt: off - precipitation = xr.DataArray( - np.array( - [ - [[1.0, 1.0, 1.0], - [nan, nan, nan], - [1.0, 1.0, 1.0]], - - [[2.0, 2.0, 1.0], - [nan, nan, nan], - [1.0, 2.0, 1.0]], - ] - ), - dims=("time", "y", "x"), - coords={"time": time, "y": y, "x": x, "dx": dx, "dy": dy} - ) +def test_meteo_no_time_grid(meteo_grids): + precipitation, _ = meteo_grids evapotranspiration = 3.0 # fmt: on @@ -114,8 +62,8 @@ def test_meteo_no_time_grid(): MeteoGrid(precipitation, evapotranspiration) -def test_regrid_meteo(simple_2d_grid_with_subunits): - meteo = setup_meteo_grid() +def test_regrid_meteo(meteo_grids, simple_2d_grid_with_subunits): + meteo = MeteoGrid(*meteo_grids) new_grid = simple_2d_grid_with_subunits regrid_context = RegridderWeightsCache() @@ -126,9 +74,9 @@ def test_regrid_meteo(simple_2d_grid_with_subunits): assert np.all(regridded_ponding.dataset["y"].values == new_grid["y"].values) -def test_meteogridcopy_write(): +def test_meteogridcopy_write(meteo_grids): # Arrange - meteo_grid = setup_meteo_grid() + meteo_grid = MeteoGrid(*meteo_grids) with tempfile.TemporaryDirectory() as output_dir: grid_dir = Path(output_dir) / "grid" @@ -144,8 +92,8 @@ def test_meteogridcopy_write(): assert filecmp.cmp(grid_dir / "mete_grid.inp", copy_dir / "mete_grid.inp") -def test_meteogridcopy_from_imod5(): - meteo_grid = setup_meteo_grid() +def test_meteogridcopy_from_imod5(meteo_grids): + meteo_grid = MeteoGrid(*meteo_grids) with tempfile.TemporaryDirectory() as output_dir: grid_dir = Path(output_dir) / "grid" From bb984b64937fa1cf7144231f91f90d76ba33b625 Mon Sep 17 00:00:00 2001 From: JoerivanEngelen Date: Fri, 6 Dec 2024 12:07:32 +0100 Subject: [PATCH 23/50] Add iMOD5 cap fixture --- imod/tests/fixtures/msw_imod5_cap_fixture.py | 227 +++++++++++++++++++ 1 file changed, 227 insertions(+) create mode 100644 imod/tests/fixtures/msw_imod5_cap_fixture.py diff --git a/imod/tests/fixtures/msw_imod5_cap_fixture.py b/imod/tests/fixtures/msw_imod5_cap_fixture.py new file mode 100644 index 000000000..bf201b0b3 --- /dev/null +++ b/imod/tests/fixtures/msw_imod5_cap_fixture.py @@ -0,0 +1,227 @@ +import numpy as np +import pytest +import xarray as xr +from numpy import nan + +from imod.typing import GridDataDict + + +@pytest.fixture(scope="function") +def imod5_cap_data() -> GridDataDict: + x = [1.0, 2.0, 3.0] + y = [3.0, 2.0, 1.0] + layer = [1] + dx = 1.0 + dy = -1.0 + + da_kwargs = {} + da_kwargs["dims"] = ("layer", "y", "x") + da_kwargs["coords"] = {"layer": layer, "y": y, "x": x, "dx": dx, "dy": dy} + + imod5_data = {} + d = {} + + # fmt: off + d["boundary"] = xr.DataArray( + np.array( + [ + [1, 1, 1], + [0, 0, 0], + [1, 1, 0], + ], + dtype=int), + **da_kwargs + ) + d["landuse"] = xr.DataArray( + np.array( + [ + [1, 2, 3], + [0, 0, 0], + [3, 2, 1], + ], + dtype=int), + **da_kwargs + ) + d["rootzone_thickness"] = xr.DataArray( + np.array( + [ + [0.1, 0.1, 0.1], + [nan, nan, nan], + [0.2, 0.2, 0.2], + ] + ), + **da_kwargs + ) + d["surface_elevation"] = xr.DataArray( + np.array( + [ + [1.1, 1.2, 1.3], + [nan, nan, nan], + [1.4, 1.5, 1.6], + ] + ), + **da_kwargs + ) + d["artificial_recharge"] = xr.DataArray( + np.array( + [ + [0.2, 0.2, 0.2], + [nan, nan, nan], + [0.3, 0.3, 0.3], + ] + ), + **da_kwargs + ) + d["artificial_recharge_layer"] = xr.DataArray( + np.array( + [ + [1, 2, 3], + [0, 0, 0], + [3, 2, 1], + ], + dtype=int), + **da_kwargs + ) + d["artificial_recharge_capacity"] = xr.DataArray( + np.array( + [ + [0.4, 0.4, 0.4], + [nan, nan, nan], + [0.6, 0.6, 0.6], + ] + ), + **da_kwargs + ) + d["wet_area"] = xr.DataArray( + np.array( + [ + [0.1, 0.1, 0.1], + [nan, nan, nan], + [0.2, 0.2, 0.2], + ] + ), + **da_kwargs + ) + d["urban_area"] = xr.DataArray( + np.array( + [ + [0.2, 0.2, 0.2], + [nan, nan, nan], + [0.3, 0.3, 0.3], + ] + ), + **da_kwargs + ) + d["urban_ponding_depth"] = xr.DataArray( + np.array( + [ + [2.2, 2.2, 2.2], + [nan, nan, nan], + [2.3, 2.3, 2.3], + ] + ), + **da_kwargs + ) + d["rural_ponding_depth"] = xr.DataArray( + np.array( + [ + [1.2, 1.2, 1.2], + [nan, nan, nan], + [1.3, 1.3, 1.3], + ] + ), + **da_kwargs + ) + d["urban_runoff_resistance"] = xr.DataArray( + np.array( + [ + [3.2, 3.2, 3.2], + [nan, nan, nan], + [3.3, 3.3, 3.3], + ] + ), + **da_kwargs + ) + d["rural_runoff_resistance"] = xr.DataArray( + np.array( + [ + [3.6, 3.6, 3.6], + [nan, nan, nan], + [3.7, 3.7, 3.7], + ] + ), + **da_kwargs + ) + d["urban_runon_resistance"] = xr.DataArray( + np.array( + [ + [5.2, 5.2, 5.2], + [nan, nan, nan], + [5.3, 5.3, 5.3], + ] + ), + **da_kwargs + ) + d["rural_runon_resistance"] = xr.DataArray( + np.array( + [ + [5.6, 5.6, 5.6], + [nan, nan, nan], + [5.7, 5.7, 5.7], + ] + ), + **da_kwargs + ) + d["urban_infiltration_capacity"] = xr.DataArray( + np.array( + [ + [10.2, 10.2, 10.2], + [nan, nan, nan], + [10.3, 10.3, 10.3], + ] + ), + **da_kwargs + ) + d["rural_infiltration_capacity"] = xr.DataArray( + np.array( + [ + [20.2, 20.2, 20.2], + [nan, nan, nan], + [20.3, 20.3, 20.3], + ] + ), + **da_kwargs + ) + d["perched_water_table_level"]= xr.DataArray( + np.array( + [ + [2.0, 2.0, 2.0], + [nan, nan, nan], + [2.0, 2.0, 2.0], + ] + ), + **da_kwargs + ) + d["soil_moisture_fraction"]= xr.DataArray( + np.array( + [ + [1.5, 1.5, 1.5], + [nan, nan, nan], + [1.5, 1.5, 1.5], + ] + ), + **da_kwargs + ) + d["conductivitiy_factor"]= xr.DataArray( + np.array( + [ + [2.5, 2.5, 2.5], + [nan, nan, nan], + [2.5, 2.5, 2.5], + ] + ), + **da_kwargs + ) + # fmt: on + imod5_data["cap"] = d + return imod5_data From 422faffb846332c60625abc8d911037f73b1de0a Mon Sep 17 00:00:00 2001 From: JoerivanEngelen Date: Fri, 6 Dec 2024 12:07:45 +0100 Subject: [PATCH 24/50] Import added fixtures --- imod/tests/conftest.py | 2 ++ 1 file changed, 2 insertions(+) diff --git a/imod/tests/conftest.py b/imod/tests/conftest.py index 0677d4b4f..a239de23d 100644 --- a/imod/tests/conftest.py +++ b/imod/tests/conftest.py @@ -92,4 +92,6 @@ well_test_data_transient, ) from .fixtures.msw_fixture import fixed_format_parser, simple_2d_grid_with_subunits +from .fixtures.msw_meteo_fixture import meteo_grids from .fixtures.msw_model_fixture import coupled_mf6_model, coupled_mf6wel, msw_model +from .fixtures.msw_imod5_cap_fixture import imod5_cap_data \ No newline at end of file From 6adc1eafc113bc4bcf3ba534200ddaadad83c374 Mon Sep 17 00:00:00 2001 From: JoerivanEngelen Date: Fri, 6 Dec 2024 13:12:58 +0100 Subject: [PATCH 25/50] Use tmp_path --- imod/tests/test_msw/test_meteo_grid.py | 64 +++++++++++++------------- 1 file changed, 31 insertions(+), 33 deletions(-) diff --git a/imod/tests/test_msw/test_meteo_grid.py b/imod/tests/test_msw/test_meteo_grid.py index 7d614f67a..66ff59ec1 100644 --- a/imod/tests/test_msw/test_meteo_grid.py +++ b/imod/tests/test_msw/test_meteo_grid.py @@ -74,44 +74,42 @@ def test_regrid_meteo(meteo_grids, simple_2d_grid_with_subunits): assert np.all(regridded_ponding.dataset["y"].values == new_grid["y"].values) -def test_meteogridcopy_write(meteo_grids): +def test_meteogridcopy_write(meteo_grids, tmp_path): # Arrange meteo_grid = MeteoGrid(*meteo_grids) - with tempfile.TemporaryDirectory() as output_dir: - grid_dir = Path(output_dir) / "grid" - grid_dir.mkdir(exist_ok=True, parents=True) - meteo_grid.write(grid_dir) + grid_dir = tmp_path / "grid" + grid_dir.mkdir(exist_ok=True, parents=True) + meteo_grid.write(grid_dir) - meteo_grid_copy = MeteoGridCopy(grid_dir / "mete_grid.inp") - copy_dir = Path(output_dir) / "copied" - copy_dir.mkdir(exist_ok=True, parents=True) - # Act - meteo_grid_copy.write(copy_dir) - # Assert - assert filecmp.cmp(grid_dir / "mete_grid.inp", copy_dir / "mete_grid.inp") + meteo_grid_copy = MeteoGridCopy(grid_dir / "mete_grid.inp") + copy_dir = tmp_path / "copied" + copy_dir.mkdir(exist_ok=True, parents=True) + # Act + meteo_grid_copy.write(copy_dir) + # Assert + assert filecmp.cmp(grid_dir / "mete_grid.inp", copy_dir / "mete_grid.inp") -def test_meteogridcopy_from_imod5(meteo_grids): +def test_meteogridcopy_from_imod5(meteo_grids, tmp_path): meteo_grid = MeteoGrid(*meteo_grids) - with tempfile.TemporaryDirectory() as output_dir: - grid_dir = Path(output_dir) / "grid" - grid_dir.mkdir(exist_ok=True, parents=True) - meteo_grid.write(grid_dir) - - imod5_data = {} - imod5_data["extra"] = {} - imod5_data["extra"]["paths"] = [ - ["foo"], - [(grid_dir / "mete_grid.inp").resolve()], - ["bar"], - ] - - meteo_grid_copy = MeteoGridCopy.from_imod5_data(imod5_data) - copy_dir = Path(output_dir) / "copied" - copy_dir.mkdir(exist_ok=True, parents=True) - # Act - meteo_grid_copy.write(copy_dir) - # Assert - assert filecmp.cmp(grid_dir / "mete_grid.inp", copy_dir / "mete_grid.inp") + grid_dir = tmp_path / "grid" + grid_dir.mkdir(exist_ok=True, parents=True) + meteo_grid.write(grid_dir) + + imod5_data = {} + imod5_data["extra"] = {} + imod5_data["extra"]["paths"] = [ + ["foo"], + [(grid_dir / "mete_grid.inp").resolve()], + ["bar"], + ] + + meteo_grid_copy = MeteoGridCopy.from_imod5_data(imod5_data) + copy_dir = tmp_path / "copied" + copy_dir.mkdir(exist_ok=True, parents=True) + # Act + meteo_grid_copy.write(copy_dir) + # Assert + assert filecmp.cmp(grid_dir / "mete_grid.inp", copy_dir / "mete_grid.inp") From 4ec678162ef873294830ef5e8527e24eb5b71fe0 Mon Sep 17 00:00:00 2001 From: JoerivanEngelen Date: Fri, 6 Dec 2024 17:36:48 +0100 Subject: [PATCH 26/50] Format --- imod/tests/conftest.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/imod/tests/conftest.py b/imod/tests/conftest.py index a239de23d..b0f400192 100644 --- a/imod/tests/conftest.py +++ b/imod/tests/conftest.py @@ -92,6 +92,6 @@ well_test_data_transient, ) from .fixtures.msw_fixture import fixed_format_parser, simple_2d_grid_with_subunits +from .fixtures.msw_imod5_cap_fixture import imod5_cap_data from .fixtures.msw_meteo_fixture import meteo_grids from .fixtures.msw_model_fixture import coupled_mf6_model, coupled_mf6wel, msw_model -from .fixtures.msw_imod5_cap_fixture import imod5_cap_data \ No newline at end of file From 5c41cc79515f1faac573d8df48d1b1b5cda6c31a Mon Sep 17 00:00:00 2001 From: JoerivanEngelen Date: Fri, 6 Dec 2024 17:37:07 +0100 Subject: [PATCH 27/50] Fixes to cap data fixture --- imod/tests/fixtures/msw_imod5_cap_fixture.py | 92 +++++++++++--------- 1 file changed, 51 insertions(+), 41 deletions(-) diff --git a/imod/tests/fixtures/msw_imod5_cap_fixture.py b/imod/tests/fixtures/msw_imod5_cap_fixture.py index bf201b0b3..c57925015 100644 --- a/imod/tests/fixtures/msw_imod5_cap_fixture.py +++ b/imod/tests/fixtures/msw_imod5_cap_fixture.py @@ -23,202 +23,212 @@ def imod5_cap_data() -> GridDataDict: # fmt: off d["boundary"] = xr.DataArray( - np.array( + np.array([ [ [1, 1, 1], [0, 0, 0], [1, 1, 0], - ], + ]], dtype=int), **da_kwargs ) d["landuse"] = xr.DataArray( - np.array( + np.array([ [ [1, 2, 3], [0, 0, 0], [3, 2, 1], - ], + ]], dtype=int), **da_kwargs ) d["rootzone_thickness"] = xr.DataArray( - np.array( + np.array([ [ [0.1, 0.1, 0.1], [nan, nan, nan], [0.2, 0.2, 0.2], - ] + ]] ), **da_kwargs ) + d["soil_physical_unit"] = xr.DataArray( + np.array([ + [ + [1, 1, 1], + [0, 0, 0], + [2, 1, 1], + ]], + dtype=int), + **da_kwargs + ) d["surface_elevation"] = xr.DataArray( - np.array( + np.array([ [ [1.1, 1.2, 1.3], [nan, nan, nan], [1.4, 1.5, 1.6], - ] + ]] ), **da_kwargs ) d["artificial_recharge"] = xr.DataArray( - np.array( + np.array([ [ [0.2, 0.2, 0.2], [nan, nan, nan], [0.3, 0.3, 0.3], - ] + ]] ), **da_kwargs ) d["artificial_recharge_layer"] = xr.DataArray( - np.array( + np.array([ [ [1, 2, 3], [0, 0, 0], [3, 2, 1], - ], + ]], dtype=int), **da_kwargs ) d["artificial_recharge_capacity"] = xr.DataArray( - np.array( + np.array([ [ [0.4, 0.4, 0.4], [nan, nan, nan], [0.6, 0.6, 0.6], - ] + ]] ), **da_kwargs ) - d["wet_area"] = xr.DataArray( - np.array( + d["wetted_area"] = xr.DataArray( + np.array([ [ [0.1, 0.1, 0.1], [nan, nan, nan], [0.2, 0.2, 0.2], - ] + ]] ), **da_kwargs ) d["urban_area"] = xr.DataArray( - np.array( + np.array([ [ [0.2, 0.2, 0.2], [nan, nan, nan], [0.3, 0.3, 0.3], - ] + ]] ), **da_kwargs ) d["urban_ponding_depth"] = xr.DataArray( - np.array( + np.array([ [ [2.2, 2.2, 2.2], [nan, nan, nan], [2.3, 2.3, 2.3], - ] + ]] ), **da_kwargs ) d["rural_ponding_depth"] = xr.DataArray( - np.array( + np.array([ [ [1.2, 1.2, 1.2], [nan, nan, nan], [1.3, 1.3, 1.3], - ] + ]] ), **da_kwargs ) d["urban_runoff_resistance"] = xr.DataArray( - np.array( + np.array([ [ [3.2, 3.2, 3.2], [nan, nan, nan], [3.3, 3.3, 3.3], - ] + ]] ), **da_kwargs ) d["rural_runoff_resistance"] = xr.DataArray( - np.array( + np.array([ [ [3.6, 3.6, 3.6], [nan, nan, nan], [3.7, 3.7, 3.7], - ] + ]] ), **da_kwargs ) d["urban_runon_resistance"] = xr.DataArray( - np.array( + np.array([ [ [5.2, 5.2, 5.2], [nan, nan, nan], [5.3, 5.3, 5.3], - ] + ]] ), **da_kwargs ) d["rural_runon_resistance"] = xr.DataArray( - np.array( + np.array([ [ [5.6, 5.6, 5.6], [nan, nan, nan], [5.7, 5.7, 5.7], - ] + ]] ), **da_kwargs ) d["urban_infiltration_capacity"] = xr.DataArray( - np.array( + np.array([ [ [10.2, 10.2, 10.2], [nan, nan, nan], [10.3, 10.3, 10.3], - ] + ]] ), **da_kwargs ) d["rural_infiltration_capacity"] = xr.DataArray( - np.array( + np.array([ [ [20.2, 20.2, 20.2], [nan, nan, nan], [20.3, 20.3, 20.3], - ] + ]] ), **da_kwargs ) d["perched_water_table_level"]= xr.DataArray( - np.array( + np.array([ [ [2.0, 2.0, 2.0], [nan, nan, nan], [2.0, 2.0, 2.0], - ] + ]] ), **da_kwargs ) d["soil_moisture_fraction"]= xr.DataArray( - np.array( + np.array([ [ [1.5, 1.5, 1.5], [nan, nan, nan], [1.5, 1.5, 1.5], - ] + ]] ), **da_kwargs ) d["conductivitiy_factor"]= xr.DataArray( - np.array( + np.array([ [ [2.5, 2.5, 2.5], [nan, nan, nan], [2.5, 2.5, 2.5], - ] + ]] ), **da_kwargs ) From aa6d67fc02c3b1593231a1943edd266d67b46a15 Mon Sep 17 00:00:00 2001 From: JoerivanEngelen Date: Fri, 6 Dec 2024 17:37:22 +0100 Subject: [PATCH 28/50] Type annotate --- imod/tests/fixtures/msw_meteo_fixture.py | 11 +++++++---- 1 file changed, 7 insertions(+), 4 deletions(-) diff --git a/imod/tests/fixtures/msw_meteo_fixture.py b/imod/tests/fixtures/msw_meteo_fixture.py index 0331ca889..ebfdee449 100644 --- a/imod/tests/fixtures/msw_meteo_fixture.py +++ b/imod/tests/fixtures/msw_meteo_fixture.py @@ -1,11 +1,14 @@ +import numpy as np import pandas as pd +import pytest import xarray as xr -import numpy as np from numpy import nan -import pytest + +from imod.typing import GridDataArray + @pytest.fixture(scope="function") -def meteo_grids(): +def meteo_grids() -> tuple[GridDataArray, GridDataArray]: x = [1.0, 2.0, 3.0] y = [1.0, 2.0, 3.0] time = pd.date_range(start="2000-01-01", end="2000-01-02", freq="D") @@ -37,4 +40,4 @@ def meteo_grids(): coords={"time": time} ) # fmt: on - return precipitation, evapotranspiration \ No newline at end of file + return precipitation, evapotranspiration From a09adfb153aa024bc45141040c0e80a25a387270 Mon Sep 17 00:00:00 2001 From: JoerivanEngelen Date: Fri, 6 Dec 2024 17:37:53 +0100 Subject: [PATCH 29/50] Move meteogrids setup outside tests --- imod/tests/test_msw/test_meteo_grid.py | 20 ++++++++------------ 1 file changed, 8 insertions(+), 12 deletions(-) diff --git a/imod/tests/test_msw/test_meteo_grid.py b/imod/tests/test_msw/test_meteo_grid.py index 66ff59ec1..66d7fe95e 100644 --- a/imod/tests/test_msw/test_meteo_grid.py +++ b/imod/tests/test_msw/test_meteo_grid.py @@ -7,8 +7,6 @@ import numpy as np import pandas as pd import pytest -import xarray as xr -from numpy import nan from numpy.testing import assert_equal from imod.mf6.utilities.regrid import RegridderWeightsCache @@ -74,14 +72,17 @@ def test_regrid_meteo(meteo_grids, simple_2d_grid_with_subunits): assert np.all(regridded_ponding.dataset["y"].values == new_grid["y"].values) -def test_meteogridcopy_write(meteo_grids, tmp_path): - # Arrange +def setup_written_meteo_grids(meteo_grids, tmp_path) -> Path: meteo_grid = MeteoGrid(*meteo_grids) - grid_dir = tmp_path / "grid" grid_dir.mkdir(exist_ok=True, parents=True) meteo_grid.write(grid_dir) + return grid_dir + +def test_meteogridcopy_write(meteo_grids, tmp_path): + # Arrange + grid_dir = setup_written_meteo_grids(meteo_grids, tmp_path) meteo_grid_copy = MeteoGridCopy(grid_dir / "mete_grid.inp") copy_dir = tmp_path / "copied" copy_dir.mkdir(exist_ok=True, parents=True) @@ -92,12 +93,8 @@ def test_meteogridcopy_write(meteo_grids, tmp_path): def test_meteogridcopy_from_imod5(meteo_grids, tmp_path): - meteo_grid = MeteoGrid(*meteo_grids) - - grid_dir = tmp_path / "grid" - grid_dir.mkdir(exist_ok=True, parents=True) - meteo_grid.write(grid_dir) - + # Arrange + grid_dir = setup_written_meteo_grids(meteo_grids, tmp_path) imod5_data = {} imod5_data["extra"] = {} imod5_data["extra"]["paths"] = [ @@ -105,7 +102,6 @@ def test_meteogridcopy_from_imod5(meteo_grids, tmp_path): [(grid_dir / "mete_grid.inp").resolve()], ["bar"], ] - meteo_grid_copy = MeteoGridCopy.from_imod5_data(imod5_data) copy_dir = tmp_path / "copied" copy_dir.mkdir(exist_ok=True, parents=True) From d296a1dfd7cbf7a8b69506c25c194a132edcea1c Mon Sep 17 00:00:00 2001 From: JoerivanEngelen Date: Fri, 6 Dec 2024 17:38:21 +0100 Subject: [PATCH 30/50] Add tests for MetaSwapModel.from_imod5_data --- imod/tests/test_msw/test_model.py | 110 ++++++++++++++++++++++++++++++ 1 file changed, 110 insertions(+) diff --git a/imod/tests/test_msw/test_model.py b/imod/tests/test_msw/test_model.py index 0e97a59c5..e804a1560 100644 --- a/imod/tests/test_msw/test_model.py +++ b/imod/tests/test_msw/test_model.py @@ -1,3 +1,4 @@ +from copy import copy from pathlib import Path import pytest @@ -6,6 +7,8 @@ from imod import mf6, msw from imod.mf6.utilities.regrid import RegridderWeightsCache +from imod.msw.model import DEFAULT_SETTINGS +from imod.typing import GridDataArray, Imod5DataDict def test_msw_model_write(msw_model, coupled_mf6_model, coupled_mf6wel, tmp_path): @@ -138,3 +141,110 @@ def test_coupled_model_regrid(msw_model, coupled_mf6_model, tmp_path): ) regridded_msw_model.write(tmp_path, mf6_discretization, regridded_mf6_wel) + + +def setup_written_meteo_grids( + meteo_grids: tuple[GridDataArray], tmp_path: Path +) -> Path: + precipitation, _ = meteo_grids + meteo_grid = msw.MeteoGrid(precipitation, precipitation) + grid_dir = tmp_path / "grid" + grid_dir.mkdir(exist_ok=True, parents=True) + meteo_grid.write(grid_dir) + return grid_dir + + +def setup_parasim_inp(directory: Path): + settings = copy(DEFAULT_SETTINGS) + settings["iybg"] = 2000 + settings["tdbg"] = 200.45 + settings["unsa_svat_path"] = str(directory) + + filename = directory / "para_sim.inp" + with open(filename, "w") as f: + rendered = msw.MetaSwapModel._template.render(settings=settings) + f.write(rendered) + return filename + + +def write_test_files(directory: Path, filenames: list[str]) -> list[Path]: + paths = [directory / filename for filename in filenames] + for p in paths: + with open(p, mode="w") as f: + f.write("test") + return paths + + +def setup_extra_files(meteo_grids: tuple[GridDataArray], directory: Path): + grid_dir = setup_written_meteo_grids(meteo_grids, directory) + setup_parasim_inp(grid_dir) + write_test_files(grid_dir, ["a.inp", "b.inp"]) + return { + "paths": [ + [str(grid_dir / fn)] + for fn in ["a.inp", "mete_grid.inp", "para_sim.inp", "b.inp"] + ] + } + + +def test_import_from_imod5( + imod5_cap_data: Imod5DataDict, + meteo_grids: tuple[GridDataArray], + coupled_mf6_model: mf6.Modflow6Simulation, + tmp_path: Path, +): + # Arrange + imod5_cap_data["extra"] = setup_extra_files(meteo_grids, tmp_path) + times = coupled_mf6_model["time_discretization"].dataset.coords["time"] + dis_pkg = coupled_mf6_model["GWF_1"]["dis"] + # Act + model = msw.MetaSwapModel.from_imod5_data(imod5_cap_data, dis_pkg, times) + # Assert + grid_packages = { + "grid", + "infiltration", + "ponding", + "sprinkling", + "idf_mapping", + } + expected_keys = grid_packages | { + "meteo_grid", + "prec_mapping", + "evt_mapping", + "coupling", + "extra_files", + "time_oc", + } + missing_keys = expected_keys - set(model.keys()) + assert not missing_keys + for pkgname in grid_packages: + # Test if all grid packages broadcasted to grid. + missing_dims = {"y", "x"} - set(model[pkgname].dataset.dims.keys()) + assert not missing_dims + assert "time" in model["time_oc"].dataset.dims.keys() + assert len(model["meteo_grid"].dataset.dims) == 0 + + +def test_import_from_imod5_and_write( + imod5_cap_data: Imod5DataDict, + meteo_grids: tuple[GridDataArray], + coupled_mf6_model: mf6.Modflow6Simulation, + tmp_path: Path, +): + # Arrange + imod5_cap_data["extra"] = setup_extra_files(meteo_grids, tmp_path) + times = coupled_mf6_model["time_discretization"].dataset.coords["time"] + dis_pkg = coupled_mf6_model["GWF_1"]["dis"] + npf_pkg = coupled_mf6_model["GWF_1"]["npf"] + active = dis_pkg["idomain"] == 1 + modeldir = tmp_path / "modeldir" + # Act + model = msw.MetaSwapModel.from_imod5_data(imod5_cap_data, dis_pkg, times) + well_pkg = mf6.LayeredWell.from_imod5_cap_data(imod5_cap_data) + mf6_wel_pkg = well_pkg.to_mf6_pkg( + active, dis_pkg["top"], dis_pkg["bottom"], npf_pkg["k"] + ) + model.write(modeldir, dis_pkg, mf6_wel_pkg, validate=False) + + # Assert + assert len(list(modeldir.rglob(r"*.inp"))) == 13 From 35d4e6f1e44102adc1460de3951c07311be6bf1b Mon Sep 17 00:00:00 2001 From: JoerivanEngelen Date: Fri, 6 Dec 2024 17:41:48 +0100 Subject: [PATCH 31/50] Correctly get all active from MetaSwapActive object --- imod/mf6/rch.py | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/imod/mf6/rch.py b/imod/mf6/rch.py index 392cc9d56..3fd892c90 100644 --- a/imod/mf6/rch.py +++ b/imod/mf6/rch.py @@ -249,7 +249,8 @@ def from_imod5_cap_data( cap_data = cast(GridDataDict, imod5_data["cap"]) msw_area = get_cell_area_from_imod5_data(cap_data) - active, _ = is_msw_active_cell(target_dis, cap_data, msw_area) + msw_active = is_msw_active_cell(target_dis, cap_data, msw_area) + active = msw_active.all data = {} layer_da = xr.full_like( From 82a3cd615caa2fed195c11948bc562e619f6c9e6 Mon Sep 17 00:00:00 2001 From: JoerivanEngelen Date: Fri, 6 Dec 2024 17:45:34 +0100 Subject: [PATCH 32/50] Fix type annotation --- imod/msw/utilities/mask.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/imod/msw/utilities/mask.py b/imod/msw/utilities/mask.py index 62e2ef87e..f28d9c348 100644 --- a/imod/msw/utilities/mask.py +++ b/imod/msw/utilities/mask.py @@ -23,7 +23,7 @@ def mask_and_broadcast_grid_data( def mask_package_data( - package: MetaSwapPackage, grid_data: GridDataDict, msw_active: MetaSwapActive + package: type[MetaSwapPackage], grid_data: GridDataDict, msw_active: MetaSwapActive ) -> GridDataDict: """ Mask and broadcast grid data, carefully mask per subunit if variable needs From 2b2793976ae20a0ef760ee43caf56757bf9a3b05 Mon Sep 17 00:00:00 2001 From: JoerivanEngelen Date: Fri, 6 Dec 2024 17:49:42 +0100 Subject: [PATCH 33/50] Fix some mypy issues --- imod/msw/model.py | 5 +++-- 1 file changed, 3 insertions(+), 2 deletions(-) diff --git a/imod/msw/model.py b/imod/msw/model.py index adae7c972..ec89ff22c 100644 --- a/imod/msw/model.py +++ b/imod/msw/model.py @@ -33,7 +33,7 @@ from imod.msw.utilities.mask import mask_and_broadcast_grid_data from imod.msw.utilities.parse import read_para_sim from imod.msw.vegetation import AnnualCropFactors -from imod.typing import Imod5DataDict, SelSettingsType +from imod.typing import GridDataDict, Imod5DataDict, SelSettingsType from imod.util.regrid_method_type import RegridderType REQUIRED_PACKAGES = ( @@ -343,10 +343,11 @@ def from_imod5_data( parasim_settings = read_para_sim(path_to_parasim) unsa_svat_path = cast(str, parasim_settings["unsa_svat_path"]) # Drop layer coord + cap_data = cast(GridDataDict, imod5_data["cap"]) imod5_cap_no_layer: Imod5DataDict = { "cap": { key: da.isel(**_DROP_LAYER_KWARGS).compute() - for key, da in imod5_data["cap"].items() + for key, da in cap_data.items() } } model = cls(unsa_svat_path, parasim_settings) From de09d3f4cd67f9d257cecb842edac2e1842c2db2 Mon Sep 17 00:00:00 2001 From: JoerivanEngelen Date: Mon, 9 Dec 2024 10:09:56 +0100 Subject: [PATCH 34/50] Fix mypy errors --- imod/msw/infiltration.py | 7 ++++--- imod/msw/model.py | 7 +++++-- imod/msw/ponding.py | 8 ++++---- 3 files changed, 13 insertions(+), 9 deletions(-) diff --git a/imod/msw/infiltration.py b/imod/msw/infiltration.py index ed4f27b23..a2ee77dd8 100644 --- a/imod/msw/infiltration.py +++ b/imod/msw/infiltration.py @@ -1,4 +1,5 @@ from textwrap import dedent +from typing import cast import xarray as xr @@ -8,7 +9,7 @@ from imod.msw.pkgbase import MetaSwapPackage from imod.msw.regrid.regrid_schemes import InfiltrationRegridMethod from imod.msw.utilities.common import concat_imod5 -from imod.typing import GridDataDict +from imod.typing import GridDataDict, Imod5DataDict from imod.typing.grid import ones_like @@ -101,8 +102,8 @@ def __init__( self._pkgcheck() @classmethod - def from_imod5_data(cls, imod5_data: dict[str, GridDataDict]) -> "Infiltration": - cap_data = imod5_data["cap"] + def from_imod5_data(cls, imod5_data: Imod5DataDict) -> "Infiltration": + cap_data = cast(GridDataDict, imod5_data["cap"]) data = {} # Use runon resistance as downward resistance, and runoff for downward # resistance diff --git a/imod/msw/model.py b/imod/msw/model.py index ec89ff22c..872b0d871 100644 --- a/imod/msw/model.py +++ b/imod/msw/model.py @@ -344,7 +344,7 @@ def from_imod5_data( unsa_svat_path = cast(str, parasim_settings["unsa_svat_path"]) # Drop layer coord cap_data = cast(GridDataDict, imod5_data["cap"]) - imod5_cap_no_layer: Imod5DataDict = { + imod5_cap_no_layer: dict[str, GridDataDict] = { "cap": { key: da.isel(**_DROP_LAYER_KWARGS).compute() for key, da in cap_data.items() @@ -357,7 +357,10 @@ def from_imod5_data( cap_data_masked = mask_and_broadcast_grid_data( imod5_cap_no_layer["cap"], msw_active ) - imod5_masked = {"cap": cap_data_masked, "extra": {"paths": extra_paths}} + imod5_masked: Imod5DataDict = { + "cap": cap_data_masked, + "extra": {"paths": extra_paths}, + } model["infiltration"] = Infiltration.from_imod5_data(imod5_masked) model["ponding"] = Ponding.from_imod5_data(imod5_masked) model["sprinkling"] = Sprinkling.from_imod5_data(imod5_masked) diff --git a/imod/msw/ponding.py b/imod/msw/ponding.py index 89ab6af44..8c8272856 100644 --- a/imod/msw/ponding.py +++ b/imod/msw/ponding.py @@ -1,4 +1,4 @@ -from typing import Any, TextIO +from typing import Any, TextIO, cast import pandas as pd import xarray as xr @@ -8,7 +8,7 @@ from imod.msw.pkgbase import DataDictType, MetaSwapPackage from imod.msw.regrid.regrid_schemes import PondingRegridMethod from imod.msw.utilities.common import concat_imod5 -from imod.typing import GridDataDict, IntArray +from imod.typing import GridDataDict, Imod5DataDict, IntArray class Ponding(MetaSwapPackage, IRegridPackage): @@ -73,11 +73,11 @@ def _render(self, file: TextIO, index: IntArray, svat: xr.DataArray, *args: Any) return self.write_dataframe_fixed_width(file, dataframe) @classmethod - def from_imod5_data(cls, imod5_data: dict[str, GridDataDict]) -> "Ponding": + def from_imod5_data(cls, imod5_data: Imod5DataDict) -> "Ponding": """ Concatenate ponding depths along subunits """ - cap_data = imod5_data["cap"] + cap_data = cast(GridDataDict, imod5_data["cap"]) data = {} for key in cls._with_subunit: data_ls = [cap_data[f"{landuse}_{key}"] for landuse in ["rural", "urban"]] From 42bb9642f61cf9c9de62ca2976d5670c3881fc75 Mon Sep 17 00:00:00 2001 From: JoerivanEngelen Date: Mon, 9 Dec 2024 10:47:30 +0100 Subject: [PATCH 35/50] Convert Imod5DataDict to TypedDict and remove casting redundancy --- imod/msw/infiltration.py | 3 +-- imod/msw/meteo_grid.py | 4 ++-- imod/msw/meteo_mapping.py | 4 ++-- imod/msw/model.py | 4 ++-- imod/msw/ponding.py | 6 +++--- imod/msw/scaling_factors.py | 6 ++---- imod/msw/sprinkling.py | 4 ++-- imod/typing/__init__.py | 6 +++++- 8 files changed, 19 insertions(+), 18 deletions(-) diff --git a/imod/msw/infiltration.py b/imod/msw/infiltration.py index a2ee77dd8..62e84a667 100644 --- a/imod/msw/infiltration.py +++ b/imod/msw/infiltration.py @@ -1,5 +1,4 @@ from textwrap import dedent -from typing import cast import xarray as xr @@ -103,7 +102,7 @@ def __init__( @classmethod def from_imod5_data(cls, imod5_data: Imod5DataDict) -> "Infiltration": - cap_data = cast(GridDataDict, imod5_data["cap"]) + cap_data = imod5_data["cap"] data = {} # Use runon resistance as downward resistance, and runoff for downward # resistance diff --git a/imod/msw/meteo_grid.py b/imod/msw/meteo_grid.py index 70d3625a5..00f1716f4 100644 --- a/imod/msw/meteo_grid.py +++ b/imod/msw/meteo_grid.py @@ -1,7 +1,7 @@ import csv from pathlib import Path from shutil import copyfile -from typing import Optional, Union, cast +from typing import Optional, Union import numpy as np import pandas as pd @@ -223,7 +223,7 @@ def write(self, directory: Path | str, *args): @classmethod def from_imod5_data(cls, imod5_data: Imod5DataDict) -> "MeteoGridCopy": - paths = cast(list[str], imod5_data["extra"]["paths"]) + paths = imod5_data["extra"]["paths"] filepath = find_in_file_list(cls._file_name, paths) return cls(filepath) diff --git a/imod/msw/meteo_mapping.py b/imod/msw/meteo_mapping.py index 1de00fa69..b728bb497 100644 --- a/imod/msw/meteo_mapping.py +++ b/imod/msw/meteo_mapping.py @@ -1,7 +1,7 @@ from copy import deepcopy from pathlib import Path from textwrap import dedent -from typing import Any, Optional, TextIO, cast +from typing import Any, Optional, TextIO import numpy as np import pandas as pd @@ -58,7 +58,7 @@ def open_first_meteo_grid(mete_grid_path: str | Path, column_nr: int) -> xr.Data def open_first_meteo_grid_from_imod5_data(imod5_data: Imod5DataDict, column_nr: int): - paths = cast(list[str], imod5_data["extra"]["paths"]) + paths = 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) diff --git a/imod/msw/model.py b/imod/msw/model.py index 872b0d871..21d3d683a 100644 --- a/imod/msw/model.py +++ b/imod/msw/model.py @@ -338,12 +338,12 @@ def from_imod5_data( MetaSwapModel MetaSWAP model imported from imod5 data. """ - extra_paths = cast(list[str], imod5_data["extra"]["paths"]) + extra_paths = imod5_data["extra"]["paths"] path_to_parasim = find_in_file_list("para_sim.inp", extra_paths) parasim_settings = read_para_sim(path_to_parasim) unsa_svat_path = cast(str, parasim_settings["unsa_svat_path"]) # Drop layer coord - cap_data = cast(GridDataDict, imod5_data["cap"]) + cap_data = imod5_data["cap"] imod5_cap_no_layer: dict[str, GridDataDict] = { "cap": { key: da.isel(**_DROP_LAYER_KWARGS).compute() diff --git a/imod/msw/ponding.py b/imod/msw/ponding.py index 8c8272856..86a7c5cff 100644 --- a/imod/msw/ponding.py +++ b/imod/msw/ponding.py @@ -1,4 +1,4 @@ -from typing import Any, TextIO, cast +from typing import Any, TextIO import pandas as pd import xarray as xr @@ -8,7 +8,7 @@ from imod.msw.pkgbase import DataDictType, MetaSwapPackage from imod.msw.regrid.regrid_schemes import PondingRegridMethod from imod.msw.utilities.common import concat_imod5 -from imod.typing import GridDataDict, Imod5DataDict, IntArray +from imod.typing import Imod5DataDict, IntArray class Ponding(MetaSwapPackage, IRegridPackage): @@ -77,7 +77,7 @@ def from_imod5_data(cls, imod5_data: Imod5DataDict) -> "Ponding": """ Concatenate ponding depths along subunits """ - cap_data = cast(GridDataDict, imod5_data["cap"]) + cap_data = imod5_data["cap"] data = {} for key in cls._with_subunit: data_ls = [cap_data[f"{landuse}_{key}"] for landuse in ["rural", "urban"]] diff --git a/imod/msw/scaling_factors.py b/imod/msw/scaling_factors.py index 87dccaf61..91984dd3d 100644 --- a/imod/msw/scaling_factors.py +++ b/imod/msw/scaling_factors.py @@ -1,11 +1,9 @@ -from typing import cast - from imod.mf6.interfaces.iregridpackage import IRegridPackage from imod.msw.fixed_format import VariableMetaData from imod.msw.pkgbase import MetaSwapPackage from imod.msw.regrid.regrid_schemes import ScalingRegridMethod from imod.msw.utilities.common import concat_imod5 -from imod.typing import GridDataDict, Imod5DataDict +from imod.typing import Imod5DataDict from imod.typing.grid import ones_like @@ -88,7 +86,7 @@ def from_imod5_data(cls, imod5_data: Imod5DataDict) -> "ScalingFactors": Import ScalingFactors from iMOD5 data. Pressure head factor is set to one for all factors, as well as all factors for urban areas. """ - cap_data = cast(GridDataDict, imod5_data["cap"]) + cap_data = imod5_data["cap"] grid_ones = ones_like(cap_data["boundary"]) data = {} diff --git a/imod/msw/sprinkling.py b/imod/msw/sprinkling.py index ea9925112..11fb3744b 100644 --- a/imod/msw/sprinkling.py +++ b/imod/msw/sprinkling.py @@ -1,4 +1,4 @@ -from typing import TextIO, cast +from typing import TextIO import numpy as np import pandas as pd @@ -208,7 +208,7 @@ def from_imod5_data(cls, imod5_data: Imod5DataDict) -> "Sprinkling": ------- Sprinkling package """ - cap_data = cast(GridDataDict, imod5_data["cap"]) + cap_data = imod5_data["cap"] if isinstance(cap_data["artificial_recharge_layer"], pd.DataFrame): data = _sprinkling_data_from_imod5_ipf(cap_data) else: diff --git a/imod/typing/__init__.py b/imod/typing/__init__.py index 0217f9782..da8bb6c19 100644 --- a/imod/typing/__init__.py +++ b/imod/typing/__init__.py @@ -12,7 +12,6 @@ GridDataArray: TypeAlias = Union[xr.DataArray, xu.UgridDataArray] GridDataset: TypeAlias = Union[xr.Dataset, xu.UgridDataset] GridDataDict: TypeAlias = dict[str, GridDataArray] -Imod5DataDict: TypeAlias = dict[str, GridDataDict | dict[str, list[str]]] ScalarAsDataArray: TypeAlias = Union[xr.DataArray, xu.UgridDataArray] ScalarAsDataset: TypeAlias = Union[xr.Dataset, xu.UgridDataset] UnstructuredData: TypeAlias = Union[xu.UgridDataset, xu.UgridDataArray] @@ -26,6 +25,11 @@ class SelSettingsType(TypedDict, total=False): missing_dims: Literal["raise", "warn", "ignore"] +class Imod5DataDict(TypedDict, total=False): + cap: GridDataDict + extra: dict[str, list[str]] + + # Types for optional dependencies. if TYPE_CHECKING: import geopandas as gpd From 173e98ec906a180be6e4095ac329b3ae1b5376ac Mon Sep 17 00:00:00 2001 From: JoerivanEngelen Date: Mon, 9 Dec 2024 11:18:14 +0100 Subject: [PATCH 36/50] Move drop layer coord logic to separate utility --- imod/msw/grid_data.py | 4 ++-- imod/msw/model.py | 17 +++-------------- imod/msw/utilities/select.py | 16 ++++++++++++++++ 3 files changed, 21 insertions(+), 16 deletions(-) create mode 100644 imod/msw/utilities/select.py diff --git a/imod/msw/grid_data.py b/imod/msw/grid_data.py index 4c202697c..b617070fe 100644 --- a/imod/msw/grid_data.py +++ b/imod/msw/grid_data.py @@ -13,7 +13,7 @@ is_msw_active_cell, ) from imod.msw.utilities.mask import MetaSwapActive, mask_package_data -from imod.typing import GridDataDict +from imod.typing import Imod5DataDict from imod.util.spatial import get_cell_area, spatial_reference @@ -124,7 +124,7 @@ def _pkgcheck(self): @classmethod def from_imod5_data( cls, - imod5_data: dict[str, GridDataDict], + imod5_data: Imod5DataDict, target_dis: StructuredDiscretization, ) -> tuple["GridData", MetaSwapActive]: # Get iMOD5 capillary zone data diff --git a/imod/msw/model.py b/imod/msw/model.py index 21d3d683a..eb2e0df19 100644 --- a/imod/msw/model.py +++ b/imod/msw/model.py @@ -32,8 +32,9 @@ from imod.msw.utilities.common import find_in_file_list from imod.msw.utilities.mask import mask_and_broadcast_grid_data from imod.msw.utilities.parse import read_para_sim +from imod.msw.utilities.select import drop_layer_dim_cap_data from imod.msw.vegetation import AnnualCropFactors -from imod.typing import GridDataDict, Imod5DataDict, SelSettingsType +from imod.typing import Imod5DataDict from imod.util.regrid_method_type import RegridderType REQUIRED_PACKAGES = ( @@ -78,12 +79,6 @@ "clocktime": 0, } -_DROP_LAYER_KWARGS: SelSettingsType = { - "layer": 0, - "drop": True, - "missing_dims": "ignore", -} - class Model(collections.UserDict): def __setitem__(self, key, value): @@ -343,13 +338,7 @@ def from_imod5_data( parasim_settings = read_para_sim(path_to_parasim) unsa_svat_path = cast(str, parasim_settings["unsa_svat_path"]) # Drop layer coord - cap_data = imod5_data["cap"] - imod5_cap_no_layer: dict[str, GridDataDict] = { - "cap": { - key: da.isel(**_DROP_LAYER_KWARGS).compute() - for key, da in cap_data.items() - } - } + imod5_cap_no_layer = drop_layer_dim_cap_data(imod5_data) model = cls(unsa_svat_path, parasim_settings) model["grid"], msw_active = GridData.from_imod5_data( imod5_cap_no_layer, target_dis diff --git a/imod/msw/utilities/select.py b/imod/msw/utilities/select.py new file mode 100644 index 000000000..abd28818f --- /dev/null +++ b/imod/msw/utilities/select.py @@ -0,0 +1,16 @@ +from imod.typing import Imod5DataDict, SelSettingsType + +_DROP_LAYER_KWARGS: SelSettingsType = { + "layer": 0, + "drop": True, + "missing_dims": "ignore", +} + + +def drop_layer_dim_cap_data(imod5_data: Imod5DataDict) -> Imod5DataDict: + cap_data = imod5_data["cap"] + return { + "cap": { + key: da.isel(**_DROP_LAYER_KWARGS).compute() for key, da in cap_data.items() + } + } From 11ad02dcdd8bce14a1c4c74ecd16288290856e57 Mon Sep 17 00:00:00 2001 From: JoerivanEngelen Date: Mon, 9 Dec 2024 15:07:33 +0100 Subject: [PATCH 37/50] Better function names --- imod/msw/grid_data.py | 4 ++-- imod/msw/model.py | 4 ++-- imod/msw/utilities/mask.py | 10 +++++----- 3 files changed, 9 insertions(+), 9 deletions(-) diff --git a/imod/msw/grid_data.py b/imod/msw/grid_data.py index b617070fe..a7ad1e7c0 100644 --- a/imod/msw/grid_data.py +++ b/imod/msw/grid_data.py @@ -12,7 +12,7 @@ get_rootzone_depth_from_imod5_data, is_msw_active_cell, ) -from imod.msw.utilities.mask import MetaSwapActive, mask_package_data +from imod.msw.utilities.mask import MetaSwapActive, mask_and_broadcast_pkg_data from imod.typing import Imod5DataDict from imod.util.spatial import get_cell_area, spatial_reference @@ -138,6 +138,6 @@ def from_imod5_data( data["soil_physical_unit"] = imod5_cap["soil_physical_unit"].astype(int) msw_active = is_msw_active_cell(target_dis, imod5_cap, data["area"]) - data_active = mask_package_data(cls, data, msw_active) + data_active = mask_and_broadcast_pkg_data(cls, data, msw_active) data_active["active"] = msw_active.all return cls(**data_active), msw_active diff --git a/imod/msw/model.py b/imod/msw/model.py index eb2e0df19..8c1c6a00b 100644 --- a/imod/msw/model.py +++ b/imod/msw/model.py @@ -30,7 +30,7 @@ from imod.msw.sprinkling import Sprinkling from imod.msw.timeutil import to_metaswap_timeformat from imod.msw.utilities.common import find_in_file_list -from imod.msw.utilities.mask import mask_and_broadcast_grid_data +from imod.msw.utilities.mask import mask_and_broadcast_cap_data from imod.msw.utilities.parse import read_para_sim from imod.msw.utilities.select import drop_layer_dim_cap_data from imod.msw.vegetation import AnnualCropFactors @@ -343,7 +343,7 @@ def from_imod5_data( model["grid"], msw_active = GridData.from_imod5_data( imod5_cap_no_layer, target_dis ) - cap_data_masked = mask_and_broadcast_grid_data( + cap_data_masked = mask_and_broadcast_cap_data( imod5_cap_no_layer["cap"], msw_active ) imod5_masked: Imod5DataDict = { diff --git a/imod/msw/utilities/mask.py b/imod/msw/utilities/mask.py index f28d9c348..7d083d6eb 100644 --- a/imod/msw/utilities/mask.py +++ b/imod/msw/utilities/mask.py @@ -11,18 +11,18 @@ class MetaSwapActive: per_subunit: GridDataArray -def mask_and_broadcast_grid_data( - grid_data: GridDataDict, msw_active: MetaSwapActive +def mask_and_broadcast_cap_data( + cap_data: GridDataDict, msw_active: MetaSwapActive ) -> GridDataDict: """ - Mask and broadcast grid data, + Mask and broadcast cap data, always mask with "all" of MetaSwapActive. """ return { - key: _mask_spatial_var(grid, msw_active.all) for key, grid in grid_data.items() + key: _mask_spatial_var(grid, msw_active.all) for key, grid in cap_data.items() } -def mask_package_data( +def mask_and_broadcast_pkg_data( package: type[MetaSwapPackage], grid_data: GridDataDict, msw_active: MetaSwapActive ) -> GridDataDict: """ From 50078c04d3479e98614f19390d6b942073a8bd20 Mon Sep 17 00:00:00 2001 From: JoerivanEngelen Date: Mon, 9 Dec 2024 15:08:21 +0100 Subject: [PATCH 38/50] Add tests for mask --- .../test_msw/test_utilities/test_mask.py | 131 ++++++++++++++++++ 1 file changed, 131 insertions(+) create mode 100644 imod/tests/test_msw/test_utilities/test_mask.py diff --git a/imod/tests/test_msw/test_utilities/test_mask.py b/imod/tests/test_msw/test_utilities/test_mask.py new file mode 100644 index 000000000..424f7dbc7 --- /dev/null +++ b/imod/tests/test_msw/test_utilities/test_mask.py @@ -0,0 +1,131 @@ +from copy import deepcopy +from typing import Callable + +import numpy as np +import pytest +import xarray as xr +from pytest_cases import parametrize_with_cases + +from imod.msw import GridData +from imod.msw.utilities.mask import ( + MetaSwapActive, + mask_and_broadcast_cap_data, + mask_and_broadcast_pkg_data, +) +from imod.typing import GridDataDict + + +@pytest.fixture(scope="function") +def coords_planar() -> dict: + x = [1.0, 2.0, 3.0] + y = [3.0, 2.0, 1.0] + dx = 1.0 + dy = 1.0 + return {"y": y, "x": x, "dx": dx, "dy": dy} + + +@pytest.fixture(scope="function") +def coords_subunit(coords_planar: dict) -> dict: + coords_subunit = deepcopy(coords_planar) + coords_subunit["subunit"] = [0, 1] + return coords_subunit + + +@pytest.fixture(scope="function") +def mask_fixture(coords_subunit: dict) -> MetaSwapActive: + mask_per_subunit = xr.DataArray( + np.array( + [ + [[0, 0, 0], [0, 1, 1], [0, 0, 0]], + [[1, 1, 1], [0, 1, 1], [0, 0, 0]], + ] + ).astype(bool), + dims=("subunit", "y", "x"), + coords=coords_subunit, + ) + mask_all = mask_per_subunit.any(dim="subunit") + + return MetaSwapActive(mask_all, mask_per_subunit) + + +@pytest.fixture(scope="function") +def data_subunit(coords_subunit: dict): + return xr.DataArray( + np.array( + [ + [[2.0, 2.0, 2.0], [2.0, 3.0, 3.0], [2.0, 2.0, 2.0]], + [[3.0, 3.0, 3.0], [2.0, 3.0, 3.0], [2.0, 2.0, 2.0]], + ] + ), + dims=("subunit", "y", "x"), + coords=coords_subunit, + ) + + +@pytest.fixture(scope="function") +def data_planar(data_subunit: xr.DataArray) -> xr.DataArray: + return data_subunit.sel(subunit=0, drop=True) + + +def integer_nan(da): + return da == 0 + + +def case_grid_float( + data_subunit: xr.DataArray, data_planar: xr.DataArray +) -> tuple[GridDataDict, Callable]: + return {"landuse": data_subunit, "elevation": data_planar}, np.isnan + + +def case_grid_integer( + data_subunit: xr.DataArray, data_planar: xr.DataArray +) -> tuple[GridDataDict, Callable]: + return { + "landuse": data_subunit.astype(int), + "elevation": data_planar.astype(int), + }, integer_nan + + +def case_constant_float(): + data = xr.DataArray(2.0) + return {"landuse": data, "elevation": data}, np.isnan + + +def case_constant_integer(): + data = xr.DataArray(2) + return {"landuse": data, "elevation": data}, integer_nan + + +@parametrize_with_cases("imod5_grid_data, isnan", cases=".") +def test_mask_and_broadcast_pkg_data( + imod5_grid_data: GridDataDict, isnan: Callable, mask_fixture: MetaSwapActive +): + # Act + masked_data = mask_and_broadcast_pkg_data(GridData, imod5_grid_data, mask_fixture) + + # Assert + assert isnan(masked_data["landuse"]).to_numpy().sum() == 11 + assert isnan(masked_data["elevation"]).to_numpy().sum() == 4 + + +@parametrize_with_cases("imod5_grid_data, isnan", cases=".") +def test_mask_and_broadcast_cap_data( + imod5_grid_data: GridDataDict, + isnan: Callable, + mask_fixture: MetaSwapActive, + coords_planar: dict, +): + # Arrange, cap data doesn't have subunit coords + imod5_grid_data["landuse"] = imod5_grid_data["landuse"].isel( + subunit=0, drop=True, missing_dims="ignore" + ) + + # Act + masked_data = mask_and_broadcast_cap_data(imod5_grid_data, mask_fixture) + + # Assert + coords_expected = xr.Coordinates(coords_planar) + xr.testing.assert_equal(masked_data["landuse"].coords, coords_expected) + xr.testing.assert_equal(masked_data["elevation"].coords, coords_expected) + assert isnan(masked_data["landuse"]).to_numpy().sum() == 4 + assert isnan(masked_data["elevation"]).to_numpy().sum() == 4 From 764c064b7982fae100aca69d1827609d0a801119 Mon Sep 17 00:00:00 2001 From: JoerivanEngelen Date: Tue, 10 Dec 2024 10:40:45 +0100 Subject: [PATCH 39/50] Add scaling factor if active. --- imod/msw/model.py | 5 ++- imod/msw/utilities/imod5_converter.py | 31 +++++++++++++ imod/tests/test_msw/test_model.py | 63 +++++++++++++++++++++++---- 3 files changed, 90 insertions(+), 9 deletions(-) diff --git a/imod/msw/model.py b/imod/msw/model.py index 8c1c6a00b..ee9dd731a 100644 --- a/imod/msw/model.py +++ b/imod/msw/model.py @@ -27,9 +27,11 @@ from imod.msw.meteo_mapping import EvapotranspirationMapping, PrecipitationMapping from imod.msw.output_control import TimeOutputControl from imod.msw.ponding import Ponding +from imod.msw.scaling_factors import ScalingFactors from imod.msw.sprinkling import Sprinkling from imod.msw.timeutil import to_metaswap_timeformat from imod.msw.utilities.common import find_in_file_list +from imod.msw.utilities.imod5_converter import has_active_scaling_factor from imod.msw.utilities.mask import mask_and_broadcast_cap_data from imod.msw.utilities.parse import read_para_sim from imod.msw.utilities.select import drop_layer_dim_cap_data @@ -356,7 +358,8 @@ def from_imod5_data( model["meteo_grid"] = MeteoGridCopy.from_imod5_data(imod5_masked) model["prec_mapping"] = PrecipitationMapping.from_imod5_data(imod5_masked) model["evt_mapping"] = EvapotranspirationMapping.from_imod5_data(imod5_masked) - # model["scaling_factor"] = ScalingFactors.from_imod5_data(imod5_masked) + if has_active_scaling_factor(imod5_cap_no_layer["cap"]): + model["scaling_factor"] = ScalingFactors.from_imod5_data(imod5_masked) area = model["grid"]["area"].isel(subunit=0, drop=True) model["idf_mapping"] = IdfMapping(area, -9999.0) model["coupling"] = CouplerMapping() diff --git a/imod/msw/utilities/imod5_converter.py b/imod/msw/utilities/imod5_converter.py index fc5066a78..2da4fe9c2 100644 --- a/imod/msw/utilities/imod5_converter.py +++ b/imod/msw/utilities/imod5_converter.py @@ -1,3 +1,5 @@ +from xarray.core.utils import is_scalar + from imod.logging import LogLevel, logger from imod.mf6 import StructuredDiscretization from imod.msw.utilities.common import concat_imod5 @@ -86,3 +88,32 @@ def is_msw_active_cell( ) active = subunit_active.any(dim="subunit") return MetaSwapActive(active, subunit_active) + + +def _is_equal_scalar_value(da, value): + """ + Helper function to guarantee that the check in ``has_active_scaling_factor`` + can shortcut after is_scalar returns False. + """ + return da.to_numpy()[()] == value + + +def has_active_scaling_factor(imod5_cap: GridDataDict): + """ + Check if scaling factor grids are active. Carefully checks if data is + provided as constant (scalar) and if it matches an inactivity value. The + function shortcuts if data is provided as constant. + """ + variable_inactive_mapping = { + "perched_water_table_level": -9999.0, + "soil_moisture_fraction": 1.0, + "conductivitiy_factor": 1.0, + } + scaling_factor_inactive = True + for var, inactive_value in variable_inactive_mapping.items(): + da = imod5_cap[var] + scaling_factor_inactive &= is_scalar(da) and _is_equal_scalar_value( + da, inactive_value + ) + + return not scaling_factor_inactive diff --git a/imod/tests/test_msw/test_model.py b/imod/tests/test_msw/test_model.py index e804a1560..f0904e96f 100644 --- a/imod/tests/test_msw/test_model.py +++ b/imod/tests/test_msw/test_model.py @@ -4,6 +4,7 @@ import pytest import xarray as xr from numpy.testing import assert_almost_equal, assert_equal +from pytest_cases import parametrize_with_cases from imod import mf6, msw from imod.mf6.utilities.regrid import RegridderWeightsCache @@ -187,18 +188,58 @@ def setup_extra_files(meteo_grids: tuple[GridDataArray], directory: Path): } +class Imod5DataCases: + def case_grid(self, imod5_cap_data: Imod5DataDict) -> tuple[Imod5DataDict, bool]: + has_scaling_factor = True + return imod5_cap_data, has_scaling_factor + + def case_no_scaling_factors( + self, imod5_cap_data: Imod5DataDict + ) -> tuple[Imod5DataDict, bool]: + has_scaling_factor = False + cap_data = imod5_cap_data["cap"] + # open_projectfile_data adds layer kwargs to constants + layer_kwargs = {"coords": {"layer": [1]}, "dims": ("layer",)} + cap_data["perched_water_table_level"] = xr.DataArray([-9999.0], **layer_kwargs) + cap_data["soil_moisture_fraction"] = xr.DataArray([1.0], **layer_kwargs) + cap_data["conductivitiy_factor"] = xr.DataArray([1.0], **layer_kwargs) + return imod5_cap_data, has_scaling_factor + + def case_constants( + self, imod5_cap_data: Imod5DataDict + ) -> tuple[Imod5DataDict, bool]: + has_scaling_factor = False + cap_data = imod5_cap_data["cap"] + # open_projectfile_data adds layer kwargs to constants + layer_kwargs = {"coords": {"layer": [1]}, "dims": ("layer",)} + cap_data["perched_water_table_level"] = xr.DataArray([-9999.0], **layer_kwargs) + cap_data["soil_moisture_fraction"] = xr.DataArray([1.0], **layer_kwargs) + cap_data["conductivitiy_factor"] = xr.DataArray([1.0], **layer_kwargs) + cap_data["urban_ponding_depth"] = xr.DataArray([1.0], **layer_kwargs) + cap_data["rural_ponding_depth"] = xr.DataArray([1.0], **layer_kwargs) + cap_data["urban_runoff_resistance"] = xr.DataArray([1.0], **layer_kwargs) + cap_data["rural_runoff_resistance"] = xr.DataArray([1.0], **layer_kwargs) + cap_data["urban_runon_resistance"] = xr.DataArray([1.0], **layer_kwargs) + cap_data["rural_runon_resistance"] = xr.DataArray([1.0], **layer_kwargs) + cap_data["urban_infiltration_capacity"] = xr.DataArray([1.0], **layer_kwargs) + cap_data["rural_infiltration_capacity"] = xr.DataArray([1.0], **layer_kwargs) + return imod5_cap_data, has_scaling_factor + + +@parametrize_with_cases("imod5_data, has_scaling_factor", cases=Imod5DataCases) def test_import_from_imod5( - imod5_cap_data: Imod5DataDict, + imod5_data: Imod5DataDict, + has_scaling_factor: bool, meteo_grids: tuple[GridDataArray], coupled_mf6_model: mf6.Modflow6Simulation, tmp_path: Path, ): # Arrange - imod5_cap_data["extra"] = setup_extra_files(meteo_grids, tmp_path) + imod5_data["extra"] = setup_extra_files(meteo_grids, tmp_path) times = coupled_mf6_model["time_discretization"].dataset.coords["time"] dis_pkg = coupled_mf6_model["GWF_1"]["dis"] # Act - model = msw.MetaSwapModel.from_imod5_data(imod5_cap_data, dis_pkg, times) + model = msw.MetaSwapModel.from_imod5_data(imod5_data, dis_pkg, times) # Assert grid_packages = { "grid", @@ -223,28 +264,34 @@ def test_import_from_imod5( assert not missing_dims assert "time" in model["time_oc"].dataset.dims.keys() assert len(model["meteo_grid"].dataset.dims) == 0 + assert ("scaling_factor" in model.keys()) == has_scaling_factor +@parametrize_with_cases("imod5_data, has_scaling_factor", cases=Imod5DataCases) def test_import_from_imod5_and_write( - imod5_cap_data: Imod5DataDict, + imod5_data: Imod5DataDict, + has_scaling_factor: bool, meteo_grids: tuple[GridDataArray], coupled_mf6_model: mf6.Modflow6Simulation, tmp_path: Path, ): # Arrange - imod5_cap_data["extra"] = setup_extra_files(meteo_grids, tmp_path) + imod5_data["extra"] = setup_extra_files(meteo_grids, tmp_path) times = coupled_mf6_model["time_discretization"].dataset.coords["time"] dis_pkg = coupled_mf6_model["GWF_1"]["dis"] npf_pkg = coupled_mf6_model["GWF_1"]["npf"] active = dis_pkg["idomain"] == 1 modeldir = tmp_path / "modeldir" # Act - model = msw.MetaSwapModel.from_imod5_data(imod5_cap_data, dis_pkg, times) - well_pkg = mf6.LayeredWell.from_imod5_cap_data(imod5_cap_data) + model = msw.MetaSwapModel.from_imod5_data(imod5_data, dis_pkg, times) + well_pkg = mf6.LayeredWell.from_imod5_cap_data(imod5_data) mf6_wel_pkg = well_pkg.to_mf6_pkg( active, dis_pkg["top"], dis_pkg["bottom"], npf_pkg["k"] ) model.write(modeldir, dis_pkg, mf6_wel_pkg, validate=False) # Assert - assert len(list(modeldir.rglob(r"*.inp"))) == 13 + expected_n_files = 13 + if has_scaling_factor: + expected_n_files += 1 + assert len(list(modeldir.rglob(r"*.inp"))) == expected_n_files From c8f9ec7b8144dcc5b07751f166adc1e5beecdd35 Mon Sep 17 00:00:00 2001 From: JoerivanEngelen Date: Tue, 10 Dec 2024 10:58:51 +0100 Subject: [PATCH 40/50] Add from_imod5_data method to public API --- docs/api/changelog.rst | 2 ++ docs/api/msw.rst | 10 ++++++++++ 2 files changed, 12 insertions(+) diff --git a/docs/api/changelog.rst b/docs/api/changelog.rst index 1f073881c..fae50cb27 100644 --- a/docs/api/changelog.rst +++ b/docs/api/changelog.rst @@ -21,6 +21,8 @@ Added MetaSWAP). Currently only griddata (IDF) is supported. - :meth:`imod.mf6.Recharge.from_imod5_cap_data` to construct a recharge package for coupling a MODFLOW6 model to MetaSWAP. +- :meth:`imod.msw.MetaSwapModel.from_imod5_data` to construct a MetaSWAP model + from data in an iMOD5 projectfile. Fixed ~~~~~ diff --git a/docs/api/msw.rst b/docs/api/msw.rst index f2260ed0f..f96b64f17 100644 --- a/docs/api/msw.rst +++ b/docs/api/msw.rst @@ -7,10 +7,15 @@ MetaSWAP :toctree: generated/msw GridData + GridData.from_imod5_data Infiltration + Infiltration.from_imod5_data Ponding + Ponding.from_imod5_data ScalingFactors + ScalingFactors.from_imod5_data Sprinkling + Sprinkling.from_imod5_data IdfMapping TimeOutputControl @@ -25,9 +30,14 @@ MetaSWAP AnnualCropFactors MeteoGrid + MeteoGridCopy + MeteoGridCopy.from_imod5_data EvapotranspirationMapping + EvapotranspirationMapping.from_imod5_data PrecipitationMapping + PrecipitationMapping.from_imod5_data CouplerMapping MetaSwapModel + MetaSwapModel.from_imod5_data From ef7bbf0513617bf7e24ce8746e088989c186ae92 Mon Sep 17 00:00:00 2001 From: JoerivanEngelen Date: Tue, 10 Dec 2024 13:37:25 +0100 Subject: [PATCH 41/50] Update API docs --- docs/api/changelog.rst | 4 ++++ docs/api/msw.rst | 1 + 2 files changed, 5 insertions(+) diff --git a/docs/api/changelog.rst b/docs/api/changelog.rst index fae50cb27..f3927430e 100644 --- a/docs/api/changelog.rst +++ b/docs/api/changelog.rst @@ -23,6 +23,10 @@ Added for coupling a MODFLOW6 model to MetaSWAP. - :meth:`imod.msw.MetaSwapModel.from_imod5_data` to construct a MetaSWAP model from data in an iMOD5 projectfile. +- :meth:`imod.msw.MetaSwapModel.write` has a ``validate`` argument, which can be + used to turn off validation upon writing, use at your own risk! +- :class:`imod.msw.MetaSwapModel` got ``settings`` argument to set simulation + settings. Fixed ~~~~~ diff --git a/docs/api/msw.rst b/docs/api/msw.rst index f96b64f17..067710445 100644 --- a/docs/api/msw.rst +++ b/docs/api/msw.rst @@ -40,4 +40,5 @@ MetaSWAP CouplerMapping MetaSwapModel + MetaSwapModel.write MetaSwapModel.from_imod5_data From 4ef85f0489eba210f110de99da6cb3bc764bf263 Mon Sep 17 00:00:00 2001 From: JoerivanEngelen Date: Tue, 10 Dec 2024 14:01:48 +0100 Subject: [PATCH 42/50] Fix failing unittests --- imod/msw/model.py | 1 - imod/tests/test_msw/test_grid_data.py | 2 +- imod/tests/test_msw/test_infiltration.py | 10 +++++----- 3 files changed, 6 insertions(+), 7 deletions(-) diff --git a/imod/msw/model.py b/imod/msw/model.py index ee9dd731a..20378e52f 100644 --- a/imod/msw/model.py +++ b/imod/msw/model.py @@ -273,7 +273,6 @@ def write( # write package contents for pkgname in self: - print(pkgname) self[pkgname].write(directory, index, svat, mf6_dis, mf6_wel) def regrid_like( diff --git a/imod/tests/test_msw/test_grid_data.py b/imod/tests/test_msw/test_grid_data.py index 0d0fbf0c3..a775a93a9 100644 --- a/imod/tests/test_msw/test_grid_data.py +++ b/imod/tests/test_msw/test_grid_data.py @@ -461,7 +461,7 @@ def test_from_imod5_data(grid_data_dict: dict[str, xr.DataArray]): # Dis only needed for idomain dis = StructuredDiscretization(top, top - 0.1, idomain, validate=False) - griddata = GridData.from_imod5_data(imod5_data, target_dis=dis) + griddata, _ = GridData.from_imod5_data(imod5_data, target_dis=dis) expected_rootzone_depth = cap_data["rootzone_thickness"] * 0.01 xr.testing.assert_allclose( expected_rootzone_depth, griddata["rootzone_depth"].sel(subunit=0, drop=True) diff --git a/imod/tests/test_msw/test_infiltration.py b/imod/tests/test_msw/test_infiltration.py index b87940a88..47990e3a6 100644 --- a/imod/tests/test_msw/test_infiltration.py +++ b/imod/tests/test_msw/test_infiltration.py @@ -283,11 +283,11 @@ def test_from_imod5_data(data_infiltration): imod5_data = {"cap": cap_data} actual_pkg = Infiltration.from_imod5_data(imod5_data) - ones_vars = ["bottom_resistance", "extra_storage_coefficient"] - expected = expected_pkg.dataset.drop_vars(ones_vars) - actual = actual_pkg.dataset.drop_vars(ones_vars) + hardcoded_vars = {"bottom_resistance": -9999.0, "extra_storage_coefficient": 1.0} + expected = expected_pkg.dataset.drop_vars(hardcoded_vars.keys()) + actual = actual_pkg.dataset.drop_vars(hardcoded_vars.keys()) xr.testing.assert_equal(actual, expected) - for var in ones_vars: - assert (actual_pkg.dataset[var] == 1.0).all() + for var, value in hardcoded_vars.items(): + assert (actual_pkg.dataset[var] == value).all() From 628fab443bad9eda7168532cc7d3671dd182b095 Mon Sep 17 00:00:00 2001 From: JoerivanEngelen Date: Tue, 10 Dec 2024 14:18:43 +0100 Subject: [PATCH 43/50] Add docstring --- imod/msw/grid_data.py | 29 ++++++++++++++++++++++++++++- 1 file changed, 28 insertions(+), 1 deletion(-) diff --git a/imod/msw/grid_data.py b/imod/msw/grid_data.py index a7ad1e7c0..11e24e38f 100644 --- a/imod/msw/grid_data.py +++ b/imod/msw/grid_data.py @@ -127,7 +127,34 @@ def from_imod5_data( imod5_data: Imod5DataDict, target_dis: StructuredDiscretization, ) -> tuple["GridData", MetaSwapActive]: - # Get iMOD5 capillary zone data + """ + Construct a MetaSWAP GridData package from iMOD5 data in the CAP + package, loaded with the :func:`imod.formats.prj.open_projectfile_data` + function. + + Method does the following things: + - Computes rural area from the wetted area and urban area grids. + - Sets a second subunit for urban landuse (== 18). + - Rootzone depth is converted from cm's to m's. + - Determines where MetaSWAP should be active based on area grids, + boundary grid, and MODFLOW6 idomain. + + Parameters + ---------- + imod5_data: Imod5DataDict + iMOD5 data as returned by + :func:`imod.formats.prj.open_projectfile_data` + target_dis: imod.mf6.StructuredDiscretization + MODFLOW6 discretization to couple the MetaSWAP model to. + + Returns + ------- + imod.msw.GridData + GridData package + MetaSwapActive + Dataclass containing where MetaSwAP is active, per subunit, as well + as aggregated over subunits. + """ imod5_cap = imod5_data["cap"] data = {} From 4c673cc9a331c3a55cb386577be961f079bd5f00 Mon Sep 17 00:00:00 2001 From: JoerivanEngelen Date: Tue, 10 Dec 2024 14:18:54 +0100 Subject: [PATCH 44/50] Format --- imod/msw/grid_data.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/imod/msw/grid_data.py b/imod/msw/grid_data.py index 11e24e38f..26237059f 100644 --- a/imod/msw/grid_data.py +++ b/imod/msw/grid_data.py @@ -130,11 +130,11 @@ def from_imod5_data( """ Construct a MetaSWAP GridData package from iMOD5 data in the CAP package, loaded with the :func:`imod.formats.prj.open_projectfile_data` - function. + function. Method does the following things: - Computes rural area from the wetted area and urban area grids. - - Sets a second subunit for urban landuse (== 18). + - Sets a second subunit for urban landuse (== 18). - Rootzone depth is converted from cm's to m's. - Determines where MetaSWAP should be active based on area grids, boundary grid, and MODFLOW6 idomain. From 3ae0c8d603afcb7b18908fbb4c642aad1d19ffc8 Mon Sep 17 00:00:00 2001 From: JoerivanEngelen Date: Thu, 12 Dec 2024 08:08:13 +0100 Subject: [PATCH 45/50] Add MaskValues dataclass to store sentinel values --- imod/msw/infiltration.py | 5 +++-- imod/msw/meteo_grid.py | 5 ++++- imod/msw/model.py | 4 ++-- imod/msw/utilities/imod5_converter.py | 4 ++-- imod/msw/utilities/mask.py | 11 ++++++++++- 5 files changed, 21 insertions(+), 8 deletions(-) diff --git a/imod/msw/infiltration.py b/imod/msw/infiltration.py index 62e84a667..7e5472c7b 100644 --- a/imod/msw/infiltration.py +++ b/imod/msw/infiltration.py @@ -8,6 +8,7 @@ from imod.msw.pkgbase import MetaSwapPackage from imod.msw.regrid.regrid_schemes import InfiltrationRegridMethod from imod.msw.utilities.common import concat_imod5 +from imod.msw.utilities.mask import MaskValues from imod.typing import GridDataDict, Imod5DataDict from imod.typing.grid import ones_like @@ -28,7 +29,7 @@ def deactivate_small_resistances_in_data(data: GridDataDict): message=message.format(var=var), additional_depth=1, ) - data[var] = data[var].where(~to_deactivate, -9999.0) + data[var] = data[var].where(~to_deactivate, MaskValues.default) return data @@ -120,7 +121,7 @@ def from_imod5_data(cls, imod5_data: Imod5DataDict) -> "Infiltration": data = deactivate_small_resistances_in_data(data) like = data["downward_resistance"].isel(subunit=0, drop=True) - data["bottom_resistance"] = ones_like(like) * -9999.0 + data["bottom_resistance"] = ones_like(like) * MaskValues.default data["extra_storage_coefficient"] = ones_like(like) return cls(**data) diff --git a/imod/msw/meteo_grid.py b/imod/msw/meteo_grid.py index 00f1716f4..3974eb838 100644 --- a/imod/msw/meteo_grid.py +++ b/imod/msw/meteo_grid.py @@ -13,6 +13,7 @@ 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.msw.utilities.mask import MaskValues from imod.typing import Imod5DataDict from imod.util.regrid_method_type import EmptyRegridMethod, RegridMethodType @@ -176,7 +177,9 @@ def write(self, directory: Union[str, Path], *args): path = (directory / self._meteo_dirname / str(varname)).with_suffix( ".asc" ) - imod.rasterio.save(path, self.dataset[str(varname)], nodata=-9999.0) + imod.rasterio.save( + path, self.dataset[str(varname)], nodata=MaskValues.default + ) def _pkgcheck(self): for varname in self.dataset.data_vars: diff --git a/imod/msw/model.py b/imod/msw/model.py index 20378e52f..6075ca59a 100644 --- a/imod/msw/model.py +++ b/imod/msw/model.py @@ -32,7 +32,7 @@ from imod.msw.timeutil import to_metaswap_timeformat from imod.msw.utilities.common import find_in_file_list from imod.msw.utilities.imod5_converter import has_active_scaling_factor -from imod.msw.utilities.mask import mask_and_broadcast_cap_data +from imod.msw.utilities.mask import MaskValues, mask_and_broadcast_cap_data from imod.msw.utilities.parse import read_para_sim from imod.msw.utilities.select import drop_layer_dim_cap_data from imod.msw.vegetation import AnnualCropFactors @@ -360,7 +360,7 @@ def from_imod5_data( if has_active_scaling_factor(imod5_cap_no_layer["cap"]): model["scaling_factor"] = ScalingFactors.from_imod5_data(imod5_masked) area = model["grid"]["area"].isel(subunit=0, drop=True) - model["idf_mapping"] = IdfMapping(area, -9999.0) + model["idf_mapping"] = IdfMapping(area, MaskValues.default) model["coupling"] = CouplerMapping() model["extra_files"] = FileCopier.from_imod5_data(imod5_masked) diff --git a/imod/msw/utilities/imod5_converter.py b/imod/msw/utilities/imod5_converter.py index 2da4fe9c2..0c4d21604 100644 --- a/imod/msw/utilities/imod5_converter.py +++ b/imod/msw/utilities/imod5_converter.py @@ -3,7 +3,7 @@ from imod.logging import LogLevel, logger from imod.mf6 import StructuredDiscretization from imod.msw.utilities.common import concat_imod5 -from imod.msw.utilities.mask import MetaSwapActive +from imod.msw.utilities.mask import MaskValues, MetaSwapActive from imod.typing import GridDataArray, GridDataDict from imod.typing.grid import ones_like from imod.util.spatial import get_cell_area @@ -105,7 +105,7 @@ def has_active_scaling_factor(imod5_cap: GridDataDict): function shortcuts if data is provided as constant. """ variable_inactive_mapping = { - "perched_water_table_level": -9999.0, + "perched_water_table_level": MaskValues.default, "soil_moisture_fraction": 1.0, "conductivitiy_factor": 1.0, } diff --git a/imod/msw/utilities/mask.py b/imod/msw/utilities/mask.py index 7d083d6eb..6a26c9b9d 100644 --- a/imod/msw/utilities/mask.py +++ b/imod/msw/utilities/mask.py @@ -11,6 +11,15 @@ class MetaSwapActive: per_subunit: GridDataArray +@dataclass +class MaskValues: + """Stores sentinel values for nodata. Most cases use -9999.0, but exceptions + can be added here.""" + + default = -9999.0 + integer = 0 + + def mask_and_broadcast_cap_data( cap_data: GridDataDict, msw_active: MetaSwapActive ) -> GridDataDict: @@ -42,7 +51,7 @@ def mask_and_broadcast_pkg_data( def _mask_spatial_var(da: GridDataArray, active: GridDataArray) -> GridDataArray: if issubclass(da.dtype.type, numbers.Integral): - return da.where(active, other=0) + return da.where(active, other=MaskValues.integer) elif issubclass(da.dtype.type, numbers.Real): return da.where(active) else: From ea86e4da27e8886e77fb92da6c408d7d95c38f20 Mon Sep 17 00:00:00 2001 From: JoerivanEngelen Date: Thu, 12 Dec 2024 08:16:53 +0100 Subject: [PATCH 46/50] Avoid nested try except --- imod/msw/utilities/parse.py | 13 +++++++------ 1 file changed, 7 insertions(+), 6 deletions(-) diff --git a/imod/msw/utilities/parse.py b/imod/msw/utilities/parse.py index 6177775e6..3b67ec429 100644 --- a/imod/msw/utilities/parse.py +++ b/imod/msw/utilities/parse.py @@ -13,13 +13,14 @@ def _try_parsing_string_to_number(s: str) -> ScalarType: "a" -> "a" """ try: - value: ScalarType = int(s) + return int(s) except ValueError: - try: - value = float(s) - except ValueError: - value = s - return value + pass + try: + return float(s) + except ValueError: + pass + return s def read_para_sim(file: Path | str) -> dict[str, ScalarType]: From d85c1127357f3aed49281f6191dec73337023cdd Mon Sep 17 00:00:00 2001 From: JoerivanEngelen Date: Thu, 12 Dec 2024 08:21:36 +0100 Subject: [PATCH 47/50] Better varname --- imod/msw/utilities/parse.py | 10 ++++++---- 1 file changed, 6 insertions(+), 4 deletions(-) diff --git a/imod/msw/utilities/parse.py b/imod/msw/utilities/parse.py index 3b67ec429..efe7f4428 100644 --- a/imod/msw/utilities/parse.py +++ b/imod/msw/utilities/parse.py @@ -28,9 +28,11 @@ def read_para_sim(file: Path | str) -> dict[str, ScalarType]: lines = f.readlines() out = {} for line in lines: - ll = line[0 : line.find("!")].split("=") - if len(ll) > 1: - key = ll[0].strip() - value = _try_parsing_string_to_number(ll[1].strip()) + # Strip comments starting with "!" and split keys from values at the + # equals sign. + key_values = line[0 : line.find("!")].split("=") + if len(key_values) > 1: + key = key_values[0].strip() + value = _try_parsing_string_to_number(key_values[1].strip()) out[key] = value return out From 78d1e7a6e2c1a6a76ebf135fd1fde9147ed74b45 Mon Sep 17 00:00:00 2001 From: JoerivanEngelen Date: Thu, 12 Dec 2024 16:03:53 +0100 Subject: [PATCH 48/50] Add newlines in docstrings for indented lists --- imod/mf6/model_gwf.py | 8 +++++--- imod/mf6/simulation.py | 2 ++ imod/mf6/wel.py | 2 ++ imod/msw/grid_data.py | 1 + imod/msw/sprinkling.py | 2 ++ 5 files changed, 12 insertions(+), 3 deletions(-) diff --git a/imod/mf6/model_gwf.py b/imod/mf6/model_gwf.py index e83ee5e1c..ebc878876 100644 --- a/imod/mf6/model_gwf.py +++ b/imod/mf6/model_gwf.py @@ -218,9 +218,11 @@ def from_imod5_data( times: list[datetime] Time discretization of the simulation. These times are used for the following: - * Times of wells with associated timeseries are resampled to these times - * Start- and end times in the list are used to repeat the stresses - of periodic data (e.g. river stages in iMOD5 for "summer", "winter") + + * Times of wells with associated timeseries are resampled to these times + * Start- and end times in the list are used to repeat the stresses + of periodic data (e.g. river stages in iMOD5 for "summer", "winter") + allocation_options: SimulationAllocationOptions object containing the allocation options per package type. If you want a package to have a different allocation option, diff --git a/imod/mf6/simulation.py b/imod/mf6/simulation.py index 5bc2469f1..4f324f31b 100644 --- a/imod/mf6/simulation.py +++ b/imod/mf6/simulation.py @@ -1361,11 +1361,13 @@ def from_imod5_data( times: list[datetime] time discretization of the model to be imported. This is used for the following: + * Times of wells with associated timeseries are resampled to these times * Start and end times in the list are used to repeat the stresses of periodic data (e.g. river stages in iMOD5 for "summer", "winter") * The simulation is discretized to these times, using :meth:`imod.mf6.Modflow6Simulation.create_time_discretization` + allocation_options: SimulationAllocationOptions, optional object containing the allocation options per package type. If you want a specific package to have a different allocation option, then diff --git a/imod/mf6/wel.py b/imod/mf6/wel.py index d30b58d30..11334b2cf 100644 --- a/imod/mf6/wel.py +++ b/imod/mf6/wel.py @@ -1312,9 +1312,11 @@ def from_imod5_cap_data(cls, imod5_data: Imod5DataDict): (IDF), this grid defines in which layer a groundwater abstraction well should be placed. The ``"artificial_recharge"`` grid contains types which point to the type of abstraction: + * 0: no abstraction * 1: groundwater abstraction * 2: surfacewater abstraction + The ``"artificial_recharge_capacity"`` grid/constant defines the capacity of each groundwater or surfacewater abstraction. This is an ``1:1`` mapping: Each grid cell maps to a separate well. diff --git a/imod/msw/grid_data.py b/imod/msw/grid_data.py index 26237059f..d4c387bd3 100644 --- a/imod/msw/grid_data.py +++ b/imod/msw/grid_data.py @@ -133,6 +133,7 @@ def from_imod5_data( function. Method does the following things: + - Computes rural area from the wetted area and urban area grids. - Sets a second subunit for urban landuse (== 18). - Rootzone depth is converted from cm's to m's. diff --git a/imod/msw/sprinkling.py b/imod/msw/sprinkling.py index 11fb3744b..b5fa8cc69 100644 --- a/imod/msw/sprinkling.py +++ b/imod/msw/sprinkling.py @@ -180,9 +180,11 @@ def from_imod5_data(cls, imod5_data: Imod5DataDict) -> "Sprinkling": (IDF), this grid defines in which layer a groundwater abstraction well should be placed. The ``"artificial_recharge"`` grid contains types which point to the type of abstraction: + * 0: no abstraction * 1: groundwater abstraction * 2: surfacewater abstraction + The ``"artificial_recharge_capacity"`` grid/constant defines the capacity of each groundwater or surfacewater abstraction. This is an ``1:1`` mapping: Each grid cell maps to a separate well. From 439aba67b46f26e59d2bbe2c8de9701ce9b9e1ea Mon Sep 17 00:00:00 2001 From: JoerivanEngelen Date: Thu, 12 Dec 2024 16:04:24 +0100 Subject: [PATCH 49/50] Import Simulation option dataclasses for public API --- imod/prepare/__init__.py | 2 ++ 1 file changed, 2 insertions(+) diff --git a/imod/prepare/__init__.py b/imod/prepare/__init__.py index 03342c87e..3e297c448 100644 --- a/imod/prepare/__init__.py +++ b/imod/prepare/__init__.py @@ -52,6 +52,8 @@ distribute_drn_conductance, distribute_ghb_conductance, distribute_riv_conductance, + SimulationAllocationOptions, + SimulationDistributingOptions, ) from imod.prepare.voxelize import Voxelizer from imod.prepare.wells import assign_wells From 298ad707cddcf939f5d179a1a40f18e8a8bd11d3 Mon Sep 17 00:00:00 2001 From: JoerivanEngelen Date: Thu, 12 Dec 2024 16:05:14 +0100 Subject: [PATCH 50/50] Format --- imod/prepare/__init__.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/imod/prepare/__init__.py b/imod/prepare/__init__.py index 3e297c448..41e16f46d 100644 --- a/imod/prepare/__init__.py +++ b/imod/prepare/__init__.py @@ -43,6 +43,8 @@ from imod.prepare.topsystem import ( ALLOCATION_OPTION, DISTRIBUTING_OPTION, + SimulationAllocationOptions, + SimulationDistributingOptions, allocate_drn_cells, allocate_ghb_cells, allocate_rch_cells, @@ -52,8 +54,6 @@ distribute_drn_conductance, distribute_ghb_conductance, distribute_riv_conductance, - SimulationAllocationOptions, - SimulationDistributingOptions, ) from imod.prepare.voxelize import Voxelizer from imod.prepare.wells import assign_wells