# Code by Alastair Potts
# 9 October 2023
# 

setwd("D:\\Google Drive\\0_SpekboomResearchGroup2.0\\3_Experiments\\RegionalTransect21(10)_ParentPlantMaterialExperiment\\manuscript/")


library(groundhog)
set.groundhog.folder("C:\\Groundhog_R\\R4.2.3_2023-04-01")
GroundhogDay <- '2023-04-01'

#groundhog.library("tidyverse",GroundhogDay)
groundhog.library("slider",GroundhogDay)
groundhog.library("glue",GroundhogDay)
groundhog.library("lubridate",GroundhogDay)
groundhog.library("drought",GroundhogDay)
groundhog.library("dplyr",GroundhogDay)
groundhog.library("tidyr",GroundhogDay)
groundhog.library("ggplot2",GroundhogDay)
groundhog.library("patchwork",GroundhogDay)

qCols <- rev(c("#3B9AB2", "#78B7C5", "#74A089", "#E1AF00", "#F21A00"))

read.csv("data_manuscript/chirps_data.csv")->chirps.dat

### Converting daily values into monthly values
chirps.dat%>%
  as_tibble()%>%
  select(Site,date,chirps)%>%
  mutate(year=year(date),
         month=month(date))%>%
  group_by(Site,year,month)%>%
  summarise(ppt=sum(chirps))%>%
  ungroup()->
  chirps.dat_monthly

# Generating 12 month spi values
chirps.dat_monthly%>%
  group_by(Site)%>%
  #filter(year<2022)%>%
  group_by(Site)%>%
  mutate(spi=
           SDI(ppt,12)%>%as.data.frame()%>%
           pivot_longer(cols=1:12,names_to = "Temp",values_to = "spi3")%>%
           select("spi3")%>%unlist())%>%
  mutate(Date=(ym(paste(year,month))))->
  forGG

forGG%>%
  group_by(Site,month)%>%
  summarise(u_ppt=mean(ppt))->
  meanPPT
###

### Figure S1
forGG%>%
  ungroup()%>%
  left_join(meanPPT)%>%
  group_by(Site,month)%>%
  mutate(posneg=ifelse(ppt>u_ppt,"darkgreen","red"))%>%
  ungroup()%>%
  filter(year>=2020,year<2022)%>%
  mutate(Date=ymd(paste(year,month,14)))%>%
  mutate(Site=factor(paste("Site:",Site),
                     levels=paste("Site:",1:10)))%>%
  
  #filter(Site==7)%>%
  ggplot()+
  #geom_point(aes(Date,u_ppt))+
  geom_line(aes(Date,u_ppt),linewidth=1.5,alpha=0.5)+
  #geom_line(aes(Date,q33))+
  geom_segment(aes(x=Date,xend=Date,y=u_ppt,
                   yend=ppt,
                   colour=as.factor(posneg)),
               linewidth=3.5)+
  scale_colour_manual(values=qCols[c(5,1)],labels=c("+","-"))+
  geom_point(aes(Date,u_ppt),size=3)+
  #geom_line(aes(Date,ppt))+
  geom_vline(xintercept=ymd(c("2021-01-01","2022-01-1")))+
  geom_vline(xintercept=ymd("2021-10-13"),linetype='dashed',lwd=2)+
  theme_bw()+
  theme(legend.position="none")+
  xlab("")+
  ylab("Monthly precipitation (mm)")+
  scale_x_date(date_breaks = "1 month", date_labels =  "%b %Y")+
  theme(axis.text.x=element_text(angle=60, hjust=1))+
  facet_wrap(~Site,ncol=2)+
  coord_cartesian(ylim=c(0,160))

ggsave("FinalFiguresV2/Fig_S1.jpg",width=7,height=8)