Search for products by tile
Contents
Search for products by tile#
Built-in search by tile#
Many providers already support search by tile (Sentinel 2 MGRS tiling grid): peps
, theia
, mundi
, onda
, creodias
, cop_dataspace
, planetary_computer
, earth_search
.
For these providers, you can use tileIdentifier
as search parameter using EODAG:
[1]:
from pprint import pprint
from eodag import EODataAccessGateway, setup_logging
dag = EODataAccessGateway()
dag.set_preferred_provider("peps")
products, _ = dag.search(
productType="S2_MSI_L1C",
start="2018-06-01",
end="2018-06-15",
tileIdentifier="31TFK"
)
pprint([(product, product.properties["tileIdentifier"]) for product in products])
[(EOProduct(id=S2A_MSIL1C_20180614T103021_N0206_R108_T31TFK_20180614T124154, provider=peps),
'31TFK'),
(EOProduct(id=S2B_MSIL1C_20180612T104019_N0206_R008_T31TFK_20180612T124717, provider=peps),
'31TFK'),
(EOProduct(id=S2B_MSIL1C_20180609T103019_N0206_R108_T31TFK_20180609T131037, provider=peps),
'31TFK'),
(EOProduct(id=S2A_MSIL1C_20180607T104021_N0206_R008_T31TFK_20180607T124742, provider=peps),
'31TFK'),
(EOProduct(id=S2A_MSIL1C_20180604T103021_N0206_R108_T31TFK_20180604T141551, provider=peps),
'31TFK'),
(EOProduct(id=S2B_MSIL1C_20180602T104019_N0206_R008_T31TFK_20180602T132118, provider=peps),
'31TFK')]
Search by tile using custom locations#
eodag
allows to search for products by geometric features that match a location query, see the API user guide for an introduction to this concept.
In this tutorial we will use a shapefile that represents the Sentinel 2 tiling grid to search for Sentinel 2 Level-1C products with PEPS at a specific tile. In this shapefile each tile is defined by its centroid and a tile_id
attribute (e.g. 29PMT). This shapefile was created by downloading first the Sentinel 2 tiling grid (MGRS) provided by ESA as a KML file. It
was then converted as a shapefile and processed to compute the centroids. We use the tile’s centroid here as eodag
returns products that intersects the user defined search area. Since tiles overlap with each other, using the polygons instead of the centroids would return more tiles than just the one we target.
[2]:
import os
from zipfile import ZipFile
# Interactive mapping
import folium
from folium.plugins import TimestampedGeoJson
# pyshp: to read shapefiles
import shapefile
Setup#
A workspace directory is created to store the files that will be generated.
[3]:
workspace = "eodag_workspace_locations_tiles"
if not os.path.isdir(workspace):
os.mkdir(workspace)
You should have an auxdata
folder next to this tutorial’s file. It contains a shapefile that is needed to run this tutorial correctly.
[4]:
sentinel2_grid_zip = os.path.join("auxdata", "sentinel2_tiling_grid_centroids.zip")
if not os.path.isfile(sentinel2_grid_zip):
raise FileNotFoundError("Auxdata not found, please check your configuration.")
[5]:
# We unzip the archived shapefile.
with ZipFile(sentinel2_grid_zip, "r") as fzip:
fzip.extractall("auxdata")
In this tutorial products will just be searched for, not downloaded. We don’t need to set up PEPS credentials to search for products. If you wish to download them, you should set the credentials beforehand, using these two environment variables for instance.
[6]:
# os.environ["EODAG__PEPS__AUTH__CREDENTIALS__USERNAME"] = "PLEASE_CHANGE_ME"
# os.environ["EODAG__PEPS__AUTH__CREDENTIALS__PASSWORD"] = "PLEASE_CHANGE_ME"
Logging is activated to better inspect what eodag
does internally.
[7]:
setup_logging(2) # INFO level
The default search criteria consists of a time period in June 2018 and eodag
’s product type identifier for Sentinel 2 Level-1C products.
[8]:
default_search_criteria = dict(
productType="S2_MSI_L1C",
start="2018-06-01",
end="2018-06-15"
)
Add a locations configuration#
We check and store the content of this shapefile.
[9]:
sentinel2_shp = os.path.join('auxdata', 'sentinel2_tiling_grid_centroids.shp')
with shapefile.Reader(sentinel2_shp) as shp:
print(shp, "\n")
print("fields:", shp.fields)
shaperecs = shp.shapeRecords()
shapefile Reader
56686 shapes (type 'POINT')
56686 records (2 fields)
fields: [('DeletionFlag', 'C', 1, 0), ['tile_id', 'C', 5, 0]]
It has about 57 000 tiles/polygons and a field tile_id
.
We create a YAML file to configure this new location selector, we will refer to it with s2_tile_centroid
.
[10]:
# Save the locations configuration file.
locations_yaml_content = """
shapefiles:
- name: s2_tile_centroid
path: {}
attr: tile_id
""".format(os.path.abspath(sentinel2_shp))
locations_filepath = os.path.abspath(os.path.join(workspace, "custom_locations.yml"))
with open(locations_filepath, "w") as f_yml:
f_yml.write(locations_yaml_content.strip())
An instance of an EODataAccessGateway class is created, it makes use of this location configuration file.
[11]:
dag = EODataAccessGateway(locations_conf_path=locations_filepath)
2023-04-17 10:15:54,835 eodag.config [INFO ] Loading user configuration from: /home/sylvain/.config/eodag/eodag.yml
2023-04-17 10:15:55,054 eodag.core [INFO ] Locations configuration loaded from /home/sylvain/workspace/eodag/docs/notebooks/tutos/eodag_workspace_locations_tiles/custom_locations.yml
We want to look for Sentinel 2 Level-1C products. We can check whether this product type is offered by PEPS (as configured in eodag
). If so, PEPS is set as the provider used to search for products.
[12]:
"peps" in dag.available_providers("S2_MSI_L1C")
[12]:
True
[13]:
dag.set_preferred_provider("peps")
Search#
A single tile#
Our target tile is 31TFK
and is located in the South-East of France. Its feature is retrieved from the shapefil to be displayed later on an interactive map.
[14]:
targeted_tile_name = "31TFK"
# Get the targeted tile feature
targeted_tile = [
sr
for sr in shaperecs
if sr.record["tile_id"] == "31TFK"
][0]
We search for all the products that intersect with the centroid of this tile.
[15]:
products = dag.search_all(
locations=dict(s2_tile_centroid="31TFK"),
**default_search_criteria
)
print(f"{len(products)} were found given the above search criteria")
2023-04-17 10:16:38,765 eodag.core [INFO ] Searching product type 'S2_MSI_L1C' on provider: peps
2023-04-17 10:16:38,766 eodag.core [INFO ] Iterate search over multiple pages: page #1
2023-04-17 10:16:38,778 eodag.plugins.search.qssearch [INFO ] Sending search request: https://peps.cnes.fr/resto/api/collections/S2ST/search.json?startDate=2018-06-01&completionDate=2018-06-15&geometry=POINT (4.9533 44.6422)&productType=S2MSI1C&maxRecords=500&page=1
2023-04-17 10:16:40,033 eodag.core [INFO ] Found 6 result(s) on provider 'peps'
6 were found given the above search criteria
The products found are displayed on an interactive map along with the centroid of the targeted tile. A time player allows to see when the products were sensed.
[16]:
# The GeoJSON representation has to be slightly adapted for the time slider
adapted_prods = products.as_geojson_object()
for feature in adapted_prods["features"]:
feature["properties"]["time"] = feature["properties"]["startTimeFromAscendingNode"]
# Create a map zoomed over the search area
fmap = folium.Map([44.5, 5], zoom_start=8)
# Add a layer that map the tile's centroid
folium.GeoJson(
data=targeted_tile,
tooltip = targeted_tile_name,
).add_to(fmap)
# Add layer that temporally maps the products found
TimestampedGeoJson(
adapted_prods,
transition_time=50, # Transition duration in ms
period="PT3H", # Array of times, here every 3 hours
duration="PT12H", # Feature display duragion, here 6 hours
time_slider_drag_update=True, # Update the map when the slider is dragged
auto_play=False, # Don't auto play the animation
).add_to(fmap)
fmap
[16]:
Note
Instead of using the tile’s centroid we could have directly used its extent and filter the returned products to keep only those fully contained within the tile. Check out the section dedicated to filtering products in the API user guide.
Multiple tiles#
We can search for products that overlap with several tiles using a regular expression. We use the expression "31T[CDE][MLK]"
to look for products over 9 different tiles (31TCM, 31TCL, 31TCK, 31TDM, etc.) over France.
[17]:
products = dag.search_all(
locations=dict(s2_tile_centroid="31T[CDE][MLK]"),
**default_search_criteria
)
2023-04-17 10:17:26,059 eodag.core [INFO ] Searching product type 'S2_MSI_L1C' on provider: peps
2023-04-17 10:17:26,060 eodag.core [INFO ] Iterate search over multiple pages: page #1
2023-04-17 10:17:26,061 eodag.plugins.search.qssearch [INFO ] Sending search request: https://peps.cnes.fr/resto/api/collections/S2ST/search.json?startDate=2018-06-01&completionDate=2018-06-15&geometry=MULTIPOINT (3.7147 46.4567, 3.7032 45.5565, 3.6922 44.6568, 2.4122 46.4574, 2.4216 45.5572, 1.1109 46.4433, 1.1413 45.5436, 1.1703 44.6442, 2.4306 44.6575)&productType=S2MSI1C&maxRecords=500&page=1
2023-04-17 10:17:26,745 eodag.core [INFO ] Found 32 result(s) on provider 'peps'
[18]:
print(f"{len(products)} were found given the above search criteria")
32 were found given the above search criteria
The products are displayed on an interactive map. By hovering over them you can observe that the MGRS number of the tiles match with the regular expressions we used.
[20]:
# Create a map zoomed over the search area
fmap = folium.Map([44.5, 1.5], zoom_start=6)
# Create a layer that maps the products found
folium.GeoJson(
data=products,
tooltip=folium.GeoJsonTooltip(
fields=[
"title", # The product's title
"tileIdentifier", # The tile number on the MGRS grid
]
),
).add_to(fmap)
fmap
[20]:
This example has demonstrated the possibilities offered by eodag
to easily select products from a tile grid by using regular expressions over their identifier.