Get monthly climate normals using ClimateR
Ivan Lizarazo
01.06.2022
1. Introduction
This is a simple notebook to illustrate how to obtain “normal” monthly precipitation and temperature values for the period: 1981- 2010 using the climateR library.
“Normal” precipitation does not equal “what you should expect”. “Normal” precipitation is an average of the precipitation values over a 30-year period. Actual precipitation in a given year may be either above or below the seasonal average, or “normal”.
See this notebook for additional information on the climateR library.
Setup
Clean the R memory:
rm(list=ls())Load libraries:
library(AOI)
library(climateR)
library(sf)
library(raster)
library(rasterVis)
library(dplyr)
library(leaflet)
library(leafem)
library(RColorBrewer)
library(stars)Read study area
(stder <- st_read("./68_SANTANDER/STDER_MPIO.shp"))## Reading layer `STDER_MPIO' from data source
## `/Users/ivanlizarazo/Documents/GB/68_SANTANDER/STDER_MPIO.shp'
## using driver `ESRI Shapefile'
## Simple feature collection with 87 features and 9 fields
## Geometry type: POLYGON
## Dimension: XY
## Bounding box: xmin: -74.52895 ymin: 5.707536 xmax: -72.47706 ymax: 8.14501
## Geodetic CRS: WGS 84Get monthly precipitation normals
One year spans twelve months:
oneyear = seq(1,12)Download “normal” monthly precipitation from 1981 to 2010:
# CWD
precip = getTerraClimNormals(stder, param = "prcp", period = "19812010", month=oneyear)## Spherical geometry (s2) switched off## Spherical geometry (s2) switched onLet’s check what we got:
precip## $terraclim_19812010_prcp
## class : RasterStack
## dimensions : 60, 50, 3000, 12 (nrow, ncol, ncell, nlayers)
## resolution : 0.04166667, 0.04166667 (x, y)
## extent : -74.54167, -72.45833, 5.666667, 8.166667 (xmin, xmax, ymin, ymax)
## crs : +proj=longlat +datum=WGS84 +no_defs
## names : X01, X02, X03, X04, X05, X06, X07, X08, X09, X10, X11, X12
## min values : 3.0, 16.5, 37.3, 118.8, 105.9, 95.0, 39.3, 45.9, 61.6, 102.6, 83.3, 11.8
## max values : 122.8, 157.0, 222.4, 355.7, 502.6, 415.7, 450.5, 443.0, 444.3, 525.2, 413.1, 277.3Note that the precip object has twelve layers.
capas <-names(precip)# histogram
par(mfcol=c(12,12), mfrow=c(4,3), oma=c(1,1,0,0),
mar=c(1,1,1,0), tcl=-0.1, mgp=c(0,0,0))
for (penta in capas) {
hist(precip[[penta]], ylab=NA, cex.axis=0.5, font.main=1, cex.main=0.8)
}
Then, calculate the normal precipitation for the whole year:
year_precip <- sum(precip$terraclim_19812010_prcp)What we got?
year_precip## class : RasterLayer
## dimensions : 60, 50, 3000 (nrow, ncol, ncell)
## resolution : 0.04166667, 0.04166667 (x, y)
## extent : -74.54167, -72.45833, 5.666667, 8.166667 (xmin, xmax, ymin, ymax)
## crs : +proj=longlat +datum=WGS84 +no_defs
## source : memory
## names : layer
## values : 943.8, 3601 (min, max)Create a precipitation palette:
ppal <- colorNumeric(c("orange", "yellow", "#a1d99b", "#deebf7", "#9ecae1","#3182bd", "darkblue"), values(year_precip$layer),
na.color = "transparent")Visualize the accumulated precipitation:
leaflet() %>% addTiles() %>%
addRasterImage(year_precip, colors = ppal, opacity = 0.8) %>%
addLegend(pal = ppal, values = values(year_precip$layer),
title = "Annual Precipitation")Get monthly temperature normals
# Maximum temperature
tmax = getTerraClimNormals(stder, param = "tmax", period = "19812010", month=oneyear)## Spherical geometry (s2) switched off## Spherical geometry (s2) switched ontmax## $terraclim_19812010_tmax
## class : RasterStack
## dimensions : 60, 50, 3000, 12 (nrow, ncol, ncell, nlayers)
## resolution : 0.04166667, 0.04166667 (x, y)
## extent : -74.54167, -72.45833, 5.666667, 8.166667 (xmin, xmax, ymin, ymax)
## crs : +proj=longlat +datum=WGS84 +no_defs
## names : X01, X02, X03, X04, X05, X06, X07, X08, X09, X10, X11, X12
## min values : 9.4, 9.5, 9.2, 8.6, 8.2, 7.0, 6.7, 7.3, 8.0, 7.9, 8.2, 8.9
## max values : 34.1, 34.7, 34.7, 33.9, 33.4, 33.4, 34.2, 33.7, 33.7, 32.9, 32.2, 32.8(year_tmax = max(tmax$terraclim_19812010_tmax))## class : RasterLayer
## dimensions : 60, 50, 3000 (nrow, ncol, ncell)
## resolution : 0.04166667, 0.04166667 (x, y)
## extent : -74.54167, -72.45833, 5.666667, 8.166667 (xmin, xmax, ymin, ymax)
## crs : +proj=longlat +datum=WGS84 +no_defs
## source : memory
## names : layer
## values : 9.5, 34.7 (min, max)# Minimum temperature
tmin = getTerraClimNormals(stder, param = "tmin", period = "19812010", month=oneyear)## Spherical geometry (s2) switched off## Spherical geometry (s2) switched ontmin## $terraclim_19812010_tmin
## class : RasterStack
## dimensions : 60, 50, 3000, 12 (nrow, ncol, ncell, nlayers)
## resolution : 0.04166667, 0.04166667 (x, y)
## extent : -74.54167, -72.45833, 5.666667, 8.166667 (xmin, xmax, ymin, ymax)
## crs : +proj=longlat +datum=WGS84 +no_defs
## names : X01, X02, X03, X04, X05, X06, X07, X08, X09, X10, X11, X12
## min values : -1.0, -1.0, 0.3, 1.0, 0.9, 1.0, 0.5, 0.4, 0.5, 0.9, 0.8, 0.1
## max values : 23.9, 24.0, 24.7, 24.5, 24.3, 24.2, 23.9, 24.0, 24.0, 23.9, 24.2, 24.3(year_tmin = min(tmin$terraclim_19812010_tmin))## class : RasterLayer
## dimensions : 60, 50, 3000 (nrow, ncol, ncell)
## resolution : 0.04166667, 0.04166667 (x, y)
## extent : -74.54167, -72.45833, 5.666667, 8.166667 (xmin, xmax, ymin, ymax)
## crs : +proj=longlat +datum=WGS84 +no_defs
## source : memory
## names : layer
## values : -1, 23.4 (min, max)Create a stack layer:
(tempera <- stack(year_tmin, year_tmax))## class : RasterStack
## dimensions : 60, 50, 3000, 2 (nrow, ncol, ncell, nlayers)
## resolution : 0.04166667, 0.04166667 (x, y)
## extent : -74.54167, -72.45833, 5.666667, 8.166667 (xmin, xmax, ymin, ymax)
## crs : +proj=longlat +datum=WGS84 +no_defs
## names : layer.1, layer.2
## min values : -1.0, 9.5
## max values : 23.4, 34.7Change layers’ names:
names(tempera) <- c("tmin", "tmax")Define display values:
#(valores <- seq(from=-1,to=36,by=6))
tvals <- sort(c(values(tempera$tmin),values(tempera$tmax)))Create color palette:
pal <- colorNumeric(rev(c("red", "orange", "#fff9ae",
"#77ab59", "#dcf3ff", "#a2d2df", "#257ca3")),
tvals, na.color = "transparent")Now, it is time for visualization:
capas <- names(tempera)Create the map view:
# Create leaflet widget --------------------------------------------------------
m <- leaflet() %>%
addTiles(group = "OSM (default)") %>%
addRasterImage(year_tmin, colors = pal, opacity = 0.8, group= capas[1]) %>%
addRasterImage(year_tmax, colors = pal, opacity = 0.8, group= capas[2])
#
# Additional leaflet options
m <- m %>%
# Add layers controls
addLayersControl(
overlayGroups = rev(capas),
options = layersControlOptions(collapsed = FALSE)
) %>%
# Add common legend
addLegend(pal = pal, values = tvals,
title = "Annual Temperature")Visualize it:
mReproducibility
Cite as: Lizarazo, I., 2022. Get monthly climate normals using ClimateR. Available at https://rpubs.com/ials2un/climate.
sessionInfo()## R version 4.1.0 (2021-05-18)
## Platform: x86_64-apple-darwin17.0 (64-bit)
## Running under: macOS Catalina 10.15.7
##
## Matrix products: default
## BLAS: /Library/Frameworks/R.framework/Versions/4.1/Resources/lib/libRblas.dylib
## LAPACK: /Library/Frameworks/R.framework/Versions/4.1/Resources/lib/libRlapack.dylib
##
## locale:
## [1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8
##
## attached base packages:
## [1] stats graphics grDevices utils datasets methods base
##
## other attached packages:
## [1] stars_0.5-5 abind_1.4-5 RColorBrewer_1.1-3
## [4] leafem_0.2.0 leaflet_2.0.3.9000 dplyr_1.0.7
## [7] rasterVis_0.50.3 latticeExtra_0.6-29 lattice_0.20-44
## [10] terra_1.3-4 raster_3.4-13 sp_1.4-7
## [13] sf_1.0-7 climateR_0.1.0 AOI_0.2.0.9000
##
## loaded via a namespace (and not attached):
## [1] httr_1.4.2 rnaturalearth_0.1.0 sass_0.4.0
## [4] jsonlite_1.7.2 viridisLite_0.4.0 foreach_1.5.2
## [7] bslib_0.2.5.1 shiny_1.6.0 assertthat_0.2.1
## [10] highr_0.9 yaml_2.2.1 pillar_1.7.0
## [13] glue_1.6.2 digest_0.6.29 promises_1.2.0.1
## [16] rvest_1.0.1 colorspace_2.0-3 htmltools_0.5.1.1
## [19] httpuv_1.6.1 pkgconfig_2.0.3 s2_1.0.7
## [22] purrr_0.3.4 xtable_1.8-4 scales_1.2.0
## [25] jpeg_0.1-9 later_1.2.0 tibble_3.1.7
## [28] proxy_0.4-26 leaflet.extras_1.0.0 farver_2.1.0
## [31] generics_0.1.0 ellipsis_0.3.2 hexbin_1.28.2
## [34] cli_3.3.0 magrittr_2.0.3 crayon_1.5.1
## [37] mime_0.11 evaluate_0.14 fansi_1.0.3
## [40] doParallel_1.0.16 xml2_1.3.2 lwgeom_0.2-8
## [43] class_7.3-19 tools_4.1.0 USAboundaries_0.3.1
## [46] lifecycle_1.0.1 stringr_1.4.0 munsell_0.5.0
## [49] compiler_4.1.0 jquerylib_0.1.4 e1071_1.7-8
## [52] RNetCDF_2.4-2 rlang_1.0.2 classInt_0.4-3
## [55] units_0.7-2 grid_4.1.0 iterators_1.0.14
## [58] rstudioapi_0.13 htmlwidgets_1.5.3 crosstalk_1.1.1
## [61] base64enc_0.1-3 rmarkdown_2.9 wk_0.5.0
## [64] codetools_0.2-18 DBI_1.1.1 R6_2.5.1
## [67] zoo_1.8-9 rgdal_1.5-23 knitr_1.33
## [70] fastmap_1.1.0 utf8_1.2.2 KernSmooth_2.23-20
## [73] stringi_1.7.3 parallel_4.1.0 Rcpp_1.0.7
## [76] vctrs_0.4.1 png_0.1-7 tidyselect_1.1.1
## [79] xfun_0.24