Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Added resampling method for 3 layers and updated tests #101

Open
wants to merge 3 commits into
base: main
Choose a base branch
from
Open

Added resampling method for 3 layers and updated tests #101

Show file tree
Hide file tree
Changes from 2 commits
Commits
Show all changes
3 commits
Select commit Hold shift + click to select a range
cad903b
Added resampling method for 3 layers and updated tests
kcartier-wri Dec 19, 2024
66baa49
turned off manual testing
kcartier-wri Dec 19, 2024
b8965b9
improved resampling and renamed function for clarity
kcartier-wri Dec 20, 2024
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
17 changes: 14 additions & 3 deletions city_metrix/layers/albedo.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -2,7 +2,7 @@
import xarray
from dask.diagnostics import ProgressBar

from .layer import Layer, get_utm_zone_epsg, get_image_collection
from .layer import Layer, get_utm_zone_epsg, get_image_collection, set_resampling_method


class Albedo(Layer):
Expand All @@ -11,14 +11,17 @@ class Albedo(Layer):
start_date: starting date for data retrieval
end_date: ending date for data retrieval
spatial_resolution: raster resolution in meters (see https://github.com/stac-extensions/raster)
resampling_method: interpolation method used by Google Earth Engine. Default is 'bilinear'. All options are: ('bilinear', 'bicubic', None).
threshold: threshold value for filtering the retrieval
"""

def __init__(self, start_date="2021-01-01", end_date="2022-01-01", spatial_resolution=10, threshold=None, **kwargs):
def __init__(self, start_date="2021-01-01", end_date="2022-01-01", spatial_resolution=10,
resampling_method='bilinear', threshold=None, **kwargs):
super().__init__(**kwargs)
self.start_date = start_date
self.end_date = end_date
self.spatial_resolution = spatial_resolution
self.resampling_method = resampling_method
self.threshold = threshold

def get_data(self, bbox):
Expand Down Expand Up @@ -118,7 +121,15 @@ def calc_s2_albedo(image):
## S2 MOSAIC AND ALBEDO
dataset = get_masked_s2_collection(ee.Geometry.BBox(*bbox), self.start_date, self.end_date)
s2_albedo = dataset.map(calc_s2_albedo)
albedo_mean = s2_albedo.reduce(ee.Reducer.mean())
if self.resampling_method is not None:
albedo_mean = (s2_albedo
.map(lambda x: set_resampling_method(x, self.resampling_method))
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

small nitpick to reduce code repetition, but I think because of how GEE works, you could just do this more like

resampled_albedo = s2_albedo
if self.resampling_method is not None:
            resampled_albedo = s2_albedo
                           .map(lambda x: set_resampling_method(x, self.resampling_method))

albedo_mean = (resampled_albedo.
                           .reduce(ee.Reducer.mean())
                           )

Maybe that's not any better though... your choice

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Alternatively, based on comment below, maybe you could just define the resampling in the Xee call and just pass it along to get_image_collection?

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

As far as I have determined, XEE does not yet directly support resampling. Resampling must instead be specified on the image.

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I found a way to fully implement in the called function.

.reduce(ee.Reducer.mean())
)
else:
albedo_mean = (s2_albedo
.reduce(ee.Reducer.mean())
)

albedo_mean_ic = ee.ImageCollection(albedo_mean)
data = get_image_collection(
Expand Down
26 changes: 19 additions & 7 deletions city_metrix/layers/alos_dsm.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -3,27 +3,39 @@
import xarray as xr


from .layer import Layer, get_image_collection
from .layer import Layer, get_image_collection, set_resampling_method


class AlosDSM(Layer):
"""
Attributes:
spatial_resolution: raster resolution in meters (see https://github.com/stac-extensions/raster)
resampling_method: interpolation method used by Google Earth Engine. Default is 'bilinear'. All options are: ('bilinear', 'bicubic', None).
"""

def __init__(self, spatial_resolution=30, **kwargs):
def __init__(self, spatial_resolution=30, resampling_method='bilinear', **kwargs):
super().__init__(**kwargs)
self.spatial_resolution = spatial_resolution
self.resampling_method = resampling_method

def get_data(self, bbox):
alos_dsm = ee.ImageCollection("JAXA/ALOS/AW3D30/V3_2")

alos_dsm_ic = ee.ImageCollection(alos_dsm
.filterBounds(ee.Geometry.BBox(*bbox))
.select('DSM')
.mean()
)
if self.resampling_method is not None:
alos_dsm_ic = ee.ImageCollection(
alos_dsm
.filterBounds(ee.Geometry.BBox(*bbox))
.map(lambda x: set_resampling_method(x, self.resampling_method), )
.select('DSM')
.mean()
)
else:
alos_dsm_ic = ee.ImageCollection(
alos_dsm
.filterBounds(ee.Geometry.BBox(*bbox))
.select('DSM')
.mean()
)

data = get_image_collection(
alos_dsm_ic,
Expand Down
12 changes: 12 additions & 0 deletions city_metrix/layers/layer.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -325,6 +325,18 @@ def get_stats_funcs(stats_func):
return [stats_func]


def set_resampling_method(image: ee.Image, resampling_method: str):
valid_raster_resampling_methods = ['bilinear', 'bicubic']
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Should you restrict the resampling methods on the dataset level instead of in this function? I think some datasets wouldn't want these resampling method (e.g. if they're not continuous).

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Wherever possible, I am striving to implement as much code as possible in layer functions so that the logic does not need to be separately implemented in multiple sub-classes.

As you mentioned, this function can only be used for continuous raster. I am renaming the function and adding documentation to further clarify the point.


if resampling_method not in valid_raster_resampling_methods:
raise ValueError(f"Invalid resampling method ('{resampling_method}'). "
f"Valid methods: {valid_raster_resampling_methods}")

data = image.resample(resampling_method)
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I'm not sure how this would play into, but in get_image_collection, you might accidentally trigger some sort of resampling by virtue of reading a raster at a given scale (which may be different than the native resolution of the raster). Do you need to make sure you're resampling correctly in the Xee read, or does it seem to handle it properly by calling this function here?

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Resampling needs to be specified on the images before calling XEE as far as I've found with my investigations.


return data


def get_image_collection(
image_collection: ImageCollection,
bbox: Tuple[float],
Expand Down
24 changes: 17 additions & 7 deletions city_metrix/layers/nasa_dem.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -2,27 +2,37 @@
import xee
import xarray as xr

from .layer import Layer, get_image_collection
from .layer import Layer, get_image_collection, set_resampling_method


class NasaDEM(Layer):
"""
Attributes:
spatial_resolution: raster resolution in meters (see https://github.com/stac-extensions/raster)
resampling_method: interpolation method used by Google Earth Engine. Default is 'bilinear'. All options are: ('bilinear', 'bicubic', None).
"""

def __init__(self, spatial_resolution=30, **kwargs):
def __init__(self, spatial_resolution=30, resampling_method='bilinear', **kwargs):
super().__init__(**kwargs)
self.spatial_resolution = spatial_resolution
self.resampling_method = resampling_method

def get_data(self, bbox):
nasa_dem = ee.Image("NASA/NASADEM_HGT/001")

nasa_dem_elev = (ee.ImageCollection(nasa_dem)
.filterBounds(ee.Geometry.BBox(*bbox))
.select('elevation')
.mean()
)
if self.resampling_method is not None:
nasa_dem_elev = (ee.ImageCollection(nasa_dem)
.filterBounds(ee.Geometry.BBox(*bbox))
.map(lambda x: set_resampling_method(x, self.resampling_method), )
.select('elevation')
.mean()
)
else:
nasa_dem_elev = (ee.ImageCollection(nasa_dem)
.filterBounds(ee.Geometry.BBox(*bbox))
.select('elevation')
.mean()
)

nasa_dem_elev_ic = ee.ImageCollection(nasa_dem_elev)
data = get_image_collection(
Expand Down
9 changes: 9 additions & 0 deletions tests/resources/bbox_constants.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -24,3 +24,12 @@
-38.39993,-12.93239
)

BBOX_NLD_AMSTERDAM_TEST = (
4.9012,52.3720,
4.9083,52.3752
)

BBOX_NLD_AMSTERDAM_LARGE_TEST = (
4.884629880473071,52.34146514406914,
4.914180290924863,52.359560786247165
)
0 ...r_dumps_for_br_lauro_de_freitas/README.md → tests/resources/layer_dumps/README.md
View file
Open in desktop
File renamed without changes.
0 ...dumps_for_br_lauro_de_freitas/__init__.py → tests/resources/layer_dumps/__init__.py
View file
Open in desktop
File renamed without changes.
5 changes: 4 additions & 1 deletion ...dumps_for_br_lauro_de_freitas/conftest.py → tests/resources/layer_dumps/conftest.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -4,7 +4,8 @@
import shutil
from collections import namedtuple

from tests.resources.bbox_constants import BBOX_BRA_LAURO_DE_FREITAS_1
from tests.resources.bbox_constants import BBOX_BRA_LAURO_DE_FREITAS_1, BBOX_NLD_AMSTERDAM_TEST, \
BBOX_NLD_AMSTERDAM_LARGE_TEST
from tests.tools.general_tools import create_target_folder, is_valid_path

# RUN_DUMPS is the master control for whether the writes and tests are executed
Expand All @@ -19,6 +20,8 @@
# 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
# BBOX = BBOX_NLD_AMSTERDAM_TEST
# BBOX = BBOX_NLD_AMSTERDAM_LARGE_TEST

# Specify None to write to a temporary default folder otherwise specify a valid custom target path.
CUSTOM_DUMP_DIRECTORY = None
Expand Down
0 ...lauro_de_freitas/test_write_all_layers.py → ...rces/layer_dumps/test_write_all_layers.py
View file
Open in desktop
File renamed without changes.
0 ...uro_de_freitas/test_write_layers_other.py → ...es/layer_dumps/test_write_layers_other.py
View file
Open in desktop
File renamed without changes.
12 changes: 8 additions & 4 deletions ...st_write_layers_using_fixed_resolution.py → ...st_write_layers_using_fixed_resolution.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -3,20 +3,23 @@
import pytest

from city_metrix.layers import *
from .conftest import RUN_DUMPS, prep_output_path, verify_file_is_populated, get_file_count_in_folder
from .conftest import RUN_DUMPS, prep_output_path, verify_file_is_populated

TARGET_RESOLUTION = 5
TARGET_RESAMPLING_METHOD = 'bilinear'

@pytest.mark.skipif(RUN_DUMPS == False, reason='Skipping since RUN_DUMPS set to False')
def test_write_albedo_fixed_res(target_folder, bbox_info, target_spatial_resolution_multiplier):
file_path = prep_output_path(target_folder, f'albedo_{TARGET_RESOLUTION}m.tif')
Albedo(spatial_resolution=TARGET_RESOLUTION).write(bbox_info.bounds, file_path, tile_degrees=None)
(Albedo(spatial_resolution=TARGET_RESOLUTION, resampling_method=TARGET_RESAMPLING_METHOD)
.write(bbox_info.bounds, file_path, tile_degrees=None))
assert verify_file_is_populated(file_path)

@pytest.mark.skipif(RUN_DUMPS == False, reason='Skipping since RUN_DUMPS set to False')
def test_write_alos_dsm_fixed_res(target_folder, bbox_info, target_spatial_resolution_multiplier):
file_path = prep_output_path(target_folder, f'alos_dsm_{TARGET_RESOLUTION}m.tif')
AlosDSM(spatial_resolution=TARGET_RESOLUTION).write(bbox_info.bounds, file_path, tile_degrees=None)
(AlosDSM(spatial_resolution=TARGET_RESOLUTION, resampling_method=TARGET_RESAMPLING_METHOD)
.write(bbox_info.bounds, file_path, tile_degrees=None))
assert verify_file_is_populated(file_path)

@pytest.mark.skipif(RUN_DUMPS == False, reason='Skipping since RUN_DUMPS set to False')
Expand Down Expand Up @@ -64,7 +67,8 @@ def test_write_land_surface_temperature_fixed_res(target_folder, bbox_info, targ
@pytest.mark.skipif(RUN_DUMPS == False, reason='Skipping since RUN_DUMPS set to False')
def test_write_nasa_dem_fixed_res(target_folder, bbox_info, target_spatial_resolution_multiplier):
file_path = prep_output_path(target_folder, f'nasa_dem_{TARGET_RESOLUTION}m.tif')
NasaDEM(spatial_resolution=TARGET_RESOLUTION).write(bbox_info.bounds, file_path, tile_degrees=None)
(NasaDEM(spatial_resolution=TARGET_RESOLUTION, resampling_method=TARGET_RESAMPLING_METHOD)
.write(bbox_info.bounds, file_path, tile_degrees=None))
assert verify_file_is_populated(file_path)

@pytest.mark.skipif(RUN_DUMPS == False, reason='Skipping since RUN_DUMPS set to False')
Expand Down
21 changes: 18 additions & 3 deletions tests/test_layer_metrics.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -47,8 +47,23 @@ def test_read_image_collection_scale():
pytest.approx(expected_y_size, rel=EE_IMAGE_DIMENSION_TOLERANCE) == actual_y_size
)

def test_albedo_metrics():
data = Albedo().get_data(BBOX)
def test_albedo_metrics_default_resampling():
# Default resampling_method is bilinear
data = Albedo(spatial_resolution=10).get_data(BBOX)

# Bounding values
expected_min_value = _convert_fraction_to_rounded_percent(0.03)
expected_max_value = _convert_fraction_to_rounded_percent(0.32)
actual_min_value = _convert_fraction_to_rounded_percent(data.values.min())
actual_max_value = _convert_fraction_to_rounded_percent(data.values.max())

# Value range
assert expected_min_value == actual_min_value
assert expected_max_value == actual_max_value


def test_albedo_metrics_no_resampling():
data = Albedo(spatial_resolution=10, resampling_method= None).get_data(BBOX)

# Bounding values
expected_min_value = _convert_fraction_to_rounded_percent(0.03)
Expand All @@ -62,7 +77,7 @@ def test_albedo_metrics():


def test_alos_dsm_metrics():
data = AlosDSM().get_data(BBOX)
data = AlosDSM(resampling_method=None).get_data(BBOX)

# Bounding values
expected_min_value = 16
Expand Down
Loading