Masking data using the `Fmask` cloud mask

Sign up to the DEA Sandbox to run this notebook interactively from a browser
Compatibility: Notebook currently compatible with both the NCI and DEA Sandbox environments
Products used: ga_ls8c_ard_3

Background

In the past, remote sensing researchers would reject partly cloud-affected scenes in favour of cloud-free scenes. However, multi-temporal analysis techniques increasingly make use of every quality assured pixel within a time series of observations. The capacity to automatically exclude low quality pixels (e.g. clouds, shadows or other invalid data) is essential for enabling these new remote sensing methods.

Analysis-ready satellite data from Digital Earth Australia includes pixel quality information that can be used to easily “mask” data (i.e. keep only certain pixels in an image) to obtain a time series containing only clear or cloud-free pixels.

Description

In this notebook, we show how to mask Digital Earth Australia satellite data using boolean masks. The notebook demonstrates how to:

Load in a time series of satellite data including the fmask pixel quality band
Inspect the band’s flags_definition attributes
Create clear and cloud-free masks and apply these to the satellite data
Buffer cloudy/shadowed pixels by dilating the masks by a specified number of pixels
Clean cloudy/shadowed pixels to reduce false positive cloud detection
Mask out invalid nodata values and replace them with nan

Getting started

First we import relevant packages and connect to the datacube. Then we define our example area of interest and load in a time series of satellite data.

[1]:

import scipy.ndimage
import xarray
import numpy
import datacube
from datacube.utils.masking import make_mask
from datacube.utils.masking import mask_invalid_data
from odc.algo import mask_cleanup

import sys
sys.path.insert(1, '../Tools/')
from dea_tools.plotting import rgb

Connect to the datacube

[2]:

dc = datacube.Datacube(app="Masking_data")

Cloud-free images

Creating the clear-pixel mask

We create a mask by specifying conditions that our pixels must satisfy. But we will only need the labels (e.g. fmask="valid") to create a mask.

[7]:

# Create the mask based on "valid" pixels
clear_mask = make_mask(data.fmask, fmask="valid")
clear_mask.plot(col="time", col_wrap=4)

[7]:

<xarray.plot.facetgrid.FacetGrid at 0x7fef0e8a8520>

../../../_images/notebooks_How_to_guides_Masking_data_17_1.png

Applying the clear-pixel mask

We can now get the clear images we want.

[8]:

# Apply the mask
clear = data.where(clear_mask)
rgb(clear, col="time")

../../../_images/notebooks_How_to_guides_Masking_data_19_0.png

Cloud-free images

If we look carefully, we can see that we have lost the ocean too. Sometimes we may instead want to create a mask using a combination of different fmask features. For example, below we create a mask that will preserve pixels that are flagged as either valid, water or snow (and mask out any cloud or cloud shadow pixels):

[9]:

# Identify pixels that are either "valid", "water" or "snow"
cloud_free_mask = (
    make_mask(data.fmask, fmask="valid") |
    make_mask(data.fmask, fmask="water") |
    make_mask(data.fmask, fmask="snow")
)

# Apply the mask
cloud_free = data.where(cloud_free_mask)
rgb(cloud_free, col="time")

../../../_images/notebooks_How_to_guides_Masking_data_22_0.png

Mask dilation and cleaning

Sometimes we want our cloud masks to be more conservative and mask out more than just the pixels that fmask classified as cloud or cloud shadow. That is, sometimes we want a buffer around the cloud and the shadow. We can achieve this by dilating the mask using the mask_cleanup function from odc.algo.

Because we now want to focus on cloud and cloud shadow pixels, we first do the opposite to our previous examples and create a mask which has True values if a pixel contains either cloud or cloud shadow, and False for all others (e.g. valid, snow, water). When we plot this data, cloud and cloud shadow will appear as yellow, and other pixels as purple.

[10]:

# Identify pixels that are either "cloud" or "cloud_shadow"
cloud_shadow_mask = (
    make_mask(data.fmask, fmask="cloud") |
    make_mask(data.fmask, fmask="shadow")
)

# Plot
cloud_shadow_mask.plot(col="time", col_wrap=4)

[10]:

<xarray.plot.facetgrid.FacetGrid at 0x7fef0c199be0>

../../../_images/notebooks_How_to_guides_Masking_data_25_1.png

We now apply mask_cleanup. This function allows us to modify areas of cloud and cloud shadow using image processing techniques called morphological operations. These tools can be used to expand and contract regions of our image to change the shape of specific features - in this case, clouds and cloud shadow. Four of the most useful morphological techniques include:

Dilation: Expand (i.e. “dilate”) True values outward, resulting in larger True features
Erosion: Shrink (i.e. “erode”) True values inward, resulting in smaller True features
Closing: First dilate, then erode True pixels. This is used to fill small or narrow False gaps inside or between True features.
Opening: First erode, then dilate True pixels. This is used to remove small or narrow areas of True features, but preserve larger features.

All of these operations are applied using a specific radius to control how many pixels our clouds and shadows are dilated, eroded, closed or opened. For example, we can specify that our clouds and shadows are expanded by 5 pixels in all directions as follows:

[11]:

# Dilate all cloud and cloud shadow pixels by 5 pixels in all directions
cloud_shadow_buffered = mask_cleanup(mask=cloud_shadow_mask,
                                     mask_filters=[("dilation", 5)])
cloud_shadow_buffered.plot(col="time", col_wrap=4)

[11]:

<xarray.plot.facetgrid.FacetGrid at 0x7feef9b65e50>

../../../_images/notebooks_How_to_guides_Masking_data_27_1.png

Our clouds and shadows (yellow) have now been expanded outwards (compare this to the cloud_shadow_mask data we plotted earlier).

We can now apply this dilated cloud and shadow mask to our original data (note we need to reverse our mask using ~ so that valid pixels are marked with True as in our original un-dilated example).

[12]:

# Apply the mask
buffered_cloud_free = data.where(~cloud_shadow_buffered)
rgb(buffered_cloud_free, col="time")

../../../_images/notebooks_How_to_guides_Masking_data_29_0.png

Cloud mask data from fmask can commonly include false positive clouds over features including urban areas, bright sandy beaches and salt pans. We can see an example of this in the sixth panel, where a length of beach has been erroneously mapped as cloud.

To reduce these false positives, it can often be useful to apply a morphological “opening” operation before we dilate our clouds and shadows. This operation can remove small or narrow clouds by first “shrinking” then “expanding” our cloud and shadow pixels. To apply an opening operation with a radius of 3, we can supply an extra ("opening", 3) processing step using mask_cleanup:

[13]:

# Dilate all cloud and cloud shadow pixels by 5 pixels in all directions
cloud_shadow_buffered = mask_cleanup(mask=cloud_shadow_mask,
                                     mask_filters=[("opening", 3), ("dilation", 5)])

# Apply the mask
buffered_cloud_free = data.where(~cloud_shadow_buffered)

Now when we plot our data, we can see that the sixth panel is no longer affected by false positive clouds along the bright narrow beach.

Note: Morphological operations and radiuses are application-specific; make sure to experiment with a range of options to make sure they have the effect you are looking for.

[14]:

rgb(buffered_cloud_free, col="time")

../../../_images/notebooks_How_to_guides_Masking_data_33_0.png

Additional information

License: The code in this notebook is licensed under the Apache License, Version 2.0. Digital Earth Australia data is licensed under the Creative Commons by Attribution 4.0 license.

Contact: If you need assistance, please post a question on the Open Data Cube Slack channel or on the GIS Stack Exchange using the open-data-cube tag (you can view previously asked questions here). If you would like to report an issue with this notebook, you can file one on GitHub.

Last modified: December 2023

Compatible datacube version:

[17]:

print(datacube.__version__)

1.8.6

Masking data using the `Fmask` cloud mask

Background

Description

Getting started

Connect to the datacube

Create a query and load satellite data

Cloud-free images

Creating the clear-pixel mask

Applying the clear-pixel mask

Cloud-free images

Mask dilation and cleaning

Masking out invalid data

Additional information

Tags

Masking data using the Fmask cloud mask

Background

Description

Getting started

Connect to the datacube

Create a query and load satellite data

Cloud-free images

Creating the clear-pixel mask

Applying the clear-pixel mask

Cloud-free images

Mask dilation and cleaning

Masking out invalid data

Additional information

Tags

Masking data using the `Fmask` cloud mask