Minor changes

nlmod/gwf/hfb.py

@@ -54,9 +54,9 @@ def get_hfb_spd(ds, linestrings, hydchr, depth=None, elevation=None):
     spd : List of Tuple
         Stress period data used to configure the hfb package of Flopy.
     """
-    assert sum([depth is None, elevation is None]) == 1, (
-        "Use either depth or elevation argument"
-    )
+    assert (
+        sum([depth is None, elevation is None]) == 1
+    ), "Use either depth or elevation argument"
 
     if not isinstance(ds, xr.Dataset):
         raise TypeError("Please pass a model dataset!")
@@ -414,7 +414,9 @@ def line_to_hfb_buffer(gdf, ds, buffer_distance=None, gi=None):
     return hfb_seg
 
 
-def polygon_to_hfb(gdf, ds, hydchr, column=None, gwf=None, lay=0, add_data=False):
+def polygon_to_hfb(
+    gdf, ds, hydchr, column=None, gwf=None, lay=0, add_data=False, include_nan=True
+):
     """Snap polygon exterior to grid to form a horizontal flow barrier.
 
     Parameters
@@ -435,6 +437,8 @@ def polygon_to_hfb(gdf, ds, hydchr, column=None, gwf=None, lay=0, add_data=False
         Layer number. The default is 0.
     add_data : bool, optional
         If True, add the data to the stress period data. The default is False.
+    include_nan : bool, optional
+        If true, add hfb's to and from nan-values
 
     Returns
     -------
@@ -475,6 +479,8 @@ def polygon_to_hfb(gdf, ds, hydchr, column=None, gwf=None, lay=0, add_data=False
 
         edges = []
         for icell2d in range(icvert.shape[0]):
+            if not include_nan and np.isnan(data[icell2d]):
+                continue
             for j in range(icvert.shape[1] - 1):
                 if icvert[icell2d, j + 1] == nodata:
                     break
@@ -493,7 +499,7 @@ def polygon_to_hfb(gdf, ds, hydchr, column=None, gwf=None, lay=0, add_data=False
                 icell2ds.append(edges[mask, 0])
         # icell2ds = np.array(icell2ds)
         for icell2d1, icell2d2 in icell2ds:
-            spd.append([(lay, icell2d1), (lay, icell2d2), hydchr])
+            spd.append([(lay, int(icell2d1)), (lay, int(icell2d2)), hydchr])
             if add_data:
                 spd[-1].extend([data[icell2d1], data[icell2d2]])
     if gwf is None:
@@ -502,49 +508,58 @@ def polygon_to_hfb(gdf, ds, hydchr, column=None, gwf=None, lay=0, add_data=False
         return flopy.mf6.ModflowGwfhfb(gwf, stress_period_data={0: spd})
 
 
-def plot_hfb(cellids, gwf, ax=None, color="red", **kwargs):
+def plot_hfb(cellids, modelgrid, ax=None, color="red", **kwargs):
     """Plots a horizontal flow barrier.
 
     Parameters
     ----------
     cellids : list of lists of integers or flopy.mf6.ModflowGwfhfb
         list with the ids of adjacent cells that should get a horizontal
-        flow barrier, hfb is the output of line2hfb.
-    gwf : flopy groundwater flow model
-        DESCRIPTION.
+        flow barrier, hfb is the output of line_to_hfb.
+    modelgrid : flopy modelgrid, flopy groundwater flow model, xarray Dataset or
+        The object containing the properties of the model grid. It can be a flopy
+        modelgrid-object, a flopy groundwater flow model, or an nlmod xarray dataset.
     ax : matplotlib axes
+        The Axes to plot the hfb-data in. When ax is None, a new figure is created.
+        The default is None.
+    color : str
+        The color of the hfb-lines in the plot. The default is 'red'.
 
 
     Returns
     -------
-    fig : TYPE
-        DESCRIPTION.
-    ax : TYPE
-        DESCRIPTION.
+    ax : matplotlib.Axes
+        The ax that conatins the hfb-lines.
     """
     if ax is None:
         _, ax = plt.subplots()
 
-    if gwf.modelgrid.grid_type == "structured":
+    if not isinstance(modelgrid, flopy.discretization.grid.Grid):
+        if isinstance(modelgrid, xr.Dataset):
+            modelgrid = modelgrid_from_ds(modelgrid)
+        else:
+            modelgrid = modelgrid.modelgrid
+
+    if modelgrid.grid_type == "structured":
         if isinstance(cellids, flopy.mf6.ModflowGwfhfb):
             spd = cellids.stress_period_data.data[0]
             cellids = [[row[0][1:], row[1][1:]] for row in spd]
         for line in cellids:
-            pc1 = Polygon(gwf.modelgrid.get_cell_vertices(*line[0]))
-            pc2 = Polygon(gwf.modelgrid.get_cell_vertices(*line[1]))
+            pc1 = Polygon(modelgrid.get_cell_vertices(*line[0]))
+            pc2 = Polygon(modelgrid.get_cell_vertices(*line[1]))
             x, y = pc1.intersection(pc2).xy
             ax.plot(x, y, color=color, **kwargs)
 
-    elif gwf.modelgrid.grid_type == "vertex":
+    elif modelgrid.grid_type == "vertex":
         if isinstance(cellids, flopy.mf6.ModflowGwfhfb):
             spd = cellids.stress_period_data.data[0]
             cellids = [[line[0][1], line[1][1]] for line in spd]
         for line in cellids:
-            pc1 = Polygon(gwf.modelgrid.get_cell_vertices(line[0]))
-            pc2 = Polygon(gwf.modelgrid.get_cell_vertices(line[1]))
+            pc1 = Polygon(modelgrid.get_cell_vertices(line[0]))
+            pc2 = Polygon(modelgrid.get_cell_vertices(line[1]))
             x, y = pc1.intersection(pc2).xy
             ax.plot(x, y, color=color, **kwargs)
     else:
-        raise ValueError(f"not supported gridtype -> {gwf.modelgrid.grid_type}")
+        raise ValueError(f"not supported gridtype -> {modelgrid.grid_type}")
 
     return ax

nlmod/gwf/surface_water.py

@@ -687,7 +687,7 @@ def _get_waterboard_selection(gdf=None, extent=None, config=None):
             if config[wb]["bgt_code"] in bronhouders:
                 wbs.append(wb)
     elif extent is not None:
-        wb_gdf = waterboard.get_polygons()
+        wb_gdf = waterboard.download_polygons()
         wbs = wb_gdf.index[wb_gdf.intersects(extent_to_polygon(extent))]
     return wbs
 

nlmod/modpath/modpath.py

@@ -491,8 +491,14 @@ def sim(
         None
     ref_time : TYPE, optional
         DESCRIPTION. The default is None.
+    stoptime : float, optional
+        User-specified value of tracking time at which to stop a particle tracking
+        simulation. The default is None
     stoptime : TYPE, optional
         DESCRIPTION. The default is None.
+    simulationtype : str
+        MODPATH 7 simulation type. Valid simulation types are 'endpoint',
+        'pathline', 'timeseries', or 'combined' (default is 'pathline').
 
     Returns
     -------

nlmod/read/waterboard.py

@@ -72,7 +72,14 @@ def get_configuration():
         },
         "weirs": {
             "url": "https://gisservices.aaenmaas.nl/arcgis/rest/services/EXTERN/Legger/MapServer",
-            "layer": 7,
+            "layer": 7,  # Stuw (7)
+            "index": "CODE",
+            "summer_stage": "HOOGSTEDOORSTROOMHOOGTE",
+            "winter_stage": "LAAGSTEDOORSTROOMHOOGTE",
+        },
+        "culverts": {
+            "url": "https://gisservices.aaenmaas.nl/arcgis/rest/services/EXTERN/Legger/MapServer",
+            "layer": 4,  # Duiker / Sifon (4)
             "index": "CODE",
             "summer_stage": "HOOGSTEDOORSTROOMHOOGTE",
             "winter_stage": "LAAGSTEDOORSTROOMHOOGTE",
@@ -121,6 +128,14 @@ def get_configuration():
             # "layer": 14,  # categorie B
             # "layer": 15,  # categorie C
         },
+        "weirs": {
+            "url": "https://geoservices.brabantsedelta.nl/arcgis/rest/services/EXTERN/WEB_Vastgestelde_Legger_Oppervlaktewaterlichamen/FeatureServer",
+            "layer": 4,  # Stuw (4)
+        },
+        "culverts": {
+            "url": "https://geoservices.brabantsedelta.nl/arcgis/rest/services/EXTERN/WEB_Vastgestelde_Legger_Oppervlaktewaterlichamen/FeatureServer",
+            "layer": 8,  # Duiker (8)
+        },
         "level_areas": {
             "url": "https://geoservices.brabantsedelta.nl/arcgis/rest/services/EXTERN/WEB_Peilbesluiten/MapServer",
             "layer": 0,  # Peilgebied vigerend
@@ -148,9 +163,22 @@ def get_configuration():
             # "layer": 2,  # LOW_2021_Profielpunt
             # "layer": 13,  # LOW_2021_Profiellijn
             # "index": "WS_PROFIELID",
-            "url": "https://services8.arcgis.com/dmR647kStmcYa6EN/arcgis/rest/services/A_watergang_Beheerregister/FeatureServer",
-            "layer": 0,  # A_watergang_Beheerregister (0)
-            "index": "LOKAALID",
+            # "url": "https://services8.arcgis.com/dmR647kStmcYa6EN/arcgis/rest/services/A_watergang_Beheerregister/FeatureServer",
+            # "layer": 0,  # A_watergang_Beheerregister (0)
+            # "index": "LOKAALID",
+            "url": "https://services8.arcgis.com/dmR647kStmcYa6EN/arcgis/rest/services/LWW_2023_A_water_V/FeatureServer",
+        },
+        "culverts": {
+            "url": "https://services8.arcgis.com/dmR647kStmcYa6EN/arcgis/rest/services/LWW_2023_Duiker_V/FeatureServer"
+        },
+        "weirs": {
+            "url": "https://services8.arcgis.com/dmR647kStmcYa6EN/arcgis/rest/services/LWW_2023_Stuw_V/FeatureServer"
+        },
+        "profile_lines": {
+            "url": "https://services8.arcgis.com/dmR647kStmcYa6EN/ArcGIS/rest/services/LWW_2023_A_water_profiellijn_V/FeatureServer"
+        },
+        "profile_points": {
+            "url": "https://services8.arcgis.com/dmR647kStmcYa6EN/ArcGIS/rest/services/LWW_2023_A_water_profielpunt_V/FeatureServer",
         },
     }
 
@@ -403,10 +431,20 @@ def get_configuration():
     config["Rivierenland"] = {
         "bgt_code": "W0621",
         "watercourses": {
-            "url": "https://kaarten.wsrl.nl/arcgis/rest/services/Kaarten/WatersysteemLeggerVastgesteld/MapServer",
+            # "url": "https://kaarten.wsrl.nl/arcgis/rest/services/Kaarten/WatersysteemLeggerVastgesteld/MapServer",
             # "layer": 13,  # profiellijn
-            "layer": 14,  # waterloop
+            # "layer": 14,  # waterloop
             # "index": "code",
+            "url": "https://kaarten.wsrl.nl/arcgis/rest/services/Kaarten/Watersysteem_beheerregister/MapServer",
+            "layer": 12,  # Waterlopen beheerregister (12)
+        },
+        "weirs": {
+            "url": "https://kaarten.wsrl.nl/arcgis/rest/services/Kaarten/Watersysteem_beheerregister/MapServer",
+            "layer": 8,  # Stuw beheer (9)
+        },
+        "culverts": {
+            "url": "https://kaarten.wsrl.nl/arcgis/rest/services/Kaarten/Watersysteem_beheerregister/MapServer",
+            "layer": 5,  # Duiker sifon hevel beheer (5)
         },
         "level_areas": {
             "url": "https://kaarten.wsrl.nl/arcgis/rest/services/Kaarten/Peilgebieden_praktijk/FeatureServer",