TITLE: Calculating Vapour Pressure Deficit and other atmospheric things
DATE: 2018-12-05
AUTHOR: John L. Godlee
====================================================================


I needed to calculate Vapour Pressure Deficit for a productivity
model that I’m running, using ECMWF spatial data timeseries as an
input. I wrote some functions in R that I think do the right thing:

    # Saturation vapour pressure (millibars)
    vps <- function(T_k){
      T_c <- T_k - 273.15
      return(0.6108 * exp(17.27 * T_c / (T_c + 273.15))) 
    }

    # Relative humidity (%)
    rh <- function(T_k, Td_k){
      es_calc <- vps(T_k)  # Saturation vp
      ea_calc <- vps(Td_k)  # Actual vp
      return((ea_calc / es_calc) * 100) 
    }

    # Vapour pressure of air (millibars)
    vpa <- function(T_k, Td_k){
      rh_calc <- rh(T_k, Td_k)
      vps_calc <- vps(T_k)
      return(rh_calc / 100 * vps_calc)
    }

    # Vapour pressure deficit (kPa)
    vpd <- function(T_k, Td_k){ 
      vps_calc <- vps(T_k)
      vpa_calc <- vpa(T_k, Td_k)
      return(vpa_calc - vps_calc)
    }