On Feb 29, 2012, at 2:14 PM, lidaky wrote:

index.refraction <- function(Temp,Press, RH, CO2)
{
 #define constant
 lambda  <-  531 #nm
 Co  <- 299792458 #m/s
 ww <- c(295.235,2.6422,-0.03238,0.004028)
 kk  <- c( 238.0185 ,5792105 ,57.362 ,167917)
 aa  <- c( 1.58123 * 10^(-6) , -2.9331 * 10^(-8) , 1.1043 * 10^(-10))
 bb  <- c( 5.707 * 10^(-6) , -2.051 * 10^(-8))
 cc  <-  c(1.9898 * 10^(-4) , -2.376 * 10^(-6))
 dd  <- 1.83 * 10^(-11)
 ee  <-  -0.765 * 10^(-8)
 pp <- 101325
 tt <- 288.15
 zz <- 0.9995922115
 ro <- 0.00985938
 R <- 8.314472
 M <- 0.018015

 S <- 1/(lambda*10^-3)^2
 ras <- 10^-8*((kk[2]/(kk[1]-S))+(kk[4]/(kk[3]-S)))
 rvs  <- 1.022*10^-8*(ww[1]+S*ww[2]+ww[3]*S^2+ww[4]*S^3)
 Ma  <-  0.0289635 + 1.2011 * 10^(-8) * (CO2 - 400)
 raxs  <- ras*(1+5.34*10^-7*(CO2-450))
 TT <- *Temp*+273.15

That was the first line with an error. Why are you using "*"'s in that manner?


#Saturation vapor pressure over ice
# A1 <- -13.928169
# A2 <- 34.7078238
# Theta <- (*Temp*+273.15)/273.16
# Y  <-  A1 * (1 - Theta^(-1.5)) + A2 * (1 - Theta^(-1.25))
# Psv <- 611.657*exp(Y)

 #Saturation vapor pressure over water
 Tc  <-  647.096       # K   : Temperature at the critical point
 Pc  <-  22.064 * 10^4 # hPa : Vapor pressure at the critical point
 nu  <-  (1-TT/Tc)
 a1  <-  -7.85951783
 a2  <-  1.84408259
 a3  <-  -11.7866497
 a4  <-  22.6807411
 a5  <-  -15.9618719
 a6  <-  1.80122502
 Psv  <-  Pc * exp(Tc/TT * (a1*nu + a2*nu^1.5 + a3*nu^3. + a4*nu^3.5 +
a5*nu^4. + a6*nu^7.5))

Alpha  <-  1.00062
Beta  <-  3.14 * 10^(-8)
Gamma  <-  5.60 * 10^(-7)
fpt <- Alpha+Beta*Press+Gamma**Temp*^2
Xv <- (RH/100)*fpt*Psv/Press

 Zm <- 1-(Press/TT)*(aa[1]+aa[2]**Temp*+aa[3]**Temp*^2
                     +(bb[1]+bb[2]*Temp)*Xv+(cc[1]+cc[2]**Temp*)*Xv^2)
+(Press/TT)^2*(dd+ee*Xv^2)
Roaxs <- pp*Ma/(zz*R*tt)
Rov <- Xv*Press*M/(Zm*R*TT)
Roa <- (1-Xv)*Press*Ma/(Zm*R*TT)
n <- 1+(Roa/Roaxs)*raxs+(Rov/ro)*rvs
return(n)
}

David Winsemius, MD
West Hartford, CT

______________________________________________
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.

Reply via email to