diff --git a/scripts/build_clustered_population_layouts.py b/scripts/build_clustered_population_layouts.py
index e7f4db5a..94cbd133 100644
--- a/scripts/build_clustered_population_layouts.py
+++ b/scripts/build_clustered_population_layouts.py
@@ -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)
 
diff --git a/scripts/build_heat_demand.py b/scripts/build_heat_demand.py
index 1ed57ea4..865f12bd 100644
--- a/scripts/build_heat_demand.py
+++ b/scripts/build_heat_demand.py
@@ -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)
 
 
diff --git a/scripts/build_solar_thermal_profiles.py b/scripts/build_solar_thermal_profiles.py
index a3a171d9..c26266aa 100644
--- a/scripts/build_solar_thermal_profiles.py
+++ b/scripts/build_solar_thermal_profiles.py
@@ -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)
 
 
diff --git a/scripts/build_temperature_profiles.py b/scripts/build_temperature_profiles.py
index 28a4b858..a55bd606 100644
--- a/scripts/build_temperature_profiles.py
+++ b/scripts/build_temperature_profiles.py
@@ -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)
 
 
diff --git a/scripts/helper.py b/scripts/helper.py
new file mode 100644
index 00000000..6f060db4
--- /dev/null
+++ b/scripts/helper.py
@@ -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)