Explore EOPF Zarr Sentinel-2 files using GDAL
Explore how to open, visualise and plot Sentinel-2 L1C MSI EOPF Zarr format
Explore EOPF Zarr Sentinel-2 files using GDAL¶
Table of Contents¶
Introduction¶
This notebook demonstrates accessing and visualizing Sentinel-2 satellite imagery stored in ZARR format using GDAL. Learn to connect to remote ZARR files, extract metadata, and visualize multispectral bands.
Import Libraries¶
from osgeo import gdal
import matplotlib.pyplot as pltBuilding GDAL Connection¶
Note: Always prepend ZARR:"/vsicurl/ to ZARR file paths
dsn = 'ZARR:"/vsicurl/https://objectstore.eodc.eu:2222/e05ab01a9d56408d82ac32d69a5aae2a:sample-data/tutorial_data/cpm_v253/S2B_MSIL1C_20250113T103309_N0511_R108_T32TLQ_20250113T122458.zarr"'Exploring the Zarr Dataset¶
ds = gdal.Open(dsn, gdal.GA_ReadOnly)
if ds is None:
raise RuntimeError("Failed to open the Zarr dataset with GDAL.")print("Dataset opened successfully!")
print("Metadata:", ds.GetMetadata())Dataset opened successfully!
Metadata: {}
Metadata Extraction¶
sub_ds = ds.GetSubDatasets()
for idx, s in enumerate(sub_ds[:5]): # Show first 5 subdatasets
print(f"{idx}: {s[0]}")0: ZARR:"/vsicurl/https://objectstore.eodc.eu:2222/e05ab01a9d56408d82ac32d69a5aae2a:sample-data/tutorial_data/cpm_v253/S2B_MSIL1C_20250113T103309_N0511_R108_T32TLQ_20250113T122458.zarr":/conditions/geometry/angle
1: ZARR:"/vsicurl/https://objectstore.eodc.eu:2222/e05ab01a9d56408d82ac32d69a5aae2a:sample-data/tutorial_data/cpm_v253/S2B_MSIL1C_20250113T103309_N0511_R108_T32TLQ_20250113T122458.zarr":/conditions/geometry/band
2: ZARR:"/vsicurl/https://objectstore.eodc.eu:2222/e05ab01a9d56408d82ac32d69a5aae2a:sample-data/tutorial_data/cpm_v253/S2B_MSIL1C_20250113T103309_N0511_R108_T32TLQ_20250113T122458.zarr":/conditions/geometry/detector
3: ZARR:"/vsicurl/https://objectstore.eodc.eu:2222/e05ab01a9d56408d82ac32d69a5aae2a:sample-data/tutorial_data/cpm_v253/S2B_MSIL1C_20250113T103309_N0511_R108_T32TLQ_20250113T122458.zarr":/conditions/geometry/x
4: ZARR:"/vsicurl/https://objectstore.eodc.eu:2222/e05ab01a9d56408d82ac32d69a5aae2a:sample-data/tutorial_data/cpm_v253/S2B_MSIL1C_20250113T103309_N0511_R108_T32TLQ_20250113T122458.zarr":/conditions/geometry/y
target_subds = gdal.Open(sub_ds[80][0], gdal.GA_ReadOnly)
print("Band metadata:", target_subds.GetMetadata())Band metadata: {'_eopf_attrs': '{ "add_offset": -0.1, "coordinates": [ "x", "y" ], "dimensions": [ "y", "x" ], "fill_value": 0, "scale_factor": 0.0001, "units": "digital_counts" }', 'dtype': '<u2', 'fill_value': '0', 'long_name': 'TOA reflectance from MSI acquisition at spectral band b01 442.3 nm', 'proj:bbox': '{300000,4890240,409800,5000040}', 'proj:epsg': '32632', 'proj:shape': '{1830,1830}', 'proj:transform': '{60,0,300000,0,-60,5000040,0,0,1}', 'proj:wkt2': 'PROJCS["WGS 84 / UTM zone 32N",GEOGCS["WGS 84",DATUM["WGS_1984",SPHEROID["WGS 84",6378137,298.257223563,AUTHORITY["EPSG","7030"]],AUTHORITY["EPSG","6326"]],PRIMEM["Greenwich",0,AUTHORITY["EPSG","8901"]],UNIT["degree",0.0174532925199433,AUTHORITY["EPSG","9122"]],AUTHORITY["EPSG","4326"]],PROJECTION["Transverse_Mercator"],PARAMETER["latitude_of_origin",0],PARAMETER["central_meridian",9],PARAMETER["scale_factor",0.9996],PARAMETER["false_easting",500000],PARAMETER["false_northing",0],UNIT["metre",1,AUTHORITY["EPSG","9001"]],AXIS["Easting",EAST],AXIS["Northing",NORTH],AUTHORITY["EPSG","32632"]]', 'valid_max': '65535', 'valid_min': '1'}
Data Visualization¶
band = target_subds.GetRasterBand(1)
data = band.ReadAsArray()
print(f"Data dimensions: {data.shape}")Data dimensions: (1830, 1830)
plt.figure(figsize=(10, 8))
plt.imshow(data, cmap="viridis")
plt.colorbar(label="Reflectance Value")
plt.title("Sentinel-2 Band B01 (60m Resolution)")
plt.axis("off")
plt.show()
