• API
• FAQ
• Tools
• Archive
daily pastebin goal
63%
SHARE
TWEET

# Untitled

a guest Oct 12th, 2018 52 Never
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
1. model {
2.   for(i in 2:nday){
3.     dayl = 0
4.     pi = 2*arcsin(1)
5.     Tr = .5*exp(p10*.5*(mint[i] + maxt[i]))
6.     lai = max(.1, Cf/LMA)
7.     trange = (maxt[i] - mint[i])/2
8.     gs = pow(2, 0.789798) / (0.37836)
9.     pp = lai*Nit / gs *p11*exp(0.011136*maxt[i])
10.     qq = -204.6453
11.     Caqqpp.sq = pow((Ca[i] + qq - pp), 2)
12.     ci = .5*(Ca[i] + qq - pp + pow( (Caqqpp.sq) - 4*(Ca[i]*qq - pp*4.22273), 1/2 ))
13.
14.     e0 = 7.19298*pow(lai, 2) / (pow(lai, 2) + 2.1001)
15.
16.     dec = -23.4*cos(360*(yearday+10)/365 *pi/180)*pi / 180
17.     mult = tan(dec)*tan(lat)
18.
19.     dayl = ifelse(abs(mult) < 1, 24*arccos(-mult)/pi, ifelse(mult >= 1, 24, 0) )
20.
21.     cps = e0*rad[i]*gs*(Ca[i] - ci) / (e0*rad[i] + gs*(Ca[i] - ci))
22.     G[i] = cps*(0.0156935*dayl + 0.0453194)
23.
24.     Cf.pred[i] = (1-p5)*Cf[i-1] + G[i]*(1 - p2)*p3
25.     Cw.pred[i] = (1 - p6)*Cw[i-1] + G[i]*(1-p2)*(1-p3)*(1-p4)
26.     Cr.pred[i] = (1 - p7)*Cr[i-1] + G*(1-p2)*(1-p3)*p4
27.     Clit.pred[i] = Cf[i-1]*p5 + Cr[i-1]*p7 + (1 - Tr*(p8+p1))*Clit[i-1]
28.     Csom.pred[i] = p6*Cw[i-1] + Tr*p1*Clit[i-1] + (1 - Tr*p9)*Csom[i-1]
29.
30.     Cf[i] ~ dnorm(Cf.pred[i], sd_cf)
31.     Cw[i] ~ dnorm(Cw.pred[i], sd_cw)
32.     Cr[i] ~ dnorm(Cr.pred[i], sd_cr)
33.     Clit[i] ~ dnorm(Clit.pred[i], sd_clit)
34.     Csom[i] ~ dnorm(Csom.pred[i], sd_csom)
35.   }
36.
37.   for(i in 1:nday){
38.     Cf.obs[i] ~ dnorm(Cf[i], sd_cf_obs)
39.     Cw.obs[i] ~ dnorm(Cw[i], sd_cw_obs)
40.     Cr.obs[i] ~ dnorm(Cr[i], sd_cr_obs)
41.     Clit.obs[i] ~ dnorm(Clit[i], sd_clit_obs)
42.     Csom.obs[i] ~ dnorm(Csom[i], sd_csom_obs)
43.   }
44.
45.   ## Initial conditions
46.   Cf[1] ~ dnorm(100, 2)
47.   Cw[1] ~ dnorm(9200, 25)
48.   Cr[1] ~ dnorm(100, 2)
49.   Clit[1] ~ dnorm(20, 2)
50.   Csom[1] ~ dnorm(11000, 25)
51.
52.   ## Priors on parameters
53.   p1 ~ dunif((10^(-6)), .01)
54.   p2 ~ dunif(.2, .7)
55.   p3 ~ dunif(.01, .5)
56.   p4 ~ dunif(.01, .5)
57.   p5 ~ dunif((10^(-4)), .1)
58.   p6 ~ dunif((10^(-6)), .01)
59.   p7 ~ dunif((10^(-4)), .01)
60.   p8 ~ dunif((10^(-5)), .1)
61.   p9 ~ dunif((10^(-6)), .01)
62.   p10 ~ dunif(.05, .2)
63.   p11 ~ dunif(5, 20)
64.
65.   ## Priors on process errors
66.   sd_cf ~ dnorm(0, 100)T(0,)
67.   sd_cw ~ dnorm(0, 100)T(0,)
68.   sd_cr ~ dnorm(0, 100)T(0,)
69.   sd_clit ~ dnorm(0, 100)T(0,)
70.   sd_csom ~ dnorm(0, 100)T(0,)
71.
72.   sd_cf_obs ~ dnorm(0, 100)T(0,)
73.   sd_cw_obs ~ dnorm(0, 100)T(0,)
74.   sd_cr_obs ~ dnorm(0, 100)T(0,)
75.   sd_clit_obs ~ dnorm(0, 100)T(0,)
76.   sd_csom_obs ~ dnorm(0, 100)T(0,)
77. }
RAW Paste Data
We use cookies for various purposes including analytics. By continuing to use Pastebin, you agree to our use of cookies as described in the Cookies Policy.

Top