Fix bug due to invalid cluster bus region geometries

These geometries are apparently invalid due to self-crossing or
self-touching polygons.

The geometries are created by pypsa-eur/scripts/cluster_network.py but
appear to be valid before being written to file.

They are only valid after being read back in from file.

This seems to indicate some numerical issue relating to file reading
and writing.

Now the geometries are cleaned after being read in.
This commit is contained in:
Tom Brown 2020-05-07 14:45:14 +02:00
parent e92558716e
commit 8538233bba
5 changed files with 26 additions and 10 deletions

View File

@ -3,7 +3,7 @@ import geopandas as gpd
import xarray as xr
import pandas as pd
import atlite
import helper
cutout = atlite.Cutout(snakemake.config['atlite']['cutout_name'],
@ -14,7 +14,7 @@ clustered_busregions_as_geopd = gpd.read_file(snakemake.input.regions_onshore).s
clustered_busregions = pd.Series(clustered_busregions_as_geopd.geometry, index=clustered_busregions_as_geopd.index)
helper.clean_invalid_geometries(clustered_busregions)
I = cutout.indicatormatrix(clustered_busregions)

View File

@ -4,7 +4,7 @@ import atlite
import pandas as pd
import xarray as xr
import scipy as sp
import helper
if 'snakemake' not in globals():
from vresutils import Dict
@ -26,6 +26,8 @@ clustered_busregions_as_geopd = gpd.read_file(snakemake.input.regions_onshore).s
clustered_busregions = pd.Series(clustered_busregions_as_geopd.geometry, index=clustered_busregions_as_geopd.index)
helper.clean_invalid_geometries(clustered_busregions)
I = cutout.indicatormatrix(clustered_busregions)

View File

@ -1,15 +1,10 @@
import sys
#override atlite
sys.path = ["/home/vres/data/tom/lib/atlite"] + sys.path
import geopandas as gpd
import atlite
import pandas as pd
import xarray as xr
import scipy as sp
import helper
if 'snakemake' not in globals():
from vresutils import Dict
@ -33,6 +28,8 @@ clustered_busregions_as_geopd = gpd.read_file(snakemake.input.regions_onshore).s
clustered_busregions = pd.Series(clustered_busregions_as_geopd.geometry, index=clustered_busregions_as_geopd.index)
helper.clean_invalid_geometries(clustered_busregions)
I = cutout.indicatormatrix(clustered_busregions)

View File

@ -4,7 +4,7 @@ import atlite
import pandas as pd
import xarray as xr
import scipy as sp
import helper
if 'snakemake' not in globals():
from vresutils import Dict
@ -27,6 +27,8 @@ clustered_busregions_as_geopd = gpd.read_file(snakemake.input.regions_onshore).s
clustered_busregions = pd.Series(clustered_busregions_as_geopd.geometry, index=clustered_busregions_as_geopd.index)
helper.clean_invalid_geometries(clustered_busregions)
I = cutout.indicatormatrix(clustered_busregions)

15
scripts/helper.py Normal file
View File

@ -0,0 +1,15 @@
import logging
logger = logging.getLogger(__name__)
#https://stackoverflow.com/questions/20833344/fix-invalid-polygon-in-shapely
#https://stackoverflow.com/questions/13062334/polygon-intersection-error-in-shapely-shapely-geos-topologicalerror-the-opera
#https://shapely.readthedocs.io/en/latest/manual.html#object.buffer
def clean_invalid_geometries(geometries):
"""Fix self-touching or self-crossing polygons; these seem to appear
due to numerical problems from writing and reading, since the geometries
are valid before being written in pypsa-eur/scripts/cluster_network.py"""
for i,p in geometries.items():
if not p.is_valid:
logger.warning(f'Clustered region {i} had an invalid geometry, fixing using zero buffer.')
geometries[i] = p.buffer(0)