RedFlags_Mer.Rmd

---
title: "Redflags"
author: "Gustavo Facincani Dourado"
date: "6/12/2020"
output: html_document
---

```{r}
library(dplyr)
library(lfstat)
library(reshape2)
```

```{r}
Livneh_Redflag_Mer <- read_csv("C:/Users/gusta/Desktop/PhD/CERCWET/Hydropower/Merced/Red flags/Livneh/IFR_redflag.csv", col_types = cols(node = col_date("%m/%d/%Y"), .default = col_double()))[-c(1:3), c(1,5,9)]%>%
    filter(between(node, as.Date("1980-10-01"), as.Date("2010-09-30")))
#CanESM2_RedFlag_Mer$CanESM2 <- rowSums(CanESM2_RedFlag_Mer[c(2,3)])
Livneh_Redflag_Mer$GCM <- as.factor("Livneh (Historical)")
Livneh_Redflag_Mer


CanESM2_RedFlag_Mer <- read_csv("C:/Users/gusta/Desktop/PhD/CERCWET/Hydropower/Merced/Red flags/CanESM2_rcp85/IFR_redflag.csv", col_types = cols(node = col_date("%m/%d/%Y"), .default = col_double()))[-c(1:3), c(1,5,9)]
#CanESM2_RedFlag_Mer$CanESM2 <- rowSums(CanESM2_RedFlag_Mer[c(2,3)])
CanESM2_RedFlag_Mer$GCM <- as.factor("CanESM2")
CanESM2_RedFlag_Mer

CNRM_RedFlag_Mer <- read_csv("C:/Users/gusta/Desktop/PhD/CERCWET/Hydropower/Merced/Red flags/CNRM-CM5_rcp85/IFR_redflag.csv", col_types = cols(node = col_date("%m/%d/%Y"), .default = col_double()))[-c(1:3), c(1,5,9)] 
#CNRM_RedFlag_Mer$`CNRM-CM5` <- rowSums(CNRM_RedFlag_Mer[c(2,3)])
CNRM_RedFlag_Mer$GCM <- as.factor("CNRM-CM5")
CNRM_RedFlag_Mer


HadGEM2_RedFlag_Mer <- read_csv("C:/Users/gusta/Desktop/PhD/CERCWET/Hydropower/Merced/Red flags/HadGEM2-ES_rcp85/IFR_redflag.csv", col_types = cols(node = col_date("%m/%d/%Y"), .default = col_double()))[-c(1:3), c(1,5,9)]
#HadGEM2_RedFlag_Mer$`HadGEM2-ES` <- rowSums(HadGEM2_RedFlag_Mer[c(2,3)])
HadGEM2_RedFlag_Mer$GCM <- as.factor("HadGEM2-ES")
HadGEM2_RedFlag_Mer 

MIROC5_RedFlag_Mer <- read_csv("C:/Users/gusta/Desktop/PhD/CERCWET/Hydropower/Merced/Red flags/MIROC5_rcp85/IFR_redflag.csv", col_types = cols(node = col_date("%m/%d/%Y"), .default = col_double()))[-c(1:3), c(1,5,9)] 
#MIROC5_RedFlag_Mer$MIROC5 <- rowSums(MIROC5_RedFlag_Mer[c(2,3)])
MIROC5_RedFlag_Mer$GCM <- as.factor("MIROC5")
MIROC5_RedFlag_Mer
```

```{r}
IFR_redgflag_Mer1 <- rbind(Livneh_Redflag_Mer, CanESM2_RedFlag_Mer, CNRM_RedFlag_Mer, 
                 HadGEM2_RedFlag_Mer, MIROC5_RedFlag_Mer) %>% 
  mutate(WaterYear = as.factor(water_year(node, origin= "usgs")),
         Scenario = "RCP 8.5",
         Basin = "Merced") %>% 
  group_by(WaterYear, GCM, Scenario, Basin) %>%
dplyr::summarize(`IFR bl New Exchequer Dam` = sum(`IFR bl New Exchequer Dam_3`),
                 `IFR at Shaffer Bridge` = sum(`IFR at Shaffer Bridge_3`))
IFR_redgflag_Mer1

IFR_redgflag_Mer <- melt(IFR_redgflag_Mer1, id = c("WaterYear", "GCM", "Scenario", "Basin"))
IFR_redgflag_Mer

IFR_redgflag_Mer$GCM <- factor(IFR_redgflag_Mer$GCM, levels = c("Livneh (Historical)","HadGEM2-ES", "CanESM2", "MIROC5", "CNRM-CM5"))

```


```{r}
library(ggthemes)
Merced_IFR <- ggplot(IFR_redgflag_Mer) +
     theme_bw(base_size=12, base_family='Times New Roman') +
  geom_bar(aes(x= WaterYear, y=value, fill = variable), color = "black", position = "stack", stat = "identity") +
scale_color_colorblind() +
   scale_y_continuous(expand = c(0, NA)) +
  scale_x_discrete(breaks = c("1981", "1986", "1991", "1996", "2001", "2006", 
                              "2031", "2036", "2041", "2046", "2051", "2056"),
                    expand = c(0, NA))+

  facet_wrap(~ GCM, ncol =2, scales = "free_x") +
  
  labs(title = "Merced River",
       subtitle = "RCP 8.5 Scenario",
  #  subtitle = "CanESM2",
    x = element_blank(),
       y = "Total Instream Flow Requirement Gross Deficit (days)") +
 theme(legend.title = element_blank(),
    legend.position = "bottom",
        legend.direction = "horizontal",
        legend.box.margin = margin(t = -17),
        plot.title = element_text(hjust = 0.5),
        plot.subtitle = element_text(hjust = 0.5),
        strip.placement = "outside",
        strip.background = element_blank(),
        panel.spacing.y = unit(0, "lines"))
Merced_IFR +
png("Mer_IFR_Deficit2.png", units ="in", width=6, height=6, res = 300)

```


```{r}
#Merced facet with 1 column
Merced_IFR <- ggplot(IFR_redgflag_Mer) +
     theme_bw(base_size=12, base_family='Times New Roman') +
  geom_bar(aes(x= WaterYear, y=value, fill = variable), color = "black", position = "stack", stat = "identity") +
scale_color_colorblind() +
   scale_y_continuous(expand = c(0, NA)) +
  scale_x_discrete(breaks = c("1981", "1986", "1991", "1996", "2001", "2006", "2010", 
                              "2031", "2036", "2041", "2046", "2051", "2056", "2060"),
                    expand = c(0, NA))+

  facet_wrap(~ GCM, ncol =1, scales = "free_x", strip.position = 'left') +
  
  labs(title = "Merced River",
       subtitle = "RCP 8.5 Scenario",
  #  subtitle = "CanESM2",
    x = element_blank(),
       y = "Total Instream Flow Requirement Gross Deficit (days)") +
 theme(legend.key.size = unit(0.75,"line"),
       legend.text = element_text(size = 9),
       legend.title = element_blank(),
    legend.position = "bottom",
        legend.direction = "horizontal",
        legend.box.margin = margin(t = -17),
        plot.title = element_text(hjust = 0.5),
        plot.subtitle = element_text(hjust = 0.5),
        strip.placement = "outside",
        strip.background = element_blank(),
        panel.spacing.y = unit(0, "lines"))
Merced_IFR +
png("Mer_IFR_Deficit.png", units ="in", width=4.5, height=7.5, res = 300)

```

```{r}
g1 <- Merced_IFR %+% dplyr::filter(IFR_redgflag_Mer, GCM == "Livneh (Historical)")+ theme(legend.position = "none")
g2 <- Merced_IFR %+% dplyr::filter(IFR_redgflag_Mer, GCM != "Livneh (Historical)") + theme(legend.position = "none")

gridExtra::grid.arrange(g1, g2,
                        layout_matrix = 
                          matrix(c(NA, NA, 1, 1, 1, 1,NA, NA,
                                   2, 2, 2, 2, 2, 2, 2, 2,
                                   2, 2, 2, 2, 2, 2, 2, 2),
                                 byrow = TRUE, nrow = 3))
```

```{r}
Livneh_Redflag_Tuo <- read_csv("C:/Users/gusta/Desktop/PhD/CERCWET/Hydropower/tuolumne/Redflags/Livneh/IFR_redflag.csv", col_types = cols(node = col_date("%m/%d/%Y"), .default = col_double()))[-c(1:3), c(1,9, 13, 17, 21)]%>%
    filter(between(node, as.Date("1980-10-01"), as.Date("2010-09-30")))
#CanESM2_RedFlag_Mer$CanESM2 <- rowSums(CanESM2_RedFlag_Mer[c(2,3)])
Livneh_Redflag_Tuo$GCM <- as.factor("Livneh (Historical)")
Livneh_Redflag_Tuo


CanESM2_RedFlag_Tuo <- read_csv("C:/Users/gusta/Desktop/PhD/CERCWET/Hydropower/tuolumne/Redflags/CanESM2_rcp85/IFR_redflag.csv", col_types = cols(node = col_date("%m/%d/%Y"), .default = col_double()))[-c(1:3), c(1,9, 13, 17, 21)]
#CanESM2_RedFlag_Mer$CanESM2 <- rowSums(CanESM2_RedFlag_Mer[c(2,3)])
CanESM2_RedFlag_Tuo$GCM <- as.factor("CanESM2")
CanESM2_RedFlag_Tuo

CNRM_RedFlag_Tuo <- read_csv("C:/Users/gusta/Desktop/PhD/CERCWET/Hydropower/tuolumne/Redflags/CNRM-CM5_rcp85/IFR_redflag.csv", col_types = cols(node = col_date("%m/%d/%Y"), .default = col_double()))[-c(1:3), c(1,9, 13, 17, 21)] 
#CNRM_RedFlag_Mer$`CNRM-CM5` <- rowSums(CNRM_RedFlag_Mer[c(2,3)])
CNRM_RedFlag_Tuo$GCM <- as.factor("CNRM-CM5")
CNRM_RedFlag_Tuo


HadGEM2_RedFlag_Tuo <- read_csv("C:/Users/gusta/Desktop/PhD/CERCWET/Hydropower/tuolumne/Redflags/HadGEM2-ES_rcp85/IFR_redflag.csv", col_types = cols(node = col_date("%m/%d/%Y"), .default = col_double()))[-c(1:3), c(1,9, 13, 17, 21)]
#HadGEM2_RedFlag_Mer$`HadGEM2-ES` <- rowSums(HadGEM2_RedFlag_Mer[c(2,3)])
HadGEM2_RedFlag_Tuo$GCM <- as.factor("HadGEM2-ES")
HadGEM2_RedFlag_Tuo 

MIROC5_RedFlag_Tuo <- read_csv("C:/Users/gusta/Desktop/PhD/CERCWET/Hydropower/tuolumne/Redflags/MIROC5_rcp85/IFR_redflag.csv", col_types = cols(node = col_date("%m/%d/%Y"), .default = col_double()))[-c(1:3), c(1,9, 13, 17, 21)] 
#MIROC5_RedFlag_Mer$MIROC5 <- rowSums(MIROC5_RedFlag_Mer[c(2,3)])
MIROC5_RedFlag_Tuo$GCM <- as.factor("MIROC5")
MIROC5_RedFlag_Tuo


```

```{r}
IFR_redgflag_Tuo<- rbind(Livneh_Redflag_Tuo, CanESM2_RedFlag_Tuo, CNRM_RedFlag_Tuo, 
                 HadGEM2_RedFlag_Tuo, MIROC5_RedFlag_Tuo) %>% 
  mutate(WaterYear = as.factor(water_year(node, origin= "usgs")),
         Scenario = "RCP 8.5",
         Basin = "Tuolumne") %>% 
  group_by(WaterYear, GCM, Scenario, Basin) %>%
dplyr::summarize(`IFR bl Lake Eleanor` = sum(`IFR bl Lake Eleanor_3`),
                 `IFR bl Hetch Hetchy Reservoir` = sum(`IFR bl Hetch Hetchy Reservoir_3`),
                 `IFR at La Grange` = sum(`IFR at La Grange_3`),
                 `IFR bl Cherry Lake` = sum(`IFR bl Cherry Lake_3`))
IFR_redgflag_Tuo

IFR_redgflag_Tuo <- melt(IFR_redgflag_Tuo, id = c("WaterYear", "GCM", "Scenario", "Basin"))
IFR_redgflag_Tuo

IFR_redgflag_Tuo$GCM <- factor(IFR_redgflag_Tuo$GCM, levels = c("Livneh (Historical)","HadGEM2-ES", "CanESM2", "MIROC5", "CNRM-CM5"))

```


```{r}
Tuolumne_IFR <- ggplot(IFR_redgflag_Tuo) +
     theme_bw(base_size=12, base_family='Times New Roman') +
  geom_bar(aes(x= WaterYear, y=value, fill = variable), color = "black", position = "stack", stat = "identity") +
scale_color_colorblind() +
   scale_y_continuous(expand = c(0, NA)) +
  scale_x_discrete(breaks = c("1981", "1986", "1991", "1996", "2001", "2006", 
                              "2031", "2036", "2041", "2046", "2051", "2056"),
                    expand = c(0, NA))+

  facet_wrap(~ GCM, ncol =2, scales = "free_x") +
  
  labs(title = "Merced River",
       subtitle = "RCP 8.5 Scenario",
  #  subtitle = "CanESM2",
    x = element_blank(),
       y = "Total Instream Flow Requirement Gross Deficit (days)") +
 theme(legend.title = element_blank(),
    legend.position = "bottom",
        legend.direction = "horizontal",
        legend.box.margin = margin(t = -17),
        plot.title = element_text(hjust = 0.5),
        plot.subtitle = element_text(hjust = 0.5),
        strip.placement = "outside",
        strip.background = element_blank(),
        panel.spacing.y = unit(0, "lines"))
Tuolumne_IFR +
png("Tuo_IFR_Deficit2.png", units ="in", width=6, height=6, res = 300)

```


```{r}
#Tuolumne facet with 1 column
Tuolumne_IFR2 <- ggplot(IFR_redgflag_Tuo) +
     theme_bw(base_size=12, base_family='Times New Roman') +
  geom_bar(aes(x= WaterYear, y=value, fill = variable), color = "black", position = "stack", stat = "identity") +
scale_color_colorblind() +
   scale_y_continuous(expand = c(0, NA)) +
  scale_x_discrete(breaks = c("1981", "1986", "1991", "1996", "2001", "2006", "2010", 
                              "2031", "2036", "2041", "2046", "2051", "2056", "2060"),
                    expand = c(0, NA))+

  facet_wrap(~ GCM, ncol =1, scales = "free_x", strip.position = 'left') +
  
  labs(title = "Tuolumne River",
       subtitle = "RCP 8.5 Scenario",
  #  subtitle = "CanESM2",
    x = element_blank(),
       y = "Total Instream Flow Requirement Gross Deficit (days)") +
 theme(legend.title = element_blank(),
       legend.key.size = unit(0.75,"line"),
       legend.text = element_text(size = 9),
    legend.position = "bottom",
        legend.direction = "horizontal",
        legend.box.margin = margin(t = -17),
        plot.title = element_text(hjust = 0.5),
        plot.subtitle = element_text(hjust = 0.5),
        strip.placement = "outside",
        strip.background = element_blank(),
        panel.spacing.y = unit(0, "lines"))+
  guides(fill = guide_legend(nrow = 2))
Tuolumne_IFR2 +
png("Tuo_IFR_Deficit.png", units ="in", width=4.5, height=7.5, res = 300)

```

```{r}
Livneh_Redflag_Stn <- read_csv("C:/Users/gusta/Desktop/PhD/CERCWET/Hydropower/stanislaus/Redflags/Livneh/IFR_redflag.csv", col_types = cols(node = col_date("%m/%d/%Y"), .default = col_double()))[-c(1:3), c(1, 5, 17, 21, 25,29,33,37,41,45,49,53,57,61,65,69,73,77,81)]%>%
    filter(between(node, as.Date("1980-10-01"), as.Date("2010-09-30")))
#CanESM2_RedFlag_Mer$CanESM2 <- rowSums(CanESM2_RedFlag_Mer[c(2,3)])
Livneh_Redflag_Stn$GCM <- as.factor("Livneh (Historical)")
Livneh_Redflag_Stn


CanESM2_RedFlag_Stn <- read_csv("C:/Users/gusta/Desktop/PhD/CERCWET/Hydropower/stanislaus/Redflags/CanESM2_rcp85/IFR_redflag.csv", col_types = cols(node = col_date("%m/%d/%Y"), .default = col_double()))[-c(1:3), c(1, 5, 17, 21, 25,29,33,37,41,45,49,53,57,61,65,69,73,77,81)]
#CanESM2_RedFlag_Mer$CanESM2 <- rowSums(CanESM2_RedFlag_Mer[c(2,3)])
CanESM2_RedFlag_Stn$GCM <- as.factor("CanESM2")
CanESM2_RedFlag_Stn

CNRM_RedFlag_Stn <- read_csv("C:/Users/gusta/Desktop/PhD/CERCWET/Hydropower/stanislaus/Redflags/CNRM-CM5_rcp85/IFR_redflag.csv", col_types = cols(node = col_date("%m/%d/%Y"), .default = col_double()))[-c(1:3), c(1, 5, 17, 21, 25,29,33,37,41,45,49,53,57,61,65,69,73,77,81)]
#CNRM_RedFlag_Mer$`CNRM-CM5` <- rowSums(CNRM_RedFlag_Mer[c(2,3)])
CNRM_RedFlag_Stn$GCM <- as.factor("CNRM-CM5")
CNRM_RedFlag_Stn


HadGEM2_RedFlag_Stn <- read_csv("C:/Users/gusta/Desktop/PhD/CERCWET/Hydropower/stanislaus/Redflags/HadGEM2-ES_rcp85/IFR_redflag.csv", col_types = cols(node = col_date("%m/%d/%Y"), .default = col_double()))[-c(1:3), c(1, 5, 17, 21, 25,29,33,37,41,45,49,53,57,61,65,69,73,77,81)]
#HadGEM2_RedFlag_Mer$`HadGEM2-ES` <- rowSums(HadGEM2_RedFlag_Mer[c(2,3)])
HadGEM2_RedFlag_Stn$GCM <- as.factor("HadGEM2-ES")
HadGEM2_RedFlag_Stn 

MIROC5_RedFlag_Stn <- read_csv("C:/Users/gusta/Desktop/PhD/CERCWET/Hydropower/stanislaus/Redflags/MIROC5_rcp85/IFR_redflag.csv", col_types = cols(node = col_date("%m/%d/%Y"), .default = col_double()))[-c(1:3), c(1, 5, 17, 21, 25,29,33,37,41,45,49,53,57,61,65,69,73,77,81)]
#MIROC5_RedFlag_Mer$MIROC5 <- rowSums(MIROC5_RedFlag_Mer[c(2,3)])
MIROC5_RedFlag_Stn$GCM <- as.factor("MIROC5")
MIROC5_RedFlag_Stn


```

```{r}
IFR_redgflag_Stn<- rbind(Livneh_Redflag_Stn, CanESM2_RedFlag_Stn, CNRM_RedFlag_Stn, 
                 HadGEM2_RedFlag_Stn, MIROC5_RedFlag_Stn) %>% 
  mutate(WaterYear = as.factor(water_year(node, origin= "usgs")),
         Scenario = "RCP 8.5") %>% 
  group_by(WaterYear, GCM, Scenario) %>%
dplyr::summarize(`IFR at Murphys Park` = sum(`IFR at Murphys Park_3`),
                 `IFR bl Utica Reservoir` = sum(`IFR bl Utica Reservoir_3`),
                 `IFR bl Pinecrest Lake` = sum(`IFR bl Pinecrest Lake_3`),
                 `IFR bl conf of NF Stanislaus and Beaver Creek` = sum(`IFR bl confluence of NF Stanislaus and Beaver Creek_3`),
                 `IFR bl Beaver Creek Div Dam` = sum(`IFR bl Beaver Creek Diversion Dam_3`),
                 `IFR bl Beardsley Afterbay` = sum(`IFR bl Beardsley Afterbay_3`),
                 `IFR bl Goodwin Reservoir` = sum(`IFR bl Goodwin Reservoir_3`),
                 `IFR bl NF Stanislaus Div Reservoir`= sum(`IFR bl NF Stanislaus Div Res_3`),
                 `IFR bl Sand Bar Diversion` = sum(`IFR bl Sand Bar Div_3`),
                 `IFR bl Relief Reservoir` = sum(`IFR bl Relief Reservoir_3`),
                 `IFR bl Philadelphia Div` = sum(`IFR bl Philadelphia Div_3`),
                 `IFR bl New Spicer Meadow Reservoir` = sum(`IFR bl New Spicer Meadow Reservoir_3`),
                 `IFR bl McKays Point Div` = sum(`IFR bl McKays Point Div_3`),
                 `IFR bl Lyons Reservoir` = sum(`IFR bl Lyons Res_3`),
                 `IFR bl Hunter Reservoir` = sum(`IFR bl Hunter Reservoir_3`),
                 `IFR bl Donnell Lake` = sum(`IFR bl Donnell Lake_3`),
                 `IFR bl Collierville PH discharge` = sum(`IFR bl Collierville PH discharge_3`),
                 `IFR bl Angels Div` = sum(`IFR bl Angels Div_3`))
IFR_redgflag_Stn

IFR_redgflag_Stn <- melt(IFR_redgflag_Stn, id = c("WaterYear", "GCM", "Scenario"))
IFR_redgflag_Stn

IFR_redgflag_Stn$GCM <- factor(IFR_redgflag_Stn$GCM, levels = c("Livneh (Historical)","HadGEM2-ES", "CanESM2", "MIROC5", "CNRM-CM5"))

```

```{r}
Stanislaus_IFR <- ggplot(IFR_redgflag_Stn) +
     theme_bw(base_size=13.6, base_family='Times New Roman') +
  geom_bar(aes(x= WaterYear, y=value, fill = variable), color = "black", position = "stack", stat = "identity") +
  scale_color_colorblind() +
   scale_y_continuous(expand = c(0, NA)) +
  scale_x_discrete(breaks = c("1981", "1986", "1991", "1996", "2001", "2006", 
                              "2031", "2036", "2041", "2046", "2051", "2056"),
                    expand = c(0, NA))+

  facet_wrap(~ GCM, ncol =2, scales = "free_x") +
  
  labs(title = "Stanislaus River",
       subtitle = "RCP 8.5 Scenario",
  #  subtitle = "CanESM2",
    x = element_blank(),
       y = "Total Instream Flow Requirement Gross Deficit (days)") +
 theme(legend.title = element_blank(),
    legend.position = "bottom",
        legend.direction = "horizontal",
        legend.box.margin = margin(t = -17),
        plot.title = element_text(hjust = 0.5),
        plot.subtitle = element_text(hjust = 0.5),
        strip.placement = "outside",
        strip.background = element_blank(),
        panel.spacing.y = unit(0, "lines"))+
    guides(fill = guide_legend(nrow = 6)) #number of rows of the legend
Stanislaus_IFR +
png("Stn_IFR_Deficit.png", units ="in", width=9, height=8, res = 300)

```

```{r}
#Merced facet with 1 column
Stanislaus_IFR2 <- ggplot(IFR_redgflag_Stn) +
     theme_bw(base_size=14, base_family='Times New Roman') +
  geom_bar(aes(x= WaterYear, y=value, fill = variable), color = "black", position = "stack", stat = "identity") +
scale_colour_colorblind() +
   scale_y_continuous(expand = c(0, NA)) +
  scale_x_discrete(breaks = c("1981", "1986", "1991", "1996", "2001", "2006", "2010", 
                              "2031", "2036", "2041", "2046", "2051", "2056", "2059"),
                    expand = c(0, NA))+

  facet_wrap(~ GCM, ncol =1, scales = "free_x", strip.position = 'left') +
  
  labs(title = "Stanislaus River",
       subtitle = "RCP 8.5 Scenario",
  #  subtitle = "CanESM2",
    x = element_blank(),
       y = "Total Instream Flow Requirement Gross Deficit (days)") +
 theme(legend.title = element_blank(),
       legend.key.size = unit(0.75,"line"), #change the size of icons
    legend.position = "bottom",
     legend.text = element_text(size = 9),
        legend.direction = "horizontal",
        legend.box.margin = margin(t = -17),
        plot.title = element_text(hjust = 0.5),
        plot.subtitle = element_text(hjust = 0.5),
        strip.placement = "outside",
        strip.background = element_blank(),
        panel.spacing.y = unit(0, "lines"))+
      guides(fill = guide_legend(nrow = 9)) #change number of rows in the legend
Stanislaus_IFR2 +
png("Stn_IFR_Deficit2.png", units ="in", width=6, height=11, res = 300)

```

```{r}
Livneh_Redflag_USJ <- read_csv("C:/Users/gusta/Desktop/PhD/CERCWET/Hydropower/upper_san_joaquin/historical/Livneh/IFR_redflag.csv", col_types = cols(node = col_date("%m/%d/%Y"), .default = col_double()))[-c(1:3), c(1, 5, 17, 21, 25,29,33,37,41,45,49,53,57,61,65,69,73,77,81,85,89, 93)]%>%
    filter(between(node, as.Date("1980-10-01"), as.Date("2010-09-30")))
#CanESM2_RedFlag_Mer$CanESM2 <- rowSums(CanESM2_RedFlag_Mer[c(2,3)])
Livneh_Redflag_USJ$GCM <- as.factor("Livneh (Historical)")
Livneh_Redflag_USJ


CanESM2_RedFlag_USJ <- read_csv("C:/Users/gusta/Desktop/PhD/CERCWET/Hydropower/upper_san_joaquin/gcms/CanESM2_rcp85/IFR_redflag.csv", col_types = cols(node = col_date("%m/%d/%Y"), .default = col_double()))[-c(1:3), c(1, 5, 17, 21, 25,29,33,37,41,45,49,53,57,61,65,69,73,77,81,85,89,93)]
#CanESM2_RedFlag_Mer$CanESM2 <- rowSums(CanESM2_RedFlag_Mer[c(2,3)])
CanESM2_RedFlag_USJ$GCM <- as.factor("CanESM2")
CanESM2_RedFlag_USJ

CNRM_RedFlag_USJ <- read_csv("C:/Users/gusta/Desktop/PhD/CERCWET/Hydropower/upper_san_joaquin/gcms/CNRM-CM5_rcp85/IFR_redflag.csv", col_types = cols(node = col_date("%m/%d/%Y"), .default = col_double()))[-c(1:3), c(1, 5, 17, 21, 25,29,33,37,41,45,49,53,57,61,65,69,73,77,81,85,89, 93)]
#CNRM_RedFlag_Mer$`CNRM-CM5` <- rowSums(CNRM_RedFlag_Mer[c(2,3)])
CNRM_RedFlag_USJ$GCM <- as.factor("CNRM-CM5")
CNRM_RedFlag_USJ


HadGEM2_RedFlag_USJ <- read_csv("C:/Users/gusta/Desktop/PhD/CERCWET/Hydropower/upper_san_joaquin/gcms/HadGEM2-ES_rcp85/IFR_redflag.csv", col_types = cols(node = col_date("%m/%d/%Y"), .default = col_double()))[-c(1:3), c(1, 5, 17, 21, 25,29,33,37,41,45,49,53,57,61,65,69,73,77,81,85,89, 93)]
#HadGEM2_RedFlag_Mer$`HadGEM2-ES` <- rowSums(HadGEM2_RedFlag_Mer[c(2,3)])
HadGEM2_RedFlag_USJ$GCM <- as.factor("HadGEM2-ES")
HadGEM2_RedFlag_USJ 

MIROC5_RedFlag_USJ <- read_csv("C:/Users/gusta/Desktop/PhD/CERCWET/Hydropower/upper_san_joaquin/gcms/MIROC5_rcp85/IFR_redflag.csv", col_types = cols(node = col_date("%m/%d/%Y"), .default = col_double()))[-c(1:3), c(1, 5, 17, 21, 25,29,33,37,41,45,49,53,57,61,65,69,73,77,81,85,89, 93)]
#MIROC5_RedFlag_Mer$MIROC5 <- rowSums(MIROC5_RedFlag_Mer[c(2,3)])
MIROC5_RedFlag_USJ$GCM <- as.factor("MIROC5")
MIROC5_RedFlag_USJ


```

```{r}
IFR_redgflag_USJ1<- rbind(Livneh_Redflag_USJ, CanESM2_RedFlag_USJ, CNRM_RedFlag_USJ, 
                 HadGEM2_RedFlag_USJ, MIROC5_RedFlag_USJ) %>% 
  mutate(WaterYear = as.factor(water_year(node, origin= "usgs")),
         Scenario = "RCP 8.5")

names(IFR_redgflag_USJ1) = gsub(pattern = "_3", replacement = "", x = names(IFR_redgflag_USJ1))
IFR_redgflag_USJ1 # remove all "_3" at once
```

```{r}
IFR_redgflag_USJ <- IFR_redgflag_USJ1 %>% 
  group_by(WaterYear, GCM, Scenario) %>%
dplyr::summarize_each(funs(sum), `IFR No. Fk. Stevenson Creek above Shaver Lake`:`IFR above Shakeflat Creek`) #summarizing the sum of deficit per location per year all at once
IFR_redgflag_USJ

IFR_redgflag_USJ <- melt(IFR_redgflag_USJ, id = c("WaterYear", "GCM", "Scenario"))
IFR_redgflag_USJ

IFR_redgflag_USJ$GCM <- factor(IFR_redgflag_USJ$GCM, levels = c("Livneh (Historical)","HadGEM2-ES", "CanESM2", "MIROC5", "CNRM-CM5"))

```

```{r}
#Merced facet with 1 column
USJ_IFR <- ggplot(IFR_redgflag_USJ) +
     theme_bw(base_size=14, base_family='Times New Roman') +
  geom_bar(aes(x= WaterYear, y=value, fill = variable), color = "black", position = "stack", stat = "identity") +
scale_colour_colorblind() +
   scale_y_continuous(expand = c(0, NA)) +
  scale_x_discrete(breaks = c("1981", "1986", "1991", "1996", "2001", "2006", "2010", 
                              "2031", "2036", "2041", "2046", "2051", "2056", "2059"),
                    expand = c(0, NA))+

  facet_wrap(~ GCM, ncol =1, scales = "free_x", strip.position = 'left') +
  
  labs(title = "Upper San Joaquin River",
       subtitle = "RCP 8.5 Scenario",
  #  subtitle = "CanESM2",
    x = element_blank(),
       y = "Total Instream Flow Requirement Gross Deficit (days)") +
 theme(legend.title = element_blank(),
       legend.key.size = unit(0.75,"line"), #change the size of icons
    legend.position = "bottom",
     legend.text = element_text(size = 9),
        legend.direction = "horizontal",
        legend.box.margin = margin(t = -17),
        plot.title = element_text(hjust = 0.5),
        plot.subtitle = element_text(hjust = 0.5),
        strip.placement = "outside",
        strip.background = element_blank(),
        panel.spacing.y = unit(0, "lines"))+
      guides(fill = guide_legend(nrow = 11)) #change number of rows in the legend
USJ_IFR +
png("USJ_IFR_Deficit.png", units ="in", width=6, height=11, res = 300)

```


```{r}
#Spill redflags

Livneh_Spill_Mer <- read_csv("C:/Users/gusta/Desktop/PhD/CERCWET/Hydropower/Merced/Red flags/Livneh/IFR_redflag.csv", col_types = cols(node = col_date("%m/%d/%Y"), .default = col_double()))[-c(1:3), c(1,6)]%>%
    filter(between(node, as.Date("1980-10-01"), as.Date("2010-09-30")))
#CanESM2_RedFlag_Mer$CanESM2 <- rowSums(CanESM2_RedFlag_Mer[c(2,3)])
Livneh_Spill_Mer$GCM <- as.factor("Livneh (Historical)")
Livneh_Spill_Mer


CanESM2_Spill_Mer <- read_csv("C:/Users/gusta/Desktop/PhD/CERCWET/Hydropower/Merced/Red flags/CanESM2_rcp85/IFR_redflag.csv", col_types = cols(node = col_date("%m/%d/%Y"), .default = col_double()))[-c(1:3), c(1,6)]
#CanESM2_RedFlag_Mer$CanESM2 <- rowSums(CanESM2_RedFlag_Mer[c(2,3)])
CanESM2_Spill_Mer$GCM <- as.factor("CanESM2")
CanESM2_Spill_Mer

CNRM_Spill_Mer <- read_csv("C:/Users/gusta/Desktop/PhD/CERCWET/Hydropower/Merced/Red flags/CNRM-CM5_rcp85/IFR_redflag.csv", col_types = cols(node = col_date("%m/%d/%Y"), .default = col_double()))[-c(1:3), c(1,6)] 
#CNRM_RedFlag_Mer$`CNRM-CM5` <- rowSums(CNRM_RedFlag_Mer[c(2,3)])
CNRM_Spill_Mer$GCM <- as.factor("CNRM-CM5")
CNRM_Spill_Mer


HadGEM2_Spill_Mer <- read_csv("C:/Users/gusta/Desktop/PhD/CERCWET/Hydropower/Merced/Red flags/HadGEM2-ES_rcp85/IFR_redflag.csv", col_types = cols(node = col_date("%m/%d/%Y"), .default = col_double()))[-c(1:3), c(1,6)]
#HadGEM2_RedFlag_Mer$`HadGEM2-ES` <- rowSums(HadGEM2_RedFlag_Mer[c(2,3)])
HadGEM2_Spill_Mer$GCM <- as.factor("HadGEM2-ES")
HadGEM2_Spill_Mer 

MIROC5_Spill_Mer <- read_csv("C:/Users/gusta/Desktop/PhD/CERCWET/Hydropower/Merced/Red flags/MIROC5_rcp85/IFR_redflag.csv", col_types = cols(node = col_date("%m/%d/%Y"), .default = col_double()))[-c(1:3), c(1,6)] 
#MIROC5_RedFlag_Mer$MIROC5 <- rowSums(MIROC5_RedFlag_Mer[c(2,3)])
MIROC5_Spill_Mer$GCM <- as.factor("MIROC5")
MIROC5_Spill_Mer

```


```{r}
Spill_Mer1 <- rbind(Livneh_Spill_Mer, CanESM2_Spill_Mer, CNRM_Spill_Mer, 
                 HadGEM2_Spill_Mer, MIROC5_Spill_Mer) %>% 
  mutate(WaterYear = as.factor(water_year(node, origin= "usgs")),
         Scenario = "RCP 8.5",
         Basin = "Merced")
Spill_Mer1

Spill_Mer <- melt(Spill_Mer1, id = c("node", "WaterYear", "GCM", "Scenario", "Basin"))
Spill_Mer
```

```{r}
#Maximum flow at Shaffer Bridge is 6500 cfs (561,600,000 cfs in a day) = 15.90274111328398 mcm
Spill_Mer2 <- Spill_Mer %>%
  mutate(MaxFlow = as.numeric(15.90274111328398),
         Spill = as.numeric(value - MaxFlow),
         Count = ifelse(Spill > 0, 1,0))

Spill_Mer2
```


```{r}
#Goodwin - Final point of spill at the Tuolumne River
Livneh_Spill_Tuo <- read_csv("C:/Users/gusta/Desktop/PhD/CERCWET/Hydropower/tuolumne/Redflags/Livneh/IFR_redflag.csv", col_types = cols(node = col_date("%m/%d/%Y"), .default = col_double()))[-c(1:3), c(1,14)]%>%
    filter(between(node, as.Date("1980-10-01"), as.Date("2010-09-30")))
#CanESM2_RedFlag_Mer$CanESM2 <- rowSums(CanESM2_RedFlag_Mer[c(2,3)])
Livneh_Spill_Tuo$GCM <- as.factor("Livneh (Historical)")
Livneh_Spill_Tuo


CanESM2_Spill_Tuo <- read_csv("C:/Users/gusta/Desktop/PhD/CERCWET/Hydropower/tuolumne/Redflags/CanESM2_rcp85/IFR_redflag.csv", col_types = cols(node = col_date("%m/%d/%Y"), .default = col_double()))[-c(1:3), c(1,14)]
#CanESM2_RedFlag_Mer$CanESM2 <- rowSums(CanESM2_RedFlag_Mer[c(2,3)])
CanESM2_Spill_Tuo$GCM <- as.factor("CanESM2")
CanESM2_Spill_Tuo

CNRM_Spill_Tuo <- read_csv("C:/Users/gusta/Desktop/PhD/CERCWET/Hydropower/tuolumne/Redflags/CNRM-CM5_rcp85/IFR_redflag.csv", col_types = cols(node = col_date("%m/%d/%Y"), .default = col_double()))[-c(1:3), c(1,14)] 
#CNRM_RedFlag_Mer$`CNRM-CM5` <- rowSums(CNRM_RedFlag_Mer[c(2,3)])
CNRM_Spill_Tuo$GCM <- as.factor("CNRM-CM5")
CNRM_Spill_Tuo


HadGEM2_Spill_Tuo <- read_csv("C:/Users/gusta/Desktop/PhD/CERCWET/Hydropower/tuolumne/Redflags/HadGEM2-ES_rcp85/IFR_redflag.csv", col_types = cols(node = col_date("%m/%d/%Y"), .default = col_double()))[-c(1:3), c(1,14)]
#HadGEM2_RedFlag_Mer$`HadGEM2-ES` <- rowSums(HadGEM2_RedFlag_Mer[c(2,3)])
HadGEM2_Spill_Tuo$GCM <- as.factor("HadGEM2-ES")
HadGEM2_Spill_Tuo 

MIROC5_Spill_Tuo <- read_csv("C:/Users/gusta/Desktop/PhD/CERCWET/Hydropower/tuolumne/Redflags/MIROC5_rcp85/IFR_redflag.csv", col_types = cols(node = col_date("%m/%d/%Y"), .default = col_double()))[-c(1:3), c(1,14)] 
#MIROC5_RedFlag_Mer$MIROC5 <- rowSums(MIROC5_RedFlag_Mer[c(2,3)])
MIROC5_Spill_Tuo$GCM <- as.factor("MIROC5")
MIROC5_Spill_Tuo

```

```{r}
Spill_Tuo1 <- rbind(Livneh_Spill_Tuo, CanESM2_Spill_Tuo, CNRM_Spill_Tuo, 
                 HadGEM2_Spill_Tuo, MIROC5_Spill_Tuo) %>% 
  mutate(WaterYear = as.factor(water_year(node, origin= "usgs")),
         Scenario = "RCP 8.5",
         Basin = "Tuolumne")
Spill_Tuo1

Spill_Tuo <- melt(Spill_Tuo1, id = c("node", "WaterYear", "GCM", "Scenario", "Basin"))
Spill_Tuo
```

```{r}
#Maximum flow at La Grande is 9000 cfs (777,600,000 cfs in a day) = 22.01918000300859 mcm
Spill_Tuo2 <- Spill_Tuo %>%
  mutate(MaxFlow = as.numeric(22.01918000300859),
         Spill = as.numeric(value - MaxFlow),
         Count = ifelse(Spill > 0, 1,0))

Spill_Tuo2
```

```{r}
#Spill from Don Pedro
Livneh_Spill_DonPedro <- read_csv("C:/Users/gusta/Desktop/PhD/CERCWET/Hydropower/tuolumne/Redflags/Livneh/IFR_redflag.csv", col_types = cols(node = col_date("%m/%d/%Y"), .default = col_double()))[-c(1:3), c(1,14)]%>%
    filter(between(node, as.Date("1980-10-01"), as.Date("2010-09-30")))
#CanESM2_RedFlag_Mer$CanESM2 <- rowSums(CanESM2_RedFlag_Mer[c(2,3)])
Livneh_Spill_DonPedro$GCM <- as.factor("Livneh (Historical)")
Livneh_Spill_DonPedro


CanESM2_Spill_DonPedro <- read_csv("C:/Users/gusta/Desktop/PhD/CERCWET/Hydropower/tuolumne/Redflags/CanESM2_rcp85/IFR_redflag.csv", col_types = cols(node = col_date("%m/%d/%Y"), .default = col_double()))[-c(1:3), c(1,14)]
#CanESM2_RedFlag_Mer$CanESM2 <- rowSums(CanESM2_RedFlag_Mer[c(2,3)])
CanESM2_Spill_DonPedro$GCM <- as.factor("CanESM2")
CanESM2_Spill_DonPedro

CNRM_Spill_DonPedro <- read_csv("C:/Users/gusta/Desktop/PhD/CERCWET/Hydropower/tuolumne/Redflags/CNRM-CM5_rcp85/IFR_redflag.csv", col_types = cols(node = col_date("%m/%d/%Y"), .default = col_double()))[-c(1:3), c(1,14)] 
#CNRM_RedFlag_Mer$`CNRM-CM5` <- rowSums(CNRM_RedFlag_Mer[c(2,3)])
CNRM_Spill_DonPedro$GCM <- as.factor("CNRM-CM5")
CNRM_Spill_DonPedro


HadGEM2_Spill_DonPedro <- read_csv("C:/Users/gusta/Desktop/PhD/CERCWET/Hydropower/tuolumne/Redflags/HadGEM2-ES_rcp85/IFR_redflag.csv", col_types = cols(node = col_date("%m/%d/%Y"), .default = col_double()))[-c(1:3), c(1,14)]
#HadGEM2_RedFlag_Mer$`HadGEM2-ES` <- rowSums(HadGEM2_RedFlag_Mer[c(2,3)])
HadGEM2_Spill_DonPedro$GCM <- as.factor("HadGEM2-ES")
HadGEM2_Spill_DonPedro 

MIROC5_Spill_DonPedro <- read_csv("C:/Users/gusta/Desktop/PhD/CERCWET/Hydropower/tuolumne/Redflags/MIROC5_rcp85/IFR_redflag.csv", col_types = cols(node = col_date("%m/%d/%Y"), .default = col_double()))[-c(1:3), c(1,14)] 
#MIROC5_RedFlag_Mer$MIROC5 <- rowSums(MIROC5_RedFlag_Mer[c(2,3)])
MIROC5_Spill_DonPedro$GCM <- as.factor("MIROC5")
MIROC5_Spill_DonPedro

```

```{r}
Spill_DonPedro1 <- rbind(Livneh_Spill_DonPedro, CanESM2_Spill_DonPedro, CNRM_Spill_DonPedro, 
                 HadGEM2_Spill_DonPedro, MIROC5_Spill_DonPedro) %>% 
  mutate(WaterYear = as.factor(water_year(node, origin= "usgs")),
         Scenario = "RCP 8.5",
         Basin = "Tuolumne")
Spill_DonPedro1

Spill_DonPedro <- melt(Spill_DonPedro1, id = c("node", "WaterYear", "GCM", "Scenario", "Basin")) %>%
  mutate(Count = ifelse(value > 0, 1,0))
Spill_DonPedro
```

```{r}
ggplot(Spill_DonPedro) +
     theme_bw(base_size=12, base_family='Times New Roman') +
  geom_bar(aes(x= WaterYear, y=Count), stat = "identity") +
scale_color_colorblind() +
   scale_y_continuous(expand = c(0, NA),
                      limits = c(0,366), breaks = c(0, 100, 200, 300, 366)) +
  scale_x_discrete(breaks = c("1981", "1986", "1991", "1996", "2001", "2006", 
                              "2031", "2036", "2041", "2046", "2051", "2056"),
                    expand = c(0, NA))+

  facet_wrap(~GCM, scales = "free_x", ncol = 1, strip.position = "left") +
  
  
  labs(title = "IFR at La Grange",
       subtitle = "RCP 8.5 Scenario",
       
  #  subtitle = "CanESM2",
    x = element_blank(),
       y = "Spill Occurrence (Days)") +
 theme(#legend.title = element_blank(),
    #legend.position = "bottom",
     #   legend.direction = "horizontal",
      #  legend.box.margin = margin(t = -17),
        plot.title = element_text(hjust = 0.5),
        plot.subtitle = element_text(hjust = 0.5),
        strip.placement = "outside",
        strip.background = element_blank())+
  png("DonPedro_Spill.png", units ="in", width=3.5, height=8, res = 300)
```

```{r}
Livneh_Spill_Stn <- read_csv("C:/Users/gusta/Desktop/PhD/CERCWET/Hydropower/stanislaus/Redflags/Livneh/IFR_redflag.csv", col_types = cols(node = col_date("%m/%d/%Y"), .default = col_double()))[-c(1:3), c(1,34)]%>%
    filter(between(node, as.Date("1980-10-01"), as.Date("2010-09-30")))
#CanESM2_RedFlag_Mer$CanESM2 <- rowSums(CanESM2_RedFlag_Mer[c(2,3)])
Livneh_Spill_Stn$GCM <- as.factor("Livneh (Historical)")
Livneh_Spill_Stn


CanESM2_Spill_Stn <- read_csv("C:/Users/gusta/Desktop/PhD/CERCWET/Hydropower/stanislaus/Redflags/CanESM2_rcp85/IFR_redflag.csv", col_types = cols(node = col_date("%m/%d/%Y"), .default = col_double()))[-c(1:3), c(1, 34)]
#CanESM2_RedFlag_Mer$CanESM2 <- rowSums(CanESM2_RedFlag_Mer[c(2,3)])
CanESM2_Spill_Stn$GCM <- as.factor("CanESM2")
CanESM2_Spill_Stn

CNRM_Spill_Stn <- read_csv("C:/Users/gusta/Desktop/PhD/CERCWET/Hydropower/stanislaus/Redflags/CNRM-CM5_rcp85/IFR_redflag.csv", col_types = cols(node = col_date("%m/%d/%Y"), .default = col_double()))[-c(1:3), c(1, 34)]
#CNRM_RedFlag_Mer$`CNRM-CM5` <- rowSums(CNRM_RedFlag_Mer[c(2,3)])
CNRM_Spill_Stn$GCM <- as.factor("CNRM-CM5")
CNRM_Spill_Stn


HadGEM2_Spill_Stn <- read_csv("C:/Users/gusta/Desktop/PhD/CERCWET/Hydropower/stanislaus/Redflags/HadGEM2-ES_rcp85/IFR_redflag.csv", col_types = cols(node = col_date("%m/%d/%Y"), .default = col_double()))[-c(1:3), c(1, 34)]
#HadGEM2_RedFlag_Mer$`HadGEM2-ES` <- rowSums(HadGEM2_RedFlag_Mer[c(2,3)])
HadGEM2_Spill_Stn$GCM <- as.factor("HadGEM2-ES")
HadGEM2_Spill_Stn 

MIROC5_Spill_Stn <- read_csv("C:/Users/gusta/Desktop/PhD/CERCWET/Hydropower/stanislaus/Redflags/MIROC5_rcp85/IFR_redflag.csv", col_types = cols(node = col_date("%m/%d/%Y"), .default = col_double()))[-c(1:3), c(1,34)]
#MIROC5_RedFlag_Mer$MIROC5 <- rowSums(MIROC5_RedFlag_Mer[c(2,3)])
MIROC5_Spill_Stn$GCM <- as.factor("MIROC5")
MIROC5_Spill_Stn
```

```{r}
Spill_Stn1 <- rbind(Livneh_Spill_Stn, CanESM2_Spill_Stn, CNRM_Spill_Stn, 
                 HadGEM2_Spill_Stn, MIROC5_Spill_Stn) %>% 
  mutate(WaterYear = as.factor(water_year(node, origin= "usgs")),
         Scenario = "RCP 8.5",
         Basin = "Stanislaus")
Spill_Stn1

Spill_Stn <- melt(Spill_Stn1, id = c("node", "WaterYear", "GCM", "Scenario", "Basin"))
Spill_Stn
```

```{r}
#Maximum flow at Goodwin Reservoir is 8000 cfs (691,200,000 cfs in a day) = 19.572604447118746 mcm
Spill_Stn2 <- Spill_Stn %>%
  mutate(MaxFlow = as.numeric(19.572604447118746),
         Spill = as.numeric(value - MaxFlow),
         Count = ifelse(Spill > 0, 1,0))

Spill_Stn2
```

```{r}
Livneh_Spill_NewMlones <- read_csv("C:/Users/gusta/Desktop/PhD/CERCWET/Hydropower/stanislaus/Redflags/Livneh/IFR_redflag.csv", col_types = cols(node = col_date("%m/%d/%Y"), .default = col_double()))[-c(1:3), c(1,10)]%>%
    filter(between(node, as.Date("1980-10-01"), as.Date("2010-09-30")))
#CanESM2_RedFlag_Mer$CanESM2 <- rowSums(CanESM2_RedFlag_Mer[c(2,3)])
Livneh_Spill_NewMlones$GCM <- as.factor("Livneh (Historical)")
Livneh_Spill_NewMlones


CanESM2_Spill_NewMelones <- read_csv("C:/Users/gusta/Desktop/PhD/CERCWET/Hydropower/stanislaus/Redflags/CanESM2_rcp85/IFR_redflag.csv", col_types = cols(node = col_date("%m/%d/%Y"), .default = col_double()))[-c(1:3), c(1, 10)]
#CanESM2_RedFlag_Mer$CanESM2 <- rowSums(CanESM2_RedFlag_Mer[c(2,3)])
CanESM2_Spill_NewMelones$GCM <- as.factor("CanESM2")
CanESM2_Spill_NewMelones

CNRM_Spill_NewMelones <- read_csv("C:/Users/gusta/Desktop/PhD/CERCWET/Hydropower/stanislaus/Redflags/CNRM-CM5_rcp85/IFR_redflag.csv", col_types = cols(node = col_date("%m/%d/%Y"), .default = col_double()))[-c(1:3), c(1, 10)]
#CNRM_RedFlag_Mer$`CNRM-CM5` <- rowSums(CNRM_RedFlag_Mer[c(2,3)])
CNRM_Spill_NewMelones$GCM <- as.factor("CNRM-CM5")
CNRM_Spill_NewMelones


HadGEM2_Spill_NewMelones <- read_csv("C:/Users/gusta/Desktop/PhD/CERCWET/Hydropower/stanislaus/Redflags/HadGEM2-ES_rcp85/IFR_redflag.csv", col_types = cols(node = col_date("%m/%d/%Y"), .default = col_double()))[-c(1:3), c(1, 10)]
#HadGEM2_RedFlag_Mer$`HadGEM2-ES` <- rowSums(HadGEM2_RedFlag_Mer[c(2,3)])
HadGEM2_Spill_NewMelones$GCM <- as.factor("HadGEM2-ES")
HadGEM2_Spill_NewMelones 

MIROC5_Spill_NewMelones <- read_csv("C:/Users/gusta/Desktop/PhD/CERCWET/Hydropower/stanislaus/Redflags/MIROC5_rcp85/IFR_redflag.csv", col_types = cols(node = col_date("%m/%d/%Y"), .default = col_double()))[-c(1:3), c(1,10)]
#MIROC5_RedFlag_Mer$MIROC5 <- rowSums(MIROC5_RedFlag_Mer[c(2,3)])
MIROC5_Spill_NewMelones$GCM <- as.factor("MIROC5")
MIROC5_Spill_NewMelones
```


```{r}
Spill_NewMelones1 <- rbind(Livneh_Spill_NewMlones, CanESM2_Spill_NewMelones, CNRM_Spill_NewMelones, 
                 HadGEM2_Spill_NewMelones, MIROC5_Spill_NewMelones) %>% 
  mutate(WaterYear = as.factor(water_year(node, origin= "usgs")),
         Scenario = "RCP 8.5",
         Basin = "Stanislaus")
Spill_NewMelones1

Spill_NewMelones <- melt(Spill_NewMelones1, id = c("node", "WaterYear", "GCM", "Scenario", "Basin")) %>%
  mutate(Count = ifelse(value > 0, 1,0))
Spill_NewMelones
```

```{r}
ggplot(Spill_DonPedro) +
     theme_bw(base_size=12, base_family='Times New Roman') +
  geom_bar(aes(x= WaterYear, y=Count), stat = "identity") +
scale_color_colorblind() +
   scale_y_continuous(expand = c(0, NA),
                      limits = c(0,366), breaks = c(0, 100, 200, 300, 366)) +
  scale_x_discrete(breaks = c("1981", "1986", "1991", "1996", "2001", "2006", 
                              "2031", "2036", "2041", "2046", "2051", "2056"),
                    expand = c(0, NA))+

  facet_wrap(~GCM, scales = "free_x", ncol = 1, strip.position = "left") +
  
  
  labs(title = "New Melones Lake Flood Control",
       subtitle = "RCP 8.5 Scenario",
       
  #  subtitle = "CanESM2",
    x = element_blank(),
       y = "Spill Occurrence (Days)") +
 theme(#legend.title = element_blank(),
    #legend.position = "bottom",
     #   legend.direction = "horizontal",
      #  legend.box.margin = margin(t = -17),
        plot.title = element_text(hjust = 0.5),
        plot.subtitle = element_text(hjust = 0.5),
        strip.placement = "outside",
        strip.background = element_blank())+
  png("LaGrange_Spill.png", units ="in", width=3.5, height=8, res = 300)
```

```{r}
Livneh_Spill_USJ <- read_csv("C:/Users/gusta/Box/VICE Lab/RESEARCH/PROJECTS/CERC-WET/Task7_San_Joaquin_Model/Pywr models/results/Binary IFRs x Prices/upper_san_joaquin/historical/Livneh/InstreamFlowRequirement_Flow_mcm.csv", col_types = cols(node = col_date(), .default = col_double()))[-c(1:3), c(1,30)]%>%
    filter(between(node, as.Date("1980-10-01"), as.Date("2010-09-30")))
#CanESM2_RedFlag_Mer$CanESM2 <- rowSums(CanESM2_RedFlag_Mer[c(2,3)])
Livneh_Spill_USJ$GCM <- as.factor("Livneh (Historical)")
Livneh_Spill_USJ


CanESM2_Spill_USJ <- read_csv("C:/Users/gusta/Box/VICE Lab/RESEARCH/PROJECTS/CERC-WET/Task7_San_Joaquin_Model/Pywr models/results/Binary IFRs x Prices/upper_san_joaquin/gcms/CanESM2_rcp85/InstreamFlowRequirement_Flow_mcm.csv", col_types = cols(node = col_date(), .default = col_double()))[-c(1:3), c(1, 30)]
#CanESM2_RedFlag_Mer$CanESM2 <- rowSums(CanESM2_RedFlag_Mer[c(2,3)])
CanESM2_Spill_USJ$GCM <- as.factor("CanESM2")
CanESM2_Spill_USJ

CNRM_Spill_USJ <- read_csv("C:/Users/gusta/Box/VICE Lab/RESEARCH/PROJECTS/CERC-WET/Task7_San_Joaquin_Model/Pywr models/results/Binary IFRs x Prices/upper_san_joaquin/gcms/CNRM-CM5_rcp85/InstreamFlowRequirement_Flow_mcm.csv", col_types = cols(node = col_date(), .default = col_double()))[-c(1:3), c(1, 30)]
#CNRM_RedFlag_Mer$`CNRM-CM5` <- rowSums(CNRM_RedFlag_Mer[c(2,3)])
CNRM_Spill_USJ$GCM <- as.factor("CNRM-CM5")
CNRM_Spill_USJ


HadGEM2_Spill_USJ <- read_csv("C:/Users/gusta/Box/VICE Lab/RESEARCH/PROJECTS/CERC-WET/Task7_San_Joaquin_Model/Pywr models/results/Binary IFRs x Prices/upper_san_joaquin/gcms/HadGEM2-ES_rcp85/InstreamFlowRequirement_Flow_mcm.csv", col_types = cols(node = col_date(), .default = col_double()))[-c(1:3), c(1, 30)]
#HadGEM2_RedFlag_Mer$`HadGEM2-ES` <- rowSums(HadGEM2_RedFlag_Mer[c(2,3)])
HadGEM2_Spill_USJ$GCM <- as.factor("HadGEM2-ES")
HadGEM2_Spill_USJ 

MIROC5_Spill_USJ <- read_csv("C:/Users/gusta/Box/VICE Lab/RESEARCH/PROJECTS/CERC-WET/Task7_San_Joaquin_Model/Pywr models/results/Binary IFRs x Prices/upper_san_joaquin/gcms/MIROC5_rcp85/InstreamFlowRequirement_Flow_mcm.csv", col_types = cols(node = col_date(), .default = col_double()))[-c(1:3), c(1,30)]
#MIROC5_RedFlag_Mer$MIROC5 <- rowSums(MIROC5_RedFlag_Mer[c(2,3)])
MIROC5_Spill_USJ$GCM <- as.factor("MIROC5")
MIROC5_Spill_USJ
```

```{r}
Spill_USJ1 <- rbind(Livneh_Spill_USJ, CanESM2_Spill_USJ, CNRM_Spill_USJ, 
                 HadGEM2_Spill_USJ, MIROC5_Spill_USJ) %>% 
  mutate(WaterYear = as.factor(water_year(node, origin= "usgs")),
         Scenario = "RCP 8.5",
         Basin = "Upper San Joaquin")
Spill_USJ1

Spill_USJ <- melt(Spill_USJ1, id = c("node", "WaterYear", "GCM", "Scenario", "Basin"))
Spill_USJ
```

```{r}
#Maximum flow bl Millerton Lake is 8000 cfs (691,200,000 cfs in a day) = 19.572604447118746 mcm
Spill_USJ2 <- Spill_USJ %>%
  mutate(MaxFlow = as.numeric(19.572604447118746),
         Spill = as.numeric(value - MaxFlow),
         Count = ifelse(Spill > 0, 1,0))

Spill_USJ2
```

```{r}
Total_Spill <- rbind(Spill_Mer2, Spill_Stn2, Spill_Tuo2, Spill_USJ2)
Total_Spill
```

```{r}
ggplot(Total_Spill) +
     theme_bw(base_size=12, base_family='Times New Roman') +
  geom_bar(aes(x= WaterYear, y=Count), stat = "identity") +
scale_color_colorblind() +
   scale_y_continuous(expand = c(0, NA),
                      limits = c(0,160), breaks = c(0, 40, 80, 120, 160)) +
  scale_x_discrete(breaks = c("1981", "1986", "1991", "1996", "2001", "2006", 
                              "2031", "2036", "2041", "2046", "2051", "2056"),
                    expand = c(0, NA))+

  facet_wrap(GCM~variable, scales = "free_x", ncol = 4) +
  
  labs(#title = "Storage Below One Third of Reservoir Operating Capacity",
       #subtitle = "RCP 8.5 Scenario",
       
  #  subtitle = "CanESM2",
    x = element_blank(),
       y = "Spill occurrence (Days)") +
 theme(#strip.text.x = element_blank(),
       #legend.title = element_blank(),
    #legend.position = "bottom",
     #   legend.direction = "horizontal",
      #  legend.box.margin = margin(t = -17),
        plot.title = element_text(hjust = 0.5),
        plot.subtitle = element_text(hjust = 0.5),
        strip.placement = "outside",
        strip.background = element_blank())+
       # panel.spacing.y = unit(0, "lines")) +
png("AllBasins_Spill_FacetWrap.png", units ="in", width=8, height=8, res = 300)
```


```{r}
ggplot(Total_Spill) +
     theme_bw(base_size=12, base_family='Times New Roman') +
  geom_bar(aes(x= WaterYear, y=Count), stat = "identity") +
scale_color_colorblind() +
   scale_y_continuous(expand = c(0, NA),
                      limits = c(0,160), breaks = c(0, 40, 80, 120, 160)) +
  scale_x_discrete(breaks = c("1981", "1991", "2001", 
                              "2031", "2041","2051"),
                    expand = c(0, NA))+

  facet_grid(variable~GCM, scales = "free_x") +
  
  labs(#title = "Storage Below One Third of Reservoir Operating Capacity",
       #subtitle = "RCP 8.5 Scenario",
       
  #  subtitle = "CanESM2",
    x = element_blank(),
       y = "Spill occurrence (Days)") +
 theme(#strip.text.x = element_blank(),
       #legend.title = element_blank(),
    #legend.position = "bottom",
     #   legend.direction = "horizontal",
      #  legend.box.margin = margin(t = -17),
        plot.title = element_text(hjust = 0.5),
        plot.subtitle = element_text(hjust = 0.5),
        strip.placement = "outside",
        strip.background = element_blank())+
       # panel.spacing.y = unit(0, "lines")) +
png("AllBasins_Spill_FacetGrid.png", units ="in", width=8, height=7, res = 300)
```


```{r}
#Storage Capacity for STN
Total_Spill%>%
  filter(variable == "IFR bl Millerton Lake") %>%
ggplot() +
     theme_bw(base_size=12, base_family='Times New Roman') +
  geom_bar(aes(x= WaterYear, y=Count), stat = "identity") +
scale_color_colorblind() +
   scale_y_continuous(expand = c(0, NA),
                      limits = c(0,160), breaks = c(0, 40, 80, 120, 160)) +
  scale_x_discrete(breaks = c("1981", "1986", "1991", "1996", "2001", "2006", 
                              "2031", "2036", "2041", "2046", "2051", "2056"),
                    expand = c(0, NA))+

  facet_wrap(~GCM, scales = "free_x", ncol = 1, strip.position = "left") +
  
  
  labs(title = "IFR bl Millerton Lake",
       subtitle = "RCP 8.5 Scenario",
       
  #  subtitle = "CanESM2",
    x = element_blank(),
       y = "Spill Occurrence (Days)") +
 theme(#legend.title = element_blank(),
    #legend.position = "bottom",
     #   legend.direction = "horizontal",
      #  legend.box.margin = margin(t = -17),
        plot.title = element_text(hjust = 0.5),
        plot.subtitle = element_text(hjust = 0.5),
        strip.placement = "outside",
        strip.background = element_blank())+
#        panel.spacing.y = unit(0, "lines")) +
png("USJ_Spill.png", units ="in", width=3.5, height=8, res = 300)
```