This vignette is a short introduction to the functionalities of the Bigleaf.jl package. It is directed to first-time package users who are familiar with the basic concepts of Julia. After presenting the use of several key functions of the package, some useful hints and guidelines are given at the end of the vignette.

Package scope and important conceptual considerations

Bigleaf.jl calculates physical and physiological ecosystem properties from eddy covariance data. Examples for such properties are aerodynamic and surface conductance, surface conditions(e.g. temperature, VPD), wind profile, roughness parameters, vegetation-atmosphere decoupling, potential evapotranspiration, (intrinsic) water-use efficiency, stomatal sensitivity to VPD, or intercellular CO2 concentration. All calculations in the Bigleaf.jl package assume that the ecosystem behaves like a "big-leaf", i.e. a single, homogenous plane which acts as the only source and sink of the measured fluxes. This assumption comes with the advantages that calculations are simplified considerably and that (in most cases) little ancillary information on the EC sites is required. It is important to keep in mind that these simplifications go hand in hand with critical limitations. All derived variables are bulk ecosystem characteristics and have to be interpreted as such. It is for example not possible to infer within-canopy variations of a certain property.

Please also keep in mind that the Bigleaf.jl package does NOT provide formulations for bottom-up modelling. The principle applied here is to use an inversion approach in which ecosystem properties are inferred top-down from the measured fluxes. Such an inversion can, in principle, be also be conducted with more complex models (e.g. sun/shade or canopy/soil models), but keep in mind that these approaches also require that the additional, site-specific parameters are adequately well known.

The use of more detailed models is not within the scope of the Bigleaf.jl package, but it is preferable to use such approaches when important assumptions of the "big-leaf" approach are not met. This is the case in particular when the ecosystem is sparsely covered with vegetation (low LAI, e.g. sparse crops, some savanna systems).

Preparing the data

In this tutorial, we will work with a dataset from the eddy covariance site Tharandt (DE-Tha), a spruce forest in Eastern Germany. The DataFrame DE_Tha_Jun_2014 is downloaded from the bigleaf R package repository and contains half-hourly data of meteorological and flux measurements made in June 2014. For loading the RData into Julia, see the source of this file. We give the data.frame a shorter name here and create a timestamp.

using Bigleaf
using DataFrames

tha = DE_Tha_Jun_2014
set_datetime_ydh!(tha)

And the first six rows of tha:

More information on the data (e.g. meaning of column names and units) can be found at the bigleaf R package. For more information on the site see e.g. Grünwald & Bernhofer 2007. In addition, we will need some ancillary data for this site throughout this tutorial. To ensure consistency, we define them here at the beginning:

thal = (
     LAI = 7.6,   # leaf area index
     zh  = 26.5,  # average vegetation height (m)
     zr  = 42,    # sensor height (m)
     Dl  = 0.01,  # leaf characteristic dimension (m)
)

General guidelines on package usage

There are a few general guidelines that are important to consider when using the Bigleaf.jl package.

Units

It is imperative that variables are provided in the right units, as the plausibility of the input units is not checked in most cases. The required units of the input arguments can be found in the respective help file of the function. The good news is that units do not change across functions. For example, pressure is always required in kPa, and temperature always in °C.

Function arguments

Bigleaf.jl usually provides functions in two flavours.

• providing all arguments seperately as scalars and output being a single scalar or a NamedTuple
• providing a DataFrame as first argument with columns corresponding to the inputs and ouput being the in-place modified DataFrame. Most keyword arguments accept both, vectors or scalars. The column names in the DataFrame should correspond to the argument names of the corresponding method with individual inputs.

We can demonstrate the usage with a simple example:

# explicit inputs
Tair, pressure, Rn, =  14.81, 97.71, 778.17
potential_ET(Val(:PriestleyTaylor), Tair, pressure, Rn)

# DataFrame
potential_ET!(copy(tha), Val(:PriestleyTaylor))

# DataFrame with a few columns overwritten by user values
potential_ET!(transform(tha, :Tair => x -> 25.0; renamecols=false), Val(:PriestleyTaylor))

# varying one input only: scalar form with dot-notation
Tair_vec =  10.0:1.0:20.0
DataFrame(potential_ET.(Val(:PriestleyTaylor), Tair_vec, pressure, Rn))

For functions operating only on vectors, e.g. roughness_parameters vectors are provided with the individual inputs. Sometimes, an additional non-mutating DataFrame variant is provided for convenience, however, the output value of this variant is not type-stable.

Ground heat flux and storage fluxes

Many functions require the available energy ($A$), which is defined as ($A = R_n - G - S$, all in $\text{W m}^{-2}$), where $R_n$ is the net radiation, $G$ is the ground heat flux, and $S$ is the sum of all storage fluxes of the ecosystem (see e.g. Leuning et al. 2012 for an overview). For some sites, $G$ is not available, and for most sites, only a few components of $S$ are measured.

In Bigleaf.jl it is not a problem if $G$ and/or $S$ are missing (other than the results might be (slightly) biased), but special options exist for the treatment of missing $S$ and $G$ values.

Note that the default for G and S in the dataframe variant is missing (and assumed zero), even if those columns are present in the DataFrame. You need to explictly pass those columns with the optional arguments: e.g. potential_ET(df, Val(:PriestleyTaylor); G = df.G)

Note that in difference to the bigleaf R package missing entries in an input vector are not relaced by zero by default. You need to explitly use coalesce when specifying a ground heat flux for which missings should be replaced by zero: ;G = coalesce(df.G, zero(df.G))

Function walkthrough

Data filtering

For most applications it is meaningful to filter your data. There are two main reasons why we want to filter our data before we start calculating ecosystem properties. The first one is to exclude datapoints that do not fulfill the requirements of the EC technique or that are of bad quality due to e.g. instrument failure or gap-filling with poor confidence. Note that the quality assessment of the EC data is not the purpose of the Bigleaf.jl package. This is done by other packages (e.g. REddyProc), which often provide quality control flags for the variables. These quality control flags are used here to filter out bad-quality datapoints.

A second reason for filtering our data is that some derived properties are only meaningful if certain meteorological conditions are met. For instance, if we are interested in properties related to plant gas exchange, it makes most sense to focus on time periods when plants are photosynthetically active (i.e. in the growing season and at daytime).

Bigleaf.jl provides methods that update (or create) the :valid column in a DataFrame. Records, i.e. rows, that contain non valid conditions are set to false. If the valid column was false before, it stays at false.

setinvalid_qualityflag!

One can check quality flags. By default (argument setvalmissing = true) this also replaces the non-valid values in the data-columns by missing.

thaf = copy(tha)   # keep the original tha DataFrame
# if the :valid columns does not exist yet, it is created with all values true
setinvalid_qualityflag!(thaf; vars = ["LE", "NEE"])
sum(.!thaf.valid) # 7 records marked non-valid
sum(ismissing.(thaf.NEE)) # 7 NEE values set to missing

7

In the function call above, vars lists the variables that should be filtered with respect to their quality. Optional parameter qc_suffix="_qc" denotes the extension of the variable name that identifies the column as a quality control indicator of a given variable. The variables LE and LE_qc, for example, denote the variable itself (latent heat flux), and the quality of the variable LE, respectively. The optional argument good_quality_threshold = 1.0 specifies the values of the quality column below which the quality control to be considered as acceptable quality (i.e. to not be filtered). For example with default value 1, all LE values whose LE_qc variable is larger than 1 are set to missing. The variable missing_qc_as_bad is required to decide what to do in case of missing values in the quality control variable. By default this is (conservatively) set to true, i.e. all entries where the qc variable is missing is set invalid.

setinvalid_range!

We can filter for meteorological conditions to be in acceptable ranges. For each variable to check we supply the valid minimum and valid maxium as a two-tuple as the second component of a pair. If their is no limit towards small or large values, supply -Inf or Inf as the minimum or maximum respectively.

setinvalid_range!(thaf,
     :PPFD => (200, Inf),
     :ustar => (0.2, Inf),
     :LE =>(0, Inf),
     :VPD => (0.01, Inf)
     )
sum(.!thaf.valid) # many more records marked invalid
minimum(skipmissing(thaf.PPFD)) >= 200 # values outsides range some set to missing
sum(ismissing.(thaf.PPFD))

697

About half of the data were filtered because radiation was not high enough (night-time). Another quarter was filtered because they showed negative LE values. However, most of them occured during the night:

sum(ismissing.(thaf.PPFD)) / nrow(thaf) # 0.48
sum(.!ismissing.(thaf.PPFD) .&& ismissing.(thaf.LE)) / nrow(thaf) # 0.05

setinvalid_nongrowingseason!

A third method filters periods outside the growing season:

setinvalid_nongrowingseason!(thaf, 0.4)
sum(.!thaf.valid) # tha dataset is all within growing season - no additional invalids

812

This function implements a simple growing season filter based on daily smoothed GPP time series. Arguments tGPP determines how high daily GPP has to be in relation to its peak value within the year. In this case, the value of 0.4 denotes that smoothed GPP has to be at least 40% of the 95th quantile. Argument ws controls the degree of smoothing in the timeseries and should be between 10-20 days. The purpose of which is to minimize the high variation of GPP between days, Argument min_int is a parameter that avoids that data are switching from inside the growing season and out from one day to the next. It determines the minimum number of days that a season should have. The growing season filter is applicable to all sites, with one more more growing seasons, but its advisable that site-specific parameter settings are used.

In this example, it does not really make sense to filter for growing season, since it uses only one month of data of which we know that vegetation is active at the site. The algorithm realizes that and does not mark any additional data as invalid.

setinvalid_afterprecip!

As a last step we will filter for precipitation events. This is often meaningful for ecophysiological studies because data during and shortly after rainfall events do not contain much information on the physiological activity of the vegetation because they comprise significant fractions of evaporation from the soil and plant surfaces. The purpose of such a filter is mostly to minimize the fraction of soil and interception evaporation on the total water flux. This filter simply excludes periods following a precipitation event. A precipitation event, here, is defined as any time step with a recorded precipitation higher than min_precip (in mm per timestep). The function then filters all time periods following a precipitation event. The number of subsequent time periods excluded is controlled by the argument precip_hours. Here, we exclude rainfall events and the following 24 hours. The timestamps in the DataFrame must be sorted in increasing order.

setinvalid_afterprecip!(thaf; min_precip=0.02, hours_after=24)
sum(.!thaf.valid) # some more invalids

1013

In this walkthrough we use the data as filtered above:

thas = subset(thaf, :valid)