From GIS to Remote Sensing: From Image Download to NDVI Calculation in One Move: SCP Batch

This tutorial describes the easy steps for downloading multiple satellite images (Download images), converting the DN values to reflectance, and automatically calculate NDVI (i.e. Normalized Difference Vegetation Index).

In particular, the SCP tools are intertwined through options in the plugin interface, allowing for the definition of an automatic workflow. The same phases of this tutorial can be adapted to automatically calculate any spectral index for several satellite images. Finally, some of the Batch functions are illustrated, which allow for the processing of already downloaded images.

Following the video of this tutorial.

1. Search the Images

First we need to search and select the satellite images to be downloaded and automatically processed. For more details about Sentinel-2 and Landsat images see Data Download and Data Download and Conversion.

In this tutorial we are going to use Landsat 8 images (data available from the U.S. Geological Survey), but of course the same method can be applied to Sentinel-2 or ASTER images.

Start a new QGIS project. Open the tab Download images clicking the button

in the SCP menu, or the SCP Tools, or the SCP dock. Select the tab Landsat download.

In Login https://ers.cr.usgs.gov/ you should enter the user name and password for accessing data (free registration at USGS EROS is required) in User and Password. Landsat 8 images are largely available for download from the Amazon Web Services (AWS) that does not require to login.

In Search area click the button

and draw a rectangle in the map to define the area coordinates.

Now click the button Find

and after a few seconds the image will be listed in the Image list.

Landsat search result

In the result table, click some of the items and click the button

. A preview for each image will be downloaded and displayed in the map, which is useful for assessing the quality of the image and the cloud cover. Remove the previews to be excluded from download, and leave loaded in QGIS only the images that we want to process.

Image previews

2. Download and Process Data

Now we can define the settings of the tools used for download and processing. It is worth noticing the options that link one tool to the following one, from the download, to the processing, to the Band set creation, and to the NDVI calculation.

Click the tab Download options and check only bands 4 (Red) and 5 (Near Infrared) (NDVI calculation requires Red and Near Infrared bands, therefore we are going to download only these two bands). Also, leave checked the options

Only if preview in Layers and

Preprocess images (this is the link to the tool Landsat).

Download settings

Now we define the options for the conversion to reflectance. Open the tab Landsat, check

Apply DOS1 atmospheric correction and

Create Band set and use Band set tools (this is the link to the tool Band set allowing for the creation of a Band set for every image).

Preprocessing settings

Open the tab Band set clicking the button

in the SCP menu, or the SCP Tools, or the SCP dock. Check

Band calc expression, which is the link to the tool Band calc.

Band set settings

Now we can define the expression used for calculation based on the Band set. Open the tab

Band calc, in Index calculation select NDVI, which enters the following Expression:

(( "#NIR#" - "#RED#") / ( "#NIR#" + "#RED#") @NDVI)

The “#NIR#” (Near Infrared) and “#RED#” bands are automatically identified from the Band set.

Band calc settings

In order to start the download and conversion process, open the tab Landsat download, click the button

and select the directory where converted bands are saved (e.g. Desktop). After a few minutes, converted bands are loaded and displayed (file name starts with RT_). Also, the NDVI calculation is performed for all the downloaded images and a file NDVI.tifsaved in the same directory of converted images.

NDVI calculation

3. Batch Data Processing

In case we have previously downloaded images, we can use the Batch functions to automatically perform the preprocessing, to the Band set creation, and to the NDVI calculation.

Open the tab Batch and select !working_dir! in Functions, which is the list of all available functions.

The text !working_dir!;'' is added to the text, which is the working directory; after the semicolon (following the name of the function) enter the path to the image directory such as:

!working_dir!;'!working_dir!;'/home/user/Desktop/'

The !working_dir! is used in the following functions where the relative path to the images is entered.

In Functions click landsat_conversion to enter all the options of the tool Landsat. After input_dir : enter the path to the image directory inside ' '; after apply_dos1 : enter 1 to apply the DOS1 correction (for checkboxes, value 1 means checked and value 0 means unchecked); after output_dir : enter the output directory for converted bands (not existing directory is automatically created), such as:

landsat_conversion;input_dir : '!working_dir!LC81900312014284LGN00/';mtl_file_path : '';celsius_temperature : 0;apply_dos1 : 1;use_nodata : 1;nodata_value : 0;pansharpening : 0;create_bandset : 1;output_dir : '!working_dir!/LC81900312014284LGN00/converted'

The folder converted is created inside the !working_dir!. The option create_bandset : 1 creates the Band set required for the calculation.

In Functions click band_calc; after expression : enter ( "#NIR#" - "#RED#" ) / ( "#NIR#" + "#RED#" ); after

output_raster_path:`!working_dir!/LC81900312014284LGN00/converted/NDVI.tif`; remove entirely the option ``extent_same_as_raster_name : '';

such as:

band_calc;expression : '( "#NIR#" - "#RED#" ) / ( "#NIR#" + "#RED#" )';output_raster_path : '!working_dir!/LC81900312014284LGN00/converted/LC81900312014284LGN00_NDVI.tif';align : 1;extent_intersection : 1;set_nodata : 0;nodata_value : 0

We can enter the same functions of conversion and calculation for other images. When the expression is green we can run the batch tool clicking the button

Batch expression

After the calculation the converted bands and NDVI rasters are displayed in the map.

NDVI rasters

4. Batch Data Calculation

In case we already have converted images, and we only need to perform the NDVI calculation, we can use the following function for creating the Band set (remove entirely the options ;multiplicative_factor : '';additive_factor : ''):

create_bandset;raster_path_list : '/home/user/Desktop/LC81900312014284LGN00/converted/RT_LC81900312014284LGN00_B4.TIF, /home/user/Desktop/LC81900312014284LGN00/converted/RT_LC81900312014284LGN00_B5.TIF';center_wavelength : '0.655,0.865';wavelength_unit : 1

The option raster_path_list is the list of paths to image bands (separated by ,). The option center_wavelength is required for the automatic identification of bands "#NIR#" and"#RED#". Now we can copy the same command used before for the NDVI calculation:

We can copy the same functions of Band set creation and calculation for other images. Finally, we could also clip the resulting NDVI rasters using upper left and lower right coordinates (identified in the map), adding the following function (remove entirely the options use_shapefile : 0;shapefile_path : '';):

clip_multiple_rasters;input_raster_path : '!working_dir!/LC81900312014284LGN00/converted/LC81900312014284LGN00_NDVI.tif,!working_dir!/LC81900312015079LGN00/converted/LC81900312015079LGN00_NDVI.tif ';output_dir : '!working_dir!/clip';use_shapefile : 0;ul_x : '346683';ul_y : '4622461';lr_x : '385483';lr_y : '4594931';nodata_value : 0;output_name_prefix : 'clip'

When the expression is green we can run the batch tool clicking the button

Batch expression

After the calculation the clipped NDVI rasters are displayed in the map.

NDVI clipped rasters

According to our needs, we can add other Batch functions and create an automatic workflow.

Moreover, with the tool Band calc we can create a new raster using several mathematical functions of NumPy , such as the function amax for returning the maximum value between input rasters.

We need to append np. before the NumPy function name such as:

band_calc;expression : 'np.amax(["LC81900312014284LGN00_NDVI.tif","LC81900312015079LGN00_NDVI.tif"],0)';output_raster_path : '!working_dir!/max.tif';align : 1;extent_intersection : 1;set_nodata : 0;nodata_value : 0

where np.amax is the function name, ["LC81900312014284LGN00_NDVI.tif","LC81900312015079LGN00_NDVI.tif"] is the list of input rasters, and the following ,0 is the axis of calculation (for comparing the corresponding pixels of the two input rasters). The tool Band calc will be better described in a future tutorial.

For any comment or question, join the Facebook group and the Google+ Community about the Semi-Automatic Classification Plugin.