Learn R Programming

MODIS (version 1.1.2)

repDoy: Repair MODIS "composite_day_of_the_year" SDS

Description

Currently works only for MODIS 16 days composites! In MODIS composites, the Julian dates inside the 'composite_day_of_the_year' SDS are referring always to the year they are effectively in. The problem is that the layer/SDS name from the last files from Terra and Aqua within a year can include dates from the following year and so starting again with 1. The problem comes if you want to sort values of a time series by date (e.g. for precise time series functions). This function generates a sequential vector beginning always with the earielst SDS/layer date and ending with the total sum of days of the time serie length.

Usage

repDoy(pixX, layerDate = NULL, bias = 0)

Arguments

pixX

matrix of values, usually derived from as.matrix.

layerDate

If NULL (default), try to autodetect layer dates. If you want to be sure, use the result from extractDate(..., asDate = TRUE) or orgTime.

bias

integer. Bias applied to all values in pixX.

Value

A matrix with sequential Julian dates.

Examples

Run this code
# NOT RUN {
runGdal(product="M.D13A2", begin="2010350", end="2011016", extent="Luxembourg",
        job="deleteme", SDSstring="100000000010")

ndviFiles <- preStack(path=paste(MODIS:::.getDef()$outDirPath,"deleteme",sep="/")
                      ,pattern="*_NDVI.tif$")

ndvi <- stack(ndviFiles)

doyFiles <- preStack(path=paste(MODIS:::.getDef()$outDirPath,"deleteme",sep="/")
                     ,pattern="*_composite_day_of_the_year.tif$")

doy <- stack(doyFiles)
layerDates <- extractDate(names(doy))

pixX <- 169

y <- ndvi[pixX]
x1 <- doy[pixX]
x2 <- repDoy(doy[pixX],layerDates)

x1
x2
# the plotting example is not really good. 
# To create a figurative example it would be necessary to dolwnload to much data! 
plot("",xlim=c(1,max(x1,x2)),ylim=c(0,2000),xlab="time",ylab="NDVI*10000")
lines(y=y,x=x1,col="red",lwd=3)
lines(y=y,x=x2,col="green",lwd=2)

# repDoy function is thought to be enbeded in something like that:
tr <- blockSize(ndvi)

doyOk <- brick(doy)
doyOk <- writeStart(doyOk, filename='test.tif',  overwrite=TRUE)

for (i in 1:tr$n)
{
  pixX  <- getValues(doy,tr$row[i],tr$nrows[i])
  ok    <- repDoy(pixX,layerDates)
  doyOk <- writeValues(x=doyOk,v=ok,start=tr$row[i])
}
doyOk <- writeStop(doyOk)

# unlink(filename(doyOk))
# }
# NOT RUN {
 
# }

Run the code above in your browser using DataLab