Repository: PyPSA/atlite

Issue:
Shapefile path added through `add_geometry` not read when no rasters are provided
<!---
SPDX-FileCopyrightText: 2021 The Atlite Authors

SPDX-License-Identifier: CC0-1.0
--->


<!-- Provide a general summary of the issue -->

## Description
<!-- Provide a more detailed introduction to the issue itself, and why you consider it to be a bug -->
<!-- If you can, add a minimal example which reproduces the bug -->

This is an odd one: If the exclusion container only has geometries specified as a string but no rasters through `add_raster`, the shapefiles are not read in.

## Expected Behavior
<!-- Tell us what should happen -->

`excluder.add_geometry(str)` should read file str into geopandas dataframe/series

## Actual Behavior
<!-- Tell us what goes wrong and happens instead -->

Somewhere iteration over the string characters (see error message) not list of paths.

Not ok:

```py
excluder = atlite.ExclusionContainer(crs=3035, res=250)
shapes = regions.geometry.to_crs(excluder.crs)
excluder.add_geometry(NATURA_SHP)
masked, transform = shape_availability(shapes, excluder)
```

Also not ok:

```py
excluder = atlite.ExclusionContainer(crs=3035, res=250)
shapes = regions.geometry.to_crs(excluder.crs)
excluder.add_geometry(NATURA_SHP)
excluder.add_geometry(EEZ, invert=True)
masked, transform = shape_availability(shapes, excluder)
```

Ok:

```py
excluder = atlite.ExclusionContainer(crs=3035, res=250)
shapes = regions.geometry.to_crs(excluder.crs)
excluder.add_geometry(NATURA_SHP)
excluder.add_raster(CORINE, codes=[1], invert=True, crs=3035)
masked, transform = shape_availability(shapes, excluder)
```

## Error Message
<!-- Paste any terminal output and error message you encounter here to help illustrate the problem -->

```
/home/fneum/miniconda3/envs/pypsa-eur/lib/python3.9/site-packages/rasterio/features.py:318: ShapeSkipWarning: Invalid or empty shape . at index 0 will not be rasterized.
  warnings.warn('Invalid or empty shape {} at index {} will not be rasterized.'.format(geom, index), ShapeSkipWarning)
/home/fneum/miniconda3/envs/pypsa-eur/lib/python3.9/site-packages/rasterio/features.py:318: ShapeSkipWarning: Invalid or empty shape . at index 1 will not be rasterized.
  warnings.warn('Invalid or empty shape {} at index {} will not be rasterized.'.format(geom, index), ShapeSkipWarning)
/home/fneum/miniconda3/envs/pypsa-eur/lib/python3.9/site-packages/rasterio/features.py:318: ShapeSkipWarning: Invalid or empty shape / at index 2 will not be rasterized.
  warnings.warn('Invalid or empty shape {} at index {} will not be rasterized.'.format(geom, index), ShapeSkipWarning)
/home/fneum/miniconda3/envs/pypsa-eur/lib/python3.9/site-packages/rasterio/features.py:318: ShapeSkipWarning: Invalid or empty shape . at index 3 will not be rasterized.
  warnings.warn('Invalid or empty shape {} at index {} will not be rasterized.'.format(geom, index), ShapeSkipWarning)
/home/fneum/miniconda3/envs/pypsa-eur/lib/python3.9/site-packages/rasterio/features.py:318: ShapeSkipWarning: Invalid or empty shape . at index 4 will not be rasterized.
  warnings.warn('Invalid or empty shape {} at index {} will not be rasterized.'.format(geom, index), ShapeSkipWarning)
/home/fneum/miniconda3/envs/pypsa-eur/lib/python3.9/site-packages/rasterio/features.py:318: ShapeSkipWarning: Invalid or empty shape / at index 5 will not be rasterized.
```

## Your Environment
<!-- Include relevant details about the environment you experienced the bug in -->
* The `atlite` version used:
* How you installed `atlite` (`conda`, `pip` or `github`):
* Operating System:
* My environment:
    <details>
      <summary>(output of `conda list`)</summary>
      ```
        <!-- output of `conda list` -->
      ```
    </details>




Discussion:

Patch:
diff --git a/RELEASE_NOTES.rst b/RELEASE_NOTES.rst
index b3f92d6..05df63c 100644
--- a/RELEASE_NOTES.rst
+++ b/RELEASE_NOTES.rst
@@ -14,6 +14,9 @@ Release Notes
 * Bugfix: When creating cutouts using SARAH2 data, an error was previously wrongly thrown if exactly
   the data was available as input as required. The error is now correctly thrown only if
   insufficient SARAH data is available.
+* Bugfix: When only adding geometries to an `atlite.ExclusionContainer` the geometries were previously
+  not opened and an error was thrown. The error did not occur if one or more shapes were included.
+  Error is corrected and geometry-only exclusions can now be calculated. (GH Issue #225)
 
 Version 0.2.7 
 ==============
diff --git a/atlite/gis.py b/atlite/gis.py
index 3f3dab5..23b3827 100644
--- a/atlite/gis.py
+++ b/atlite/gis.py
@@ -309,12 +309,16 @@ class ExclusionContainer:
     @property
     def all_closed(self):
         """Check whether all files in the raster container are closed."""
-        return all(isinstance(d["raster"], (str, Path)) for d in self.rasters)
+        return all(isinstance(d["raster"], (str, Path)) for d in self.rasters) and all(
+            isinstance(d["geometry"], (str, Path)) for d in self.geometries
+        )
 
     @property
     def all_open(self):
         """Check whether all files in the raster container are open."""
-        return all(isinstance(d["raster"], rio.DatasetReader) for d in self.rasters)
+        return all(
+            isinstance(d["raster"], rio.DatasetReader) for d in self.rasters
+        ) and all(isinstance(d["geometry"], gpd.GeoSeries) for d in self.geometries)
 
     def __repr__(self):
         return (


Test patch:
diff --git a/test/test_gis.py b/test/test_gis.py
index 36a0c9e..d425a93 100755
--- a/test/test_gis.py
+++ b/test/test_gis.py
@@ -44,6 +44,24 @@ def ref():
     )
 
 
+@pytest.fixture(scope="session")
+def geometry(tmp_path_factory):
+    tmp_path = tmp_path_factory.mktemp("geometries")
+    geometry = gpd.GeoSeries(
+        [box(X0 / 2 + X1 / 2, Y0 / 2 + Y1 / 2, X1, Y1)],
+        crs="EPSG:4326",
+        index=[0],
+        name="boxes",
+    )
+    geometry = geometry.to_frame().set_geometry("boxes")
+
+    path = tmp_path / "geometry.gpkg"
+
+    geometry.to_file(path, driver="GPKG")
+
+    return path
+
+
 @pytest.fixture(scope="session")
 def raster(tmp_path_factory):
     tmp_path = tmp_path_factory.mktemp("rasters")
@@ -116,6 +134,54 @@ def raster_codes(tmp_path_factory):
     return path
 
 
+def test_open_closed_checks(ref, geometry, raster):
+    """Test atlite.ExclusionContainer(...) file open/closed checks for plausibility. C.f. GH issue #225."""
+
+    res = 0.01
+    excluder = ExclusionContainer(ref.crs, res=res)
+
+    # Without raster/shapes, both should evaluate to True
+    assert excluder.all_closed and excluder.all_open
+
+    # First add geometries, than raster
+    excluder.add_geometry(geometry)
+    assert excluder.all_closed and not excluder.all_open
+
+    # Check if still works with 2nd geometry
+    excluder.add_geometry(geometry)
+    assert excluder.all_closed and not excluder.all_open
+
+    excluder.add_raster(raster)
+    assert excluder.all_closed and not excluder.all_open
+
+    # Check if still works with 2nd raster
+    excluder.add_raster(raster)
+    assert excluder.all_closed and not excluder.all_open
+
+    excluder.open_files()
+    assert not excluder.all_closed and excluder.all_open
+
+    # First add raster, then geometries
+    excluder = ExclusionContainer(ref.crs, res=res)
+
+    excluder.add_raster(raster)
+    assert excluder.all_closed and not excluder.all_open
+
+    # 2nd raster
+    excluder.add_raster(raster)
+    assert excluder.all_closed and not excluder.all_open
+
+    excluder.add_geometry(geometry)
+    assert excluder.all_closed and not excluder.all_open
+
+    # 2nd geometry
+    excluder.add_geometry(geometry)
+    assert excluder.all_closed and not excluder.all_open
+
+    excluder.open_files()
+    assert not excluder.all_closed and excluder.all_open
+
+
 def test_transform():
     """Test the affine transform. It always has to point to cell origin."""
     cutout = Cutout(