gjoseph92 / stackstac

Turn a STAC catalog into a dask-based xarray
https://stackstac.readthedocs.io
MIT License
232 stars 49 forks source link

'bounds' is messing up the x, y coordinnates #242

Open Berhinj opened 5 months ago

Berhinj commented 5 months ago

This used to not be the case but it's suddenly popping up, I don't know what happened.

Suddenly the 'bounds' parameter from stackstac.stack function messes up the original coordinates of the images. When using the 'bounds' parameter, instead of preserving the original x/y coordinnates, stackstac modifies them and they'll give you float values:

import numpy as np
import rasterio
import stackstac
import rioxarray as rxr  # noqa: F401
import json
import pystac_client

bbox = [2.3575522081018616, 46.9821947096507, 2.373435778967132, 46.98773884336863]
bbox_utm = [451141.38241652114, 5203375.974548349, 452354.2625236209, 5204001.850759215]
datetime="2022-12-01/2023-02-10"
cloud_cover=30
processing_level="l2"

URL = "https://landsatlook.usgs.gov/stac-server"

catalog = pystac_client.Client.open(URL)

params = dict(
    max_items=None,
    limit=1000,
    bbox=bbox,
    datetime=datetime,
    query={
        "eo:cloud_cover": {"lt": cloud_cover},
        "landsat:collection_category": {"eq": "T1"},
        "landsat:collection_number": {"eq": "02"},
        "platform": {"in": ["LANDSAT_8", "LANDSAT_9"]},
    },
)

xsr = catalog.search(collections=["landsat-c2l2-sr"] , **params).item_collection()

da = stackstac.stack(
    xsr[0],
    assets={"image/tiff", "image/vnd.stac.geotiff", "image/x.geotiff"},
    bounds=list(aoi.to_crs(crs_utm).bbox),
    snap_bounds=False,
    gdal_env=stackstac.DEFAULT_GDAL_ENV.updated(
        always=dict(session=rasterio.session.AWSSession(requester_pays=True))
    ),
)
da.x

Will give you

Screenshot 2024-03-01 at 14 43 04

While removing the 'bounds' will give you the

da = stackstac.stack(
    xsr[0],
    assets={"image/tiff", "image/vnd.stac.geotiff", "image/x.geotiff"},
    # bounds=list(aoi.to_crs(crs_utm).bbox),
    snap_bounds=False,
    gdal_env=stackstac.DEFAULT_GDAL_ENV.updated(
        always=dict(session=rasterio.session.AWSSession(requester_pays=True))
    ),
)
da.x

will give you the original integer like x, y coordinnates

Screenshot 2024-03-01 at 14 44 07