earthdaily / earthdaily-python-client

EarthDaily python client
https://earthdaily.github.io/earthdaily-python-client/
MIT License
11 stars 6 forks source link

Report nan pixels #11

Open robmarkcole opened 10 months ago

robmarkcole commented 10 months ago

In addition to calculating cloud pixels I find it essential to calc the nan pixels, I am using the function below which can be adapted:

def get_item_cover(item: Item, bbox: List[float]):

    """
    Attempt to get the nodata cover & cloud cover from the item within the bbox region.

    Args:
        item: A STAC item

    Returns:
        A dictionary containing the nodata cover and cloud cover
    """
    try:
        dataset = stac_load(
            [item],
            bands=("scl"),
            crs="epsg:3857",
            resolution=10,  # the resolution of the output image in metres
            groupby="id",
            bbox=bbox,
        )
        scl_array = dataset.scl.values[0, :, :]
        nan_pixels = (scl_array == 0).sum() # count the number of no data pixels in scl_array
        all_pixels = scl_array.shape[0] * scl_array.shape[1]
        nan_pct = (nan_pixels / all_pixels) * 100
        nan_pct = round(nan_pct, 1)
        non_zero_pixels = all_pixels - nan_pixels
        cloud_pixels = ((scl_array == 8) | (scl_array == 9)).sum()
        cloud_pct = (cloud_pixels / non_zero_pixels) * 100
        cloud_pct = round(cloud_pct, 1)
        return {
            "nan_pct": nan_pct,
            "cloud_pct": cloud_pct
        }
    except Exception as e:
        print(f"Failed to get cloud cover for item {item.id}")
        return None
nkarasiak commented 10 months ago

Today the clear pixels (% and number) computation is only done when you choose to mask with a cloudmask your dataset and you ask for it (_maskstatistics=True when you create the datacube). There is a specific part in earthdaily py that manages that according to the cloudmask you select (and it is generic so you can add easily new mask providers) : https://github.com/GEOSYS/earthdaily-python-client/blob/main/earthdaily/earthdatastore/mask/__init__.py and nan pixels (gerally mainly those outside the bbox/geometry) are not counted in it. It is the clear percent and number of clear pixels inside your geometry. The nan size per date is so the total available pixels in the geometry minus the x*y size (except if you have nan pixels values inside your geometry).