[R] vectorization of rolling function

[A Duranel]

Hello
I use R to run a simple model of rainfall interception by vegetation:
rainfall falls on vegetation, some is retained by the vegetation (part of
which can evaporate), the rest falls on the ground (quite crude but very
similar to those used in SWAT or MikeSHE, for the hydrologists among you).
It uses a loop on zoo time-series of rainfall and potential
evapotranspiration. Unfortunately I did not find a way to vectorize it and
it takes ages to run on long datasets. Could anybody help me to make it run
faster?

library(zoo)
set.seed(1)
# artificial potential evapotranspiration time-series
ETmax<-zoo(runif(10, min=1, max=6), c(1:10))
# artificial rainfall time-series
RR<-zoo(runif(10, min=0, max=6), c(1:10))

## create empty time-series to fill
# effective rainfall (i.e. rainfall minus intercepted rainfall)
RReff<-zoo(NA, c(1:10))
# intercepted rainfall
int<-zoo(NA, c(1:10))
# residual potential evapotranspiration (ie ETmax minus evaporation from
interception)
ETres<-zoo(NA, c(1:10))

# define maximum interception storage capacity (maximum volume of rainfall
that can be intercepted per time step, provided the interception store is
empty at start of time-step)
Smax<-3
# volume of water in interception storage at start of computation
storage<-0

for (i in 1:length(ETmax)) {
  # compute interception capacity for time step i (maximum interception
capacity minus any water intercepted but not evaporated during previous
time-step).
  int[i]<-Smax-storage
  # compute intercepted rainfall: equal to rainfall if smaller than
interception capacity, and to interception capacity if larger.
  if(RR[i]<int[i]) int[i]<-RR[i]
  # compute effective rainfall (rainfall minus intercepted rainfall).
  RReff[i]<-RR[i]-int[i]
  # update interception storage: initial interception storage + intercepted
rainfall.
  storage<-storage+coredata(int[i])
  # compute evaporation from interception storage: equal to potential
evapotranspiration if the latter is smaller than interception storage, and
to interception storage if larger. 
  if(storage>coredata(ETmax[i])) evap<-coredata(ETmax[i]) else evap<-storage
  # compute residual potentiel evapotranspiration: potential
evapotranspiration minus evaporation from interception storage.
  ETres[i]<-ETmax[i]-evap
  # update interception storage, to be carried over to next time-step:
interception storage minus evaporation from interception storage.
  storage<-storage-evap
}

Many thanks for your help!

Arnaud
UCL Department of Geography, UK



--
View this message in context: http://r.789695.n4.nabble.com/vectorization-of-rolling-function-tp4700487.html
Sent from the R help mailing list archive at Nabble.com.