Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Fix rasterization of spatial criteria relative to InVEST 3.9 #1121

Merged
merged 2 commits into from
Nov 22, 2022
Merged

Fix rasterization of spatial criteria relative to InVEST 3.9 #1121

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
2 commits
Select commit Hold shift + click to select a range
130e965
Restricting which vectors rasterized with ALL_TOUCHED=TRUE.
phargogh Nov 21, 2022
7981e9d
Updating tests to reflect new functionality.
phargogh Nov 21, 2022
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
6 changes: 6 additions & 0 deletions HISTORY.rst
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -45,6 +45,12 @@ Unreleased Changes
now reprojected to the ``lulc_cur_path`` raster. This fixes a bug where
rasters with a different SRS would appear to not intersect the
``lulc_cur_path`` even if they did. (https://github.com/natcap/invest/issues/1093)
* HRA
* Fixed a regression relative to InVEST 3.9.0 outputs where spatial
criteria vectors were being rasterized with the ``ALL_TOUCHED=TRUE``
flag, leading to a perceived buffering of spatial criteria in certain
cases. In InVEST 3.9.0, these were rasterized with ``ALL_TOUCHED=FALSE``.
https://github.com/natcap/invest/issues/1120
* Workbench
* Fixed a bug where the Workbench would become unresponsive during an
InVEST model run if the model emitted a very high volume of log messages.
Expand Down
20 changes: 18 additions & 2 deletions src/natcap/invest/hra.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -361,6 +361,7 @@ def execute(args):
alignment_dependent_tasks = []
aligned_habitat_raster_paths = {}
aligned_stressor_raster_paths = {}
habitat_stressor_vectors = set([])
for name, attributes in itertools.chain(habitats_info.items(),
stressors_info.items(),
spatial_criteria_attrs.items()):
Expand Down Expand Up @@ -415,6 +416,10 @@ 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
Expand Down Expand Up @@ -472,6 +477,7 @@ def execute(args):
'vector_path_map': alignment_source_vector_paths,
'target_pixel_size': (resolution, -resolution),
'target_srs_wkt': target_srs_wkt,
'all_touched_vectors': habitat_stressor_vectors,
},
task_name='Align raster stack',
target_path_list=list(itertools.chain(
Expand Down Expand Up @@ -1395,7 +1401,7 @@ def _write_info_json(subregion_rasters_list, subregion_ids_to_names):


def _align(raster_path_map, vector_path_map, target_pixel_size,
target_srs_wkt):
target_srs_wkt, all_touched_vectors=None):
"""Align a stack of rasters and/or vectors.

In HRA, habitats and stressors (and optionally criteria) may be defined as
Expand Down Expand Up @@ -1429,6 +1435,9 @@ def _align(raster_path_map, vector_path_map, target_pixel_size,
target_pixel_size (tuple): The pixel size of the target rasters, in
the form (x, y), expressed in meters.
target_srs_wkt (string): The target SRS of the aligned rasters.
all_touched_vectors=None (set): A set of vector paths found in
``vector_path_map`` that should be rasterized with
``ALL_TOUCHED=TRUE``.

Returns:
``None``
Expand Down Expand Up @@ -1485,6 +1494,8 @@ def _align(raster_path_map, vector_path_map, target_pixel_size,
# that align with the bounding box we determined earlier.
# This approach yields more precise rasters than resampling an
# already-rasterized vector through align_and_resize_raster_stack.
if all_touched_vectors is None:
all_touched_vectors = set([])
if vector_path_map:
LOGGER.info(f'Aligning {len(vector_path_map)} vectors')
for source_vector_path, target_raster_path in vector_path_map.items():
Expand All @@ -1495,7 +1506,12 @@ def _align(raster_path_map, vector_path_map, target_pixel_size,
raster_type = _TARGET_GDAL_TYPE_BYTE
nodata_value = _TARGET_NODATA_BYTE
burn_values = [1]
rasterize_option_list = ['ALL_TOUCHED=TRUE']
rasterize_option_list = []

# Only rasterize with ALL_TOUCHED=TRUE if we know it should be
# rasterized as such (habitats/stressors)
if source_vector_path in all_touched_vectors:
rasterize_option_list.append('ALL_TOUCHED=TRUE')

for field in layer.schema:
fieldname = field.GetName()
Expand Down
12 changes: 7 additions & 5 deletions tests/test_hra.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -608,7 +608,8 @@ def test_align(self):
self.workspace_dir, 'aligned_criterion_vector.tif'),
}

hra._align(raster_path_map, vector_path_map, (30, -30), SRS_WKT)
hra._align(raster_path_map, vector_path_map, (30, -30), SRS_WKT,
all_touched_vectors=set([habitat_vector_path]))

# Calculated by hand given the above spatial inputs and
# (30, -30) pixels. All rasters should share the same extents and
Expand Down Expand Up @@ -654,15 +655,16 @@ def test_align(self):

# The aligned criterion raster should have been rasterized from the
# rating column.
# This is an ALL_TOUCHED=FALSE rasterization.
ndta = hra._TARGET_NODATA_FLOAT32
expected_criterion_array = numpy.array([
[0.12, 0.12, 0.12, 0.12, 0.12, 0.12, ndta, ndta, ndta, ndta],
[0.12, 0.12, 0.12, 0.12, 0.12, 0.12, 0.12, ndta, ndta, ndta],
[0.12, 0.12, 0.12, 0.12, 0.12, 0.12, 0.12, ndta, ndta, ndta],
[0.12, 0.12, 0.12, 0.12, 0.12, 0.12, 0.12, ndta, ndta, ndta],
[ndta, ndta, 0.12, 0.12, 0.12, ndta, ndta, ndta, ndta, ndta],
[ndta, 0.12, 0.12, 0.12, 0.12, 0.12, ndta, ndta, ndta, ndta],
[0.12, 0.12, 0.12, 0.12, 0.12, 0.12, 0.12, ndta, ndta, ndta],
[0.12, 0.12, 0.12, 0.12, 0.12, 0.12, 0.12, ndta, ndta, ndta],
[0.12, 0.12, 0.12, 0.12, 0.12, 0.12, ndta, ndta, ndta, ndta],
[ndta, 0.12, 0.12, 0.12, 0.12, 0.12, ndta, ndta, ndta, ndta],
[ndta, ndta, 0.12, 0.12, ndta, ndta, ndta, ndta, ndta, ndta],
[ndta, ndta, ndta, ndta, ndta, ndta, ndta, ndta, ndta, ndta],
[ndta, ndta, ndta, ndta, ndta, ndta, ndta, ndta, ndta, ndta]],
dtype=numpy.float32)
Expand Down