Repository: dcs4cop/xcube-sh

Issue:
Landsat 8 L2 Surface Reflectance seems to have wrong format (int instead of float)
I am using CubeGeneratorRequest to request Landsat 8 Level 2 SR data from sentinelHub. I request multiple SR bands plus the BQA quality band. The BQA band shows the correct integer values, while the SR bands are either 0 or 1. The format of the SR bands is integer. I assume the SR bands actually should be in float format. 

The source of the issue might be the documentation on the sentinelHub Website: https://docs.sentinel-hub.com/api/latest/data/landsat-8-l2/
Here it states "Source Format: Uint16". But the typical range is given as 0 - 0.4.

Here is a config I use:

`input_config = InputConfig(
                store_id="@sentinelhub_uswest",
                data_id="LOTL2",
                open_params=dict(
                    time_tolerance="30m"
                )
            )
cube_config = CubeConfig(variable_names=["B02", "B03", "B04", "B05", "B06", "B07", "B10", "BQA"],
                 bbox=bbox,
                 spatial_res=spatial_res * 3,
                 crs=crs,
                 time_range=['2017-07-01', '2017-07-16'])

output_config = OutputConfig(
            # store_id="@eurodatacube_scratch",
            store_id="s3",
            store_params={'root': "s2gmbatchtest-bc"},
            replace=True,
            data_id=prod_name
            )`

 `request = CubeGeneratorRequest(input_config=input_config, cube_config=cube_config, output_config=output_config)`

Discussion:

Patch:
diff --git a/CHANGES.md b/CHANGES.md
index 3be9a38..a75d85a 100644
--- a/CHANGES.md
+++ b/CHANGES.md
@@ -1,5 +1,10 @@
 ## Changes in 0.9.6 (in development)
 
+* Fixed Landsat radiance/reflectance bands to have 
+  sample type `FLOAT32` by default, was `UNIT16`. (#89) 
+* Avoid having `wavelength: NaN` values in band metadata attributes.
+* The open parameter `variable_names` in the `sentinelhub` data store  
+  is now optional. 
 * Fixed some import statements to be compatible with Python 3.10.
 
 ## Changes in 0.9.5
diff --git a/xcube_sh/chunkstore.py b/xcube_sh/chunkstore.py
index 842ae1a..2eafb32 100644
--- a/xcube_sh/chunkstore.py
+++ b/xcube_sh/chunkstore.py
@@ -597,7 +597,8 @@ class SentinelHubChunkStore(RemoteStore):
             time_range=(
                 time_start.strftime(datetime_format),
                 time_end.strftime(datetime_format)
-            )
+            ),
+            bad_request_ok=False
         )
         if not features:
             # Maybe the dataset has no data in given time_range
@@ -608,6 +609,7 @@ class SentinelHubChunkStore(RemoteStore):
                 collection_name=collection_name,
                 bbox=self._cube_config.bbox,
                 crs=self._cube_config.crs,
+                bad_request_ok=True
             )
             if not features:
                 # Still no results.
diff --git a/xcube_sh/metadata.py b/xcube_sh/metadata.py
index a3eeec6..a705c3c 100644
--- a/xcube_sh/metadata.py
+++ b/xcube_sh/metadata.py
@@ -328,10 +328,29 @@ S5PL2_BAND_METADATA = {
     for i in range(len(S5PL2_BAND_NAMES))
 }
 
-L8L1C_BAND_NAMES = ['B01', 'B02', 'B03', 'B04', 'B05', 'B06', 'B07',
-                    'B08', 'B09', 'B10', 'B11', 'BQA']
-L8L1C_WAVELENGTHS = [443, 482, 561.5, 654.5, 865, 1608.5, 2200.5, 589.5,
-                     1373.5, 10895, 12005, float('nan')]
+
+def _get_band_attrs(units: str, wavelength: float, resolution: int):
+    return dict(
+        sample_type='FLOAT32',
+        units=units,
+        wavelength=wavelength,
+        resolution=resolution,
+        fill_value=0.0
+    ) if wavelength > 0 else dict(
+        sample_type='UINT16',
+        units=units,
+        resolution=resolution,
+    )
+
+
+L8L1C_BAND_NAMES = ['B01', 'B02', 'B03', 'B04',
+                    'B05', 'B06', 'B07', 'B08',
+                    'B09', 'B10', 'B11',
+                    'BQA']
+L8L1C_WAVELENGTHS = [443., 482., 561.5, 654.5, 865.,
+                     1608.5, 2200.5, 589.5,
+                     1373.5, 10895., 12005.,
+                     0]
 L8L1C_RESOLUTIONS = [30, 30, 30, 30, 30, 30, 30, 15, 30, 100, 100, 30]
 L8L1C_UNITS = [
     'reflectance', 'reflectance', 'reflectance', 'reflectance', 'reflectance',
@@ -340,22 +359,24 @@ L8L1C_UNITS = [
 ]
 
 L8L1C_BAND_METADATA = {
-    L8L1C_BAND_NAMES[i]: dict(sample_type='UINT16',
-                              units=L8L1C_UNITS[i],
-                              wavelength=L8L1C_WAVELENGTHS[i],
-                              resolution=L8L1C_RESOLUTIONS[i],
-                              # fill_value=0.0
-                              )
+    L8L1C_BAND_NAMES[i]: _get_band_attrs(
+        units=L8L1C_UNITS[i],
+        wavelength=L8L1C_WAVELENGTHS[i],
+        resolution=L8L1C_RESOLUTIONS[i]
+    )
     for i in range(len(L8L1C_BAND_NAMES))
 }
 
-LOTL1_BAND_NAMES = ['B01', 'B02', 'B03', 'B04', 'B05', 'B06', 'B07', 'B08',
-                    'B09', 'B10', 'B11', 'BQA', 'QA_RADSAT',
+LOTL1_BAND_NAMES = ['B01', 'B02', 'B03', 'B04',
+                    'B05', 'B06', 'B07', 'B08',
+                    'B09', 'B10', 'B11',
+                    'BQA', 'QA_RADSAT',
                     'VAA', 'VZA', 'SAA', 'SZA']
-LOTL1_WAVELENGTHS = [443, 482, 561.5, 654.5, 865, 1608.5, 2200.5, 589.5,
-                     1373.5, 10895, 12005, float('nan'),
-                     float('nan'), float('nan'), float('nan'), float('nan'),
-                     float('nan')]
+LOTL1_WAVELENGTHS = [443., 482., 561.5, 654.5,
+                     865., 1608.5, 2200.5, 589.5,
+                     1373.5, 10895., 12005.,
+                     0, 0,
+                     0, 0, 0, 0]
 LOTL1_RESOLUTIONS = [30, 30, 30, 30, 30, 30, 30, 15, 30, 100, 100, 30, 30,
                      30, 30, 30, 30]
 LOTL1_UNITS = [
@@ -367,24 +388,24 @@ LOTL1_UNITS = [
 ]
 
 LOTL1_BAND_METADATA = {
-    LOTL1_BAND_NAMES[i]: dict(sample_type='UINT16',
-                              units=LOTL1_UNITS[i],
-                              wavelength=LOTL1_WAVELENGTHS[i],
-                              resolution=LOTL1_RESOLUTIONS[i],
-                              # fill_value=0.0
-                              )
+    LOTL1_BAND_NAMES[i]: _get_band_attrs(
+        units=LOTL1_UNITS[i],
+        wavelength=LOTL1_WAVELENGTHS[i],
+        resolution=LOTL1_RESOLUTIONS[i],
+    )
     for i in range(len(LOTL1_BAND_NAMES))
 }
 
-LOTL2_BAND_NAMES = ['B01', 'B02', 'B03', 'B04', 'B05', 'B06', 'B07',
-                    'B10', 'BQA', 'QA_RADSAT', 'SR_QA_AEROSOL',
-                    'ST_QA', 'ST_TRAD', 'ST_URAD', 'ST_DRAD', 'ST_ATRAN',
+LOTL2_BAND_NAMES = ['B01', 'B02', 'B03', 'B04',
+                    'B05', 'B06', 'B07', 'B10',
+                    'BQA', 'QA_RADSAT', 'SR_QA_AEROSOL', 'ST_QA',
+                    'ST_TRAD', 'ST_URAD', 'ST_DRAD', 'ST_ATRAN',
                     'ST_EMIS', 'ST_EMSD', 'ST_CDIST']
-LOTL2_WAVELENGTHS = [443, 482, 561.5, 654.5, 865, 1608.5, 2200.5, 10895,
-                     float('nan'), float('nan'), float('nan'),
-                     float('nan'), float('nan'), float('nan'), float('nan'),
-                     float('nan'), float('nan'), float('nan'),
-                     float('nan')]
+LOTL2_WAVELENGTHS = [443, 482, 561.5, 654.5,
+                     865., 1608.5, 2200.5, 10895.,
+                     0, 0, 0, 0,
+                     0, 0, 0, 0,
+                     0, 0, 0]
 LOTL2_RESOLUTIONS = [30, 30, 30, 30, 30, 30, 30, 100, 30, 30, 30, 30, 30,
                      30, 30, 30, 30, 30, 30]
 LOTL2_UNITS = [
@@ -396,12 +417,10 @@ LOTL2_UNITS = [
 ]
 
 LOTL2_BAND_METADATA = {
-    LOTL2_BAND_NAMES[i]: dict(
-        sample_type='UINT16',
+    LOTL2_BAND_NAMES[i]: _get_band_attrs(
         units=LOTL2_UNITS[i],
         wavelength=LOTL2_WAVELENGTHS[i],
         resolution=LOTL2_RESOLUTIONS[i],
-        # fill_value=0.0
     )
     for i in range(len(LOTL2_BAND_NAMES))
 }
@@ -419,7 +438,7 @@ MODIS_BANDWIDTHS = [50, 35, 20, 20, 20, 24, 50]
 MODIS_RESOLUTIONS = [500, 500, 500, 500, 500, 500, 500]
 MODIS_BAND_METADATA = {
     MODIS_BAND_NAMES[i]: dict(
-        sample_type='UINT16',
+        sample_type='FLOAT32',
         units='reflectance',
         wavelength=MODIS_WAVELENGTHS[i],
         bandwidth=MODIS_BANDWIDTHS[i],
diff --git a/xcube_sh/sentinelhub.py b/xcube_sh/sentinelhub.py
index 072c269..7f07b1f 100644
--- a/xcube_sh/sentinelhub.py
+++ b/xcube_sh/sentinelhub.py
@@ -98,7 +98,8 @@ class SentinelHub:
                  session: Union["SerializableOAuth2Session", Any] = None):
         if instance_id:
             warnings.warn('instance_id has been deprecated,'
-                          ' it is no longer used')
+                          ' it is no longer used',
+                          DeprecationWarning)
         self.api_url = api_url or os.environ.get('SH_API_URL',
                                                  DEFAULT_SH_API_URL)
         self.oauth2_url = oauth2_url or os.environ.get('SH_OAUTH2_URL',
@@ -195,7 +196,8 @@ class SentinelHub:
                      collection_name: str,
                      bbox: Tuple[float, float, float, float] = None,
                      crs: str = None,
-                     time_range: Tuple[str, str] = None) \
+                     time_range: Tuple[str, str] = None,
+                     bad_request_ok: bool = False) \
             -> List[Dict[str, Any]]:
         """
         Get geometric intersections of dataset given by *collection_name*
@@ -207,6 +209,8 @@ class SentinelHub:
         :param crs: Name of a coordinate reference system of the coordinates
             given by *bbox*. Ignored if *bbox* is not given.
         :param time_range: time range
+        :param bad_request_ok: return empty list rather than raise error
+            on bad request
         :return: list of features that include a "datetime" field for
             all intersections.
         """
@@ -262,6 +266,9 @@ class SentinelHub:
                                          json=request,
                                          headers=headers)
 
+            if bad_request_ok and response.status_code == 400:
+                break
+
             SentinelHubError.maybe_raise_for_response(response)
 
             feature_collection = json.loads(response.content)
@@ -273,6 +280,8 @@ class SentinelHub:
                                        response=response)
 
             features = feature_collection['features']
+            if not features:
+                break
             all_features.extend(features)
             features_count = len(features)
             feature_offset += features_count
diff --git a/xcube_sh/store.py b/xcube_sh/store.py
index aec80aa..3b4f58f 100644
--- a/xcube_sh/store.py
+++ b/xcube_sh/store.py
@@ -223,7 +223,8 @@ class SentinelHubDataOpener(DataOpener):
                 # Note, in xcube-sh < 0.9.2 formerly used an enum here.
                 # However, that doesn't work for BYOC, so we omit that
                 # limitation.
-                items=JsonStringSchema()
+                items=JsonStringSchema(),
+                nullable=True
             ),
             variable_fill_values=JsonArraySchema(
                 items=JsonNumberSchema(nullable=True)


Test patch:
diff --git a/test/test_chunkstore.py b/test/test_chunkstore.py
index 72e992d..b80877d 100644
--- a/test/test_chunkstore.py
+++ b/test/test_chunkstore.py
@@ -24,6 +24,7 @@ import zlib
 from abc import ABCMeta
 from abc import abstractmethod
 from collections import namedtuple
+from typing import Tuple, List, Dict, Any
 
 import numpy as np
 import pandas as pd
@@ -463,7 +464,13 @@ class SentinelHubMock:
         return [dict(name=b) for b in S2_BAND_NAMES]
 
     # noinspection PyUnusedLocal
-    def get_features(self, collection_name, bbox, crs, time_range):
+    def get_features(self,
+                     collection_name: str,
+                     bbox: Tuple[float, float, float, float] = None,
+                     crs: str = None,
+                     time_range: Tuple[str, str] = None,
+                     bad_request_ok: bool = False) \
+            -> List[Dict[str, Any]]:
         """Return dummy catalog features"""
         start_time, end_time = map(pd.to_datetime, time_range)
         datetime = start_time + pd.to_timedelta('8H')
diff --git a/test/test_metadata.py b/test/test_metadata.py
index 5950bb9..6e8d41c 100644
--- a/test/test_metadata.py
+++ b/test/test_metadata.py
@@ -72,9 +72,9 @@ class SentinelHubMetadataTest(unittest.TestCase):
         self.assertEqual((21, 21), self._assert_all_wavelengths_numeric('S3OLCI'))
         self.assertEqual((11, 11), self._assert_all_wavelengths_numeric('S3SLSTR'))
         self.assertEqual((0, 0), self._assert_all_wavelengths_numeric('S5PL2'))
-        self.assertEqual((12, 0), self._assert_all_wavelengths_numeric('L8L1C'))
-        self.assertEqual((17, 0), self._assert_all_wavelengths_numeric('LOTL1'))
-        self.assertEqual((19, 0), self._assert_all_wavelengths_numeric('LOTL2'))
+        self.assertEqual((11, 0), self._assert_all_wavelengths_numeric('L8L1C'))
+        self.assertEqual((11, 0), self._assert_all_wavelengths_numeric('LOTL1'))
+        self.assertEqual((8, 0), self._assert_all_wavelengths_numeric('LOTL2'))
         self.assertEqual((7, 7), self._assert_all_wavelengths_numeric('MODIS'))
 
     def _assert_all_wavelengths_numeric(self, ds_name: str) -> Tuple[int, int]:
diff --git a/test/test_reprod_89.py b/test/test_reprod_89.py
new file mode 100644
index 0000000..b8cd437
--- /dev/null
+++ b/test/test_reprod_89.py
@@ -0,0 +1,57 @@
+# The MIT License (MIT)
+# Copyright (c) 2022 by the xcube development team and contributors
+#
+# Permission is hereby granted, free of charge, to any person obtaining a
+# copy of this software and associated documentation files (the "Software"),
+# to deal in the Software without restriction, including without limitation
+# the rights to use, copy, modify, merge, publish, distribute, sublicense,
+# and/or sell copies of the Software, and to permit persons to whom the
+# Software is furnished to do so, subject to the following conditions:
+#
+# The above copyright notice and this permission notice shall be included in
+# all copies or substantial portions of the Software.
+#
+# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
+# IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
+# FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
+# AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
+# LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
+# FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
+# DEALINGS IN THE SOFTWARE.
+
+import unittest
+
+import numpy as np
+
+from xcube_sh.config import CubeConfig
+from xcube_sh.cube import open_cube
+from xcube_sh.sentinelhub import SentinelHub
+
+
+@unittest.skip(reason="should be fixed")
+class ReproduceIssue89Test(unittest.TestCase):
+    """
+    This test reproduces issue #89.
+    """
+
+    def test_reproduce_issue_89(self):
+        sentinel_hub = SentinelHub(
+            api_url='https://services-uswest2.sentinel-hub.com'
+        )
+
+        cube_config = CubeConfig(
+            dataset_name="LOTL2",
+            band_names=["B03", "B04", "B05", "BQA"],
+            crs='EPSG:4326',
+            bbox=(12.49, 41.88, 12.53, 41.92),
+            spatial_res=0.0001,
+            time_range=('2017-07-01', '2017-07-16'),
+            time_tolerance="30m",
+        )
+
+        cube = open_cube(cube_config, sentinel_hub=sentinel_hub)
+
+        expected_data_types = [np.float32, np.float32, np.float32, np.uint16]
+        for i, b in enumerate(["B03", "B04", "B05", "BQA"]):
+            band = cube[b]
+            self.assertEqual(expected_data_types[i], band.dtype, b)
diff --git a/test/test_sentinelhub.py b/test/test_sentinelhub.py
index 3bd5e8b..a3b04f3 100644
--- a/test/test_sentinelhub.py
+++ b/test/test_sentinelhub.py
@@ -88,9 +88,9 @@ class SentinelHubCatalogSearchTest(unittest.TestCase):
                   1705367, 5857046),
             crs='EPSG:3857'
         )
-        print(json.dumps(features, indent=2))
+        # print(json.dumps(features, indent=2))
         self.assertEqual(1, len(features))
-        self.assertEqual({}, features[0])
+        self.assertEqual({'properties': {}}, features[0])
 
         features = SentinelHub().get_features(
             collection_name=user_collection_id,
@@ -253,6 +253,11 @@ class SentinelHubGetDataTest(unittest.TestCase):
                                                 dtype=np.float32),
                                        np_array[0, 511, -10:, 0])
 
+    # As of 2022-12-06 this tests produces HTTP code 500 (server error)
+
+    @unittest.skip(
+        'As of 2022-12-06 this tests produces HTTP code 500 (server error)'
+    )
     def test_get_data_single_binary_byod(self):
         with open(REQUEST_SINGLE_BYOD_JSON, 'r') as fp:
             request = json.load(fp)