Skip to content

Generating Rotated Images to Save Space

Kirill Kouzoubov edited this page Jun 12, 2022 · 2 revisions

Image below shows NIR band of Lansat8 captured on 2021-07-01 displayed in QGIS.

screenshot-qgis-ls8-nir-rotated-30m

This image was loaded using STAC metadata with odc-stac==0.3.0. Exact query shown below:

catalog = pystac_client.Client.open("https://planetarycomputer.microsoft.com/api/stac/v1")
items = catalog.search(
    collections=["landsat-8-c2-l2"],
    datetime="2021-07-01T08:00:00Z/2021-07-01T09:00:00Z",
    bbox=(-180, -50, 180, 50))

Compute Polygon of the pass in EPSG:3857:

ls8_pass = geom.unary_union(
    geom.Geometry(item.geometry, "epsg:4326").to_crs("epsg:3857") for item in items
)

Constructing rotated GeoBox that fits the data:

  • rotate geometry by appropriate amount for the data (about 13 degrees)
  • construct axis aligned geobox in rotated space with desired resolution
  • rotate geobox the other way
gbox = Affine.rotation(-a) * GeoBox.from_geopolygon(
    ls8_pass.transform(Affine.rotation(a)),
    resolution=res,
)

Full code of the script is in docs/samples/save-cog-from-stac.py. You can run it on Planetary Computer or "at home", although I'd recommend reducing resolution if running outside of Azure cloud.

pip install odc-stac==0.3.0
wget https://raw.githubusercontent.com/opendatacube/odc-stac/c0cb7096de38632c333caabe59da7a9e0cf33c9e/docs/samples/save-cog-from-stac.py
python save-cog-from-stac.py
Clone this wiki locally