include function find_buses_area for clustering

etrago/cluster/electrical.py

@@ -34,7 +34,6 @@
         get_clustering_from_busmap,
     )
     from six import iteritems
-    import geopandas as gpd
     import numpy as np
     import pandas as pd
     import pypsa.io as io
@@ -46,9 +45,9 @@
         kmedoids_dijkstra_clustering,
         strategies_generators,
         strategies_one_ports,
+        find_buses_area
     )
     from etrago.tools.utilities import *
-    from shapely.geometry import Point
 
     logger = logging.getLogger(__name__)
 
@@ -508,42 +507,9 @@ def select_elec_network(etrago):
     settings = etrago.args["network_clustering"]
 
     # Exclude buses in the area that should not be clustered
-    if settings["exclusion_area"]:
-        con = etrago.engine
-        query = "SELECT gen, geometry FROM boundaries.vg250_krs"
-
-        de_areas = gpd.read_postgis(query, con, geom_col="geometry")
-        de_areas = de_areas[de_areas["gen"].isin(settings["exclusion_area"])]
-
-        try:
-            buses_area = gpd.GeoDataFrame(
-                etrago.network.buses, geometry="geom", crs=4326
-            )
-        except:
-            buses_area = etrago.network.buses[
-                [
-                    "x",
-                    "y",
-                ]
-            ]
-            buses_area["geom"] = buses_area.apply(
-                lambda x: Point(x["x"], x["y"]), axis=1
-            )
-            buses_area = gpd.GeoDataFrame(
-                buses_area, geometry="geom", crs=4326
-            )
-
-        buses_area = gpd.clip(buses_area, de_areas)
-        elec_network.buses = elec_network.buses[
-            ~elec_network.buses.index.isin(buses_area.index)
-        ]
-
-        busmap_area = pd.Series(
-            buses_area.index.rename("bus_area"),
-            index=buses_area.index.rename("bus"),
-        )
-    else:
-        busmap_area = pd.DataFrame()
+    busmap_area = find_buses_area(etrago, "AC")
+    elec_network.buses = elec_network.buses[
+        ~elec_network.buses.index.isin(busmap_area.index)]
 
     # Exclude foreign buses when it is set to don't include them in the clustering
     if settings["cluster_foreign_AC"]:

etrago/cluster/gas.py

@@ -39,6 +39,7 @@
         group_links,
         kmedoids_dijkstra_clustering,
         sum_with_inf,
+        find_buses_area
     )
     from etrago.tools.utilities import *
 
@@ -101,7 +102,8 @@ def preprocessing(etrago):
         )
     ]
 
-    # select buses dependent on whether they should be clustered in (only DE or DE+foreign)
+    # select buses dependent on whether they should be clustered in (only DE
+    # or DE+foreign DE except specified area)
     if not settings["cluster_foreign_gas"]:
         network_ch4.buses = network_ch4.buses.loc[
             ch4_filter & (network_ch4.buses["country"].values == "DE")
@@ -122,6 +124,11 @@ def preprocessing(etrago):
         network_ch4.buses = network_ch4.buses.loc[ch4_filter]
         n_clusters = settings["n_clusters_gas"]
 
+    # Exclude buses in the area that should not be clustered
+    busmap_area = find_buses_area(etrago, "CH4")
+    network_ch4.buses = network_ch4.buses[
+        ~network_ch4.buses.index.isin(busmap_area.index)]
+
     def weighting_for_scenario(ch4_buses, save=None):
         """
         Calculate CH4-bus weightings dependant on the connected
@@ -215,7 +222,7 @@ def weighting_for_scenario(ch4_buses, save=None):
         weight_ch4.loc[loaded_weights.index] = loaded_weights
     else:
         weight_ch4 = weighting_for_scenario(network_ch4.buses, save=False)
-    return network_ch4, weight_ch4.squeeze(), n_clusters
+    return network_ch4, weight_ch4.squeeze(), n_clusters, busmap_area
 
 
 def kmean_clustering_gas(etrago, network_ch4, weight, n_clusters):
@@ -293,7 +300,7 @@ def get_h2_clusters(etrago, busmap_ch4):
     return busmap
 
 
-def gas_postprocessing(etrago, busmap, medoid_idx=None):
+def gas_postprocessing(etrago, busmap, medoid_idx=None, busmap_area=pd.DataFrame()):
     """
     Performs the postprocessing for the gas grid clustering based on the
     provided busmap
@@ -345,7 +352,11 @@ def gas_postprocessing(etrago, busmap, medoid_idx=None):
                 + str(settings["n_clusters_gas"])
                 + "_result.csv"
             )
-
+    breakpoint()
+    ###########################################################################
+    #include busmap_area to busmap
+    ###########################################################################
+
     busmap = get_h2_clusters(etrago, busmap)
 
     # Add all other buses to busmap
@@ -947,7 +958,7 @@ def run_spatial_clustering_gas(self):
             method = settings["method_gas"]
             logger.info(f"Start {method} clustering GAS")
 
-            gas_network, weight, n_clusters = preprocessing(self)
+            gas_network, weight, n_clusters, busmap_area = preprocessing(self)
 
             if method == "kmeans":
                 if settings["k_gas_busmap"]:
@@ -991,7 +1002,7 @@ def run_spatial_clustering_gas(self):
                     "spatial clustering method for the gas network"
                 )
                 raise ValueError(msg)
-            self.network, busmap = gas_postprocessing(self, busmap, medoid_idx)
+            self.network, busmap = gas_postprocessing(self, busmap, medoid_idx, busmap_area)
 
             self.update_busmap(busmap)
 

etrago/cluster/spatial.py

@@ -38,9 +38,11 @@
     )
     from sklearn.cluster import KMeans
     from threadpoolctl import threadpool_limits
+    from shapely.geometry import Point
     import networkx as nx
     import numpy as np
     import pandas as pd
+    import geopandas as gpd
 
     from etrago.tools.utilities import *
 
@@ -797,3 +799,47 @@ def kmedoids_dijkstra_clustering(
         busmap.index.name = "bus_id"
 
     return busmap, medoid_idx
+
+def find_buses_area(etrago, carrier):
+    """
+    Find buses of a specified carrier in a defined area. Usually used to
+    findout the buses that sould not be clustered.
+    """
+    settings = etrago.args["network_clustering"]
+
+    if settings["exclusion_area"]:
+        con = etrago.engine
+        query = "SELECT gen, geometry FROM boundaries.vg250_krs"
+
+        de_areas = gpd.read_postgis(query, con, geom_col="geometry")
+        de_areas = de_areas[de_areas["gen"].isin(settings["exclusion_area"])]
+
+        try:
+            buses_area = gpd.GeoDataFrame(
+                etrago.network.buses, geometry="geom", crs=4326
+            )
+        except:
+            buses_area = etrago.network.buses[
+                [
+                    "x",
+                    "y",
+                    "carrier"
+                ]
+            ]
+            buses_area["geom"] = buses_area.apply(
+                lambda x: Point(x["x"], x["y"]), axis=1
+            )
+            buses_area = gpd.GeoDataFrame(
+                buses_area, geometry="geom", crs=4326
+            )
+
+        buses_area = gpd.clip(buses_area, de_areas)
+        buses_area = buses_area[buses_area.carrier == carrier]
+        busmap_area = pd.Series(
+            buses_area.index.rename("bus_area"),
+            index=buses_area.index.rename("bus"),
+        )
+    else:
+        busmap_area = pd.DataFrame()
+
+    return busmap_area