modified to avoid returning null values and upgraded cdsapi version

city_metrix/layers/era_5_hottest_day.py

@@ -80,46 +80,47 @@ def hourly_mean_temperature(image):
         fc_list = []
         c = cdsapi.Client()
         for i in range(len(utc_dates)):
+            target_file = f'download_{i}.grib'
             c.retrieve(
                 'reanalysis-era5-single-levels',
                 {
                     'product_type': 'reanalysis',
                     'variable': [
-                        '10m_u_component_of_wind', 
-                        '10m_v_component_of_wind', 
+                        '10m_u_component_of_wind',
+                        '10m_v_component_of_wind',
                         '2m_dewpoint_temperature',
-                        '2m_temperature', 
-                        'clear_sky_direct_solar_radiation_at_surface', 
+                        '2m_temperature',
+                        'clear_sky_direct_solar_radiation_at_surface',
                         'mean_surface_direct_short_wave_radiation_flux_clear_sky',
-                        'mean_surface_downward_long_wave_radiation_flux_clear_sky', 
-                        'sea_surface_temperature', 
+                        'mean_surface_downward_long_wave_radiation_flux_clear_sky',
+                        'sea_surface_temperature',
                         'total_precipitation',
                     ],
                     'year': utc_dates[i].year,
                     'month': utc_dates[i].month,
                     'day': utc_dates[i].day,
                     'time': [
-                        '00:00', '01:00', '02:00', '03:00', '04:00', '05:00', 
-                        '06:00', '07:00', '08:00', '09:00', '10:00', '11:00', 
-                        '12:00', '13:00', '14:00', '15:00', '16:00', '17:00', 
+                        '00:00', '01:00', '02:00', '03:00', '04:00', '05:00',
+                        '06:00', '07:00', '08:00', '09:00', '10:00', '11:00',
+                        '12:00', '13:00', '14:00', '15:00', '16:00', '17:00',
                         '18:00', '19:00', '20:00', '21:00', '22:00', '23:00'
                     ],
                     'area': [max_lat, min_lon, min_lat, max_lon],
                     'data_format': 'grib',
                     'download_format': 'unarchived'
                 },
-                f'download_{i}.grib')
+                target_file)
 
             # {"dataType": "an"(analysis)/"fc"(forecast)/"pf"(perturbed forecast)}
-            with xr.open_dataset(f'download_{i}.grib', backend_kwargs={"filter_by_keys": {"dataType": "an"}}) as ds:
+            with xr.open_dataset(target_file, backend_kwargs={"filter_by_keys": {"dataType": "an"}}) as ds:
                 # Subset times for the day
                 times = [time.astype('datetime64[s]').astype(datetime).replace(tzinfo=pytz.UTC) for time in ds['time'].values]
                 indices = [i for i, value in enumerate(times) if value in utc_times]
                 subset_ds = ds.isel(time=indices).load()
 
             an_list.append(subset_ds)
 
-            with xr.open_dataset(f'download_{i}.grib', backend_kwargs={"filter_by_keys": {"dataType": "fc"}}) as ds:
+            with xr.open_dataset(target_file, backend_kwargs={"filter_by_keys": {"dataType": "fc"}}) as ds:
                 # reduce dimension
                 ds = ds.assign_coords(datetime=ds.time + ds.step)
                 ds = ds.stack(new_time=("time", "step"))
@@ -136,8 +137,9 @@ def hourly_mean_temperature(image):
             for file in glob.glob(f'download_{i}.grib*'):
                 os.remove(file)
 
-        an_data = xr.concat(an_list, dim='time').dropna(dim='time')
-        fc_data = xr.concat(fc_list, dim='time').dropna(dim='time')
+        an_data = xr.concat(an_list, dim='time')
+        fc_data = xr.combine_nested(fc_list, concat_dim='time').dropna(dim='time')
+
         fc_data = fc_data.sel(time=~fc_data.indexes['time'].duplicated())
         fc_data = fc_data.transpose(*an_data.dims)
 

environment.yml

@@ -23,7 +23,7 @@ dependencies:
   - geemap=0.32.0
   - pip=23.3.1
   - boto3=1.34.124
-  - cdsapi=0.7.3
+  - cdsapi=0.7.5
   - timezonefinder=6.5.2
   - scikit-image=0.24.0
   - exactextract=0.2.0

tests/conftest.py

@@ -30,9 +30,10 @@ def _create_gdf_from_coords(xmin, ymin, xmax, ymax):
 
 
 # Test zones of a regular 0.01x0.01 grid over a 0.1x0.1 extent by degrees
-ZONES = test_create_fishnet_grid(106.7, -6.3, 106.8, -6.2, 0.01).reset_index()
-LARGE_ZONES = test_create_fishnet_grid(106, -7, 107, -6, 0.1).reset_index()
+IDN_JAKARTA_TILED_ZONES = test_create_fishnet_grid(106.7, -6.3, 106.8, -6.2, 0.01).reset_index()
+LARGE_IDN_JAKARTA_TILED_ZONES = test_create_fishnet_grid(106, -7, 107, -6, 0.1).reset_index()
 OR_PORTLAND_NO_TILE_ZONE = _create_gdf_from_coords(-122.7037,45.51995,-122.6923117,45.5232773)
+NLD_AMSTERDAM_NO_TILE_ZONE = _create_gdf_from_coords(4.9012, 52.372, 4.9062057, 52.3735242)
 
 
 class MockLayer(Layer):
@@ -41,7 +42,7 @@ class MockLayer(Layer):
     """
     def get_data(self, bbox):
         arr = make_geocube(
-            vector_data=ZONES,
+            vector_data=IDN_JAKARTA_TILED_ZONES,
             measurements=['index'],
             resolution=(0.001, 0.001),
             output_crs=4326,
@@ -92,7 +93,7 @@ class MockLargeLayer(Layer):
     """
     def get_data(self, bbox):
         arr = make_geocube(
-            vector_data=LARGE_ZONES,
+            vector_data=LARGE_IDN_JAKARTA_TILED_ZONES,
             measurements=['index'],
             resolution=(0.01, 0.01),
             output_crs=4326,

tests/test_methods.py

@@ -2,8 +2,8 @@
 import pytest
 from city_metrix.layers.layer import create_fishnet_grid, offset_meters_to_geographic_degrees
 from .conftest import (
-    LARGE_ZONES,
-    ZONES,
+    LARGE_IDN_JAKARTA_TILED_ZONES,
+    IDN_JAKARTA_TILED_ZONES,
     MockGroupByLayer,
     MockLargeGroupByLayer,
     MockLargeLayer,
@@ -12,24 +12,24 @@
 )
 
 def test_count():
-    counts = MockLayer().groupby(ZONES).count()
+    counts = MockLayer().groupby(IDN_JAKARTA_TILED_ZONES).count()
     assert counts.size == 100
     assert all([count == 100 for count in counts])
 
 def test_mean():
-    means = MockLayer().groupby(ZONES).mean()
+    means = MockLayer().groupby(IDN_JAKARTA_TILED_ZONES).mean()
     assert means.size == 100
     assert all([mean == i for i, mean in enumerate(means)])
 
 
 def test_fishnetted_count():
-    counts = MockLargeLayer().groupby(LARGE_ZONES).count()
+    counts = MockLargeLayer().groupby(LARGE_IDN_JAKARTA_TILED_ZONES).count()
     assert counts.size == 100
     assert all([count == 100 for count in counts])
 
 
 def test_fishnetted_mean():
-    means = MockLargeLayer().groupby(LARGE_ZONES).mean()
+    means = MockLargeLayer().groupby(LARGE_IDN_JAKARTA_TILED_ZONES).mean()
     assert means.size == 100
     assert all([mean == i for i, mean in enumerate(means)])
 
@@ -58,7 +58,7 @@ def test_meters_to_offset_degrees():
 
 
 def test_masks():
-    counts = MockLayer().mask(MockMaskLayer()).groupby(ZONES).count()
+    counts = MockLayer().mask(MockMaskLayer()).groupby(IDN_JAKARTA_TILED_ZONES).count()
     assert counts.size == 100
     for i, count in enumerate(counts):
         if i % 2 == 0:
@@ -68,13 +68,13 @@ def test_masks():
 
 
 def test_group_by_layer():
-    counts = MockLayer().groupby(ZONES, layer=MockGroupByLayer()).count()
+    counts = MockLayer().groupby(IDN_JAKARTA_TILED_ZONES, layer=MockGroupByLayer()).count()
     assert all([count == {1: 50.0, 2: 50.0} for count in counts])
 
 
 def test_group_by_large_layer():
     counts = (
-        MockLargeLayer().groupby(LARGE_ZONES, layer=MockLargeGroupByLayer()).count()
+        MockLargeLayer().groupby(LARGE_IDN_JAKARTA_TILED_ZONES, layer=MockLargeGroupByLayer()).count()
     )
     assert all([count == {1: 50.0, 2: 50.0} for count in counts])
 

tests/test_metrics.py

@@ -1,87 +1,98 @@
 from city_metrix import *
-from .conftest import ZONES, EXECUTE_IGNORED_TESTS, OR_PORTLAND_NO_TILE_ZONE
+from .conftest import IDN_JAKARTA_TILED_ZONES, EXECUTE_IGNORED_TESTS, OR_PORTLAND_NO_TILE_ZONE, NLD_AMSTERDAM_NO_TILE_ZONE
 import pytest
 
 
 def test_built_land_with_high_lst():
-    indicator = built_land_with_high_land_surface_temperature(ZONES)
-    expected_zone_size = ZONES.geometry.size
+    indicator = built_land_with_high_land_surface_temperature(IDN_JAKARTA_TILED_ZONES)
+    expected_zone_size = IDN_JAKARTA_TILED_ZONES.geometry.size
     actual_indicator_size = indicator.size
     assert expected_zone_size == actual_indicator_size
 
 
 def test_built_land_with_low_surface_reflectivity():
-    indicator = built_land_with_low_surface_reflectivity(ZONES)
-    expected_zone_size = ZONES.geometry.size
+    indicator = built_land_with_low_surface_reflectivity(IDN_JAKARTA_TILED_ZONES)
+    expected_zone_size = IDN_JAKARTA_TILED_ZONES.geometry.size
     actual_indicator_size = indicator.size
     assert expected_zone_size == actual_indicator_size
 
 
 def test_built_land_without_tree_cover():
-    indicator = built_land_without_tree_cover(ZONES)
-    expected_zone_size = ZONES.geometry.size
+    indicator = built_land_without_tree_cover(IDN_JAKARTA_TILED_ZONES)
+    expected_zone_size = IDN_JAKARTA_TILED_ZONES.geometry.size
     actual_indicator_size = indicator.size
     assert expected_zone_size == actual_indicator_size
 
 
 @pytest.mark.skipif(EXECUTE_IGNORED_TESTS == False, reason="CDS API needs personal access token file to run")
 def test_era_5_met_preprocess_portland():
+    # Useful site: https://projects.oregonlive.com/weather/temps/
     indicator = era_5_met_preprocessing(OR_PORTLAND_NO_TILE_ZONE)
-    has_nan_values = indicator.direct_rad.isna().any()
-    assert has_nan_values == False
+    non_nullable_variables = ['temp','rh','global_rad','direct_rad','diffuse_rad','wind','vpd']
+    has_empty_required_cells = indicator[non_nullable_variables].isnull().any().any()
+    # p1= indicator[non_nullable_variables].isnull().any()
+    # p2 = indicator['global_rad'].values
+    # p3 = indicator['temp'].values
+    assert has_empty_required_cells == False
     assert len(indicator) == 24
 
 
 @pytest.mark.skipif(EXECUTE_IGNORED_TESTS == False, reason="CDS API needs personal access token file to run")
-def test_era_5_met_preprocess():
-    indicator = era_5_met_preprocessing(ZONES)
+def test_era_5_met_preprocess_amsterdam():
+    indicator = era_5_met_preprocessing(NLD_AMSTERDAM_NO_TILE_ZONE)
+    non_nullable_variables = ['temp','rh','global_rad','direct_rad','diffuse_rad','wind','vpd']
+    has_empty_required_cells = indicator[non_nullable_variables].isnull().any().any()
+    # p1= indicator[non_nullable_variables].isnull().any()
+    # p2 = indicator['global_rad'].values
+    # p3 = indicator['temp'].values
+    assert has_empty_required_cells == False
     assert len(indicator) == 24
 
 
 def test_mean_tree_cover():
-    indicator = mean_tree_cover(ZONES)
-    expected_zone_size = ZONES.geometry.size
+    indicator = mean_tree_cover(IDN_JAKARTA_TILED_ZONES)
+    expected_zone_size = IDN_JAKARTA_TILED_ZONES.geometry.size
     actual_indicator_size = indicator.size
     assert expected_zone_size == actual_indicator_size
 
 
 def test_natural_areas():
-    indicator = natural_areas(ZONES)
-    expected_zone_size = ZONES.geometry.size
+    indicator = natural_areas(IDN_JAKARTA_TILED_ZONES)
+    expected_zone_size = IDN_JAKARTA_TILED_ZONES.geometry.size
     actual_indicator_size = indicator.size
     assert expected_zone_size == actual_indicator_size
 
 
 def test_recreational_space_per_capita():
-    indicator = recreational_space_per_capita(ZONES)
-    expected_zone_size = ZONES.geometry.size
+    indicator = recreational_space_per_capita(IDN_JAKARTA_TILED_ZONES)
+    expected_zone_size = IDN_JAKARTA_TILED_ZONES.geometry.size
     actual_indicator_size = indicator.size
     assert expected_zone_size == actual_indicator_size
 
 
 def test_urban_open_space():
-    indicator = urban_open_space(ZONES)
-    expected_zone_size = ZONES.geometry.size
+    indicator = urban_open_space(IDN_JAKARTA_TILED_ZONES)
+    expected_zone_size = IDN_JAKARTA_TILED_ZONES.geometry.size
     actual_indicator_size = indicator.size
     assert expected_zone_size == actual_indicator_size
 
 
 def test_vegetation_water_change_gain_area():
-    indicator = vegetation_water_change_gain_area(ZONES)
-    expected_zone_size = ZONES.geometry.size
+    indicator = vegetation_water_change_gain_area(IDN_JAKARTA_TILED_ZONES)
+    expected_zone_size = IDN_JAKARTA_TILED_ZONES.geometry.size
     actual_indicator_size = indicator.size
     assert expected_zone_size == actual_indicator_size
 
 
 def test_vegetation_water_change_loss_area():
-    indicator = vegetation_water_change_loss_area(ZONES)
-    expected_zone_size = ZONES.geometry.size
+    indicator = vegetation_water_change_loss_area(IDN_JAKARTA_TILED_ZONES)
+    expected_zone_size = IDN_JAKARTA_TILED_ZONES.geometry.size
     actual_indicator_size = indicator.size
     assert expected_zone_size == actual_indicator_size
 
 
 def test_vegetation_water_change_gain_loss_ratio():
-    indicator = vegetation_water_change_gain_loss_ratio(ZONES)
-    expected_zone_size = ZONES.geometry.size
+    indicator = vegetation_water_change_gain_loss_ratio(IDN_JAKARTA_TILED_ZONES)
+    expected_zone_size = IDN_JAKARTA_TILED_ZONES.geometry.size
     actual_indicator_size = indicator.size
     assert expected_zone_size == actual_indicator_size