Using AppEEARS Quality Service to Extract Information from MODIS Quality Layers


Quality information for LP DAAC (Land Processes Distributed Active Archive Center) data products is important to consider when determining the usability and usefulness of a data product for a particular science application. However, this information has been notoriously difficult for users to access. Quality information is often stored as an integer value that requires the user to decode into binary strings. In order to interpret the binary string users must map the unique combinations of bits found in separate subsets of the binary string (i.e. bit-fields) to quality tables that characterize the particular quality attribute associated with each bit-field. This process is difficult for most users and is often passed over.

The LP DAAC has added a quality web service to its AppEEARS service API to alleviate the difficult process of decoding and interpreting quality layers for tiled MODIS and WELD products. The quality web service will decode an integer value for a data product/layer and return the quality attribute information associated with that particular integer value. Additionally, users can determine if a quality layer exists for a particular data product/layer as well as get details about individual quality layers.

NOTE

  • Users of this tutorial will need internet access in order to successfully execute the python code.
  • This is a python specific tutorial and has been tested in Python 3.4 and 2.7.

Let's start by importing the required packages.

In [1]:
import requests
import pandas as pd

Next, assign the AppEEARS services API URL to a variable

In [2]:
SERVICES_URL = 'https://lpdaacsvc.cr.usgs.gov/services/appeears-api/'

NOTE: Placing https://lpdaacsvc.cr.usgs.gov/services/appeears-api/ in your browser will take you to the AppEEARS services API page. The page contains all of the AppEEARS services accompanied by their documentation.

Get all of the data layers that contain associated quality layers


Say we wish to know what layers have a quality layer associated with them. We can find that out! To start, let's print out a count of the number of layers with associated quality layers.

In [3]:
qualityLayers = requests.get('{}/quality?format=json'.format(SERVICES_URL)).json()
#print(qualityLayers)
print('Quality information was found for {0} data layers!'.format(len(qualityLayers)))
Quality information was found for 454 data layers!

The request returns a list of dictionaries that can be indexed to filter the information returned

In [4]:
qualityLayers[218:222]
Out[4]:
[{'Layer': '_500m_16_days_EVI',
  'ProductAndVersion': 'MOD13A1.005',
  'QualityLayers': ['_500m_16_days_VI_Quality'],
  'QualityProductAndVersion': 'MOD13A1.005'},
 {'Layer': '_500m_16_days_NDVI',
  'ProductAndVersion': 'MOD13A1.005',
  'QualityLayers': ['_500m_16_days_VI_Quality'],
  'QualityProductAndVersion': 'MOD13A1.005'},
 {'Layer': '_500m_16_days_EVI',
  'ProductAndVersion': 'MOD13A1.006',
  'QualityLayers': ['_500m_16_days_VI_Quality'],
  'QualityProductAndVersion': 'MOD13A1.006'},
 {'Layer': '_500m_16_days_NDVI',
  'ProductAndVersion': 'MOD13A1.006',
  'QualityLayers': ['_500m_16_days_VI_Quality'],
  'QualityProductAndVersion': 'MOD13A1.006'}]

The list of products contained in qualityLayers is long, and trying to find the product/layer of interest via indexing can be inefficient. Instead, let's use a python expression, known as a list comprehension, to get all of the layers in MYD13A2.005 that have quality layers associated with them.

In [5]:
productQuality_A = [q for q in qualityLayers if q['ProductAndVersion'] == 'MYD13A2.005']
productQuality_A
Out[5]:
[{'Layer': '_1_km_16_days_EVI',
  'ProductAndVersion': 'MYD13A2.005',
  'QualityLayers': ['_1_km_16_days_VI_Quality'],
  'QualityProductAndVersion': 'MYD13A2.005'},
 {'Layer': '_1_km_16_days_NDVI',
  'ProductAndVersion': 'MYD13A2.005',
  'QualityLayers': ['_1_km_16_days_VI_Quality'],
  'QualityProductAndVersion': 'MYD13A2.005'}]

If you don't want to use indexing or list comprehensions (the previous two steps), just use the product parameter in the quality service URL to get all of the layers in MYD13A2.005 that have associated quality layers.

In [6]:
productQuality_B = requests.get('{}/quality/MYD13A2.005?format=json'.format(SERVICES_URL)).json()
productQuality_B
Out[6]:
[{'Layer': '_1_km_16_days_EVI',
  'ProductAndVersion': 'MYD13A2.005',
  'QualityLayers': ['_1_km_16_days_VI_Quality'],
  'QualityProductAndVersion': 'MYD13A2.005'},
 {'Layer': '_1_km_16_days_NDVI',
  'ProductAndVersion': 'MYD13A2.005',
  'QualityLayers': ['_1_km_16_days_VI_Quality'],
  'QualityProductAndVersion': 'MYD13A2.005'}]

NOTE: In both cases above, we can see that EVI (_1_km_16_days_EVI) and NDVI (_1_km_16_days_NDVI) have the same quality layer (_1_km_16_days_VI_Quality). This is the case for many of the MODIS products. However, there are occasions when the layers within products have their own associated quality layers (e.g., MOD/MYD11 products).

Get the quality layer description for MYD13A2.005


Information within quality layers is often stored as bit-packed integers. This means that in order to get the actual quality information, one must convert the integer to binary and parse the binary string into bit-fields. Bit-fields are distinct combinations of bits representing quality categories that give an indication as to the usability and usefulness of the data product. Let's get the quality layer description for MYD13A2.005 by adding the quality layer parameter to the quality service URL.

In [7]:
qualityLayerInfo = requests.get('{}/quality/MYD13A2.005/_1_km_16_days_VI_Quality?format=json'.format(SERVICES_URL)).json()
#qualityLayerInfo
pd.DataFrame(qualityLayerInfo)
Out[7]:
Acceptable Description Name ProductAndVersion QualityLayer Value
0 True VI produced, good quality MODLAND MYD13A2.005 _1_km_16_days_VI_Quality 0
1 False VI produced, but check other QA MODLAND MYD13A2.005 _1_km_16_days_VI_Quality 1
2 False Pixel produced, but most probably cloudy MODLAND MYD13A2.005 _1_km_16_days_VI_Quality 2
3 False Pixel not produced due to other reasons than c... MODLAND MYD13A2.005 _1_km_16_days_VI_Quality 3
4 None Higest quality VI Usefulness MYD13A2.005 _1_km_16_days_VI_Quality 0
5 None Lower quality VI Usefulness MYD13A2.005 _1_km_16_days_VI_Quality 1
6 None Decreasing quality VI Usefulness MYD13A2.005 _1_km_16_days_VI_Quality 2
7 None Decreasing quality VI Usefulness MYD13A2.005 _1_km_16_days_VI_Quality 3
8 None Decreasing quality VI Usefulness MYD13A2.005 _1_km_16_days_VI_Quality 4
9 None Decreasing quality VI Usefulness MYD13A2.005 _1_km_16_days_VI_Quality 5
10 None Decreasing quality VI Usefulness MYD13A2.005 _1_km_16_days_VI_Quality 6
11 None Decreasing quality VI Usefulness MYD13A2.005 _1_km_16_days_VI_Quality 7
12 None Decreasing quality VI Usefulness MYD13A2.005 _1_km_16_days_VI_Quality 8
13 None Decreasing quality VI Usefulness MYD13A2.005 _1_km_16_days_VI_Quality 9
14 None Decreasing quality VI Usefulness MYD13A2.005 _1_km_16_days_VI_Quality 10
15 None Decreasing quality VI Usefulness MYD13A2.005 _1_km_16_days_VI_Quality 11
16 None Lowest quality VI Usefulness MYD13A2.005 _1_km_16_days_VI_Quality 12
17 None Quality so low that it is not useful VI Usefulness MYD13A2.005 _1_km_16_days_VI_Quality 13
18 None L1B data faulty VI Usefulness MYD13A2.005 _1_km_16_days_VI_Quality 14
19 None Not useful for any other reason/not processed VI Usefulness MYD13A2.005 _1_km_16_days_VI_Quality 15
20 None Climatology Aerosol Quantity MYD13A2.005 _1_km_16_days_VI_Quality 0
21 None Low Aerosol Quantity MYD13A2.005 _1_km_16_days_VI_Quality 1
22 None Average Aerosol Quantity MYD13A2.005 _1_km_16_days_VI_Quality 2
23 None High Aerosol Quantity MYD13A2.005 _1_km_16_days_VI_Quality 3
24 None No Adjacent cloud detected MYD13A2.005 _1_km_16_days_VI_Quality 0
25 None Yes Adjacent cloud detected MYD13A2.005 _1_km_16_days_VI_Quality 1
26 None No Atmosphere BRDF Correction MYD13A2.005 _1_km_16_days_VI_Quality 0
27 None Yes Atmosphere BRDF Correction MYD13A2.005 _1_km_16_days_VI_Quality 1
28 None No Mixed Clouds MYD13A2.005 _1_km_16_days_VI_Quality 0
29 None Yes Mixed Clouds MYD13A2.005 _1_km_16_days_VI_Quality 1
30 None Shallow ocean Land/Water Flag MYD13A2.005 _1_km_16_days_VI_Quality 0
31 None Land (Nothing else but land) Land/Water Flag MYD13A2.005 _1_km_16_days_VI_Quality 1
32 None Ocean coastlines and lake shorelines Land/Water Flag MYD13A2.005 _1_km_16_days_VI_Quality 2
33 None Shallow inland water Land/Water Flag MYD13A2.005 _1_km_16_days_VI_Quality 3
34 None Ephemeral water Land/Water Flag MYD13A2.005 _1_km_16_days_VI_Quality 4
35 None Deep inland water Land/Water Flag MYD13A2.005 _1_km_16_days_VI_Quality 5
36 None Moderate or continental ocean Land/Water Flag MYD13A2.005 _1_km_16_days_VI_Quality 6
37 None Deep ocean Land/Water Flag MYD13A2.005 _1_km_16_days_VI_Quality 7
38 None No Possible snow/ice MYD13A2.005 _1_km_16_days_VI_Quality 0
39 None Yes Possible snow/ice MYD13A2.005 _1_km_16_days_VI_Quality 1
40 None No Possible shadow MYD13A2.005 _1_km_16_days_VI_Quality 0
41 None Yes Possible shadow MYD13A2.005 _1_km_16_days_VI_Quality 1

NOTE: The table above contains the description for each bit-field in MYD13A2.005's quality layer (i.e., _1_km_16_days_VI_Quality). For more information of quality layer descriptions, please consult with the product pages on the LP DAAC website.

Decode Quality Value


Say that we have extracted a pixel value from the EVI layer (i.e., _1_km_16_days_EVI) of product MYD13A2.005. We can use the quality service to decode the associated quality value for the pixel by inserting the pixel value from the quality layer (_1_km_16_days_VI_Quality) in the quality service URL. We see that our pixel value from the quality layer equals 4160 and that our EVI value equals 0.3091.

In [8]:
qualityIntDecoder = requests.get('{}/quality/MYD13A2.005/_1_km_16_days_VI_Quality/4160?format=json'.format(SERVICES_URL)).json()
#qualityIntDecoder
pd.DataFrame(qualityIntDecoder)
Out[8]:
Adjacent cloud detected Aerosol Quantity Atmosphere BRDF Correction Binary Representation Land/Water Flag MODLAND Mixed Clouds Possible shadow Possible snow/ice VI Usefulness
bits 0b0 0b01 0b0 0b0001000001000000 0b010 0b00 0b0 0b0 0b0 0b0000
description No Low No 0b0001000001000000 Ocean coastlines and lake shorelines VI produced, good quality No No No Higest quality

The above code snippet returns the decoded quality information for value 4160, however, the order in which the quality categories are displayed does not match the order in which the bit-fields are assembled in the binary string. To correct this, the same service call can be used, but the manner in which we bring the request into python changes.

In [9]:
import json
from collections import OrderedDict
qualityIntDecoder = json.loads(
    requests.get('{}/quality/MYD13A2.005/_1_km_16_days_VI_Quality/4160?format=json'.format(SERVICES_URL)).text, 
    object_pairs_hook=OrderedDict
)

Now let's print out the decoded quality value in its proper order.

In [10]:
#qualityIntDecoder
pd.DataFrame(qualityIntDecoder)
Out[10]:
Binary Representation MODLAND VI Usefulness Aerosol Quantity Adjacent cloud detected Atmosphere BRDF Correction Mixed Clouds Land/Water Flag Possible snow/ice Possible shadow
bits 0b0001000001000000 0b00 0b0000 0b01 0b0 0b0 0b0 0b010 0b0 0b0
description 0b0001000001000000 VI produced, good quality Higest quality Low No No No Ocean coastlines and lake shorelines No No

We can see in the MODLAND column that this is a good quality pixel and the Land/Water Flag column indicates that this pixel is likely near a lake shoreline (given that the pixel extracted was from tile h10v04).

Finally, we can take the decoded python object (a python dictionary) and filter it based on the bit-field description.

In [11]:
qualityIntDecoder['MODLAND']
Out[11]:
OrderedDict([('bits', '0b00'), ('description', 'VI produced, good quality')])