Hello everyone,
I am working for a few days already on a basic algorithm, very common in
applied agronomy, that aims to determine the degree-days necessary for a
given individual to reach a given growth stade. The algorithm (and context)
is explained here: http://www.oardc.ohio-state.edu/gdd/glossary.htm , and
so I implemented my function in R as follows:
DD <- function(Tmin, Tmax, Tseuil, meanT, method = "DDsin")
### function that calculates the degree-days based on
### minimum and maximum recorded temperatures and the
### minimal threshold temperature (lower growth temperature)
{
### method arcsin
if(method == "DDsin"){
cond1 <- (Tmax <= Tseuil)
cond2 <- (Tmin >= Tseuil)
amp <- ((Tmax - Tmin) / 2)
print((Tseuil-meanT)/amp)
alpha <- asin((Tseuil - meanT) / amp)
DD_ifelse3 <- ((1 / pi) * ((meanT - Tseuil) * ((pi/2) - alpha))
+
amp*cos(alpha))
DD <- ifelse(cond1, 0, ifelse(cond2, (meanT - Tseuil),
DD_ifelse3))
}
### method (Tmin + Tmax) / 2
else if(method == "DDt2"){
cond1 <- (meanT > Tseuil)
DD <- ifelse(cond1,(meanT - Tseuil),0)
}
else{
stop("\nMethod name is invalid.\nMethods available = DDsin
(sinus) or DDt2
(mean)\n")
}
return(DD)
}
BUT! When I try to process random data:
library(reshape2)
library(plyr)
station <- rep(c("station1","station2","station3"), 20)
values_min <- sample(-5:20, size = 60, replace = T)
values_max <- sample(20:40, size = 60, replace = T)
meanT <- ((values_min+values_max)/2)
d <- data.frame(station,values_min,values_max,meanT)
names(d) <- c("station", "values_min","values_max","meanT")
x<-ddply(d, .(station), transform, t1 =
cumsum(DD(values_min,values_max,0,meanT)))
I get a warning on my alpha calculation (NaN produced); indeed, the values I
give as argument to asin() are out of the range [-1:1], as the print()
reveals. I can't figure out how to solve this issue, because the same
algorithm works in Excel (visual basic).
It is very annoying, especially because it seems that no occurence of such
error using that algorithm can be found on Internet.
Any help is welcome :) Thanks for your time
P.
--
View this message in context:
http://r.789695.n4.nabble.com/Issue-with-asin-tp4484462p4484462.html
Sent from the R help mailing list archive at Nabble.com.
______________________________________________
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
