---
title:  "Introduction to CropWaterBalance"
output: rmarkdown::html_vignette
vignette: >
  %\VignetteIndexEntry{Introduction to CropWaterBalance}
  %\VignetteEngine{knitr::rmarkdown}
  %\VignetteEncoding{UTF-8}
  bibliography: bibliography.bib
---

```{r, include = FALSE}
knitr::opts_chunk$set(
  collapse = TRUE,
  comment = "#>"
)
```
# Introduction

Irrigated areas have increased throughout the globe to support the growing global population and to cope with climate change impacts (@Siyal2023). In this context, the {CropWaterBalance} allows users to keep track of the soil water deficit in the root zone through the crop water balance accounting (@Andales2012). The goal of the package is to assist users in making decisions about when and how much to irrigate.
The most important function of the package is the `CWB()`, which calculates several parameters of the crop water balance, including crop evapotranspiration (ETc), actual crop evapotranspiration, stored water in the root zone and soil water deficit (D). The function also suggests when irrigate, considering the management allowed depletion (MAD) provided by the users.
Although the FAO-Penman and Monteith equation is recognized as the standard method for estimating daily amounts of reference evapotranspiration (ETO) (@Allen1998), the high number of variables required for calculating this model limits its operational use in several regions in the world. Thus, the package includes the function `ETO_PM()`, which estimate daily ETO amounts through the FAO-Penman and Monteith equation, and the functions `ETO_PT()` and `ETO_HS()`), which calculates this agrometeorological parameter through other two alternative (and simpler) methods: Priestley-Taylor (@Priestley1972) and Hargreaves-Samani (@Hargreaves1985), respectively. 
Additionally, the {CropWaterBalance} has other two functions (`Compare()` and `Descriptive()`). The first calculates measures of accuracy and agreement between two data samples and the second calculates descriptive statistics for these samples. Therefore, these two functions assist users in selecting suitable ET0 estimating methods for a particular region or season.

# Getting Started
Load the library in your R session.

```{r setup}
library(CropWaterBalance)
```
## Using ETO_PM() to calculate daily amounts of reference evapotranspiration 

Reference evapotranspiration (ET0) is the combined process of evaporation and transpiration that occurs from well-fertilized and disease-free hypothetical grass reference crop, grown in large fields under no soil water restriction and achieving full production (@Allen1998). 
The `ETO_PM()` function calculates daily ETO amounts using the FAO-Penman and Monteith equation (in millimetres).

## Step 1: Applying `ETO_PM()` in Campinas-SP, Brazil


```{r CpsET0PM}
Tavg <- DataForCWB[, 2]
Tmax <- DataForCWB[, 3]
Tmin <- DataForCWB[, 4]
Rn <- DataForCWB[, 6]
WS <- DataForCWB[, 7]
RH <- DataForCWB[, 8]
G <- DataForCWB[, 9]
CpsET0PM <-
  ET0_PM(
    Tavg = Tavg,
    Tmax = Tmax,
    Tmin = Tmin,
    Rn = Rn,
    RH = RH,
    WS = WS,
    G = G,
    Alt = 700
  )
head(CpsET0PM)
```

By analyzing the `ET0_PM()` function, the user verify that the FAO-Penman and Monteith equation requires several inputs, which includes the soil heat flux (G). This latter variable is rarely measure. This is the reason why the {CropWaterBalance} has an auxiliary function (`Soil_Heat_Flux()`) that estimates G as function of daily average air temperature values. The users may use this auxiliary function to estimate G and then apply the `ET0_PM()` function. Alternatively, the users may simply run `ET0_PM()` without the G argument. In this case, this latter function will automatically use `Soil_Heat_Flux()` and then estimate ETO. See the example below.

```{r CpsET0PM estimating G}
Tavg <- DataForCWB[, 2]
Tmax <- DataForCWB[, 3]
Tmin <- DataForCWB[, 4]
Rn <- DataForCWB[, 6]
WS <- DataForCWB[, 7]
RH <- DataForCWB[, 8]
G <- Soil_Heat_Flux(Tavg)
CpsET0PM_WithG <-
  ET0_PM(
    Tavg = Tavg,
    Tmax = Tmax,
    Tmin = Tmin,
    Rn = Rn,
    RH = RH,
    WS = WS,
    G = G,
    Alt = 700
  )
CpsET0PM_WithoutG <-
  ET0_PM(
    Tavg = Tavg,
    Tmax = Tmax,
    Tmin = Tmin,
    Rn = Rn,
    RH = RH,
    WS = WS,
    Alt = 700
  )
head(cbind(CpsET0PM_WithG, CpsET0PM_WithoutG))
```

## Step 2: Applying `ETO_PT()` and `ET0_HS()` in Campinas-SP, Brazil

As described in the introduction, the users may need to estimate ETO using non-standard methods, which require less input than the FAO-Penman and Monteith. The {CropWaterBalance} allows users to estimate this agrometeorological parameter using the Priestley-Taylor (`ETO_PT()`) and Hargreaves-Samani (`ET0_HS()`) methods. Note that the same considerations made for G in respect to `ETO_PM()` are valid for `ETO_PT()`.

```{r CpsET0PT and CpsET0HS}
Tavg <- DataForCWB[, 2]
Tmax <- DataForCWB[, 3]
Tmin <- DataForCWB[, 4]
Ra <- DataForCWB[, 5]
Rn <- DataForCWB[, 6]
G <- DataForCWB[, 9]
CpsET0PT <- ET0_PT(Tavg = Tavg, Rn = Rn, G = G)
CpsET0HS <- ET0_HS(
  Ra = Ra,
  Tavg = Tavg,
  Tmax = Tmax,
  Tmin = Tmin
)
head(cbind(CpsET0PT, CpsET0HS))
```

The estimation of CpsET0PT and/or CpsET0HS raises the following question: Can these alternative methods really replace the FAO-Penman and Monteith model? Although the answer to this question may be regarded as a complex function involving local weather conditions and users' subjective choices, the (`Compare()`) and (`Descriptive()`) functions provide statistical information, which may assist users in such decision. In this context, the `Compare()` and `Descriptive()` functions may be used to verify how well an alternative ET0 estimating method approaches ET0_PM.

```{r Evaluating alternative methods for ETO}
Tavg <- DataForCWB[, 2]
Tmax <- DataForCWB[, 3]
Tmin <- DataForCWB[, 4]
Ra <- DataForCWB[, 5]
Rn <- DataForCWB[, 6]
WS <- DataForCWB[, 7]
RH <- DataForCWB[, 8]
G <- DataForCWB[, 9]
CpsET0PM <-
  ET0_PM(
    Tavg = Tavg,
    Tmax = Tmax,
    Tmin = Tmin,
    Rn = Rn,
    RH = RH,
    WS = WS,
    G = G,
    Alt = 700
  )
CpsET0PT <- ET0_PT(Tavg = Tavg, Rn = Rn, G = G)
CpsET0HS <- ET0_HS(
  Ra = Ra,
  Tavg = Tavg,
  Tmax = Tmax,
  Tmin = Tmin
)
PM_PT <- Compare(Sample1 = CpsET0PM, Sample2 = CpsET0PT)
PM_PT
Descrp_PM_PT <-
  cbind(Descriptive(Sample = CpsET0PM), Descriptive(Sample = CpsET0PT))
Descrp_PM_PT
PM_HS <- Compare(Sample1 = CpsET0PM, Sample2 = CpsET0HS)
PM_HS
Descrp_PM_HS <-
  cbind(Descriptive(Sample = CpsET0PM), Descriptive(Sample = CpsET0HS))
Descrp_PM_HS
```

The results provided by the `Compare()` and `Descriptive()` functions, indicate that none of the two alternative methods can be used to replace the for calculating daily ET0 amounts. For instance, the values of the modified index of agreement (dmod) (@Willmott1985) remained below 0.4 for both comparisons  (PM vs PT) and (PM vs HS). Additionally, the corresponding absolute mean errors (AME; 1.69222 and 1.468722) represents, approximately, 50% and 43% of the average value of the ET0_PM (3.367377).  

## Step 3: The Crop Water Balance Accounting. Applying `InitialD()` and `CWB()` in Campinas-SP, Brazil

Considering the previous results, we applied the `CWB()` using daily values of rainfall and ET0_PM obtained/estimated from daily data of the weather station of Campinas. We included this meteorological data in the package (DataForCWB), along with parameters required for calculating the crop water balance: depth of the root zone (Drz), available water capacity (AWC; amount of water between field capacity and permanent wilting point), management allowed depletion (MAD), and crop coefficient (Kc). The example below loaded DataForCWB to apply `CWB()` in Campinas-SP. Only for the sake of simplicity, Kc values were set to 1 for the entire period. The initial D value (Dinitial) required for initiating the crop water balance account was obtained  using function `DInitial()` with teta_FC and AWC values obtained from two data sets included in the package (DataForSWC.rda and DataForAWC). These data sets provide soil water content values (m3/m3) for the effective root zone at the field capacity and at permanent wilting point, and AWC values (mm/m), respectively (teta_obs was measured in the field). We applied the `CWB()` function considering two scenarios. Scenario 1 in which no irrigation was applied and scenario 2 in which the package’s recommendation about when and how much to irrigate was met. 

```{r `CWB()` and `DInitial()` in Campinas-SP, Brazil: scenario 1}
Tavg <- DataForCWB[, 2]
Tmax <- DataForCWB[, 3]
Tmin <- DataForCWB[, 4]
Rn <- DataForCWB[, 6]
WS <- DataForCWB[, 7]
RH <- DataForCWB[, 8]
G <- DataForCWB[, 9]
ET0 <- ET0_PM(Tavg, Tmax, Tmin, Rn, RH, WS, G, Alt = 700)
Dinitial <-
  Dinitial(teta_FC = 0.30,
           teta_Obs = 0.17,
           Drz = DataForCWB[1, 11])
Rain <- DataForCWB[, 10]
Drz <- DataForCWB[, 11]
AWC <- DataForCWB[, 12]
MAD <- DataForCWB[, 13]
Kc <- DataForCWB[, 14]
Irrig <- DataForCWB[, 15]
Scenario1 <- CWB(
  Rain = Rain,
  ET0 = ET0,
  AWC = AWC,
  Drz = Drz,
  Kc = Kc,
  Irrig = Irrig,
  MAD = MAD,
  InitialD = Dinitial,
  start.date = "2011-11-23"
)
Scenario1[1:8, ]
Irrig[7] <- 16
Scenario2 <- CWB(
  Rain = Rain,
  ET0 = ET0,
  AWC = AWC,
  Drz = Drz,
  Kc = Kc,
  Irrig = Irrig,
  MAD = MAD,
  InitialD = Dinitial,
  start.date = "2011-11-23"
)
Scenario2[1:8, ]
```

In both scenarios, the package suggested irrigating on 29/11/2011 because the soil water deficit had become larger than dmad. Since in scenario 1 we applied no irrigation, the Ks coefficient started to assume values smaller than 1 and the NonStandardCropEvap becomes smaller than ETc, which may prevent the crop from achieving its potential yield. In scenario 2, we followed the package's recommendation and applied irrigation on day 29 (Irrig = 16 mm). In this case, no crop evapotranspiration deficit was observed, indicating no water shortage in the root zone. 

## Step 3 (Alternative): The Crop Water Balance Accounting. Applying `InitialD()` and `CWB_fixedSchedule` in Campinas-SP, Brazil

Until now, we have assumed that the decision of when and how much to irrigate was solely dependent on soil-plant-atmosphere factors. However, it is well-known that this decision may also be significantly influenced by the design and operation of the irrigation system, as well as the availability of labor and water. This is why we included in the package the `CWB_fixedSchedule` fuction that, as the `CWB()`, also calculates the crop water balance in the root zone. However, it allows users to specify the number of days between consecutive irrigations.This is accompliched by setting the parameter Scheduling to the number of days defining this fixed interval between two consecutive irrigations. The estimation of how much to irrigate is performed on this specific days. This is why we included in the package the `CWB_fixedSchedule()` function that, as the `CWB`, also calculates the crop water balance in the root zone. However, it allows users to specify the number of days between consecutive irrigations.This is accomplished by setting the parameter Scheduling to the number of days defining this fixed interval between two consecutive irrigations. The estimation of how much to irrigate is performed on this specific days. 

```{r `CWB_fixedSchedule()` and `DInitial()` in Campinas-SP, Brazil: scenario 1}
Tavg <- DataForCWB[, 2]
Tmax <- DataForCWB[, 3]
Tmin <- DataForCWB[, 4]
Rn <- DataForCWB[, 6]
WS <- DataForCWB[, 7]
RH <- DataForCWB[, 8]
G <- DataForCWB[, 9]
ET0 <- ET0_PM(Tavg, Tmax, Tmin, Rn, RH, WS, G, Alt = 700)
Dinitial <-
  Dinitial(teta_FC = 0.30,
           teta_Obs = 0.17,
           Drz = DataForCWB[1, 11])
Rain <- DataForCWB[, 10]
Drz <- DataForCWB[, 11]
AWC <- DataForCWB[, 12]
MAD <- DataForCWB[, 13]
Kc <- DataForCWB[, 14]
Scheduling <- 5
Irrig <- DataForCWB[, 15]
Scenario_FixedSchedule <-
  CWB_FixedSchedule(
    Rain = Rain,
    ET0 = ET0,
    AWC = AWC,
    Drz = Drz,
    Kc = Kc,
    Irrig = Irrig,
    MAD = MAD,
    InitialD = Dinitial,
    Scheduling = Scheduling,
    start.date = "2011-11-23"
  )
Scenario_FixedSchedule[1:15, ]
```

# References
