From a2fd22d2a5336d9da74ff7bb36b1149f9fd7f373 Mon Sep 17 00:00:00 2001 From: Daniel Tollenaar Date: Mon, 9 Dec 2024 13:24:49 +0100 Subject: [PATCH] Fix vlakken (#196) - misc updates op fix-scripts voor NZV en DOD - samenvoegen voor LHM, inclusief patch RWS-HWS en reset van static tables AAM - uploaden feedbackformulieren --------- Co-authored-by: ngoorden <49673619+ngoorden@users.noreply.github.com> Co-authored-by: Martijn Visser --- notebooks/aa_en_maas/01b_fix_basin_area.py | 29 +-- .../00_review_model_network.py | 61 ++++++ .../01_fix_model_network.py | 15 +- .../noorderzijlvest/01b_fix_basin_area.py | 182 ++++++++++++++++-- .../rijkswaterstaat/8c_patch_model_202470 | 18 ++ notebooks/samenvoegen_modellen.py | 107 +++------- notebooks/upload_feedback_formulieren.py | 72 +++++++ notebooks/uploaden_modellen.py | 40 ++++ pixi.lock | 105 +++++++++- pixi.toml | 1 + src/ribasim_nl/ribasim_nl/model.py | 113 ++++++++--- src/ribasim_nl/ribasim_nl/network.py | 17 +- src/ribasim_nl/ribasim_nl/reset_index.py | 53 +++-- .../ribasim_nl/reset_static_tables.py | 107 ++++++++++ 14 files changed, 740 insertions(+), 180 deletions(-) create mode 100644 notebooks/drents_overijsselse_delta/00_review_model_network.py create mode 100644 notebooks/rijkswaterstaat/8c_patch_model_202470 create mode 100644 notebooks/upload_feedback_formulieren.py create mode 100644 notebooks/uploaden_modellen.py create mode 100644 src/ribasim_nl/ribasim_nl/reset_static_tables.py diff --git a/notebooks/aa_en_maas/01b_fix_basin_area.py b/notebooks/aa_en_maas/01b_fix_basin_area.py index f27f28c1..88c52c4f 100644 --- a/notebooks/aa_en_maas/01b_fix_basin_area.py +++ b/notebooks/aa_en_maas/01b_fix_basin_area.py @@ -1,10 +1,10 @@ # %% Import Libraries and Initialize Variables import geopandas as gpd -import numpy as np import pandas as pd from ribasim.nodes import level_boundary, outlet from ribasim_nl import CloudStorage, Model, NetworkValidator +from ribasim_nl.reset_static_tables import reset_static_tables # Initialize cloud storage and set authority/model parameters cloud_storage = CloudStorage() @@ -209,32 +209,7 @@ # %% -# basin-profielen/state updaten -df = pd.DataFrame( - { - "node_id": np.repeat(model.basin.node.df.index.to_numpy(), 2), - "level": [0.0, 1.0] * len(model.basin.node.df), - "area": [0.01, 1000.0] * len(model.basin.node.df), - } -) -df.index.name = "fid" -model.basin.profile.df = df - -df = model.basin.profile.df.groupby("node_id")[["level"]].max().reset_index() -df.index.name = "fid" -model.basin.state.df = df - - -# tabulated_rating_curves updaten -df = pd.DataFrame( - { - "node_id": np.repeat(model.tabulated_rating_curve.node.df.index.to_numpy(), 2), - "level": [0.0, 5] * len(model.tabulated_rating_curve.node.df), - "flow_rate": [0, 0.1] * len(model.tabulated_rating_curve.node.df), - } -) -df.index.name = "fid" -model.tabulated_rating_curve.static.df = df +model = reset_static_tables(model) model.write(ribasim_model_dir.with_stem("AaenMaas_fix_model_area") / "aam.toml") diff --git a/notebooks/drents_overijsselse_delta/00_review_model_network.py b/notebooks/drents_overijsselse_delta/00_review_model_network.py new file mode 100644 index 00000000..49153a08 --- /dev/null +++ b/notebooks/drents_overijsselse_delta/00_review_model_network.py @@ -0,0 +1,61 @@ +# %% + +import pandas as pd + +from ribasim_nl import CloudStorage, Model, NetworkValidator + +cloud = CloudStorage() + +authority = "DrentsOverijsselseDelta" +short_name = "dod" + +ribasim_toml = cloud.joinpath(authority, "modellen", f"{authority}_fix_model_network", f"{short_name}.toml") +database_gpkg = ribasim_toml.with_name("database.gpkg") + +model = Model.read(ribasim_toml) +network_validator = NetworkValidator(model) + +verwerkt_dir = cloud.joinpath(authority, "verwerkt") +verwerkt_dir.mkdir(exist_ok=True) + +modelfouten_gpkg = cloud.joinpath(authority, "verwerkt", "modelfouten.gpkg") + +# %% verwijderen duplicated edges + +duplicated_edges = len(model.edge.df[model.edge.df.duplicated()]) +model.edge.df.drop_duplicates(inplace=True) + +model.fix_unassigned_basin_area() + +# %% wegschrijven fouten + +# niet-bestaande fouten +mask = model.edge.df.to_node_id.isin(model.node_table().df.index) & model.edge.df.from_node_id.isin( + model.node_table().df.index +) + +edge_mist_node_df = model.edge.df[~mask] +model.edge.df = model.edge.df[mask] + +mask = model.edge.df.geometry.length == 0 +model.edge.df[mask].centroid.to_file(modelfouten_gpkg, layer="edge_zonder_lengte") +model.edge.df = model.edge.df[~mask] + +# niet-gekoppelde areas +model.basin.area.df[~model.basin.area.df.node_id.isin(model.basin.node.df.index)].to_file( + modelfouten_gpkg, layer="area_niet_een_basin" +) + +model.basin.node.df[~model.basin.node.df.index.isin(model.basin.area.df.node_id)].to_file( + modelfouten_gpkg, layer="basin_zonder_area" +) + +# ontbrekende basins +network_validator.node_invalid_connectivity().to_file(modelfouten_gpkg, layer="node_mist") +pd.concat([network_validator.edge_incorrect_connectivity(), edge_mist_node_df]).to_file( + modelfouten_gpkg, layer="ege_mist_node" +) + +# nodes met verkeerde richting + +model.invalid_topology_at_node().to_file(modelfouten_gpkg, layer="node_met_verkeerde_instroom_uitstroom_egde") diff --git a/notebooks/drents_overijsselse_delta/01_fix_model_network.py b/notebooks/drents_overijsselse_delta/01_fix_model_network.py index 7779a3cf..77e1420c 100644 --- a/notebooks/drents_overijsselse_delta/01_fix_model_network.py +++ b/notebooks/drents_overijsselse_delta/01_fix_model_network.py @@ -160,7 +160,7 @@ # Aansluiten NW boezem op Fryslan # basin /area 1681 op te knippen nabij basin 1717 (rode lijn) -model.split_basin(split_line_gdf.at[14, "geometry"]) +model.split_basin(geometry=split_line_gdf.at[14, "geometry"]) model.basin.area.df = model.basin.area.df[model.basin.area.df.node_id != 1717] # basin 1682 te veranderen in een LevelBoundary @@ -398,17 +398,22 @@ # %% +model.explode_basin_area() # all multipolygons to singles # update from layers actions = [ + "remove_basin_area", + "split_basin", + "merge_basins", + "update_node", + "add_basin_area", "add_basin", + "update_basin_area", "redirect_edge", - "merge_basins", "reverse_edge", - "update_node", "deactivate_node", - "remove_node", "move_node", + "remove_node", ] actions = [i for i in actions if i in gpd.list_layers(model_edits_path).name.to_list()] @@ -418,6 +423,8 @@ method = getattr(model, action) keywords = inspect.getfullargspec(method).args df = gpd.read_file(model_edits_path, layer=action, fid_as_index=True) + if "order" in df.columns: + df.sort_values("order", inplace=True) for row in df.itertuples(): # filter kwargs by keywords kwargs = {k: v for k, v in row._asdict().items() if k in keywords} diff --git a/notebooks/noorderzijlvest/01b_fix_basin_area.py b/notebooks/noorderzijlvest/01b_fix_basin_area.py index fc09a2d8..a958908f 100644 --- a/notebooks/noorderzijlvest/01b_fix_basin_area.py +++ b/notebooks/noorderzijlvest/01b_fix_basin_area.py @@ -1,9 +1,12 @@ # %% Import Libraries and Initialize Variables +import inspect import geopandas as gpd import pandas as pd +from networkx import all_shortest_paths, shortest_path +from shapely.geometry import MultiLineString -from ribasim_nl import CloudStorage, Model +from ribasim_nl import CloudStorage, Model, Network # Initialize cloud storage and set authority/model parameters cloud_storage = CloudStorage() @@ -14,19 +17,84 @@ ribasim_model_dir = cloud_storage.joinpath(authority_name, "modellen", f"{authority_name}_fix_model_network") ribasim_model_path = ribasim_model_dir / f"{model_short_name}.toml" model = Model.read(ribasim_model_path) -node_df = model.node_table().df + # read hydrologische eenheden he_df = gpd.read_file( - cloud_storage.joinpath(authority_name, "verwerkt", "1_ontvangen_data", "20241113", "HydrologischeEenheden_v45.shp") + cloud_storage.joinpath(authority_name, "verwerkt", "1_ontvangen_data", "20241113", "HydrologischeEenheden_v45.shp"), ) he_df.loc[:, "node_id"] = pd.Series() he_snap_df = gpd.read_file( - cloud_storage.joinpath(authority_name, "verwerkt", "1_ontvangen_data", "20241113", "HE_v45_snappingpoints.shp") + cloud_storage.joinpath(authority_name, "verwerkt", "1_ontvangen_data", "20241113", "HE_v45_snappingpoints.shp"), +) + +lines_gdf = gpd.read_file( + cloud_storage.joinpath( + authority_name, "verwerkt", "5_D_HYDRO_export", "hydroobjecten", "Noorderzijlvest_hydroobjecten.shp" + ), ) +lines_gdf.crs = 28992 +network = Network(lines_gdf) +network.to_file(cloud_storage.joinpath(authority_name, "verwerkt", "network.gpkg")) + +# %% add snap_point to he_df + +# Purpose is to define 1 point for every 1 he polygon +# We can use that later in finding he over the network + +# drop duplicated kwk points +he_snap_df = he_snap_df[~(he_snap_df["Kunstwerk"].duplicated() & he_snap_df["Kunstwerk"].notna())] + +# strip GPG from Kunstwerk column +mask = he_snap_df["Kunstwerk"].notna() & he_snap_df["Kunstwerk"].str.startswith("GPG") +he_snap_df.loc[mask, "Kunstwerk"] = he_snap_df[mask].Kunstwerk.str[3:] + +# we define 1 outlet for every he +he_outlet_df = gpd.GeoDataFrame(geometry=gpd.GeoSeries(), crs=he_snap_df.crs, index=he_df.index) + +mask = he_df["KWKuit"].notna() +he_outlet_df.loc[mask, "geometry"] = he_df[mask]["KWKuit"].apply( + lambda x: he_snap_df.set_index("Kunstwerk").at[x, "geometry"] +) +he_outlet_df.loc[~mask, "geometry"] = he_df[~mask]["HEIDENT"].apply( + lambda x: he_snap_df.set_index("HEIDENT").at[x, "geometry"] +) + +he_outlet_df.loc[:, "HEIDENT"] = he_df["HEIDENT"] +he_outlet_df.set_index("HEIDENT", inplace=True) +he_df.set_index("HEIDENT", inplace=True) +he_outlet_df.to_file(cloud_storage.joinpath(authority_name, "verwerkt", "HydrologischeEenheden_v45_outlets.gpkg")) -# %% assign basin-id if we can find KWKuit in node_table +# %% Edit network + +# We modify the network: +# merge basins in Lauwersmeer + +# Load node edit data +model_edits_url = cloud_storage.joinurl(authority_name, "verwerkt", "model_edits.gpkg") +model_edits_path = cloud_storage.joinpath(authority_name, "verwerkt", "model_edits.gpkg") +if not model_edits_path.exists(): + cloud_storage.download_file(model_edits_url) + + +for action in gpd.list_layers(model_edits_path).name.to_list(): + print(action) + # get method and args + method = getattr(model, action) + keywords = inspect.getfullargspec(method).args + df = gpd.read_file(model_edits_path, layer=action, fid_as_index=True) + for row in df.itertuples(): + # filter kwargs by keywords + kwargs = {k: v for k, v in row._asdict().items() if k in keywords} + method(**kwargs) + +# %% assign Basin / Area using KWKuit + +node_df = model.node_table().df + + +# we find Basin area if we kan find KWKuit in the model def find_basin_id(kwk_code): kwk_node_id = node_df[node_df.node_type != "Basin"].reset_index().set_index("name").at[kwk_code, "node_id"] basin_node_id = model.upstream_node_id(kwk_node_id) @@ -42,24 +110,108 @@ def find_basin_id(kwk_code): he_df.loc[mask, "node_id"] = he_df[mask]["KWKuit"].apply(lambda x: find_basin_id(x)) -# %% assign basin-id if we can find KWKuit in node_table -point, df = next( - i - for i in iter(he_snap_df[he_snap_df.duplicated("geometry", keep=False)].groupby("geometry")) - if i[1].Kunstwerk.notna().any() -) +# %% find he on network within basin + + +# We find all hydrologische eenheden using outlets between basin and it's connector-nodes +def get_network_node(point): + node = network.move_node(point, max_distance=1, align_distance=10) + if node is None: + node = network.add_node(point, max_distance=1) + if node is None: + node = network.nodes.distance(point).idxmin() + return node + + +for node_id in model.basin.node.df.index: + print(node_id) + # get basin_node_id + network_basin_node = get_network_node(node_df.at[node_id, "geometry"]) + network._graph.nodes[network_basin_node]["node_id"] = node_id + + upstream_node_ids = model.upstream_node_id(node_id) + if isinstance(upstream_node_ids, pd.Series): + upstream_node_ids = upstream_node_ids.to_list() + else: + upstream_node_ids = [upstream_node_ids] + + downstream_node_ids = model.downstream_node_id(node_id) + if isinstance(downstream_node_ids, pd.Series): + downstream_node_ids = downstream_node_ids.to_list() + else: + downstream_node_ids = [downstream_node_ids] + + upstream_nodes = [get_network_node(node_df.at[i, "geometry"]) for i in upstream_node_ids if i is not None] + downstream_nodes = [get_network_node(node_df.at[i, "geometry"]) for i in downstream_node_ids] + + # empty list of LineStrings + data = [] + + # draw edges from upstream nodes + for idx, network_node in enumerate(upstream_nodes): + all_paths = list(all_shortest_paths(network.graph_undirected, source=network_node, target=network_basin_node)) + if len(all_paths) > 1: + all_nodes = [i for i in upstream_nodes + downstream_nodes if i != network_node] + all_paths = [i for i in all_paths if not any(_node_id in all_nodes for _node_id in i)] + if len(all_paths) != 1: + all_paths = [shortest_path(network.graph_undirected, source=network_node, target=network_basin_node)] + else: + edge = network.path_to_line(all_paths[0]) + if edge.length > 0: + data += [edge] + + mask = (model.edge.df["from_node_id"] == upstream_node_ids[idx]) & ( + model.edge.df["to_node_id"] == node_id + ) + model.edge.df.loc[mask, ["geometry"]] = edge + + # draw edges to downstream nodes + for idx, network_node in enumerate(downstream_nodes): + all_paths = list(all_shortest_paths(network.graph_undirected, target=network_node, source=network_basin_node)) + if len(all_paths) > 1: + all_nodes = [i for i in upstream_nodes + downstream_nodes if i != network_node] + all_paths = [i for i in all_paths if not any(_node_id in all_nodes for _node_id in i)] + if len(all_paths) != 1: + all_paths = [shortest_path(network.graph_undirected, target=network_node, source=network_basin_node)] + else: + edge = network.path_to_line(all_paths[0]) + if edge.length > 0: + data += [edge] + + mask = (model.edge.df["to_node_id"] == downstream_node_ids[idx]) & ( + model.edge.df["from_node_id"] == node_id + ) + model.edge.df.loc[mask, ["geometry"]] = edge + + mask = he_df.node_id.isna() & (he_outlet_df.distance(MultiLineString(data)) < 1) + he_df.loc[mask, ["node_id"]] = node_id + +# %% add last missings + +# We add last missing hydrologische eenheden on downstream basin +for row in he_df[he_df["node_id"].isna()].itertuples(): + # row = next(i for i in he_df.itertuples() if i.Index == "GFE04712") + print(row.Index) + point = he_outlet_df.at[row.Index, "geometry"] + + network_node = get_network_node(point) + + basin_node_id = network.find_downstream(network_node, attribute="node_id") + he_df.loc[row.Index, ["node_id"]] = basin_node_id + + +# %% renew model Basin Area + +# based on actions above we update the Basin / Area table: +# basins we found on KWKuit defined in the he-table -# %% dissolve geometry and take min summer target level as streefpeil data = [] for node_id, df in he_df[he_df["node_id"].notna()].groupby("node_id"): geometry = df.union_all() streefpeil = df["OPVAFWZP"].min() data += [{"node_id": node_id, "meta_streefpeil": streefpeil, "geometry": geometry}] - -# %% -# add new basins to model df = gpd.GeoDataFrame(data, crs=model.crs) df.loc[:, "geometry"] = df.buffer(0.1).buffer(-0.1) df.index.name = "fid" diff --git a/notebooks/rijkswaterstaat/8c_patch_model_202470 b/notebooks/rijkswaterstaat/8c_patch_model_202470 new file mode 100644 index 00000000..0c8431a1 --- /dev/null +++ b/notebooks/rijkswaterstaat/8c_patch_model_202470 @@ -0,0 +1,18 @@ +# %% +from ribasim_nl import CloudStorage, Model + +cloud = CloudStorage() + +ribasim_toml = cloud.joinpath("Rijkswaterstaat", "modellen", "hws_2024_7_0", "hws.toml") +model = Model.read(ribasim_toml) + +# %%update stuw driel +node_id = 8413 +min_level = model.upstream_profile(node_id).level.min() + 0.1 +mask = model.tabulated_rating_curve.static.df.node_id == node_id +mask = mask & (model.tabulated_rating_curve.static.df.level < min_level) +model.tabulated_rating_curve.static.df.loc[mask, ["level"]] = min_level + + +model.write(cloud.joinpath("Rijkswaterstaat", "modellen", "hws", "hws.toml")) +# %% diff --git a/notebooks/samenvoegen_modellen.py b/notebooks/samenvoegen_modellen.py index 18134265..424c1c2b 100644 --- a/notebooks/samenvoegen_modellen.py +++ b/notebooks/samenvoegen_modellen.py @@ -1,12 +1,12 @@ # %% from datetime import datetime -import numpy as np import ribasim -from ribasim_nl import CloudStorage +from ribasim_nl import CloudStorage, Model from ribasim_nl.case_conversions import pascal_to_snake_case from ribasim_nl.concat import concat +from ribasim_nl.reset_static_tables import reset_static_tables # %% cloud = CloudStorage() @@ -22,138 +22,118 @@ download_latest_model = True upload_model = False +RESET_TABLES = ["AaenMaas", "HunzeenAas", "DrentsOverijsselseDelta", "Vechtstromen"] + models = [ { "authority": "Rijkswaterstaat", "model": "hws", "find_toml": False, - "zoom_level": 0, }, { "authority": "AmstelGooienVecht", - "model": "AmstelGooienVecht_boezemmodel", + "model": "AmstelGooienVecht_parametrized", "find_toml": True, - "zoom_level": 3, }, { "authority": "Delfland", - "model": "Delfland_boezemmodel", + "model": "Delfland_parametrized", "find_toml": True, - "zoom_level": 3, }, { "authority": "HollandseDelta", - "model": "HollandseDelta_boezemmodel", + "model": "HollandseDelta_parametrized", "find_toml": True, - "zoom_level": 3, }, { "authority": "HollandsNoorderkwartier", - "model": "HollandsNoorderkwartier_boezemmodel", + "model": "HollandsNoorderkwartier_parametrized", "find_toml": True, - "zoom_level": 3, }, { "authority": "Rijnland", - "model": "Rijnland_boezemmodel", + "model": "Rijnland_parametrized", "find_toml": True, - "zoom_level": 3, }, { "authority": "Rivierenland", - "model": "Rivierenland_boezemmodel", + "model": "Rivierenland_parametrized", "find_toml": True, - "zoom_level": 3, }, { "authority": "Scheldestromen", - "model": "Scheldestromen_boezemmodel", + "model": "Scheldestromen_parametrized", "find_toml": True, - "zoom_level": 3, }, { "authority": "SchielandendeKrimpenerwaard", - "model": "SchielandendeKrimpenerwaard_boezemmodel", + "model": "SchielandendeKrimpenerwaard_parametrized", "find_toml": True, - "zoom_level": 3, }, { "authority": "WetterskipFryslan", - "model": "WetterskipFryslan_boezemmodel", + "model": "WetterskipFryslan_parametrized", "find_toml": True, - "zoom_level": 3, }, { "authority": "Zuiderzeeland", - "model": "Zuiderzeeland_boezemmodel", + "model": "Zuiderzeeland_parametrized", "find_toml": True, - "zoom_level": 3, }, { "authority": "AaenMaas", "model": "AaenMaas", "find_toml": True, - "zoom_level": 3, }, { "authority": "BrabantseDelta", "model": "BrabantseDelta", "find_toml": True, - "zoom_level": 3, }, { "authority": "DeDommel", "model": "DeDommel", "find_toml": True, - "zoom_level": 3, }, { "authority": "DrentsOverijsselseDelta", "model": "DrentsOverijsselseDelta", "find_toml": True, - "zoom_level": 3, }, { "authority": "HunzeenAas", "model": "HunzeenAas", "find_toml": True, - "zoom_level": 3, }, { "authority": "Limburg", "model": "Limburg", "find_toml": True, - "zoom_level": 3, }, { "authority": "Noorderzijlvest", "model": "Noorderzijlvest", "find_toml": True, - "zoom_level": 3, }, { "authority": "RijnenIJssel", "model": "RijnenIJssel", "find_toml": True, - "zoom_level": 3, }, { "authority": "StichtseRijnlanden", "model": "StichtseRijnlanden", "find_toml": True, - "zoom_level": 3, }, { "authority": "ValleienVeluwe", "model": "ValleienVeluwe", "find_toml": True, - "zoom_level": 3, }, { "authority": "Vechtstromen", "model": "Vechtstromen", "find_toml": True, - "zoom_level": 3, }, ] @@ -203,68 +183,25 @@ def get_model_path(model, model_version): model_path = model_path.joinpath(f"{model["model"]}.toml") # read model - ribasim_model = ribasim.Model.read(model_path) + ribasim_model = Model.read(model_path) - # zoom_level to model + # add meta_waterbeheerder for node_type in ribasim_model.node_table().df.node_type.unique(): ribasim_node = getattr(ribasim_model, pascal_to_snake_case(node_type)) - ribasim_node.node.df.loc[:, "meta_zoom_level"] = model["zoom_level"] ribasim_node.node.df.loc[:, "meta_waterbeheerder"] = model["authority"] - ribasim_model.edge.df.loc[:, "meta_zoom_level"] = model["zoom_level"] - ribasim_model.edge.df.loc[:, "meta_waterbeheerder"] = model["authority"] - if idx == 0: lhm_model = ribasim_model else: - if "meta_categorie" in ribasim_model.edge.df.columns: - # tric to use higher zoom-level for primary waters - - ribasim_model.edge.df.loc[ - ribasim_model.edge.df["meta_categorie"] == "hoofdwater", - "meta_zoom_level", - ] = 2 - df = ribasim_model.edge.df[ribasim_model.edge.df["meta_zoom_level"] == 2] - if not df.empty: - node_ids = set( - np.concatenate( - ribasim_model.edge.df[ribasim_model.edge.df["meta_zoom_level"] == 2][ - ["from_node_id", "to_node_id"] - ].to_numpy() - ) - ) - # zoom_level to model - for node_type in ribasim_model.node_table().df.node_type.unique(): - ribasim_node = getattr(ribasim_model, pascal_to_snake_case(node_type)) - ribasim_node.node.df.loc[ - ribasim_node.node.df["node_id"].isin(node_ids), - "meta_zoom_level", - ] = 2 - - cols = [i for i in lhm_model.edge.df.columns if i != "meta_index"] - lhm_model.edge.df = lhm_model.edge.df[cols] - try: - lhm_model = concat([lhm_model, ribasim_model]) - except KeyError: - print(f"model {model['authority']} is waarschijnlijk niet aan de haak") - pass + # TODO: make sure this isn't needed next round! + if model["authority"] in RESET_TABLES: + ribasim_model.remove_unassigned_basin_area() + ribasim_model = reset_static_tables(ribasim_model) + lhm_model = concat([lhm_model, ribasim_model]) readme += f""" - **{model["authority"]}**: {model["model"]} ({model_version.version})""" - -# %% color -# color_cycle = itertools.cycle(Category10[10]) -# lhm_model.basin.area.df.loc[ -# lhm_model.basin.area.df["meta_streefpeil"] == "Onbekend streefpeil", -# "meta_streefpeil", -# ] = None -# lhm_model.basin.area.df.loc[:, "meta_streefpeil"] = lhm_model.basin.area.df[ -# "meta_streefpeil" -# ].astype(float) +**{model["authority"]}**: {model["model"]} ({model_version.version})""" -# lhm_model.basin.area.df.loc[:, "meta_color"] = [ -# next(color_cycle) for _ in range(len(lhm_model.basin.area.df)) -# ] # %% print("write lhm model") diff --git a/notebooks/upload_feedback_formulieren.py b/notebooks/upload_feedback_formulieren.py new file mode 100644 index 00000000..a25af365 --- /dev/null +++ b/notebooks/upload_feedback_formulieren.py @@ -0,0 +1,72 @@ +# %% +from datetime import date + +import xlwings as xw + +from ribasim_nl import CloudStorage, Model + +cloud = CloudStorage() + +WATER_AUTHORITIES = [ + "AaenMaas", + "BrabantseDelta", + "DeDommel", + "DrentsOverijsselseDelta", + "HunzeenAas", + "Limburg", + "Noorderzijlvest", + "RijnenIJssel", + "StichtseRijnlanden", + "ValleienVeluwe", + "Vechtstromen", +] + +FEEDBACK_XLS = cloud.joinpath("Basisgegevens", "feedbackformulier", "Feedback Formulier.xlsx") + +for authority in WATER_AUTHORITIES: + print(authority) + model_dir = cloud.joinpath(authority, "modellen", authority) + toml_file = next(model_dir.glob("*.toml")) + model = Model.read(toml_file) + + # load feedback xlsx + app = xw.App(visible=False) + workbook = xw.Book(FEEDBACK_XLS) + + # add authority to Feedback_Formulier + sheet = workbook.sheets["Feedback_Formulier"] + sheet.range("B1").value = authority + sheet.range("B3").value = date.today() + + # add Node_Data + sheet = workbook.sheets["Node_Data"] + df = model.node_table().df.reset_index()[["node_id", "name", "node_type", "subnetwork_id"]] + sheet.range("A2").value = df.to_numpy() + + # add Edge_Data + sheet = workbook.sheets["Edge_Data"] + df = model.edge.df + df.loc[:, "from_node_type"] = model.edge_from_node_type + df.loc[:, "to_node_type"] = model.edge_to_node_type + df = model.edge.df.reset_index()[ + [ + "edge_id", + "name", + "from_node_type", + "from_node_id", + "to_node_type", + "to_node_id", + "edge_type", + "subnetwork_id", + ] + ] + sheet.range("A2").value = df.to_numpy() + + # write copy for authority + excel_file = cloud.joinpath(authority, "verwerkt", FEEDBACK_XLS.name) + workbook.save(excel_file) + workbook.close() + app.quit() + + # upload file to cloud + cloud.upload_file(excel_file) diff --git a/notebooks/uploaden_modellen.py b/notebooks/uploaden_modellen.py new file mode 100644 index 00000000..711c9be3 --- /dev/null +++ b/notebooks/uploaden_modellen.py @@ -0,0 +1,40 @@ +# %% +from ribasim_nl import CloudStorage, Model + +cloud = CloudStorage() + +FIND_POST_FIXES = ["bergend", "fix_model_area", "fix_model_network"] +SELECTION = [] + + +def get_model_dir(authority, post_fix): + return cloud.joinpath(authority, "modellen", f"{authority}_{post_fix}") + + +if len(SELECTION) == 0: + authorities = cloud.water_authorities +else: + authorities = SELECTION + +for authority in authorities: + # find model directory + model_dir = next( + ( + get_model_dir(authority, post_fix) + for post_fix in FIND_POST_FIXES + if get_model_dir(authority, post_fix).exists() + ), + None, + ) + if model_dir is not None: + print(f"uploading copy for {authority}") + # read model + toml_file = next(model_dir.glob("*.toml")) + model = Model.read(toml_file) + + # write a local copy in root + toml_file = toml_file.parents[1].joinpath(authority, toml_file.name) + model.write(toml_file) + + # create version and upload + cloud.upload_model(authority=authority, model=authority) diff --git a/pixi.lock b/pixi.lock index e04dd7a5..c7ed3833 100644 --- a/pixi.lock +++ b/pixi.lock @@ -442,6 +442,7 @@ environments: - conda: https://conda.anaconda.org/conda-forge/linux-64/xcb-util-wm-0.4.2-hb711507_0.conda - conda: https://conda.anaconda.org/conda-forge/linux-64/xerces-c-3.2.5-h988505b_2.conda - conda: https://conda.anaconda.org/conda-forge/linux-64/xkeyboard-config-2.43-hb9d3cd8_0.conda + - conda: https://conda.anaconda.org/conda-forge/linux-64/xlwings-0.33.4-py312h7900ff3_0.conda - conda: https://conda.anaconda.org/conda-forge/linux-64/xorg-libice-1.1.1-hb9d3cd8_1.conda - conda: https://conda.anaconda.org/conda-forge/linux-64/xorg-libsm-1.2.4-he73a12e_1.conda - conda: https://conda.anaconda.org/conda-forge/linux-64/xorg-libx11-1.8.10-h4f16b4b_0.conda @@ -481,6 +482,7 @@ environments: - conda: https://conda.anaconda.org/conda-forge/osx-64/aom-3.9.1-hf036a51_0.conda - conda: https://conda.anaconda.org/conda-forge/noarch/appdirs-1.4.4-pyh9f0ad1d_0.tar.bz2 - conda: https://conda.anaconda.org/conda-forge/noarch/appnope-0.1.4-pyhd8ed1ab_0.conda + - conda: https://conda.anaconda.org/conda-forge/osx-64/appscript-1.3.0-py312hde3c1db_0.conda - conda: https://conda.anaconda.org/conda-forge/noarch/argon2-cffi-23.1.0-pyhd8ed1ab_0.conda - conda: https://conda.anaconda.org/conda-forge/osx-64/argon2-cffi-bindings-21.2.0-py312hb553811_5.conda - conda: https://conda.anaconda.org/conda-forge/noarch/arrow-1.3.0-pyhd8ed1ab_0.conda @@ -704,6 +706,7 @@ environments: - conda: https://conda.anaconda.org/conda-forge/osx-64/llvm-openmp-19.1.4-ha54dae1_0.conda - conda: https://conda.anaconda.org/conda-forge/osx-64/llvmlite-0.43.0-py312hcc8fd36_1.conda - conda: https://conda.anaconda.org/conda-forge/noarch/locket-1.0.0-pyhd8ed1ab_0.tar.bz2 + - conda: https://conda.anaconda.org/conda-forge/osx-64/lxml-5.3.0-py312h4feaf87_2.conda - conda: https://conda.anaconda.org/conda-forge/osx-64/lz4-4.3.3-py312h83408cd_1.conda - conda: https://conda.anaconda.org/conda-forge/osx-64/lz4-c-1.9.4-hf0c8a7f_0.conda - conda: https://conda.anaconda.org/conda-forge/osx-64/lzo-2.10-h10d778d_1001.conda @@ -872,6 +875,7 @@ environments: - conda: https://conda.anaconda.org/conda-forge/osx-64/x265-3.5-hbb4e6a2_3.tar.bz2 - conda: https://conda.anaconda.org/conda-forge/noarch/xarray-2024.11.0-pyhd8ed1ab_0.conda - conda: https://conda.anaconda.org/conda-forge/osx-64/xerces-c-3.2.5-h197e74d_2.conda + - conda: https://conda.anaconda.org/conda-forge/osx-64/xlwings-0.33.4-py312hb401068_0.conda - conda: https://conda.anaconda.org/conda-forge/osx-64/xorg-libxau-1.0.11-h00291cd_1.conda - conda: https://conda.anaconda.org/conda-forge/osx-64/xorg-libxdmcp-1.1.5-h00291cd_0.conda - conda: https://conda.anaconda.org/conda-forge/noarch/xugrid-0.12.1-pyhd8ed1ab_0.conda @@ -1285,6 +1289,7 @@ environments: - conda: https://conda.anaconda.org/conda-forge/win-64/x265-3.5-h2d74725_3.tar.bz2 - conda: https://conda.anaconda.org/conda-forge/noarch/xarray-2024.11.0-pyhd8ed1ab_0.conda - conda: https://conda.anaconda.org/conda-forge/win-64/xerces-c-3.2.5-he0c23c2_2.conda + - conda: https://conda.anaconda.org/conda-forge/win-64/xlwings-0.33.4-py312h2e8e312_0.conda - conda: https://conda.anaconda.org/conda-forge/win-64/xorg-libxau-1.0.11-h0e40799_1.conda - conda: https://conda.anaconda.org/conda-forge/win-64/xorg-libxdmcp-1.1.5-h0e40799_0.conda - conda: https://conda.anaconda.org/conda-forge/noarch/xugrid-0.12.1-pyhd8ed1ab_0.conda @@ -1744,6 +1749,7 @@ environments: - conda: https://conda.anaconda.org/conda-forge/linux-64/xcb-util-wm-0.4.2-hb711507_0.conda - conda: https://conda.anaconda.org/conda-forge/linux-64/xerces-c-3.2.5-h988505b_2.conda - conda: https://conda.anaconda.org/conda-forge/linux-64/xkeyboard-config-2.43-hb9d3cd8_0.conda + - conda: https://conda.anaconda.org/conda-forge/linux-64/xlwings-0.33.4-py312h7900ff3_0.conda - conda: https://conda.anaconda.org/conda-forge/linux-64/xorg-libice-1.1.1-hb9d3cd8_1.conda - conda: https://conda.anaconda.org/conda-forge/linux-64/xorg-libsm-1.2.4-he73a12e_1.conda - conda: https://conda.anaconda.org/conda-forge/linux-64/xorg-libx11-1.8.10-h4f16b4b_0.conda @@ -1784,6 +1790,7 @@ environments: - conda: https://conda.anaconda.org/conda-forge/osx-64/aom-3.9.1-hf036a51_0.conda - conda: https://conda.anaconda.org/conda-forge/noarch/appdirs-1.4.4-pyh9f0ad1d_0.tar.bz2 - conda: https://conda.anaconda.org/conda-forge/noarch/appnope-0.1.4-pyhd8ed1ab_0.conda + - conda: https://conda.anaconda.org/conda-forge/osx-64/appscript-1.3.0-py312hde3c1db_0.conda - conda: https://conda.anaconda.org/conda-forge/noarch/argon2-cffi-23.1.0-pyhd8ed1ab_0.conda - conda: https://conda.anaconda.org/conda-forge/osx-64/argon2-cffi-bindings-21.2.0-py312hb553811_5.conda - conda: https://conda.anaconda.org/conda-forge/noarch/arrow-1.3.0-pyhd8ed1ab_0.conda @@ -2007,6 +2014,7 @@ environments: - conda: https://conda.anaconda.org/conda-forge/osx-64/llvm-openmp-19.1.4-ha54dae1_0.conda - conda: https://conda.anaconda.org/conda-forge/osx-64/llvmlite-0.43.0-py312hcc8fd36_1.conda - conda: https://conda.anaconda.org/conda-forge/noarch/locket-1.0.0-pyhd8ed1ab_0.tar.bz2 + - conda: https://conda.anaconda.org/conda-forge/osx-64/lxml-5.3.0-py312h4feaf87_2.conda - conda: https://conda.anaconda.org/conda-forge/osx-64/lz4-4.3.3-py312h83408cd_1.conda - conda: https://conda.anaconda.org/conda-forge/osx-64/lz4-c-1.9.4-hf0c8a7f_0.conda - conda: https://conda.anaconda.org/conda-forge/osx-64/lzo-2.10-h10d778d_1001.conda @@ -2174,6 +2182,7 @@ environments: - conda: https://conda.anaconda.org/conda-forge/osx-64/x265-3.5-hbb4e6a2_3.tar.bz2 - conda: https://conda.anaconda.org/conda-forge/noarch/xarray-2024.11.0-pyhd8ed1ab_0.conda - conda: https://conda.anaconda.org/conda-forge/osx-64/xerces-c-3.2.5-h197e74d_2.conda + - conda: https://conda.anaconda.org/conda-forge/osx-64/xlwings-0.33.4-py312hb401068_0.conda - conda: https://conda.anaconda.org/conda-forge/osx-64/xorg-libxau-1.0.11-h00291cd_1.conda - conda: https://conda.anaconda.org/conda-forge/osx-64/xorg-libxdmcp-1.1.5-h00291cd_0.conda - conda: https://conda.anaconda.org/conda-forge/noarch/xugrid-0.12.1-pyhd8ed1ab_0.conda @@ -2587,6 +2596,7 @@ environments: - conda: https://conda.anaconda.org/conda-forge/win-64/x265-3.5-h2d74725_3.tar.bz2 - conda: https://conda.anaconda.org/conda-forge/noarch/xarray-2024.11.0-pyhd8ed1ab_0.conda - conda: https://conda.anaconda.org/conda-forge/win-64/xerces-c-3.2.5-he0c23c2_2.conda + - conda: https://conda.anaconda.org/conda-forge/win-64/xlwings-0.33.4-py312h2e8e312_0.conda - conda: https://conda.anaconda.org/conda-forge/win-64/xorg-libxau-1.0.11-h0e40799_1.conda - conda: https://conda.anaconda.org/conda-forge/win-64/xorg-libxdmcp-1.1.5-h0e40799_0.conda - conda: https://conda.anaconda.org/conda-forge/noarch/xugrid-0.12.1-pyhd8ed1ab_0.conda @@ -2828,6 +2838,24 @@ packages: - pkg:pypi/appnope?source=hash-mapping size: 10241 timestamp: 1707233195627 +- kind: conda + name: appscript + version: 1.3.0 + build: py312hde3c1db_0 + subdir: osx-64 + url: https://conda.anaconda.org/conda-forge/osx-64/appscript-1.3.0-py312hde3c1db_0.conda + sha256: 4b730025a22510de399d412c86b63dd40a61e60e07227b083f9370c2d1243653 + md5: 16b1cec30fe32451f1a2c72976d40ab2 + depends: + - __osx >=10.9 + - lxml >=4.7.1 + - python >=3.12,<3.13.0a0 + - python_abi 3.12.* *_cp312 + license: Public Domain + purls: + - pkg:pypi/appscript?source=hash-mapping + size: 166574 + timestamp: 1729580107578 - kind: conda name: argon2-cffi version: 23.1.0 @@ -13176,6 +13204,27 @@ packages: - pkg:pypi/locket?source=hash-mapping size: 8250 timestamp: 1650660473123 +- kind: conda + name: lxml + version: 5.3.0 + build: py312h4feaf87_2 + build_number: 2 + subdir: osx-64 + url: https://conda.anaconda.org/conda-forge/osx-64/lxml-5.3.0-py312h4feaf87_2.conda + sha256: 4e6c79df93b27024435763d6216cc11d423b28bd01c676246f86a69c1c691846 + md5: d1d2583ce9c06e63bf90a5c523e30ab0 + depends: + - __osx >=10.13 + - libxml2 >=2.12.7,<3.0a0 + - libxslt >=1.1.39,<2.0a0 + - libzlib >=1.3.1,<2.0a0 + - python >=3.12,<3.13.0a0 + - python_abi 3.12.* *_cp312 + license: BSD-3-Clause and MIT-CMU + purls: + - pkg:pypi/lxml?source=hash-mapping + size: 1194243 + timestamp: 1730194714394 - kind: conda name: lz4 version: 4.3.3 @@ -14103,7 +14152,7 @@ packages: sha256: df5d4365b724cf81b8c6a7312509d0c22386097011ad1abe274afd5e9d3bbc5f requires_dist: - numpy>=1.24 ; extra == 'default' - - scipy>=1.10,!=1.11.0,!=1.11.1 ; extra == 'default' + - scipy!=1.11.0,!=1.11.1,>=1.10 ; extra == 'default' - matplotlib>=3.7 ; extra == 'default' - pandas>=2.0 ; extra == 'default' - changelist==0.5 ; extra == 'developer' @@ -19880,6 +19929,60 @@ packages: purls: [] size: 389475 timestamp: 1727840188958 +- kind: conda + name: xlwings + version: 0.33.4 + build: py312h2e8e312_0 + subdir: win-64 + url: https://conda.anaconda.org/conda-forge/win-64/xlwings-0.33.4-py312h2e8e312_0.conda + sha256: 47842133ecd7a17d671f5a06bb4724ab8cfd02f88b5f606a17eea19ef3a9d327 + md5: 78e4937784b3033a92805744f68a5ad6 + depends: + - python >=3.12,<3.13.0a0 + - python_abi 3.12.* *_cp312 + - pywin32 + license: BSD-3-Clause + license_family: BSD + purls: + - pkg:pypi/xlwings?source=hash-mapping + size: 1268684 + timestamp: 1732365098367 +- kind: conda + name: xlwings + version: 0.33.4 + build: py312h7900ff3_0 + subdir: linux-64 + url: https://conda.anaconda.org/conda-forge/linux-64/xlwings-0.33.4-py312h7900ff3_0.conda + sha256: 0235c6d5c8f8f65bb0842e09472df67eea1484beef97369ff971ec690767a902 + md5: a2792e91957f6f9a62b14d7e75ee77ff + depends: + - python >=3.12,<3.13.0a0 + - python_abi 3.12.* *_cp312 + license: BSD-3-Clause + license_family: BSD + purls: + - pkg:pypi/xlwings?source=hash-mapping + size: 939450 + timestamp: 1732364869337 +- kind: conda + name: xlwings + version: 0.33.4 + build: py312hb401068_0 + subdir: osx-64 + url: https://conda.anaconda.org/conda-forge/osx-64/xlwings-0.33.4-py312hb401068_0.conda + sha256: f8e86a0ba6ebedf2168caace974f5e1abbc9804887fd1d3ca7c3692ee6798390 + md5: 7f1cae62da4fceb51131be562cf7d57d + depends: + - appscript >=1.0.1 + - psutil >=2.0 + - python >=3.12,<3.13.0a0 + - python_abi 3.12.* *_cp312 + license: BSD-3-Clause + license_family: BSD + purls: + - pkg:pypi/xlwings?source=hash-mapping + size: 942090 + timestamp: 1732364881987 - kind: conda name: xorg-libice version: 1.1.1 diff --git a/pixi.toml b/pixi.toml index 7bc700a1..667c12e8 100644 --- a/pixi.toml +++ b/pixi.toml @@ -66,6 +66,7 @@ shapely = ">=2" tomli = "*" tomli-w = "*" types-requests = "*" +xlwings = "*" xugrid = "*" [feature.common.pypi-dependencies] diff --git a/src/ribasim_nl/ribasim_nl/model.py b/src/ribasim_nl/ribasim_nl/model.py index 4da12781..bce94e4a 100644 --- a/src/ribasim_nl/ribasim_nl/model.py +++ b/src/ribasim_nl/ribasim_nl/model.py @@ -132,7 +132,9 @@ def basin_node_without_area(self): def upstream_node_id(self, node_id: int): """Get upstream node_id(s)""" - return self.edge.df.set_index("to_node_id").loc[node_id].from_node_id + _df = self.edge.df.set_index("to_node_id") + if node_id in _df.index: + return _df.loc[node_id].from_node_id def upstream_profile(self, node_id: int): """Get upstream basin-profile""" @@ -411,10 +413,28 @@ def reverse_direction_at_node(self, node_id): ].index: self.reverse_edge(edge_id=edge_id) - def update_basin_area(self, node_id: int, basin_area_fid: int | None = None): - if basin_area_fid is None: - basin_area_fid = self.basin.area.df.index.max() + 1 - self.basin.area.df.loc[basin_area_fid, ["node_id"]] = node_id + def select_basin_area(self, geometry): + geometry = shapely.force_2d(geometry) + if isinstance(geometry, MultiPolygon): + polygons = list(geometry.geoms) + elif isinstance(geometry, Polygon): + polygons = [geometry] + else: + raise TypeError("geometry cannot be used for selection, is not a (Multi)Polygon") + + mask = self.basin.area.df.geometry.apply(lambda x: any(x.equals(i) for i in polygons)) + + if not mask.any(): + raise ValueError("Not any basin area equals input geometry") + return mask + + def update_basin_area(self, node_id: int, geometry: Polygon | MultiPolygon, basin_area_fid: int | None = None): + if pd.isna(basin_area_fid): + mask = self.select_basin_area(geometry) + else: + mask = self.basin.area.df.index == basin_area_fid + + self.basin.area.df.loc[mask, ["node_id"]] = node_id def add_basin_area(self, geometry: MultiPolygon, node_id: int | None = None): # if node_id is None, get an available node_id @@ -553,19 +573,45 @@ def edge_to_node_type(self): node_df = self.node_table().df return self.edge.df.to_node_id.apply(lambda x: node_df.at[x, "node_type"] if x in node_df.index else None) - def split_basin(self, line: LineString, basin_id: int | None = None): + def split_basin( + self, + line: LineString | None = None, + basin_id: int | None = None, + geometry: LineString | None = None, + assign_unique_node: bool = True, + ): + if geometry is None: + if line is None: + raise ValueError("geometry cannot be None") + else: + DeprecationWarning("value `line` in split_basin funtion is deprecated. Use `geometry` in stead.") + line = geometry + if self.basin.area.df is None: raise ValueError("provide basin / area table first") - line_centre = line.interpolate(0.5, normalized=True) + line_centre = geometry.interpolate(0.5, normalized=True) + + # if basin_id is supplied, we select by that first if basin_id is not None: basin_area_df = self.basin.area.df.loc[self.basin.area.df.node_id == basin_id] - basin_area_df = self.basin.area.df[self.basin.area.df.contains(line_centre)] + basin_area_df = basin_area_df[basin_area_df.intersects(geometry)] + if len(basin_area_df) > 1: + mask = ~( + basin_area_df.contains(geometry.boundary.geoms[0]) + | basin_area_df.contains(geometry.boundary.geoms[1]) + ) + basin_area_df = basin_area_df[mask] + + else: + basin_area_df = self.basin.area.df[self.basin.area.df.contains(geometry.interpolate(0.5, normalized=True))] - if len(basin_area_df) != 1: - raise ValueError("Overlapping basin-areas at cut_line") + if len(basin_area_df) == 0: + raise ValueError("No basin-areas intersecting cut_line") + elif len(basin_area_df) > 1: + raise ValueError("Multiple Overlapping basin-areas intersecting cut_line") - # get all we need + # get all we need and remove area from basin.area.df basin_fid = int(basin_area_df.iloc[0].name) basin_geometry = basin_area_df.iloc[0].geometry self.basin.area.df = self.basin.area.df[self.basin.area.df.index != basin_fid] @@ -575,7 +621,7 @@ def split_basin(self, line: LineString, basin_id: int | None = None): cut_idx = next(idx for idx, i in enumerate(basin_geoms) if i.contains(line_centre)) # split it - right_basin_poly, left_basin_poly = split_basin(basin_geoms[cut_idx], line).geoms + right_basin_poly, left_basin_poly = split_basin(basin_geoms[cut_idx], geometry).geoms # concat left-over polygons to the right-side right_basin_poly = [right_basin_poly] @@ -591,21 +637,22 @@ def split_basin(self, line: LineString, basin_id: int | None = None): right_basin_poly = MultiPolygon(right_basin_poly) left_basin_poly = MultiPolygon(left_basin_poly) + # add polygons to area for poly in [right_basin_poly, left_basin_poly]: - if self.basin.node.df.geometry.within(poly).any(): - node_ids = self.basin.node.df[self.basin.node.df.geometry.within(poly)].index.to_list() - if len(node_ids) == 1: - if node_ids[0] in self.basin.area.df.node_id.to_numpy(): - self.basin.area.df.loc[self.basin.area.df.node_id.isin(node_ids), ["geometry"]] = poly - else: - self.basin.area.df.loc[self.basin.area.df.index.max() + 1] = { - "node_id": node_ids[0], - "geometry": poly, - } - else: - self.basin.area.df.loc[self.basin.area.df.index.max() + 1, ["geometry"]] = poly - else: - self.basin.area.df.loc[self.basin.area.df.index.max() + 1, ["geometry"]] = poly + # by default we assign provided basin_id + kwargs = { + "node_id": basin_id, + "geometry": poly, + } + + # we override node_id to a unique basin-node within area if that is specified + if assign_unique_node: + if self.basin.node.df.geometry.within(poly).any(): + node_ids = self.basin.node.df[self.basin.node.df.geometry.within(poly)].index.to_list() + if node_ids[0] not in self.basin.area.df.node_id.to_numpy(): + kwargs["node_id"] = node_ids[0] + + self.basin.area.df.loc[self.basin.area.df.index.max() + 1] = kwargs if self.basin.area.df.crs is None: self.basin.area.df.crs = self.crs @@ -631,6 +678,20 @@ def remove_unassigned_basin_area(self): df.index.name = "fid" self.basin.area.df = df + def explode_basin_area(self, remove_z=True): + df = self.basin.area.df.explode().reset_index(drop=True) + df.index.name = "fid" + self.basin.area.df = df + + if remove_z: + self.basin.area.df.loc[:, "geometry"] = gpd.GeoSeries( + shapely.force_2d(self.basin.area.df.geometry.array), crs=self.basin.area.df.crs + ) + + def remove_basin_area(self, geometry): + mask = self.select_basin_area(geometry) + self.basin.area.df = self.basin.area.df[~mask] + def merge_basins( self, basin_id: int | None = None, diff --git a/src/ribasim_nl/ribasim_nl/network.py b/src/ribasim_nl/ribasim_nl/network.py index 6e628091..138df43c 100644 --- a/src/ribasim_nl/ribasim_nl/network.py +++ b/src/ribasim_nl/ribasim_nl/network.py @@ -9,6 +9,7 @@ import pandas as pd from geopandas import GeoDataFrame, GeoSeries from networkx import DiGraph, Graph, NetworkXNoPath, shortest_path, traversal +from shapely import force_2d from shapely.geometry import LineString, Point, box from shapely.ops import snap, split @@ -91,6 +92,11 @@ def validate_inputs(self): elif "MultiLineString" in geom_types: self.lines_gdf = self.lines_gdf.explode(index_parts=False) + # remove z-coordinates + self.lines_gdf.loc[:, "geometry"] = gpd.GeoSeries( + force_2d(self.lines_gdf.geometry.array), crs=self.lines_gdf.crs + ) + @classmethod def from_lines_gpkg(cls, gpkg_file: str | Path, layer: str | None = None, **kwargs): """Instantiate class from a lines_gpkg""" @@ -494,9 +500,16 @@ def get_path(self, node_from, node_to, directed=True, weight="length"): return shortest_path(self.graph_undirected, node_from, node_to, weight=weight) def get_links(self, node_from, node_to, directed=True, weight="length"): - path = self.get_path(node_from, node_to, directed, weight) + # get path and edges on path + path = self.get_path(node_from, node_to, directed=directed, weight=weight) edges_on_path = list(zip(path[:-1], path[1:])) - return self.links.set_index(["node_from", "node_to"]).loc[edges_on_path] + + try: + return self.links.set_index(["node_from", "node_to"]).loc[edges_on_path] + except KeyError: # if path only undirected we need to fix edges_on_path + idx = self.links.set_index(["node_from", "node_to"]).index + edges_on_path = [i if i in idx else (i[1], i[0]) for i in edges_on_path] + return self.links.set_index(["node_from", "node_to"]).loc[edges_on_path] def subset_links(self, nodes_from, nodes_to): gdf = pd.concat([self.get_links(node_from, node_to) for node_from, node_to in product(nodes_from, nodes_to)]) diff --git a/src/ribasim_nl/ribasim_nl/reset_index.py b/src/ribasim_nl/ribasim_nl/reset_index.py index 8c21426e..a474e1d9 100644 --- a/src/ribasim_nl/ribasim_nl/reset_index.py +++ b/src/ribasim_nl/ribasim_nl/reset_index.py @@ -5,24 +5,22 @@ from ribasim_nl.case_conversions import pascal_to_snake_case -def reset_index(model: Model, node_start=1): - # only reset if we have to - node_id_min = model.node_table().df.index.min() - node_id_max = model.node_table().df.index.max() +def reset_index(model: Model, node_start=1, edge_start=1, original_index_postfix: str | None = "waterbeheerder"): + # only reset nodes if we have to + node_ids = model.node_table().df.index + node_id_min = node_ids.min() + node_id_max = node_ids.max() expected_length = node_id_max - node_id_min + 1 if not ((node_start == node_id_min) and (expected_length == len(model.node_table().df))): - # make sure column node_id == index - node_ids = model.node_table().df.index + # create a re-index for nodes - # create a new index for re-indexing all tables - index = pd.Series(data=[i + node_start for i in range(len(node_ids))], index=node_ids).astype("int32") + node_index = pd.Series(data=[i + node_start for i in range(len(node_ids))], index=node_ids).astype("int32") - # renumber edges - model.edge.df.loc[:, ["from_node_id"]] = model.edge.df["from_node_id"].apply(lambda x: index[x]) + # re-number from_node_id and to_node_id + model.edge.df.loc[:, ["from_node_id"]] = model.edge.df["from_node_id"].apply(lambda x: node_index[x]) + model.edge.df.loc[:, ["to_node_id"]] = model.edge.df["to_node_id"].apply(lambda x: node_index[x]) - model.edge.df.loc[:, ["to_node_id"]] = model.edge.df["to_node_id"].apply(lambda x: index[x]) - - # renumber tables + # renumber all node-tables (node, static, area, ...) for node_type in model.node_table().df.node_type.unique(): ribasim_node = getattr(model, pascal_to_snake_case(node_type)) for attr in ribasim_node.model_fields.keys(): @@ -30,17 +28,32 @@ def reset_index(model: Model, node_start=1): try: if table.df is not None: if "node_id" in table.df.columns: - table.df.loc[:, "node_id"] = table.df["node_id"].apply(lambda x: index[x]) + table.df.loc[:, "node_id"] = table.df["node_id"].apply(lambda x: node_index[x]) table.df.index += 1 if "listen_node_id" in table.df.columns: - table.df.loc[:, "listen_node_id"] = table.df["listen_node_id"].apply(lambda x: index[x]) + table.df.loc[:, "listen_node_id"] = table.df["listen_node_id"].apply( + lambda x: node_index[x] + ) if table.df.index.name == "node_id": - table.df.index = table.df.reset_index("node_id")["node_id"].apply(lambda x: index[x]) - - # table.df.index.name = "node_id" + table.df.loc[:, f"meta_node_id_{original_index_postfix}"] = table.df.index + table.df.index = pd.Index( + table.df.reset_index("node_id")["node_id"].apply(lambda x: node_index[x]).to_list(), + name="node_id", + ) except KeyError as e: raise KeyError(f"node_id {e} not in {node_type} / {attr} not in node-table") - return model + # only reset nodes if we have to + edge_ids = model.edge.df.index + edge_id_min = edge_ids.min() + edge_id_max = edge_ids.max() + expected_length = edge_id_max - edge_id_min + 1 -# %% + if not ((edge_start == edge_id_min) and (expected_length == len(model.edge.df))): + # create a re-index for edges + model.edge.df.index = pd.Index([i + node_start for i in range(len(edge_ids))], name="edge_id") + + # keep original index if + if original_index_postfix is not None: + model.edge.df.loc[:, f"meta_edge_id_{original_index_postfix}"] = edge_ids + return model diff --git a/src/ribasim_nl/ribasim_nl/reset_static_tables.py b/src/ribasim_nl/ribasim_nl/reset_static_tables.py new file mode 100644 index 00000000..d0a77a72 --- /dev/null +++ b/src/ribasim_nl/ribasim_nl/reset_static_tables.py @@ -0,0 +1,107 @@ +# %% +import numpy as np +import pandas as pd + + +def reset_static_tables(model): + # basin / profile + if model.basin.node.df is not None: + df = pd.DataFrame( + { + "node_id": np.repeat(model.basin.node.df.index.to_numpy(), 2), + "level": [0.0, 1.0] * len(model.basin.node.df), + "area": [0.01, 1000.0] * len(model.basin.node.df), + } + ) + df.index.name = "fid" + model.basin.profile.df = df + + # basin / state + df = model.basin.profile.df.groupby("node_id")[["level"]].max().reset_index() + df.index.name = "fid" + model.basin.state.df = df + + # basin / static + df = pd.DataFrame( + { + "node_id": model.basin.node.df.index.to_numpy(), + "precipitation": [0] * model.basin.node.df.index.to_numpy(), + "potential_evaporation": [0] * model.basin.node.df.index.to_numpy(), + "drainage": [0] * model.basin.node.df.index.to_numpy(), + "infiltration": [0] * model.basin.node.df.index.to_numpy(), + } + ) + df.index.name = "fid" + model.basin.static.df = df + + # tabulated_rating_curves / static + if model.tabulated_rating_curve.static.df is not None: + df = pd.DataFrame( + { + "node_id": np.repeat(model.tabulated_rating_curve.node.df.index.to_numpy(), 2), + "level": [0.0, 5] * len(model.tabulated_rating_curve.node.df), + "flow_rate": [0, 0.1] * len(model.tabulated_rating_curve.node.df), + } + ) + df.index.name = "fid" + model.tabulated_rating_curve.static.df = df + + # pump / static + if model.pump.static.df is not None: + df = pd.DataFrame( + { + "node_id": model.pump.node.df.index.to_numpy(), + "flow_rate": [0.1] * len(model.pump.node.df), + } + ) + df.index.name = "fid" + model.pump.static.df = df + + # outlet / static + if model.outlet.static.df is not None: + df = pd.DataFrame( + { + "node_id": model.outlet.node.df.index.to_numpy(), + "flow_rate": [0.1] * len(model.outlet.node.df), + } + ) + df.index.name = "fid" + model.outlet.static.df = df + + # flow_boundary / static + if model.flow_boundary.static.df is not None: + df = pd.DataFrame( + { + "node_id": model.flow_boundary.node.df.index.to_numpy(), + "flow_rate": [0.1] * len(model.flow_boundary.node.df), + } + ) + df.index.name = "fid" + model.flow_boundary.static.df = df + + # level_boundary / static + if model.level_boundary.static.df is not None: + df = pd.DataFrame( + { + "node_id": model.level_boundary.node.df.index.to_numpy(), + "level": [0.1] * len(model.level_boundary.node.df), + } + ) + df.index.name = "fid" + model.level_boundary.static.df = df + + # manning_resistance / static + if model.manning_resistance.static.df is not None: + df = pd.DataFrame( + { + "node_id": model.manning_resistance.node.df.index.to_numpy(), + "length": [100] * len(model.manning_resistance.node.df), + "manning_n": [0.04] * len(model.manning_resistance.node.df), + "profile_width": [10] * len(model.manning_resistance.node.df), + "profile_slope": [1] * len(model.manning_resistance.node.df), + } + ) + df.index.name = "fid" + model.manning_resistance.static.df = df + + return model