Black stripes on mosaicking about stackstac HOT 3 CLOSED

Thanks for the good report @lhcorralo. It looks like the underlying issue is rasterio/rasterio#2916.

Short explanation: the xarray is basically composed of multiple tiles, spatially. It's also reading from multiple items, which cover different areas. When a tile doesn't overlap with an item at all, stackstac has an optimization to return all NaNs for the tile. But when the tile overlaps with an item, we use rasterio to read it. It looks like in this particular dataset, the GeoTIFFs don't have a nodata value set, so GDAL is filling in 0s for the pixels outside the bounds of the item, but not actually telling us that they should be masked out. However, the mosaic operation is only treating NaNs as nodata, so the 0s are left in, resulting in the stripes. That's also why max works: any actual data value is greater than 0.

I'll push up a fix for this, but for now, you can probably work around this with:

stack = stackstac.stack(items_pc, fill_value=0)
mosaic = stackstac.mosaic(stack, nodata=0)

That will make stackstac use 0s for empty tiles, matching up with GDAL's nodata value. Then you mosaic out all the 0s.

As a sidenote, this is yet another reason to move away from VRTs, which already needs to happen: #196.

from stackstac.

Thanks!

Yes, it works as a workaround: it is easy to get a fill_value/nodata for elevation values (I have seen -9999 sometimes), so no problem.

from stackstac.

FYI @lhcorralo, the fix is now released in 0.5.0, so you should be able to go back to using mosaic. Thanks for the great bug report!

from stackstac.