Coverage for C:\src\imod-python\imod\msw\grid_data.py: 36%
39 statements
« prev ^ index » next coverage.py v7.4.4, created at 2024-04-08 10:26 +0200
« prev ^ index » next coverage.py v7.4.4, created at 2024-04-08 10:26 +0200
1import numpy as np
2import xarray as xr
4from imod.msw.fixed_format import VariableMetaData
5from imod.msw.pkgbase import MetaSwapPackage
6from imod.util.spatial import spatial_reference
9class GridData(MetaSwapPackage):
10 """
11 This contains the grid data of MetaSWAP.
13 This class is responsible for the file `area_svat.inp`
15 Parameters
16 ----------
17 area: array of floats (xr.DataArray)
18 Describes the area of SVAT units. This array must have a subunit coordinate
19 to describe different landuses.
20 landuse: array of integers (xr.DataArray)
21 Describes the landuse type of SVAT units.
22 This array must have a subunit coordinate.
23 rootzone_depth: array of floats (xr.DataArray)
24 Describes the rootzone depth of SVAT units.
25 This array must have a subunit coordinate to describe different landuses.
26 surface_elevation: array of floats (xr.DataArray)
27 Describes the surface elevation of SVAT units.
28 This array must not have a subunit coordinate.
29 soil_physical_unit: array of integers (xr.DataArray)
30 Describes the physical parameters of SVAT units.
31 These parameters will be looked up in a table according to the given integers.
32 This array must not have a subunit coordinate.
33 active: array of bools (xr.DataArray)
34 Describes whether SVAT units are active or not.
35 This array must not have a subunit coordinate.
36 """
38 _file_name = "area_svat.inp"
39 _metadata_dict = {
40 "svat": VariableMetaData(10, 1, 99999999, int),
41 "area": VariableMetaData(10, 0.0, 999999.0, float),
42 "surface_elevation": VariableMetaData(8, -9999.0, 9999.0, float),
43 "temp": VariableMetaData(8, None, None, str),
44 "soil_physical_unit": VariableMetaData(6, 1, 999999, int),
45 "soil_physical_unit_string": VariableMetaData(16, None, None, str),
46 "landuse": VariableMetaData(6, 1, 999999, int),
47 "rootzone_depth": VariableMetaData(8, 0.0, 10.0, float),
48 }
49 _with_subunit = ("area", "landuse", "rootzone_depth")
50 _without_subunit = ("surface_elevation", "soil_physical_unit")
51 _to_fill = ("soil_physical_unit_string", "temp")
53 def __init__(
54 self,
55 area: xr.DataArray,
56 landuse: xr.DataArray,
57 rootzone_depth: xr.DataArray,
58 surface_elevation: xr.DataArray,
59 soil_physical_unit: xr.DataArray,
60 active: xr.DataArray,
61 ):
62 super().__init__()
63 self.dataset["area"] = area
64 self.dataset["landuse"] = landuse
65 self.dataset["rootzone_depth"] = rootzone_depth
66 self.dataset["surface_elevation"] = surface_elevation
67 self.dataset["soil_physical_unit"] = soil_physical_unit
68 self.dataset["active"] = active
70 self._pkgcheck()
72 def generate_index_array(self):
73 """
74 Generate index arrays to be used on other packages
75 """
76 area = self.dataset["area"]
77 active = self.dataset["active"]
79 isactive = area.where(active).notnull()
81 index = isactive.values.ravel()
83 svat = xr.full_like(area, fill_value=0, dtype=np.int64).rename("svat")
84 svat.values[isactive.values] = np.arange(1, index.sum() + 1)
86 return index, svat
88 def _pkgcheck(self):
89 super()._pkgcheck()
91 dx, _, _, dy, _, _ = spatial_reference(self.dataset)
93 if (not np.isscalar(dx)) or (not np.isscalar(dy)):
94 raise ValueError("MetaSWAP only supports equidistant grids")
96 active = self.dataset["active"]
98 cell_area = active.astype(float) * dx * abs(dy)
99 total_area = self.dataset["area"].sum(dim="subunit")
101 # Apparently all regional models intentionally provided area grids
102 # smaller than cell area, to allow surface waters as workaround.
103 unequal_area = (total_area > cell_area).values[active.values]
105 if np.any(unequal_area):
106 raise ValueError(
107 "Provided area grid with total areas larger than cell area"
108 )