Calculate NDVI

Learn to calculate a vegetation index to assess the quality of vegetation within an area of interest

Normalized Difference Vegetation Index (NDVI), developed by a NASA scientist named Compton Tucker in 1977, is commonly used to assess whether an area contains live green vegetation or not. It can show the difference between water and plants, bare soil and grass, whether plants are under stress, and what lifecycle stage a crop is in.

It compares how much more near-infrared light is reflected by chlorophyll vs visible red light:

In this guide we will perform a basic NDVI calculation in python, using Planet's 4-Band imagery.

You will need,

To export your Planet API Key to your environment:

export PLANET_API_KEY=a3a64774d30c4749826b6be445489d3b # (not a real key)

You can find Numpy installation instructions here.

To install the Python libraries, do:

pip install rasterio
pip install matplotlib
pip install retrying
pip install requests

Once installation is complete, you're ready to start with Step 1 below.

A note about Planet Imagery Product Compatability: this tutorial is compatable with only analytic and basic_analytic asset types of the following ItemTypes:

  • PSOrthoTile
  • REOrthoTile
  • PSScene4Band

Download a 4-Band image

First, download a 4-band Planetscope image of agricultural land in Toledo, Spain (id: 20161218_101700_0e0d). You can do this using the Planet API, or with Planet Explorer, by filtering for '4 Band PlanetScope scene' (PSScene4Band) or 'Planetscope ortho tile' (PSOrthoTile), and downloading an 'analytic' asset.

Using the Planet CLI:

planet data download --item-type PSScene4Band --asset-type analytic,analytic_xml --string-in id 20161218_101700_0e0d

You now have a file called '20161228_101647_0e26_3B_AnalyticMS.tif' in your directory.


Thumbnail preview of the 4-band scene

Extract the Visible Red and NIR bands

In a python script, open the image and extract just the Red and Near-infrared bands, numbers 3 and 4.

import rasterio
import numpy

image_file = "20161228_101647_0e26_3B_AnalyticMS.tif"

# Load red and NIR bands - note all PlanetScope 4-band images have band order BGRN
with rasterio.open(image_file) as src:
    band_red = src.read(3)

with rasterio.open(image_file) as src:
    band_nir = src.read(4)

Normalize to Top of Atmosphere Reflectance

Converting the pixel values to TOA Reflectance makes the analysis more accurate, and comparable with other scenes. Load the TOA Reflectance coefficients from the metadata XML asset.

from xml.dom import minidom

xmldoc = minidom.parse("20161218_101700_0e0d_3B_AnalyticMS_metadata.xml")
nodes = xmldoc.getElementsByTagName("ps:bandSpecificMetadata")

# XML parser refers to bands by numbers 1-4
coeffs = {}
for node in nodes:
    bn = node.getElementsByTagName("ps:bandNumber")[0].firstChild.data
    if bn in ['1', '2', '3', '4']:
        i = int(bn)
        value = node.getElementsByTagName("ps:reflectanceCoefficient")[0].firstChild.data
        coeffs[i] = float(value)

Multiply the band values by the TOA Reflectance coefficients.

# Multiply by corresponding coefficients
band_red = band_red * coeffs[3]
band_nir = band_nir * coeffs[4]

Perform the NDVI calculation

Next we perform the NDVI calculation through subtraction and division of the pixel values.

# Allow division by zero
numpy.seterr(divide='ignore', invalid='ignore')

# Calculate NDVI
ndvi = (band_nir.astype(float) - band_red.astype(float)) / (band_nir + band_red)

Save the NDVI image

Finally we output these new pixel values to a new image file, making sure we mirror the GeoTIFF spatial metadata:

# Set spatial characteristics of the output object to mirror the input
kwargs = src.meta
kwargs.update(
    dtype=rasterio.float32,
    count = 1)

# Create the file
with rasterio.open('ndvi.tif', 'w', **kwargs) as dst:
        dst.write_band(1, ndvi.astype(rasterio.float32))

Apply a color map

Applying a color map can help to visually distinguish vegetation. For more color map options check out the docs for matplotlib.

import matplotlib.pyplot as plt
plt.imsave("ndvi_cmap.png", ndvi, cmap=plt.cm.summer)

Questions or comments about this guide? Join the conversation at Planet Community.