Stackstac.stack rescaling and landsat-c2-l2

[Jesse Anderson (jessjaco)]

The stackstac.stack rescale option doesn't seem to work for landsat-c2-l2, e.g.

import planetary_computer
import pystac_client
import stackstac

catalog = pystac_client.Client.open(
        "https://planetarycomputer.microsoft.com/api/stac/v1",
        modifier=planetary_computer.sign_inplace,
    )

items = catalog.search(collections=["landsat-c2-l2"],
                datetime="2021-01-01T00:00:00Z/2021-01-08T00:00:00Z",
                bbox=[-111.1, 35, -111, 35.1]).item_collection()


stackstac.stack(items, rescale=True).isel(time=0).sel(band="blue").min().values

Gives a value of 3205. I see the correct scale and offset values in the collection metadata. Do I need to perform the calculations myself?

I can ask on the stackstac issues page too, just not sure whose issue this is.

[Tom Augspurger (TomAugspurger)]

Mmm, I'm not 100% what rescale=True means in stackstac. I see https://github.com/gjoseph92/stackstac/blob/9106708bcc20daec4c5975e9c1240de38c38f2f1/stackstac/stack.py#L204-L208

IIUC, that's actually getting the scales and offsets from actual data assets (COGs in this case) not the STAC metadata. That's based off https://github.com/gjoseph92/stackstac/blob/9106708bcc20daec4c5975e9c1240de38c38f2f1/stackstac/rio_reader.py#L98.

import rasterio
import pystac
import planetary_computer

item = pystac.read_file("https://planetarycomputer.microsoft.com/api/stac/v1/collections/landsat-c2-l2/items/LC08_L2SR_085082_20220323_02_T2")
ds = rasterio.open(planetary_computer.sign(item.assets["blue"].href))

print("STAC:", item.assets["blue"].extra_fields["raster:bands"][0])
print(" COG:", ds.scales, ds.offsets)

prints out

STAC: {'scale': 2.75e-05, 'nodata': 0, 'offset': -0.2, 'data_type': 'uint16', 'spatial_resolution': 30}
 COG: (1.0,) (0.0,)

IIUC, these values are coming from the stactools source: https://github.com/stactools-packages/landsat/blob/5cfff89ab39030aef6c1dc5c818453469963cd76/src/stactools/landsat/fragments/collections/landsat-c2-l2.json#L372-L380

So I would have two questions:

1. Is there a problem with the COGs generated by USGS not including scale / offset appropriately? (the landsat-c2-l2 COGs are unmodified from what USGS produces)
2. Should stackstac look for scale / offset information in the raster:bands? Maybe... I'm not sure.

[Jesse Anderson (jessjaco)]

Ok, looks like the USGS COGs are the issue. I will scale manually.

[Gabe Joseph (gjoseph92)]

FYI, stackstac also wasn't always handling the scale and offset fields set on COGs correctly, but that's been fixed in gjoseph92/stackstac#183 (but not yet released).

Also, stackstac doesn't look at raster:bands (where the scale/offset metadata is set in STAC) at all. I'd like to do this at some point, but haven't had the time: gjoseph92/stackstac#91.

So short answer is that unless scale and offset are set correctly on the COG, you'll have to rescale manually for now.