||
全国580+ 站点,20分钟全搞定。
详情邮件联系:liyj AT msu.edu
部分代码如下:
mutate(Jd = dayid,
P = 101.3*(((293-0.0065*h)/293)^5.26), # P: atmospheric pressure [kPa]
# h: elevation above sea level [m]
lat_rad = deg2rad(lat),
cp = 1.013*10^(-3), # cp: specific heat at constant pressure, 1.013 10-3 [MJ kg-1 °C-1]
e = 0.622, # e: ratio molecular weight of water vapour/dry air = 0.622
lam = 2.45, # lam: latent heat of vaporization, 2.45 [MJ kg-1],
r = cp*P/e/lam, # psychrometric constant [kPa °C-1]
Tmean = (tmax+tmin)/2, # mean daily air temperature at 2 m height [°C]
e0T = 0.6108*exp(17.27*Tmean/(Tmean+237.3)), # e°(T) saturation vapour pressure at the air temperature T [kPa]
e0Tmax = 0.6108*exp(17.27*tmax/(tmax+237.3)),
e0Tmin = 0.6108*exp(17.27*tmin/(tmin+237.3)),
es = (e0Tmax+e0Tmin)/2, # saturation vapour pressure [kPa]
ea = (rh/100)*e0T, # actual vapour pressure [kPa]
# es - ea saturation vapour pressure deficit [kPa]
delta_low = 0.409*sin(2*3.1415927*Jd/365-1.39),
ws = acos(-tan(lat_rad)*tan(delta_low)),
n = tsun, # n actual duration of sunshine [hour]
N = (24/3.1415927)*ws, # N maximum possible duration of sunshine or daylight hours [hour]
dr = 1+0.033*cos(2*3.1415927*Jd/365),
Ra = (24*60*0.082*dr/3.1415926)*(ws*sin(lat_rad)*sin(delta_low)+cos(lat_rad)*cos(delta_low)*sin(ws)),
# Ra extraterrestrial radiation [MJ m-2 day-1]
a = 0.25, # empirical coefficient 0.12~0.29
b = 0.50, # empirical coefficient 0.45~0.73
Rs = (a + b*n/N)*Ra, # Rs solar or shortwave radiation [MJ m-2 day-1]
Rns = (1-0.23)*Rs,
Rso = (0.75+2*(10)^(-5)*h)*Ra,
Tkm = (tmax+273.2),
Tkn = (tmin+273.2),
Rnl = 4.903*(10^(-9))*(0.1+0.9*n/N)*(0.34-0.14*(ea^(0.5)))*((Tkm^4+Tkn^4)/2),
Rn = Rns - Rnl, # net radiation at the crop surface [MJ m-2 day-1]
Delta = 4098*es/((Tmean+237.3)^2), # slope vapour pressure curve [kPa °C-1]
uz = wind * 0.72, # uz measured wind speed at z m above ground surface [m s-1]
# should multiply 0.72 (empirical coefficient)to correct
z = 10, # z height of measurement above ground surface [m]
# measurement at 10 m (Changming Liu)
u2 = 4.87*uz/log(67.8*z-5.42), # u2 wind speed at 2 m height [m s-1]
G = 0, # soil heat flux density [MJ m-2 day-1]
# For day and ten-day period: Gday = 0
# For hourly or shorter periods: Ghr = 0.1*Rn (daylight periods); Ghr = 0.5*Rn (during nighttime period)
ET0 = (0.408*Delta*(Rn-G) + 900*r*u2*(es-ea)/(Tmean+273)) / (Delta+r*(1+0.34*u2)) # reference evapotranspiration [mm day-1]
)
Archiver|手机版|科学网 ( 京ICP备07017567号-12 )
GMT+8, 2024-4-25 04:23
Powered by ScienceNet.cn
Copyright © 2007- 中国科学报社