Site index plotting code

[Dave Larsen]

Hello everyone -

I know most of the people on the list. Aaron invited me to join. I
have quite abit of R/Splus code for forestry problems that I have
accumulated over the last 20 years.

I thought I would offer abit of code as a way to say Hello. Here is
one that will help Eastern US people more that Western US people.

Based on:

Carmean, W. H., J. T. Hahn and R. D. Jacobs. 1989. Site index curves
for forest tree species in the eastern United States. NCEFS GTR
NC-128.

http://www.treesearch.fs.fed.us/pubs/10192

Here is R code that plots site index any of the 110 site index curves
in the publication.


----------------------------------------------------------------------------------------
"plot.siteidx"<-function( param=mo.oak )
{
#
# Routine to plot a standard site index functions
# As presented in "Site index curves for forest tree species in
# the eastern United States" by W. Carmean, J. t. Hahn and R. D.
Jacobs
# NCFES General Technical Report NC-128
#
# by David Larsen, Copyright May 12, 1995
#

oldpar <- par()
par(pty = "s", tck = 1, lab = c(16, 16, 18),lty=2)
x <- seq(param$age[1],param$age[2])
yy <- numeric()
y <- vector()
s <- seq(param$site[1],param$site[2],5)
for( i in s ){
tmp <- param$bht + (param$b[1]*(i^param$b[2])*
((1-exp(param$b[3]*x))^(param$b[4]*(i^param$b[5]))))
y_cbind(y,tmp)
yy_c(yy,tmp[length(x)])
}
matplot(x=x,y=y,type="l",lty=c(1,3),col=1,xlim=c(0,(param
$age[2]+10)),
xlab="Age (Years)", ylab = "Height (ft)", main=param$title )
par(lty=1)
box()
tmp <- param$bht + (param$b[1]*(param$site[2]^param$b[2])*
((1-exp(param$b[3]*(param$age[2]-2)))^(param$b[4]*
(param$site[2]^param$b[5]))))
text( param$age[1], tmp, labels=param$sub )
text( param$age[2]+3, yy, labels=as.character(s) )
par(oldpar)
invisible()
}
---------------------------------------------------------------------------------------------------------

Here is the input data structure required by the R-code

---------------------------------------------------------------------------------------------------------

"mo.oak" <-
list("title" = "Black, Scarlet, and White Oaks"
, "sub" = "McQuilkin, 1978"
, "bht" = 0.
, "b" = c(4.9676, 0.7459, -0.0154, 7.364, -0.5144)
, "site" = c(40., 80.)
, "age" = c(10., 80.)
)

"nash.pine" <-
list("title" = "Shortleaf Pine"
, "sub" = "Nash, 1963"
, "bht" = 0.
, "b" = c(0.8977, 1.0624, -0.0398, 0.1419, 0.4709)
, "site" = c(30., 85.)
, "age" = c(10., 70.)
)

"red.cedar" <-
list("title" = "Eastern Red Cedar"
, "sub" = "Hampf, 1965"
, "bht" = 0.
, "b" = c(0.9276, 1.0591, -0.0424, 0.3529, 0.3114)
, "site" = c(20., 50.)
, "age" = c(10., 90.)
)

"schnur" <-
list("title" = "Upland Oaks"
, "sub" = "Schnur, 1937"
, "bht" = 0.
, "b" = c(2.1037, 0.914, -0.0275, 3.7962, -0.253)
, "site" = c(40., 80.)
, "age" = c(10., 100.)
)
----------------------------------------------------------------------------------------------

I use this code to create Site Index graphs to hand out to students.


Enjoy

Dave Larsen