| Title: | Tools for hydrological model |
|---|---|
| Description: | What the package does (one paragraph). |
| Authors: | Dongdong Kong [aut, cre] (ORCID: <https://orcid.org/0000-0003-1836-8172>) |
| Maintainer: | Dongdong Kong <[email protected]> |
| License: | MIT + file LICENSE |
| Version: | 0.1.8 |
| Built: | 2026-05-14 09:20:50 UTC |
| Source: | https://github.com/CUG-hydro/hydroTools |
Budyko Curve
budyko(PET.P, par = 2, method = "Fu1981") ET_budyko(PET, P, par = 2, method = "Fu1981", ...) budyko_goal(AET.P, PET.P, par = 2, method = "Fu1981", ...) budyko_fit(AET.P, PET.P, par = 2, method = "Fu1981", ...)budyko(PET.P, par = 2, method = "Fu1981") ET_budyko(PET, P, par = 2, method = "Fu1981", ...) budyko_goal(AET.P, PET.P, par = 2, method = "Fu1981", ...) budyko_fit(AET.P, PET.P, par = 2, method = "Fu1981", ...)
PET.P |
PET divide P,
|
par |
scalar value, the parameter of Budyko curve |
method |
one of |
PET |
potential evapotranspiration (mm) |
P |
precipitation (mm) |
... |
ignored |
AET.P |
|
Budyko, M. I., Climate and Life, 508 pp., Academic, San Diego, Calif., 1974.
Fu, Baw-Puh., 1981. On the Calculation of the Evaporation from Land Surface. Chinese Journal of Atmospheric Sciences, 5(1), 23-31.
Zhang, L., Dawes, W. R., & Walker, G. R. (2001). Response of mean annual evapotranspiration to vegetation changes at catchment scale. Water Resources Research, 37(3), 701–708. https://doi.org/10.1029/2000WR900325
Yang, H., Yang, D., Lei, Z., & Sun, F. (2008). New analytical derivation of the mean annual water-energy balance equation. Water Resources Research. https://doi.org/10.1029/2007WR006135
Rotation angle of two speed vector
cal_angle(v1, v2, clockwise = -1) cal_azimuth(v, refer = c(0, 1))cal_angle(v1, v2, clockwise = -1) cal_azimuth(v, refer = c(0, 1))
v1, v2
|
speed vector, in deg, (lon, lat).
|
clockwise |
one of |
v2 -> v1:
rotation angle (in deg, [-180, 180])
counterclockwise is positive
clockwise is negative
azimuth angle (in deg, [-180, 180])
the angle with north (c(0, 1))
cal_angle(c(1, 0), c(2, 0)) # == 0 cal_angle(c(0, 1), c(0, 2)) # == 0 p1 <- c(1, 0) p2 <- c(0, 1) cal_angle(p1, p2) # == 90, counterclockwise cal_angle(p2, p1) # == -90, clockwisecal_angle(c(1, 0), c(2, 0)) # == 0 cal_angle(c(0, 1), c(0, 2)) # == 0 p1 <- c(1, 0) p2 <- c(0, 1) cal_angle(p1, p2) # == 90, counterclockwise cal_angle(p2, p1) # == -90, clockwise
black body radiation
cal_Rl_out(Ta, emiss = 0.97)cal_Rl_out(Ta, emiss = 0.97)
Ta |
temperature in C deg |
emiss |
Emissivity coefficient |
radiation in W m-2
Estimate Incoming longwave radiation. Not included in FAO56 paper but added for convinience.
cal_Rli( temp, ea = NULL, s = 1, method = c("MAR", "SWI", "IJ", "BRU", "SAT", "KON") )cal_Rli( temp, ea = NULL, s = 1, method = c("MAR", "SWI", "IJ", "BRU", "SAT", "KON") )
temp |
Near surface air temperature |
ea |
Near surface actual vapour pressure |
s |
Cloud air emissivity, the ratio between acutal incomming shortwave radiation and clear sky incoming shortwave radiation. |
method |
The method to esteimate the air emissivity. Must be one of 'MAR', 'SWI', 'IJ', 'BRU', 'T', and 'KON'. |
Incomming longwave radiation [W /m2].
Satterlund, D. R. (1979), An improved equation for estimating long-wave radiation from the atmosphere, Water Resour. Res., 15( 6), 1649– 1650, doi:10.1029/WR015i006p01649.
Sedlar, J., & Hock, R. (2009). Testing longwave radiation parameterizations under clear and overcast skies at Storglaciären, Sweden. The Cryosphere, 3(1), 75-84.
Net outgoing longwave radiation.
cal_Rln(tmax, tmin, ea, Rsi = NULL, Rsi_o = NULL, cld = NULL) cal_Rln_yang2019( Ts, Rsi, Rsi_toa, lat = 30, emiss = 0.96, n1 = 2.52, n2 = 2.37, n3 = 0.035 )cal_Rln(tmax, tmin, ea, Rsi = NULL, Rsi_o = NULL, cld = NULL) cal_Rln_yang2019( Ts, Rsi, Rsi_toa, lat = 30, emiss = 0.96, n1 = 2.52, n2 = 2.37, n3 = 0.035 )
tmax |
Daily maximum air temperature at 2m height |
tmin |
Daily minimum air temperature at 2m height |
ea |
Actual vapor pressure |
Rsi |
Surface incoming short-wave radiation (Rsi) |
Rsi_o |
Clear sky incoming shortwave radiation, i. e. extraterrestrial
radiation multiply by clear sky transmissivity (i. e. a + b,
a and b are coefficients of Angstrom formula. Normally 0.75)
|
cld |
Cloud cover |
Ts |
land surface temperature |
Rsi_toa |
incoming short-wave radiation at the top of the atmosphere |
lat |
latitude (in deg) |
emiss |
Emissivity |
n1, n2, n3
|
parameter for |
Net outgoing longwave radiation [MJ m-2 day-1]
Yang, Y., & Roderick, M. L. (2019). Radiation, surface temperature and evaporation over wet surfaces. Quarterly Journal of the Royal Meteorological Society, 145(720), 1118–1129. https://doi.org/10.1002/qj.3481
Tu, Z., & Yang, Y. (2022). On the Estimation of Potential Evaporation Under Wet and Dry Conditions. Water Resources Research, 58(4). https://doi.org/10.1029/2021wr031486
Surface net radiation
cal_Rn( lat, J, Tmin, Tmax, ea = NULL, RH, ssd = NULL, cld = NULL, Rsi = NULL, albedo = 0.23, Z = 0, ... )cal_Rn( lat, J, Tmin, Tmax, ea = NULL, RH, ssd = NULL, cld = NULL, Rsi = NULL, albedo = 0.23, Z = 0, ... )
lat |
float, latitude |
J |
integer, day of year |
Tmin, Tmax
|
numeric vector, min and max 2m-air temperature (degC) |
ea |
Actual vapor pressure (kPa) |
RH |
Relative humidity (%). If |
ssd |
numeric vector, sun shine duration (hour) |
cld |
(optional) cloud coverage (0-1). At least one of |
Rsi |
(optional) Surface downward shortwave radiation (MJ m-2 d-1).
|
albedo |
(optional), |
Z |
(optional) elevation (m), for the calculation of |
... |
ignored |
radiation in [MJ d-1]
Rn : Surface net radiation
Rln : Surface outward net longwave radiation (negative means outgoing)
Rsn : Surface downward net shortwave radiation
Rsi : Surface downward shortwave radiation (Rsi)
Rsi_o : Clear-sky surface downward shortwave radiation
Rsi_toa : Extraterrestrial radiation (top of atmosphere, Rsi_toa) #'
Rn might <= 0. Users need to constrain the min value by pmax(Rn, 0).
Xie YuXuan and Kong Dongdong
cal_Rn(lat = 30, J = 1, RH = 70, Tmin = 20, Tmax = 30, ssd = 10)cal_Rn(lat = 30, J = 1, RH = 70, Tmin = 20, Tmax = 30, ssd = 10)
Daily inward shortwave solar radiation at crop surface [MJ m-2 day-1] by
providing sunshine duration (SSD) in hours or cloud cover in fraction.
cal_Rsi(lat, J, ssd = NULL, cld = NULL, Z = 0, a = 0.25, b = 0.5)cal_Rsi(lat, J, ssd = NULL, cld = NULL, Z = 0, a = 0.25, b = 0.5)
lat |
Latitude |
J |
day of year |
ssd |
sunshine duration lubridate::hours. If |
cld |
Cloud cover |
Z |
elevation (m) |
a |
Coefficient of the Angstrom formula. Determine the relationship between ssd and radiation. Default 0.25. |
b |
Coefficient of the Angstrom formula. Default 0.50. |
A data.table, solar radiation at crop surface [MJ m-2 day-1].
Rsi_toa: Clear-sky surface downward shortwave radiation
Rsi: Surface downward shortwave radiation
Martinezlozano J A, Tena F, Onrubia J E, et al. The historical evolution of the Angstrom formula and its modifications: review and bibliography.data.table::J. Agricultural & Forest Meteorology, 1985, 33(2):109-128.
https://github.com/sbegueria/SPEI/blob/master/R/penman.R
Estimate daily extraterrestrial radiation [MJ m-2 day-1].
cal_Rsi_toa(lat, J)cal_Rsi_toa(lat, J)
lat |
Latitude |
J |
|
extraterrestrial radiation [MJ m-2 day-1].
cal_Rsi_toa(20, 1) # [1] 25.84874cal_Rsi_toa(20, 1) # [1] 25.84874
Calculating sunset hour angle (ws) according to Allen Eq. 25.
cal_sunset_angle(lat, J) ws2ssd(ws) cal_ssd(lat, J)cal_sunset_angle(lat, J) ws2ssd(ws) cal_ssd(lat, J)
lat |
Latitude |
J |
|
ws |
sunset hour angle |
Evaporative surface temperature (land surface, Tland; or leaf surface Tleaf). Land surface temperature infered by Monteith 1965 Equation.
cal_Ts( Rn, Tair, D, U2, Pa = atm, rH = NULL, rs = 0, method = c("simple", "full", "ma2021"), ... )cal_Ts( Rn, Tair, D, U2, Pa = atm, rH = NULL, rs = 0, method = c("simple", "full", "ma2021"), ... )
Rn |
land surface net radiation, W m-2 |
Tair |
2m air temperature (degC) |
D |
vapor pressure deficit (kPa) |
U2 |
wind speed at 2m |
Pa |
surface air pressure (kPa) |
rH |
conductance for heat |
rs |
If |
method |
|
... |
ignored |
rH: aerodynamic resistance of heat
rs: stamotal resistance of water
Monteith, J. P. (1965). An introduction to environmental physics.
Yang, Y., & Roderick, M. L. (2019). Radiation, surface temperature and evaporation over wet surfaces. Quarterly Journal of the Royal Meteorological Society, 145(720), 1118–1129. https://doi.org/10.1002/qj.3481
Ma, N., Szilagyi, J., & Zhang, Y. (2021). Calibration-Free Complementary Relationship Estimates Terrestrial Evapotranspiration Globally. Water Resources Research, 57(9), 1–27. https://doi.org/10.1029/2021WR029691
library(dplyr) # The relationship of Rn, Ts Rn = 0:200 cal_Ts(Rn, 25, D = 1, U2 = 2) cal_Ts(200, 25, D = 1, U2 = 2) # plot(Rn, dat$Ts, type = "l", main = "(a) Ts ~ Rn") # plot(Rn, dat$Eeq, type = "l", main = "(b) Eeq ~ Rn") # # plot(Rn, dat$Evp, type = "l") # a constant value # dat %<>% mutate(Rn = Rn, bowen = ET0 / (Rn * 0.086400 / lambda - ET0)) # plot(bowen ~ Rn, dat, type = "l", main = "(b) Eeq ~ Rn")library(dplyr) # The relationship of Rn, Ts Rn = 0:200 cal_Ts(Rn, 25, D = 1, U2 = 2) cal_Ts(200, 25, D = 1, U2 = 2) # plot(Rn, dat$Ts, type = "l", main = "(a) Ts ~ Rn") # plot(Rn, dat$Eeq, type = "l", main = "(b) Eeq ~ Rn") # # plot(Rn, dat$Evp, type = "l") # a constant value # dat %<>% mutate(Rn = Rn, bowen = ET0 / (Rn * 0.086400 / lambda - ET0)) # plot(bowen ~ Rn, dat, type = "l", main = "(b) Eeq ~ Rn")
wetbulb temperature
cal_wetbulb(Tair, Tdew, Pa = atm) cal_Tw(ea, Tair, Pa = atm)cal_wetbulb(Tair, Tdew, Pa = atm) cal_Tw(ea, Tair, Pa = atm)
Tair |
2m air temperature (degC) |
Pa |
surface air pressure (kPa) |
ea |
Actual vapor pressure (kPa) |
https://github.com/wrf-model/WRF/blob/master/phys/module_diag_functions.F#L1154
T0: Absolute zero in
Kelvin (K)
Es.T0: is the saturation vapor pressure at the
absolute zero .
L: Latent heat of water
vapor
atm: standard pressure at sea surface, 101.325 kPa
Cp: specific heat at constant pressure, 1.013 * 1e-3 MJ kg-1 degC-1
T0 K0 L Es.T0 atmT0 K0 L Es.T0 atm
An object of class numeric of length 1.
An object of class numeric of length 1.
An object of class numeric of length 1.
An object of class numeric of length 1.
An object of class numeric of length 1.
This function allows you to convert a vector of temperature values between Fahrenheit, Celsius, and degrees Kelvin.
convert_temperature(temperature, old_metric, new_metric) celsius.to.fahrenheit(T.celsius) celsius.to.kelvin(T.celsius) fahrenheit.to.celsius(T.fahrenheit) fahrenheit.to.kelvin(T.fahrenheit) kelvin.to.celsius(T.kelvin) kelvin.to.fahrenheit(T.kelvin)convert_temperature(temperature, old_metric, new_metric) celsius.to.fahrenheit(T.celsius) celsius.to.kelvin(T.celsius) fahrenheit.to.celsius(T.fahrenheit) fahrenheit.to.kelvin(T.fahrenheit) kelvin.to.celsius(T.kelvin) kelvin.to.fahrenheit(T.kelvin)
temperature |
A numeric vector of temperatures to be converted. |
old_metric |
The metric from which you want to convert. Possible options are:
|
new_metric |
The metric to which you want to convert. The same options
are possible as for |
T.kelvin, T.fahrenheit, T.celsius
|
Numeric vector of temperatures in Kelvin, Fahrenheit and Celsius. |
A numeric vector with temperature converted to the metric specified
by the argument new_metric.
A numeric vector of temperature values in the new unit.
Joshua Ferreri [email protected], Brooke Anderson [email protected], Roger Peng [email protected]
correct_basin_net_shapefile
correct_basin_net_shapefile( ID, shp, tree, outdir = "basins/subbasin", prefix = "shed_珠江", overwrite = FALSE )correct_basin_net_shapefile( ID, shp, tree, outdir = "basins/subbasin", prefix = "shed_珠江", overwrite = FALSE )
ID |
the ID of basins |
shp |
basin shapefile, sf polygon object |
tree |
returned by |
outdir |
outdir of shapefile |
prefix |
prefix of shapefile |
detect_flood_events
flood_divide
detect_groups(df, inds, gap_max = 5, extend = ddays(5)) detect_flood_events( date, Q, Q_min = 2, Q_peak = 10, gap_max = 5, extend = ddays(5) ) flood_divide(df, ...) merge_flood(df, info_flood, format = "%Y.%m")detect_groups(df, inds, gap_max = 5, extend = ddays(5)) detect_flood_events( date, Q, Q_min = 2, Q_peak = 10, gap_max = 5, extend = ddays(5) ) flood_divide(df, ...) merge_flood(df, info_flood, format = "%Y.%m")
df |
A data.table with date and Q columns |
gap_max |
Default
|
extend |
Default |
Q_min |
minimum discharge to detect flood events |
Q_peak |
peak discharge to detect flood events |
... |
parameters passed to |
ET complementary
ET_CR_Xiao2020( Rn, Tair, D, U2, Pa = atm, method = c("simple", "full", "ma2021"), rs = 0 ) ET_CR_Ma2021( Rn, Tair, D, U2, Pa = atm, method = c("simple", "full", "ma2021"), rs = 0 )ET_CR_Xiao2020( Rn, Tair, D, U2, Pa = atm, method = c("simple", "full", "ma2021"), rs = 0 ) ET_CR_Ma2021( Rn, Tair, D, U2, Pa = atm, method = c("simple", "full", "ma2021"), rs = 0 )
Rn |
land surface net radiation, W m-2 |
Tair |
2m air temperature (degC) |
D |
vapor pressure deficit (kPa) |
U2 |
2m wind speed (m/s) |
Pa |
surface air pressure (kPa) |
method |
method for |
rs |
stamotal conductance for water on capony scale |
Xiao, M., Yu, Z., Kong, D., Gu, X., Mammarella, I., Montagnani, L., … Gioli, B. (2020). Stomatal response to decreased relative humidity constrains the acceleration of terrestrial evapotranspiration. Environmental Research Letters, 15(9). doi:10.1088/1748-9326/ab9967
Ma, N., Szilagyi, J., & Zhang, Y. (2021). Calibration-Free Complementary Relationship Estimates Terrestrial Evapotranspiration Globally. Water Resources Research, 57(9), 1–27. https://doi.org/10.1029/2021WR029691
ET_CR_Xiao2020(250, 25, D = 1, U2 = 2) ET_CR_Ma2021(250, 25, D = 1, U2 = 2)ET_CR_Xiao2020(250, 25, D = 1, U2 = 2) ET_CR_Ma2021(250, 25, D = 1, U2 = 2)
lambda: latent heat of vaporization, about [2.5 MJ kg-1].
slope: The slope of the saturation vapour pressure curve at certain air
temperature Tair, [kPa degC-1].
gamma: psychrometric constant ([kPa degC-1]), Cp*Pa/(epsilon*lambda).
U2: 10m wind speed (m/s). According to wind profile relationship, convert U_z to U_2.
es: saturation vapor pressure (kPa)
ea: actual vapor pressure (kPa)
RH_mean | Tmin, Tmax: (es(Tmax) + es(Tmin)) * RH_mean/200
Tmin : es(Tmin)
cal_TvK: vitual temperature (K)
Tair + q
Tair + ea/Pa
1.01 * (Tair + 273) , FAO56 Eq. 3-7
cal_rho_a: air density (kg/m^3)
rH: resistance of aerodynamic (s/m), about 100 s/m.
cal_U2(Uz, z.wind = 10) cal_ea(Tmin, Tmax = NULL, RH_mean = NULL) cal_VPD(Tair, Tdew = NULL) cal_lambda(Tair) cal_slope(Tair) cal_gamma(Tair, Pa = atm) cal_bowen(Tair, Pa = atm) cal_Pa(z = NULL) cal_rH(U2, h = 0.12) cal_rH2(U2, Tair, Pa = atm) cal_TvK(Tair, q = NULL, ea = NULL, Pa = atm) cal_rho_a(Tair, Pa = atm, q = NULL, ea = NULL) cal_es(Tair)cal_U2(Uz, z.wind = 10) cal_ea(Tmin, Tmax = NULL, RH_mean = NULL) cal_VPD(Tair, Tdew = NULL) cal_lambda(Tair) cal_slope(Tair) cal_gamma(Tair, Pa = atm) cal_bowen(Tair, Pa = atm) cal_Pa(z = NULL) cal_rH(U2, h = 0.12) cal_rH2(U2, Tair, Pa = atm) cal_TvK(Tair, q = NULL, ea = NULL, Pa = atm) cal_rho_a(Tair, Pa = atm, q = NULL, ea = NULL) cal_es(Tair)
Uz |
wind speed at height |
z.wind |
Height where measure the wind speed |
Tmin |
Daily minimum air temperature at 2m height |
Tmax |
Daily maximum air temperature at 2m height |
RH_mean |
daily mean relative humidity |
Tair |
2m air temperature, ( |
Tdew |
dew temperature, ( |
Pa |
land surface Air pressure |
z |
Elevation above sea level |
U2 |
wind speed at 2m |
h |
canopy height (m), for the calculation of the following intermediate variables:
|
q |
specific humidity (g/g) |
ea |
actual vapor pressure (kPa) |
Allen, R. G., & Luis S. Pereira. (1998). Crop evapotranspiration-Guidelines for computing crop water requirements-FAO Irrigation and drainage paper 56. European Journal of Agronomy, 34(3), 144–152. doi:10.1016/j.eja.2010.12.001
cal_VPD(10, 5) par(mfrow = c(2, 2), mar = c(3, 1.8, 2, 1), mgp = c(2, 0.6, 0)) Tair <- -10:50 plot(cal_es(-10:50), main = "es (kPa)", xlab = "Tair") plot(cal_bowen(-10:50), main = "Bowen ratio", xlab = "Tair") plot(cal_slope(-10:50), main = "slope (kPa/degC)", xlab = "Tair") plot(cal_gamma(-10:50), main = "gamma (kPa/degC)", xlab = "Tair")cal_VPD(10, 5) par(mfrow = c(2, 2), mar = c(3, 1.8, 2, 1), mgp = c(2, 0.6, 0)) Tair <- -10:50 plot(cal_es(-10:50), main = "es (kPa)", xlab = "Tair") plot(cal_bowen(-10:50), main = "Bowen ratio", xlab = "Tair") plot(cal_slope(-10:50), main = "slope (kPa/degC)", xlab = "Tair") plot(cal_gamma(-10:50), main = "gamma (kPa/degC)", xlab = "Tair")
ET0_eq : Equilibrium evaporation, slope / (slope + gamma) * Rn
ET0_Penman48 : Penman 1948 equation simplified by Shuttleworth 1993
ET0_Monteith65: Penman-Monteith 1965
ET0_PT72 : Priestley Taylor 1972, ET0_PT72 = ET0_eq * 1.26
ET0_FAO98 : Penman-Monteith reference crop evapotranspiration, FAO56
ET0_eq(Rn, Tair, Pa = atm, ...) ET0_Penman48(Rn, Tair, Pa = atm, D, wind, z.wind = 2) ET0_Monteith65( Rn, Tair, Pa = atm, D, wind, z.wind = 2, rs = 70, rH = NULL, ... ) ET0_PT72(Rn, Tair, Pa = atm, alpha = 1.26, ...) ET0_FAO98(Rn, Tair, Pa = atm, D, wind, z.wind = 2, tall_crop = FALSE)ET0_eq(Rn, Tair, Pa = atm, ...) ET0_Penman48(Rn, Tair, Pa = atm, D, wind, z.wind = 2) ET0_Monteith65( Rn, Tair, Pa = atm, D, wind, z.wind = 2, rs = 70, rH = NULL, ... ) ET0_PT72(Rn, Tair, Pa = atm, alpha = 1.26, ...) ET0_FAO98(Rn, Tair, Pa = atm, D, wind, z.wind = 2, tall_crop = FALSE)
Rn |
land surface net radiation, W m-2 |
Tair |
2m air temperature (degC) |
Pa |
surface air pressure (kPa) |
D |
vapor pressure deficit (kPa) |
wind |
wind speed at the height of |
z.wind |
Height where measure the wind speed |
alpha |
albedo |
tall_crop |
tall crop or not, see FAO1998 for details |
lambda: Latent heat of vaporization, about [2.5 MJ kg-1].
slope : The slope of the saturation vapour pressure curve at certain air
temperature Tair, [kPa degC-1].
gamma : Psychrometric constant ([kPa degC-1]), Cp*Pa/(epsilon*lambda)
Eeq : Equilibrium (radiative) evaporation, mm
Evp : Aerodynamic evaporation, mm
ET0 : Potential evapotranspiration, mm
Penman equation. Wikipedia 2021. https://en.wikipedia.org/w/index.php?title=Penman_equation
Shuttleworth, W. J. (1993) Evaporation. In: Handbook of Hydrology (ed. by D. Maidment). McGraw-Hill, New York.
Allen, R. G., & Luis S. Pereira. (1998). Crop evapotranspiration-Guidelines for computing crop water requirements-FAO Irrigation and drainage paper 56. European Journal of Agronomy, 34(3), 144-152. https://doi.org/10.1016/j.eja.2010.12.001
Allen, R. G., Pereira, L. S., Raes, D., & Smith, M. (1998). FAO Irrigation and drainage paper No. 56. Rome: Food and Agriculture Organization of the United Nations, 56(97), e156.
Chapter 2 - FAO Penman-Monteith equation, https://www.fao.org/3/x0490E/x0490e06.htm
ET0_PT72(250, 25, D = 1) ET0_Penman48(250, 25, D = 1, wind = 2) ET0_FAO98(250, 25, D = 1, wind = 2) ET0_Monteith65(250, 25, D = 1, wind = 2)ET0_PT72(250, 25, D = 1) ET0_Penman48(250, 25, D = 1, wind = 2) ET0_FAO98(250, 25, D = 1, wind = 2) ET0_Monteith65(250, 25, D = 1, wind = 2)
Get local time
get_localtime(time, lon, timeZone = 8)get_localtime(time, lon, timeZone = 8)
time |
POSIXct object |
lon |
longitude in deg |
heat_index converts a numeric scalar of temperature
(in Fahrenheit) and a numeric scalar of relative humidity (in %)
to heat index (in Fahrenheit). This function is not meant to be used
outside of the heat_index_vec() function.
heat_index(t, rh) .heat.index(t = NA, rh = NA) heat_index_vec(t = NA, rh = NA)heat_index(t, rh) .heat.index(t = NA, rh = NA) heat_index_vec(t = NA, rh = NA)
t |
Numeric scalar of air temperature, in celsius. |
rh |
Numeric scalar of relative humidity, in %. |
If an impossible value of relative humidity is given
(below 0% or above 100%), heat index is returned as NA.
A numeric scalar of heat index, in celsius
Equations are from the source code for the US National Weather Service's online heat index calculator.
Brooke Anderson [email protected], Roger Peng [email protected]
Anderson GB, Bell ML, Peng RD. 2013. Methods to calculate the heat index as an exposure Metric in environmental health research. Environmental Health Perspectives 121(10):1111-1119.
National Weather Service Hydrometeorological Prediction Center Web Team. Heat Index Calculator. 30 Jan 2015. http://www.wpc.ncep.noaa.gov/html/heatindex.shtml. Accessed 18 Dec 2015.
Rothfusz L. 1990. The heat index (or, more than you ever wanted to know about heat index) (Technical Attachment SR 90-23). Fort Worth: Scientific Services Division, National Weather Service.
R. Steadman, 1979. The assessment of sultriness. Part I: A temperature-humidity index based on human physiology and clothing science. Journal of Applied Meteorology, 18(7):861–873.
n <- 1e3 t <- seq(-50, 100, length.out = n) rh <- seq(0, 100, length.out = length(t)) ## bug exist i <- 88 heat_index(t[i], rh[i]) heat_index_vec(t[i], rh[i]) # heat_index_julia(t[i], rh[i]) # julia_setup() r <- heat_index(t, rh) r_vec <- heat_index_vec(t, rh) # r_jl <- heat_index_julia(t, rh) # add tests for julia NA values # r1 <- heat_index_julia(t*NA, rh*NA) all.equal(r, r_vec) # all.equal(r_jl, r_vec) # microbenchmark::microbenchmark( # r2 <- heat_index(t, rh), # r1 <- heat_index_julia(t, rh) # # r_vec <- heat_index_vec(t, rh) # )n <- 1e3 t <- seq(-50, 100, length.out = n) rh <- seq(0, 100, length.out = length(t)) ## bug exist i <- 88 heat_index(t[i], rh[i]) heat_index_vec(t[i], rh[i]) # heat_index_julia(t[i], rh[i]) # julia_setup() r <- heat_index(t, rh) r_vec <- heat_index_vec(t, rh) # r_jl <- heat_index_julia(t, rh) # add tests for julia NA values # r1 <- heat_index_julia(t*NA, rh*NA) all.equal(r, r_vec) # all.equal(r_jl, r_vec) # microbenchmark::microbenchmark( # r2 <- heat_index(t, rh), # r1 <- heat_index_julia(t, rh) # # r_vec <- heat_index_vec(t, rh) # )
LCL : Bolton 1980, Eq. 15
LCL_RH: Bolton 1980, Eq. 22
LCL(P0, T0, Td) LCL_RH(T, RH)LCL(P0, T0, Td) LCL_RH(T, RH)
P0 |
pressure at surface (hPa) |
T0 |
temperature at surface (C) |
Td |
dew point temperature (C) |
Mw: Molecular weight of dry air
Md: Molecular weight of water vapor
epsilon: the ratio of Mw to Md
Mw Md epsilonMw Md epsilon
An object of class numeric of length 1.
An object of class numeric of length 1.
An object of class numeric of length 1.
Good of fitting
NSE(yobs, ysim, w, ...) KGE(yobs, ysim, ...) GOF(yobs, ysim, w, include.cv = FALSE, include.r = TRUE)NSE(yobs, ysim, w, ...) KGE(yobs, ysim, ...) GOF(yobs, ysim, w, include.cv = FALSE, include.r = TRUE)
yobs |
Numeric vector, observations |
ysim |
Numeric vector, corresponding simulated values |
w |
Numeric vector, weights of every points. If w included, when calculating mean, Bias, MAE, RMSE and NSE, w will be taken into considered. |
... |
ignored |
include.cv |
If true, cv will be included. |
include.r |
If true, r and R2 will be included. |
RMSE root mean square error
NSE NASH coefficient
MAE mean absolute error
AI Agreement index (only good points (w == 1)) participate to
calculate. See details in Zhang et al., (2015).
Bias bias
Bias_perc bias percentage
n_sim number of valid obs
cv Coefficient of variation
R2 correlation of determination
R pearson correlation
pvalue pvalue of R
https://en.wikipedia.org/wiki/Coefficient_of_determination
https://en.wikipedia.org/wiki/Explained_sum_of_squares
https://en.wikipedia.org/wiki/Nash%E2%80%93Sutcliffe_model_efficiency_coefficient
Zhang Xiaoyang (2015), http://dx.doi.org/10.1016/j.rse.2014.10.012
yobs = rnorm(100) ysim = yobs + rnorm(100)/4 GOF(yobs, ysim)yobs = rnorm(100) ysim = yobs + rnorm(100)/4 GOF(yobs, ysim)
Estimate reference evapotranspiration (ET0) from a hypothetical short grass reference surface using the the Hargreaves equation.
PET_hg(Tmax, Tmin, Tavg = NULL, Rsi_toa = NULL, lat = NULL, dates = NULL)PET_hg(Tmax, Tmin, Tavg = NULL, Rsi_toa = NULL, lat = NULL, dates = NULL)
Tmax |
Daily maximum air temperature at 2m height |
Tmin |
Daily minimum air temperature at 2m height |
Tavg |
Daily mean air temperature at 2m height |
Rsi_toa |
Clear sky incoming shortwave radiation, i. e. extraterrestrial
radiation multiply by clear sky transmissivity (i. e. a + b, a and b are
coefficients of Angstrom formula. Normally 0.75) |
lat |
Latitude |
dates |
A R Date type of a vector of Date type. If not provided, it will Regard the ssd series is begin on the first day of a year. |
Reference evapotranspiration ET0 from a hypothetical grass reference
surface [mm day-1].
plot_calib
plot_calib(data, ...)plot_calib(data, ...)
data |
with the columns of |
... |
others to plot precipitation |
plot_runoff
plot_runoff( df_q, df_prcp = NULL, xlim = NULL, ylim2 = c(50, 0), legend.position = c(1, 1), legend.justification = c(1, 1), ... )plot_runoff( df_q, df_prcp = NULL, xlim = NULL, ylim2 = c(50, 0), legend.position = c(1, 1), legend.justification = c(1, 1), ... )
df_q |
A data.frame with the columns of |
df_prcp |
A data.frame with the columns of |
xlim |
limit of x axis |
ylim2 |
limit of second y axis (for precipitation) |
... |
other parameters passed to |
plot_tree
plot_StreamNet(info, title = "", fout = "tree.pdf", show = FALSE, root = -1)plot_StreamNet(info, title = "", fout = "tree.pdf", show = FALSE, root = -1)
info |
A data.table with the column of |
title |
The title of the tree |
R: gas constant, J/(mol K)
Rw: Specific gas constant
of water vapor .
Rd: Specific gas constant of dry air (J/(kg K)).
R Rw Rd CpR Rw Rd Cp
An object of class numeric of length 1.
An object of class numeric of length 1.
An object of class numeric of length 1.
An object of class numeric of length 1.
air temperature under solar radiation
radiation_Ta(Rs, Rl = 0, albedo = 0.3, emiss = 0.97) cal_radiation_Ta( Ta, Rs, Rl = 0, dt = 3600, method = c("approx", "exact"), ..., albedo = 0.3, emiss = 0.97 )radiation_Ta(Rs, Rl = 0, albedo = 0.3, emiss = 0.97) cal_radiation_Ta( Ta, Rs, Rl = 0, dt = 3600, method = c("approx", "exact"), ..., albedo = 0.3, emiss = 0.97 )
Rs, Rl
|
numeric vector, inward shortwave and longwave radiation, in |
Ta |
scalar, the initial air temperature, in |
dt |
delta |
soil texture class, based on USDA triangle
soil_class(clay, sand)soil_class(clay, sand)
clay |
Vector, percent of clay |
sand |
Vector, percent of sand |
Character vector of soil texture class, one of
c("Cl", "SiCl", "SaCl", "ClLo", "SiClLo", "SaClLo", "Lo", "SiLo", "SaLo", "Si", "LoSa", "Sa").
Julien Moeys (2016). soiltexture: Functions for Soil Texture Plot, Classification and Transformation. R package version 1.5.1. https://CRAN.R-project.org/package=soiltexture
clay = c(05, 60, 15, 05, 25, 05, 25, 45, 65, 75, 13, 47) sand = c(90, 32, 70, 70, 20, 10, 10, 10, 20, 10, 70, 10) soil_class(clay, sand) # c("Sa", "Cl", "SaLo", "SaLo", "SiLo", "Si", "SiLo", "SiCl", "Cl", "Cl", "SaLo", "SiCl")clay = c(05, 60, 15, 05, 25, 05, 25, 45, 65, 75, 13, 47) sand = c(90, 32, 70, 70, 20, 10, 10, 10, 20, 10, 70, 10) soil_class(clay, sand) # c("Sa", "Cl", "SaLo", "SaLo", "SiLo", "Si", "SiLo", "SiCl", "Cl", "Cl", "SaLo", "SiCl")
Estimate monthly soil heat flux (G) from the mean air temperature of the previous, current or next month assuming as grass crop.
soil_heat_flux(t.p, t.n = NULL, t.c = NULL)soil_heat_flux(t.p, t.n = NULL, t.c = NULL)
t.p |
Mean air temperature of the previous month |
t.n |
Mean air temperature of the next month |
t.c |
Mean air temperature of the current month |
Soil heat flux [MJ m-2 day-1].
Wilting point, field capacity and saturated moisture
wilting point: 1500 kPa moisture, %v, (Eq. 1)
field capacity: 33 kPa moisture, %v, (Eq. 2)
saturated moisture: Saturation (0 kPa moisture), %v, (Eq. 5)
wilting_point(S, C, OM) field_capacity(S, C, OM) saturated_mois(S, C, OM) wilting_point_salinity(S = 0.2, C = 0.2, OM = 2.5, EC = 3)wilting_point(S, C, OM) field_capacity(S, C, OM) saturated_mois(S, C, OM) wilting_point_salinity(S = 0.2, C = 0.2, OM = 2.5, EC = 3)
S |
the weight ratio of Sand, (weight ratio, 0-1) |
C |
the weight ratio of Clay, (weight ratio, 0-1) |
OM |
the percentage of Organic Matter, (%w, 0-100) |
EC |
Electrical conductance of a saturated soil extract (dS/m = mili-mho/cm) |
%v
The unit of OM is %w, which is different from S and C.
Saxton, K. E., & Rawls, W. J. (2006). Soil water characteristic estimates by texture and organic matter for hydrologic solutions. Soil Science Society of America Journal, 70(5), 1569-1578. https://doi.org/10.2136/sssaj2005.0117
S = C = 0.20; OM = 2.5 wilting_point(S, C, OM) # 13.8% field_capacity(S, C, OM) # 32.1% saturated_mois(S, C, OM) # 48.2%S = C = 0.20; OM = 2.5 wilting_point(S, C, OM) # 13.8% field_capacity(S, C, OM) # 32.1% saturated_mois(S, C, OM) # 48.2%
The default location is ZuoLing.
suncalc( time, lon = 114.6053, lat = 30.49694, ..., year = year(Sys.time()), verbose = TRUE ) sunrise(date = Sys.Date(), lon = 80, lat = 30, timeZone = 8)suncalc( time, lon = 114.6053, lat = 30.49694, ..., year = year(Sys.time()), verbose = TRUE ) sunrise(date = Sys.Date(), lon = 80, lat = 30, timeZone = 8)
suncalc(Sys.Date()) sunrise(Sys.Date())suncalc(Sys.Date()) sunrise(Sys.Date())
equivalent potential temperature
theta(P0, T0) theta_wet(P0, T0, Td) theta_wet_bolton(P0, T0, Td)theta(P0, T0) theta_wet(P0, T0, Td) theta_wet_bolton(P0, T0, Td)
P0 |
pressure at surface (hPa) |
T0 |
temperature at surface (Cdeg) |
Td |
dew point temperature (Cdeg) |
https://unidata.github.io/MetPy/latest/api/generated/metpy.calc.equivalent_potential_temperature.html
https://github.com/wrf-model/WRF/blob/master/phys/module_diag_functions.F#L122
https://github.com/Unidata/MetPy/blob/e0e24d51702787943fc3c0481fa9a6632abe9d20/src/metpy/calc/thermo.py#L1512
theta_wet(850, 20, 18) theta_wet_bolton(850, 20, 18)theta_wet(850, 20, 18) theta_wet_bolton(850, 20, 18)
R2Q: Q(m^3/s) = R(mm) * area (km^2) * 1000/(3600*24)
W2mm: 1 MJ /m2/day =0.408 mm /day, 1 Watt /m2 = 0.0864 MJ /m2/day
R2Q(R, area = dt * 3.6, dt = 24) Q2R(Q, area = dt * 3.6, dt = 24) MJ_2W(x) MJ_2mm(x) W2_MJ(x) W2mm(x, tmean = 0) K2T(x) T2K(x) deg2rad(deg) rad2deg(rad)R2Q(R, area = dt * 3.6, dt = 24) Q2R(Q, area = dt * 3.6, dt = 24) MJ_2W(x) MJ_2mm(x) W2_MJ(x) W2mm(x, tmean = 0) K2T(x) T2K(x) deg2rad(deg) rad2deg(rad)
R |
runoff (mm per dt) |
area |
basin area (km^2), if area = , |
dt |
time duration (hour) |
Q |
runoff flow (m^3/s) |
x |
scalar or numeric vector |
tmean |
daily mean temperature |
deg |
angle in degree |
rad |
angle in radian |
http://www.fao.org/3/X0490E/x0490e0i.htm
USDA soil class classification
usda_sfusda_sf
An object of class sf (inherits from data.frame) with 12 rows and 2 columns.
https://en.wikipedia.org/wiki/Soil_texture
Helper functions for vapour pressure
w2q(w) q2w(q, Pa = atm) q2ea(q, Pa = atm) w2ea(w, Pa = atm) ea2VPD(ea, RH) vapour_press(q, Pa = atm) q2RH(q, Tair, Pa = atm) RH2q(RH, Tair, Pa = atm) q_from_RH(RH, Tair, Pa = atm) cal_qs(Tair, Pa = atm) cal_ws(Tair, Pa = atm) ea2w(ea, Pa = atm) ea2q(ea, Pa = atm) RH2q(RH, Tair, Pa = atm) Tdew2q(Tdew, Pa = atm) Tdew2w(Tdew, Pa = atm) Tdew2RH(Tdew, Tair) Tdew_from_q(q, Pa = atm) Tdew_from_w(w, Pa = atm)w2q(w) q2w(q, Pa = atm) q2ea(q, Pa = atm) w2ea(w, Pa = atm) ea2VPD(ea, RH) vapour_press(q, Pa = atm) q2RH(q, Tair, Pa = atm) RH2q(RH, Tair, Pa = atm) q_from_RH(RH, Tair, Pa = atm) cal_qs(Tair, Pa = atm) cal_ws(Tair, Pa = atm) ea2w(ea, Pa = atm) ea2q(ea, Pa = atm) RH2q(RH, Tair, Pa = atm) Tdew2q(Tdew, Pa = atm) Tdew2w(Tdew, Pa = atm) Tdew2RH(Tdew, Tair) Tdew_from_q(q, Pa = atm) Tdew_from_w(w, Pa = atm)
w |
mix ratio, m_w / m_d |
q |
specific humidity in kg/kg or g/g |
Pa |
surface air pressure |
ea |
actual vapor pressure (kPa) |
RH |
relative humidity, in % |
Tair |
air temperature, in degC |
Tdew |
dew temperature (in degC) |
e in the same unit as Pa
https://www.eol.ucar.edu/projects/ceop/dm/documents/refdata_report/eqns.html, Eq-17
Bolton, David. “The Computation of Equivalent Potential Temperature.” Monthly Weather Review 108, no. 7 (July 1980): 1046–53. <https://doi.org/10.1175/1520-0493(1980)108<1046:TCOEPT>2.0.CO;2>.
https://earthscience.stackexchange.com/questions/2360/how-do-i-convert-specific-humidity-to-relative-humidity
RH = 90 Pa = atm Tair = 30 q = RH2q(RH, Pa, Tair) RH2 = q2RH(q, Pa, Tair) e = q2ea(q, Pa) es = cal_es(Tair)RH = 90 Pa = atm Tair = 30 q = RH2q(RH, Pa, Tair) RH2 = q2RH(q, Pa, Tair) e = q2ea(q, Pa) es = cal_es(Tair)
Soil parameters for VIC model (calculated from HWSD database)
expt : 3+2b (Eq. 17)
Kstat: Saturated conductivity (mm/hr) (Eq. 16)
bubble: Bubbling pressure (cm) (Eq. 4)
Wcr_FT : Field Capacity, Fractional soil moisture content at the
critical point (~70% of field capacity) (fraction of maximum moisture)
Wpwp_FT : Wilting Point, Fractional soil moisture content at the
wilting point (fraction of maximum moisture)
expt(S, C, OM) Ksat(S, C, OM) bubble(S, C, OM) Wcr_FT(S, C, OM) Wpwp_FT(S, C, OM)expt(S, C, OM) Ksat(S, C, OM) bubble(S, C, OM) Wcr_FT(S, C, OM) Wpwp_FT(S, C, OM)
S |
the weight ratio of Sand, (weight ratio, 0-1) |
C |
the weight ratio of Clay, (weight ratio, 0-1) |
OM |
the percentage of Organic Matter, (%w, 0-100) |
Kstat: unit has changed from mm/h to mm/day
Saxton, K. E., & Rawls, W. J. (2006). Soil water characteristic estimates by texture and organic matter for hydrologic solutions. Soil Science Society of America Journal, 70(5), 1569-1578. https://doi.org/10.2136/sssaj2005.0117
S = C = 0.20; OM = 2.5 expt(S, C, OM) Ksat(S, C, OM) # 12.19 bubble(S, C, OM) field_capacity(S, C, OM) # 32.1% wilting_point(S, C, OM) # 13.8% # Fractional soil moisture Wcr_FT(S, C, OM) Wpwp_FT(S, C, OM)S = C = 0.20; OM = 2.5 expt(S, C, OM) Ksat(S, C, OM) # 12.19 bubble(S, C, OM) field_capacity(S, C, OM) # 32.1% wilting_point(S, C, OM) # 13.8% # Fractional soil moisture Wcr_FT(S, C, OM) Wpwp_FT(S, C, OM)
XAJ model Parameter calibration
XAJ_calib( Qobs, prcp, ET0, area, dt = 24, date = NULL, maxn = 1000, index = "KGE", seed = 1, ... )XAJ_calib( Qobs, prcp, ET0, area, dt = 24, date = NULL, maxn = 1000, index = "KGE", seed = 1, ... )
Qobs |
Observed Total runoff, (m^3/s) |
prcp |
Precipitation (mm/d) |
ET0 |
Pan evaporation or potential evapotranspiration (mm/d) |
area |
basin area (km^2). |
dt |
time step (hour) |
date |
(optional) corresponding date of |
index |
KGE or NSE |
seed |
(can be ignored) starting number of random number generator,
see |
... |
ignored |
XAJ_predict
XAJ_predict(model, newdata)XAJ_predict(model, newdata)
model |
XAJ model returned by |
newdata |
A data.frame, with the column of |