Resolve merge conflits

.github/requirements.txt

@@ -15,4 +15,7 @@ s3fs==2024.5.0
 geemap==0.32.0
 pip==23.3.1
 boto3==1.34.124
-cartoframes==1.2.5
+cartoframes==1.2.5
+scikit-learn==1.5.0
+overturemaps==0.6.0
+git+https://github.com/isciences/exactextract

.github/workflows/build-image.yml

@@ -0,0 +1,45 @@
+name: build-image
+on:
+  push:
+    branches:
+      - main
+    paths:
+      - container/Containerfile
+      - container/VERSION
+      - environment.yml
+jobs:
+  build-image:
+    name: build-image
+    runs-on: ubuntu-22.04
+    steps:
+    - name: Clean up Ubuntu
+      run: |
+        sudo rm -rf /usr/share/dotnet
+        sudo rm -rf /opt/ghc
+        sudo rm -rf "/usr/local/share/boost"
+        sudo rm -rf "$AGENT_TOOLSDIRECTORY"
+    - name: Checkout
+      uses: actions/checkout@v4
+      with:
+        fetch-depth: 0
+    - name: Get VERSION
+      run: echo "VERSION=$(cat VERSION)" >> $GITHUB_ENV
+    - name: Build Image
+      id: build-image
+      uses: redhat-actions/buildah-build@v2
+      with:
+        image: wri-cities-cif-environment
+        tags: latest ${{ env.VERSION }} ${{ github.sha }}
+        containerfiles: |
+          ./Containerfile
+    - name: Push image to container registry
+      id: push-image-to-registry
+      uses: redhat-actions/push-to-registry@v2
+      with:
+        image: ${{ steps.build-image.outputs.image }}
+        tags: ${{ steps.build-image.outputs.tags }}
+        registry: ghcr.io/wri
+        username: ${{ secrets.REGISTRY_USER }}
+        password: ${{ secrets.REGISTRY_PASSWORD }}
+    - name: Print image url
+      run: echo "Image pushed to ${{ steps.push-image-to-registry.outputs.registry-paths }}"

.gitignore

@@ -163,3 +163,6 @@ cython_debug/
 /keys
 keys/
 wri-gee-358d958ce7c6.json
+
+# data
+/data

README.md

@@ -24,13 +24,23 @@ The Cities Indicator Framework (CIF) is a set of Python tools to make it easier
 
 ## Dependencies
 
+There are 2 ways to install dependencies. Choose one...
+
 ### Conda
 
 `conda env create -f environment.yml`
 
+### Setuptools
+
+`python setup.py`
+NOTE: If you are using this method you may want to use something like pyenv to manage Python environments
+
 ## Credentials
 
-To run the module, you need access to Google Earth Engine.
+To run the module,
+
+  1. You need access to Google Earth Engine
+  2. Install <https://cloud.google.com/sdk/docs/install>
 
 ### Interactive development
 

city_metrix/layers/__init__.py

@@ -3,7 +3,7 @@
 from .land_surface_temperature import LandSurfaceTemperature
 from .tree_cover import TreeCover
 from .high_land_surface_temperature import HighLandSurfaceTemperature
-from .smart_cities_lulc import SmartCitiesLULC
+from .smart_surface_lulc import SmartSurfaceLULC
 from .open_street_map import OpenStreetMap, OpenStreetMapClass
 from .urban_land_use import UrbanLandUse
 from .natural_areas import NaturalAreas
@@ -14,5 +14,7 @@
 from .built_up_height import BuiltUpHeight
 from .average_net_building_height import AverageNetBuildingHeight
 from .open_buildings import OpenBuildings
-from .tree_canopy_hight import TreeCanopyHeight
+from .tree_canopy_height import TreeCanopyHeight
 from .alos_dsm import AlosDSM
+from .overture_buildings import OvertureBuildings
+from .nasa_dem import NasaDEM

city_metrix/layers/nasa_dem.py

@@ -0,0 +1,21 @@
+import ee
+import xee
+import xarray as xr
+
+from .layer import Layer, get_image_collection
+
+
+class NasaDEM(Layer):
+    def __init__(self, **kwargs):
+        super().__init__(**kwargs)
+
+    def get_data(self, bbox):
+        dataset = ee.Image("NASA/NASADEM_HGT/001")
+        nasa_dem = ee.ImageCollection(ee.ImageCollection(dataset)
+                                      .filterBounds(ee.Geometry.BBox(*bbox))
+                                      .select('elevation')
+                                      .mean()
+                                      )
+        data = get_image_collection(nasa_dem, bbox, 30, "NASA DEM").elevation
+
+        return data

city_metrix/layers/overture_buildings.py

@@ -0,0 +1,30 @@
+import geopandas as gpd
+import subprocess
+from io import StringIO
+
+from .layer import Layer
+
+
+class OvertureBuildings(Layer):
+    def __init__(self, **kwargs):
+        super().__init__(**kwargs)
+
+    def get_data(self, bbox):
+        bbox_str = ','.join(map(str, bbox))
+
+        command = [
+            "overturemaps", "download",
+            "--bbox="+bbox_str,
+            "-f", "geojson",
+            "--type=building"
+        ]
+
+        result = subprocess.run(command, capture_output=True, text=True)
+
+        if result.returncode == 0:
+            geojson_data = result.stdout
+            overture_buildings = gpd.read_file(StringIO(geojson_data))
+        else:
+            print("Error occurred:", result.stderr)
+
+        return overture_buildings

city_metrix/layers/smart_surface_lulc.py

@@ -0,0 +1,131 @@
+import xarray as xr
+import numpy as np
+import pandas as pd
+import geopandas as gpd
+from shapely.geometry import CAP_STYLE, JOIN_STYLE
+from shapely.geometry import box
+import psutil
+from exactextract import exact_extract
+import pickle
+import importlib.resources as pkg_resources
+import warnings
+warnings.filterwarnings('ignore',category=UserWarning)
+
+from .layer import Layer, get_utm_zone_epsg, create_fishnet_grid, MAX_TILE_SIZE
+from .open_street_map import OpenStreetMap, OpenStreetMapClass
+from ..models.building_classifier.building_classifier import BuildingClassifier
+
+
+class SmartSurfaceLULC(Layer):
+    def __init__(self, land_cover_class=None, **kwargs):
+        super().__init__(**kwargs)
+        self.land_cover_class = land_cover_class
+
+    def get_data(self, bbox):
+        crs = get_utm_zone_epsg(bbox)
+
+        # load building roof slope classifier
+        with pkg_resources.files('city_metrix.models.building_classifier').joinpath('building_classifier.pkl').open('rb') as f:
+            clf = pickle.load(f)
+
+        # ESA world cover
+        esa_1m = BuildingClassifier().get_data_esa_reclass(bbox, crs)
+
+        # Open space
+        open_space_osm = OpenStreetMap(osm_class=OpenStreetMapClass.OPEN_SPACE_HEAT).get_data(bbox).to_crs(crs).reset_index()
+        open_space_osm['Value'] = np.int8(10)
+
+
+        # Water
+        water_osm = OpenStreetMap(osm_class=OpenStreetMapClass.WATER).get_data(bbox).to_crs(crs).reset_index()
+        water_osm['Value'] = np.int8(20)
+
+
+        # Roads
+        roads_osm = OpenStreetMap(osm_class=OpenStreetMapClass.ROAD).get_data(bbox).to_crs(crs).reset_index()
+        if len(roads_osm) > 0:
+            roads_osm['lanes'] = pd.to_numeric(roads_osm['lanes'], errors='coerce')
+            # Get the average number of lanes per highway class
+            lanes = (roads_osm.drop(columns='geometry')
+                    .groupby('highway')
+                    # Calculate average and round up
+                    .agg(avg_lanes=('lanes', lambda x: np.ceil(np.nanmean(x)) if not np.isnan(x).all() else np.NaN))
+                    )
+            # Handle NaN values in avg_lanes
+            lanes['avg_lanes'] = lanes['avg_lanes'].fillna(2)
+
+            # Fill lanes with avg lane value when missing
+            roads_osm = roads_osm.merge(lanes, on='highway', how='left')
+            roads_osm['lanes'] = roads_osm['lanes'].fillna(roads_osm['avg_lanes'])
+
+            # Add value field (30)
+            roads_osm['Value'] = np.int8(30)
+
+            # Buffer roads by lanes * 10 ft (3.048 m)
+            # https://nacto.org/publication/urban-street-design-guide/street-design-elements/lane-width/#:~:text=wider%20lane%20widths.-,Lane%20widths%20of%2010%20feet%20are%20appropriate%20in%20urban%20areas,be%20used%20in%20each%20direction
+            # 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,
+                cap_style=CAP_STYLE.flat,
+                join_style=JOIN_STYLE.mitre),
+                axis=1
+            )
+        else:
+            # Add value field (30)
+            roads_osm['Value'] = np.int8(30)
+
+
+        # Building
+        ulu_lulc_1m = BuildingClassifier().get_data_ulu(bbox, crs, esa_1m)
+        anbh_1m = BuildingClassifier().get_data_anbh(bbox, esa_1m)
+        # get building features
+        buildings = BuildingClassifier().get_data_buildings(bbox, crs)
+        # extract ULU, ANBH, and Area_m
+        buildings['ULU'] = exact_extract(ulu_lulc_1m, buildings, ["majority"], output='pandas')['majority']
+        buildings['ANBH'] = exact_extract(anbh_1m, buildings, ["mean"], output='pandas')['mean']
+        buildings['Area_m'] = buildings.geometry.area
+        # classify buildings
+        unclassed_buildings = buildings[buildings['ULU'] == 0]
+        classed_buildings = buildings[buildings['ULU'] != 0]
+
+        if len(classed_buildings) > 0:
+            classed_buildings['Value'] = clf.predict(classed_buildings[['ULU', 'ANBH', 'Area_m']])
+            # Define conditions and choices
+            case_when_class = [
+                # "residential" & "high"
+                (classed_buildings['Value'] == 40) & (classed_buildings['ULU'] == 2),
+                # "non-residential" & "high"
+                (classed_buildings['Value'] == 40) & (classed_buildings['ULU'] == 1),
+                # "residential" & "low"
+                (classed_buildings['Value'] == 42) & (classed_buildings['ULU'] == 2),
+                # "non-residential" & "low"
+                (classed_buildings['Value'] == 42) & (classed_buildings['ULU'] == 1)
+            ]
+            case_when_value = [40, 41, 42, 43]
+            classed_buildings['Value'] = np.select(case_when_class, case_when_value, default=44)
+            unclassed_buildings['Value'] = 44
+            buildings = pd.concat([classed_buildings, unclassed_buildings])
+        else:
+            buildings['Value'] = 44
+
+
+        # Parking
+        parking_osm = OpenStreetMap(osm_class=OpenStreetMapClass.PARKING).get_data(bbox).to_crs(crs).reset_index()
+        parking_osm['Value'] = np.int8(50)
+
+
+        # combine features: open space, water, road, building, parking
+        feature_df = pd.concat([open_space_osm[['geometry','Value']], water_osm[['geometry','Value']], roads_osm[['geometry','Value']], buildings[['geometry','Value']], parking_osm[['geometry','Value']]], axis=0)
+        feature_1m = BuildingClassifier().rasterize_polygon(feature_df, esa_1m)
+
+        # Combine rasters
+        datasets = [esa_1m, feature_1m]
+        # not all raster has 'time', concatenate without 'time' dimension
+        aligned_datasets = [ds.drop_vars('time', errors='ignore') for ds in datasets]
+        # use chunk 512x512
+        aligned_datasets = [ds.chunk({'x': 512, 'y': 512}) for ds in aligned_datasets]
+        lulc = xr.concat(aligned_datasets, dim='Value').max(dim='Value')
+        lulc = lulc.compute()
+
+        return lulc

city_metrix/layers/tree_canopy_hight.py → city_metrix/layers/tree_canopy_height.py

@@ -7,9 +7,7 @@
 
 
 class TreeCanopyHeight(Layer):
-
-    name = "tree_canopy_hight"
-
+    name = "tree_canopy_height"
     NO_DATA_VALUE = 0
 
     def __init__(self, **kwargs):
@@ -20,10 +18,6 @@ def get_data(self, bbox):
         # aggregate time series into a single image
         canopy_ht = canopy_ht.reduce(ee.Reducer.mean()).rename("cover_code")
 
-
-
-
-        data = get_image_collection(ee.ImageCollection(canopy_ht), bbox, 1, "tree canopy hight")
+        data = get_image_collection(ee.ImageCollection(canopy_ht), bbox, 1, "tree canopy height")
 
         return data.cover_code
-