Package 'bigleaf'

Aerodynamic Conductance

Description

Bulk aerodynamic conductance, including options for the boundary layer conductance formulation and stability correction functions.

Usage

aerodynamic.conductance(
  data,
  Tair = "Tair",
  pressure = "pressure",
  wind = "wind",
  ustar = "ustar",
  H = "H",
  zr,
  zh,
  d,
  z0m = NULL,
  Dl,
  N = 2,
  fc = NULL,
  LAI,
  Cd = 0.2,
  hs = 0.01,
  wind_profile = FALSE,
  stab_correction = TRUE,
  stab_formulation = c("Dyer_1970", "Businger_1971"),
  Rb_model = c("Thom_1972", "Choudhury_1988", "Su_2001", "constant_kB-1"),
  kB_h = NULL,
  Sc = NULL,
  Sc_name = NULL,
  constants = bigleaf.constants()
)

Arguments

data	Data.frame or matrix containing all required variables
Tair	Air temperature (degC)
pressure	Atmospheric pressure (kPa)
wind	Wind speed (m s$^{-1}$)
ustar	Friction velocity (m s$^{-1}$)
H	Sensible heat flux (W m$^{-2}$)
zr	Instrument (reference) height (m)
zh	Canopy height (m)
d	Zero-plane displacement height (m)
z0m	Roughness length for momentum (m), optional; if not provided, it is estimated from roughness.parameters (method="wind_profile"). Only used if wind_profile = TRUE and/or Rb_model = "Su_2001" or "Choudhury_1988".
Dl	Characteristic leaf dimension (m) (if Rb_model = "Su_2001") or leaf width (if Rb_model = "Choudhury_1988"); ignored otherwise.
N	Number of leaf sides participating in heat exchange (1 or 2); only used if Rb_model = "Su_2001". Defaults to 2.
fc	Fractional vegetation cover (-); only used if Rb_model = "Su_2001". See Details.
LAI	One-sided leaf area index (m$^2$ m$^{-2}$); only used if Rb_model = "Choudhury_1988" or "Su_2001".
Cd	Foliage drag coefficient (-); only used if Rb_model = "Su_2001".
hs	Roughness length of bare soil (m); only used if Rb_model = "Su_2001".
wind_profile	Should Ga for momentum be calculated based on the logarithmic wind profile equation? Defaults to FALSE.
stab_correction	Should stability correction be applied? Defaults to TRUE. Ignored if wind_profile = FALSE.
stab_formulation	Stability correction function. Either "Dyer_1970" (default) or "Businger_1971". Ignored if wind_profile = FALSE or if stab_correction = FALSE.
Rb_model	Boundary layer resistance formulation. One of "Thom_1972","Choudhury_1988","Su_2001","constant_kB-1".
kB_h	$kB^{-1}$ value for heat transfer; only used if Rb_model = "constant_kB-1"
Sc	Optional: Schmidt number of additional quantities to be calculated
Sc_name	Optional: Name of the additional quantities, has to be of same length than Sc_name
constants	k - von Karman constant cp - specific heat of air for constant pressure (J K$^{-1}$ kg$^{-1}$) Kelvin - conversion degree Celsius to Kelvin g - gravitational acceleration (m s$^{-2}$) pressure0 - reference atmospheric pressure at sea level (Pa) Tair0 - reference air temperature (K) Sc_CO2 - Schmidt number for CO$_2$ Pr - Prandtl number (if Sc is provided)

Details

Aerodynamic conductance for heat (Ga_h) is calculated as:

$Ga_h = 1 / (Ra_m + Rb_h)$

where Ra_m is the aerodynamic resistance for momentum and Rb the (quasi-laminar) canopy boundary layer resistance ('excess resistance').

The aerodynamic resistance for momentum Ra_m is given by:

$Ra_m = u/ustar^2$

Note that this formulation accounts for changes in atmospheric stability, and does not require an additional stability correction function.

An alternative method to calculate Ra_m is provided (calculated if wind_profile = TRUE):

$Ra_m = (ln((zr - d)/z0m) - psi_h) / (k \: ustar)$

If the roughness parameters z0m and d are unknown, they can be estimated using roughness.parameters. The argument stab_formulation determines the stability correction function used to account for the effect of atmospheric stability on Ra_m (Ra_m is lower for unstable and higher for stable stratification). Stratification is based on a stability parameter zeta (z-d/L), where z = reference height, d the zero-plane displacement height, and L the Monin-Obukhov length, calculated with Monin.Obukhov.length The stability correction function is chosen by the argument stab_formulation. Options are "Dyer_1970" and "Businger_1971".

The model used to determine the canopy boundary layer resistance for heat (Rb_h) is specified by the argument Rb_model. The following options are implemented: "Thom_1972" is an empirical formulation based on the friction velocity (ustar) (Thom 1972):

$Rb_h = 6.2ustar^{-0.667}$

The model by Choudhury & Monteith 1988 (Rb_model = "Choudhury_1988"), calculates Rb_h based on leaf width, LAI and ustar (Note that function argument Dl represents leaf width (w) and not characteristic leaf dimension (Dl) if Rb_model = "Choudhury_1988"):

$Gb_h = LAI((0.02/\alpha)*\sqrt{u(zh)/w}*(1-exp(-\alpha/2)))$

where $\alpha$ is a canopy attenuation coefficient modeled in dependence on LAI, u(zh) is wind speed at canopy height (calculated from wind.profile), and w is leaf width (m). See Gb.Choudhury for further details.

The option Rb_model = "Su_2001" calculates Rb_h based on the physically-based Rb model by Su et al. 2001, a simplification of the model developed by Massman 1999:

$kB_h = (k \: Cd \: fc^2) / (4Ct \: ustar/u(zh)) + kBs^{-1}(1 - fc)^2$

where Cd is a foliage drag coefficient (defaults to 0.2), fc is fractional vegetation cover, $Bs^{-1}$ is the inverse Stanton number for bare soil surface, and Ct is a heat transfer coefficient. See Gb.Su for details on the model.

The models calculate the parameter kB$^{-1}$ (in the code referred to as kB_h), which is related to Rb_h:

$kB_h = Rb_h * (k * ustar)$

From version 0.7.6 onwards, the roughness length for heat (z0h) is added to the output if z0m is available (i.e. provided as input or calculated within this function). z0h is calculated from roughness.length.heat:

$z0h = z0m / exp(kB_h)$

Rb (and Gb) for water vapor and heat are assumed to be equal in this package. Gb for other quantities x is calculated as (Hicks et al. 1987):

$Gb_x = Gb / (Sc_x / Pr)^{0.67}$

where $Sc_x$ is the Schmidt number of quantity x, and Pr is the Prandtl number (0.71).

Value

a data.frame with the following columns:

Ga_m	Aerodynamic conductance for momentum transfer (m s$^{−1}$)
Ra_m	Aerodynamic resistance for momentum transfer (s m$^{-1}$)
Ga_h	Aerodynamic conductance for heat transfer (m s$^{-1}$)
Ra_h	Aerodynamic resistance for heat transfer (s m$^{-1}$)
Gb_h	Canopy boundary layer conductance for heat transfer (m s$^{-1}$)
Rb_h	Canopy boundary layer resistance for heat transfer (s m$^{-1}$)
kB_h	$kB^{-1}$ parameter for heat transfer
z0h	Roughness length for heat (m) (NA if not input z0m not provided as input or not estimated in this function)
zeta	Stability parameter 'zeta' (NA if wind_profile = FALSE)
psi_h	Integrated stability correction function (NA if wind_profile = FALSE)
Ra_CO2	Aerodynamic resistance for CO$_2$ transfer (s m$^{-1}$)
Ga_CO2	Aerodynamic conductance for CO$_2$ transfer (m s$^{-1}$)
Gb_CO2	Canopy boundary layer conductance for CO$_2$ transfer (m s$^{-1}$)
Ga_Sc_name	Aerodynamic conductance for Sc_name (m s$^{-1}$). Only added if Sc_name and the respective Sc are provided
Gb_Sc_name	Boundary layer conductance for Sc_name (m s$^{-1}$). Only added if Sc_name and the respective Sc are provided

Note

Input variables such as LAI, Dl, or zh can be either constants, or vary with time (i.e. vectors of the same length as data).

Note that boundary layer conductance to water vapor transfer ($Gb_w$) is often assumed to equal $Gb_h$. This assumption is also made in this R package, for example in the function surface.conductance.

If the roughness length for momentum (z0m) is not provided as input, it is estimated from the function roughness.parameters within wind.profile if wind_profile = TRUE and/or Rb_model = "Su_2001" or "Choudhury_1988" The roughness.parameters function estimates a single z0m value for the entire time period! If a varying z0m value (e.g. across seasons or years) is required, z0m should be provided as input argument.

References

Verma, S., 1989: Aerodynamic resistances to transfers of heat, mass and momentum. In: Estimation of areal evapotranspiration, IAHS Pub, 177, 13-20.

Verhoef, A., De Bruin, H., Van Den Hurk, B., 1997: Some practical notes on the parameter kB-1 for sparse vegetation. Journal of Applied Meteorology, 36, 560-572.

Hicks, B.B., Baldocchi, D.D., Meyers, T.P., Hosker, J.R., Matt, D.R., 1987: A preliminary multiple resistance routine for deriving dry deposition velocities from measured quantities. Water, Air, and Soil Pollution 36, 311-330.

Monteith, J.L., Unsworth, M.H., 2008: Principles of environmental physics. Third Edition. Elsevier Academic Press, Burlington, USA.

See Also

Gb.Thom, Gb.Choudhury, Gb.Su for calculations of Rb / Gb only

Examples

df <- data.frame(Tair=25,pressure=100,wind=c(3,4,5),ustar=c(0.5,0.6,0.65),H=c(200,230,250))   
 
# simple calculation of Ga  
aerodynamic.conductance(df,Rb_model="Thom_1972") 

# calculation of Ga using a model derived from the logarithmic wind profile
aerodynamic.conductance(df,Rb_model="Thom_1972",zr=40,zh=25,d=17.5,z0m=2,wind_profile=TRUE) 

# simple calculation of Ga_m, but a physically based canopy boundary layer model
aerodynamic.conductance(df,Rb_model="Su_2001",zr=40,zh=25,d=17.5,Dl=0.05,N=2,fc=0.8)

---

Air Density

Description

Air density of moist air from air temperature and pressure.

Usage

air.density(Tair, pressure, constants = bigleaf.constants())

Arguments

Tair	Air temperature (degC)
pressure	Atmospheric pressure (kPa)
constants	Kelvin - conversion degC to Kelvin Rd - gas constant of dry air (J kg$^{-1}$ K$^{-1}$) kPa2Pa - conversion kilopascal (kPa) to pascal (Pa)

Details

Air density ($\rho$) is calculated as:

$\rho = pressure / (Rd * Tair)$

Value

$\rho$ - air density (kg m$^{-3}$)

References

Foken, T, 2008: Micrometeorology. Springer, Berlin, Germany.

Examples

# air density at 25degC and standard pressure (101.325kPa)
air.density(25,101.325)

---

(Modified) Arrhenius Temperature Response Function

Description

(Modified) Arrhenius function describing the temperature response of biochemical parameters.

Usage

Arrhenius.temp.response(
  param,
  Temp,
  Ha,
  Hd,
  dS,
  constants = bigleaf.constants()
)

Arguments

param	Parameter measured at measurement temperature (umol m$^{-2}$ s$^{-1}$)
Temp	Measurement temperature (degC)
Ha	Activation energy for param (kJ mol$^{-1}$)
Hd	Deactivation energy for param (kJ mol$^{-1}$)
dS	Entropy term for param (kJ mol$^{-1}$ K$^{-1}$)
constants	Kelvin - conversion degree Celsius to Kelvin Rgas - universal gas constant (J mol$^{-1}$ K$^{-1}$) kJ2J - conversion kilojoule (kJ) to joule (J)

Details

The function returns the biochemical rate at a reference temperature of 25degC given a predefined temperature response function. This temperature response is given by a modified form of the Arrhenius function:

$param25 = param / ( exp(Ha * (Temp - Tref) / (Tref*Rgas*Temp)) * (1 + exp((Tref*dS - Hd) / (Tref * Rgas))) / (1 + exp((Temp*dS - Hd) / (Temp * Rgas))) )$

where param is the value/rate of the parameter at measurement temperature, Temp is temperature in K, Tref is reference temperature (298.15K), and Rgas is the universal gas constant (8.314 J K$^{-1}$ mol$^{-1}$). Ha is the activation energy (kJ mol$^{-1}$), Hd is the deactivation energy (kJ mol$^{-1}$), and dS the entropy term (kJ mol$^{-1}$ K$^{-1}$) of the respective parameter.

If either Hd or dS or both are not provided, the equation above reduces to the first term (i.e. the common Arrhenius equation without the deactivation term.)

Value

param25 - value of the input parameter at the reference temperature of 25degC (umol m$^{-2}$ s$^{-1}$)

References

Johnson F.H., Eyring H., Williams R.W. 1942: The nature of enzyme inhibitions in bacterial luminescence: sulfanilamide, urethane, temperature and pressure. Journal of cellular and comparative physiology 20, 247-268.

Kattge J., Knorr W., 2007: Temperature acclimation in a biochemical model of photosynthesis: a reanalysis of data from 36 species. Plant, Cell and Environment 30, 1176-1190.

---

Eddy Covariance Data of AT-Neu (Neustift)

Description

Halfhourly eddy covariance Data of the site AT-Neu, a mountain meadow in Austria. (https://sites.fluxdata.org/AT-Neu/). Data are from July 2010.

Usage

AT_Neu_Jul_2010

Format

A data frame with 1488 observations and 31 columns:

year

year of measurement

month

month of measurement

doy

day of year

hour

hour (0 - 23.5)

Tair

Air temperature (degC) [TA_F]

Tair_qc

Quality control of Tair [TA_F_QC]

PPFD

Photosynthetic photon flux density (umol m$^{-2}$ s$^{-1}$) [PPFD_IN]

PPFD_qc

Quality control of PPFD [PPFD_IN_QC]

VPD

Vapor pressure deficit (kPa) [VPD_F]

VPD_qc

Quality control of VPD [VPD_F_QC]

pressure

Atmospheric pressure (kPa) [PA_F]

precip

precipitation (mm) [P_F]

precip_qc

Quality control of precip [P_F_QC]

ustar

friction velocity (m s$^{-1}$) [USTAR]

wind

horizontal wind velocity (m s$^{-1}$) [WS_F]

wind_qc

Quality control of wind [WS_F_QC]

Ca

Atmospheric CO$_2$ concentration (ppm) [CO2_F_MDS]

Ca_qc

Quality control of Ca [CO2_F_MDS_QC]

LW_up

upward longwave radiation (W m$^{-2}$) [LW_OUT]

Rn

Net radiation (W m$^{-2}$) [NETRAD]

LE

Latent heat flux (W m$^{-2}$) [LE_F_MDS]

LE_qc

Quality control of LE [LE_F_MDS_QC]

H

Sensible heat flux (W m$^{-2}$) [H_F_MDS]

H_qc

Quality control of H [H_F_MDS_QC]

G

Ground heat flux (W m$^{-2}$) [G_F_MDS]

G_qc

Quality control of G [G_F_MDS_QC]

NEE

Net ecosystem exchange (umol m$^{-2}$ s$^{-1}$) [NEE_VUT_USTAR50]

NEE_qc

Quality control of NEE [NEE_VUT_USTAR50_QC]

GPP

Gross primary productivity from nighttime partitioning (umol m$^{-2}$ s$^{-1}$) [GPP_NT_VUT_USTAR50]

GPP_qc

Quality control of GPP [NEE_VUT_USTAR50_QC]

Reco

Ecosystem respiration from nighttime partitioning (umol m$^{-2}$ s$^{-1}$) [RECO_NT_VUT_USTAR50]

Note

The original variable names as provided by the FLUXNET2015 dataset are given in squared brackets. Note that variable units have been converted in some cases (e.g. VPD from hPa to kPa).

Source

original data were downloaded from https://fluxnet.org/ (accessed 09 November 2016)

---

Constants Used in the bigleaf Package

Description

This function defines the following constants:

Usage

bigleaf.constants(
  cp = 1004.834,
  Rgas = 8.31451,
  Rv = 461.5,
  Rd = 287.0586,
  Md = 0.0289645,
  Mw = 0.0180153,
  eps = 0.622,
  g = 9.81,
  solar_constant = 1366.1,
  pressure0 = 101325,
  Tair0 = 273.15,
  k = 0.41,
  Cmol = 0.012011,
  Omol = 0.0159994,
  H2Omol = 0.01801528,
  sigma = 5.670367e-08,
  Pr = 0.71,
  Sc_CO2 = 1.07,
  Le067 = 0.93,
  Kelvin = 273.15,
  DwDc = 1.6,
  days2seconds = 86400,
  kPa2Pa = 1000,
  Pa2kPa = 0.001,
  umol2mol = 1e-06,
  mol2umol = 1e+06,
  kg2g = 1000,
  g2kg = 0.001,
  kJ2J = 1000,
  J2kJ = 0.001,
  se_median = 1.253,
  frac2percent = 100
)

Arguments

cp	Specific heat of air for constant pressure (J K$^{-1}$ kg$^{-1}$)
Rgas	Universal gas constant (J mol$^{-1}$ K$^{-1}$)
Rv	Gas constant of water vapor (J kg$^{-1}$ K$^{-1}$) (Stull 1988 p.641)
Rd	Gas constant of dry air (J kg$^{-1}$ K$^{-1}$) (Foken p. 245)
Md	Molar mass of dry air (kg mol$^{-1}$)
Mw	Molar mass of water vapor (kg mol$^{-1}$)
eps	Ratio of the molecular weight of water vapor to dry air (=Mw/Md)
g	Gravitational acceleration (m s$^{-2}$)
solar_constant	Solar constant (W m$^{-2}$)
pressure0	Reference atmospheric pressure at sea level (Pa)
Tair0	Reference air temperature (K)
k	von Karman constant
Cmol	Molar mass of carbon (kg mol$^{-1}$)
Omol	Molar mass of oxygen (kg mol$^{-1}$)
H2Omol	Molar mass of water (kg mol$^{-1}$)
sigma	Stefan-Boltzmann constant (W m$^{-2}$ K$^{-4}$)
Pr	Prandtl number
Sc_CO2	Schmidt number for CO$_2$
Le067	Lewis number for water vapor to the power of 0.67
Kelvin	Conversion degree Celsius to Kelvin
DwDc	Ratio of the molecular diffusivities for water vapor and CO$_2$
days2seconds	Seconds per day
kPa2Pa	Conversion kilopascal (kPa) to pascal (Pa)
Pa2kPa	Conversion pascal (Pa) to kilopascal (kPa)
umol2mol	Conversion micromole (umol) to mole (mol)
mol2umol	Conversion mole (mol) to micromole (umol)
kg2g	Conversion kilogram (kg) to gram (g)
g2kg	Conversion gram (g) to kilogram (kg)
kJ2J	Conversion kilojoule (kJ) to joule (J)
J2kJ	Conversion joule (J) to kilojoule (kJ)
se_median	Conversion standard error (SE) of the mean to SE of the median
frac2percent	Conversion between fraction and percent

Details

This function is passed as an argument to every function that uses one or more constants. Individual constants passed to a function can be easily altered. E.g. the following command will change the value of the von Karman constant from 0.41 to 0.4:

bigleaf.constants(k=0.4)

the value of a constant can be returned by calling:

bigleaf.constants()$*name_of_constant*

To permanently change the constants contained within this function (which makes sense for some of them, e.g. for the von Karman constant), the command fixInNamespace can be used. E.g.

fixInNamespace(bigleaf.constants,ns="bigleaf")

Note that this has to be repeated every time the package is newly installed/loaded.

---

Biochemical Energy

Description

Radiant energy absorbed in photosynthesis or heat release by respiration calculated from net ecosystem exchange of CO$_2$ (NEE).

Usage

biochemical.energy(NEE, alpha = 0.422)

Arguments

NEE	Net ecosystem exchange (umol CO$_2$ m$^{-2}$ s$^{-1}$)
alpha	Energy taken up/released by photosynthesis/respiration per mol CO$_2$ fixed/respired (J umol$^{-1}$)

Details

The following sign convention is employed: NEE is negative when carbon is taken up by the ecosystem. Positive values of the resulting biochemical energy mean that energy (heat) is taken up by the ecosystem, negative ones that heat is released. The value of alpha is taken from Nobel 1974 (see Meyers & Hollinger 2004), but other values have been used (e.g. Blanken et al., 1997)

Value

Sp -	biochemical energy (W m$^{-2}$)

References

Meyers, T.P., Hollinger, S.E. 2004: An assessment of storage terms in the surface energy balance of maize and soybean. Agricultural and Forest Meteorology 125, 105-115.

Nobel, P.S., 1974: Introduction to Biophysical Plant Physiology. Freeman, New York.

Blanken, P.D. et al., 1997: Energy balance and canopy conductance of a boreal aspen forest: Partitioning overstory and understory components. Journal of Geophysical Research 102, 28915-28927.

Examples

# Calculate biochemical energy taken up by the ecosystem with 
# a measured NEE of -30umol CO2 m-2 s-1             
biochemical.energy(NEE=-30)

---

Eddy Covariance Data of DE-Tha (Tharandt)

Description

Halfhourly eddy covariance Data of the site DE-Tha, a spruce forest in Eastern Germany (https://sites.fluxdata.org/DE-Tha/). Data are from June 2014.

Usage

DE_Tha_Jun_2014

Format

A data frame with 1440 observations and 32 columns:

year

year of measurement

month

month of measurement

doy

day of year

hour

hour (0 - 23.5)

Tair

Air temperature (degC) [TA_F]

Tair_qc

Quality control of Tair [TA_F_QC]

PPFD

Photosynthetic photon flux density (umol m$^{-2}$ s$^{-1}$) [PPFD_IN]

PPFD_qc

Quality control of PPFD [PPFD_IN_QC]

VPD

Vapor pressure deficit (kPa) [VPD_F]

VPD_qc

Quality control of VPD [VPD_F_QC]

pressure

Atmospheric pressure (kPa) [PA_F]

precip

precipitation (mm) [P_F]

precip_qc

Quality control of precip [P_F_QC]

ustar

friction velocity (m s$^{-1}$) [USTAR]

wind

horizontal wind velocity (m s$^{-1}$) [WS_F]

wind_qc

Quality control of wind [WS_F_QC]

Ca

Atmospheric CO$_2$ concentration (ppm) [CO2_F_MDS]

Ca_qc

Quality control of Ca [CO2_F_MDS_QC]

LW_up

upward longwave radiation (W m$^{-2}$) [LW_OUT]

LW_down

downward longwave radiation (W m$^{-2}$) [LW_IN_F]

Rn

Net radiation (W m$^{-2}$) [NETRAD]

LE

Latent heat flux (W m$^{-2}$) [LE_F_MDS]

LE_qc

Quality control of LE [LE_F_MDS_QC]

H

Sensible heat flux (W m$^{-2}$) [H_F_MDS]

H_qc

Quality control of H [H_F_MDS_QC]

G

Ground heat flux (W m$^{-2}$) [G_F_MDS]

G_qc

Quality control of G [G_F_MDS_QC]

NEE

Net ecosystem exchange (umol m$^{-2}$ s$^{-1}$) [NEE_VUT_USTAR50]

NEE_qc

Quality control of NEE [NEE_VUT_USTAR50_QC]

GPP

Gross primary productivity from nighttime partitioning (umol m$^{-2}$ s$^{-1}$) [GPP_NT_VUT_USTAR50]

GPP_qc

Quality control of GPP [NEE_VUT_USTAR50_QC]

Reco

Ecosystem respiration from nighttime partitioning (umol m$^{-2}$ s$^{-1}$) [RECO_NT_VUT_USTAR50]

Note

The original variable names as provided by the FLUXNET2015 dataset are given in squared brackets. Note that variable units have been converted in some cases (e.g. VPD from hPa to kPa).

Source

original data were downloaded from https://fluxnet.org/ (accessed 09 November 2016)

---

Canopy-Atmosphere Decoupling Coefficient

Description

The canopy-atmosphere decoupling coefficient 'Omega'.

Usage

decoupling(
  data,
  Tair = "Tair",
  pressure = "pressure",
  Ga = "Ga_h",
  Gs = "Gs_ms",
  approach = c("Jarvis&McNaughton_1986", "Martin_1989"),
  LAI,
  Esat.formula = c("Sonntag_1990", "Alduchov_1996", "Allen_1998"),
  constants = bigleaf.constants()
)

Arguments

data	Data.frame or matrix containing all required input variables
Tair	Air temperature (deg C)
pressure	Atmospheric pressure (kPa)
Ga	Aerodynamic conductance to heat/water vapor (m s$^{-1}$)
Gs	Surface conductance (m s$^{-1}$)
approach	Approach used to calculate omega. Either "Jarvis&McNaughton_1986" (default) or "Martin_1989".
LAI	Leaf area index (m$^{2}$ m$^{-2}$), only used if approach = "Martin_1989".
Esat.formula	Optional: formula to be used for the calculation of esat and the slope of esat. One of "Sonntag_1990" (Default), "Alduchov_1996", or "Allen_1998". See Esat.slope.
constants	Kelvin - conversion degree Celsius to Kelvin cp - specific heat of air for constant pressure (J K$^{-1}$ kg$^{-1}$) eps - ratio of the molecular weight of water vapor to dry air (-) sigma - Stefan-Boltzmann constant (W m$^{-2}$ K$^{-4}$) Pa2kPa - conversion pascal (Pa) to kilopascal (kPa)

Details

The decoupling coefficient Omega ranges from 0 to 1 and quantifies the linkage of the conditions (foremost humidity and temperature) at the canopy surface to the ambient air. Values close to 0 indicate well coupled conditions characterized by high physiological (i.e. stomatal) control on transpiration and similar conditions at the canopy surface compared to the atmosphere above the canopy. Values close to 1 indicate the opposite, i.e. decoupled conditions and a low stomatal control on transpiration (Jarvis & McNaughton 1986).
The "Jarvis&McNaughton_1986" approach (default option) is the original formulation for the decoupling coefficient, given by (for an amphistomatous canopy):

$\Omega = \frac{\epsilon + 1}{\epsilon + 1 + \frac{Ga}{Gc}}$

where $\epsilon = \frac{s}{\gamma}$ is a dimensionless coefficient with s being the slope of the saturation vapor pressure curve (Pa K$^{-1}$), and $\gamma$ the psychrometric constant (Pa K$^{-1}$).

The approach "Martin_1989" by Martin 1989 additionally takes radiative coupling into account:

$\Omega = \frac{\epsilon + 1 + \frac{Gr}{Ga}}{\epsilon + (1 + \frac{Ga}{Gs}) (1 + \frac{Gr}{Ga})}$

Value

$\Omega$ - the decoupling coefficient Omega (-)

References

Jarvis P.G., McNaughton K.G., 1986: Stomatal control of transpiration: scaling up from leaf to region. Advances in Ecological Research 15, 1-49.

Martin P., 1989: The significance of radiative coupling between vegetation and the atmosphere. Agricultural and Forest Meteorology 49, 45-53.

See Also

aerodynamic.conductance, surface.conductance, equilibrium.imposed.ET

Examples

# Omega calculated following Jarvis & McNaughton 1986
set.seed(3)
df <- data.frame(Tair=rnorm(20,25,1),pressure=100,Ga_h=rnorm(20,0.06,0.01),
                 Gs_ms=rnorm(20,0.005,0.001))
decoupling(df,approach="Jarvis&McNaughton_1986")

# Omega calculated following Martin 1989 (requires LAI)
decoupling(df,approach="Martin_1989",LAI=4)

---

Dew Point

Description

calculates the dew point, the temperature to which air must be cooled to become saturated (i.e. e = Esat(Td))

Usage

dew.point(
  Tair,
  VPD,
  accuracy = 0.001,
  Esat.formula = c("Sonntag_1990", "Alduchov_1996", "Allen_1998"),
  constants = bigleaf.constants()
)

Arguments

Tair	Air temperature (degC)
VPD	Vapor pressure deficit (kPa)
accuracy	Accuracy of the result (deg C)
Esat.formula	Optional: formula to be used for the calculation of esat and the slope of esat. One of "Sonntag_1990" (Default), "Alduchov_1996", or "Allen_1998". See Esat.slope.
constants	Pa2kPa - conversion pascal (Pa) to kilopascal (kPa)

Details

Dew point temperature (Td) is defined by:

$e = Esat(Td)$

where e is vapor pressure of the air and Esat is the vapor pressure deficit. This equation is solved for Td using optimize.

Value

Td -	dew point temperature (degC)

References

Monteith J.L., Unsworth M.H., 2008: Principles of Environmental Physics. 3rd edition. Academic Press, London.

Examples

dew.point(c(25,30),1.5)

---

Energy Balance Closure

Description

Calculates the degree of the energy balance non-closure for the entire time span based on the ratio of two sums (energy balance ratio), and ordinary least squares (OLS).

Usage

energy.closure(
  data,
  Rn = "Rn",
  G = NULL,
  S = NULL,
  LE = "LE",
  H = "H",
  instantaneous = FALSE,
  missing.G.as.NA = FALSE,
  missing.S.as.NA = FALSE
)

Arguments

data	Data.frame or matrix containing all required variables.
Rn	Net radiation (W m$^{-2}$)
G	Ground heat flux (W m$^{-2}$); optional
S	Sum of all storage fluxes (W m$^{-2}$); optional
LE	Latent heat flux (W m$^{-2}$)
H	Sensible heat flux (W m$^{-2}$)
instantaneous	should the energy balance be calculated at the time step of the observations (TRUE), or over the entire time period provided as input (FALSE)
missing.G.as.NA	if TRUE, missing G are treated as NAs ,otherwise set to 0.
missing.S.as.NA	if TRUE, missing S are treated as NAs, otherwise set to 0.

Details

The energy balance ratio (EBR) is calculated as:

$EBR = sum(LE + H)/sum(Rn - G - S)$

the sum is taken for all time steps with complete observations (i.e. where all energy balance terms are available).

Value

a named vector containing:

n	number of complete (all energy balance terms available) observations
intercept	intercept of the OLS regression
slope	slope of the OLS regression
r_squared	r^2 of the OLS regression
EBR	energy balance ratio

if instantaneous = TRUE, only EBR is returned.

References

Wilson K., et al. 2002: Energy balance closure at FLUXNET sites. Agricultural and Forest Meteorology 113, 223-243.

Examples

## characterize energy balance closure for DE-Tha in June 2014
energy.closure(DE_Tha_Jun_2014,instantaneous=FALSE)

## look at half-hourly closure 
EBR_inst <- energy.closure(DE_Tha_Jun_2014,instantaneous=TRUE)
summary(EBR_inst)

---

Energy-Use Efficiency (EUE)

Description

Fraction of net radiation fixed by primary productivity.

Usage

energy.use.efficiency(GPP, alpha = 0.422, Rn)

Arguments

GPP	Gross primary productivity exchange (umol CO$_{2}$ m$^{-2}$ s$^{-1}$)
alpha	Energy taken up/released by photosynthesis/respiration (J umol$^{-1}$)
Rn	Net radiation (W m$^{-2}$)

Details

Energy use efficiency is calculated as:

$EUE = sum(GPP)/sum(Rn)$

where the sums are calculated for complete cases of GPP and Rn over the entire time period.

Value

EUE -

Energy use efficiency (-)

See Also

light.use.efficiency

Examples

energy.use.efficiency(GPP=20,Rn=500)

---

Equilibrium and Imposed Evapotranspiration

Description

Evapotranspiration (ET) split up into imposed ET and equilibrium ET.

Usage

equilibrium.imposed.ET(
  data,
  Tair = "Tair",
  pressure = "pressure",
  VPD = "VPD",
  Gs = "Gs_ms",
  Rn = "Rn",
  G = NULL,
  S = NULL,
  missing.G.as.NA = FALSE,
  missing.S.as.NA = FALSE,
  Esat.formula = c("Sonntag_1990", "Alduchov_1996", "Allen_1998"),
  constants = bigleaf.constants()
)

Arguments

data	Data.frame or matrix containing all required input variables
Tair	Air temperature (deg C)
pressure	Atmospheric pressure (kPa)
VPD	Air vapor pressure deficit (kPa)
Gs	surface conductance to water vapor (m s$^{-1}$)
Rn	Net radiation (W m$^{-2}$)
G	Ground heat flux (W m$^{-2}$); optional
S	Sum of all storage fluxes (W m$^{-2}$); optional
missing.G.as.NA	if TRUE, missing G are treated as NAs, otherwise set to 0.
missing.S.as.NA	if TRUE, missing S are treated as NAs, otherwise set to 0.
Esat.formula	Optional: formula to be used for the calculation of esat and the slope of esat. One of "Sonntag_1990" (Default), "Alduchov_1996", or "Allen_1998". See Esat.slope.
constants	cp - specific heat of air for constant pressure (J K$^{-1}$ kg$^{-1}$) eps - ratio of the molecular weight of water vapor to dry air (-) Pa2kPa - conversion pascal (Pa) to kilopascal (kPa)

Details

Total evapotranspiration can be written in the form (Jarvis & McNaughton 1986):

$ET = \Omega ET_eq + (1 - \Omega)ET_imp$

where $\Omega$ is the decoupling coefficient as calculated from decoupling. ET_eq is the equilibrium evapotranspiration rate, the ET rate that would occur under uncoupled conditions, where the heat budget is dominated by radiation (when Ga -> 0):

$ET_eq = (\Delta * (Rn - G - S) * \lambda) / (\Delta + \gamma)$

where $\Delta$ is the slope of the saturation vapor pressure curve (kPa K$^{-1}$), $\lambda$ is the latent heat of vaporization (J kg$^{-1}$), and $\gamma$ is the psychrometric constant (kPa K$^{-1}$). ET_imp is the imposed evapotranspiration rate, the ET rate that would occur under fully coupled conditions (when Ga -> inf):

$ET_imp = (\rho * cp * VPD * Gs * \lambda) / \gamma$

where $\rho$ is the air density (kg m$^{-3}$).

Value

A data.frame with the following columns:

ET_eq	Equilibrium ET (kg m$^{-2}$ s$^{-1}$)
ET_imp	Imposed ET (kg m$^{-2}$ s$^{-1}$)
LE_eq	Equilibrium LE (W m$^{-2}$)
LE_imp	Imposed LE (W m$^{-2}$)

Note

Surface conductance (Gs) can be calculated with surface.conductance. Aerodynamic conductance (Ga) can be calculated using aerodynamic.conductance.

References

Jarvis, P.G., McNaughton, K.G., 1986: Stomatal control of transpiration: scaling up from leaf to region. Advances in Ecological Research 15, 1-49.

Monteith, J.L., Unsworth, M.H., 2008: Principles of Environmental Physics. 3rd edition. Academic Press, London.

See Also

decoupling

Examples

df <- data.frame(Tair=20,pressure=100,VPD=seq(0.5,4,0.5),
                 Gs_ms=seq(0.01,0.002,length.out=8),Rn=seq(50,400,50))            
equilibrium.imposed.ET(df)

---

Saturation Vapor Pressure (Esat) and Slope of the Esat Curve

Description

Calculates saturation vapor pressure (Esat) over water and the corresponding slope of the saturation vapor pressure curve.

Usage

Esat.slope(
  Tair,
  formula = c("Sonntag_1990", "Alduchov_1996", "Allen_1998"),
  constants = bigleaf.constants()
)

Arguments

Tair	Air temperature (degC)
formula	Formula to be used. Either "Sonntag_1990" (Default), "Alduchov_1996", or "Allen_1998".
constants	Pa2kPa - conversion pascal (Pa) to kilopascal (kPa)

Details

Esat (kPa) is calculated using the Magnus equation:

$Esat = a * exp((b * Tair) / (c + Tair)) / 1000$

where the coefficients a, b, c take different values depending on the formula used. The default values are from Sonntag 1990 (a=611.2, b=17.62, c=243.12). This version of the Magnus equation is recommended by the WMO (WMO 2008; p1.4-29). Alternatively, parameter values determined by Alduchov & Eskridge 1996 or Allen et al. 1998 can be used (see references). The slope of the Esat curve ($\Delta$) is calculated as the first derivative of the function:

$\Delta = dEsat / dTair$

which is solved using D.

Value

A dataframe with the following columns:

Esat	Saturation vapor pressure (kPa)
Delta	Slope of the saturation vapor pressure curve (kPa K$^{-1}$)

References

Sonntag D. 1990: Important new values of the physical constants of 1986, vapor pressure formulations based on the ITS-90 and psychrometric formulae. Zeitschrift fuer Meteorologie 70, 340-344.

World Meteorological Organization 2008: Guide to Meteorological Instruments and Methods of Observation (WMO-No.8). World Meteorological Organization, Geneva. 7th Edition.

Alduchov, O. A. & Eskridge, R. E., 1996: Improved Magnus form approximation of saturation vapor pressure. Journal of Applied Meteorology, 35, 601-609

Allen, R.G., Pereira, L.S., Raes, D., Smith, M., 1998: Crop evapotranspiration - Guidelines for computing crop water requirements - FAO irrigation and drainage paper 56, FAO, Rome.

Examples

Esat.slope(seq(0,45,5))[,"Esat"]  # Esat in kPa
Esat.slope(seq(0,45,5))[,"Delta"] # the corresponding slope of the Esat curve (Delta) in kPa K-1

---

Extraterrestrial solar radiation

Description

Compute the extraterrestrial solar radiation with the

Usage

extraterrestrial.radiation(doy, constants = bigleaf.constants())

Arguments

doy	integer vector with day of year (DoY)
constants	solar_constant - solar constant (W m$^{-2}$)

Details

Computation follows Lanini, 2010 (Master thesis, Bern University)

Value

numeric vector of extraterrestrial radiation (W m$^{-2}$)

Examples

plot(1:365, extraterrestrial.radiation(1:365), type = "l"
  , ylab = "radiation (W m-2)", xlab = "day of year")

---

Basic Eddy Covariance Data Filtering

Description

Filters time series of EC data for high-quality values and specified meteorological conditions.

Usage

filter.data(
  data,
  quality.control = TRUE,
  filter.growseas = FALSE,
  filter.precip = FALSE,
  filter.vars = NULL,
  filter.vals.min,
  filter.vals.max,
  NA.as.invalid = TRUE,
  vars.qc = NULL,
  quality.ext = "_qc",
  good.quality = c(0, 1),
  missing.qc.as.bad = TRUE,
  GPP = "GPP",
  doy = "doy",
  year = "year",
  tGPP = 0.5,
  ws = 15,
  min.int = 5,
  precip = "precip",
  tprecip = 0.01,
  precip.hours = 24,
  records.per.hour = 2,
  filtered.data.to.NA = TRUE,
  constants = bigleaf.constants()
)

Arguments

data	Data.frame or matrix containing all required input variables in half-hourly or hourly resolution. Including year, month, day information
quality.control	Should quality control be applied? Defaults to TRUE.
filter.growseas	Should data be filtered for growing season? Defaults to FALSE.
filter.precip	Should precipitation filtering be applied? Defaults to FALSE.
filter.vars	Additional variables to be filtered. Vector of type character.
filter.vals.min	Minimum values of the variables to be filtered. Numeric vector of the same length than filter.vars. Set to NA to be ignored.
filter.vals.max	Maximum values of the variables to be filtered. Numeric vector of the same length than filter.vars. Set to NA to be ignored.
NA.as.invalid	If TRUE (the default) missing data are filtered out (applies to all variables).
vars.qc	Character vector indicating the variables for which quality filter should be applied. Ignored if quality.control = FALSE.
quality.ext	The extension to the variables' names that marks them as quality control variables. Ignored if quality.control = FALSE.
good.quality	Which values indicate good quality (i.e. not to be filtered) in the quality control (qc) variables? Ignored if quality.control = FALSE.
missing.qc.as.bad	If quality control variable is NA, should the corresponding data point be treated as bad quality? Defaults to TRUE. Ignored if quality.control = FALSE.
GPP	Gross primary productivity (umol m$^{-2}$ s$^{-1}$); Ignored if filter.growseas = FALSE.
doy	Day of year; Ignored if filter.growseas = FALSE.
year	Year; Ignored if filter.growseas = FALSE.
tGPP	GPP threshold (fraction of 95th percentile of the GPP time series). Must be between 0 and 1. Ignored if filter.growseas is FALSE.
ws	Window size used for GPP time series smoothing. Ignored if filter.growseas = FALSE.
min.int	Minimum time interval in days for a given state of growing season. Ignored if filter.growseas = FALSE.
precip	Precipitation (mm time$^{-1}$)
tprecip	Precipitation threshold used to identify a precipitation event (mm). Ignored if filter.precip = FALSE.
precip.hours	Number of hours removed following a precipitation event (h). Ignored if filter.precip = FALSE.
records.per.hour	Number of observations per hour. I.e. 2 for half-hourly data.
filtered.data.to.NA	Logical. If TRUE (the default), all variables in the input data.frame/matrix are set to NA for the time step where ANY of the filter.vars were beyond their acceptable range (as determined by filter.vals.min and filter.vals.max). If FALSE, values are not filtered, and an additional column 'valid' is added to the data.frame/matrix, indicating if any value of a row did (1) or did not fulfill the filter criteria (0).
constants	frac2percent - conversion between fraction and percent

Details

This routine consists of two parts:

1) Quality control: All variables included in vars.qc are filtered for good quality data. For these variables a corresponding quality variable with the same name as the variable plus the extension as specified in quality.ext must be provided. For time steps where the value of the quality indicator is not included in the argument good.quality, i.e. the quality is not considered as 'good', its value is set to NA.

2) Meteorological filtering. Under certain conditions (e.g. low ustar), the assumptions of the EC method are not fulfilled. Further, some data analysis require certain meteorological conditions, such as periods without rainfall, or active vegetation (growing season, daytime). The filter applied in this second step serves to exclude time periods that do not fulfill the criteria specified in the arguments. More specifically, time periods where one of the variables is higher or lower than the specified thresholds (filter.vals.min and filter.vals.max) are set to NA for all variables. If a threshold is set to NA, it will be ignored.

Value

If filtered.data.to.NA = TRUE (default), the input data.frame/matrix with observations which did not fulfill the filter criteria set to NA. If filtered.data.to.NA = FALSE, the input data.frame/matrix with an additional column "valid", which indicates whether all the data of a time step fulfill the filtering criteria (1) or not (0).

Note

The thresholds set with filter.vals.min and filter.vals.max filter all data that are smaller than ("<"), or greater than (">") the specified thresholds. That means if a variable has exactly the same value as the threshold, it will not be filtered. Likewise, tprecip filters all data that are greater than tprecip.

Variables considered of bad quality (as specified by the corresponding quality control variables) will be set to NA by this routine. Data that do not fulfill the filtering criteria are set to NA if filtered.data.to.NA = TRUE. Note that with this option *all* variables of the same time step are set to NA. Alternatively, if filtered.data.to.NA = FALSE data are not set to NA, and a new column "valid" is added to the data.frame/matrix, indicating if any value of a row did (1) or did not fulfill the filter criteria (0).

Examples

# Example of data filtering; data are for a month within the growing season,
# hence growing season is not filtered.
# If filtered.data.to.NA=TRUE, all values of a row are set to NA if one filter
# variable is beyond its bounds. 
DE_Tha_Jun_2014_2 <- filter.data(DE_Tha_Jun_2014,quality.control=FALSE,
                                 vars.qc=c("Tair","precip","H","LE"),
                                 filter.growseas=FALSE,filter.precip=TRUE,
                                 filter.vars=c("Tair","PPFD","ustar"),
                                 filter.vals.min=c(5,200,0.2),
                                 filter.vals.max=c(NA,NA,NA),NA.as.invalid=TRUE,
                                 quality.ext="_qc",good.quality=c(0,1),
                                 missing.qc.as.bad=TRUE,GPP="GPP",doy="doy",
                                 year="year",tGPP=0.5,ws=15,min.int=5,precip="precip",
                                 tprecip=0.1,precip.hours=24,records.per.hour=2,
                                 filtered.data.to.NA=TRUE)

 ## same, but with filtered.data.to.NA=FALSE
 DE_Tha_Jun_2014_3 <- filter.data(DE_Tha_Jun_2014,quality.control=FALSE,
                                 vars.qc=c("Tair","precip","H","LE"),
                                 filter.growseas=FALSE,filter.precip=TRUE,
                                 filter.vars=c("Tair","PPFD","ustar"),
                                 filter.vals.min=c(5,200,0.2),
                                 filter.vals.max=c(NA,NA,NA),NA.as.invalid=TRUE,
                                 quality.ext="_qc",good.quality=c(0,1),
                                 missing.qc.as.bad=TRUE,GPP="GPP",doy="doy",
                                 year="year",tGPP=0.5,ws=15,min.int=5,precip="precip",
                                 tprecip=0.1,precip.hours=24,records.per.hour=2,
                                 filtered.data.to.NA=FALSE)
                                 
 # note the additional column 'valid' in DE_Tha_Jun_2014_3.
 # To remove time steps marked as filtered out (i.e. 0 values in column 'valid'):
 DE_Tha_Jun_2014_3[DE_Tha_Jun_2014_3["valid"] == 0,] <- NA

---

GPP-based Growing Season Filter

Description

Filters annual time series for growing season based on smoothed daily GPP data.

Usage

filter.growing.season(GPPd, tGPP, ws = 15, min.int = 5)

Arguments

GPPd	daily GPP (any unit)
tGPP	GPP threshold (fraction of 95th percentile of the GPP time series). Takes values between 0 and 1.
ws	window size used for GPP time series smoothing
min.int	minimum time interval in days for a given state of growing season

Details

The basic idea behind the growing season filter is that vegetation is considered to be active when its carbon uptake (GPP) is above a specified threshold, which is defined relative to the peak GPP (95th percentile) observed in the year. The GPP-threshold is calculated as:

$GPP_threshold = quantile(GPPd,0.95)*tGPP$

GPPd time series are smoothed with a moving average to avoid fluctuations in the delineation of the growing season. The window size defaults to 15 days, but depending on the ecosystem, other values can be appropriate.

The argument min.int serves to avoid short fluctuations in the status growing season vs. no growing season by defining a minimum length of the status. If a time interval shorter than min.int is labeled as growing season or non-growing season, it is changed to the status of the neighboring values.

Value

a vector of type integer of the same length as the input GPPd in which 0 indicate no growing season (dormant season) and 1 indicate growing season.

---

Eddy Covariance Data of FR-Pue (Puechabon)

Description

Halfhourly eddy covariance Data of the site FR-Pue, a Mediterranean evergreen oak forest in Southern France (https://sites.fluxdata.org/FR-Pue/). Data are from May 2012.

Usage

FR_Pue_May_2012

Format

A data frame with 1488 observations and 29 columns:

year

year of measurement

month

month of measurement

doy

day of year

hour

hour (0 - 23.5)

Tair

Air temperature (degC) [TA_F]

Tair_qc

Quality control of Tair [TA_F_QC]

PPFD

Photosynthetic photon flux density (umol m$^{-2}$ s$^{-1}$) [PPFD_IN]

PPFD_qc

Quality control of PPFD [PPFD_IN_QC]

VPD

Vapor pressure deficit (kPa) [VPD_F]

VPD_qc

Quality control of VPD [VPD_F_QC]

pressure

Atmospheric pressure (kPa) [PA_F]

precip

precipitation (mm) [P_F]

precip_qc

Quality control of precip [P_F_QC]

ustar

friction velocity (m s$^{-1}$) [USTAR]

wind

horizontal wind velocity (m s$^{-1}$) [WS_F]

wind_qc

Quality control of wind [WS_F_QC]

Ca

Atmospheric CO$_{2}$ concentration (ppm) [CO2_F_MDS]

Ca_qc

Quality control of Ca [CO2_F_MDS_QC]

LW_up

upward longwave radiation (W m$^{-2}$) [LW_OUT]

Rn

Net radiation (W m$^{-2}$) [NETRAD]

LE

Latent heat flux (W m$^{-2}$) [LE_F_MDS]

LE_qc

Quality control of LE [LE_F_MDS_QC]

H

Sensible heat flux (W m$^{-2}$) [H_F_MDS]

H_qc

Quality control of H [H_F_MDS_QC]

NEE

Net ecosystem exchange (umol m$^{-2}$ s$^{-1}$) [NEE_VUT_USTAR50]

NEE_qc

Quality control of NEE [NEE_VUT_USTAR50_QC]

GPP

Gross primary productivity from nighttime partitioning (umol m$^{-2}$ s$^{-1}$) [GPP_NT_VUT_USTAR50]

GPP_qc

Quality control of GPP [NEE_VUT_USTAR50_QC]

Reco

Ecosystem respiration from nighttime partitioning (umol m$^{-2}$ s$^{-1}$) [RECO_NT_VUT_USTAR50]

Note

The original variable names as provided by the FLUXNET2015 dataset are given in squared brackets. Note that variable units have been converted in some cases (e.g. VPD from hPa to kPa).

Source

original data were downloaded from https://fluxnet.org/ (accessed 09 November 2016)

---

Boundary Layer Conductance according to Choudhury & Monteith 1988

Description

A formulation for the canopy boundary layer conductance for heat transfer according to Choudhury & Monteith 1988.

Usage

Gb.Choudhury(
  data,
  Tair = "Tair",
  pressure = "pressure",
  wind = "wind",
  ustar = "ustar",
  H = "H",
  leafwidth,
  LAI,
  zh,
  zr,
  d,
  z0m = NULL,
  stab_formulation = c("Dyer_1970", "Businger_1971"),
  Sc = NULL,
  Sc_name = NULL,
  constants = bigleaf.constants()
)

Arguments

data	Data.frame or matrix containing all required variables
Tair	Air temperature (degC)
pressure	Atmospheric pressure (kPa)
wind	Wind speed at sensor height (m s$^{-1}$)
ustar	Friction velocity (m s$^{-1}$)
H	Sensible heat flux (W m$^{-2}$)
leafwidth	Leaf width (m)
LAI	One-sided leaf area index
zh	Canopy height (m)
zr	Instrument (reference) height (m)
d	Zero-plane displacement height (-), can be calculated using roughness.parameters
z0m	Roughness length for momentum (m). If not provided, calculated from roughness.parameters within wind.profile
stab_formulation	Stability correction function used (If stab_correction = TRUE). Either "Dyer_1970" or "Businger_1971".
Sc	Optional: Schmidt number of additional quantities to be calculated
Sc_name	Optional: Name of the additonal quantities, has to be of same length than Sc_name
constants	k - von-Karman constant Sc_CO2 - Schmidt number for CO$_2$ Pr - Prandtl number (if Sc is provided)

Details

Boundary layer conductance according to Choudhury & Monteith 1988 is given by:

$Gb_h = LAI((2a/\alpha)*sqrt(u(h)/w)*(1-exp(-\alpha/2)))$

where u(zh) is the wind speed at the canopy surface, approximated from measured wind speed at sensor height zr and a wind extinction coefficient $\alpha$:

$u(zh) = u(zr) / (exp(\alpha(zr/zh -1)))$

.

$\alpha$ is modeled as an empirical relation to LAI (McNaughton & van den Hurk 1995):

$\alpha = 4.39 - 3.97*exp(-0.258*LAI)$

Gb (=1/Rb) for water vapor and heat are assumed to be equal in this package. Gb for other quantities x is calculated as (Hicks et al. 1987):

$Gb_x = Gb / (Sc_x / Pr)^0.67$

where Sc_x is the Schmidt number of quantity x, and Pr is the Prandtl number (0.71).

Value

A data frame with the following columns:

Gb_h	Boundary layer conductance for heat transfer (m s$^{-1}$)
Rb_h	Boundary layer resistance for heat transfer (s m$^{-1}$)
kB_h	kB$^{-1}$ parameter for heat transfer
Gb_Sc_name	Boundary layer conductance for Sc_name (m s$^{-1}$). Only added if Sc_name and Sc_name are provided

Note

If the roughness length for momentum (z0m) is not provided as input, it is estimated from the function roughness.parameters within wind.profile. This function estimates a single z0m value for the entire time period! If a varying z0m value (e.g. across seasons or years) is required, z0m should be provided as input argument.

References

Choudhury, B. J., Monteith J.L., 1988: A four-layer model for the heat budget of homogeneous land surfaces. Q. J. R. Meteorol. Soc. 114, 373-398.

McNaughton, K. G., Van den Hurk, B.J.J.M., 1995: A 'Lagrangian' revision of the resistors in the two-layer model for calculating the energy budget of a plant canopy. Boundary-Layer Meteorology 74, 261-288.

Hicks, B.B., Baldocchi, D.D., Meyers, T.P., Hosker, J.R., Matt, D.R., 1987: A preliminary multiple resistance routine for deriving dry deposition velocities from measured quantities. Water, Air, and Soil Pollution 36, 311-330.

See Also

Gb.Thom, Gb.Su, aerodynamic.conductance

Examples

## bulk canopy boundary layer resistance for a closed canopy (LAI=5) 
## with large leaves (leafwdith=0.1)            
df <- data.frame(Tair=25,pressure=100,wind=c(3,4,5),ustar=c(0.5,0.6,0.65),H=c(200,230,250))    
Gb.Choudhury(data=df,leafwidth=0.1,LAI=5,zh=25,d=17.5,zr=40)

## same conditions, but smaller leaves (leafwidth=0.01)
Gb.Choudhury(data=df,leafwidth=0.01,LAI=5,zh=25,d=17.5,zr=40)

---

Boundary Layer Conductance according to Su et al. 2001

Description

A physically based formulation for the canopy boundary layer conductance to heat transfer according to Su et al. 2001.

Usage

Gb.Su(
  data,
  Tair = "Tair",
  pressure = "pressure",
  ustar = "ustar",
  wind = "wind",
  H = "H",
  zh,
  zr,
  d,
  z0m = NULL,
  Dl,
  fc = NULL,
  LAI = NULL,
  N = 2,
  Cd = 0.2,
  hs = 0.01,
  stab_formulation = c("Dyer_1970", "Businger_1971"),
  Sc = NULL,
  Sc_name = NULL,
  constants = bigleaf.constants()
)

Arguments

data	Data.frame or matrix containing all required variables
Tair	Air temperature (degC)
pressure	Atmospheric pressure (kPa)
ustar	Friction velocity (m s$^{-1}$)
wind	Wind speed (m s$^{-1}$)
H	Sensible heat flux (W m$^{-2}$)
zh	Canopy height (m)
zr	Reference height (m)
d	Zero-plane displacement height (-), can be calculated using roughness.parameters
z0m	Roughness length for momentum (m). If not provided, calculated from roughness.parameters within wind.profile
Dl	Leaf characteristic dimension (m)
fc	Fractional vegetation cover [0-1] (if not provided, calculated from LAI)
LAI	One-sided leaf area index (-)
N	Number of leaf sides participating in heat exchange (defaults to 2)
Cd	Foliage drag coefficient (-)
hs	Roughness height of the soil (m)
stab_formulation	Stability correction function used (If stab_correction = TRUE). Either "Dyer_1970" or "Businger_1971".
Sc	Optional: Schmidt number of additional quantities to be calculated
Sc_name	Optional: Name of the additional quantities, has to be of same length than Sc_name
constants	Kelvin - conversion degree Celsius to Kelvin pressure0 - reference atmospheric pressure at sea level (Pa) Tair0 - reference air temperature (K) Sc_CO2 - Schmidt number for CO$_2$ Pr - Prandtl number (if Sc is provided)

Details

The formulation is based on the kB$^{-1}$ model developed by Massman 1999. Su et al. 2001 derived the following approximation:

$kB^{-1} = (k Cd fc^2) / (4Ct ustar/u(zh)) + kBs-1(1 - fc)^2$

If fc (fractional vegetation cover) is missing, it is estimated from LAI:

$fc = 1 - exp(-LAI/2)$

The wind speed at the top of the canopy is calculated using function wind.profile.

Ct is the heat transfer coefficient of the leaf (Massman 1999):

$Ct = Pr^{-2/3} Reh^-1/2 N$

where Pr is the Prandtl number (set to 0.71), and Reh is the Reynolds number for leaves:

$Reh = Dl wind(zh) / v$

kBs$^{-1}$, the kB$^{-1}$ value for bare soil surface, is calculated according to Su et al. 2001:

$kBs^{-1} = 2.46(Re)^0.25 - ln(7.4)$

Gb (=1/Rb) for water vapor and heat are assumed to be equal in this package. Gb for other quantities x is calculated as (Hicks et al. 1987):

$Gb_x = Gb / (Sc_x / Pr)^{0.67}$

where Sc_x is the Schmidt number of quantity x, and Pr is the Prandtl number (0.71).

Value

A data.frame with the following columns:

Gb_h	Boundary layer conductance for heat transfer (m s$^{-1}$)
Rb_h	Boundary layer resistance for heat transfer (s m$^{-1}$)
kB_h	kB$^{-1}$ parameter for heat transfer
Gb_Sc_name	Boundary layer conductance for Sc_name (m s$^{-1}$). Only added if Sc_name and Sc_name are provided

Note

If the roughness length for momentum (z0m) is not provided as input, it is estimated from the function roughness.parameters within wind.profile. This function estimates a single z0m value for the entire time period! If a varying z0m value (e.g. across seasons or years) is required, z0m should be provided as input argument.

References

Su, Z., Schmugge, T., Kustas, W. & Massman, W., 2001: An evaluation of two models for estimation of the roughness height for heat transfer between the land surface and the atmosphere. Journal of Applied Meteorology 40, 1933-1951.

Massman, W., 1999: A model study of kB H- 1 for vegetated surfaces using 'localized near-field' Lagrangian theory. Journal of Hydrology 223, 27-43.

Hicks, B.B., Baldocchi, D.D., Meyers, T.P., Hosker, J.R., Matt, D.R., 1987: A preliminary multiple resistance routine for deriving dry deposition velocities from measured quantities. Water, Air, and Soil Pollution 36, 311-330.

See Also

Gb.Thom, Gb.Choudhury, aerodynamic.conductance

Examples

# Canopy boundary layer resistance (and kB-1 parameter) for a set of meteorological conditions,
# a leaf characteristic dimension of 1cm, and an LAI of 5
df <- data.frame(Tair=25,pressure=100,wind=c(3,4,5),ustar=c(0.5,0.6,0.65),H=c(200,230,250)) 
Gb.Su(data=df,zh=25,zr=40,d=17.5,Dl=0.01,LAI=5)

# the same meteorological conditions, but larger leaves
Gb.Su(data=df,zh=25,zr=40,d=17.5,Dl=0.1,LAI=5)

# same conditions, large leaves, and sparse canopy cover (LAI = 1.5)
Gb.Su(data=df,zh=25,zr=40,d=17.5,Dl=0.1,LAI=1.5)

---

Boundary Layer Conductance according to Thom 1972

Description

An empirical formulation for the canopy boundary layer conductance for heat transfer based on a simple ustar dependency.

Usage

Gb.Thom(ustar, Sc = NULL, Sc_name = NULL, constants = bigleaf.constants())

Arguments

ustar	Friction velocity (m s$^{-1}$)
Sc	Optional: Schmidt number of additional quantities to be calculated
Sc_name	Optional: Name of the additional quantities, has to be of same length than Sc_name
constants	k - von-Karman constant Sc_CO2 - Schmidt number for CO$_{2}$ Pr - Prandtl number (if Sc is provided)

Details

The empirical equation for Rb suggested by Thom 1972 is:

$Rb = 6.2ustar^-0.67$

Gb (=1/Rb) for water vapor and heat are assumed to be equal in this package. Gb for other quantities x is calculated as (Hicks et al. 1987):

$Gb_x = Gb / (Sc_x / Pr)^0.67$

where Sc_x is the Schmidt number of quantity x, and Pr is the Prandtl number (0.71).

Value

a data.frame with the following columns:

Gb_h	Boundary layer conductance for heat transfer (m s$^{-1}$)
Rb_h	Boundary layer resistance for heat transfer (s m$^{-1}$)
kB_h	kB$^{-1}$ parameter for heat transfer
Gb_Sc_name	Boundary layer conductance for Sc_name (m s$^{-1}$). Only added if Sc_name and Sc_name are provided

References

Thom, A., 1972: Momentum, mass and heat exchange of vegetation. Quarterly Journal of the Royal Meteorological Society 98, 124-134.

Hicks, B.B., Baldocchi, D.D., Meyers, T.P., Hosker, J.R., Matt, D.R., 1987: A preliminary multiple resistance routine for deriving dry deposition velocities from measured quantities. Water, Air, and Soil Pollution 36, 311-330.

See Also

Gb.Choudhury, Gb.Su, aerodynamic.conductance

Examples

Gb.Thom(seq(0.1,1.4,0.1))

## calculate Gb for SO2 as well
Gb.Thom(seq(0.1,1.4,0.1),Sc=1.25,Sc_name="SO2")

---

Bulk Intercellular CO$_{2}$ Concentration

Description

Bulk canopy intercellular CO$_{2}$ concentration (Ci) calculated based on Fick's law given surface conductance (Gs), gross primary productivity (GPP) and atmospheric CO$_{2}$ concentration (Ca).

Usage

intercellular.CO2(
  data,
  Ca = "Ca",
  GPP = "GPP",
  Gs = "Gs_mol",
  Rleaf = NULL,
  missing.Rleaf.as.NA = FALSE,
  constants = bigleaf.constants()
)

Arguments

data	Data.Frame or matrix with all required columns
Ca	Atmospheric or surface CO$_{2}$ concentration (umol mol$^{-1}$)
GPP	Gross primary productivity (umol CO$_{2}$ m$^{-2}$ s$^{-1}$)
Gs	Surface conductance to water vapor (mol m$^{-2}$ s$^{-1}$)
Rleaf	Ecosystem respiration stemming from leaves (umol CO$_{2}$ m$^{-2}$ s$^{-1}$); defaults to 0
missing.Rleaf.as.NA	if Rleaf is provided, should missing values be treated as NA (TRUE) or set to 0 (FALSE, the default)?
constants	DwDc - Ratio of the molecular diffusivities for water vapor and CO$_{2}$ (-)

Details

Bulk intercellular CO$_{2}$ concentration (Ci) is given by:

$Ci = Ca - (GPP - Rleaf)/(Gs/1.6)$

where Gs/1.6 (mol m$^{-2}$ s$^{-1}$) represents the surface conductance to CO$_{2}$. Note that Gs is required in mol m$^{-2}$ s$^{-1}$ for water vapor. Gs is converted to its value for CO$_{2}$ internally. Ca can either be atmospheric CO$_{2}$ concentration (as measured), or surface CO$_{2}$ concentration as calculated from surface.CO2.

Value

Ci -	Bulk canopy intercellular CO$_{2}$ concentration (umol mol$^{-1}$)

Note

The equation is based on Fick's law of diffusion and is equivalent to the often used equation at leaf level (ci = ca - An/gs). Note that GPP and Gs have a different interpretation than An and gs. Gs comprises non-physiological contributions (i.e. physical evaporation) and is confounded by physical factors (e.g. energy balance non-closure). GPP does not account for dark respiration and is further subject to uncertainties in the NEE partitioning algorithm used. Leaf respiration can be provided, but it is usually not known at ecosystem level (as a consequence, Ci is likely to be slightly underestimated) This function should be used with care and the resulting Ci might not be readily comparable to its leaf-level analogue and/or physiological meaningful.

References

Kosugi Y. et al., 2013: Determination of the gas exchange phenology in an evergreen coniferous forest from 7 years of eddy covariance flux data using an extended big-leaf analysis. Ecol Res 28, 373-385.

Keenan T., Sabate S., Gracia C., 2010: The importance of mesophyll conductance in regulating forest ecosystem productivity during drought periods. Global Change Biology 16, 1019-1034.

Examples

# calculate bulk canopy Ci of a productive ecosystem
intercellular.CO2(Ca=400,GPP=40,Gs=0.7)
 
# note the sign convention for NEE

---

Isothermal Net Radiation

Description

Calculates the isothermal net radiation, i.e. the net radiation that the surface would receive if it had the same temperature than the air.

Usage

isothermal.Rn(
  data,
  Rn = "Rn",
  Tair = "Tair",
  Tsurf = "Tsurf",
  emissivity,
  constants = bigleaf.constants()
)

Arguments

data	Data.frame or matrix containing all required variables
Rn	Net radiation (W m$^{-2}$)
Tair	Air temperature (degC)
Tsurf	Surface temperature (degC)
emissivity	Emissivity of the surface (-)
constants	sigma - Stefan-Boltzmann constant (W m$^{-2}$ K$^{-4}$) Kelvin - conversion degree Celsius to Kelvin

Details

The isothermal net radiation (Rni) is given by:

$Rni = Rn + \epsilon * \sigma * (Tsurf^4 - Tair^4)$

where $\epsilon$ is the emissivity of the surface. Tsurf and Tair are in Kelvin.

Value

Rni -

isothermal net radiation (W m$^{-2}$)

References

Jones, H. 2014: Plants and Microclimate. 3rd edition, Cambridge University Press.

Examples

# calculate isothermal net radiation of a surface that is 2degC warmer than the air.
isothermal.Rn(Rn=400,Tair=25,Tsurf=27,emissivity=0.98)

---

Conversion between Mass and Molar Units

Description

Converts mass units of a substance to the corresponding molar units and vice versa.

Usage

kg.to.mol(mass, molarMass = bigleaf.constants()$H2Omol)

Arguments

mass	Numeric vector of mass in kg
molarMass	Numeric vector of molar mass of the substance (kg mol$^{-1}$) e.g. as provided by bigleaf.constants()$H2Omol Default is molar mass of Water.

Value

Numeric vector of amount of substance in mol.

---

Kinematic Viscosity of Air

Description

calculates the kinematic viscosity of air.

Usage

kinematic.viscosity(Tair, pressure, constants = bigleaf.constants())

Arguments

Tair	Air temperature (degC)
pressure	Atmospheric pressure (kPa)
constants	Kelvin - conversion degree Celsius to Kelvin pressure0 - reference atmospheric pressure at sea level (Pa) Tair0 - reference air temperature (K) kPa2Pa - conversion kilopascal (kPa) to pascal (Pa)

Details

where v is the kinematic viscosity of the air (m$^2$ s$^{-1}$), given by (Massman 1999b):

$v = 1.327 * 10^{-5}(pressure0/pressure)(Tair/Tair0)^{1.81}$

Value

v -	kinematic viscosity of air (m$^2$ s$^{-1}$)

References

Massman, W.J., 1999b: Molecular diffusivities of Hg vapor in air, O2 and N2 near STP and the kinematic viscosity and thermal diffusivity of air near STP. Atmospheric Environment 33, 453-457.

Examples

kinematic.viscosity(25,100)

---

Latent Heat of Vaporization

Description

Latent heat of vaporization as a function of air temperature.

Usage

latent.heat.vaporization(Tair)

Arguments

Tair	Air temperature (degC)

Details

The following formula is used:

$\lambda = (2.501 - 0.00237*Tair)10^6$

Value

$\lambda$ - Latent heat of vaporization (J kg$^{-1}$)

References

Stull, B., 1988: An Introduction to Boundary Layer Meteorology (p.641) Kluwer Academic Publishers, Dordrecht, Netherlands

Foken, T, 2008: Micrometeorology. Springer, Berlin, Germany.

Examples

latent.heat.vaporization(seq(5,45,5))

---