Merge branch 'main' into add_impervious_surface_layer

.github/requirements.txt

@@ -17,4 +17,4 @@ pip==23.3.1
 boto3==1.34.124
 scikit-learn==1.5.0
 overturemaps==0.6.0
-git+https://github.com/isciences/exactextract
+exactextract==0.2.0.dev252

.github/workflows/dev_ci_cd.yml

@@ -2,6 +2,7 @@ name: Dev CIF API CI/CD
 
 on:
   pull_request:
+  workflow_dispatch:
 
 permissions:
   contents: read
@@ -14,9 +15,9 @@ jobs:
         python-version: ["3.10"]
 
     steps:
-    - uses: actions/checkout@v3
+    - uses: actions/checkout@v4
     - name: Set up Python ${{ matrix.python-version }}
-      uses: actions/setup-python@v3
+      uses: actions/setup-python@v5
       with:
         python-version: ${{ matrix.python-version }}
     - name: Install Linux dependencies

README.md

@@ -11,6 +11,8 @@ The Cities Indicator Framework (CIF) is a set of Python tools to make it easier
 * `pip install git+https://github.com/wri/cities-cif/releases/latest` gives you the latest stable release.
 * `pip install git+https://github.com/wri/cities-cif` gives you the main branch with is not stable.
 
+NOTE: If you have already installed the package and want to update to the latest code you may need to add the `--force-reinstall` flag
+
 ## PR Review
 
 0. Prerequisites
@@ -28,7 +30,7 @@ There are 2 ways to install dependencies. Choose one...
 
 ### Conda
 
-`conda env create -f environment.yml`
+`conda env create -f environment.yml` or `conda env update -f environment.yml`
 
 ### Setuptools
 

city_metrix/layers/__init__.py

@@ -1,4 +1,5 @@
 from .albedo import Albedo
+from .ndvi_sentinel2_gee import NdviSentinel2
 from .esa_world_cover import EsaWorldCover, EsaWorldCoverClass
 from .land_surface_temperature import LandSurfaceTemperature
 from .tree_cover import TreeCover

city_metrix/layers/albedo.py

@@ -13,7 +13,7 @@ def __init__(self, start_date="2021-01-01", end_date="2022-01-01", threshold=Non
         self.threshold = threshold
 
     def get_data(self, bbox):
-        S2 = ee.ImageCollection("COPERNICUS/S2_SR")
+        S2 = ee.ImageCollection("COPERNICUS/S2_SR_HARMONIZED")
         S2C = ee.ImageCollection("COPERNICUS/S2_CLOUD_PROBABILITY")
 
         MAX_CLOUD_PROB = 30

city_metrix/layers/high_land_surface_temperature.py

@@ -54,8 +54,3 @@ def addDate(image):
 
         # convert to date object
         return datetime.datetime.strptime(hottest_date, "%Y%m%d").date()
-
-    def write(self, output_path):
-        self.data.rio.to_raster(output_path)
-
-

city_metrix/layers/land_surface_temperature.py

@@ -33,8 +33,3 @@ def apply_scale_factors(image):
 
         data = get_image_collection(ee.ImageCollection(l8_st), bbox, 30, "LST").ST_B10_mean
         return data
-
-    def write(self, output_path):
-        self.data.rio.to_raster(output_path)
-
-

city_metrix/layers/landsat_collection_2.py

@@ -29,8 +29,7 @@ def get_data(self, bbox):
             fail_on_error=False,
         )
 
+        # TODO: Determine how to output xarray
+
         qa_lst = lc2.where((lc2.qa_pixel & 24) == 0)
         return qa_lst.drop_vars("qa_pixel")
-
-
-

city_metrix/layers/layer.py

@@ -18,10 +18,8 @@
 import shapely.geometry as geometry
 import pandas as pd
 
-
 MAX_TILE_SIZE = 0.5
 
-
 class Layer:
     def __init__(self, aggregate=None, masks=[]):
         self.aggregate = aggregate
@@ -56,7 +54,7 @@ def groupby(self, zones, layer=None):
         """
         return LayerGroupBy(self.aggregate, zones, layer, self.masks)
 
-    def write(self, bbox, output_path, tile_degrees=None):
+    def write(self, bbox, output_path, tile_degrees=None, **kwargs):
         """
         Write the layer to a path. Does not apply masks.
 
@@ -301,21 +299,23 @@ def get_image_collection(
 
     return data
 
-
 def write_layer(path, data):
     if isinstance(data, xr.DataArray):
-        # for rasters, need to write to locally first then copy to cloud storage
-        if path.startswith("s3://"):
-            tmp_path = f"{uuid4()}.tif"
-            data.rio.to_raster(raster_path=tmp_path, driver="COG")
-
-            s3 = boto3.client('s3')
-            s3.upload_file(tmp_path, path.split('/')[2], '/'.join(path.split('/')[3:]))
-
-            os.remove(tmp_path)
-        else:
-            data.rio.to_raster(raster_path=path, driver="COG")
+        write_dataarray(path, data)
     elif isinstance(data, gpd.GeoDataFrame):
         data.to_file(path, driver="GeoJSON")
     else:
-        raise NotImplementedError("Can only write DataArray or GeoDataFrame")
+        raise NotImplementedError("Can only write DataArray, Dataset, or GeoDataFrame")
+
+def write_dataarray(path, data):
+    # for rasters, need to write to locally first then copy to cloud storage
+    if path.startswith("s3://"):
+        tmp_path = f"{uuid4()}.tif"
+        data.rio.to_raster(raster_path=tmp_path, driver="COG")
+
+        s3 = boto3.client('s3')
+        s3.upload_file(tmp_path, path.split('/')[2], '/'.join(path.split('/')[3:]))
+
+        os.remove(tmp_path)
+    else:
+        data.rio.to_raster(raster_path=path, driver="COG")

city_metrix/layers/ndvi_sentinel2_gee.py

@@ -0,0 +1,46 @@
+import ee
+from .layer import Layer, get_image_collection
+
+class NdviSentinel2(Layer):
+    """"
+    NDVI = Sentinel-2 Normalized Difference Vegetation Index
+    param: year: The satellite imaging year.
+    return: a rioxarray-format DataArray
+    Author of associated Jupyter notebook: [email protected]
+    Notebook: https://github.com/wri/cities-cities4forests-indicators/blob/dev-eric/scripts/extract-VegetationCover.ipynb
+    Reference: https://en.wikipedia.org/wiki/Normalized_difference_vegetation_index
+    """
+    def __init__(self, year=None, **kwargs):
+        super().__init__(**kwargs)
+        self.year = year
+
+    def get_data(self, bbox):
+        if self.year is None:
+            raise Exception('NdviSentinel2.get_data() requires a year value')
+
+        start_date = "%s-01-01" % self.year
+        end_date = "%s-12-31" % self.year
+
+        # Compute NDVI for each image
+        def calculate_ndvi(image):
+            ndvi = (image
+                    .normalizedDifference(['B8', 'B4'])
+                    .rename('NDVI'))
+            return image.addBands(ndvi)
+
+        s2 = ee.ImageCollection("COPERNICUS/S2_HARMONIZED")
+        ndvi = (s2
+                .filterBounds(ee.Geometry.BBox(*bbox))
+                .filterDate(start_date, end_date)
+                .map(calculate_ndvi)
+                .select('NDVI')
+                 )
+
+        ndvi_mosaic = ndvi.qualityMosaic('NDVI')
+
+        ic = ee.ImageCollection(ndvi_mosaic)
+        ndvi_data = get_image_collection(ic, bbox, 10, "NDVI")
+
+        xdata = ndvi_data.to_dataarray()
+
+        return xdata

city_metrix/layers/open_street_map.py

@@ -19,6 +19,13 @@ class OpenStreetMapClass(Enum):
     BUILDING = {'building': True}
     PARKING = {'amenity': ['parking'],
                'parking': True}
+    ECONOMIC_OPPORTUNITY = {'landuse': ['commercial', 'industrial', 'retail', 'institutional', 'education'],
+							'building': ['office', 'commercial', 'industrial', 'retail', 'supermarket'],
+							'shop': True}
+    SCHOOLS = {'building': ['school',],
+				'amenity': ['school', 'kindergarten']}
+    HIGHER_EDUCATION = {'amenity': ['college', 'university'],
+						'building': ['college', 'university']}
 
 
 class OpenStreetMap(Layer):
@@ -54,3 +61,10 @@ def get_data(self, bbox):
         osm_feature = osm_feature.reset_index()[keep_col]
 
         return osm_feature
+
+    def write(self, output_path):
+        self.data['bbox'] = str(self.data.total_bounds)
+        self.data['osm_class'] = str(self.osm_class.value)
+
+        # Write to a GeoJSON file
+        self.data.to_file(output_path, driver='GeoJSON')

city_metrix/layers/sentinel_2_level_2.py

@@ -50,4 +50,6 @@ def get_data(self, bbox):
         cloud_masked = s2.where(s2 != 0).where(s2.scl != 3).where(s2.scl != 8).where(s2.scl != 9).where(
             s2.scl != 10)
 
+        # TODO: Determine how to output as an xarray
+
         return cloud_masked.drop_vars("scl")

city_metrix/layers/smart_surface_lulc.py

@@ -66,7 +66,7 @@ def get_data(self, bbox):
             # cap is flat to the terminus of the road
             # join style is mitred so intersections are squared
             roads_osm['geometry'] = roads_osm.apply(lambda row: row['geometry'].buffer(
-                row['lanes'] * 3.048,
+                row['lanes'] * 3.048 / 2,
                 cap_style=CAP_STYLE.flat,
                 join_style=JOIN_STYLE.mitre),
                 axis=1

docs/developer.md

@@ -132,7 +132,5 @@ Ensure the columns in the `GeoDataFrame` align with the [boundaries table](https
 You can run the tests by setting the credentials above and running the following:
 
 ```
-cd ./tests
-pytest layers.py
-pytest metrics.py
+pytest
 ```

environment.yml

@@ -22,6 +22,6 @@ dependencies:
   - pip=23.3.1
   - boto3=1.34.124
   - scikit-learn=1.5.0
+  - exactextract=0.2.0.dev252
   - pip:
-      - git+https://github.com/isciences/exactextract
       - overturemaps==0.6.0

setup.py

@@ -28,7 +28,7 @@
         "s3fs",
         "dask>=2023.11.0",
         "boto3",
-        "exactextract",
+        "exactextract<=0.2.0.dev252",
         "overturemaps",
         "scikit-learn>=1.5.0",
     ],

tests/resources/bbox_constants.py

@@ -0,0 +1,22 @@
+# File defines bboxes using in the test code
+
+
+BBOX_BRA_LAURO_DE_FREITAS_1 = (
+    -38.35530428121955,
+    -12.821710300686393,
+    -38.33813814352424,
+    -12.80363249765361,
+)
+
+BBOX_BRA_SALVADOR_ADM4 = (
+    -38.647320153390055,
+    -13.01748678217598787,
+    -38.3041637148564007,
+    -12.75607703449720631
+)
+
+BBOX_SMALL_TEST = (
+    -38.43864,-12.97987,
+    -38.39993,-12.93239
+)
+

tests/resources/layer_dumps_for_br_lauro_de_freitas/README.md

@@ -0,0 +1,4 @@
+# QGIS manual analysis for Lauro de Freitas, Brazil
+Folder contains:
+1. Test code that can be set to output the layers as geotiff files. Execution is controlled by a "master switch"
+1. A QGIS file used for manually inspecting the generated geotiff files

tests/resources/layer_dumps_for_br_lauro_de_freitas/conftest.py

@@ -0,0 +1,67 @@
+import tempfile
+import pytest
+import os
+import shutil
+from collections import namedtuple
+
+from tests.resources.bbox_constants import BBOX_BRA_LAURO_DE_FREITAS_1
+from tools.general_tools import create_target_folder, is_valid_path
+
+# RUN_DUMPS is the master control for whether the writes and tests are executed
+# Setting RUN_DUMPS to True turns on code execution.
+# Values should normally be set to False in order to avoid unnecessary execution.
+RUN_DUMPS = False
+
+# Specify None to write to a temporary default folder otherwise specify a valid custom target path.
+CUSTOM_DUMP_DIRECTORY = None
+
+# Both the tests and QGIS file are implemented for the same bounding box in Brazil.
+COUNTRY_CODE_FOR_BBOX = 'BRA'
+BBOX = BBOX_BRA_LAURO_DE_FREITAS_1
+
+def pytest_configure(config):
+    qgis_project_file = 'layers_for_br_lauro_de_freitas.qgz'
+
+    source_folder = os.path.dirname(__file__)
+    target_folder = get_target_folder_path()
+    create_target_folder(target_folder, True)
+
+    source_qgis_file = os.path.join(source_folder, qgis_project_file)
+    target_qgis_file = os.path.join(target_folder, qgis_project_file)
+    shutil.copyfile(source_qgis_file, target_qgis_file)
+
+    print("\n\033[93m QGIS project file and layer files written to folder %s.\033[0m\n" % target_folder)
+
+@pytest.fixture
+def target_folder():
+    return get_target_folder_path()
+
+@pytest.fixture
+def bbox_info():
+    bbox = namedtuple('bbox', ['bounds', 'country'])
+    bbox_instance = bbox(bounds=BBOX, country=COUNTRY_CODE_FOR_BBOX)
+    return bbox_instance
+
+def get_target_folder_path():
+    if CUSTOM_DUMP_DIRECTORY is not None:
+        if is_valid_path(CUSTOM_DUMP_DIRECTORY) is False:
+            raise ValueError(f"The custom path '%s' is not valid. Stopping." % CUSTOM_DUMP_DIRECTORY)
+        else:
+            output_dir = CUSTOM_DUMP_DIRECTORY
+    else:
+        sub_directory_name = 'test_result_tif_files'
+        scratch_dir_name = tempfile.TemporaryDirectory(ignore_cleanup_errors=True).name
+        dir_path = os.path.dirname(scratch_dir_name)
+        output_dir = os.path.join(dir_path, sub_directory_name)
+
+    return output_dir
+
+def prep_output_path(output_folder, file_name):
+    file_path = os.path.join(output_folder, file_name)
+    if os.path.isfile(file_path):
+        os.remove(file_path)
+    return file_path
+
+def verify_file_is_populated(file_path):
+    is_populated = True if os.path.getsize(file_path) > 0 else False
+    return is_populated