rm(list=ls())
library(AOI)
library(climateR)
library(sf)
## Linking to GEOS 3.8.0, GDAL 3.0.4, PROJ 6.3.1; sf_use_s2() is TRUE
library(raster)
## Loading required package: sp
library(rasterVis)
## Loading required package: lattice
library(dplyr)
##
## Attaching package: 'dplyr'
## The following objects are masked from 'package:raster':
##
## intersect, select, union
## The following objects are masked from 'package:stats':
##
## filter, lag
## The following objects are masked from 'package:base':
##
## intersect, setdiff, setequal, union
library(sf)
(narino <- sf::st_read("/cloud/project/MunicipiosNarino.shp"))
## Reading layer `MunicipiosNarino' from data source
## `/cloud/project/MunicipiosNarino.shp' using driver `ESRI Shapefile'
## Simple feature collection with 66 features and 11 fields
## Geometry type: POLYGON
## Dimension: XY
## Bounding box: xmin: -79.01021 ymin: 0.3613481 xmax: -76.83368 ymax: 2.683898
## Geodetic CRS: WGS 84
## Simple feature collection with 66 features and 11 fields
## Geometry type: POLYGON
## Dimension: XY
## Bounding box: xmin: -79.01021 ymin: 0.3613481 xmax: -76.83368 ymax: 2.683898
## Geodetic CRS: WGS 84
## First 10 features:
## DPTO_CCDGO MPIO_CCDGO MPIO_CNMBR
## 1 52 083 BELÉN
## 2 52 110 BUESACO
## 3 52 203 COLÓN
## 4 52 480 NARIÑO
## 5 52 506 OSPINA
## 6 52 720 SAPUYES
## 7 52 786 TAMINANGO
## 8 52 788 TANGUA
## 9 52 240 CHACHAGÜÍ
## 10 52 254 EL PEÑOL
## MPIO_CRSLC MPIO_NAREA
## 1 Ordenanza 53 Noviembre 29 de 1985 41.84541
## 2 1899 635.96083
## 3 Ordenanza 37 de 1921 61.75053
## 4 Ordenanza 027 de Noviembre. 29 de 1999. Decreto 0312 del 24 25.31281
## 5 Ordenanza 50 de 1865 64.84321
## 6 1849 115.54851
## 7 1834 234.65783
## 8 Ordenanza 103 de 1874 217.95977
## 9 Ordenanza 20 Noviembre 24 de 1992 146.27176
## 10 Ordenanza 036 de Diciembre 7 de 1998 119.85744
## MPIO_CCNCT MPIO_NANO DPTO_CNMBR SHAPE_AREA SHAPE_LEN ORIG_FID
## 1 52083 2020 NARIÑO 0.003391678 0.3732840 0
## 2 52110 2020 NARIÑO 0.051533090 1.2292312 1
## 3 52203 2020 NARIÑO 0.005005108 0.4592866 2
## 4 52480 2020 NARIÑO 0.002050175 0.2642048 3
## 5 52506 2020 NARIÑO 0.005249269 0.3371496 4
## 6 52720 2020 NARIÑO 0.009351438 0.6599792 5
## 7 52786 2020 NARIÑO 0.019009395 0.6601636 6
## 8 52788 2020 NARIÑO 0.017652117 0.7808421 7
## 9 52240 2020 NARIÑO 0.011849554 0.7373282 8
## 10 52254 2020 NARIÑO 0.009707107 0.5550189 9
## geometry
## 1 POLYGON ((-77.07227 1.63422...
## 2 POLYGON ((-77.23516 1.45240...
## 3 POLYGON ((-77.04473 1.67173...
## 4 POLYGON ((-77.34282 1.31465...
## 5 POLYGON ((-77.55776 1.07006...
## 6 POLYGON ((-77.71499 1.0915,...
## 7 POLYGON ((-77.32644 1.67981...
## 8 POLYGON ((-77.36152 1.19568...
## 9 POLYGON ((-77.30295 1.51777...
## 10 POLYGON ((-77.39239 1.60127...
narino %>% dplyr::select( MPIO_CCNCT, MPIO_CNMBR) -> narino
narino
## Simple feature collection with 66 features and 2 fields
## Geometry type: POLYGON
## Dimension: XY
## Bounding box: xmin: -79.01021 ymin: 0.3613481 xmax: -76.83368 ymax: 2.683898
## Geodetic CRS: WGS 84
## First 10 features:
## MPIO_CCNCT MPIO_CNMBR geometry
## 1 52083 BELÉN POLYGON ((-77.07227 1.63422...
## 2 52110 BUESACO POLYGON ((-77.23516 1.45240...
## 3 52203 COLÓN POLYGON ((-77.04473 1.67173...
## 4 52480 NARIÑO POLYGON ((-77.34282 1.31465...
## 5 52506 OSPINA POLYGON ((-77.55776 1.07006...
## 6 52720 SAPUYES POLYGON ((-77.71499 1.0915,...
## 7 52786 TAMINANGO POLYGON ((-77.32644 1.67981...
## 8 52788 TANGUA POLYGON ((-77.36152 1.19568...
## 9 52240 CHACHAGÜÍ POLYGON ((-77.30295 1.51777...
## 10 52254 EL PEÑOL POLYGON ((-77.39239 1.60127...
rename(narino, MPIO_CCDGO = MPIO_CCNCT)
## Simple feature collection with 66 features and 2 fields
## Geometry type: POLYGON
## Dimension: XY
## Bounding box: xmin: -79.01021 ymin: 0.3613481 xmax: -76.83368 ymax: 2.683898
## Geodetic CRS: WGS 84
## First 10 features:
## MPIO_CCDGO MPIO_CNMBR geometry
## 1 52083 BELÉN POLYGON ((-77.07227 1.63422...
## 2 52110 BUESACO POLYGON ((-77.23516 1.45240...
## 3 52203 COLÓN POLYGON ((-77.04473 1.67173...
## 4 52480 NARIÑO POLYGON ((-77.34282 1.31465...
## 5 52506 OSPINA POLYGON ((-77.55776 1.07006...
## 6 52720 SAPUYES POLYGON ((-77.71499 1.0915,...
## 7 52786 TAMINANGO POLYGON ((-77.32644 1.67981...
## 8 52788 TANGUA POLYGON ((-77.36152 1.19568...
## 9 52240 CHACHAGÜÍ POLYGON ((-77.30295 1.51777...
## 10 52254 EL PEÑOL POLYGON ((-77.39239 1.60127...
library(climateR)
tc_prcp = getTerraClim(narino, param = "prcp", startDate = "2019-02-01")
## Spherical geometry (s2) switched off
## Spherical geometry (s2) switched on
tc_tmp <- tc_prcp[[1]]
tc_tmp
## class : RasterStack
## dimensions : 57, 53, 3021, 1 (nrow, ncol, ncell, nlayers)
## resolution : 0.04166667, 0.04166667 (x, y)
## extent : -79.04167, -76.83333, 0.3333333, 2.708333 (xmin, xmax, ymin, ymax)
## crs : +proj=longlat +datum=WGS84 +no_defs
## names : X2019.02
## min values : 61.8
## max values : 580.5
library(leaflet)
library(RColorBrewer)
library(raster)
pal <- colorNumeric(c("red", "orange", "#fcc000","yellow", "cyan", "blue", "#3240cd"), values(tc_tmp$X2019.02),
na.color = "transparent")
leaflet() %>% addTiles() %>%
addRasterImage(tc_tmp$X2019.02 , colors = pal, opacity = 0.8) %>%
addLegend(pal = pal, values = values(tc_tmp$X2019.02),
title = "Rainfall-Feb.2019 [mm]")
## Warning in showSRID(uprojargs, format = "PROJ", multiline = "NO", prefer_proj =
## prefer_proj): Discarded ellps WGS 84 in Proj4 definition: +proj=merc +a=6378137
## +b=6378137 +lat_ts=0 +lon_0=0 +x_0=0 +y_0=0 +k=1 +units=m +nadgrids=@null
## +wktext +no_defs +type=crs
## Warning in showSRID(uprojargs, format = "PROJ", multiline = "NO", prefer_proj =
## prefer_proj): Discarded datum World Geodetic System 1984 in Proj4 definition
## Warning in colors(.): Some values were outside the color scale and will be
## treated as NA
(param_meta$terraclim)
## common.name call description units
## 1 aet aet Actual Evapotranspiration mm
## 2 water_deficit def Climatic Water Deficit mm
## 3 palmer PDSI Palmer Drought Severity Index
## 4 pet pet Reference Evapotranspiration mm
## 5 prcp ppt Accumulated Precipitation mm
## 6 q q Runoff mm
## 7 soilm soil Soil Moisture at End of Month mm
## 8 srad srad Downward Shortwave Radiation Flux at the Surface W/m2
## 9 swe swe Snow Water Equivalent at End of Month mm
## 10 tmax tmax Maximum 2-m Temperature C
## 11 tmin tmin Minimum 2-m Temperature C
## 12 vp vap 2-m Vapor Pressure kPa
## 13 vpd vpd Mean Vapor Pressure Deficit kPa
## 14 wind ws Wind Speed at 10-m m/s
library(climateR)
library(s2)
tc_prcp = getTerraClim(narino, param = "prcp", startDate = "2019-02-01")
## Spherical geometry (s2) switched off
## Spherical geometry (s2) switched on
tc_tmp <- tc_prcp[[1]]
tc_tmp
## class : RasterStack
## dimensions : 57, 53, 3021, 1 (nrow, ncol, ncell, nlayers)
## resolution : 0.04166667, 0.04166667 (x, y)
## extent : -79.04167, -76.83333, 0.3333333, 2.708333 (xmin, xmax, ymin, ymax)
## crs : +proj=longlat +datum=WGS84 +no_defs
## names : X2019.02
## min values : 61.8
## max values : 580.5
library(leaflet)
pal <- colorNumeric(c("#fc7300","orange", "yellow","#9acd32", "green"), values(tc_tmp$X2019.02),
na.color = "transparent")
leaflet() %>% addTiles() %>%
addRasterImage(tc_tmp$X2019.02, colors = pal, opacity = 0.8) %>%
addLegend(pal = pal, values = values(tc_tmp$X2019.02),
title = "PDSI-Feb.2019")
## Warning in colors(.): Some values were outside the color scale and will be
## treated as NA
library(climateR)
wat_def = getTerraClimNormals(narino, param = "water_deficit", period = "19812010", month=2)
## Spherical geometry (s2) switched off
## Spherical geometry (s2) switched on
wat_def
## $terraclim_19812010_water_deficit
## class : RasterStack
## dimensions : 57, 53, 3021, 1 (nrow, ncol, ncell, nlayers)
## resolution : 0.04166667, 0.04166667 (x, y)
## extent : -79.04167, -76.83333, 0.3333333, 2.708333 (xmin, xmax, ymin, ymax)
## crs : +proj=longlat +datum=WGS84 +no_defs
## names : X02
## min values : 0
## max values : 31.8
tc_tmp <- wat_def[[1]]
pal <- colorNumeric(c("green", "#9acd32","yellow", "orange",
"#fc7300"), values(tc_tmp$X02),
na.color = "transparent")
leaflet() %>% addTiles() %>%
addRasterImage(tc_tmp$X02, colors = pal, opacity = 0.8) %>%
addLegend(pal = pal, values = values(tc_tmp$X02),
title = "WaterDeficit-February")
## Warning in colors(.): Some values were outside the color scale and will be
## treated as NA
3.2 Obtaining climatic variables
Precipitacion
oneyear = seq(1,12)
library(climateR)
library(AOI)
precip = getTerraClimNormals(narino, param = "prcp", period = "19812010", month=oneyear)
## Spherical geometry (s2) switched off
## Spherical geometry (s2) switched on
precip
## $terraclim_19812010_prcp
## class : RasterStack
## dimensions : 57, 53, 3021, 12 (nrow, ncol, ncell, nlayers)
## resolution : 0.04166667, 0.04166667 (x, y)
## extent : -79.04167, -76.83333, 0.3333333, 2.708333 (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 : 40.0, 53.5, 66.7, 89.4, 44.5, 23.3, 12.3, 12.8, 27.4, 59.0, 68.6, 53.5
## max values : 443.8, 404.0, 408.1, 533.4, 691.7, 548.0, 478.5, 439.5, 552.8, 559.4, 513.1, 522.9
capas <-names(precip)
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)
}
year_precip <- sum(precip$terraclim_19812010_prcp)
year_precip
## class : RasterLayer
## dimensions : 57, 53, 3021 (nrow, ncol, ncell)
## resolution : 0.04166667, 0.04166667 (x, y)
## extent : -79.04167, -76.83333, 0.3333333, 2.708333 (xmin, xmax, ymin, ymax)
## crs : +proj=longlat +datum=WGS84 +no_defs
## source : memory
## names : layer
## values : 608.7, 5597.1 (min, max)
ppal <- colorNumeric(c("orange", "yellow", "#a1d99b", "#deebf7", "#9ecae1","#3182bd", "darkblue"), values(year_precip$layer),
na.color = "transparent")
leaflet() %>% addTiles() %>%
addRasterImage(year_precip, colors = ppal, opacity = 0.8) %>%
addLegend(pal = ppal, values = values(year_precip$layer),
title = "Annual Precipitation")
## Warning in colors(.): Some values were outside the color scale and will be
## treated as NA
writeRaster(tc_tmp, filename="precip", format="GTiff", overwrite=TRUE)
Bibliografia Lizarazo, I., 2021. A tutorial on getting global TerraClimate and CHIRPS data from R. Available at https://rpubs.com/ials2un/trrclmt.