Using Overloaded Operators to Compute Snow Cover (Sentinel-2, MOD10A1)

Tutorial created by **David Montero Loaiza**: GitHub | Twitter

GitHub Repo: https://github.com/davemlz/eemont
PyPI link: https://pypi.org/project/eemont/
Conda-forge: https://anaconda.org/conda-forge/eemont
Documentation: https://eemont.readthedocs.io/
More tutorials: https://github.com/davemlz/eemont/tree/master/docs/tutorials

Let’s start!

If required, please uncomment:

[1]:

#!pip install eemont
#!pip install geemap

Import the required packages.

[2]:

import ee, eemont, geemap

Authenticate and Initialize Earth Engine and geemap.

[3]:

Map = geemap.Map()

Point of interest.

[4]:

point = ee.Geometry.Point([-76.0269,2.92846])

Get, filter, mask clouds and scale the image collection.

[5]:

S2 = (ee.ImageCollection('COPERNICUS/S2_SR')
      .filterBounds(point)
      .sort('CLOUDY_PIXEL_PERCENTAGE')
      .first()
      .maskClouds()
      .scaleAndOffset()
      .spectralIndices('NDSI')) # Let's compute the NDSI, we'll need it!

Let’s select the required bands:

[6]:

NDSI = S2.select('NDSI')
N = S2.select('B8')
G = S2.select('B3')

Overloaded Operators

eemont has overloaded the binary operators, rich comparisons and unary operators in the following list for the ee.Image class:

(+, -, *, /, //, %, **, <<, >>, &, |, <, <=, ==, !=, >, >=, -, ~)

Therefore, you can now use them for image operations!

The following line computes snow cover according to (Hall et al., 2001):

[7]:

snowPixels = (NDSI > 0.4) & (N >= 0.1) & (G > 0.11) # Overloaded operators used here: >, >=, &.

Now, update the mask of the NDSI.

[8]:

NDSI = NDSI.updateMask(snowPixels)

Let’s save the date of the image to get the closest MOD10A1 image for comparison:

[9]:

dateOfInterest = ee.Date(S2.get('system:time_start'))

Get, filter and scale the MOD10A1 product:

[10]:

MOD10A1 = (ee.ImageCollection('MODIS/006/MOD10A1')
           .filterBounds(point)
           .closest(dateOfInterest)
           .scaleAndOffset() # NEW! Note that the scaleAndOffset() method supports the MOD10A1 product!
           .first())

This product already has the snow cover, therefore, we just need to compute the NDSI operation for snow cover pixels (Let’s use overloaded operators!):

[11]:

NDSI_MODIS = MOD10A1.select('NDSI')
NDSI_MODIS = NDSI_MODIS.updateMask(NDSI_MODIS > 0.4) # The overloaded operator > is used here!

Visualization

Let’s define the NDSI visualization parameters:

[12]:

visNDSI = {
    'min':0.4,
    'max':1,
    'palette': ['000000', '0dffff', '0524ff', 'ffffff']
}

Let’s define the RGB visualization parameters:

[13]:

visRGB = {
    'min':0,
    'max':0.3,
    'bands':['B4', 'B3', 'B2']
}

Use geemap to display results:

[14]:

Map.addLayer(S2,visRGB,'Sentinel-2 RGB')
Map.addLayer(NDSI_MODIS,visNDSI,'MODIS NDSI')
Map.addLayer(NDSI,visNDSI,'Sentinel-2 NDSI')
Map.add_colorbar(visNDSI['palette'],caption = 'NDSI')
Map.centerObject(point,13)
Map

Using Overloaded Operators to Compute the EVI (Sentinel-2)

Overloaded Operators for the ee.Number Object Class