Repository: natcap/invest
Pull request #1810 (2025-02-25 19:27:16)

Merge urban cooling changes on main into release/3.15
## Description

Merge changes from PR #1791 into release 3.15 branch and resolve a merge conflict (as per auto PR #1809).

This involved updating the values in the regression test, as prior values (from changes due to (1) #1737 aligning rasters to the grid of the LULC raster, rather than the ET0 raster and (2) #1791 incorporating CC into the CC_park calculation) were incorrect. 


## Checklist
- [ ] Updated HISTORY.rst and link to any relevant issue (if these changes are user-facing)
- [ ] Updated the user's guide (if needed)
- [ ] Tested the Workbench UI (if relevant)



Patch:
diff --git a/HISTORY.rst b/HISTORY.rst
index d3483fe4a..d3a985a22 100644
--- a/HISTORY.rst
+++ b/HISTORY.rst
@@ -65,8 +65,8 @@ Unreleased Changes
     * Updated styling of the HTML report generated by the carbon model, for
       visual consistency with the Workbench (`InVEST #1732
       <https://github.com/natcap/invest/issues/1732>`_).
-    * Raster outputs that previously contained per-pixel values (e.g., t/pixel)
-      now contain per-hectare values (e.g., t/ha). (`InVEST #1270
+    * Raster outputs that previously contained per-pixel values (e.g., Mg/px)
+      now contain per-hectare values (e.g., Mg/ha). (`InVEST #1270
       <https://github.com/natcap/invest/issues/1270>`_).
 * Coastal Blue Carbon
     * The ``code`` column in the model's biophysical table input, as well as
@@ -75,31 +75,28 @@ Unreleased Changes
       for consistency with other InVEST models (`InVEST #1249
       <https://github.com/natcap/invest/issues/1249>`_).
 * Crop Production
-    * Raster outputs that previously contained per-pixel values (e.g., t/pixel)
-      now contain per-hectare values (e.g., t/ha). This change affects both
+    * Raster outputs that previously contained per-pixel values (e.g., Mg/px)
+      now contain per-hectare values (e.g., Mg/ha). This change affects both
       the Percentile and Regression models (`InVEST #1270
       <https://github.com/natcap/invest/issues/1270>`_).
 * Forest Carbon Edge Effect
-    * Raster outputs that previously contained per-pixel values (e.g., t/pixel)
-      now contain per-hectare values (e.g., t/ha). (`InVEST #1270
+    * Raster outputs that previously contained per-pixel values (e.g., Mg/px)
+      now contain per-hectare values (e.g., Mg/ha). (`InVEST #1270
       <https://github.com/natcap/invest/issues/1270>`_).
 * Habitat Quality
     * The ``lulc`` column in the sensitivity table input, and the ``lulc_code``
       column in the rarity table outputs, have been renamed ``lucode``, for
       consistency with other InVEST models (`InVEST #1249
       <https://github.com/natcap/invest/issues/1249>`_).
-    * The model now expects the maximum threat distance (``max_dist`` in the
-      threats table) to be specified in ``m`` instead of ``km`` (`InVEST #1252
-      <https://github.com/natcap/invest/issues/1252>`_).
 * NDR
     * Align rasters to the grid of the DEM raster
       (`#1488 <https://github.com/natcap/invest/issues/1488>`_).
-    * Raster outputs that previously contained per-pixel values (e.g., kg/pixel)
+    * Raster outputs that previously contained per-pixel values (e.g., kg/px)
       now contain per-hectare values (e.g., kg/ha). (`InVEST #1270
       <https://github.com/natcap/invest/issues/1270>`_).
 * SDR
-    * Raster outputs that previously contained per-pixel values (e.g., t/pixel)
-      now contain per-hectare values (e.g., t/ha). (`InVEST #1270
+    * Raster outputs that previously contained per-pixel values (e.g., Mg/px)
+      now contain per-hectare values (e.g., Mg/ha). (`InVEST #1270
       <https://github.com/natcap/invest/issues/1270>`_).
 * Urban Cooling
     * Align rasters to the grid of the LULC raster, rather than the ET0 raster
@@ -107,6 +104,9 @@ Unreleased Changes
     * Updated the documentation for the ``mean_t_air`` attribute of the
       ``buildings_with_stats.shp`` output to clarify how the value is
       calculated.  https://github.com/natcap/invest/issues/1746
+    * Fixed bug in the calculation of Cooling Capacity (CC) provided by parks,
+      where the CC Index was not being properly incorporated.
+      https://github.com/natcap/invest/issues/1726
 * Wind Energy
     * Fixed a bug that could cause the Workbench to crash when running the Wind
       Energy model with ``Taskgraph`` logging set to ``DEBUG`` (`InVEST #1497
diff --git a/src/natcap/invest/__init__.py b/src/natcap/invest/__init__.py
index 3d73a9089..2ba01b07c 100644
--- a/src/natcap/invest/__init__.py
+++ b/src/natcap/invest/__init__.py
@@ -16,8 +16,12 @@ __all__ = ['local_dir', ]
 try:
     __version__ = importlib.metadata.version('natcap.invest')
 except importlib.metadata.PackageNotFoundError:
-    # package is not installed.  Log the exception for debugging.
-    LOGGER.exception('Could not load natcap.invest version information')
+    try:
+        __version__ = importlib.metadata.version('natcap_invest')
+    except importlib.metadata.PackageNotFoundError:
+        # package is not installed.  Log the exception for debugging.
+        LOGGER.exception(
+            'Could not load natcap.invest (or natcap_invest) version information')
 
 # location of our translation message catalog directory
 LOCALE_DIR = os.path.join(
diff --git a/src/natcap/invest/habitat_quality.py b/src/natcap/invest/habitat_quality.py
index 813e94f8c..a2bcc761e 100644
--- a/src/natcap/invest/habitat_quality.py
+++ b/src/natcap/invest/habitat_quality.py
@@ -84,7 +84,7 @@ MODEL_SPEC = {
                         "corresponding column in the Sensitivity table.")},
                 "max_dist": {
                     "type": "number",
-                    "units": u.meter,
+                    "units": u.kilometer,
                     "about": gettext(
                         "The maximum distance over which each threat affects "
                         "habitat quality. The impact of each degradation "
@@ -1056,7 +1056,7 @@ def _decay_distance(dist_raster_path, max_dist, decay_type, target_path):
         dist_raster_path (string): a filepath for the raster to decay.
             The raster is expected to be a euclidean distance transform with
             values measuring distance in pixels.
-        max_dist (float): max distance of threat in meters.
+        max_dist (float): max distance of threat in KM.
         decay_type (string): a string defining which decay method to use.
             Options include: 'linear' | 'exponential'.
         target_path (string): a filepath for a float output raster.
@@ -1068,9 +1068,12 @@ def _decay_distance(dist_raster_path, max_dist, decay_type, target_path):
     threat_pixel_size = pygeoprocessing.get_raster_info(
         dist_raster_path)['pixel_size']
 
+    # convert max distance (given in KM) to meters
+    max_dist_m = max_dist * 1000
+
     # convert max distance from meters to the number of pixels that
     # represents on the raster
-    max_dist_pixel = max_dist / abs(threat_pixel_size[0])
+    max_dist_pixel = max_dist_m / abs(threat_pixel_size[0])
     LOGGER.debug(f'Max distance in pixels: {max_dist_pixel}')
 
     def linear_op(dist):
diff --git a/src/natcap/invest/urban_cooling_model.py b/src/natcap/invest/urban_cooling_model.py
index a45909cff..8c1458a7a 100644
--- a/src/natcap/invest/urban_cooling_model.py
+++ b/src/natcap/invest/urban_cooling_model.py
@@ -318,6 +318,10 @@ MODEL_SPEC = {
                     "about": "Map of cooling capacity",
                     "bands": {1: {"type": "ratio"}}
                 },
+                "cc_masked_green_areas.tif": {
+                    "about": "Cooling capacity map masked by non-green areas",
+                    "bands": {1: {"type": "ratio"}}
+                },
                 "T_air.tif": {
                     "about": "Map of air temperature with air mixing.",
                     "bands": {1: {"type": "number", "units": u.degree_Celsius}}
@@ -347,7 +351,47 @@ MODEL_SPEC = {
                         "reference of the LULC."),
                     "geometries": spec_utils.POLYGONS,
                     "fields": {}
-                }
+                },
+                "albedo.tif": {
+                    "about": "Map of albedo.",
+                    "bands": {1: {"type": "ratio"}}
+                },
+                "area_kernel.tif": {
+                    "about": "Area kernel for green area convolution.",
+                    "bands": {1: {"type": "integer"}}
+                },
+                "green_area_sum.tif": {
+                    "about": (
+                        "Map of green area pixels within a search distance "
+                        "around each pixel."
+                    ),
+                    "bands": {1: {"type": "number", "units": u.none}}
+                },
+                "kc.tif": {
+                    "about": "Map of crop coefficient.",
+                    "bands": {1: {"type": "number", "units": u.none}}
+                },
+                "lulc.tif": {
+                    "about": "Map of land use/land cover.",
+                    "bands": {1: {"type": "integer"}}
+                },
+                "ref_eto.tif": {
+                    "about": (
+                        "Map of reference evapotranspiration reprojected and "
+                        "aligned to the intersection of the AOI, ET0, and LULC."
+                        ),
+                    "bands": {1: {"type": "number", "units": u.millimeter}}
+                },
+                "shade.tif": {
+                    "about": "Map of shade.",
+                    "bands": {1: {"type": "ratio"}}
+                },
+                "cc_ref_aoi_stats.pickle": {
+                    "about": "Cooling capacity zonal statistics for aoi.",
+                },
+                "t_air_aoi_stats.pickle": {
+                    "about": "Air temperature zonal statistics for aoi.",
+                },
             }
         },
         "taskgraph_cache": spec_utils.TASKGRAPH_DIR
@@ -517,50 +561,6 @@ def execute(args):
             task_name=f'reclassify to {prop}')
         task_path_prop_map[prop] = (prop_task, prop_raster_path)
 
-    green_area_decay_kernel_distance = int(numpy.round(
-        float(args['green_area_cooling_distance']) / cell_size))
-    cc_park_raster_path = os.path.join(
-        intermediate_dir, f'cc_park{file_suffix}.tif')
-    cc_park_task = task_graph.add_task(
-        func=convolve_2d_by_exponential,
-        args=(
-            green_area_decay_kernel_distance,
-            task_path_prop_map['green_area'][1],
-            cc_park_raster_path),
-        target_path_list=[cc_park_raster_path],
-        dependent_task_list=[
-            task_path_prop_map['green_area'][0]],
-        task_name='calculate T air')
-
-    # Calculate the area of greenspace within a search radius of each pixel.
-    area_kernel_path = os.path.join(
-        intermediate_dir, f'area_kernel{file_suffix}.tif')
-    area_kernel_task = task_graph.add_task(
-        func=pygeoprocessing.kernels.dichotomous_kernel,
-        kwargs=dict(
-            target_kernel_path=area_kernel_path,
-            max_distance=green_area_decay_kernel_distance,
-            normalize=False),
-        target_path_list=[area_kernel_path],
-        task_name='area kernel')
-
-    green_area_sum_raster_path = os.path.join(
-        intermediate_dir, f'green_area_sum{file_suffix}.tif')
-    green_area_sum_task = task_graph.add_task(
-        func=pygeoprocessing.convolve_2d,
-        args=(
-            (task_path_prop_map['green_area'][1], 1),  # green area path
-            (area_kernel_path, 1),
-            green_area_sum_raster_path),
-        kwargs={
-            'working_dir': intermediate_dir,
-            'ignore_nodata_and_edges': True},
-        target_path_list=[green_area_sum_raster_path],
-        dependent_task_list=[
-            task_path_prop_map['green_area'][0],  # reclassed green area task
-            area_kernel_task],
-        task_name='calculate green area')
-
     align_task.join()
 
     cc_raster_path = os.path.join(intermediate_dir, f'cc{file_suffix}.tif')
@@ -623,6 +623,63 @@ def execute(args):
                 task_path_prop_map['building_intensity'][0]],
             task_name='calculate cc index (intensity)')
 
+    cc_masked_green_area_raster_path = os.path.join(
+        intermediate_dir, f'cc_masked_green_areas{file_suffix}.tif')
+
+    green_area_cc_task = task_graph.add_task(
+        func=pygeoprocessing.raster_calculator,
+        args=(
+            [(task_path_prop_map['green_area'][1], 1),
+             (cc_raster_path, 1)],
+            mask_cc_green_areas_op,
+            cc_masked_green_area_raster_path, gdal.GDT_Float32, TARGET_NODATA),
+        target_path_list=[cc_masked_green_area_raster_path],
+        dependent_task_list=[task_path_prop_map['green_area'][0], cc_task],
+        task_name='Compute green area cooling effect')
+
+    green_area_decay_kernel_distance = int(numpy.round(
+        float(args['green_area_cooling_distance']) / cell_size))
+    cc_park_raster_path = os.path.join(
+        intermediate_dir, f'cc_park{file_suffix}.tif')
+    cc_park_task = task_graph.add_task(
+        func=convolve_2d_by_exponential,
+        args=(
+            green_area_decay_kernel_distance,
+            cc_masked_green_area_raster_path,
+            cc_park_raster_path),
+        target_path_list=[cc_park_raster_path],
+        dependent_task_list=[green_area_cc_task],
+        task_name='calculate CC park')
+
+    # Calculate the area of greenspace within a search radius of each pixel.
+    area_kernel_path = os.path.join(
+        intermediate_dir, f'area_kernel{file_suffix}.tif')
+    area_kernel_task = task_graph.add_task(
+        func=pygeoprocessing.kernels.dichotomous_kernel,
+        kwargs=dict(
+            target_kernel_path=area_kernel_path,
+            max_distance=green_area_decay_kernel_distance,
+            normalize=False),
+        target_path_list=[area_kernel_path],
+        task_name='area kernel')
+
+    green_area_sum_raster_path = os.path.join(
+        intermediate_dir, f'green_area_sum{file_suffix}.tif')
+    green_area_sum_task = task_graph.add_task(
+        func=pygeoprocessing.convolve_2d,
+        args=(
+            (task_path_prop_map['green_area'][1], 1),  # green area path
+            (area_kernel_path, 1),
+            green_area_sum_raster_path),
+        kwargs={
+            'working_dir': intermediate_dir,
+            'ignore_nodata_and_edges': True},
+        target_path_list=[green_area_sum_raster_path],
+        dependent_task_list=[
+            task_path_prop_map['green_area'][0],  # reclassed green area task
+            area_kernel_task],
+        task_name='calculate green area')
+
     # Compute Heat Mitigation (HM) index.
     #
     # convert 2 hectares to number of pixels
@@ -1313,6 +1370,42 @@ def hm_op(cc_array, green_area_sum, cc_park_array, green_area_threshold):
     return result
 
 
+def mask_cc_green_areas_op(green_area_array, cc_array):
+    """
+    Mask out non-green areas from the cooling capacity (CC) raster.
+
+    Both `green_area_array` and `cc_array` originate from the LULC raster:
+        - `green_area_array` is a reclassification of LULC.
+        - `cc_array` is derived from raster calculations using reclassified
+           LULC arrays and ETI.
+    This function sets `cc_array` to:
+        - Its original cc value in green areas.
+        - 0 in non-green areas.
+        - Nodata if either `cc_array` or `green_area_array` are nodata
+
+    Args:
+        cc_array (numpy.ndarray): this is the raw cooling index mapped from
+            landcover values.
+        green_area_array (numpy.ndarray): this is the boolean array of green
+            areas where 1 corresponds with lulc classes that are green areas
+            and 0 represent areas that are not green.
+
+    Returns:
+        A modified `cc_array` where only green areas retain their original values, 
+        non-green areas are set to 0, and nodata values are preserved.
+
+    """
+    result = numpy.empty(cc_array.shape, dtype=numpy.float32)
+    result[:] = TARGET_NODATA
+    valid_mask = ~(pygeoprocessing.array_equals_nodata(cc_array, TARGET_NODATA) |
+                   pygeoprocessing.array_equals_nodata(green_area_array, TARGET_NODATA))
+    green_area_mask = green_area_array.astype(bool)
+    result[green_area_mask & valid_mask] = cc_array[green_area_mask & valid_mask]
+    result[~green_area_mask & valid_mask] = 0
+
+    return result
+
+
 def map_work_loss(
         work_temp_threshold_array, temperature_raster_path,
         work_loss_raster_path):
@@ -1386,7 +1479,7 @@ def convolve_2d_by_exponential(
     Args:
         decay_kernel_distance (float): radius of 1/e cutoff of decay kernel
             raster in pixels.
-        signal_rater_path (str): path to single band signal raster.
+        signal_raster_path (str): path to single band signal raster.
         target_convolve_raster_path (str): path to convolved raster.
 
     Returns:


Test patch:
diff --git a/.github/workflows/build-and-test.yml b/.github/workflows/build-and-test.yml
index 2719b2eac..376624d4d 100644
--- a/.github/workflows/build-and-test.yml
+++ b/.github/workflows/build-and-test.yml
@@ -121,7 +121,7 @@ jobs:
           python -m pip uninstall -y natcap.invest
           python -m build --wheel
           ls -la dist
-          pip install $(find dist -name "natcap.invest*.whl")
+          pip install $(find dist -name "natcap[._-]invest*.whl")
 
       - name: Run model tests
         run: make test
diff --git a/tests/test_habitat_quality.py b/tests/test_habitat_quality.py
index 46d8e38a8..fd3179797 100644
--- a/tests/test_habitat_quality.py
+++ b/tests/test_habitat_quality.py
@@ -226,9 +226,9 @@ class HabitatQualityTests(unittest.TestCase):
             open_table.write(
                 'MAX_DIST,WEIGHT,THREAT,DECAY,BASE_PATH,CUR_PATH,FUT_PATH\n')
             open_table.write(
-                '40,0.7,threat_1,linear,,threat_1_c.tif,threat_1_f.tif\n')
+                '0.04,0.7,threat_1,linear,,threat_1_c.tif,threat_1_f.tif\n')
             open_table.write(
-                '70,1.0,threat_2,exponential,,threat_2_c.tif,'
+                '0.07,1.0,threat_2,exponential,,threat_2_c.tif,'
                 'threat_2_f.tif\n')
 
         habitat_quality.execute(args)
@@ -334,9 +334,9 @@ class HabitatQualityTests(unittest.TestCase):
             open_table.write(
                 'MAX_DIST,WEIGHT,THREAT,DECAY,BASE_PATH,CUR_PATH,FUT_PATH\n')
             open_table.write(
-                '40,0.7,threat_1,linear,,threat_reprojected_1_c.tif,threat_1_f.tif\n')
+                '0.04,0.7,threat_1,linear,,threat_reprojected_1_c.tif,threat_1_f.tif\n')
             open_table.write(
-                '70,1.0,threat_2,exponential,,threat_2_c.tif,'
+                '0.07,1.0,threat_2,exponential,,threat_2_c.tif,'
                 'threat_2_f.tif\n')
 
         habitat_quality.execute(args)
@@ -395,9 +395,9 @@ class HabitatQualityTests(unittest.TestCase):
             open_table.write(
                 'MAX_DIST,WEIGHT,THREAT,DECAY,BASE_PATH,CUR_PATH,FUT_PATH\n')
             open_table.write(
-                '40,0.7,threat_1,linear,,threat_1_c.tif,threat_1_f.tif\n')
+                '0.04,0.7,threat_1,linear,,threat_1_c.tif,threat_1_f.tif\n')
             open_table.write(
-                '70,1.0,threat_2,exponential,,threat_2_c.tif,'
+                '0.07,1.0,threat_2,exponential,,threat_2_c.tif,'
                 'threat_2_f.tif\n')
 
         habitat_quality.execute(args)
@@ -457,9 +457,9 @@ class HabitatQualityTests(unittest.TestCase):
             open_table.write(
                 'MAX_DIST,WEIGHT,THREAT,DECAY,BASE_PATH,CUR_PATH,FUT_PATH\n')
             open_table.write(
-                '40,0.7,threat_1,linear,,threat_1_c.tif,threat_1_f.tif\n')
+                '0.04,0.7,threat_1,linear,,threat_1_c.tif,threat_1_f.tif\n')
             open_table.write(
-                '70,1.0,threat_2,exponential,,threat_2_c.tif,'
+                '0.07,1.0,threat_2,exponential,,threat_2_c.tif,'
                 'threat_2_f.tif\n')
 
         habitat_quality.execute(args)
@@ -517,9 +517,9 @@ class HabitatQualityTests(unittest.TestCase):
             open_table.write(
                 'MAX_DIST,WEIGHT,THREAT,DECAY,BASE_PATH,CUR_PATH,FUT_PATH\n')
             open_table.write(
-                '40,0.7,threat_1,linear,,threat_1_c.tif,threat_1_f.tif\n')
+                '0.04,0.7,threat_1,linear,,threat_1_c.tif,threat_1_f.tif\n')
             open_table.write(
-                '70,1.0,threat_2,exponential,,threat_2_c.tif,'
+                '0.07,1.0,threat_2,exponential,,threat_2_c.tif,'
                 'threat_2_f.tif\n')
 
         habitat_quality.execute(args)
@@ -571,9 +571,9 @@ class HabitatQualityTests(unittest.TestCase):
             open_table.write(
                 'MAX_DIST,WEIGHT,THREAT,DECAY,BASE_PATH,CUR_PATH,FUT_PATH\n')
             open_table.write(
-                '40,0.7,%s,linear,,1111_c.tif,1111_f.tif\n' % threatnames[0])
+                '0.04,0.7,%s,linear,,1111_c.tif,1111_f.tif\n' % threatnames[0])
             open_table.write(
-                '70,1.0,%s,exponential,,2222_c.tif,2222_f.tif\n'
+                '0.07,1.0,%s,exponential,,2222_c.tif,2222_f.tif\n'
                 % threatnames[1])
 
         args = {
@@ -655,9 +655,9 @@ class HabitatQualityTests(unittest.TestCase):
             open_table.write(
                 'MAX_DIST,WEIGHT,THREAT,DECAY,BASE_PATH,CUR_PATH,FUT_PATH\n')
             open_table.write(
-                '40,0.7,threat_1,linear,,threat_1_c.tif,threat_1_f.tif\n')
+                '0.04,0.7,threat_1,linear,,threat_1_c.tif,threat_1_f.tif\n')
             open_table.write(
-                '70,1.0,threat_2,exponential,,threat_2_c.tif,'
+                '0.07,1.0,threat_2,exponential,,threat_2_c.tif,'
                 'threat_2_f.tif\n')
 
         with self.assertRaises(KeyError):
@@ -699,7 +699,7 @@ class HabitatQualityTests(unittest.TestCase):
             open_table.write(
                 'MAX_DIST,WEIGHT,THREAT,DECAY,BASE_PATH,CUR_PATH,FUT_PATH\n')
             open_table.write(
-                '70,0.8,missing_threat,linear,,missing_threat_c.tif,'
+                '0.07,0.8,missing_threat,linear,,missing_threat_c.tif,'
                 'missing_threat_f.tif\n')
 
         with self.assertRaises(ValueError):
@@ -745,9 +745,9 @@ class HabitatQualityTests(unittest.TestCase):
             open_table.write(
                 'MAX_DIST,WEIGHT,THREAT,DECAY,BASE_PATH,CUR_PATH,FUT_PATH\n')
             open_table.write(
-                '40,0.7,threat_1,linear,,threat_1_c.tif,threat_1_f.tif\n')
+                '0.04,0.7,threat_1,linear,,threat_1_c.tif,threat_1_f.tif\n')
             open_table.write(
-                '70,1.0,threat_2,exponential,,threat_2_c.tif,'
+                '0.07,1.0,threat_2,exponential,,threat_2_c.tif,'
                 'threat_2_f.tif\n')
 
         with self.assertRaises(ValueError):
@@ -793,7 +793,7 @@ class HabitatQualityTests(unittest.TestCase):
             open_table.write(
                 '0.0,0.7,threat_1,linear,,threat_1_c.tif,threat_1_f.tif\n')
             open_table.write(
-                '70,1.0,threat_2,exponential,,threat_2_c.tif,'
+                '0.07,1.0,threat_2,exponential,,threat_2_c.tif,'
                 'threat_2_f.tif\n')
 
         with self.assertRaises(ValueError) as cm:
@@ -838,7 +838,7 @@ class HabitatQualityTests(unittest.TestCase):
             open_table.write(
                 'MAX_DIST,WEIGHT,THREAT,DECAY,BASE_PATH,CUR_PATH,FUT_PATH\n')
             open_table.write(
-                '40,0.7,threat_1,invalid,,threat_1_c.tif,threat_1_f.tif\n')
+                '0.04,0.7,threat_1,invalid,,threat_1_c.tif,threat_1_f.tif\n')
 
         with self.assertRaises(ValueError):
             habitat_quality.execute(args)
@@ -882,9 +882,9 @@ class HabitatQualityTests(unittest.TestCase):
             open_table.write(
                 'MAX_DIST,WEIGHT,THREAT,DECAY,CUR_PATH,FUT_PATH\n')
             open_table.write(
-                '40,0.7,threat_1,linear,threat_1_c.tif,threat_1_f.tif\n')
+                '0.04,0.7,threat_1,linear,threat_1_c.tif,threat_1_f.tif\n')
             open_table.write(
-                '70,1.0,threat_2,exponential,threat_2_c.tif,'
+                '0.07,1.0,threat_2,exponential,threat_2_c.tif,'
                 'threat_2_f.tif\n')
 
         try:
@@ -933,9 +933,9 @@ class HabitatQualityTests(unittest.TestCase):
             open_table.write(
                 'MAX_DIST,WEIGHT,THREAT,DECAY,BASE_PATH,CUR_PATH\n')
             open_table.write(
-                '40,0.7,threat_1,linear,,threat_1_c.tif\n')
+                '0.04,0.7,threat_1,linear,,threat_1_c.tif\n')
             open_table.write(
-                '70,1.0,threat_2,exponential,,threat_2_c.tif\n')
+                '0.07,1.0,threat_2,exponential,,threat_2_c.tif\n')
 
         try:
             habitat_quality.execute(args)
@@ -977,9 +977,9 @@ class HabitatQualityTests(unittest.TestCase):
             open_table.write(
                 'MAX_DIST,WEIGHT,THREAT,DECAY,BASE_PATH,CUR_PATH,FUT_PATH\n')
             open_table.write(
-                '40,0.7,threat_1,linear,,threat_1_c.tif,threat_1_f.tif\n')
+                '0.04,0.7,threat_1,linear,,threat_1_c.tif,threat_1_f.tif\n')
             open_table.write(
-                '70,1.0,threat_2,exponential,,threat_2_c.tif,'
+                '0.07,1.0,threat_2,exponential,,threat_2_c.tif,'
                 'threat_2_f.tif\n')
 
         with self.assertRaises(ValueError) as cm:
@@ -1021,9 +1021,9 @@ class HabitatQualityTests(unittest.TestCase):
             open_table.write(
                 'MAX_DIST,WEIGHT,THREAT,DECAY,BASE_PATH,CUR_PATH,FUT_PATH\n')
             open_table.write(
-                '40,0.7,threat_1,linear,,threat_1_c.tif,threat_1_f.tif\n')
+                '0.04,0.7,threat_1,linear,,threat_1_c.tif,threat_1_f.tif\n')
             open_table.write(
-                '70,1.0,threat_2,exponential,,threat_2_c.tif,'
+                '0.07,1.0,threat_2,exponential,,threat_2_c.tif,'
                 'threat_2_f.tif\n')
 
         habitat_quality.execute(args)
@@ -1069,9 +1069,9 @@ class HabitatQualityTests(unittest.TestCase):
             open_table.write(
                 'Max_Dist,Weight,threat,Decay,BASE_PATH,cur_PATH,fut_path\n')
             open_table.write(
-                '40,0.7,threat_1,linear,,threat_1_c.tif,threat_1_f.tif\n')
+                '0.04,0.7,threat_1,linear,,threat_1_c.tif,threat_1_f.tif\n')
             open_table.write(
-                '70,1.0,threat_2,exponential,,threat_2_c.tif,'
+                '0.07,1.0,threat_2,exponential,,threat_2_c.tif,'
                 'threat_2_f.tif\n')
 
         habitat_quality.execute(args)
@@ -1114,9 +1114,9 @@ class HabitatQualityTests(unittest.TestCase):
             open_table.write(
                 'MAX_DIST,WEIGHT,THREAT,DECAY,BASE_PATH,CUR_PATH,FUT_PATH\n')
             open_table.write(
-                '40,0.7,threat_1,linear,,threat_1_c.tif,threat_1_f.tif\n')
+                '0.04,0.7,threat_1,linear,,threat_1_c.tif,threat_1_f.tif\n')
             open_table.write(
-                '70,1.0,threat_2,exponential,,threat_2_c.tif,'
+                '0.07,1.0,threat_2,exponential,,threat_2_c.tif,'
                 'threat_2_f.tif\n')
 
         habitat_quality.execute(args)
@@ -1180,9 +1180,9 @@ class HabitatQualityTests(unittest.TestCase):
             open_table.write(
                 'MAX_DIST,WEIGHT,THREAT,DECAY,BASE_PATH,CUR_PATH,FUT_PATH\n')
             open_table.write(
-                '40,0.7,threat_1,linear,,threat_1_c.tif,threat_1_f.tif\n')
+                '0.04,0.7,threat_1,linear,,threat_1_c.tif,threat_1_f.tif\n')
             open_table.write(
-                '70,1.0,threat_2,exponential,,threat_2_c.tif,'
+                '0.07,1.0,threat_2,exponential,,threat_2_c.tif,'
                 'threat_2_f.tif\n')
 
         with self.assertRaises(ValueError) as cm:
@@ -1252,9 +1252,9 @@ class HabitatQualityTests(unittest.TestCase):
             open_table.write(
                 'MAX_DIST,WEIGHT,THREAT,DECAY,BASE_PATH,CUR_PATH,FUT_PATH\n')
             open_table.write(
-                '40,0.7,threat_1,linear,,threat_1_c.tif,threat_1_f.tif\n')
+                '0.04,0.7,threat_1,linear,,threat_1_c.tif,threat_1_f.tif\n')
             open_table.write(
-                '70,1.0,threat_2,exponential,,threat_2_c.tif,'
+                '0.07,1.0,threat_2,exponential,,threat_2_c.tif,'
                 'threat_2_f.tif\n')
 
         validate_result = habitat_quality.validate(args, limit_to=None)
@@ -1302,9 +1302,9 @@ class HabitatQualityTests(unittest.TestCase):
             open_table.write(
                 'MAX_DIST,WEIGHT,THREAT,DECAY,BASE_PATH,CUR_PATH,FUT_PATH\n')
             open_table.write(
-                '40,0.7,threat_1,linear,,threat_1_c.tif,threat_1_f.tif\n')
+                '0.04,0.7,threat_1,linear,,threat_1_c.tif,threat_1_f.tif\n')
             open_table.write(
-                '70,1.0,threat_2,exponential,,threat_2_c.tif,'
+                '0.07,1.0,threat_2,exponential,,threat_2_c.tif,'
                 'threat_2_f.tif\n')
 
         args['lulc_cur_path'], args['access_vector_path'] = (
@@ -1360,9 +1360,9 @@ class HabitatQualityTests(unittest.TestCase):
             open_table.write(
                 'MAX_DIST,WEIGHT,THREAT,DECAY,BASE_PATH,CUR_PATH,FUT_PATH\n')
             open_table.write(
-                '40,0.7,threat_1,linear,,threat_1_c.tif,threat_1_f.tif\n')
+                '0.04,0.7,threat_1,linear,,threat_1_c.tif,threat_1_f.tif\n')
             open_table.write(
-                '70,1.0,threat_3,exponential,,threat_2_c.tif,'
+                '0.07,1.0,threat_3,exponential,,threat_2_c.tif,'
                 'threat_2_f.tif\n')
 
         # At least one threat header is expected, so there should be a message
@@ -1411,9 +1411,9 @@ class HabitatQualityTests(unittest.TestCase):
             open_table.write(
                 'MAX_DIST,WEIGHT,THREAT,DECAY,BASE_PATH,CUR_PATH,FUT_PATH\n')
             open_table.write(
-                '40,0.7,threat_1,linear,,threat_1_c.tif,threat_1_f.tif\n')
+                '0.04,0.7,threat_1,linear,,threat_1_c.tif,threat_1_f.tif\n')
             open_table.write(
-                '70,1.0,threat_2,exponential,,threat_2_c.tif,'
+                '0.07,1.0,threat_2,exponential,,threat_2_c.tif,'
                 'threat_2_f.tif\n')
 
         validate_result = habitat_quality.validate(args, limit_to=None)
@@ -1461,9 +1461,9 @@ class HabitatQualityTests(unittest.TestCase):
             open_table.write(
                 'MAX_DIST,WEIGHT,THREAT,DECAY,BASE_PATH,CUR_PATH,FUT_PATH\n')
             open_table.write(
-                '40,0.7,threat_1,linear,,,threat_1_f.tif\n')
+                '0.04,0.7,threat_1,linear,,,threat_1_f.tif\n')
             open_table.write(
-                '70,1.0,threat_2,exponential,,threat_2_c.tif,'
+                '0.07,1.0,threat_2,exponential,,threat_2_c.tif,'
                 'threat_2_f.tif\n')
 
         with self.assertRaises(ValueError) as cm:
@@ -1514,9 +1514,9 @@ class HabitatQualityTests(unittest.TestCase):
             open_table.write(
                 'MAX_DIST,WEIGHT,THREAT,DECAY,BASE_PATH,CUR_PATH,FUT_PATH\n')
             open_table.write(
-                '40,0.7,threat_1,linear,,threat_1_c.tif,\n')
+                '0.04,0.7,threat_1,linear,,threat_1_c.tif,\n')
             open_table.write(
-                '70,1.0,threat_2,exponential,,threat_2_c.tif,'
+                '0.07,1.0,threat_2,exponential,,threat_2_c.tif,'
                 'threat_2_f.tif\n')
 
         with self.assertRaises(ValueError) as cm:
@@ -1567,9 +1567,9 @@ class HabitatQualityTests(unittest.TestCase):
             open_table.write(
                 'MAX_DIST,WEIGHT,THREAT,DECAY,BASE_PATH,CUR_PATH,FUT_PATH\n')
             open_table.write(
-                '40,0.7,threat_1,linear,,threat_1_cur.tif,threat_1_c.tif\n')
+                '0.04,0.7,threat_1,linear,,threat_1_cur.tif,threat_1_c.tif\n')
             open_table.write(
-                '70,1.0,threat_2,exponential,,threat_2_c.tif,'
+                '0.07,1.0,threat_2,exponential,,threat_2_c.tif,'
                 'threat_2_f.tif\n')
 
         with self.assertRaises(ValueError) as cm:
@@ -1617,9 +1617,9 @@ class HabitatQualityTests(unittest.TestCase):
             open_table.write(
                 'MAX_DIST,WEIGHT,THREAT,DECAY,BASE_PATH,CUR_PATH,FUT_PATH\n')
             open_table.write(
-                '40,0.7,threat_1,linear,,,threat_1_f.tif\n')
+                '0.04,0.7,threat_1,linear,,,threat_1_f.tif\n')
             open_table.write(
-                '70,1.0,threat_2,exponential,,threat_2_c.tif,'
+                '0.07,1.0,threat_2,exponential,,threat_2_c.tif,'
                 'threat_2_f.tif\n')
 
         validate_result = habitat_quality.validate(args, limit_to=None)
@@ -1670,9 +1670,9 @@ class HabitatQualityTests(unittest.TestCase):
             open_table.write(
                 'MAX_DIST,WEIGHT,THREAT,DECAY,BASE_PATH,CUR_PATH,FUT_PATH\n')
             open_table.write(
-                '40,0.7,threat_1,linear,,threat_1_c.tif,\n')
+                '0.04,0.7,threat_1,linear,,threat_1_c.tif,\n')
             open_table.write(
-                '70,1.0,threat_2,exponential,,threat_2_c.tif,'
+                '0.07,1.0,threat_2,exponential,,threat_2_c.tif,'
                 'threat_2_f.tif\n')
 
         validate_result = habitat_quality.validate(args, limit_to=None)
@@ -1736,9 +1736,9 @@ class HabitatQualityTests(unittest.TestCase):
             open_table.write(
                 'MAX_DIST,WEIGHT,THREAT,DECAY,BASE_PATH,CUR_PATH,FUT_PATH\n')
             open_table.write(
-                '40,0.7,threat_1,linear,,threat_1_c.tif,threat_1_c.tif\n')
+                '0.04,0.7,threat_1,linear,,threat_1_c.tif,threat_1_c.tif\n')
             open_table.write(
-                '70,1.0,threat_2,exponential,,threat_2_c.tif,'
+                '0.07,1.0,threat_2,exponential,,threat_2_c.tif,'
                 'threat_2_f.tif\n')
 
         validate_result = habitat_quality.validate(args, limit_to=None)
@@ -1801,9 +1801,9 @@ class HabitatQualityTests(unittest.TestCase):
             open_table.write(
                 'MAX_DIST,WEIGHT,THREAT,DECAY,BASE_PATH,CUR_PATH,FUT_PATH\n')
             open_table.write(
-                '40,0.7,threat_1,linear,,threat_1_c.tif,threat_1_c.tif\n')
+                '0.04,0.7,threat_1,linear,,threat_1_c.tif,threat_1_c.tif\n')
             open_table.write(
-                '70,1.0,threat_2,exponential,,threat_2_c.tif,'
+                '0.07,1.0,threat_2,exponential,,threat_2_c.tif,'
                 'threat_2_f.tif\n')
 
         with self.assertRaises(ValueError) as cm:
@@ -1891,9 +1891,9 @@ class HabitatQualityTests(unittest.TestCase):
             open_table.write(
                 'MAX_DIST,WEIGHT,THREAT,DECAY,BASE_PATH,CUR_PATH,FUT_PATH\n')
             open_table.write(
-                '40,0.7,threat_1,linear,,threat_1_c.tif,threat_1_f.tif\n')
+                '0.04,0.7,threat_1,linear,,threat_1_c.tif,threat_1_f.tif\n')
             open_table.write(
-                '70,1.0,threat_2,exponential,,threat_2_c.tif,'
+                '0.07,1.0,threat_2,exponential,,threat_2_c.tif,'
                 'threat_2_f.tif\n')
 
         validate_result = habitat_quality.validate(args)
@@ -1970,9 +1970,9 @@ class HabitatQualityTests(unittest.TestCase):
             open_table.write(
                 'MAX_DIST,WEIGHT,THREAT,DECAY,BASE_PATH,CUR_PATH,FUT_PATH\n')
             open_table.write(
-                '40,0.7,threat_1,linear,,threat_1_c.tif,threat_1_f.tif\n')
+                '0.04,0.7,threat_1,linear,,threat_1_c.tif,threat_1_f.tif\n')
             open_table.write(
-                '70,1.0,threat_2,exponential,,threat_2_c.tif,'
+                '0.07,1.0,threat_2,exponential,,threat_2_c.tif,'
                 'threat_2_f.tif\n')
 
         validate_result = habitat_quality.validate(args)
@@ -2019,9 +2019,9 @@ class HabitatQualityTests(unittest.TestCase):
             open_table.write(
                 'MAX_DIST,WEIGHT,THREAT,CUR_PATH,FUT_PATH\n')
             open_table.write(
-                '40,0.7,threat_1,threat_1_cur.tif,threat_1_c.tif\n')
+                '0.04,0.7,threat_1,threat_1_cur.tif,threat_1_c.tif\n')
             open_table.write(
-                '70,1.0,threat_2,threat_2_c.tif,threat_2_f.tif\n')
+                '0.07,1.0,threat_2,threat_2_c.tif,threat_2_f.tif\n')
 
         validate_result = habitat_quality.validate(args, limit_to=None)
         expected = [(
@@ -2070,9 +2070,9 @@ class HabitatQualityTests(unittest.TestCase):
             open_table.write(
                 'MAX_DIST,WEIGHT,THREAT,DECAY,CUR_PATH,FUT_PATH\n')
             open_table.write(
-                '40,0.7,threat_1,linear,threat_1_c.tif,threat_1_f.tif\n')
+                '0.04,0.7,threat_1,linear,threat_1_c.tif,threat_1_f.tif\n')
             open_table.write(
-                '70,1.0,threat_2,exponential,threat_2_c.tif,'
+                '0.07,1.0,threat_2,exponential,threat_2_c.tif,'
                 'threat_2_f.tif\n')
 
         validate_result = habitat_quality.validate(args, limit_to=None)
@@ -2123,9 +2123,9 @@ class HabitatQualityTests(unittest.TestCase):
             open_table.write(
                 'MAX_DIST,WEIGHT,THREAT,DECAY,BASE_PATH,CUR_PATH\n')
             open_table.write(
-                '40,0.7,threat_1,linear,,threat_1_c.tif\n')
+                '0.04,0.7,threat_1,linear,,threat_1_c.tif\n')
             open_table.write(
-                '70,1.0,threat_2,exponential,,threat_2_c.tif')
+                '0.07,1.0,threat_2,exponential,,threat_2_c.tif')
 
         validate_result = habitat_quality.validate(args, limit_to=None)
         expected = [(
@@ -2172,9 +2172,9 @@ class HabitatQualityTests(unittest.TestCase):
             open_table.write(
                 'MAX_DIST,WEIGHT,THREAT,DECAY,BASE_PATH,CUR_PATH,FUT_PATH\n')
             open_table.write(
-                '40,0.7,threat_1,linear,,threat_1_c.tif,threat_1_f.tif\n')
+                '0.04,0.7,threat_1,linear,,threat_1_c.tif,threat_1_f.tif\n')
             open_table.write(
-                '70,1.0,threat_2,exponential,,threat_2_c.tif,'
+                '0.07,1.0,threat_2,exponential,,threat_2_c.tif,'
                 'threat_2_f.tif\n')
 
         with self.assertRaises(TypeError):
diff --git a/tests/test_ucm.py b/tests/test_ucm.py
index 21995d6a8..5d267da3b 100644
--- a/tests/test_ucm.py
+++ b/tests/test_ucm.py
@@ -6,13 +6,65 @@ import unittest
 
 import numpy
 import pandas
-from osgeo import gdal
+from osgeo import gdal, osr, ogr
+import pygeoprocessing
+from shapely import Polygon
 
 gdal.UseExceptions()
 REGRESSION_DATA = os.path.join(
     os.path.dirname(__file__), '..', 'data', 'invest-test-data', 'ucm')
 
 
+def make_simple_vector(path_to_shp):
+    """
+    Generate shapefile with one rectangular polygon
+    Args:
+        path_to_shp (str): path to target shapefile
+    Returns:
+        None
+    """
+    # (xmin, ymin), (xmax, ymin), (xmax, ymax), (xmin, ymax), (xmin, ymin)
+    shapely_geometry_list = [
+        Polygon([(461251, 4923195), (461501, 4923195),
+                 (461501, 4923445), (461251, 4923445),
+                 (461251, 4923195)])
+    ]
+
+    srs = osr.SpatialReference()
+    srs.ImportFromEPSG(26910)
+    projection_wkt = srs.ExportToWkt()
+
+    vector_format = "ESRI Shapefile"
+    fields = {"id": ogr.OFTReal}
+    attribute_list = [{"id": 0}]
+
+    pygeoprocessing.shapely_geometry_to_vector(shapely_geometry_list,
+                                               path_to_shp, projection_wkt,
+                                               vector_format, fields,
+                                               attribute_list)
+
+
+def make_simple_raster(base_raster_path, array):
+    """Create a raster on designated path with arbitrary values.
+    Args:
+        base_raster_path (str): the raster path for making the new raster.
+    Returns:
+        None.
+    """
+    # UTM Zone 10N
+    srs = osr.SpatialReference()
+    srs.ImportFromEPSG(26910)
+    projection_wkt = srs.ExportToWkt()
+
+    origin = (461251, 4923445)
+    pixel_size = (30, -30)
+    no_data = -1
+
+    pygeoprocessing.numpy_array_to_raster(
+        array, no_data, pixel_size, origin, projection_wkt,
+        base_raster_path)
+
+
 class UCMTests(unittest.TestCase):
     """Regression tests for InVEST Urban Cooling Model."""
 
@@ -68,10 +120,10 @@ class UCMTests(unittest.TestCase):
 
         expected_results = {
             'avg_cc': 0.221991,
-            'avg_tmp_v': 37.321395,
-            'avg_tmp_an': 2.321395,
-            'avd_eng_cn': 3525589.694141,
-            'avg_wbgt_v': 32.600393,
+            'avg_tmp_v': 37.303789,
+            'avg_tmp_an': 2.303789,
+            'avd_eng_cn': 3602851.784639,
+            'avg_wbgt_v': 32.585935,
             'avg_ltls_v': 75.000000000000000,
             'avg_hvls_v': 75.000000000000000,
         }
@@ -89,7 +141,7 @@ class UCMTests(unittest.TestCase):
 
         # Assert that the decimal value of the energy savings value is what we
         # expect.
-        expected_energy_sav = 3562738.321307
+        expected_energy_sav = 3641030.461044
 
         energy_sav = 0.0
         n_nonetype = 0
@@ -110,7 +162,7 @@ class UCMTests(unittest.TestCase):
             # Expected energy savings is an accumulated value and may differ
             # past about 4 decimal places.
             numpy.testing.assert_allclose(energy_sav, expected_energy_sav, rtol=1e-4)
-            self.assertEqual(n_nonetype, 134)
+            self.assertEqual(n_nonetype, 136)
         finally:
             buildings_layer = None
             buildings_vector = None
@@ -149,7 +201,7 @@ class UCMTests(unittest.TestCase):
             # and may differ past about 4 decimal places.
             numpy.testing.assert_allclose(energy_sav, expected_energy_sav,
                                           rtol=1e-4)
-            self.assertEqual(n_nonetype, 134)
+            self.assertEqual(n_nonetype, 136)
         finally:
             buildings_layer = None
             buildings_vector = None
@@ -405,7 +457,7 @@ class UCMTests(unittest.TestCase):
                 header='column', header_name='green_area'))]
         self.assertEqual(result, expected)
 
-    def test_do_energy_valuation_option(self):
+    def test_do_energy_valuation_option(self): 
         """UCM: test separate valuation options."""
         import natcap.invest.urban_cooling_model
         args = {
@@ -534,3 +586,109 @@ class UCMTests(unittest.TestCase):
         for path in [intermediate_building_vector_path, t_air_stats_pickle_path,
                      energy_consumption_vector_path]:
             self.assertFalse(os.path.exists(path))
+
+    def test_cc_rasters(self):
+        """Test that `execute` creates correct cooling coefficient rasters with
+        synthetic data"""
+        from natcap.invest import urban_cooling_model
+
+        args = {}
+        args['workspace_dir'] = self.workspace_dir
+        args['results_suffix'] = "_01"
+        args['t_ref'] = 23.5
+        args['lulc_raster_path'] = os.path.join(self.workspace_dir,
+                                                "lulc.tif")
+        args['ref_eto_raster_path'] = os.path.join(self.workspace_dir,
+                                                   "evapotranspiration.tif")
+        args['aoi_vector_path'] = os.path.join(self.workspace_dir, "aoi.shp")
+        args['biophysical_table_path'] = os.path.join(self.workspace_dir,
+                                                      "biophysical_table.csv")
+        args['green_area_cooling_distance'] = 90
+        args['t_air_average_radius'] = 300
+        args['uhi_max'] = 2.05
+        args['do_energy_valuation'] = False
+        args['do_productivity_valuation'] = False
+        args['avg_rel_humidity'] = ''
+        args['building_vector_path'] = os.path.join(self.workspace_dir,
+                                                    "buildings.shp")
+        args['energy_consumption_table_path'] = os.path.join(self.workspace_dir,
+                                                             "ucm_energy.csv")
+        args['cc_method'] = "factors"
+        args['cc_weight_shade'] = ''  # 0.6
+        args['cc_weight_albedo'] = ''  # 0.2
+        args['cc_weight_eti'] = ''  # 0.2
+
+        def _make_input_data(args):
+            """ Create aoi shapefile, biophysical table csv, lulc tif,
+            and evapotranspiration tif"""
+
+            make_simple_vector(args['aoi_vector_path'])
+
+            biophysical_table = pandas.DataFrame({
+                "lucode": [1, 2, 3, 4, 5],
+                "lu_desc": ["water", "forest", "grassland", "urban", "barren"],
+                "green_area": [0, 1, 1, 0, 0],
+                "kc": [1, 1.1, .9, .3, .2],
+                "albedo": [.05, .1, .2, .3, .4],
+                "shade": [0, 1, 0, 0.2, 0]
+            })
+
+            biophysical_csv_path = args['biophysical_table_path']
+            biophysical_table.to_csv(biophysical_csv_path, index=False)
+
+            lulc_array = numpy.array([
+                [2, 3, 1, 5, 5, 5],
+                [3, 3, 1, 1, 4, 5],
+                [5, 5, 4, 2, 3, 1],
+                [4, 1, 4, 2, 2, 1],
+                [1, 5, 4, 1, 1, 2]
+            ], dtype=numpy.float32)
+            make_simple_raster(args['lulc_raster_path'], lulc_array)
+
+            et_array = numpy.array([
+                [800, 799, 567, 234, 422, 422],
+                [765, 867, 765, 654, 456, 677],
+                [556, 443, 456, 265, 876, 890],
+                [433, 266, 677, 776, 900, 687],
+                [456, 832, 234, 234, 234, 554]
+            ], dtype=numpy.float32)
+            make_simple_raster(args['ref_eto_raster_path'], et_array)
+
+        _make_input_data(args)
+        urban_cooling_model.execute(args)
+
+        # This array was generated by manually running through calculations
+        # Equation: cc = 0.6 * shade + 0.2 * albedo + 0.2 * eti
+        cc_array = numpy.array(
+            [[0.815556, 0.1998, 0.136, 0.0904, 0.098756, 0.098756],
+             [0.193, 0.2134, 0.18, 0.155333, 0.2104, 0.110089],
+             [0.10471112, 0.0996889, 0.2104, 0.6847778, 0.2152, 0.20777778],
+             [0.20886667, 0.06911112, 0.22513334, 0.8096889,  0.84, 0.1626668],
+             [0.11133333, 0.11697778, 0.1956, 0.062, 0.062, 0.75542223]]
+            )
+
+        cc_tif = gdal.Open(os.path.join(args["workspace_dir"], "intermediate",
+                                        f"cc{args['results_suffix']}.tif"))
+        band_cc = cc_tif.GetRasterBand(1)
+        actual_cc = band_cc.ReadAsArray()
+
+        numpy.testing.assert_allclose(actual_cc, cc_array, atol=1e-6)
+
+        # Check CC_park
+        cc_park_tif = gdal.Open(
+            os.path.join(args["workspace_dir"], "intermediate",
+                         f"cc_park{args['results_suffix']}.tif"))
+        band = cc_park_tif.GetRasterBand(1)
+        actual_cc_park = band.ReadAsArray()
+
+        # This array was created by running `convolve_2d_by_exponential`
+        # using a manually-calculated signal raster (cc * green_area)
+        expected_cc_park = numpy.array(
+            [[0.183541, 0.160865, 0.144615, 0.137877, 0.134337, 0.133191],
+             [0.15858, 0.152832, 0.147821, 0.150366, 0.14931, 0.146379],
+             [0.138345, 0.141573, 0.153495, 0.176796, 0.176097, 0.167245],
+             [0.12765204, 0.13550324, 0.15505734, 0.18876731, 0.20118835, 0.18913657],
+             [0.12391507, 0.13292464, 0.1511357, 0.17420742, 0.1924037, 0.21190241]])
+
+        numpy.testing.assert_allclose(
+            actual_cc_park, expected_cc_park, atol=1e-6)