Merge pull request #1382 from phargogh/feature/1180-una-proportion-of…

…-urban-nature

Redefine Urban Nature pixels to have a proportion of Urban Nature

HISTORY.rst

@@ -186,6 +186,15 @@ Unreleased Changes
         * When defining search radii for population groups, one new output
           raster is created for each population group.  These files are named
           ``accessible_urban_nature_to_[POP_GROUP].tif``.
+
+    * Urban nature classes can now be defined to occupy a proportion of a
+      pixel, such as a park that is semi-developed. This proportion is
+      provided through user input as a proportion (0-1) in the
+      ``urban_nature`` column of the LULC Attribute Table.  A value of ``0``
+      indicates that there is no urban nature in this class, ``0.333``
+      indicates that a third of the area of this LULC class is urban nature,
+      and ``1`` would indicate that the entire LULC class's area is urban
+      nature.  https://github.com/natcap/invest/issues/1180
 * Visitation: Recreation and Tourism
     * Fixed a bug where overlapping predictor polygons would be double-counted
       in ``polygon_area_coverage`` and ``polygon_percent_coverage`` calculations.

src/natcap/invest/urban_nature_access.py

@@ -1,5 +1,4 @@
 import collections
-import functools
 import logging
 import math
 import os
@@ -80,11 +79,13 @@
             'columns': {
                 'lucode': spec_utils.LULC_TABLE_COLUMN,
                 'urban_nature': {
-                    'type': 'number',
-                    'units': u.none,
+                    'type': 'ratio',
                     'about': (
-                        "Binary code indicating whether the LULC type is "
-                        "(1) or is not (0) an urban nature type."
+                        "The proportion (0-1) indicating how much of the land "
+                        "in this LULC type is urban nature. "
+                        "0 indicates no area this LULC type is urban nature, "
+                        "1 indicates that this LULC type is entirely urban "
+                        "nature."
                     ),
                 },
                 'search_radius_m': {
@@ -648,9 +649,10 @@ def execute(args):
             CSV with the following columns:
 
             * ``lucode``: (required) the integer landcover code represented.
-            * ``urban_nature``: (required) ``0`` or ``1`` indicating whether
-              this landcover code is (``1``) or is not (``0``) an urban nature
-              pixel.
+            * ``urban_nature``: (required) a proportion (0-1) representing
+              how much of this landcover type is urban nature.  ``0``
+              indicates none of this type's area is urban nature, ``1``
+              indicates all of this type's area is urban nature.
             * ``search_radius_m``: (conditionally required) the search radius
               for this urban nature LULC class in meters. Required for all
               urban nature LULC codes if ``args['search_radius_mode'] ==
@@ -939,7 +941,7 @@ def execute(args):
     if args['search_radius_mode'] == RADIUS_OPT_UNIFORM:
         search_radii = set([float(args['search_radius'])])
     elif args['search_radius_mode'] == RADIUS_OPT_URBAN_NATURE:
-        urban_nature_attrs = attr_table[attr_table['urban_nature'] == 1]
+        urban_nature_attrs = attr_table[attr_table['urban_nature'] > 0]
         try:
             search_radii = set(urban_nature_attrs['search_radius_m'].unique())
         except KeyError as missing_key:
@@ -1802,13 +1804,16 @@ def _reclassify_urban_nature_area(
     """Reclassify LULC pixels into the urban nature area they represent.
 
     After execution, urban nature pixels will have values representing the
-    pixel's area, while pixels that are not urban nature will have a pixel
-    value of 0.  Nodata values will propagate to the output raster.
+    pixel's area of urban nature (pixel area * proportion of urban nature),
+    while pixels that are not urban nature will have a pixel value of 0.
+    Nodata values will propagate to the output raster.
 
     Args:
         lulc_raster_path (string): The path to a land-use/land-cover raster.
         lulc_attribute_table (string): The path to a CSV table representing
             LULC attributes.  Must have "lucode" and "urban_nature" columns.
+            The "urban_nature" column represents a proportion 0-1 of how much
+            of the pixel's area represents urban nature.
         target_raster_path (string): Where the reclassified urban nature raster
             should be written.
         only_these_urban_nature_codes=None (iterable or None): If ``None``, all
@@ -1830,13 +1835,15 @@ def _reclassify_urban_nature_area(
         valid_urban_nature_codes = set(only_these_urban_nature_codes)
     else:
         valid_urban_nature_codes = set(
-            lulc_attribute_df[lulc_attribute_df['urban_nature'] == 1].index)
+            lulc_attribute_df[lulc_attribute_df['urban_nature'] > 0].index)
 
     urban_nature_area_map = {}
-    for lucode in lulc_attribute_df.index:
+    for row in lulc_attribute_df[['urban_nature']].itertuples():
+        lucode = row.Index
+        urban_nature_proportion = row.urban_nature
         urban_nature_area = 0
         if lucode in valid_urban_nature_codes:
-            urban_nature_area = squared_pixel_area
+            urban_nature_area = squared_pixel_area * urban_nature_proportion
         urban_nature_area_map[lucode] = urban_nature_area
 
     lulc_raster_info = pygeoprocessing.get_raster_info(lulc_raster_path)

tests/test_urban_nature_access.py

@@ -85,7 +85,8 @@ def _build_model_args(workspace):
             6,0,100
             7,1,100
             8,0,100
-            9,1,100"""))
+            9,1,100
+            """))
 
     admin_geom = [
         shapely.geometry.box(
@@ -956,9 +957,76 @@ def test_weighted_sum(self):
         numpy.testing.assert_allclose(
             numpy.sum(weighted_sum_array[~nodata_pixels]), 1122.5)
 
+    def test_urban_nature_proportion(self):
+        """UNA: Run the model with urban nature proportion."""
+        from natcap.invest import urban_nature_access
+
+        args = _build_model_args(self.workspace_dir)
+        args['search_radius_mode'] = urban_nature_access.RADIUS_OPT_UNIFORM
+        args['search_radius'] = 1000
+        with open(args['lulc_attribute_table'], 'a') as attr_table:
+            attr_table.write("10,0.5,100\n")
+
+        # make sure our inputs validate
+        validation_results = urban_nature_access.validate(args)
+        self.assertEqual(validation_results, [])
+
+        urban_nature_access.execute(args)
+
+    def test_reclassify_urban_nature(self):
+        """UNA: Test for urban nature area reclassification."""
+        from natcap.invest import urban_nature_access
+        args = _build_model_args(self.workspace_dir)
+
+        # Rewrite the lulc attribute table to use proportions of urban nature.
+        with open(args['lulc_attribute_table'], 'w') as attr_table:
+            attr_table.write(textwrap.dedent(
+                """\
+                lucode,urban_nature,search_radius_m
+                0,0,100
+                1,0.1,100
+                2,0,100
+                3,0.3,100
+                4,0,100
+                5,0.5,100
+                6,0,100
+                7,0.7,100
+                8,0,100
+                9,0.9,100
+                """))
+
+        urban_nature_area_path = os.path.join(
+            self.workspace_dir, 'urban_nature_area.tif')
+
+        for limit_to_lucodes in (None, set([1, 3])):
+            urban_nature_access._reclassify_urban_nature_area(
+                args['lulc_raster_path'], args['lulc_attribute_table'],
+                urban_nature_area_path,
+                only_these_urban_nature_codes=limit_to_lucodes)
+
+            # The source lulc is randomized, so need to programmatically build
+            # up the expected array.
+            source_lulc_array = pygeoprocessing.raster_to_numpy_array(
+                args['lulc_raster_path'])
+            pixel_area = abs(_DEFAULT_PIXEL_SIZE[0] * _DEFAULT_PIXEL_SIZE[1])
+            expected_array = numpy.zeros(source_lulc_array.shape,
+                                         dtype=numpy.float32)
+            for i in range(1, 10, 2):
+                if limit_to_lucodes is not None:
+                    if i not in limit_to_lucodes:
+                        continue
+                factor = float(f"0.{i}")
+                expected_array[source_lulc_array == i] = factor * pixel_area
+
+            reclassified_array = pygeoprocessing.raster_to_numpy_array(
+                urban_nature_area_path)
+            numpy.testing.assert_array_almost_equal(
+                reclassified_array, expected_array)
+
     def test_validate(self):
         """UNA: Basic test for validation."""
         from natcap.invest import urban_nature_access
         args = _build_model_args(self.workspace_dir)
-        args['search_radius_mode'] = urban_nature_access.RADIUS_OPT_URBAN_NATURE
+        args['search_radius_mode'] = (
+            urban_nature_access.RADIUS_OPT_URBAN_NATURE)
         self.assertEqual(urban_nature_access.validate(args), [])