Skip to content
Open
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
1 change: 1 addition & 0 deletions config/config.default.yaml
Original file line number Diff line number Diff line change
Expand Up @@ -1100,6 +1100,7 @@ clustering:
remove_islands: false
aggregate_to_tyndp: false
simplify_network:
to_380: true
to_substations: false
remove_stubs: true
remove_stubs_across_borders: false
Expand Down
24 changes: 17 additions & 7 deletions scripts/build_powerplants.py
Original file line number Diff line number Diff line change
Expand Up @@ -175,12 +175,16 @@ def map_to_country_bus(
unmatched = []

for country, plants in ppl.groupby("Country"):
country_regions = regions[regions.index.str[:2] == country]
joined = (
plants.sjoin(country_regions)
.rename(columns={"name": "bus"})
.reindex(plants.index)
if "country" in regions.columns:
country_regions = regions[regions["country"] == country]
else:
country_regions = regions[regions.index.str[:2] == country]
joined = plants.sjoin(country_regions[["geometry"]]).rename(
columns={"name": "bus"}
)
# Drop duplicate matches (plant in overlapping onshore/offshore regions)
joined = joined[~joined.index.duplicated(keep="first")]
joined = joined.reindex(plants.index)
assigned.append(joined.dropna(subset=["bus"]))
missing = joined[joined["bus"].isna()]
if not missing.empty:
Expand All @@ -189,10 +193,16 @@ def map_to_country_bus(
if unmatched:
unmatched = pd.concat(unmatched)
for country, plants in unmatched.groupby("Country"):
country_regions = regions[regions.index.str[:2] == country]
if "country" in regions.columns:
country_regions = regions[regions["country"] == country]
else:
country_regions = regions[regions.index.str[:2] == country]
nearest = (
plants.to_crs(3035)
.sjoin_nearest(country_regions.to_crs(3035), max_distance=max_distance)
.sjoin_nearest(
country_regions[["geometry"]].to_crs(3035),
max_distance=max_distance,
)
.rename(columns={"name": "bus"})
.to_crs(4326)
)
Expand Down
4 changes: 4 additions & 0 deletions scripts/build_renewable_profiles.py
Original file line number Diff line number Diff line change
Expand Up @@ -157,6 +157,7 @@
if snakemake.wildcards.technology.startswith("offwind"):
# for offshore regions, the shortest distance to the shoreline is used
offshore_regions = availability.coords["bus"].values
offshore_regions = np.intersect1d(offshore_regions, regions.index)
Comment thread
fneum marked this conversation as resolved.
regions = regions.loc[offshore_regions]
regions = regions.map(lambda g: _simplify_polys(g, minarea=1)).set_crs(
regions.crs
Expand Down Expand Up @@ -294,6 +295,9 @@
average_distance = []
bus_bins = layoutmatrix.indexes["bus_bin"]
for bus, bin in bus_bins:
if bus not in regions.index:
average_distance.append(0.0)
Comment thread
fneum marked this conversation as resolved.
continue
row = layoutmatrix.sel(bus=bus, bin=bin).data
nz_b = row != 0
row = row[nz_b]
Expand Down
6 changes: 4 additions & 2 deletions scripts/cluster_network.py
Original file line number Diff line number Diff line change
Expand Up @@ -425,7 +425,7 @@ def clustering_for_n_clusters(
busmap,
bus_strategies=bus_strategies,
line_strategies=line_strategies,
custom_line_groupers=["build_year"],
custom_line_groupers=["build_year", "v_nom"],
)

return clustering
Expand Down Expand Up @@ -703,7 +703,9 @@ def update_bus_coordinates(
# nc.shapes = n.shapes.copy()
for which in ["regions_onshore", "regions_offshore"]:
regions = gpd.read_file(snakemake.input[which])
clustered_regions = cluster_regions((clustering.busmap,), regions)
clustered_regions = cluster_regions(
(clustering.busmap,), regions, with_country=True
)
clustered_regions.to_file(snakemake.output[which])
# append_bus_shapes(nc, clustered_regions, type=which.split("_")[1])

Expand Down
58 changes: 42 additions & 16 deletions scripts/determine_availability_matrix_MD_UA.py
Original file line number Diff line number Diff line change
Expand Up @@ -15,6 +15,7 @@
import fiona
import geopandas as gpd
import numpy as np
import xarray as xr

from scripts._helpers import configure_logging, load_cutout, set_scenario_config

Expand Down Expand Up @@ -57,6 +58,8 @@ def get_wdpa_layer_name(wdpa_fn, layer_substring):
excluder = atlite.ExclusionContainer(crs=3035, res=100)

corine = config.get("corine", {})
if not isinstance(corine, dict):
corine = {}
if "grid_codes" in corine:
# Land cover codes to emulate CORINE results
if snakemake.wildcards.technology == "solar":
Expand Down Expand Up @@ -85,18 +88,35 @@ def get_wdpa_layer_name(wdpa_fn, layer_substring):
snakemake.input.copernicus, codes=codes, buffer=buffer, crs="EPSG:4326"
)

has_valid_geometry = not regions.empty and not regions.geometry.isna().all()

if not has_valid_geometry:
logger.info("No valid MD/UA regions; writing empty availability matrix.")
availability = xr.DataArray(
np.zeros((len(buses), len(cutout.data.y), len(cutout.data.x))),
dims=["bus", "y", "x"],
coords={"bus": buses, "y": cutout.data.y, "x": cutout.data.x},
)
availability.to_netcdf(snakemake.output.availability_matrix)
raise SystemExit(0)

if config["natura"]:
wdpa_fn = (
snakemake.input.wdpa_marine
if "offwind" in snakemake.wildcards.technology
else snakemake.input.wdpa
)
layer = get_wdpa_layer_name(wdpa_fn, "polygons")
wdpa = gpd.read_file(
wdpa_fn,
bbox=regions.geometry,
layer=layer,
).to_crs(3035)

if has_valid_geometry:
layer = get_wdpa_layer_name(wdpa_fn, "polygons")
wdpa = gpd.read_file(
wdpa_fn,
bbox=regions.geometry,
layer=layer,
).to_crs(3035)
else:
logger.info("No valid region geometries for MD/UA; skipping WDPA polygons.")
wdpa = gpd.GeoDataFrame()

# temporary file needed for parallelization
with NamedTemporaryFile(suffix=".geojson", delete=False) as f:
Expand All @@ -107,17 +127,23 @@ def get_wdpa_layer_name(wdpa_fn, layer_substring):
time.sleep(1)
excluder.add_geometry(plg_tmp_fn)

layer = get_wdpa_layer_name(wdpa_fn, "points")
wdpa_pts = gpd.read_file(
wdpa_fn,
bbox=regions.geometry,
layer=layer,
).to_crs(3035)
if has_valid_geometry:
layer = get_wdpa_layer_name(wdpa_fn, "points")
wdpa_pts = gpd.read_file(
wdpa_fn,
bbox=regions.geometry,
layer=layer,
).to_crs(3035)
else:
logger.info("No valid region geometries for MD/UA; skipping WDPA points.")
wdpa_pts = gpd.GeoDataFrame(columns=["REP_AREA", "geometry"])

wdpa_pts = wdpa_pts[wdpa_pts["REP_AREA"] > 1]
wdpa_pts["buffer_radius"] = np.sqrt(wdpa_pts["REP_AREA"] / np.pi) * 1000
wdpa_pts = wdpa_pts.set_geometry(
wdpa_pts["geometry"].buffer(wdpa_pts["buffer_radius"])
)
if not wdpa_pts.empty:
wdpa_pts["buffer_radius"] = np.sqrt(wdpa_pts["REP_AREA"] / np.pi) * 1000
wdpa_pts = wdpa_pts.set_geometry(
wdpa_pts["geometry"].buffer(wdpa_pts["buffer_radius"])
)

# temporary file needed for parallelization
with NamedTemporaryFile(suffix=".geojson", delete=False) as f:
Expand Down
11 changes: 9 additions & 2 deletions scripts/simplify_network.py
Original file line number Diff line number Diff line change
Expand Up @@ -432,8 +432,15 @@ def remove_converters(n: pypsa.Network) -> pypsa.Network:
Nyears = n.snapshot_weightings.objective.sum() / 8760
buses_prev, lines_prev, links_prev = len(n.buses), len(n.lines), len(n.links)

linetype_380 = snakemake.config["lines"]["types"][380]
n, trafo_map = simplify_network_to_380(n, linetype_380)
if params.simplify_network.get("to_380", True):
linetype_380 = snakemake.config["lines"]["types"][380]
n, trafo_map = simplify_network_to_380(n, linetype_380)
else:
logger.info(
"Skipping 380kV simplification — preserving original voltage levels"
" and transformers"
)
trafo_map = n.buses.index.to_series()
busmaps = [trafo_map]

n, converter_map = remove_converters(n)
Expand Down
Loading