Preprocessed MODIS data can be retrieved within R
either by accessing the single-date raster files, or by loading the saved RasterStack objects.
Any single-date image can be accessed by simply opening it with a raster
command:
library(raster)
modistsp_file <- "/my_path/my_moditsp_folder/MOD13Q1_2005_137_EVI.tif"
my_raster <- raster(modistsp_file)
rasterStack
time series containing all the processed data for a given parameter (saved in the “Time Series” subfolder - see here for details) can be opened by:
in_virtual_file <- "/my_path/my_moditsp_folder/Time_Series/MOD13Q1_MYD13Q1_NDVI_49_2000_353_2015_RData.RData"
indata <- get(load(in_virtual_file))
This second option allows accessing the complete data stack and analyzing it using the functionalities for raster/raster time series analysis, extraction and plotting provided for example by the raster
or rasterVis
packages.
MODIStsp
provides an efficient function (MODIStsp\_extract
) for extracting time series data at specific locations. The function takes as input a RasterStack virtual object created by MODIStsp
(see above), the starting and ending dates for the extraction and a standard _Sp*_ object (or an ESRI shapefile name) specifying the locations (points, lines or polygons) of interest, and provides as output a xts
object or data.frame
containing time series data for those locations.
If the input is of class SpatialPoints, the output object contains one column for each point specified, and one row for each date. If it is of class SpatialPolygons (or SpatialLines), it contains one column for each polygon (or each line), with values obtained applying the function specified as the “FUN” argument (e.g., mean, standard deviation, etc.) on pixels belonging to the polygon (or touched by the line), and one row for each date.
As an example the following code:
#Set the input paths to raster and shape file
infile <- 'in_path/MOD13Q1_MYD13Q1_NDVI_49_2000_353_2015_RData.RData'
shpname <- 'path_to_file/rois.shp'
#Set the start/end dates for extraction
startdate <- as.Date("2010-01-01")
enddate <- as.Date("2014-12-31")
#Load the RasterStack
inrts <- get(load(infile))
# Compute average and St.dev
dataavg <- MODIStsp_extract(inrts, shpname, startdate, enddate, FUN = 'mean', na.rm = T)
datasd <- MODIStsp_extract (inrts, shpname, startdate, enddate, FUN = 'sd', na.rm = T)
# Plot average time series for the polygons
plot.xts(dataavg)
loads a RasterStack object containing 8-days 250 m resolution time series for the 2000-2015 period and extracts time series of average and standard deviation values over the different polygons of a user’s selected shapefile on the 2010-2014 period.