update to return one data layer at a time

city_metrix/layers/vegetation_water_map.py

@@ -10,13 +10,15 @@ class VegetationWaterMap(Layer):
     Attributes:
         start_date: starting date for data retrieval
         end_date: ending date for data retrieval
+        greenwater_layer: select returned layer from 'startgreenwaterIndex'/'endgreenwaterIndex'/'lossgreenwaterSlope'/'gaingreenwaterSlope'
         spatial_resolution: raster resolution in meters (see https://github.com/stac-extensions/raster)
     """
 
-    def __init__(self, start_date="2018-01-01", end_date="2022-12-31", spatial_resolution=10, **kwargs):
+    def __init__(self, start_date="2016-01-01", end_date="2022-12-31", greenwater_layer='startgreenwaterIndex', spatial_resolution=10, **kwargs):
         super().__init__(**kwargs)
         self.start_date = start_date
         self.end_date = end_date
+        self.greenwater_layer = greenwater_layer
         self.spatial_resolution = spatial_resolution
 
     def get_data(self, bbox):
@@ -115,9 +117,9 @@ def senumwrap(Im):
             wP = ee.Image(1).subtract(eeCdf(wT.abs()))
             return wP
 
-        # function to generate vegetation and water trend and change maps
 
-        def get_map_vegwaterchange(IC):
+        # function to generate vegetation and water trend and change maps
+        def get_map_vegwaterchange(IC, greenwater_layer):
             gwic2019 = AnnualIC(IC, '2019')
             gwic2020 = AnnualIC(IC, '2020')
             gwic2021 = AnnualIC(IC, '2021')
@@ -206,16 +208,20 @@ def get_map_vegwaterchange(IC):
             gorwendmask = gendmask.select('NDVI').blend(wendmask.select('NDWI'))
             greenorwaterLimitLoss = glfLimitanyyeargreenLoss.blend(wlfLimitanyyearwaterLoss)
             greenorwaterLimitGain = glfLimitanyyeargreenGain.blend(wlfLimitanyyearwaterGain)
-            combinedStartLossGain = gorwstartmask.rename('startgreenwaterIndex').addBands(gorwendmask.rename('endgreenwaterIndex')).addBands(
-                greenorwaterLimitLoss.rename('lossgreenwaterSlope')).addBands(greenorwaterLimitGain.rename('gaingreenwaterSlope'))
 
-            return combinedStartLossGain
+            if greenwater_layer == 'startgreenwaterIndex':
+                return gorwstartmask.rename('greenwater_layer')
+            elif greenwater_layer == 'endgreenwaterIndex':
+                return gorwendmask.rename('greenwater_layer')
+            elif greenwater_layer == 'lossgreenwaterSlope':
+                return greenorwaterLimitLoss.rename('greenwater_layer')
+            elif greenwater_layer == 'gaingreenwaterSlope':
+                return greenorwaterLimitGain.rename('greenwater_layer')
+
 
         s2cloudmasked = Albedo().get_masked_s2_collection(ee.Geometry.BBox(*bbox), self.start_date, self.end_date)
-        vegwatermap = get_map_vegwaterchange(s2cloudmasked)
+        vegwatermap = get_map_vegwaterchange(s2cloudmasked, self.greenwater_layer)
 
         data = get_image_collection(ee.ImageCollection(vegwatermap), bbox, self.spatial_resolution, "vegetation water map")
-        # xarray.dataset to xarray.dataarray
-        data = data.to_array()
 
-        return data
+        return data.greenwater_layer

city_metrix/metrics/__init__.py

@@ -5,3 +5,4 @@
 from .urban_open_space import urban_open_space
 from .natural_areas import natural_areas
 from .era_5_met_preprocessing import era_5_met_preprocessing
+from .vegetation_water_change import vegetation_water_change

city_metrix/metrics/vegetation_water_change.py

@@ -4,11 +4,10 @@
 
 
 def vegetation_water_change(zones: GeoDataFrame) -> GeoSeries:
-    for i in range(len(zones)):
-        veg_water_map = VegetationWaterMap().get_data(zones.iloc[[i]].total_bounds)
+    start_counts = VegetationWaterMap(greenwater_layer='startgreenwaterIndex').groupby(zones).count()
+    loss_counts = VegetationWaterMap(greenwater_layer='lossgreenwaterSlope').groupby(zones).count()
+    gain_counts = VegetationWaterMap(greenwater_layer='gaingreenwaterSlope').groupby(zones).count()
 
-    counts = vegwaterImg.select('startgreenwaterIndex').reduceRegions(collection=boundary,reducer=ee.Reducer.count().setOutputs(['greenorwater2018']),scale=30)#,tileScale=10)
-    counts = vegwaterImg.select('lossgreenwaterSlope').reduceRegions(collection=counts,reducer=ee.Reducer.count().setOutputs(['greenorwaterLoss']),scale=30)#,tileScale=10)
-    counts = vegwaterImg.select('gaingreenwaterSlope').reduceRegions(collection=counts,reducer=ee.Reducer.count().setOutputs(['greenorwaterGain']),scale=30)#,tileScale=10)
+    # TODO: layer generation and zonal stats use different spatial resolutions
 
-    return 
+    return (gain_counts - loss_counts) / start_counts

tests/test_metrics.py

@@ -49,3 +49,10 @@ def test_urban_open_space():
     expected_zone_size = ZONES.geometry.size
     actual_indicator_size = indicator.size
     assert expected_zone_size == actual_indicator_size
+
+
+def test_vegetation_water_change():
+    indicator = vegetation_water_change(ZONES)
+    expected_zone_size = ZONES.geometry.size
+    actual_indicator_size = indicator.size
+    assert expected_zone_size == actual_indicator_size