La carga parcial de datos del cubo de datos abiertos no funciona

Question

He desplegado el cubo en una caja docker en nuestro entorno local y he conseguido indexar un par de datos de nivel 2 de Sentinel-2. Sin embargo, cuando he intentado cargar mi AOI desde esos conjuntos de datos recibo una lista vacía. El código utilizado se copia a continuación:

# Create a query object
lat, lon = 22.821, 28.518
buffer = 0.05

query = {
    'time': ('2020-12-11', '2020-12-13'),
    'x': (lon - buffer, lon + buffer),
    'y': (lat + buffer, lat - buffer),
    'output_crs': 'epsg:3857',
    'resolution':(-10,10),
}

# Load an image, I am loading a Senitnel2 Level1 data from the datacube
bands = ['B02_10m', 'B03_10m', 'B04_10m']

ds = dc.load(product='s2a_sen2cor_granule',
             measurements=bands,
             dask_chunks={},
             **query)

print(ds)

Si ignoro la longitud y latitud entonces funciona bien

# Create a query object
lat, lon = 22.821, 28.518
buffer = 0.05

query = {
    'time': ('2020-12-11', '2020-12-13'),
    'output_crs': 'epsg:3857',
    'resolution':(-10,10),
}

#load data
ds = dc.load(product='s2_l2a',
             dask_chunks={},
             **query)

print(ds)

Y la salida que recibo es la siguiente:

<xarray.Dataset>
Dimensions:      (time: 2, x: 484, y: 591)
Coordinates:
  * time         (time) datetime64[ns] 2021-07-10T08:53:55 2021-07-10T08:54:10
  * y            (y) float64 2.814e+06 2.814e+06 ... 2.803e+06 2.803e+06
  * x            (x) float64 2.765e+06 2.765e+06 ... 2.775e+06 2.775e+06
    spatial_ref  int32 6933
Data variables:
    B01          (time, y, x) uint16 dask.array<chunksize=(1, 591, 484), meta=np.ndarray>
    B02          (time, y, x) uint16 dask.array<chunksize=(1, 591, 484), meta=np.ndarray>
    B03          (time, y, x) uint16 dask.array<chunksize=(1, 591, 484), meta=np.ndarray>
    B04          (time, y, x) uint16 dask.array<chunksize=(1, 591, 484), meta=np.ndarray>
    B05          (time, y, x) uint16 dask.array<chunksize=(1, 591, 484), meta=np.ndarray>
    B06          (time, y, x) uint16 dask.array<chunksize=(1, 591, 484), meta=np.ndarray>
    B07          (time, y, x) uint16 dask.array<chunksize=(1, 591, 484), meta=np.ndarray>
    B08          (time, y, x) uint16 dask.array<chunksize=(1, 591, 484), meta=np.ndarray>
    B8A          (time, y, x) uint16 dask.array<chunksize=(1, 591, 484), meta=np.ndarray>
    B09          (time, y, x) uint16 dask.array<chunksize=(1, 591, 484), meta=np.ndarray>
    B11          (time, y, x) uint16 dask.array<chunksize=(1, 591, 484), meta=np.ndarray>
    B12          (time, y, x) uint16 dask.array<chunksize=(1, 591, 484), meta=np.ndarray>
    SCL          (time, y, x) uint8 dask.array<chunksize=(1, 591, 484), meta=np.ndarray>
    AOT          (time, y, x) uint16 dask.array<chunksize=(1, 591, 484), meta=np.ndarray>
    WVP          (time, y, x) uint16 dask.array<chunksize=(1, 591, 484), meta=np.ndarray>
Attributes:
    crs:           epsg:6933
    grid_mapping:  spatial_ref

¿Cuál puede ser el motivo?

Actualización para @Robbi Bishop-Taylor: Hola, he probado la consulta que has publicado a continuación, pero no ha habido suerte. Así que, borré la base de datos ODC y empecé desde el principio con un nuevo conjunto de datos, por favor, vea la imagen de abajo en el mapa de Google y la marca roja es la coordenada que he utilizado en la consulta

También puedo confirmar que los datos se indexan correctamente en datacube-explorer

Los nuevos parámetros de consulta son los siguientes:

lat, lon = 22.821, 28.518
buffer = 0.05

query = {
    'time': ('2021-07-24', '2021-07-27'),
    'x': (lon - buffer, lon + buffer),
    'y': (lat + buffer, lat - buffer),
    'output_crs': 'epsg:6933',
    'resolution':(-10,10),
}

Todavía no he conseguido cargar los datos AOI/parciales.

NOTA: Los datos de Sentinel-2 que he utilizado son S2A_MSIL2A_20210725T083601_N0301_R064_T35QPF_20210725T115425.SAFE

Answer 1

Por lo que veo, su consulta está configurada correctamente - es totalmente posible consultar por coordenadas geográficas (lat-lon) y devolver datos en un sistema de referencia de coordenadas diferente (por ejemplo. EPSG:3857 con unidades de metro como se hace arriba, aunque @whatnick tiene razón en que este CRS no se suele utilizar para el análisis). Si no se especifica ninguna extensión espacial, simplemente se obtienen todos los datos, independientemente de su ubicación; en este caso se cargan en el CRS nativo de los datos, que es EPSG:6933 - se trata de una cuadrícula global (WGS 84 / NSIDC EASE-Grid 2.0 Global) por lo que creo que también debería estar bien.

He trazado tus coordenadas geográficas de entrada con su buffer de 0,05, y los límites de los datos que has cargado en EPSG:6933 sin ninguna consulta espacial. Parece que los datos indexados sólo cubren una pequeña zona que no se solapa con los límites de la consulta:

Para descartar esto, recomendaría volver a ejecutar su código utilizando la siguiente entrada de consulta lat-lon en su lugar:

# Create a query object
lat, lon = 22.58318610, 28.70269922
buffer = 0.5

query = {
    'time': ('2020-12-11', '2020-12-13'),
    'x': (lon - buffer, lon + buffer),
    'y': (lat + buffer, lat - buffer),
    'output_crs': 'epsg:3857',
    'resolution':(-10, 10),
}

Answer 2

He conseguido solucionar este problema. La consulta estaba absolutamente bien; sin embargo, la preparación script que he utilizado del repositorio Github del cubo de datos necesitaba un ajuste. Para tu información, usé el sen2cor_prepare.py para indexar los datos de L2A descargados del hub de Copernicus, este script generó el archivo YAML con las coordenadas erróneas; simplemente añadí la siguiente línea en el archivo get_coords para que funcione. Gracias a mi amigo que me ayudó a resolver esto.

spatial_ref.SetAxisMappingStrategy(osr.OAMS_TRADITIONAL_GIS_ORDER)

NOTA: Es debido a la versión de GDAL que estoy utilizando, puede que sólo funcione con la versión >3

Answer 3

Parece que sus datos están almacenados en 'epsg:6933' (vea los atributos en el conjunto de datos que cargó sin carga parcial de datos), que es una proyección de área igual de Albers para África (las coordenadas están en metros), pero usted está tratando de usar coordenadas lat/lon para recortar y cargar los datos. Lo mejor es que averigües qué es lo que está fallando en tus sistemas de coordenadas. Además, las medidas que pasas a dc.load (B04_10m, etc.) no coinciden con las de los datos (B04.), y los nombres de tus productos no coinciden: 's2_l2a' frente a 's2a_sen2cor_granule'.

Answer 4

Hay todo tipo de cosas que suceden con este conjunto particular de variables. EPSG:3857 utiliza un sistema de coordenadas en metros en lugar de grados. Además, sólo se utiliza realmente en mapas deslizantes, para el análisis sería mejor una proyección de área igual.

Como la cónica de Albers sintonizada con su área : https://en.m.wikipedia.org/wiki/Albers_projection