## BioCro/R/weach.R by Fernando Ezequiel Miguez Copyright (C) 2007-2009
##
## This program is free software; you can redistribute it and/or modify it
## under the terms of the GNU General Public License as published by the Free
## Software Foundation; either version 2 or 3 of the License (at your option).
##
## This program is distributed in the hope that it will be useful, but WITHOUT
## ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
## FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for
## more details.
##
## A copy of the GNU General Public License is available at
## http://www.r-project.org/Licenses/
weach366 <- function(X, lati, ts = 1, temp.units = c("Farenheit", "Celsius"), rh.units = c("percent",
"fraction"), ws.units = c("mph", "mps"), pp.units = c("in", "mm"), ...) {
if (missing(lati))
stop("latitude is missing")
if ((ts < 1) || (24%%ts != 0))
stop("ts should be a divisor of 24 (e.g. 1,2,3,4,6,etc.)")
if (dim(X)[2] != 11)
stop("X should have 11 columns")
MPHTOMPERSEC <- 0.447222222222222
temp.units <- match.arg(temp.units)
rh.units <- match.arg(rh.units)
ws.units <- match.arg(ws.units)
pp.units <- match.arg(pp.units)
year <- X[, 1]
DOYm <- X[, 2]
solar <- X[, 3]
maxTemp <- X[, 4]
minTemp <- X[, 5]
avgTemp <- X[, 6]
maxRH <- X[, 7]
minRH <- X[, 8]
avgRH <- X[, 9]
WindSpeed <- X[, 10]
precip <- X[, 11]
tint <- 24/ts
tseq <- seq(0, 23, ts)
## Solar radiation
solarR <- (0.12 * solar) * 2.07 * 10^6/3600
solarR <- rep(solarR, each = tint)
ltseq <- length(tseq)
resC2 <- numeric(ltseq * 366)
for (i in 1:366) {
res <- lightME(DOY = i, t.d = tseq, lat = lati, ...)
Itot <- res$I.dir + res$I.diff
indx <- 1:ltseq + (i - 1) * ltseq
resC2[indx] <- (Itot - min(Itot))/max(Itot)
}
SolarR <- solarR * resC2
## Temperature
if (temp.units == "Farenheit") {
minTemp <- (minTemp - 32) * (5/9)
minTemp <- rep(minTemp, each = tint)
maxTemp <- (maxTemp - 32) * (5/9)
maxTemp <- rep(maxTemp, each = tint)
rangeTemp <- maxTemp - minTemp
# avgTemp <- rep((avgTemp - 32)*(5/9),each=tint)
} else {
minTemp <- rep(minTemp, each = tint)
maxTemp <- rep(maxTemp, each = tint)
rangeTemp <- maxTemp - minTemp
}
xx <- rep(tseq, 366)
temp1 <- sin(2 * pi * (xx - 10)/tint)
temp1 <- (temp1 + 1)/2
Temp <- minTemp + temp1 * rangeTemp
## Relative humidity avgRH <- rep(avgRH,each=tint)
minRH <- rep(minRH, each = tint)
maxRH <- rep(maxRH, each = tint)
temp2 <- cos(2 * pi * (xx - 10)/tint)
temp2 <- (temp2 + 1)/2
if (rh.units == "percent") {
RH <- (minRH + temp2 * (maxRH - minRH))/100
} else {
RH <- (minRH + temp2 * (maxRH - minRH))
}
## Wind Speed
if (ws.units == "mph") {
WS <- rep(WindSpeed, each = tint) * MPHTOMPERSEC
} else {
WS <- rep(WindSpeed, each = tint)
}
## Precipitation
if (pp.units == "in") {
precip <- rep(I((precip * 2.54 * 10)/tint), each = tint)
} else {
precip <- rep(I(precip/tint), each = tint)
}
hour <- rep(tseq, 366)
DOY <- 1:366
doy <- rep(DOY, each = tint)
ans <- cbind(year, doy, hour, SolarR, Temp, RH, WS, precip)
ans
}
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.