From 37023a0d3108975c50af4b682081bd6acee94e9a Mon Sep 17 00:00:00 2001 From: JoerivanEngelen Date: Thu, 23 Apr 2026 17:04:54 +0200 Subject: [PATCH 01/19] Fix typo --- imod/formats/prj/prj.py | 2 +- imod/msw/scaling_factors.py | 2 +- imod/msw/utilities/imod5_converter.py | 2 +- imod/tests/fixtures/msw_imod5_cap_fixture.py | 2 +- imod/tests/test_msw/test_model.py | 4 ++-- imod/tests/test_msw/test_scaling_factors.py | 2 +- 6 files changed, 7 insertions(+), 7 deletions(-) diff --git a/imod/formats/prj/prj.py b/imod/formats/prj/prj.py index 56e582347..dd9859a92 100644 --- a/imod/formats/prj/prj.py +++ b/imod/formats/prj/prj.py @@ -98,7 +98,7 @@ "rural_infiltration_capacity", "perched_water_table_level", "soil_moisture_fraction", - "conductivitiy_factor", + "conductivity_factor", "plot_number", "steering_location", "plot_drainage_level", diff --git a/imod/msw/scaling_factors.py b/imod/msw/scaling_factors.py index 85eb37324..e4dd10dae 100644 --- a/imod/msw/scaling_factors.py +++ b/imod/msw/scaling_factors.py @@ -109,7 +109,7 @@ def from_imod5_data(cls, imod5_data: Imod5DataDict) -> "ScalingFactors": cap_data["soil_moisture_fraction"], grid_ones ) data["scale_hydraulic_conductivity"] = concat_imod5( - cap_data["conductivitiy_factor"], grid_ones + cap_data["conductivity_factor"], grid_ones ) data["scale_pressure_head"] = concat_imod5(grid_ones, grid_ones) data["depth_perched_water_table"] = cap_data["perched_water_table_level"] diff --git a/imod/msw/utilities/imod5_converter.py b/imod/msw/utilities/imod5_converter.py index 58fcce04c..c6a3cd32d 100644 --- a/imod/msw/utilities/imod5_converter.py +++ b/imod/msw/utilities/imod5_converter.py @@ -106,7 +106,7 @@ def has_active_scaling_factor(imod5_cap: GridDataDict): variable_inactive_mapping = { "perched_water_table_level": MaskValues.msw_default, "soil_moisture_fraction": 1.0, - "conductivitiy_factor": 1.0, + "conductivity_factor": 1.0, } scaling_factor_inactive = True for var, inactive_value in variable_inactive_mapping.items(): diff --git a/imod/tests/fixtures/msw_imod5_cap_fixture.py b/imod/tests/fixtures/msw_imod5_cap_fixture.py index c57925015..28c53335b 100644 --- a/imod/tests/fixtures/msw_imod5_cap_fixture.py +++ b/imod/tests/fixtures/msw_imod5_cap_fixture.py @@ -222,7 +222,7 @@ def imod5_cap_data() -> GridDataDict: ), **da_kwargs ) - d["conductivitiy_factor"]= xr.DataArray( + d["conductivity_factor"]= xr.DataArray( np.array([ [ [2.5, 2.5, 2.5], diff --git a/imod/tests/test_msw/test_model.py b/imod/tests/test_msw/test_model.py index 80f136466..1601a3f89 100644 --- a/imod/tests/test_msw/test_model.py +++ b/imod/tests/test_msw/test_model.py @@ -320,7 +320,7 @@ def case_no_scaling_factors( 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["conductivity_factor"] = xr.DataArray([1.0], **layer_kwargs) return imod5_cap_data, has_scaling_factor def case_constants( @@ -332,7 +332,7 @@ def case_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["conductivity_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) diff --git a/imod/tests/test_msw/test_scaling_factors.py b/imod/tests/test_msw/test_scaling_factors.py index d2610577e..e61e5b1f9 100644 --- a/imod/tests/test_msw/test_scaling_factors.py +++ b/imod/tests/test_msw/test_scaling_factors.py @@ -139,7 +139,7 @@ def test_from_imod5_data(fixed_format_parser): scale_rural = scale.sel(subunit=0, drop=True) imod5_data["cap"]["boundary"] = ones_like(scale_rural) imod5_data["cap"]["soil_moisture_fraction"] = scale_rural - imod5_data["cap"]["conductivitiy_factor"] = scale_rural + imod5_data["cap"]["conductivity_factor"] = scale_rural imod5_data["cap"]["perched_water_table_level"] = depth_perched_water_table scaling_factors = ScalingFactors.from_imod5_data(imod5_data) From 0f5fc4bedcf5ef083ce0e3fab68f9832bf1a29ec Mon Sep 17 00:00:00 2001 From: JoerivanEngelen Date: Thu, 23 Apr 2026 17:57:38 +0200 Subject: [PATCH 02/19] Regrid iMOD5 cap data to MODFLOW6 dis Co-authored-by: Copilot --- imod/msw/model.py | 28 +++++++----- imod/msw/regrid/regrid_schemes.py | 65 +++++++++++++++++++++++++++ imod/msw/utilities/imod5_converter.py | 35 ++++++++++++++- 3 files changed, 117 insertions(+), 11 deletions(-) diff --git a/imod/msw/model.py b/imod/msw/model.py index eaf81fdbf..01a432e68 100644 --- a/imod/msw/model.py +++ b/imod/msw/model.py @@ -13,6 +13,7 @@ from imod.common.constants import MaskValues from imod.common.utilities.clip import clip_by_grid from imod.common.utilities.partitioninfo import create_partition_info +from imod.common.utilities.regrid import _regrid_package_data from imod.common.utilities.value_filters import enforce_scalar from imod.common.utilities.version import prepend_content_with_version_info from imod.mf6.dis import StructuredDiscretization @@ -38,11 +39,12 @@ from imod.msw.output_control import TimeOutputControl from imod.msw.pkgbase import MetaSwapPackage from imod.msw.ponding import Ponding +from imod.msw.regrid.regrid_schemes import CapDataRegridMethod 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.imod5_converter import has_active_scaling_factor, regrid_imod5_data from imod.msw.utilities.mask import mask_and_broadcast_cap_data from imod.msw.utilities.parse import read_para_sim from imod.msw.vegetation import AnnualCropFactors @@ -605,7 +607,8 @@ def from_imod5_data( imod5_data: Imod5DataDict, target_dis: StructuredDiscretization, times: list[datetime], - ): + regridder_types: Optional[CapDataRegridMethod] = None, + ) -> "MetaSwapModel": """ Construct a MetaSWAP model from iMOD5 data in the CAP package, loaded with the :func:`imod.formats.prj.open_projectfile_data` function. @@ -616,29 +619,34 @@ def from_imod5_data( 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. + Target discretization, iMOD5 CAP data will be regridded to this + 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. - + regridder_types: Optional[CapDataRegridMethod], default None + Custom regrid method for CAP data. + Returns ------- MetaSwapModel MetaSWAP model imported from imod5 data. """ + # Path and settings management 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 - imod5_cap_no_layer = drop_layer_dim_cap_data(imod5_data) + # Regrid iMOD5 CAP data to target discretization. + imod5_regridded = regrid_imod5_data(imod5_data, target_dis, regridder_types) + # Setup model model = cls(unsa_svat_path, parasim_settings) model["grid"], msw_active = GridData.from_imod5_data( - imod5_cap_no_layer, target_dis + imod5_regridded, target_dis ) cap_data_masked = mask_and_broadcast_cap_data( - imod5_cap_no_layer["cap"], msw_active + imod5_regridded["cap"], msw_active ) imod5_masked: Imod5DataDict = { "cap": cap_data_masked, @@ -650,7 +658,7 @@ 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) - if has_active_scaling_factor(imod5_cap_no_layer["cap"]): + if has_active_scaling_factor(imod5_masked["cap"]): model["scaling_factor"] = ScalingFactors.from_imod5_data(imod5_masked) area = model["grid"]["area"].isel(subunit=0, drop=True) model["idf_mapping"] = IdfMapping(area, MaskValues.msw_default) diff --git a/imod/msw/regrid/regrid_schemes.py b/imod/msw/regrid/regrid_schemes.py index 9caedc6dd..cf82018b1 100644 --- a/imod/msw/regrid/regrid_schemes.py +++ b/imod/msw/regrid/regrid_schemes.py @@ -10,6 +10,71 @@ ) +@dataclass(config=_CONFIG) +class CapDataRegridMethod(DataclassType): + """ + Object containing regridder methods for the iMOD5 CAP data. This can be + provided to the ``regrid_like`` method to regrid with custom settings. This + is used when constructing a MetaSWAP model from iMOD5 data. + + Parameters + ---------- + boundary: tuple, default (RegridderType.OVERLAP, "mode") + landuse: tuple, default (RegridderType.OVERLAP, "mode") + rootzone_thickness: tuple, default (RegridderType.OVERLAP, "mean") + soil_physical_unit: tuple, default (RegridderType.OVERLAP, "mode") + meteo_station_number: tuple, default (RegridderType.OVERLAP, "mode") + surface_elevation: tuple, default (RegridderType.OVERLAP, "mean") + artificial_recharge: tuple, default (RegridderType.OVERLAP, "mean") + artificial_recharge_layer: tuple, default (RegridderType.OVERLAP, "mode") + artificial_recharge_capacity: tuple, default (RegridderType.OVERLAP, "mean") + wetted_area: tuple, default (RegridderType.RELATIVEOVERLAP, "conductance") + urban_area: tuple, default (RegridderType.RELATIVEOVERLAP, "conductance") + urban_ponding_depth: tuple, default (RegridderType.OVERLAP, "mean") + rural_ponding_depth: tuple, default (RegridderType.OVERLAP, "mean") + urban_runoff_resistance: tuple, default (RegridderType.OVERLAP, "mean") + rural_runoff_resistance: tuple, default (RegridderType.OVERLAP, "mean") + urban_runon_resistance: tuple, default (RegridderType.OVERLAP, "mean") + rural_runon_resistance: tuple, default (RegridderType.OVERLAP, "mean") + urban_infiltration_capacity: tuple, default (RegridderType.OVERLAP, "mean") + rural_infiltration_capacity: tuple, default (RegridderType.OVERLAP, "mean") + perched_water_table_level: tuple, default (RegridderType.OVERLAP, "mean") + soil_moisture_fraction: tuple, default (RegridderType.OVERLAP, "mean") + conductivity_factor: tuple, default (RegridderType.OVERLAP, "mean") + plot_number: tuple, default (RegridderType.OVERLAP, "mode") + steering_location: tuple, default (RegridderType.OVERLAP, "mode") + plot_drainage_level: tuple, default (RegridderType.OVERLAP, "mean") + plot_drainage_resistance: tuple, default (RegridderType.OVERLAP, "mean") + """ + + boundary: RegridVarType = (RegridderType.OVERLAP, "mode") + landuse: RegridVarType = (RegridderType.OVERLAP, "mode") + rootzone_thickness: RegridVarType = (RegridderType.OVERLAP, "mean") + soil_physical_unit: RegridVarType = (RegridderType.OVERLAP, "mode") + meteo_station_number: RegridVarType = (RegridderType.OVERLAP, "mode") + surface_elevation: RegridVarType = (RegridderType.OVERLAP, "mean") + artificial_recharge: RegridVarType = (RegridderType.OVERLAP, "mean") + artificial_recharge_layer: RegridVarType = (RegridderType.OVERLAP, "mode") + artificial_recharge_capacity: RegridVarType = (RegridderType.OVERLAP, "mean") + wetted_area: RegridVarType = (RegridderType.RELATIVEOVERLAP, "conductance") + urban_area: RegridVarType = (RegridderType.RELATIVEOVERLAP, "conductance") + urban_ponding_depth: RegridVarType = (RegridderType.OVERLAP, "mean") + rural_ponding_depth: RegridVarType = (RegridderType.OVERLAP, "mean") + urban_runoff_resistance: RegridVarType = (RegridderType.OVERLAP, "mean") + rural_runoff_resistance: RegridVarType = (RegridderType.OVERLAP, "mean") + urban_runon_resistance: RegridVarType = (RegridderType.OVERLAP, "mean") + rural_runon_resistance: RegridVarType = (RegridderType.OVERLAP, "mean") + urban_infiltration_capacity: RegridVarType = (RegridderType.OVERLAP, "mean") + rural_infiltration_capacity: RegridVarType = (RegridderType.OVERLAP, "mean") + perched_water_table_level: RegridVarType = (RegridderType.OVERLAP, "mean") + soil_moisture_fraction: RegridVarType = (RegridderType.OVERLAP, "mean") + conductivity_factor: RegridVarType = (RegridderType.OVERLAP, "mean") + plot_number: RegridVarType = (RegridderType.OVERLAP, "mode") + steering_location: RegridVarType = (RegridderType.OVERLAP, "mode") + plot_drainage_level: RegridVarType = (RegridderType.OVERLAP, "mean") + plot_drainage_resistance: RegridVarType = (RegridderType.OVERLAP, "mean") + + @dataclass(config=_CONFIG) class SprinklingRegridMethod(DataclassType): """ diff --git a/imod/msw/utilities/imod5_converter.py b/imod/msw/utilities/imod5_converter.py index c6a3cd32d..bca107889 100644 --- a/imod/msw/utilities/imod5_converter.py +++ b/imod/msw/utilities/imod5_converter.py @@ -1,12 +1,18 @@ +from typing import Optional + from xarray.core.utils import is_scalar from imod.common.constants import MaskValues +from imod.common.utilities.regrid import _regrid_package_data from imod.logging import LogLevel, logger from imod.mf6 import StructuredDiscretization +from imod.msw.regrid.regrid_schemes import CapDataRegridMethod from imod.msw.utilities.common import concat_imod5 from imod.msw.utilities.mask import MetaSwapActive -from imod.typing import GridDataArray, GridDataDict +from imod.typing import GridDataArray, GridDataDict, Imod5DataDict from imod.typing.grid import ones_like +from imod.util.dims import drop_layer_dim_cap_data +from imod.util.regrid import RegridderWeightsCache from imod.util.spatial import get_cell_area @@ -116,3 +122,30 @@ def has_active_scaling_factor(imod5_cap: GridDataDict): ) return not scaling_factor_inactive + +def regrid_imod5_data( + imod5_data: Imod5DataDict, target_dis: StructuredDiscretization, regridder_types: Optional[CapDataRegridMethod] = None + ) -> Imod5DataDict: + """ + Regrid iMOD5 CAP data to consistent grid. This is necessary to be able to + use iMOD5 data in MetaSWAP, as the grid of the iMOD5 CAP data is not + necessarily the same as the grid of the target MODFLOW6 discretization. The + regridding process ensures consistency between the iMOD5 CAP data and the + target MODFLOW6 grid. + """ + # Drop layer coords + imod5_cap_no_layer = drop_layer_dim_cap_data(imod5_data) + target_grid = target_dis.dataset["idomain"].isel(layer=0, drop=True) + # Regrid the input data + if regridder_types is not None: + regridder_types = CapDataRegridMethod() + regrid_cache = RegridderWeightsCache() + cap_data_regridded = _regrid_package_data( + imod5_cap_no_layer, target_grid, regridder_types, regrid_cache + ) + extra_paths = imod5_data["extra"]["paths"] + imod5_regridded: Imod5DataDict = { + "cap": cap_data_regridded, + "extra": {"paths": extra_paths}, + } + return imod5_regridded \ No newline at end of file From ead10e54bf3d4a7e2741e58fb0e9b1bed80039e3 Mon Sep 17 00:00:00 2001 From: JoerivanEngelen Date: Tue, 28 Apr 2026 11:22:56 +0200 Subject: [PATCH 03/19] format --- imod/msw/model.py | 9 +++++---- imod/msw/utilities/imod5_converter.py | 9 ++++++--- 2 files changed, 11 insertions(+), 7 deletions(-) diff --git a/imod/msw/model.py b/imod/msw/model.py index 01a432e68..52099698d 100644 --- a/imod/msw/model.py +++ b/imod/msw/model.py @@ -13,7 +13,6 @@ from imod.common.constants import MaskValues from imod.common.utilities.clip import clip_by_grid from imod.common.utilities.partitioninfo import create_partition_info -from imod.common.utilities.regrid import _regrid_package_data from imod.common.utilities.value_filters import enforce_scalar from imod.common.utilities.version import prepend_content_with_version_info from imod.mf6.dis import StructuredDiscretization @@ -44,12 +43,14 @@ 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, regrid_imod5_data +from imod.msw.utilities.imod5_converter import ( + has_active_scaling_factor, + regrid_imod5_data, +) from imod.msw.utilities.mask import mask_and_broadcast_cap_data from imod.msw.utilities.parse import read_para_sim from imod.msw.vegetation import AnnualCropFactors from imod.typing import GridDataArray, Imod5DataDict -from imod.util.dims import drop_layer_dim_cap_data from imod.util.regrid import RegridderWeightsCache from imod.util.time import to_datetime_internal @@ -627,7 +628,7 @@ def from_imod5_data( Is also used to infer the starttime of the simulation. regridder_types: Optional[CapDataRegridMethod], default None Custom regrid method for CAP data. - + Returns ------- MetaSwapModel diff --git a/imod/msw/utilities/imod5_converter.py b/imod/msw/utilities/imod5_converter.py index bca107889..5e6ad8c09 100644 --- a/imod/msw/utilities/imod5_converter.py +++ b/imod/msw/utilities/imod5_converter.py @@ -123,9 +123,12 @@ def has_active_scaling_factor(imod5_cap: GridDataDict): return not scaling_factor_inactive + def regrid_imod5_data( - imod5_data: Imod5DataDict, target_dis: StructuredDiscretization, regridder_types: Optional[CapDataRegridMethod] = None - ) -> Imod5DataDict: + imod5_data: Imod5DataDict, + target_dis: StructuredDiscretization, + regridder_types: Optional[CapDataRegridMethod] = None, +) -> Imod5DataDict: """ Regrid iMOD5 CAP data to consistent grid. This is necessary to be able to use iMOD5 data in MetaSWAP, as the grid of the iMOD5 CAP data is not @@ -148,4 +151,4 @@ def regrid_imod5_data( "cap": cap_data_regridded, "extra": {"paths": extra_paths}, } - return imod5_regridded \ No newline at end of file + return imod5_regridded From b209d003ebb16968a263ee71f403b483214454f9 Mon Sep 17 00:00:00 2001 From: JoerivanEngelen Date: Tue, 28 Apr 2026 15:59:49 +0200 Subject: [PATCH 04/19] Squeeze if scalar to drop dimensions, instead of converting to numpy. The preservation of the dtype has the advantage that it allows treating it as a DataArray for masking. --- imod/common/utilities/regrid.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/imod/common/utilities/regrid.py b/imod/common/utilities/regrid.py index 4fbea0a31..e0fcb64c8 100644 --- a/imod/common/utilities/regrid.py +++ b/imod/common/utilities/regrid.py @@ -77,7 +77,7 @@ def _regrid_array( # the dataarray might be a scalar. If it is, then it does not need regridding. if scalar_da: - return da.values[()] + return da.squeeze(drop=True) if isinstance(da, xr.DataArray): coords = da.coords From 818ace688cd39afb7ed7a4680a782dd0a3961dfd Mon Sep 17 00:00:00 2001 From: JoerivanEngelen Date: Tue, 28 Apr 2026 16:00:41 +0200 Subject: [PATCH 05/19] Carefully assess active_scaling_factor with the right data. --- imod/msw/model.py | 6 +++++- 1 file changed, 5 insertions(+), 1 deletion(-) diff --git a/imod/msw/model.py b/imod/msw/model.py index 52099698d..836bc7cf2 100644 --- a/imod/msw/model.py +++ b/imod/msw/model.py @@ -641,6 +641,10 @@ def from_imod5_data( unsa_svat_path = cast(str, parasim_settings["unsa_svat_path"]) # Regrid iMOD5 CAP data to target discretization. imod5_regridded = regrid_imod5_data(imod5_data, target_dis, regridder_types) + # Test with regridded data instead of masked, as masking broadcasts + # scalars to grids, which causes the is_scalar check in + # has_active_scaling_factor to always return False. + active_scaling_factor = has_active_scaling_factor(imod5_regridded["cap"]) # Setup model model = cls(unsa_svat_path, parasim_settings) model["grid"], msw_active = GridData.from_imod5_data( @@ -659,7 +663,7 @@ 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) - if has_active_scaling_factor(imod5_masked["cap"]): + if active_scaling_factor: model["scaling_factor"] = ScalingFactors.from_imod5_data(imod5_masked) area = model["grid"]["area"].isel(subunit=0, drop=True) model["idf_mapping"] = IdfMapping(area, MaskValues.msw_default) From 6dbebcd3457a82c2aaae3e3e027741a184540a1e Mon Sep 17 00:00:00 2001 From: JoerivanEngelen Date: Tue, 28 Apr 2026 16:01:14 +0200 Subject: [PATCH 06/19] Fix is None check and unpack iMOD5DataDict --- imod/msw/utilities/imod5_converter.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/imod/msw/utilities/imod5_converter.py b/imod/msw/utilities/imod5_converter.py index 5e6ad8c09..5c9bac9a7 100644 --- a/imod/msw/utilities/imod5_converter.py +++ b/imod/msw/utilities/imod5_converter.py @@ -140,11 +140,11 @@ def regrid_imod5_data( imod5_cap_no_layer = drop_layer_dim_cap_data(imod5_data) target_grid = target_dis.dataset["idomain"].isel(layer=0, drop=True) # Regrid the input data - if regridder_types is not None: + if regridder_types is None: regridder_types = CapDataRegridMethod() regrid_cache = RegridderWeightsCache() cap_data_regridded = _regrid_package_data( - imod5_cap_no_layer, target_grid, regridder_types, regrid_cache + imod5_cap_no_layer["cap"], target_grid, regridder_types, regrid_cache ) extra_paths = imod5_data["extra"]["paths"] imod5_regridded: Imod5DataDict = { From db3e71c4ccee5c39992e241f299f04f824d20bbf Mon Sep 17 00:00:00 2001 From: JoerivanEngelen Date: Tue, 28 Apr 2026 16:01:43 +0200 Subject: [PATCH 07/19] Add cases for a target_dis with different grid. Co-authored-by: Copilot --- imod/tests/test_msw/test_model.py | 38 ++++++++++++++++++++++++++++--- 1 file changed, 35 insertions(+), 3 deletions(-) diff --git a/imod/tests/test_msw/test_model.py b/imod/tests/test_msw/test_model.py index 1601a3f89..c7f40fd1d 100644 --- a/imod/tests/test_msw/test_model.py +++ b/imod/tests/test_msw/test_model.py @@ -1,5 +1,6 @@ from copy import copy from pathlib import Path +from typing import cast import pytest import xarray as xr @@ -12,9 +13,11 @@ from imod.msw.meteo_mapping import MeteoMapping from imod.msw.model import DEFAULT_SETTINGS from imod.msw.utilities.parse import read_para_sim +from imod.tests.fixtures.msw_model_fixture import coupled_mf6_model from imod.typing import GridDataArray, Imod5DataDict from imod.typing.grid import zeros_like from imod.util.regrid import RegridderWeightsCache +from imod.util.spatial import empty_2d def test_msw_model_write(msw_model, coupled_mf6_model, coupled_mf6wel, tmp_path): @@ -344,19 +347,46 @@ def case_constants( return imod5_cap_data, has_scaling_factor +class TargetDisCases: + def case_dis_same_grid( + self, coupled_mf6_model: mf6.Modflow6Simulation + ) -> tuple[mf6.StructuredDiscretization :, xr.DataArray]: + """ + Modflow6 discretization on the same grid as the CAP data, thus no regridding needed. + """ + dis_pkg = cast(mf6.StructuredDiscretization, coupled_mf6_model["GWF_1"]["dis"]) + times = coupled_mf6_model["time_discretization"].dataset.coords["time"] + return dis_pkg, times + + def case_dis_different_grid( + self, coupled_mf6_model: mf6.Modflow6Simulation + ) -> tuple[mf6.StructuredDiscretization :, xr.DataArray]: + """ + Modflow6 discretization on different grid as the CAP data, thus regridding needed. + """ + dis_pkg = cast(mf6.StructuredDiscretization, coupled_mf6_model["GWF_1"]["dis"]) + target_grid = empty_2d( + dx=0.5, xmin=0.75, xmax=3.25, dy=-0.5, ymin=0.75, ymax=3.25 + ) + times = coupled_mf6_model["time_discretization"].dataset.coords["time"] + regridder_types = RegridderWeightsCache() + dis_pkg_regridded = dis_pkg.regrid_like(target_grid, regridder_types) + return dis_pkg_regridded, times + + @pytest.mark.unittest_jit @parametrize_with_cases("imod5_data, has_scaling_factor", cases=Imod5DataCases) +@parametrize_with_cases("dis_pkg, times", cases=TargetDisCases) def test_import_from_imod5( imod5_data: Imod5DataDict, has_scaling_factor: bool, meteo_grids: tuple[GridDataArray], - coupled_mf6_model: mf6.Modflow6Simulation, + dis_pkg: mf6.StructuredDiscretization, + times: xr.DataArray, tmp_path: Path, ): # Arrange 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_data, dis_pkg, times) # Assert @@ -384,6 +414,8 @@ def test_import_from_imod5( assert "time" in model["time_oc"].dataset.sizes.keys() assert len(model["meteo_grid"].dataset.dims) == 0 assert ("scaling_factor" in model.keys()) == has_scaling_factor + assert dis_pkg.dataset["x"].equals(model["grid"].dataset["x"]) + assert dis_pkg.dataset["y"].equals(model["grid"].dataset["y"]) @pytest.mark.unittest_jit From 12acec1dec086cf0ef2f7c2edcd437313953a875 Mon Sep 17 00:00:00 2001 From: JoerivanEngelen Date: Tue, 28 Apr 2026 16:06:51 +0200 Subject: [PATCH 08/19] format --- imod/tests/test_msw/test_model.py | 1 - 1 file changed, 1 deletion(-) diff --git a/imod/tests/test_msw/test_model.py b/imod/tests/test_msw/test_model.py index c7f40fd1d..8e584bf2e 100644 --- a/imod/tests/test_msw/test_model.py +++ b/imod/tests/test_msw/test_model.py @@ -13,7 +13,6 @@ from imod.msw.meteo_mapping import MeteoMapping from imod.msw.model import DEFAULT_SETTINGS from imod.msw.utilities.parse import read_para_sim -from imod.tests.fixtures.msw_model_fixture import coupled_mf6_model from imod.typing import GridDataArray, Imod5DataDict from imod.typing.grid import zeros_like from imod.util.regrid import RegridderWeightsCache From b0033569dfe6dc925c8148911facd26d30e2ccf6 Mon Sep 17 00:00:00 2001 From: JoerivanEngelen Date: Tue, 28 Apr 2026 16:11:54 +0200 Subject: [PATCH 09/19] Update changelog --- docs/api/changelog.rst | 2 ++ 1 file changed, 2 insertions(+) diff --git a/docs/api/changelog.rst b/docs/api/changelog.rst index a69de959a..37f424888 100644 --- a/docs/api/changelog.rst +++ b/docs/api/changelog.rst @@ -64,6 +64,8 @@ Fixed - :meth:`imod.mf6.NodePropertyFlow.regrid_like` now regrids ``k33`` using the correct method, namely ``mean`` instead of ``harmonic_mean``. As this is the appropriate method for horizontal regridding of ``k33``. +- :meth:`imod.msw.MetaSwapModel.from_imod5_data` now regrids the iMOD5 data to + the MODFLOW6 target discretization. Changed ~~~~~~~ From b29dc2adb632da51fb9b943cd904d00592a16827 Mon Sep 17 00:00:00 2001 From: JoerivanEngelen Date: Wed, 29 Apr 2026 11:52:18 +0200 Subject: [PATCH 10/19] Make CapDataRegridMethod not Optional anymore, change to a default value. Co-authored-by: Copilot --- imod/msw/model.py | 4 ++-- imod/msw/utilities/imod5_converter.py | 6 +----- 2 files changed, 3 insertions(+), 7 deletions(-) diff --git a/imod/msw/model.py b/imod/msw/model.py index 836bc7cf2..7d71c683c 100644 --- a/imod/msw/model.py +++ b/imod/msw/model.py @@ -608,7 +608,7 @@ def from_imod5_data( imod5_data: Imod5DataDict, target_dis: StructuredDiscretization, times: list[datetime], - regridder_types: Optional[CapDataRegridMethod] = None, + regridder_types: CapDataRegridMethod = CapDataRegridMethod(), ) -> "MetaSwapModel": """ Construct a MetaSWAP model from iMOD5 data in the CAP package, loaded @@ -626,7 +626,7 @@ def from_imod5_data( 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. - regridder_types: Optional[CapDataRegridMethod], default None + regridder_types: CapDataRegridMethod, default CapDataRegridMethod() Custom regrid method for CAP data. Returns diff --git a/imod/msw/utilities/imod5_converter.py b/imod/msw/utilities/imod5_converter.py index 5c9bac9a7..e7d1c139e 100644 --- a/imod/msw/utilities/imod5_converter.py +++ b/imod/msw/utilities/imod5_converter.py @@ -1,5 +1,3 @@ -from typing import Optional - from xarray.core.utils import is_scalar from imod.common.constants import MaskValues @@ -127,7 +125,7 @@ def has_active_scaling_factor(imod5_cap: GridDataDict): def regrid_imod5_data( imod5_data: Imod5DataDict, target_dis: StructuredDiscretization, - regridder_types: Optional[CapDataRegridMethod] = None, + regridder_types: CapDataRegridMethod, ) -> Imod5DataDict: """ Regrid iMOD5 CAP data to consistent grid. This is necessary to be able to @@ -140,8 +138,6 @@ def regrid_imod5_data( imod5_cap_no_layer = drop_layer_dim_cap_data(imod5_data) target_grid = target_dis.dataset["idomain"].isel(layer=0, drop=True) # Regrid the input data - if regridder_types is None: - regridder_types = CapDataRegridMethod() regrid_cache = RegridderWeightsCache() cap_data_regridded = _regrid_package_data( imod5_cap_no_layer["cap"], target_grid, regridder_types, regrid_cache From 4425bcda45c32fc52536ddb3f98d1f6ad4c358af Mon Sep 17 00:00:00 2001 From: JoerivanEngelen Date: Wed, 29 Apr 2026 13:16:27 +0200 Subject: [PATCH 11/19] Add regrid_cache argument to ``MetaSwapModel.from_imod5_data``. Co-authored-by: Copilot --- imod/msw/model.py | 8 +++++++- imod/msw/utilities/imod5_converter.py | 6 +++--- 2 files changed, 10 insertions(+), 4 deletions(-) diff --git a/imod/msw/model.py b/imod/msw/model.py index 7d71c683c..5d7b4e495 100644 --- a/imod/msw/model.py +++ b/imod/msw/model.py @@ -609,6 +609,7 @@ def from_imod5_data( target_dis: StructuredDiscretization, times: list[datetime], regridder_types: CapDataRegridMethod = CapDataRegridMethod(), + regrid_cache: RegridderWeightsCache = RegridderWeightsCache(), ) -> "MetaSwapModel": """ Construct a MetaSWAP model from iMOD5 data in the CAP package, loaded @@ -628,6 +629,9 @@ def from_imod5_data( Is also used to infer the starttime of the simulation. regridder_types: CapDataRegridMethod, default CapDataRegridMethod() Custom regrid method for CAP data. + regrid_cache: RegridderWeightsCache, default RegridderWeightsCache() + Cache for regridder weights, can be used to speed up regridding if the + same regridders are used multiple times. Returns ------- @@ -640,7 +644,9 @@ def from_imod5_data( parasim_settings = read_para_sim(path_to_parasim) unsa_svat_path = cast(str, parasim_settings["unsa_svat_path"]) # Regrid iMOD5 CAP data to target discretization. - imod5_regridded = regrid_imod5_data(imod5_data, target_dis, regridder_types) + imod5_regridded = regrid_imod5_data( + imod5_data, target_dis, regridder_types, regrid_cache + ) # Test with regridded data instead of masked, as masking broadcasts # scalars to grids, which causes the is_scalar check in # has_active_scaling_factor to always return False. diff --git a/imod/msw/utilities/imod5_converter.py b/imod/msw/utilities/imod5_converter.py index e7d1c139e..85b609b7b 100644 --- a/imod/msw/utilities/imod5_converter.py +++ b/imod/msw/utilities/imod5_converter.py @@ -1,10 +1,10 @@ from xarray.core.utils import is_scalar from imod.common.constants import MaskValues +from imod.common.utilities.dataclass_type import DataclassType from imod.common.utilities.regrid import _regrid_package_data from imod.logging import LogLevel, logger from imod.mf6 import StructuredDiscretization -from imod.msw.regrid.regrid_schemes import CapDataRegridMethod from imod.msw.utilities.common import concat_imod5 from imod.msw.utilities.mask import MetaSwapActive from imod.typing import GridDataArray, GridDataDict, Imod5DataDict @@ -125,7 +125,8 @@ def has_active_scaling_factor(imod5_cap: GridDataDict): def regrid_imod5_data( imod5_data: Imod5DataDict, target_dis: StructuredDiscretization, - regridder_types: CapDataRegridMethod, + regridder_types: DataclassType, + regrid_cache: RegridderWeightsCache, ) -> Imod5DataDict: """ Regrid iMOD5 CAP data to consistent grid. This is necessary to be able to @@ -138,7 +139,6 @@ def regrid_imod5_data( imod5_cap_no_layer = drop_layer_dim_cap_data(imod5_data) target_grid = target_dis.dataset["idomain"].isel(layer=0, drop=True) # Regrid the input data - regrid_cache = RegridderWeightsCache() cap_data_regridded = _regrid_package_data( imod5_cap_no_layer["cap"], target_grid, regridder_types, regrid_cache ) From 56b5a4d4e569e4a77ee9a3ad5bd6da672c1590cf Mon Sep 17 00:00:00 2001 From: JoerivanEngelen Date: Wed, 29 Apr 2026 13:42:00 +0200 Subject: [PATCH 12/19] Regrid imod5 cap data also for rch.from_imod5_data Co-authored-by: Copilot --- imod/common/utilities/regrid.py | 33 +++++++++++++++++++++++++- imod/mf6/rch.py | 13 +++++++--- imod/mf6/regrid/regrid_schemes.py | 13 ++++++++++ imod/msw/model.py | 4 ++-- imod/msw/utilities/imod5_converter.py | 34 +-------------------------- 5 files changed, 58 insertions(+), 39 deletions(-) diff --git a/imod/common/utilities/regrid.py b/imod/common/utilities/regrid.py index e0fcb64c8..70ee6fc9b 100644 --- a/imod/common/utilities/regrid.py +++ b/imod/common/utilities/regrid.py @@ -20,13 +20,15 @@ from imod.common.utilities.dataclass_type import DataclassType, EmptyRegridMethod from imod.common.utilities.dtype import is_integer from imod.common.utilities.value_filters import is_valid +from imod.mf6 import StructuredDiscretization +from imod.typing import Imod5DataDict from imod.typing.grid import ( GridDataArray, GridDataset, is_unstructured, ones_like, ) -from imod.util.dims import enforced_dim_order +from imod.util.dims import drop_layer_dim_cap_data, enforced_dim_order from imod.util.regrid import ( RegridderType, RegridderWeightsCache, @@ -459,3 +461,32 @@ def _get_regridding_domain( ).astype(int) return new_idomain + + +def regrid_imod5_cap_data( + imod5_data: Imod5DataDict, + target_dis: StructuredDiscretization, + regridder_types: DataclassType, + regrid_cache: RegridderWeightsCache, +) -> Imod5DataDict: + """ + Regrid iMOD5 CAP data to consistent grid. This is necessary to be able to + use iMOD5 data in MetaSWAP, as the grid of the iMOD5 CAP data is not + necessarily the same as the grid of the target MODFLOW6 discretization. The + regridding process ensures consistency between the iMOD5 CAP data and the + target MODFLOW6 grid. Used in both ``imod.msw.MetaSwapModel.from_imod5_data`` + and ``imod.mf6.Recharge.from_imod5_cap_data``. + """ + # Drop layer coords + imod5_cap_no_layer = drop_layer_dim_cap_data(imod5_data) + target_grid = target_dis.dataset["idomain"].isel(layer=0, drop=True) + # Regrid the input data + cap_data_regridded = _regrid_package_data( + imod5_cap_no_layer["cap"], target_grid, regridder_types, regrid_cache + ) + extra_paths = imod5_data["extra"]["paths"] + imod5_regridded: Imod5DataDict = { + "cap": cap_data_regridded, + "extra": {"paths": extra_paths}, + } + return imod5_regridded diff --git a/imod/mf6/rch.py b/imod/mf6/rch.py index 642c7b637..1933f480e 100644 --- a/imod/mf6/rch.py +++ b/imod/mf6/rch.py @@ -6,10 +6,14 @@ from imod.common.interfaces.iregridpackage import IRegridPackage from imod.common.utilities.dataclass_type import DataclassType +from imod.common.utilities.regrid import regrid_imod5_cap_data from imod.logging import init_log_decorator from imod.mf6.aggregate.aggregate_schemes import RechargeAggregationMethod from imod.mf6.dis import StructuredDiscretization, VerticesDiscretization -from imod.mf6.regrid.regrid_schemes import RechargeRegridMethod +from imod.mf6.regrid.regrid_schemes import ( + CapDataRechargeRegridMethod, + RechargeRegridMethod, +) from imod.mf6.topsystem import TopSystemBoundaryCondition from imod.mf6.utilities.imod5_converter import ( convert_unit_rch_rate, @@ -39,7 +43,6 @@ enforce_dim_order, is_planar_grid, ) -from imod.util.dims import drop_layer_dim_cap_data from imod.util.regrid import RegridderWeightsCache @@ -286,6 +289,8 @@ def from_imod5_cap_data( cls, imod5_data: Imod5DataDict, target_dis: StructuredDiscretization, + regridder_types: CapDataRechargeRegridMethod = CapDataRechargeRegridMethod(), + regrid_cache: RegridderWeightsCache = RegridderWeightsCache(), ) -> "Recharge": """ Construct an rch-package from iMOD5 data in the CAP package, loaded with @@ -293,7 +298,9 @@ def from_imod5_cap_data( used to couple MODFLOW6 to MetaSWAP models. Active cells will have a recharge rate of 0.0. """ - cap_data = drop_layer_dim_cap_data(imod5_data)["cap"] + cap_data = regrid_imod5_cap_data( + imod5_data, target_dis, regridder_types, regrid_cache + )["cap"] msw_area = get_cell_area_from_imod5_data(cap_data) msw_active = is_msw_active_cell(target_dis, cap_data, msw_area) diff --git a/imod/mf6/regrid/regrid_schemes.py b/imod/mf6/regrid/regrid_schemes.py index 0af7267af..a852eeda9 100644 --- a/imod/mf6/regrid/regrid_schemes.py +++ b/imod/mf6/regrid/regrid_schemes.py @@ -450,3 +450,16 @@ class StorageCoefficientRegridMethod(DataclassType): convertible: RegridVarType = (RegridderType.OVERLAP, "mode") storage_coefficient: RegridVarType = (RegridderType.OVERLAP, "mean") specific_yield: RegridVarType = (RegridderType.OVERLAP, "mean") + + +@dataclass(config=_CONFIG) +class CapDataRechargeRegridMethod(DataclassType): + """ + Object containing regridder methods for CAP data for the + :class:`imod.mf6.Recharge.from_imod5_data` method. This contains regridder + methods for only the relevant CAP variables for the recharge package. + """ + + boundary: RegridVarType = (RegridderType.OVERLAP, "mode") + wetted_area: RegridVarType = (RegridderType.RELATIVEOVERLAP, "conductance") + urban_area: RegridVarType = (RegridderType.RELATIVEOVERLAP, "conductance") diff --git a/imod/msw/model.py b/imod/msw/model.py index 5d7b4e495..1c604c352 100644 --- a/imod/msw/model.py +++ b/imod/msw/model.py @@ -13,6 +13,7 @@ from imod.common.constants import MaskValues from imod.common.utilities.clip import clip_by_grid from imod.common.utilities.partitioninfo import create_partition_info +from imod.common.utilities.regrid import regrid_imod5_cap_data from imod.common.utilities.value_filters import enforce_scalar from imod.common.utilities.version import prepend_content_with_version_info from imod.mf6.dis import StructuredDiscretization @@ -45,7 +46,6 @@ from imod.msw.utilities.common import find_in_file_list from imod.msw.utilities.imod5_converter import ( has_active_scaling_factor, - regrid_imod5_data, ) from imod.msw.utilities.mask import mask_and_broadcast_cap_data from imod.msw.utilities.parse import read_para_sim @@ -644,7 +644,7 @@ def from_imod5_data( parasim_settings = read_para_sim(path_to_parasim) unsa_svat_path = cast(str, parasim_settings["unsa_svat_path"]) # Regrid iMOD5 CAP data to target discretization. - imod5_regridded = regrid_imod5_data( + imod5_regridded = regrid_imod5_cap_data( imod5_data, target_dis, regridder_types, regrid_cache ) # Test with regridded data instead of masked, as masking broadcasts diff --git a/imod/msw/utilities/imod5_converter.py b/imod/msw/utilities/imod5_converter.py index 85b609b7b..c6a3cd32d 100644 --- a/imod/msw/utilities/imod5_converter.py +++ b/imod/msw/utilities/imod5_converter.py @@ -1,16 +1,12 @@ from xarray.core.utils import is_scalar from imod.common.constants import MaskValues -from imod.common.utilities.dataclass_type import DataclassType -from imod.common.utilities.regrid import _regrid_package_data 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, Imod5DataDict +from imod.typing import GridDataArray, GridDataDict from imod.typing.grid import ones_like -from imod.util.dims import drop_layer_dim_cap_data -from imod.util.regrid import RegridderWeightsCache from imod.util.spatial import get_cell_area @@ -120,31 +116,3 @@ def has_active_scaling_factor(imod5_cap: GridDataDict): ) return not scaling_factor_inactive - - -def regrid_imod5_data( - imod5_data: Imod5DataDict, - target_dis: StructuredDiscretization, - regridder_types: DataclassType, - regrid_cache: RegridderWeightsCache, -) -> Imod5DataDict: - """ - Regrid iMOD5 CAP data to consistent grid. This is necessary to be able to - use iMOD5 data in MetaSWAP, as the grid of the iMOD5 CAP data is not - necessarily the same as the grid of the target MODFLOW6 discretization. The - regridding process ensures consistency between the iMOD5 CAP data and the - target MODFLOW6 grid. - """ - # Drop layer coords - imod5_cap_no_layer = drop_layer_dim_cap_data(imod5_data) - target_grid = target_dis.dataset["idomain"].isel(layer=0, drop=True) - # Regrid the input data - cap_data_regridded = _regrid_package_data( - imod5_cap_no_layer["cap"], target_grid, regridder_types, regrid_cache - ) - extra_paths = imod5_data["extra"]["paths"] - imod5_regridded: Imod5DataDict = { - "cap": cap_data_regridded, - "extra": {"paths": extra_paths}, - } - return imod5_regridded From c752741c3543629b1af8c12197d57cef9af81ec4 Mon Sep 17 00:00:00 2001 From: JoerivanEngelen Date: Wed, 29 Apr 2026 13:50:14 +0200 Subject: [PATCH 13/19] Remove type ignore, which is not necessary anymore apparently --- imod/mf6/riv.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/imod/mf6/riv.py b/imod/mf6/riv.py index 15f82e2e3..3cf85667f 100644 --- a/imod/mf6/riv.py +++ b/imod/mf6/riv.py @@ -141,7 +141,7 @@ def _create_drain_from_leftover_riv_imod5_data( else: drain_leftover_data = infiltration_drn_data - return Drainage(**drain_leftover_data) # type: ignore[arg-type] + return Drainage(**drain_leftover_data) class River(TopSystemBoundaryCondition, IRegridPackage): From 9aa622ef561a85da0ce6170048a14d1f58acde75 Mon Sep 17 00:00:00 2001 From: JoerivanEngelen Date: Wed, 29 Apr 2026 14:28:47 +0200 Subject: [PATCH 14/19] Regrid Cap data for well data as well. Co-authored-by: Copilot --- imod/mf6/model_gwf.py | 8 ++++++-- imod/mf6/regrid/regrid_schemes.py | 14 +++++++++++++- imod/mf6/utilities/imod5_converter.py | 21 +++++++++++++++------ imod/mf6/wel.py | 14 ++++++++++++-- 4 files changed, 46 insertions(+), 11 deletions(-) diff --git a/imod/mf6/model_gwf.py b/imod/mf6/model_gwf.py index 5e7b544af..40cf867b2 100644 --- a/imod/mf6/model_gwf.py +++ b/imod/mf6/model_gwf.py @@ -354,8 +354,12 @@ def from_imod5_data( ) if "cap" in imod5_keys: - result["msw-sprinkling"] = LayeredWell.from_imod5_cap_data(imod5_data) # type: ignore - result["msw-rch"] = Recharge.from_imod5_cap_data(imod5_data, dis_pkg) # type: ignore + result["msw-sprinkling"] = LayeredWell.from_imod5_cap_data( + imod5_data, dis_pkg, regrid_cache=regrid_cache + ) # type: ignore + result["msw-rch"] = Recharge.from_imod5_cap_data( + imod5_data, dis_pkg, regrid_cache=regrid_cache + ) # type: ignore # import ghb's ghb_keys = [key for key in imod5_keys if key[0:3] == "ghb"] diff --git a/imod/mf6/regrid/regrid_schemes.py b/imod/mf6/regrid/regrid_schemes.py index a852eeda9..f2a957316 100644 --- a/imod/mf6/regrid/regrid_schemes.py +++ b/imod/mf6/regrid/regrid_schemes.py @@ -456,10 +456,22 @@ class StorageCoefficientRegridMethod(DataclassType): class CapDataRechargeRegridMethod(DataclassType): """ Object containing regridder methods for CAP data for the - :class:`imod.mf6.Recharge.from_imod5_data` method. This contains regridder + :class:`imod.mf6.Recharge.from_imod5_cap_data` method. This contains regridder methods for only the relevant CAP variables for the recharge package. """ boundary: RegridVarType = (RegridderType.OVERLAP, "mode") wetted_area: RegridVarType = (RegridderType.RELATIVEOVERLAP, "conductance") urban_area: RegridVarType = (RegridderType.RELATIVEOVERLAP, "conductance") + + +@dataclass(config=_CONFIG) +class CapDataWellRegridMethod(DataclassType): + """ + Object containing regridder methods for CAP data for the + :class:`imod.mf6.LayeredWell.from_imod5_cap_data` method. This contains + regridder methods for only the relevant CAP variables for the well package. + """ + + artificial_recharge: RegridVarType = (RegridderType.OVERLAP, "mean") + artificial_recharge_layer: RegridVarType = (RegridderType.OVERLAP, "mode") diff --git a/imod/mf6/utilities/imod5_converter.py b/imod/mf6/utilities/imod5_converter.py index 5044fd79f..b62834e21 100644 --- a/imod/mf6/utilities/imod5_converter.py +++ b/imod/mf6/utilities/imod5_converter.py @@ -5,11 +5,11 @@ import xarray as xr from imod.common.utilities.dataclass_type import DataclassType -from imod.common.utilities.regrid import _regrid_package_data +from imod.common.utilities.regrid import _regrid_package_data, regrid_imod5_cap_data +from imod.mf6.dis import StructuredDiscretization from imod.mf6.package import Package from imod.typing import GridDataArray, GridDataDict, Imod5DataDict from imod.typing.grid import full_like -from imod.util.dims import drop_layer_dim_cap_data from imod.util.regrid import RegridderWeightsCache @@ -83,7 +83,12 @@ def _well_from_imod5_cap_grid_data(cap_data: GridDataDict) -> dict[str, np.ndarr return data -def well_from_imod5_cap_data(imod5_data: Imod5DataDict) -> dict[str, np.ndarray]: +def well_from_imod5_cap_data( + imod5_data: Imod5DataDict, + target_dis: StructuredDiscretization, + regridder_types: DataclassType, + regrid_cache: RegridderWeightsCache, +) -> dict[str, np.ndarray]: """ Abstraction data for sprinkling is defined in iMOD5 either with grids (IDF) or points (IPF) combined with a grid. Depending on the type, the function does @@ -110,12 +115,16 @@ def well_from_imod5_cap_data(imod5_data: Imod5DataDict) -> dict[str, np.ndarray] capacity is already defined in the point data. This is an ``n:1`` mapping: multiple grid cells can map to one well. """ - cap_data = drop_layer_dim_cap_data(imod5_data)["cap"] + cap_data = imod5_data["cap"] + has_ipf_well = isinstance(cap_data["artificial_recharge_layer"], pd.DataFrame) - if isinstance(cap_data["artificial_recharge_layer"], pd.DataFrame): + if has_ipf_well: return _well_from_imod5_cap_point_data(cap_data) else: - return _well_from_imod5_cap_grid_data(cap_data) + cap_data_regridded = regrid_imod5_cap_data( + imod5_data, target_dis, regridder_types, regrid_cache + )["cap"] + return _well_from_imod5_cap_grid_data(cap_data_regridded) def regrid_imod5_pkg_data( diff --git a/imod/mf6/wel.py b/imod/mf6/wel.py index ab4a62fa8..83eca04bd 100644 --- a/imod/mf6/wel.py +++ b/imod/mf6/wel.py @@ -28,6 +28,7 @@ ) from imod.mf6.mf6_wel_adapter import Mf6Wel, concat_indices_to_cellid from imod.mf6.package import Package +from imod.mf6.regrid.regrid_schemes import CapDataWellRegridMethod from imod.mf6.utilities.dataset import remove_inactive from imod.mf6.utilities.imod5_converter import well_from_imod5_cap_data from imod.mf6.validation_settings import ValidationSettings @@ -45,6 +46,7 @@ from imod.typing import GridDataArray, Imod5DataDict, StressPeriodTimesType from imod.typing.grid import is_spatial_grid, ones_like from imod.util.expand_repetitions import average_timeseries, resample_timeseries +from imod.util.regrid import RegridderWeightsCache from imod.util.structured import values_within_range ABSTRACT_METH_ERROR_MSG = "Method in abstract base class called" @@ -1443,7 +1445,13 @@ def _validate_imod5_depth_information( raise ValueError(log_msg) @classmethod - def from_imod5_cap_data(cls, imod5_data: Imod5DataDict): + def from_imod5_cap_data( + cls, + imod5_data: Imod5DataDict, + target_dis: StructuredDiscretization, + regridder_types: CapDataWellRegridMethod = CapDataWellRegridMethod(), + regrid_cache: RegridderWeightsCache = RegridderWeightsCache(), + ): """ Create LayeredWell from imod5_data in "cap" package. Abstraction data for sprinkling is defined in iMOD5 either with grids (IDF) or points @@ -1481,7 +1489,9 @@ def from_imod5_cap_data(cls, imod5_data: Imod5DataDict): belongs, as returned by :func:`imod.formats.prj.open_projectfile_data`. """ - data = well_from_imod5_cap_data(imod5_data) + data = well_from_imod5_cap_data( + imod5_data, target_dis, regridder_types, regrid_cache + ) return cls(**data) # type: ignore @standard_log_decorator() From 406e5ade6cb86946e6cc00c7fd90cfe90cde4121 Mon Sep 17 00:00:00 2001 From: JoerivanEngelen Date: Wed, 29 Apr 2026 14:29:14 +0200 Subject: [PATCH 15/19] fix mypy error --- imod/mf6/model_gwf.py | 11 ++++++----- 1 file changed, 6 insertions(+), 5 deletions(-) diff --git a/imod/mf6/model_gwf.py b/imod/mf6/model_gwf.py index 40cf867b2..f834aa1a7 100644 --- a/imod/mf6/model_gwf.py +++ b/imod/mf6/model_gwf.py @@ -45,7 +45,7 @@ SimulationDistributingOptions, ) from imod.schemata import TypeSchema -from imod.typing import GridDataArray, StressPeriodTimesType +from imod.typing import GridDataArray, Imod5DataDict, StressPeriodTimesType from imod.typing.grid import zeros_like from imod.util.regrid import RegridderWeightsCache @@ -354,12 +354,13 @@ def from_imod5_data( ) if "cap" in imod5_keys: + imod5_cap_data_dict = cast(Imod5DataDict, imod5_data) result["msw-sprinkling"] = LayeredWell.from_imod5_cap_data( - imod5_data, dis_pkg, regrid_cache=regrid_cache - ) # type: ignore + imod5_cap_data_dict, dis_pkg, regrid_cache=regrid_cache + ) result["msw-rch"] = Recharge.from_imod5_cap_data( - imod5_data, dis_pkg, regrid_cache=regrid_cache - ) # type: ignore + imod5_cap_data_dict, dis_pkg, regrid_cache=regrid_cache + ) # import ghb's ghb_keys = [key for key in imod5_keys if key[0:3] == "ghb"] From 3af5ee287dfc1c3e0acb1792ef5ff0bc16338242 Mon Sep 17 00:00:00 2001 From: JoerivanEngelen Date: Wed, 29 Apr 2026 14:36:50 +0200 Subject: [PATCH 16/19] Fix mypy errors Co-authored-by: Copilot --- imod/common/utilities/regrid.py | 3 +-- imod/mf6/utilities/imod5_converter.py | 4 ++-- 2 files changed, 3 insertions(+), 4 deletions(-) diff --git a/imod/common/utilities/regrid.py b/imod/common/utilities/regrid.py index 70ee6fc9b..bffb13cb4 100644 --- a/imod/common/utilities/regrid.py +++ b/imod/common/utilities/regrid.py @@ -20,7 +20,6 @@ from imod.common.utilities.dataclass_type import DataclassType, EmptyRegridMethod from imod.common.utilities.dtype import is_integer from imod.common.utilities.value_filters import is_valid -from imod.mf6 import StructuredDiscretization from imod.typing import Imod5DataDict from imod.typing.grid import ( GridDataArray, @@ -465,7 +464,7 @@ def _get_regridding_domain( def regrid_imod5_cap_data( imod5_data: Imod5DataDict, - target_dis: StructuredDiscretization, + target_dis: IRegridPackage, regridder_types: DataclassType, regrid_cache: RegridderWeightsCache, ) -> Imod5DataDict: diff --git a/imod/mf6/utilities/imod5_converter.py b/imod/mf6/utilities/imod5_converter.py index b62834e21..adee28a6e 100644 --- a/imod/mf6/utilities/imod5_converter.py +++ b/imod/mf6/utilities/imod5_converter.py @@ -4,9 +4,9 @@ import pandas as pd import xarray as xr +from imod.common.interfaces.iregridpackage import IRegridPackage from imod.common.utilities.dataclass_type import DataclassType from imod.common.utilities.regrid import _regrid_package_data, regrid_imod5_cap_data -from imod.mf6.dis import StructuredDiscretization from imod.mf6.package import Package from imod.typing import GridDataArray, GridDataDict, Imod5DataDict from imod.typing.grid import full_like @@ -85,7 +85,7 @@ def _well_from_imod5_cap_grid_data(cap_data: GridDataDict) -> dict[str, np.ndarr def well_from_imod5_cap_data( imod5_data: Imod5DataDict, - target_dis: StructuredDiscretization, + target_dis: IRegridPackage, regridder_types: DataclassType, regrid_cache: RegridderWeightsCache, ) -> dict[str, np.ndarray]: From ca4a6930fbd9dcc8242a9a47a46a283543a40f76 Mon Sep 17 00:00:00 2001 From: JoerivanEngelen Date: Wed, 29 Apr 2026 14:37:11 +0200 Subject: [PATCH 17/19] Reimplement type ignore, as it turned out to be necessary Co-authored-by: Copilot --- imod/mf6/riv.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/imod/mf6/riv.py b/imod/mf6/riv.py index 3cf85667f..0957ebbad 100644 --- a/imod/mf6/riv.py +++ b/imod/mf6/riv.py @@ -141,7 +141,7 @@ def _create_drain_from_leftover_riv_imod5_data( else: drain_leftover_data = infiltration_drn_data - return Drainage(**drain_leftover_data) + return Drainage(**drain_leftover_data) # type: ignore class River(TopSystemBoundaryCondition, IRegridPackage): From 62c4090c7f175b0fe5bf77da1e5f0889ff5fcb48 Mon Sep 17 00:00:00 2001 From: JoerivanEngelen Date: Wed, 29 Apr 2026 14:54:47 +0200 Subject: [PATCH 18/19] Update changelog --- docs/api/changelog.rst | 6 ++++-- 1 file changed, 4 insertions(+), 2 deletions(-) diff --git a/docs/api/changelog.rst b/docs/api/changelog.rst index 37f424888..e71ba80b4 100644 --- a/docs/api/changelog.rst +++ b/docs/api/changelog.rst @@ -64,8 +64,10 @@ Fixed - :meth:`imod.mf6.NodePropertyFlow.regrid_like` now regrids ``k33`` using the correct method, namely ``mean`` instead of ``harmonic_mean``. As this is the appropriate method for horizontal regridding of ``k33``. -- :meth:`imod.msw.MetaSwapModel.from_imod5_data` now regrids the iMOD5 data to - the MODFLOW6 target discretization. +- :meth:`imod.msw.MetaSwapModel.from_imod5_data`, + :meth:`imod.mf6.Recharge.from_imod5_cap_data`, + :meth:`imod.mf6.LayeredWell.from_imod5_cap_data` now regrids the iMOD5 CAP + data to the MODFLOW6 target discretization. Changed ~~~~~~~ From 7c7e828b161eb9462281eb98b55101803e96e64d Mon Sep 17 00:00:00 2001 From: JoerivanEngelen Date: Thu, 30 Apr 2026 11:26:52 +0200 Subject: [PATCH 19/19] Update tests --- imod/tests/conftest.py | 2 + imod/tests/fixtures/imod5_cap_data.py | 45 ++++++++++++++++++++-- imod/tests/test_mf6/test_mf6_rch.py | 3 ++ imod/tests/test_mf6/test_mf6_simulation.py | 1 + imod/tests/test_mf6/test_mf6_wel.py | 20 +++++++--- imod/tests/test_msw/test_model.py | 2 +- 6 files changed, 62 insertions(+), 11 deletions(-) diff --git a/imod/tests/conftest.py b/imod/tests/conftest.py index 9bc606476..d1f66a799 100644 --- a/imod/tests/conftest.py +++ b/imod/tests/conftest.py @@ -23,6 +23,8 @@ ) from .fixtures.flow_transport_simulation_fixture import flow_transport_simulation from .fixtures.imod5_cap_data import ( + cap_coupled_dis_grid, + cap_coupled_dis_grid__big, cap_data_sprinkling_grid, cap_data_sprinkling_grid__big, cap_data_sprinkling_points, diff --git a/imod/tests/fixtures/imod5_cap_data.py b/imod/tests/fixtures/imod5_cap_data.py index e74f1664c..88a2a5d7d 100644 --- a/imod/tests/fixtures/imod5_cap_data.py +++ b/imod/tests/fixtures/imod5_cap_data.py @@ -3,6 +3,7 @@ import pytest import xarray as xr +from imod.mf6 import StructuredDiscretization from imod.typing import Imod5DataDict @@ -24,6 +25,42 @@ def zeros_dask_grid(n): return da.chunk({"x": n, "y": n}) +@pytest.fixture(scope="function") +def cap_coupled_dis_grid() -> StructuredDiscretization: + n = 3 + layer_coord = np.arange(n) + 1 + layer_template = xr.DataArray( + layer_coord, dims="layer", coords={"layer": layer_coord} + ) + template_grid = layer_template * zeros_grid(n) + + idomain = template_grid.astype(int) + 1 + top = 10.0 + bottom = xr.DataArray( + [0.0, -10.0, -20.0], dims="layer", coords={"layer": layer_coord} + ) + + return StructuredDiscretization(idomain=idomain, top=top, bottom=bottom) + + +@pytest.fixture(scope="function") +def cap_coupled_dis_grid__big() -> StructuredDiscretization: + n = 1000 + layer_coord = np.arange(3) + 1 + layer_template = xr.DataArray( + layer_coord, dims="layer", coords={"layer": layer_coord} + ) + template_grid = layer_template * zeros_dask_grid(n) + + idomain = template_grid.astype(int) + 1 + top = 10.0 + bottom = xr.DataArray( + [0.0, -10.0, -20.0], dims="layer", coords={"layer": layer_coord} + ) + + return StructuredDiscretization(idomain=idomain, top=top, bottom=bottom) + + @pytest.fixture(scope="function") def cap_data_sprinkling_grid() -> Imod5DataDict: n = 3 @@ -46,7 +83,7 @@ def cap_data_sprinkling_grid() -> Imod5DataDict: "artificial_recharge_capacity": xr.DataArray(25.0), } - return {"cap": cap_data} + return {"cap": cap_data, "extra": {"paths": ["path1", "path2"]}} @pytest.fixture(scope="function") @@ -69,7 +106,7 @@ def cap_data_sprinkling_grid__big() -> Imod5DataDict: "artificial_recharge_capacity": xr.DataArray(25.0), } - return {"cap": cap_data} + return {"cap": cap_data, "extra": {"paths": ["path1", "path2"]}} @pytest.fixture(scope="function") @@ -94,7 +131,7 @@ def cap_data_sprinkling_grid_with_layer() -> Imod5DataDict: "artificial_recharge_capacity": xr.DataArray(25.0), } - return {"cap": cap_data} + return {"cap": cap_data, "extra": {"paths": ["path1", "path2"]}} @pytest.fixture(scope="function") @@ -126,4 +163,4 @@ def cap_data_sprinkling_points() -> Imod5DataDict: "artificial_recharge_capacity": xr.DataArray(25.0), } - return {"cap": cap_data} + return {"cap": cap_data, "extra": {"paths": ["path1", "path2"]}} diff --git a/imod/tests/test_mf6/test_mf6_rch.py b/imod/tests/test_mf6/test_mf6_rch.py index 16fef6bb9..5867234f0 100644 --- a/imod/tests/test_mf6/test_mf6_rch.py +++ b/imod/tests/test_mf6/test_mf6_rch.py @@ -578,6 +578,7 @@ def test_from_imod5_cap_data(imod5_dataset): # Arrange data = deepcopy(imod5_dataset[0]) target_discretization = StructuredDiscretization.from_imod5_data(data) + data["extra"] = {"paths": ["path1", "path2"]} data["cap"] = {} msw_bound = data["bnd"]["ibound"].isel(layer=0, drop=False) data["cap"]["boundary"] = msw_bound @@ -613,6 +614,7 @@ def test_from_imod5_cap_data__regrid(imod5_dataset): # Arrange data = deepcopy(imod5_dataset[0]) target_discretization = StructuredDiscretization.from_imod5_data(data) + data["extra"] = {"paths": ["path1", "path2"]} data["cap"] = {} msw_bound = data["bnd"]["ibound"].isel(layer=0, drop=False) data["cap"]["boundary"] = msw_bound @@ -637,6 +639,7 @@ def test_from_imod5_cap_data__clip_box(imod5_dataset): # Arrange data = deepcopy(imod5_dataset[0]) target_discretization = StructuredDiscretization.from_imod5_data(data) + data["extra"] = {"paths": ["path1", "path2"]} data["cap"] = {} msw_bound = data["bnd"]["ibound"].isel(layer=0, drop=False) data["cap"]["boundary"] = msw_bound diff --git a/imod/tests/test_mf6/test_mf6_simulation.py b/imod/tests/test_mf6/test_mf6_simulation.py index ef7386a01..35406a52c 100644 --- a/imod/tests/test_mf6/test_mf6_simulation.py +++ b/imod/tests/test_mf6/test_mf6_simulation.py @@ -615,6 +615,7 @@ def test_from_imod5__has_cap_data(imod5_dataset): imod5_data = deepcopy(imod5_dataset[0]) period_data = imod5_dataset[1] + imod5_data["extra"] = {"paths": ["path1", "path2"]} imod5_data["cap"] = {} msw_bound = imod5_data["bnd"]["ibound"].isel(layer=0, drop=True) imod5_data["cap"]["boundary"] = msw_bound diff --git a/imod/tests/test_mf6/test_mf6_wel.py b/imod/tests/test_mf6/test_mf6_wel.py index c37f71e7d..94d240f87 100644 --- a/imod/tests/test_mf6/test_mf6_wel.py +++ b/imod/tests/test_mf6/test_mf6_wel.py @@ -1124,14 +1124,16 @@ def test_logmessage_for_missing_filter_settings( @pytest.mark.unittest_jit -def test_from_imod5_cap_data__grid(cap_data_sprinkling_grid): +def test_from_imod5_cap_data__grid(cap_data_sprinkling_grid, cap_coupled_dis_grid): # Arrange expected_layer = np.array([2, 1, 2, 1, 2, 1]) expected_y = np.array([3.0, 3.0, 2.0, 2.0, 1.0, 1.0]) expected_x = np.array([2.0, 3.0, 2.0, 3.0, 2.0, 3.0]) # Act - well = LayeredWell.from_imod5_cap_data(cap_data_sprinkling_grid) + well = LayeredWell.from_imod5_cap_data( + cap_data_sprinkling_grid, cap_coupled_dis_grid + ) # Assert ds = well.dataset @@ -1143,20 +1145,26 @@ def test_from_imod5_cap_data__grid(cap_data_sprinkling_grid): @pytest.mark.unittest_jit @pytest.mark.timeout(300) -def test_from_imod5_cap_data__big_grid(cap_data_sprinkling_grid__big): +def test_from_imod5_cap_data__big_grid( + cap_data_sprinkling_grid__big, cap_coupled_dis_grid__big +): """Test if performance is acceptable for large grids.""" # Arrange bnd_2d = cap_data_sprinkling_grid__big["cap"]["boundary"] layer = xr.DataArray([1, 1], coords={"layer": [1, 2]}, dims=["layer"]) bnd = layer * bnd_2d # Act - well = LayeredWell.from_imod5_cap_data(cap_data_sprinkling_grid__big) + well = LayeredWell.from_imod5_cap_data( + cap_data_sprinkling_grid__big, cap_coupled_dis_grid__big + ) mf6_wel = well.to_mf6_pkg(bnd.astype(bool), bnd, bnd, bnd) # Assert assert mf6_wel.dataset.sizes["ncellid"] == bnd_2d.size @pytest.mark.unittest_jit -def test_from_imod5_cap_data__points(cap_data_sprinkling_points): +def test_from_imod5_cap_data__points(cap_data_sprinkling_points, cap_coupled_dis_grid): with pytest.raises(NotImplementedError): - LayeredWell.from_imod5_cap_data(cap_data_sprinkling_points) + LayeredWell.from_imod5_cap_data( + cap_data_sprinkling_points, cap_coupled_dis_grid + ) diff --git a/imod/tests/test_msw/test_model.py b/imod/tests/test_msw/test_model.py index 8e584bf2e..aeaaf887c 100644 --- a/imod/tests/test_msw/test_model.py +++ b/imod/tests/test_msw/test_model.py @@ -435,7 +435,7 @@ def test_import_from_imod5_and_write( modeldir = tmp_path / "modeldir" # Act model = msw.MetaSwapModel.from_imod5_data(imod5_data, dis_pkg, times) - well_pkg = mf6.LayeredWell.from_imod5_cap_data(imod5_data) + well_pkg = mf6.LayeredWell.from_imod5_cap_data(imod5_data, dis_pkg) mf6_wel_pkg = well_pkg.to_mf6_pkg( active, dis_pkg["top"], dis_pkg["bottom"], npf_pkg["k"] )