Skip to content

Instantly share code, notes, and snippets.

@khufkens
Created May 10, 2018 09:04
  • Star 0 You must be signed in to star a gist
  • Fork 0 You must be signed in to fork a gist
Star You must be signed in to star a gist
Save khufkens/0cbd66f06829014c62a335b65509873d to your computer and use it in GitHub Desktop.
CHIRP drought / precipitation differences between the Amazon and Congo Basin during the 2015 El Nino
# map El Nino drought in the Amazon and Africa
library(raster)
library(tidyverse)
library(maps)
library(RColorBrewer)
# set colour theme
colours = brewer.pal(11,"RdBu")
# download CHIRP data
server = "ftp://ftp.chg.ucsb.edu/pub/org/chg/products/CHIRP/monthly/"
years = 2012:2016
months = 1:12
lapply(years,function(year){
lapply(months,function(month){
# set paths
url = sprintf("%sCHIRP.%s.%02d.tif",server,year,month)
filename = file.path(tempdir(),
sprintf("CHIRP.%s.%02d.tif",year,month))
# download data
download.file(url = url,
destfile = filename,
quiet = TRUE)
})
})
# read files
files = list.files(tempdir(),"*.tif", full.names = TRUE)
# get vector for names
year = lapply(files,function(file){as.numeric(strsplit(file,"\\.")[[1]][2])})
month = lapply(files,function(file){as.numeric(strsplit(file,"\\.")[[1]][3])})
# do a stackApply()
r_pre_el_nino = stack(files[year < 2015])
r_el_nino = stack(files[year >= 2015])
# calculate rasters with mean precip before and during
# the 2015 el nino event
r_pre_el_nino = stackApply(r_pre_el_nino,
indices = month[year < 2015],
fun = mean)
r_el_nino = stackApply(r_el_nino,
indices = month[year >= 2015],
fun = mean)
# calculate difference between precip regimes
el_nino_diff = r_el_nino - r_pre_el_nino
mean_precip = mean(el_nino_diff)
# count dry months
s_el_nino = sum(r_el_nino < 100, na.rm = TRUE)
s_pre_el_nino = sum(r_pre_el_nino < 100, na.rm = TRUE)
threshold_precip = s_el_nino - s_pre_el_nino
pdf("~/mean_temp_diff.pdf",10,10.5)
par(mfrow = c(2, 1))
# Africa
par(mar = c(1,5,5,5))
plot(mean_precip,
xlim = c(0, 55),
ylim = c(-15, 15),
zlim = c(-100,100),
col = colours,
ylab = expression("Longitude (" * degree * ")"),
legend = FALSE)
maps::map("world", add = TRUE)
mtext(text = expression(Delta * "precip. (mm) (El Nino - pre-El Nino)"),
line = 1.5,
cex = 2)
# Amazon
par(mar = c(5,5,1,5))
plot(mean_precip,
xlim = c(-85, -30),
ylim = c(-15, 15),
zlim = c(-100,100),
col = colours,
xlab = expression("Latitude (" * degree * ")"),
ylab = expression("Longitude (" * degree * ")"))
maps::map("world", add = TRUE)
dev.off()
pdf("~/dry_monhts.pdf",10,10.5)
par(mfrow = c(2, 1))
# Africa
par(mar = c(1,5,5,5))
plot(threshold_precip,
xlim = c(0, 55),
ylim = c(-15, 15),
zlim = c(-5,5),
col = colours,
ylab = expression("Longitude (" * degree * ")"),
legend = FALSE)
maps::map("world", add = TRUE)
mtext(text = expression(Delta * " months < 100 mm (El Nino - pre-El Nino)"),
line = 1.5,
cex = 2)
# Amazon
par(mar = c(5,5,1,5))
plot(threshold_precip,
xlim = c(-85, -30),
ylim = c(-15, 15),
zlim = c(-5,5),
col = colours,
xlab = expression("Latitude (" * degree * ")"),
ylab = expression("Longitude (" * degree * ")"))
maps::map("world", add = TRUE)
dev.off()
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment