Description Usage Arguments Details Value Note Author(s) References See Also Examples
Interpolate Raster* times series to new time stamps
1 2 3 4 5 |
x |
A |
w |
Optional, a |
t |
Optional, a |
timeInfo |
object returned by |
method |
character |
lambda |
Yearly lambda value. Default is 5000. See |
nIter |
numeric |
outlierThreshold |
Numeric. If provided, |
df |
Yearly degree of freedom for 'spline' |
cl |
cluster object for parallel processing (do not use it if beginCluster() has been called. Default is |
filename |
name of the file where the interpolated raster is saved. Default is created through |
... |
arguments passed to |
'linear' and 'constant' method
use the approx
interpolation functions while 'spline' use
splinefun
.
A Raster-class
object with values interpolated to the output time series (as defined in timeInfo
)
Most of the code is borrowed from smooth.spline.raster
and whittaker.raster
.
Compared to the MODIS
version, code is simplified and made general for reflectance values
Antoine Stevens & Matteo Mattiuzzi
Atzberger, C., and Eilers, P.H.C. (2011). A time series for monitoring vegetation activity and phenology at 10-daily time steps covering large parts of South America. International Journal of Digital Earth 4, 365-386.
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 | ## Not run:
# Get and extract NDVI MOD13A2 for 2009
runGdal(product="MOD13A2", begin = "2009.01.01", end = "2009.12.31",
tileH = 19, tileV = 5 ,
SDSstring = "100000000000",job="H19V5",
gdalPath = "C:/OSGeo4W64/bin")
# process and convert to raster
path <- "~/MODIS_DATA/PROCESSED/H19V5"
r <- convert_modis(path = path,pattern = "16_days_NDVI",
begin = "2009.01.01", end = "2009.12.31",extractAll=T)
str(r)
r_NDVI <- r$raster_NDVI
# Smooth by splines
r_smooth <- interpolate_raster(x = r_NDVI, method = "spline")
# plot
x <- r$ti_NDVI$inputLayerDates
y1 <- r_NDVI[1000]
y2 <- r_smooth[1000]
plot(x,y1)
lines(x,y2,col = "red")
# now, using weights and composite day of the year to smooth
# first extract data
runGdal(product="MOD13A2", begin = "2009.01.01", end = "2009.12.31",
tileH = 19, tileV = 5 ,
SDSstring = "001000000010", # extract VI_usefulness and composite_day_of_the_year
job="H19V5",
gdalPath = "C:/OSGeo4W64/bin")
cdoy <- convert_modis(path = path ,pattern = "16_days_composite",
begin = "2009.01.01", end = "2009.12.31")
QF <- convert_modis(path = path ,pattern = "16_days_VI_Quality",
begin = "2009.01.01", end = "2009.12.31", convertDN = FALSE)
# Convert QF to 'usefulness'
usefulness <- convert_qf_modis(r = QF,qf = "VI_usefulness")
# Convert usefulness to weights (see reference section)
w <- (15-usefulness)/15
r_smooth2 <- interpolate_raster(x = r_NDVI, w = w, t = cdoy, method = "spline")
y3 <- r_smooth2[1000]
lines(x,y3,col="green")
# One could aslo interpolate to new time stamps
ti <- orgTime(r_NDVI,nDays = 10,
begin = "2009.01.01", end = "2009.12.31") # resample every 10 days
r_smooth_resampled <- interpolate_raster(x = r$raster_NDVI, w = w, t = cdoy,
timeInfo = ti, method = "spline")
y4 <- r_smooth_resampled[1000]
points(ti$outputLayerDates,y4,pch=3)
## End(Not run)
|
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.