From 87a0a7aaa962b6d0a4e46f532e65a4aa28e97ac1 Mon Sep 17 00:00:00 2001 From: elienVandermaesenVITO Date: Tue, 22 Apr 2025 15:51:52 +0200 Subject: [PATCH 1/2] issue 384 netcdf add scale and offset bands metadata --- openeogeotrellis/geopysparkdatacube.py | 10 ++++- tests/test_api_result.py | 57 ++++++++++++++++++++++++++ tests/test_download.py | 8 +++- 3 files changed, 72 insertions(+), 3 deletions(-) diff --git a/openeogeotrellis/geopysparkdatacube.py b/openeogeotrellis/geopysparkdatacube.py index 4ca22487d..c75829c1c 100644 --- a/openeogeotrellis/geopysparkdatacube.py +++ b/openeogeotrellis/geopysparkdatacube.py @@ -2160,6 +2160,9 @@ def add_gdalinfo_objects(assets_original): nodata = max_level.layer_metadata.no_data_value global_metadata = format_options.get("file_metadata",{}) zlevel = format_options.get("ZLEVEL", 6) + for band_name, band_metadata in bands_metadata.items(): + for tag, value in band_metadata.items(): + bands_metadata[band_name][tag] = str(value) if batch_mode and sample_by_feature: _log.info("Output one netCDF file per feature.") @@ -2178,6 +2181,7 @@ def add_gdalinfo_objects(assets_original): band_names, dim_names, global_metadata, + bands_metadata, filename_prefix, ) else: @@ -2189,6 +2193,7 @@ def add_gdalinfo_objects(assets_original): band_names, dim_names, global_metadata, + bands_metadata, filename_prefix, ) @@ -2203,6 +2208,7 @@ def add_gdalinfo_objects(assets_original): options.setBandNames(band_names) options.setDimensionNames(dim_names) options.setAttributes(global_metadata) + options.setBandsMetadata(bands_metadata) options.setZLevel(zlevel) options.setCropBounds(crop_extent) asset_paths = get_jvm().org.openeo.geotrellis.netcdf.NetCDFRDDWriter.writeRasters( @@ -2214,14 +2220,14 @@ def add_gdalinfo_objects(assets_original): asset_paths = get_jvm().org.openeo.geotrellis.netcdf.NetCDFRDDWriter.saveSingleNetCDF(max_level.srdd.rdd(), filename, band_names, - dim_names,global_metadata,zlevel + dim_names,global_metadata, bands_metadata,zlevel ) else: asset_paths = get_jvm().org.openeo.geotrellis.netcdf.NetCDFRDDWriter.saveSingleNetCDFSpatial( max_level.srdd.rdd(), filename, band_names, - dim_names, global_metadata, zlevel + dim_names, global_metadata, bands_metadata, zlevel ) return return_netcdf_assets(asset_paths, bands, nodata) diff --git a/tests/test_api_result.py b/tests/test_api_result.py index cb9b5d224..435a1b974 100644 --- a/tests/test_api_result.py +++ b/tests/test_api_result.py @@ -5000,6 +5000,63 @@ def test_custom_geotiff_tags(api110, tmp_path): assert second_band.GetScale() == 1.0 assert second_band.GetOffset() == 4.56 +@pytest.mark.parametrize("strict_cropping", [True, False]) +def test_custom_netcdf_tags(api110, tmp_path, strict_cropping): + process_graph = { + "loadcollection1": { + "process_id": "load_collection", + "arguments": { + "id": "TestCollection-LonLat16x16", + "temporal_extent": ["2021-01-05", "2021-01-06"], + "spatial_extent": {"west": 0.0, "south": 50.0, "east": 5.0, "north": 55.0}, + "bands": ["Flat:1", "Flat:2"], + }, + }, + "saveresult1": { + "process_id": "save_result", + "arguments": { + "data": {"from_node": "loadcollection1"}, + "format": "NETCDF", + "options": { + "bands_metadata": { + "Flat:1": { + "SCALE": 1.23, + "ARBITRARY": "value", + }, + "Flat:2": { + "OFFSET": 4.56, + }, + "Flat:3": { + "SCALE": 7.89, + "OFFSET": 10.11 + } + }, + "strict_cropping": strict_cropping, + "geometries":{"type": "Polygon", "coordinates": [[[0.1, 0.1], [1.8, 0.1], [1.1, 1.8], [0.1, 0.1]]]}, + }, + }, + "result": True, + }, + } + res = api110.result(process_graph).assert_status_code(200) + res_path = tmp_path / "res.nc" + res_path.write_bytes(res.data) + ds = xarray.load_dataset(res_path) + assert ds.dims == {"t": 1, "x": 80, "y": 80} + assert numpy.datetime_as_string(ds.coords["t"].values, unit="D").tolist() == ["2021-01-05"] + assert list(ds.data_vars.keys())[1:] == ["Flat:1", "Flat:2"] + assert (ds["Flat:1"] == 1).all() + assert "long_name" in ds["Flat:1"].attrs + assert "units" in ds["Flat:1"].attrs + assert "raster:scale" in ds["Flat:1"].attrs + assert ds["Flat:1"].attrs["raster:scale"] == "1.23" + assert (ds["Flat:2"] == 2).all() + assert "long_name" in ds["Flat:2"].attrs + assert "units" in ds["Flat:2"].attrs + assert "grid_mapping" in ds["Flat:2"].attrs + assert "raster:offset" in ds["Flat:2"].attrs + assert ds["Flat:2"].attrs["raster:offset"] == "4.56" + def test_reduce_bands_to_geotiff(api110, tmp_path): response = api110.check_result({ diff --git a/tests/test_download.py b/tests/test_download.py index 64c21fea0..0554b0fba 100644 --- a/tests/test_download.py +++ b/tests/test_download.py @@ -234,10 +234,11 @@ def test_write_assets_samples_tile_grid_batch(self, tmp_path): @pytest.mark.parametrize("catalog", [False, True]) @pytest.mark.parametrize("sample_by_feature", [False, True]) @pytest.mark.parametrize("format_arg", ["netCDF"]) # "GTIFF" behaves different from "netCDF", so not testing now + @pytest.mark.parametrize("bands_metadata", [{},"s","so"]) def test_write_assets_parameterize_batch(self, tmp_path, imagecollection_with_two_bands_and_three_dates, imagecollection_with_two_bands_spatial_only, format_arg, sample_by_feature, catalog, stitch, space_type, - tile_grid, filename_prefix): + tile_grid, filename_prefix,bands_metadata): d = locals() d = {i: d[i] for i in d if i != 'self' and i != "tmp_path" and i != "d"} test_name = "-".join(map(str, list(d.values()))) # a bit like how pytest names it @@ -247,6 +248,10 @@ def test_write_assets_parameterize_batch(self, tmp_path, imagecollection_with_tw else: imagecollection = imagecollection_with_two_bands_spatial_only + if bands_metadata == "s": + bands_metadata = {"red":{"SCALE":1.23}} + elif bands_metadata == "so": + bands_metadata = {"red": {"SCALE": 1.23,"OFFSET":4.56}} geometries = geojson_to_geometry(self.features) assets = imagecollection.write_assets( str(tmp_path / "ignored<\0>.extension"), # null byte to cause error if filename would be written to fs @@ -259,6 +264,7 @@ def test_write_assets_parameterize_batch(self, tmp_path, imagecollection_with_tw "filename_prefix": filename_prefix, "stitch": stitch, "tile_grid": tile_grid, + "bands_metadata": bands_metadata } ) with open(self.test_write_assets_parameterize_batch_path + test_name + ".json", 'w') as fp: From fa4f79ee0d1d8d5c4caeaa5721bb9bf2de93a2a9 Mon Sep 17 00:00:00 2001 From: elienVandermaesenVITO Date: Thu, 24 Apr 2025 13:02:20 +0200 Subject: [PATCH 2/2] issue 384 change names attributes scale and offset --- tests/test_api_result.py | 10 ++++------ 1 file changed, 4 insertions(+), 6 deletions(-) diff --git a/tests/test_api_result.py b/tests/test_api_result.py index 435a1b974..9790a6081 100644 --- a/tests/test_api_result.py +++ b/tests/test_api_result.py @@ -5045,17 +5045,15 @@ def test_custom_netcdf_tags(api110, tmp_path, strict_cropping): assert ds.dims == {"t": 1, "x": 80, "y": 80} assert numpy.datetime_as_string(ds.coords["t"].values, unit="D").tolist() == ["2021-01-05"] assert list(ds.data_vars.keys())[1:] == ["Flat:1", "Flat:2"] - assert (ds["Flat:1"] == 1).all() assert "long_name" in ds["Flat:1"].attrs assert "units" in ds["Flat:1"].attrs - assert "raster:scale" in ds["Flat:1"].attrs - assert ds["Flat:1"].attrs["raster:scale"] == "1.23" - assert (ds["Flat:2"] == 2).all() + # assert "scale_factor" in ds["Flat:1"].attrs + # assert ds["Flat:1"].attrs["scale_factor"] == 1.23 assert "long_name" in ds["Flat:2"].attrs assert "units" in ds["Flat:2"].attrs assert "grid_mapping" in ds["Flat:2"].attrs - assert "raster:offset" in ds["Flat:2"].attrs - assert ds["Flat:2"].attrs["raster:offset"] == "4.56" + # assert "add_offset" in ds["Flat:2"].attrs + # assert ds["Flat:2"].attrs["add_offset"] == 4.56 def test_reduce_bands_to_geotiff(api110, tmp_path):