Nothing
sitree <- function(tree.df,
stand.df,
functions,
n.periods = 5,
period.length,
mng.options = NA,
print.comments = FALSE,
... ,
ext.mod.first = TRUE){
this.is.the.call <- match.call()
others <- list(...)
## check that functions is a list
if (!is.list(functions)) stop ('functions should be a list')
if (!is.data.frame(tree.df)) stop ('tree.df should be a data.frame')
if ( !(is.data.frame(stand.df) | is.list(stand.df) ) ) {
stop ('stand.df should be a data.frame or a list')
}
## check that all functions are defined
if (!all(c('fn.growth', 'fn.mort', 'fn.recr', 'fn.management', 'fn.modif',
'fn.tree.removal', 'fn.prep.common.vars') %in%
names(functions))){
stop ('There are some missing definitions for the functions. Please read sitree man')
}
## check that tree.df and stand.df have all neccessary information
## No missing dbh or heights.
names.tr <- c('plot.id', 'treeid', 'dbh', 'height', 'tree.sp')
if (!all(names.tr %in%
names(tree.df))){
stop ('There are some missing or incorrecly named columns in tree.df')
}
## tree list should contain at least dbh, height, treeid, plot.id, and species
if (any(!is.finite(tree.df$dbh))) stop('Some dbh are NA, please correct it')
if (any(!is.finite(tree.df$height))) stop('Some heights are NA, please correct it')
if (any(is.na(tree.df$treeid))) stop('Some treeid are missing, please correct it')
if (any(is.na(tree.df$tree.sp))) stop('Some tree.sp are missing, please correct it')
if (any(is.na(tree.df$plot.id))) stop('Some plot.id in tree.df are missing')
## stand
if (!all(c('plot.id') %in% names(stand.df))) stop ('plot.id is missing from stand.df')
if (any(is.na(stand.df$plot.id))) stop('Some plot.id in stand.df are missing')
## BUILD THE OBJECTS tr and fl
### tree data
foo <- function(mvec, n.periods) {
mdt <- matrix(0, nrow = length(mvec), ncol = (n.periods +1))
mdt[,1] <- mvec
colnames(mdt) <- paste("t", 0:n.periods, sep = "")
return(as.data.frame(mdt))
}
## make sure plot.id in tr are defined in stands
if (any(!tree.df$plot.id %in% stand.df$plot.id)) {
stop ('plot.id missing in stand.df')
}
tr.rest <- as.list(tree.df[, !names(tree.df) %in% names.tr])
## tree species
## if already a factor keep it like that, in case it has extra levels
## that might be needed after
if (is.factor(tree.df$tree.sp)) {
tree.sp <- tree.df$tree.sp} else tree.sp <- factor(tree.df$tree.sp)
trl <- list(
plot.id = tree.df$plot.id,
treeid = tree.df$treeid,
dbh.mm = foo(as.integer(round(tree.df$dbh)), n.periods = n.periods ),
height.dm = foo(as.integer(tree.df$height), n.periods = n.periods ),
yrs.sim = rep(0, nrow(tree.df)),
tree.sp = tree.sp
)
rm(tree.sp)
trl <- c(trl, tr.rest)
tr <- trList$new(data = trl, nperiods = as.integer(n.periods))
## clean up
rm(trl, tr.rest)
### stand data, a data.frame or list
fl <- as.list(stand.df)
## create the management data frame if it does not exist
if (!"management" %in% names(fl)){
fl$management <- data.frame(matrix(NA, ncol = n.periods,
nrow = length(fl$plot.id)))
names(fl$management) <- paste0("t", 1:n.periods)
}
######################################
## VARS REQUIRED
######################################
## Calculate the variables required by the functions
## This is meant to avoid recalculating common variables for each function
vars.required <- fn.vars.required(my.functions = functions, ...)
management <- list(management = rep(NA, length(fl$plot.id)))
## a loop over all periods
for(i.period in (0: (tr$nperiods - 1))){
if (print.comments) print(paste0("period: ", i.period))
this.period <- paste0("t", i.period )
next.period <- paste0("t", i.period+1)
## PREP shared variables
if (i.period == 0) {
existing.common.vars <- "NULL"
} else {
existing.common.vars <- common.vars
}
if (!ext.mod.first){
prep.common.vars <- do.call(functions$fn.prep.common.vars,
args = list(
tr = tr,
fl = fl,
this.period = this.period,
i.period = i.period,
common.vars = existing.common.vars,
mng.options = mng.options,
vars.required = vars.required,
period.length = period.length,
n.periods = n.periods,
print.comments = print.comments,
last.management = management,
...
)
)
common.vars <- prep.common.vars$res
fl <- prep.common.vars$fl
if (print.comments) print('Passed common vars')
## No hgt.inc but hgt function
if (is.null(functions$fn.hgt.inc) & !is.null(functions$fn.hgt)){
tr$data$height.dm[, this.period] <-
do.call(functions$fn.hgt,
args = list(
tr = tr,
fl = fl,
this.period = this.period,
i.period = i.period,
common.vars = common.vars,
vars.required = vars.required,
period.length = period.length,
print.comments = print.comments,
...
)
)
}
if (print.comments) print('Passed heights')
}
## EXTERNAL MODIFIERS
## External modifiers other than management, e.g. climate change
## external modifieres should be a list with the elements
## to be replaced in fl
if (!is.null(functions$fn.modif) ){
ext.modif <- do.call(functions$fn.modif,
args = list(
tr = tr,
fl = fl,
this.period = this.period,
i.period = i.period,
common.vars = existing.common.vars,
vars.required = vars.required,
period.length = period.length,
print.comments = print.comments,
...
)
)
if (!is.null(ext.modif)) fl[names(ext.modif)] <- ext.modif
if (print.comments) print('Passed Exter mod')
}
## apply external modifiers, they might change some variables
if (ext.mod.first ) {
## recalculate prep.common.vars, it probably change things
prep.common.vars <- do.call(functions$fn.prep.common.vars,
args = list(
tr = tr,
fl = fl,
this.period = this.period,
i.period = i.period,
common.vars = existing.common.vars,
mng.options = mng.options,
vars.required = vars.required,
period.length = period.length,
n.periods = n.periods,
print.comments = print.comments,
last.management = management,
...
)
)
common.vars <- prep.common.vars$res
fl <- prep.common.vars$fl
if (print.comments) print('Passed Prep common variables')
## No hgt.inc but hgt function
if (is.null(functions$fn.hgt.inc) & !is.null(functions$fn.hgt)){
tr$data$height.dm[, this.period] <-
do.call(functions$fn.hgt,
args = list(
tr = tr,
fl = fl,
this.period = this.period,
i.period = i.period,
common.vars = common.vars,
vars.required = vars.required,
period.length = period.length,
print.comments = print.comments,
...
)
)
}
if (print.comments) print('Passed heights')
}
## Management - return the management of the different plots
if ( !is.null(functions$fn.management)){
management <- do.call(functions$fn.management,#.scen,
args = list(
tr = tr,
fl = fl,
common.vars = common.vars,
this.period = this.period,
next.period = next.period,
i.period = i.period,
mng.options = mng.options,
print.comments = print.comments,
...)
)
if (print.comments) print('Passed management')
## save that certain stands have been managed in the management table,
## not sure it has to be done here, but if possible we should move it down
## this was done before inside the fn.management function--- does not work
## when working on parallel
## CHECK IF POSSIBLE
fl$management[, next.period] <- management$management
removed <- do.call(
functions$fn.tree.removal,
args = list(
tr = tr,
fl = fl,
common.vars = common.vars,
management = management,
this.period = this.period,
next.period = next.period,
i.period = i.period,
mng.options = mng.options,
print.comments = print.comments,
...)
)
if (print.comments) print(paste0('sum removed ', sum(removed)))
if (print.comments) print('Passed removed')
} else {
## we don't remove trees if there is no management
management <- list(management = rep(NA, length(fl$plot.id)))
fl$management[, next.period] <- management$management
removed <- rep(FALSE, nrow(tr$data$dbh.mm))
}
## GROWTH
growth <- do.call(
functions$fn.growth,
args = list(
tr = tr,
fl = fl,
common.vars = common.vars,
this.period = this.period,
next.period = next.period,
i.period = i.period,
functions = functions,
print.comments = print.comments,
...)
)
if (print.comments) print('Passed growth')
## DEATH - return a FALSE and TRUE(dead) vector
mort <- do.call(
functions$fn.mort,
args = list(
tr = tr,
fl = fl,
common.vars = common.vars,
this.period = this.period,
next.period = next.period,
i.period = i.period,
functions = functions,
print.comments = print.comments,
...)
)
if (print.comments) print('Passed mortality')
## INGROWTH -natural and artificial on all plots
## returns a data.frame with the new trees
max.tree.id <-
ifelse (exists('dead.trees') & exists('removed.trees'),
max(c(tr$data$treeid, dead.trees$data$treeid,
removed.trees$data$treeid)),
ifelse(exists('dead.trees'),
max(c(tr$data$treeid, dead.trees$data$treeid)),
ifelse(exists('removed.trees'),
max(c(tr$data$treeid, removed.trees$data$treeid) ),
max(tr$data$treeid)
)))
ingrowth <- do.call(
functions$fn.recr,
args = list(
tr = tr,
fl = fl,
common.vars = common.vars,
i.period = i.period,
this.period = this.period,
next.period = next.period,
management = management,
period.length = period.length,
vars.required = vars.required,
mng.options = mng.options,
print.comments = print.comments,
max.tree.id = max.tree.id,
...)
)
if (print.comments) print('Passed ingrowth')
## book keeping
i.removed <- tr$data$treeid[removed]
## YEARS IN SIMULATION
tr$data$yrs.sim <- period.length + tr$data$yrs.sim
## GROWTH
## grow dbh
tr$data$dbh.mm[ next.period] <-
tr$data$dbh.mm[, this.period] + as.integer(round(growth[, "dbh.inc.mm"]))
##i.small.dbh <- tr$data$dbh.mm[ next.period] < 50
##tr$data$dbh.mm[ next.period] [i.small.dbh] <- 50
## grow heigth
tr$data$height.dm[, next.period] <-
tr$data$height.dm[, this.period] + as.integer(growth[, "hgt.inc.dm"])
if (print.comments) print("Growth applied")
## DEAD TREES
i.dead.trees <- mort & !removed ## not sure if this is neccessary
## extract
if (sum(i.dead.trees) > 0) {
new.dead.trees <- tr$extractTrees(which(i.dead.trees))
## create dead trees class object
new.dead.trees <- trListDead$new(
data = new.dead.trees,
last.measurement = cbind(
do.call("dead.trees.growth"
, args = list(
dt = new.dead.trees,
growth = growth,
mort = i.dead.trees,
this.period = this.period)
),
found.dead = next.period
),
nperiods = tr$nperiods
)
## remove last measurement --- we could probably skip this step if
## we change the code slightly, if only alived trees are grown
new.dead.trees$remove.next.period(next.period = next.period)
## If this is the first period create the dead.trees object, if not
## add the new dead trees to the dead.trees object
if (i.period == 0 | !exists('dead.trees')){
dead.trees <- new.dead.trees
} else {
dead.trees$addTrees(new.dead.trees)
}
}
if (print.comments) print("Dead trees applied")
## REMOVED TREES
## extract
if (sum(removed) > 0){
new.removed.trees <- tr$extractTrees(
which(tr$data$treeid %in% i.removed))
## create removed trees class object
new.removed.trees <- trListDead$new(
data = new.removed.trees,
last.measurement = cbind(
do.call("dead.trees.growth"
, args = list(
dt = new.removed.trees,
growth = growth,
mort = removed,
this.period = this.period)
),
found.removed = next.period
),
nperiods = tr$nperiods
)
## remove last measurement --- we could probably skip this step if
## we change the code slightly, if only alived trees are grown
new.removed.trees$remove.next.period(next.period = next.period)
## If this is the first period create the removed.trees object, if not
## add the new removed trees to the removed.trees object
if (!exists("removed.trees")){
removed.trees <- new.removed.trees
} else {
removed.trees$addTrees(new.removed.trees)
}
}
if (print.comments) print("removed trees applied")
## Apply ingrowth
if (length(ingrowth$treeid) > 0) tr$addTrees(ingrowth)
## remove objects
rm(growth, mort, new.dead.trees, removed)
## End of the period-loop
}
if (!exists('removed.trees')) removed.trees <- NULL
if (!exists('dead.trees')) dead.trees <- NULL
to.return <- list(live = recover.state(tr = tr,
dead.trees = dead.trees,
removed.trees = removed.trees
),
dead = dead.trees,
removed = removed.trees,
plot.data = fl)
attr(to.return, 'Call') <- this.is.the.call
attr(to.return, 'par') <- others
invisible( to.return )
}
Any scripts or data that you put into this service are public.
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.