Introduction

MODIStsp is a novel “R” package allowing to automatize the creation of time series of rasters derived from MODIS Land Products data. It allows to perform several preprocessing steps on MODIS data available within a given time period.

Development of MODIStsp started from modifications of the ModisDownload “R” script by Thomas Hengl (2010), and successive adaptations by Babak Naimi (2014). The basic functionalities for download and preprocessing of MODIS datasets provided by these scripts were gradually incremented with the aim of:

All processing parameters can be easily set with a user-friendly GUI, although non-interactive execution exploiting a previously created Options File is possible. Stand-alone execution outside an “R” environment is also possible, allowing to use scheduled execution of MODIStsp to automatically update time series related to a MODIS product and extent whenever a new image is available.

Installation

MODIStsp was developed completely in the “R” Language and Environment for Statistical Computing (R Development Core Team 2014) v. 3.1.3 and is distributed as open source software under the GNU-GPL 3.0 License. Source code can be downloaded at the GITHUB repository https://github.com/lbusett/MODIStsp

To Install MODIStsp:

  1. Within “R”, install the “gWidgetsRGtk2” package

    install.packages('gWidgetsRGtk2')

    when asked, request to install GTK

  2. Install the package from GitHub. (You’ll need to have the “devtools” package installed and loaded)

    install.packages('devtools')
    require('devtools')
    install_github('lbusett/MODIStsp', ref = 'master')

Dependencies

MODIStsp exploits functionalities of several additional “R” packages. A complete list of required packages is shown in Table I.

“R” package Version Authors
gWidgetsRGtk2 >= 0.0-54 (Lawrence and Verzani 2013)
rgdal >= 0.8-16 (Bivand and Rundel 2014)
plyr >= 1.8.1 (Wickham 2011)
gdalUtils >= 0.3.1 (Greenberg and Mattiuzzi 2014)
XML >= 3.98-1.1 (Lang 2013a)
RCurl >= 1.95-4.1 (Lang 2013b)
Stringr >= 0.6.2 (Wickham 2012)
hash >= 2.26 (Brown 2013)
raster >= 2.3.40 (Hijmans 2014)
rgeos >= 0.3-8 (Bivand, Keitt, and Rowlingson 2014)

Additionaly, MODIStsp requires availability of the GDAL - Geospatial Data Abstraction Library v. >1.10, with support for HDF4 raster format.

Running the tool

Interactive mode

To run the tool in interactive mode, load the package and launch the MODIS_tsp function, with no parameters

library(MODIStsp)
MODIStsp()

This opens a GUI from which processing options can be specified and eventually saved (or loaded). A description of the different processing parameters to be selected is reported in the “Selection of Processing parameters” subsection at the end of this vignette

Non-Interactive mode (Using a previously saved options file)

MODIStsp can be also launched in non-interactive mode by setting the optional “GUI” parameter to FALSE, and the “Options_File” parameter to the path of a previously saved Options file. This allows to exploit MODIStsp functionalities within generic “R” processing scripts

library(MODIStsp) 
# --> Specify the path to a valid options file saved in advance
options_file = "X:/yourpath/youroptions.RData"   from the GUI
MODIStsp(gui = FALSE, options_File = options_file)

Specifying also the “spatial_file_path” parameter overrides the output extent of the selected Options File. This allows to perform the same preprocessing on different extents using a single Options File, by looping on an array of spatial files representing the desired output extents.

For example:

# Create a character array containing a list of shapefiles (or other spatial files)
extent_list = list.files("X:/path/containing/some/shapefiles/", "\\.shp$")  

# loop on the list of spatial files and run MODIStsp using each of them to automatically 
# define the output extent (A separate output folder is created for each input spatial file).

for (single_shape in extent_list) {
  MODIStsp(gui = FALSE, options_File = "X:/yourpath/youroptions.RData", 
          spatial_file_path = single_shape )
}

Selection of processing parameters through the GUI

When running MODIStsp in Interactive Mode, a user-friendly GUI allows selection of all processing options required for the creation of the desired MODIS time series. The main available processing options are described briefly in the following.

Figure 1: The MODIStsp main GUI

MODIS Product, Satellites and Layers

Allows to select the MODIS product of interest from a drop-down menu. The user can also select which MODIS sensors should be considered for download and creation of the time series (“Terra”, “Aqua” or both). After selecting the product, pushing the “Select Layers” button opens the “Select Processing Layers” GUI panel (Figure 2), from which the user must select which MODIS original and/or derived QI and SI layers should be processed:

Figure 2: The “Select Processing Layers” GUI
  1. The left-hand frame allows to select which original MODIS layers should be processed
  2. The central frame allows to select which Quality Indicators should be extracted from the original MODIS Quality Assurance layers.
  3. For MODIS products containing surface reflectance data, the right-hand frame allows to select which additional Spectral Indexes should be computed.1

Some of the most commonly used Spectral Indexes are available for computation by default (Table II).

Index Acronym Index name and reference
NDVI Normalized Difference Vegetation Index (Rouse et al. 1973)
EVI Enhanced Vegetation Index (Huete et al. 2002)
SR Simple Ratio(Tucker 1979) ]
NDFI Normalized Difference Flood Index (Boschetti et al. 2014)
NDII7 (NDWI) Normalized Difference Infrared Index – Band 7 (Hunt and Rock 1989)
SAVI Soil Adjusted Vegetation Index (Huete 1988)
NDSI* Normalized Difference Snow Index ((Hall et al. 2002)
NDII6* Normalized Difference Infrared Index – band 6 (Hunt and Rock 1989)
GNDVI* Green Normalized Difference Vegetation Index (Gitelson and Merzlyak 1998)
RGRI* Red Green Ration index (Gamon and Surfus 1999)
GRVI* Green-red ratio vegetation index (Tucker 1979)

Users can however specify other SIs to be computed without modifying MODIStsp source code by clicking on the “Add Custom Index” button, which allow to provide info related to the new desired SI using a simple GUI (Figure 3) interface.

Figure 3: The GUI for insertion of additional Spectral Indexes

Provided information (e.g., correct bandnames, computable formula, etc…) is automatically checked, and the new index added in the list of available ones for all products allowing its computation at the next MODIStsp execution.

Processing Period

Allows to specify the starting and ending dates (dd/mm/yyyy) to be considered for the creation of the time series.

Spatial Extent

Allows to define the area of interest for the processing. Two main options are possible: .

  1. Full Tiles Extent: the user must specify which MODIS tiles he would like to process using the “Start” and “End” horizontal and vertical sliders in the “Required MODIS Tiles” frame. During processing, data from the different tiles is mosaiced, and a single file covering the total area is produced for each acquisition date (Note: pressing the “show map” button, a representation of the MODIS tiles grid is shown to facilitate the selection).

  2. Resized: the user can specify the spatial extent of the desired outputs either;

    1. Manually inserting the coordinates of the Upper Left and Lower Right corners of the area of interest in the “Bounding Box” frame. Coordinates of the corners must be provided in the coordinate system of the selected output projection.

    2. pressing the “Load Extent from a Spatial File” and selecting a raster or vector spatial file. In this case, the bounding box of the selected file is retrieved, converted in the selected output projection, and shown in the “Bounding Box” frame.Required input MODIS tiles are also automatically retrieved from the output extent, and the tiles selection sliders modified accordingly.

Reprojection and Resize

Allows to specify the options to be used for reprojecting and resizing the MODIS images. In particular:

  1. The “Output Projection” menu allows to either select one of the pre-defined output projections or specify a user-defined one by selecting “User Defined” and then inserting a valid “Proj4” string in the pop-up window. Validity of the Proj4 string is automatically checked, and error messages issued if the check fails.

  2. The “Output Resolution”, “Pixel Size” and “Reprojection Method” menus allow to specify whether output images should inherit their spatial resolution from the original MODIS files, or be resampled to a user-defined resolution. In the latter case, output spatial resolution must be specified in the measure units of the selected output projection. Resampling method can be selected among those supported by the “gdalwarp” routine (http://www.gdal.org/gdalwarp.html).

Processing Options

Allows first of all to specify the format desired for the output images. Two of the most commonly formats used in remote sensing applications are available at the moment: ENVI binary and GeoTiff. If GeoTiff is selected, the type of file compression can be also specified among “None”, “PACKBITS”, “LZW” and “DEFLATE”.

The user can then specify if virtual multitemporal files should be created. These virtual files allow access to the entire time series of images as a single file without the need of creating large multitemporal raster images. Available virtual files formats are ENVI metafiles and GDAL “vrt” files.

Finally, users can select if the NoData values of MODIS layers should be kept at their original values, or changed to those specified within the “MODIStsp_Products_Opts” XML file. By selecting “Yes” in the “Change Original NODATA values” checkbox, NoData of outputs are set to the largest integer value possible for the data type of the processed layer (e.g., for 8-bit unsigned integer layers, NoData is set always to 255, for 16-bit signed integer layers to 32767, and for 16-bit unsigned integer layers to 65535). .

Main Output Folder for Time Series Storage

Allows to specify the main folder where the pre-processed time series data will be saved stored. The “Reprocess Existing Data” checkbox allows to specify if images already available should be reprocessed if a new run of MODIStsp is launched with the same output folder. If set to “No”, MODIStsp skips dates for which output files following the MODIStsp naming conventions are already present in the output folder. This allows to incrementally extend MODIS time series without reprocessing already available dates.

Output Folder for Original HDF Storage

Allows to specify the folder where downloaded original MODIS HDF files are stored. The “delete original HDF files” checkbox allows to specify if the downloaded images should be deleted from the file system at the end of the processing. To avoid accidental file deletion, this is always set to “No” by default, and a warning is issued before execution whenever the selection is changed to “Yes”.

Processing

Upon pressing the “Start” button, the main processing routine is launched (MODIStsp_process). This routine performs the following main tasks:

  1. Retrieve the processing options from the GUI (or the saved RData file in case of non-interactive execution

  2. Connect to the lpdaac http MODIS distribution archive, and retrieve the list of HDF images of the selected MODIS product available for the time period selected by the user, for each tile required to cover the selected study area.

  3. For each identified date of acquisition;

    1. Download all required hdf images.

    2. For each original hdf layer i) selected by the user, or ii) required to compute a selected QI or SI layer, extract the data from the original MODIS images and resize and reproject it to the selected output projection, extent and resolution. If more than one tile is needed to cover the output extent, all required tiles are automatically mosaicked before resizing (gdalbuildvrt functionalities are used to avoid creating large temporary raster files). All the main spatial processing tasks are performed using standard GDAL routines, exploiting the wrappers for “R” provided in the “gdalUtils” package (Greenberg and Mattiuzzi 2014). Results are saved as raster GeoTiff or ENVI files using MODIStsp naming conventions (See Section 5).

    3. Starting from files created at point b), compute required QI and/or SI layers and save the results as GeoTiff or ENVI files. Quality Indicators are computed from original QA layers using bitwise operators available in the “bitOps” package (Dutky and Maechler 2013), using a generalization of the “modis.qc.R” script by Yann Chemin (Chemin 2008). Computation of SI layers exploits on-the-fly parsing of the indexes’ formulae to identify the required input raster files and perform the computation.

    4. Delete raster files created at point b) used to compute the QI or SI layers but that correspond to original HDF layers not originally selected for processing.

  4. When all dates have been processed, create the virtual time series files if required.

Output format and naming conventions

Output raster files are saved in specific subfolders of the main output folder. A separate subfolder is created for each processed original MODIS layer, Quality Indicator or Spectral Index. Each subfolder contains one image for each processed date, created according to the following naming conventions:

“Product_Code”_“Layer”_“YYYY”_“DOY”.“ext” (e.g.,MOD13Q1_NDVI_2000_065.dat)

Product_Code is the code name of the MODIS product from which the image was derived (e.g., MOD13Q1), Layer is a short name describing the dataset (e.g., b1_Red, NDII, UI), YYYY and DOY corresponds to the year and DOY (Day of the Year) of acquisition of the original MODIS image, and ext is the file extension (.tif for GTiff outputs, or .dat for ENVI outputs).

If requested, ENVI and/or GDAL virtual files are stored in the “Time_Series” subfolder. Naming convention for virtual files is:

“Product_Code”_“Layer”_“StartYYYY”_“StartDOY”_“EndYYYY”_“EndDOY”.“ext”

(e.g., MOD13Q1_NDVI_49_2000_17_2015.dat)

StartYYYY, StartDOY, EndYYYY and EndDOY indicate the starting and ending years and DOYS of the time series.

Standalone execution and scheduled processing

MODIStsp can be executed as a standalone application using the MODIStsp.bat (for Windows) or MODIStsp.sh (for Linux) batch execution scripts available in the “MODIStsp/ExtData/launcher” subfolder of the package installation. Double-clicking the files or launching them from a shell without parameters launches MODIStsp in interactive mode.

Non-interactive mode is triggered by adding the “-g” argument to the call, and specifying the path to a valid Options File as “-s” argument

(launch yourpath_to_MODIStsp_sh/MODIStsp.sh -h for details).

(launch yourpath_to_MODIStsp_bat\MODIStsp.bat -h for details).

Standalone non-interactive execution easily allows to automatically update the time series of a selected product over a given study area whenever a new MODIS image is available. To do that, the user must simply:

  1. Open the MODIStsp GUI, define the parameters of the processing specifying a date in the future as the “Ending Date” and save the processing options. Then quit the program

  2. Schedule non-interactive execution of MODIStsp.bat (or MODIStsp.sh) as windows scheduled task (or linux “cron” job) according to a specified time schedule, specifying the path of a previously saved Options file as additional argument:

    • In Linux: edit your crontab by opening a terminal and typing

      crontab -e

      Then add an entry for the MODIStsp.bsh For example, if you want to run the tool every day at 23.00, add the following row:

      0 23 * * * /bin/bash /yourpath_to_MODIStsp_sh/MODIStsp.sh -g -s "/yourpath/youroptions.RData"
    • In Windows: create a Task following these instructions; add the path of the MODIStsp.bat launcher as Action (point 6), and specify -g -s "X:/yourpath/youroptions.RData" as argument.


  1. The lists of original MODIS layers, QIs and Sis available for the selected product are automatically retrieved from the “MODIStsp_Products_Opts” XML file distributed with the package in /ExtData subfolder.