From b7b827920537658ac5616b837eb0e75839318341 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?R=C3=A9mi=20BRAUN?= Date: Fri, 23 Jul 2021 18:53:58 +0200 Subject: [PATCH] - ENH: Allowing `stack` to take single band in input instead of a list - FIX: Fixing a regression loading optical bands which have been previously cleaned (Landsat, Theia, possibly PlanetScope) - FIX: `load` and `stack` always returns `float32` arrays --- CHANGES.md | 8 ++++--- CI/SCRIPTS/test_satellites.py | 22 ++++++++++++++----- eoreader/__init__.py | 2 +- eoreader/products/optical/pla_product.py | 8 ++++++- eoreader/products/optical/s2_theia_product.py | 12 +++++++--- eoreader/products/optical/s3_product.py | 2 +- 6 files changed, 39 insertions(+), 15 deletions(-) diff --git a/CHANGES.md b/CHANGES.md index d655f7e4..a9cac4fb 100644 --- a/CHANGES.md +++ b/CHANGES.md @@ -2,9 +2,11 @@ ## X.Y.Z (YYYY-MM-DD) -## 0.4.7.post1 (2021-07-23) -- Allowing `stack` to take single band in input instead of a list -- Fixing a regression loading Landsat bands which have been previously cleaned (again) +## 0.4.8 (2021-07-23) +- ENH: Allowing `stack` to take single band in input instead of a list +- FIX: Fixing a regression loading optical bands which have been previously cleaned (Landsat, Theia, possibly PlanetScope) +- FIX: `load` and `stack` always returns `float32` arrays +- CI: Testing if loading 2 times a band gives the same result ## 0.4.7.post0 (2021-07-23) diff --git a/CI/SCRIPTS/test_satellites.py b/CI/SCRIPTS/test_satellites.py index e17061f7..26c19444 100644 --- a/CI/SCRIPTS/test_satellites.py +++ b/CI/SCRIPTS/test_satellites.py @@ -3,6 +3,7 @@ import os import tempfile +import numpy as np import xarray as xr from cloudpathlib import AnyPath from geopandas import gpd @@ -213,13 +214,22 @@ def _test_core(pattern: str, prod_dir: str, possible_bands: list, debug=False): # Remove DEM tifs if existing remove_dem(prod) - # Get stack bands + # BAND TESTS LOGGER.info("Checking load and stack") # DO NOT RECOMPUTE BANDS WITH SNAP --> WAY TOO SLOW stack_bands = [ band for band in possible_bands if prod.has_band(band) ] + first_band = stack_bands[0] + + # Check that band loaded 2 times gives the same results (disregarding float uncertainties) + band_arr1 = prod.load(first_band, resolution=res)[first_band] + band_arr2 = prod.load(first_band, resolution=res)[first_band] + np.testing.assert_array_almost_equal(band_arr1, band_arr2) + assert band_arr1.dtype == np.float32 + assert band_arr2.dtype == np.float32 + # Get stack bands # Stack data ci_stack = get_ci_data_dir().joinpath( prod.condensed_name, f"{prod.condensed_name}_stack.tif" @@ -231,6 +241,7 @@ def _test_core(pattern: str, prod_dir: str, possible_bands: list, debug=False): stack = prod.stack( stack_bands, resolution=res, stack_path=curr_path ) + assert stack.dtype == np.float32 # Write to path if needed if not ci_stack.exists(): @@ -241,20 +252,19 @@ def _test_core(pattern: str, prod_dir: str, possible_bands: list, debug=False): # Load a band with the size option LOGGER.info("Checking load with size keyword") - band = stack_bands[0] ci_band = get_ci_data_dir().joinpath( prod.condensed_name, - f"{prod.condensed_name}_{band.name}_test.tif", + f"{prod.condensed_name}_{first_band.name}_test.tif", ) curr_path_band = os.path.join( - tmp_dir, f"{prod.condensed_name}_{band.name}_test.tif" + tmp_dir, f"{prod.condensed_name}_{first_band.name}_test.tif" ) if not ci_band.exists(): ci_band = curr_path_band band_arr = prod.load( - band, size=(stack.rio.width, stack.rio.height) - )[band] + first_band, size=(stack.rio.width, stack.rio.height) + )[first_band] rasters.write(band_arr, curr_path_band) assert_raster_almost_equal(curr_path_band, ci_band, decimal=4) diff --git a/eoreader/__init__.py b/eoreader/__init__.py index 97dd0c99..ce1d53c6 100644 --- a/eoreader/__init__.py +++ b/eoreader/__init__.py @@ -17,7 +17,7 @@ """ **EOReader** library """ -__version__ = "0.4.7-1" +__version__ = "0.4.8" __title__ = "eoreader" __description__ = ( "Remote-sensing opensource python library reading optical and SAR sensors, " diff --git a/eoreader/products/optical/pla_product.py b/eoreader/products/optical/pla_product.py index 7199cc79..33bc8d47 100644 --- a/eoreader/products/optical/pla_product.py +++ b/eoreader/products/optical/pla_product.py @@ -414,7 +414,13 @@ def _read_band( ) # Compute the correct radiometry of the band - band_xda = band_xda / 10000.0 + original_dtype = band_xda.encoding.get("dtype", band_xda.dtype) + if original_dtype == "uint16": + band_xda /= 10000.0 + + # Convert type if needed + if band_xda.dtype != np.float32: + band_xda = band_xda.astype(np.float32) return band_xda diff --git a/eoreader/products/optical/s2_theia_product.py b/eoreader/products/optical/s2_theia_product.py index dfa34e26..6f02ed38 100644 --- a/eoreader/products/optical/s2_theia_product.py +++ b/eoreader/products/optical/s2_theia_product.py @@ -239,10 +239,16 @@ def _read_band( # Read band band_xda = rasters.read( path, resolution=resolution, size=size, resampling=Resampling.bilinear - ).astype(np.float32) + ) + + # Compute the correct radiometry of the band (Theia product are stored into int16 bits) + original_dtype = band_xda.encoding.get("dtype", band_xda.dtype) + if original_dtype == "int16": + band_xda /= 10000.0 - # Compute the correct radiometry of the band - band_xda = band_xda / 10000.0 + # Convert type if needed + if band_xda.dtype != np.float32: + band_xda = band_xda.astype(np.float32) return band_xda diff --git a/eoreader/products/optical/s3_product.py b/eoreader/products/optical/s3_product.py index b2079828..cc15fa18 100644 --- a/eoreader/products/optical/s3_product.py +++ b/eoreader/products/optical/s3_product.py @@ -433,7 +433,7 @@ def _read_band( # Read band return rasters.read( path, resolution=resolution, size=size, resampling=Resampling.bilinear - ) + ).astype(np.float32) # pylint: disable=R0913 # R0913: Too many arguments (6/5) (too-many-arguments)