nasa_hls icon indicating copy to clipboard operation
nasa_hls copied to clipboard

Implement QA-to-binary-mask functionality

Open benmack opened this issue 5 years ago • 0 comments

What we want:

Flag pixels as invalid if the QA layer defines them as being any of the following:

  • cirrus
  • cloud
  • adjecent cloud
  • cloud shadow
  • snow/ice

Understanding the QA Layer

The example from the User Guide:

Suppose we get a decimal QA value 100, which translates into binary 01100100, indicating that the aerosol level is low (bits 6-7), it is water (bit 5), and adjacent to cloud (bit 2).

This is how we find out which decimal QA values are valid (see also Table 10 of HLS v1.4 UserGuide):

i = 100  # decimal QA value
binary_str = "{0:08b}".format(i)
lut_qa.loc[i, "binary"] = binary_str
if all(np.array(list(binary_str[-5:])) == valid):
    lut_qa.loc[i, "valid"] = True
print(binary_str)
print(lut_qa.loc[i, "valid"])
01100100
False
lut_qa = pd.DataFrame({"integer":list(range(256)), 
                       "binary": [""] * 256,
                       "valid": [False]*256})
# if all of the last 5 binary digits are 0 then the pixel is valid   
# https://hls.gsfc.nasa.gov/wp-content/uploads/2018/10/HLS.v1.4.UserGuide_draft_ver3.0_clean.pdf
valid = np.array(list("00000"))
for i in range(256):
    binary_str = "{0:08b}".format(i)
    lut_qa.loc[i, "binary"] = binary_str
    if all(np.array(list(binary_str[-5:])) == valid): 
        lut_qa.loc[i, "valid"] = True

lut_valid = lut_qa[lut_qa["valid"]]
qa_valid_values = lut_valid["integer"].values
display(lut_valid)
print(qa_valid_values)
integer binary valid
0 00000000 True
32 00100000 True
64 01000000 True
96 01100000 True
128 10000000 True
160 10100000 True
192 11000000 True
224 11100000 True

ARE THESE REALLY THE CORRECT VALUES ??? - if so :

Function to read a QA layer, create a binary mask and write the result back as raster:

def hls_qa_to_clearsky_mask(qa_path, clear_path, overwrite=False):
    clear_sky_int = [0, 32, 64, 96, 160, 192, 224]

    with rasterio.open(qa_path) as qa:
        meta = qa.meta
        qa_array = qa.read(1)

    if not clear_path.exists() or overwrite:
        with rasterio.open(qa_path) as qa:
            meta = qa.meta
            qa_array = qa.read()

        meta.update(compress="lzw")
        clear_array = np.zeros_like(qa_array)
        for num in qa_valid_values:
            clear_array[qa_array==num] = 1
        with rasterio.open(clear_path, "w", **meta) as dst:
            dst.write(clear_array)

benmack avatar Mar 05 '19 18:03 benmack