diff --git a/docs/source/config-file.rst b/docs/source/config-file.rst index 1e21a716..53460f05 100644 --- a/docs/source/config-file.rst +++ b/docs/source/config-file.rst @@ -143,8 +143,25 @@ Some additional notes on YAML * ``[12,1.2]`` parses to ``[12,1.2]`` * ``(12,1.2)`` parses to ``"(12,1.2)"`` * ``{12,1.2}`` parses to ``{12: None, 1.2: None}`` +* dictionaries can be represented with indentation, but spaces are needed after the colon(s): + .. code-block:: yaml + items: + 0:1 + 1:2 + + parses to ``'0:1 1:2'`` + + .. code-block:: yaml + + items: + 0: 1 + 1: 2 + + parses to ``{0:1, 1:2}`` + +Using a YAML-aware text editor such as VS Code can help with these issues, for example by changing the highlighting color to indicate a string in the first dictionary example above and an interpreted python data structure in the second dictionary example. .. _JSON: https://www.json.org/json-en.html .. _pyyaml: https://pyyaml.org/ diff --git a/mfsetup/grid.py b/mfsetup/grid.py index b005037c..487e5101 100644 --- a/mfsetup/grid.py +++ b/mfsetup/grid.py @@ -634,7 +634,7 @@ def write_bbox_shapefile(modelgrid, outshp): def rasterize(feature, grid, id_column=None, - include_ids=None, names_column=None, + include_ids=None, exclude_ids=None, names_column=None, crs=None, **kwargs): """Rasterize a feature onto the model grid, using the rasterio.features.rasterize method. Features are intersected @@ -649,6 +649,8 @@ def rasterize(feature, grid, id_column=None, from this column will be assigned to the output raster. include_ids : sequence Subset of IDs in id_column to include + exclude_ids : sequence + Subset of IDs in id_column to exclude names_column : str, optional By default, the IDs in id_column, or sequential integers are returned. This option allows another column of strings @@ -738,7 +740,8 @@ def rasterize(feature, grid, id_column=None, # subset to include_ids if id_column is not None and include_ids is not None: df = df.loc[df[id_column].isin(include_ids)].copy() - + if id_column is not None and exclude_ids is not None: + df = df.loc[~df[id_column].isin(exclude_ids)].copy() # create list of GeoJSON features, with unique value for each feature if id_column is None: numbers = list(range(1, len(df)+1)) @@ -986,14 +989,14 @@ def setup_structured_grid(xoff=None, yoff=None, xul=None, yul=None, if parent_model is not None and snap_to_parent: to_grid_units_parent = convert_length_units(get_model_length_units(parent_model), grid_units) # parent model grid spacing in meters - parent_delr_grid = np.round(parent_model.dis.delr.array[0] * to_grid_units_parent, 4) - if not parent_delr_grid % delr_grid % parent_delr_grid == 0: - raise ValueError('inset delr spacing of {} must be factor of parent spacing of {}'.format(delr_grid, - parent_delr_grid)) - parent_delc_grid = np.round(parent_model.dis.delc.array[0] * to_grid_units_parent, 4) - if not parent_delc_grid % delc_grid % parent_delc_grid == 0: - raise ValueError('inset delc spacing of {} must be factor of parent spacing of {}'.format(delc_grid, - parent_delc_grid)) + #parent_delr_grid = np.round(parent_model.dis.delr.array[0] * to_grid_units_parent, 4) + #if not parent_delr_grid % delr_grid % parent_delr_grid == 0: + # raise ValueError('inset delr spacing of {} must be factor of parent spacing of {}'.format(delr_grid, + # parent_delr_grid)) + #parent_delc_grid = np.round(parent_model.dis.delc.array[0] * to_grid_units_parent, 4) + #if not parent_delc_grid % delc_grid % parent_delc_grid == 0: + # raise ValueError('inset delc spacing of {} must be factor of parent spacing of {}'.format(delc_grid, + # parent_delc_grid)) # option 1: make grid from xoff, yoff and specified dimensions if xoff is not None and yoff is not None: @@ -1057,6 +1060,10 @@ def setup_structured_grid(xoff=None, yoff=None, xul=None, yul=None, df.to_crs(crs, inplace=True) # use all features by default features = df.geometry.tolist() + elif features is None and features_shapefile is not None: + raise ValueError( + "setup_grid: need one of xoff/yoff, xul/yul, features_shapefile or " + "features inputs") # alternatively, accept features as an argument # convert multiple features to a MultiPolygon if isinstance(features, list): diff --git a/mfsetup/mf6model.py b/mfsetup/mf6model.py index 11edd56e..c480f17c 100644 --- a/mfsetup/mf6model.py +++ b/mfsetup/mf6model.py @@ -405,7 +405,7 @@ def create_lgr_models(self): last_refined_layer = max(np.where(is_refined > 0)[0]) consecutive = all(np.diff(is_refined)[:last_refined_layer] == 0) if (is_refined[0] != 1) | (not consecutive): - raise ValueError("Configuration input: layer_refinement must" + raise ValueError("Configuration input: layer_refinement must " "include consecutive sequence of layers, " "starting with the top layer.") # check the specified DIS package input is consistent diff --git a/mfsetup/mfmodel.py b/mfsetup/mfmodel.py index 3829c88d..4d4bd8cc 100644 --- a/mfsetup/mfmodel.py +++ b/mfsetup/mfmodel.py @@ -1417,6 +1417,12 @@ def setup_grid(self): self._modelgrid.cfg = self.cfg['grid'] else: kwargs = get_input_arguments(cfg, setup_structured_grid) + if not set(kwargs.keys()).intersection({ + 'features_shapefile', 'features', 'xoff', 'yoff', 'xul', 'yul'}): + raise ValueError( + "No features_shapefile or xoff, yoff supplied " + "to setup_grid: block. Check configuration file input, " + "including for accidental indentation of the setup_grid: block.") self._modelgrid = setup_structured_grid(**kwargs) self.cfg['grid'] = self._modelgrid.cfg # update DIS package configuration diff --git a/mfsetup/tests/data/shellmound.yml b/mfsetup/tests/data/shellmound.yml index 3917a438..6ebae823 100644 --- a/mfsetup/tests/data/shellmound.yml +++ b/mfsetup/tests/data/shellmound.yml @@ -317,10 +317,10 @@ drn: shapefile: filename: '../../../mfsetup/tests/data/shellmound/shps/waterbodies.shp' id_column: 'COMID' - # Note: to include all features in a shapefile, - # simply omit the include_ids: key - # id_column: can also be omitted - include_ids: [18047154, 18046236] + # Include all features in the above shapefile, + # except those associated with these COMIDs + exclude_ids: [18046230, 18046226, 18046238, 17953939, 18046140, 18046162] + boundname_column: 'COMID' elev: filename: 'shellmound/rasters/meras_100m_dem.tif' elevation_units: 'feet' diff --git a/mfsetup/tests/test_mf6_shellmound.py b/mfsetup/tests/test_mf6_shellmound.py index 6391d32d..a3095cea 100644 --- a/mfsetup/tests/test_mf6_shellmound.py +++ b/mfsetup/tests/test_mf6_shellmound.py @@ -650,6 +650,11 @@ def test_basic_stress_package_setup(shellmound_model_with_dis, pckg_abbrv, # heads in CHD package should match the CSV input, # after conversion to meters assert np.allclose(in_df['head'].values, spd_heads[1:]/.3048) + if pckg_abbrv == 'drn': + exclude_ids = list(map(str, m.cfg[pckg_abbrv]['source_data']['shapefile']['exclude_ids'])) + comids = [s.replace('feature-','') for s in np.unique(m.drn.stress_period_data.data[0]['boundname'])] + assert not set(exclude_ids).intersection(comids) + assert comids == ['18046236', '18047154'] # more advanced transient input case with csv # and conductance values specified via a raster if pckg_abbrv == 'ghb':