Move GDAL constraint into constraints file and update to GDAL 3.8

.github/workflows/build-and-test.yml

@@ -79,7 +79,6 @@ jobs:
         with:
           fetch-depth: 0 # Fetch complete history for accurate versioning
 
-
       # NOTE: It takes twice as long to save the sample data cache
       # as it does to do a fresh clone (almost 5 minutes vs. 2.5 minutes)
       # Test data is way, way faster by contrast (on the order of a few
@@ -94,7 +93,7 @@ jobs:
         uses: ./.github/actions/setup_env
         with:
           python-version: ${{ matrix.python-version }}
-          requirements-files: requirements.txt requirements-dev.txt
+          requirements-files: requirements.txt requirements-dev.txt constraints_tests.txt
           requirements: ${{ env.CONDA_DEFAULT_DEPENDENCIES }}
 
       - name: Download previous conda environment.yml

.readthedocs_environment.yml

@@ -11,7 +11,7 @@ channels:
 - nodefaults
 dependencies:
 - python=3.11
-- gdal>=3.4.2,<3.6.0
+- gdal>=3.4.2
 - pip
 - pip:
   - -r requirements.txt

constraints_tests.txt

@@ -0,0 +1,5 @@
+# This file contains package constraints needed to run the invest test suite
+
+# A gdal bug caused our test suite to fail, but this issue is unlikely to
+# occur with regular use of invest. https://github.com/OSGeo/gdal/issues/8497
+GDAL!=3.6.*,!=3.7.*

pyproject.toml

@@ -78,4 +78,4 @@ where = ["src"]
 
 [tool.pytest.ini_options]
 # raise warnings to errors, except for deprecation warnings
-filterwarnings = ["error", "default::DeprecationWarning"]
+filterwarnings = ["error", "default::DeprecationWarning", "default::FutureWarning"]

requirements.txt

@@ -10,7 +10,7 @@
 # scripts/convert-requirements-to-conda-yml.py as though it can only be found
 # on pip.
 
-GDAL>=3.4.2,<3.6.0
+GDAL>=3.4.2
 Pyro4==4.77  # pip-only
 pandas>=1.2.1
 numpy>=1.11.0,!=1.16.0

src/natcap/invest/hra.py

@@ -633,10 +633,6 @@ def execute(args):
         # If the input is a vector, reproject to the AOI SRS and simplify.
         # Rasterization happens in the alignment step.
         elif gis_type == pygeoprocessing.VECTOR_TYPE:
-            # Habitats and stressors are rasterized with ALL_TOUCHED=TRUE
-            if name in habitats_info or name in stressors_info:
-                habitat_stressor_vectors.add(source_filepath)
-
             # Using Shapefile here because its driver appears to not raise a
             # warning if a MultiPolygon geometry is inserted into a Polygon
             # layer, which was happening on a real-world sample dataset while
@@ -680,6 +676,10 @@ def execute(args):
                 dependent_task_list=[reprojected_vector_task]
             ))
 
+            # Habitats and stressors are rasterized with ALL_TOUCHED=TRUE
+            if name in habitats_info or name in stressors_info:
+                habitat_stressor_vectors.add(target_simplified_vector)
+
         # Later operations make use of the habitats rasters or the stressors
         # rasters, so it's useful to collect those here now.
         if name in habitats_info:
@@ -1648,7 +1648,6 @@ def _simplify(source_vector_path, tolerance, target_vector_path,
     for source_feature in source_layer:
         target_feature = ogr.Feature(target_layer_defn)
         source_geom = source_feature.GetGeometryRef()
-
         simplified_geom = source_geom.SimplifyPreserveTopology(tolerance)
         if simplified_geom is not None:
             target_geom = simplified_geom

src/natcap/invest/ndr/ndr_core.pyx

@@ -415,9 +415,9 @@ def ndr_eff_calculation(
 
     # create direction raster in bytes
     def _mfd_to_flow_dir_op(mfd_array):
-        result = numpy.zeros(mfd_array.shape, dtype=numpy.int8)
+        result = numpy.zeros(mfd_array.shape, dtype=numpy.uint8)
         for i in range(8):
-            result[:] |= (((mfd_array >> (i*4)) & 0xF) > 0) << i
+            result[:] |= ((((mfd_array >> (i*4)) & 0xF) > 0) << i).astype(numpy.uint8)
         return result
 
     pygeoprocessing.raster_calculator(

tests/test_validation.py

@@ -2002,6 +2002,7 @@ def test_spatial_overlap_error(self):
         layer = vector.CreateLayer('layer', vector_srs, ogr.wkbPoint)
         new_feature = ogr.Feature(layer.GetLayerDefn())
         new_feature.SetGeometry(ogr.CreateGeometryFromWkt('POINT (1 1)'))
+        layer.CreateFeature(new_feature)
 
         new_feature = None
         layer = None