CHESS baseline: compute soil moisture deficit and surplus

Load packages, read data and source custom scripts

rm(list = ls())
library(stars)

#> Loading required package: abind

#> Loading required package: sf

#> Linking to GEOS 3.8.0, GDAL 3.0.4, PROJ 6.3.1

path_proj <- day2day::git_path()
path_data <- file.path(path_proj, "data")
path_cleaned <- file.path(path_data, "cleaned")
path_processed <- file.path(path_data, "processed")

source(file.path(path_proj, "src", "32-soil-moisture.R"))
source(file.path(path_proj, "src", "61-vis-stars.R"))

chess_meta <- readr::read_csv(file.path(path_cleaned, "chess-metadata.csv"))

#> Parsed with column specification:
#> cols(
#>   file = col_character(),
#>   rcp = col_logical(),
#>   variable = col_character(),
#>   summary = col_character(),
#>   temp_resolution = col_character(),
#>   from = col_double(),
#>   to = col_double(),
#>   from_year = col_double(),
#>   from_month = col_character(),
#>   to_year = col_double(),
#>   to_month = col_double(),
#>   description = col_character()
#> )

awc <- read_stars(file.path(path_processed, "ukcp18_uk_1km_soilparams_awc_bc.tif"))
uk <- st_read(file.path(path_cleaned, "uk_simple.gpkg"), "union")

#> Reading layer `union' from data source `/home/rstudio/Documents/Projects/land-suitability/data/cleaned/uk_simple.gpkg' using driver `GPKG'
#> Simple feature collection with 1 feature and 0 fields
#> geometry type:  MULTIPOLYGON
#> dimension:      XY
#> bbox:           xmin: -116.1923 ymin: 7054.099 xmax: 655644.8 ymax: 1218625
#> projected CRS:  OSGB 1936 / British National Grid

Compute soil moisture deficit and surplus

out <- soil_moisture_metadata(
    chess_meta, 1991, 2011, awc, vars = c("precip", "peti"),
    path_in = path_cleaned, path_out = path_processed, return_vars = TRUE
)

Visualize available water-holding capacity (AWC)

plot_awc(out$awc, uk, main = "")

Total monthly precipitation

The data provided is mean monthly precipitation (per second), it need to be transformed to total. It is achieved by multiplying by the number seconds in a month.

plot_prec(out$precip, uk, main = "")

Total monthly potential evapotranspiration

The data provided is mean monthly potential evapotranspiration (per day), it need to be transformed to total. It is achieved by multiplying by the number of days in a month.

plot_pet(out$pet, uk, main = "")

Visualize soil moisture (deficit and surplus)

sm_files <- list.files(path_processed, "^chess-met_sm")
sm <- lapply(sm_files, function(x) read_stars(file.path(path_processed, x)))

plot_sms(sm[[2]], uk, nbreaks = 20, pal_col = "YlGnBu")

plot_smd(sm[[1]], uk, nbreaks = 20, pal_col = "YlOrRd")

Time to execute the task

Only useful when executed with Rscript.

proc.time()

#>    user  system elapsed 
#>  48.579   2.756  51.316