STAC

STAC#

Let’s use EOReader to create SpatioTemporal Asset Catalog (STAC) items.

Note: This is experimental for now, use it at your own risk !

Warning: You will need to install pystac[validation], folium and eodag (version != 2.6.0) to run this notebook

Imports#

# Imports
import os

import pystac
import geopandas as gpd
from tempfile import TemporaryDirectory
from shapely.geometry import mapping

from eodag import setup_logging
from eodag.api.core import EODataAccessGateway

from eoreader.reader import Reader

/home/docs/checkouts/readthedocs.org/user_builds/eoreader/envs/stable/lib/python3.9/site-packages/dask/dataframe/__init__.py:42: FutureWarning: 
Dask dataframe query planning is disabled because dask-expr is not installed.

You can install it with `pip install dask[dataframe]` or `conda install dask`.
This will raise in a future version.

  warnings.warn(msg, FutureWarning)

Create logger#

# Create logger
import logging
from sertit import logs

logger = logging.getLogger("eoreader")
logs.init_logger(logger)

Linking some data#

Let’s take 3 products covering approximately the same area (over DAX city in France):

  • One Landsat-8 OLI-TIRS collection 2

  • One Landsat-5 TM collection 2

  • One Sentinel-2 L1C

prod_folder = os.path.join("/home", "prods")
paths = [
    # Landsat-8 OLI-TIRS collection 2
    os.path.join(prod_folder, "LANDSATS_COL2", "LC08_L1TP_200030_20201220_20210310_02_T1.tar"),
    # Landsat-5 TM collection 2    
    os.path.join(prod_folder, "LANDSATS_COL2", "LT05_L1TP_200030_20111110_20200820_02_T1.tar"),
    # Sentinel-2 L2A
    os.path.join(prod_folder, "S2", "PB 02.07+", "S2A_MSIL1C_20191215T110441_N0208_R094_T30TXP_20191215T114155.SAFE"),
]

Create STAC catalog#

Create a STAC catalog and add 3 STAC items to it.

# Create the reader
reader = Reader()

# Work in a temporary directory
tmp = TemporaryDirectory()

# Create STAC catalog
catalog_path = os.path.join(tmp.name, "catalog.json")
catalog = pystac.Catalog(
    id='SERTIT_101',
    description="SERTIT's Catalog",
    title='SERTIT Catalog',
    href=catalog_path
)

# Add all the products into the STAC catalog
for path in paths:
    logger.info(f"*** {os.path.basename(path)} ***")

    # Open the product
    prod = reader.open(path, remove_tmp=True)

    # Get item
    item = prod.stac.create_item()

    # Add item to catalogue
    catalog.add_item(item)

2025-05-05 11:29:39,067 - [INFO] - *** LC08_L1TP_200030_20201220_20210310_02_T1.tar ***
2025-05-05 11:29:39,071 - [WARNING] - There is no existing products in EOReader corresponding to /home/prods/LANDSATS_COL2/LC08_L1TP_200030_20201220_20210310_02_T1.tar.
2025-05-05 11:29:39,072 - [INFO] - Your given path may not be a satellite image. If it is, maybe the product isn't handled by EOReader. If you are sure this product is handled, it is either corrupted or you may need to go deeper in the filetree to find the correct path to give.
2025-05-05 11:29:39,073 - [DEBUG] - Please look at what folder you should give to EOReader by accessing the documentation: https://eoreader.readthedocs.io/latest/main_features.html#recognized-paths

---------------------------------------------------------------------------
AttributeError                            Traceback (most recent call last)
Cell In[6], line 9
      6 prod = reader.open(path, remove_tmp=True)
      8 # Get item
----> 9 item = prod.stac.create_item()
     11 # Add item to catalogue
     12 catalog.add_item(item)

AttributeError: 'NoneType' object has no attribute 'stac'

# Save catalog
catalog.describe()
catalog.normalize_and_save(tmp.name, catalog_type=pystac.CatalogType.SELF_CONTAINED)

* <Catalog id=SERTIT_101>
  * <Item id=20201220T104856_L8_200030_OLI_TIRS>
  * <Item id=20111110T103612_L5_200030_TM>
  * <Item id=20191215T110441_S2_T30TXP_L1C_114155>

list(catalog.get_items())[0]

Query the catalog#

EODAG is an opensource python library that implements STAC and allows you to query your local STAC catalog.
Look at here for a detailed tutorial.

# Create an EODAG custom STAC provider
dag = EODataAccessGateway()

# Set EODAG logging level to WARNING
setup_logging(verbose=1)

# Add the custom STAC provider, exactly like in the tutorial mentioned above
dag.update_providers_config("""
stac_http_provider:
    search:
        type: StaticStacSearch
        api_endpoint: %s
    products:
        GENERIC_PRODUCT_TYPE:
            productType: '{productType}'
    download:
        type: HTTPDownload
        base_uri: %s
        flatten_top_dirs: True
        outputs_prefix: %s
""" % (catalog_path, tmp.name, tmp.name))

# Set the custom STAC provider as preferred
dag.set_preferred_provider("stac_http_provider")

# Query every product from inside the catalog
all_products, _ = dag.search()

---------------------------------------------------------------------------
KeyError                                  Traceback (most recent call last)
File /opt/conda/lib/python3.10/site-packages/eodag/plugins/manager.py:156, in PluginManager.get_search_plugins(self, product_type, provider)
    155 try:
--> 156     config = self.providers_config[provider]
    157 except KeyError:

KeyError: 'stac_http_provider'

During handling of the above exception, another exception occurred:

UnsupportedProvider                       Traceback (most recent call last)
Cell In[10], line 2
      1 # Query every product from inside the catalog
----> 2 all_products, _ = dag.search()

File /opt/conda/lib/python3.10/site-packages/eodag/api/core.py:887, in EODataAccessGateway.search(self, page, items_per_page, raise_errors, start, end, geom, locations, **kwargs)
    825 def search(
    826     self,
    827     page=DEFAULT_PAGE,
   (...)
    834     **kwargs,
    835 ):
    836     """Look for products matching criteria on known providers.
    837 
    838     The default behaviour is to look for products on the provider with the
   (...)
    885         enforced here.
    886     """
--> 887     search_kwargs = self._prepare_search(
    888         start=start, end=end, geom=geom, locations=locations, **kwargs
    889     )
    890     search_plugin = search_kwargs.pop("search_plugin", None)
    891     if search_kwargs.get("id"):
    892         # adds minimal pagination to be able to check only 1 product is returned

File /opt/conda/lib/python3.10/site-packages/eodag/api/core.py:1306, in EODataAccessGateway._prepare_search(self, start, end, geom, locations, **kwargs)
   1299 if (
   1300     product_type
   1301     not in self._plugins_manager.product_type_to_provider_config_map.keys()
   1302 ):
   1303     logger.debug(
   1304         f"Fetching external product types sources to find {product_type} product type"
   1305     )
-> 1306     self.fetch_product_types_list()
   1308 search_plugin = next(
   1309     self._plugins_manager.get_search_plugins(product_type=product_type)
   1310 )
   1311 if search_plugin.provider != self.get_preferred_provider()[0]:

File /opt/conda/lib/python3.10/site-packages/eodag/api/core.py:613, in EODataAccessGateway.fetch_product_types_list(self, provider)
    608         continue
    609     # providers not skipped here should be user-modified
    610     # or not in ext_product_types_conf (if eodag system conf != eodag conf used for ext_product_types_conf)
    611 
    612 # discover product types for user configured provider
--> 613 provider_ext_product_types_conf = self.discover_product_types(
    614     provider=provider
    615 )
    617 # update eodag product types list with new conf
    618 self.update_product_types_list(provider_ext_product_types_conf)

File /opt/conda/lib/python3.10/site-packages/eodag/api/core.py:648, in EODataAccessGateway.discover_product_types(self, provider)
    646     return
    647 if getattr(search_plugin_config, "discover_product_types", None):
--> 648     search_plugin = next(
    649         self._plugins_manager.get_search_plugins(provider=provider)
    650     )
    651     # append auth to search plugin if needed
    652     if getattr(search_plugin.config, "need_auth", False):

File /opt/conda/lib/python3.10/site-packages/eodag/plugins/manager.py:158, in PluginManager.get_search_plugins(self, product_type, provider)
    156     config = self.providers_config[provider]
    157 except KeyError:
--> 158     raise UnsupportedProvider
    159 yield get_plugin()
    160 # Signal the end of iteration as we already have what we wanted (see PEP-479)

UnsupportedProvider: 

# Load an AOI
aoi_path = os.path.join("/home", "aois", "DAX.geojson")
aoi = gpd.read_file(aoi_path)
aoi_geojson = mapping(aoi.geometry.values[0])

# Query spatially with the AOI and temporally with a time period
query_args = {"start": "2020-05-01", "end": "2022-05-06", "geom": aoi.geometry.values[0]}
query_products, _ = dag.search(**query_args)

query_products[0]

query_products[0].assets

Display the results#

We can use folium to display the results geometry over a map.

import folium

# Create a map zoomed over the search area
fmap = folium.Map((43.2, -1.05), zoom_start=7)

# Add a layer green layer for the query over the AOI
folium.GeoJson(
    data=all_products.as_geojson_object(),
    tooltip = "All products stored in the catalog",
    style_function=lambda x: {'color': 'green'}
).add_to(fmap)

# Add a layer green layer for the query over the AOI
folium.GeoJson(
    data=query_products.as_geojson_object(),
    tooltip = "Retrieved products with the query",
    style_function=lambda x: {'color': 'red'}
).add_to(fmap)

# Add a layer blue layer for the AOI
folium.GeoJson(
    data=aoi_geojson,
    tooltip = "DAX AOI",
    style_function=lambda x: {'color': 'blue'}
).add_to(fmap)

fmap

# Clean the tmp directory
tmp.cleanup()
Contents