sertit / eosets

This library aims to simplify any process working with different sets of EO data handled by EOReader
https://eosets.readthedocs.io/en/latest/
Apache License 2.0
9 stars 0 forks source link

ERROR with Series().stack() #6

Closed bcerripromethee closed 2 months ago

bcerripromethee commented 4 months ago

Description

Hello, I am actually trying to create a series of two granule of Sentinel-2 L2A ( _S2A_MSIL2A_20220713T103041_N0400_R108_T31TFJ20220713T164557.SAFE and _S2A_MSIL2A_20220716T103641_N0400_R008_T31TFJ20220716T183414.SAFE).

After creating my Series, this latter fails to stack them, do you have an idea of the issue ?

Code To Reproduce

from eoreader.bands import RED, GREEN, SWIR_1
from eosets import Series
import geopandas as gpd

bands = [RED, GREEN, SWIR_1]
aoi_gdf = gpd.read_file("path/to/aoi.geojson")
products_paths = ["path/to/S2A_MSIL2A_20220713T103041_N0400_R108_T31TFJ_20220713T164557.SAFE.zip", "path/to/S2A_MSIL2A_20220716T103641_N0400_R008_T31TFJ_20220716T183414.SAFE.zip"]
multi_reader = Series(paths=products_paths)
products = multi_reader.stack(bands, window=aoi_gdf).astype(np.float32)

Output

products = multi_reader.stack(bands, window=aoi_gdf).astype(np.float32)
               ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
  File "C:\ProgramData\miniconda3\envs\ff_env_toml\Lib\site-packages\eosets\series.py", line 426, in stack
    band_ds = self.load(bands, pixel_size=pixel_size, **kwargs)
              ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
  File "C:\ProgramData\miniconda3\envs\ff_env_toml\Lib\site-packages\eosets\series.py", line 335, in load
    reference_ds = self.reference_mosaic.load(
                   ^^^^^^^^^^^^^^^^^^^^^^^^^^^
  File "C:\ProgramData\miniconda3\envs\ff_env_toml\Lib\site-packages\eosets\mosaic.py", line 408, in load
    merge_fct(prod_path, output_path)
  File "C:\ProgramData\miniconda3\envs\ff_env_toml\Lib\site-packages\sertit\rasters.py", line 1212, in merge_vrt
    return rasters_rio.merge_vrt(crs_paths, crs_merged_path, abs_path, **kwargs)
           ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
  File "C:\ProgramData\miniconda3\envs\ff_env_toml\Lib\site-packages\sertit\rasters_rio.py", line 1384, in merge_vrt
    with rasterio.open(str(crs_path)) as src:
         ^^^^^^^^^^^^^^^^^^^^^^^^^^^^
  File "C:\ProgramData\miniconda3\envs\ff_env_toml\Lib\site-packages\rasterio\env.py", line 451, in wrapper
    return f(*args, **kwds)
           ^^^^^^^^^^^^^^^^
  File "C:\ProgramData\miniconda3\envs\ff_env_toml\Lib\site-packages\rasterio\__init__.py", line 317, in open
    dataset = DatasetReader(path, driver=driver, sharing=sharing, **kwargs)
              ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
  File "rasterio\\_base.pyx", line 312, in rasterio._base.DatasetBase.__init__
rasterio.errors.RasterioIOError: C:/Users/BASTIE~1/AppData/Local/Temp/tmpmsfu4pnc/tmp/T31TFJ_20220713T103041_B11_20m.jp2: No such file or directory

Environment:

remi-braun commented 4 months ago

Hi @bcerripromethee!

Unfortunately I'm out of office till the first of July, so I won't investiguate this issue until then. In the meantime, could you try with unzipped products?

I think there is a copy issue (maybe when creating the VRT?)

remi-braun commented 2 months ago

Hi @bcerripromethee !

Sorry for the delay, I've only been able to investigate this now 😢 This is a bug due to the fact the S-2 SWIR band wasn't looked for ion the right place (because of different native resolution).

I'll fix this in next commit :)

remi-braun commented 2 months ago

This should work now :) (install latest version of sertit[full] though)

from eosets import Series
import geopandas as gpd
# Create logger
import logging
from sertit import AnyPath, logs
logs.init_logger(logging.getLogger("eoreader"))
logs.init_logger(logging.getLogger("eopairs"))
# Paths
data_path = AnyPath(r"D:\_EXTRACTEO\DS3\CI\eosets\DATA")
s2_paths = [data_path / "S2A_MSIL1C_20200824T110631_N0209_R137_T29TQE_20200824T150432.zip", data_path / "S2B_MSIL1C_20200908T110619_N0209_R137_T29TQE_20200908T132324.zip"]
aoi = data_path / "Fire_Spain.geojson"
# Create pair
series = Series(paths=s2_paths, output_path=r"D:\_EXTRACTEO\OUTPUT\eosets\Series2")
# Stack
stack = series.stack(bands=[RED, GREEN, SWIR_2], window=aoi, stack_path=series.output / "stack.tif")
>>> stack
stack
<xarray.DataArray '2020-08-24T11:06:31_RED_2020-09-08T11:06:19_RED_2020-08-24T11:06:31_GREEN_2020-09-08T11:06:19_GREEN_2020-08-24T11:06:31_SWIR_2_2020-09-08T11:06:19_SWIR_2' (
                                                                                                                                                                               bands: 6,
                                                                                                                                                                               y: 1625,
                                                                                                                                                                               x: 1634)> Size: 64MB
array([[[0.0876, 0.0774, 0.1084, ..., 0.1805, 0.1812, 0.1905],
        [0.1046, 0.0822, 0.0826, ..., 0.1835, 0.1854, 0.1841],
        [0.1221, 0.0942, 0.0755, ..., 0.1756, 0.1826, 0.1858],
        ...,
        [0.1456, 0.1516, 0.1559, ..., 0.0901, 0.1007, 0.0983],
        [0.1527, 0.1519, 0.1576, ..., 0.1017, 0.106 , 0.1178],
        [0.1554, 0.1519, 0.1568, ..., 0.0836, 0.1028, 0.116 ]],
       [[0.078 , 0.0964, 0.1215, ..., 0.1809, 0.1895, 0.2032],
        [0.0847, 0.0763, 0.111 , ..., 0.1814, 0.1809, 0.1971],
        [0.1038, 0.0773, 0.0824, ..., 0.1856, 0.1842, 0.1852],
        ...,
        [0.1394, 0.1537, 0.1573, ..., 0.0713, 0.0901, 0.0959],
        [0.1563, 0.161 , 0.1673, ..., 0.0816, 0.1045, 0.0871],
        [0.1652, 0.1629, 0.1667, ..., 0.1017, 0.1043, 0.1193]],
       [[0.0939, 0.094 , 0.1061, ..., 0.1538, 0.1541, 0.1581],
        [0.1008, 0.0917, 0.0914, ..., 0.153 , 0.1585, 0.1541],
        [0.1093, 0.0974, 0.0917, ..., 0.1476, 0.1536, 0.1531],
        ...,
...
        [0.1152, 0.1252, 0.1338, ..., 0.08  , 0.0868, 0.0922],
        [0.1269, 0.1286, 0.1443, ..., 0.0904, 0.0955, 0.0907],
        [0.1295, 0.1344, 0.1463, ..., 0.0969, 0.101 , 0.1017]],
       [[0.1472, 0.1468, 0.1459, ..., 0.2412, 0.2431, 0.244 ],
        [0.1537, 0.1515, 0.1473, ..., 0.2408, 0.2427, 0.2437],
        [0.1666, 0.161 , 0.15  , ..., 0.2398, 0.2421, 0.2432],
        ...,
        [0.2329, 0.241 , 0.2572, ..., 0.1185, 0.1304, 0.1364],
        [0.248 , 0.2544, 0.2673, ..., 0.1166, 0.1283, 0.1341],
        [0.2602, 0.2635, 0.27  , ..., 0.1153, 0.1264, 0.1319]],
       [[0.1518, 0.1569, 0.1671, ..., 0.2638, 0.2651, 0.2657],
        [0.1516, 0.1553, 0.1627, ..., 0.2635, 0.2649, 0.2655],
        [0.1511, 0.152 , 0.1538, ..., 0.2629, 0.2644, 0.2652],
        ...,
        [0.2272, 0.2368, 0.2559, ..., 0.122 , 0.1326, 0.1379],
        [0.2498, 0.2582, 0.2751, ..., 0.1151, 0.1279, 0.1343],
        [0.2691, 0.2761, 0.2901, ..., 0.1144, 0.1285, 0.1356]]],
      dtype=float32)
Coordinates:
  * x            (x) float64 13kB 7.705e+05 7.705e+05 ... 7.868e+05 7.868e+05
  * y            (y) float64 13kB 4.458e+06 4.458e+06 ... 4.442e+06 4.442e+06
    spatial_ref  int32 4B 0
  * bands        (bands) object 48B MultiIndex
  * variable     (bands) object 48B SpectralBandNames.RED ... SpectralBandNam...
  * band         (bands) int32 24B 1 1 1 1 1 1
  * time         (bands) datetime64[ns] 48B 2020-08-24T11:06:31 ... 2020-09-0...
Attributes:
    long_name:       2020-08-24T11:06:31_RED 2020-09-08T11:06:19_RED 2020-08-...
    condensed_name:  20200824T110631_S2_T29TQE_L1C_150432_20200908T110619_S2_...