My-Code/Getting data for BC in the USJ.Rmd at main · GustavoFacincani/My-Code · GitHub

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
---
title: "Bias Correction USJ"
author: "Gustavo Facincani Dourado"
date: "9/17/2020"
output: html_document
---

The Big Creek system is composed of subbasins 31-38, which discharge into Mammoth Pool Reservoir (subbasin 21), together with subbasins 22 and 23 (both gauged). Subbasins 32-35 have a very poor record of USGS gauge data. There are many mismatches on the time that the gauges in the region present historical records. So, here we're collecting data as of an uncontrolled system, according to the data available for each subbasin, to bias correct the Livneh data for these regions per subbasin individually, as the VIC runoff is not well represented in this part of the Upper San Joaquin river system.


```{r}
library(readr)
library(tidyverse)
library(hyfo)
library(lubridate)
```

```{r}
#Reading Livneh data, that is going to be used for bias correction and filtering them to have objects with the same length of the observed data, to be used as the hindcast vs observed
path1 <- "C:/Users/gusta/Box/VICE Lab/RESEARCH/PROJECTS/CERC-WET/Task7_San_Joaquin_Model/pywr_models/data/Upper San Joaquin River/hydrology/historical/Livneh/runoff/"

Sub22_Livneh <- read_csv(paste(path1,"tot_runoff_sb22_mcm.csv",sep="")) %>%
  rename(Inflow = tot_runoff_sb22)
Sub22_PreBC <- Sub22_Livneh %>% #Granite Creek
   filter(between(Date, as.Date("1965-10-01"), as.Date("1986-09-29")))

Sub23_Livneh <- read_csv(paste(path1,"tot_runoff_sb23_mcm.csv",sep=""))%>%
  rename(Inflow = tot_runoff_sb23)
Sub23_PreBC <- Sub23_Livneh%>% #San Joaquin River
    filter(between(Date, as.Date("1951-10-01"), as.Date("1991-09-29")))


Sub31_Livneh <- read_csv(paste(path1,"tot_runoff_sb31_mcm.csv",sep=""))%>%
  rename(Inflow = tot_runoff_sb31)
Sub31_PreBC <-  Sub31_Livneh %>% #Thomas Edison Lake
  filter(between(Date, as.Date("1954-10-12"), as.Date("2013-12-31")))


Sub36_Livneh <- read_csv(paste(path1,"tot_runoff_sb36_mcm.csv",sep=""))%>%
  rename(Inflow = tot_runoff_sb36)
Sub36_PreBC <- Sub36_Livneh# %>% #Bear Creek
   #filter(between(Date, as.Date("1950-10-01"), as.Date("2013-12-31")))

Sub37_Livneh <- read_csv(paste(path1,"tot_runoff_sb37_mcm.csv",sep=""))%>%
  rename(Inflow = tot_runoff_sb37)
Sub37_PreBC <-  Sub37_Livneh %>% #SF San Joaquin River
  filter(between(Date, as.Date("2009-10-01"), as.Date("2013-12-31")))


Sub38_Livneh <- read_csv(paste(path1,"tot_runoff_sb38_mcm.csv",sep=""))%>%
  rename(Inflow = tot_runoff_sb38)
Sub38_PreBC <- Sub38_Livneh %>% #Florence Lake
  filter(between(Date, as.Date("1950-01-01"), as.Date("1980-09-29")))


```

```{r}
#Read the gauge data for the same period of the data we have for Livneh
USJpath <- "C:/Users/gusta/Box/VICE Lab/RESEARCH/PROJECTS/CERC-WET/Task7_San_Joaquin_Model/pywr_models/data/Upper San Joaquin River/"

#Subbasins 31 and 38 have Florence Lake, Thomas A. Edison Lake and Mammoth Pool Reservoir respectively
Reservoir_stor <- read_csv(paste0(USJpath, "gauges/storage_mcm.csv", sep = ""), col_types = cols(Date = col_date(), .default = col_double()))#%>% #reading data column as data, others as numeric
#filter(between(Date, as.Date("1950-01-01"), as.Date("2013-12-31"))) #selecting period of available gauged data
Reservoir_stor

#The other subbasins have only streamflow
Reservoir_flow <- read_csv(paste0(USJpath, "gauges/streamflow_cfs.csv", sep = ""), col_types = cols(Date = col_date(), .default = col_double()))#%>% #reading data column as data, others as numeric
#filter(between(Date, as.Date("1950-01-01"), as.Date("2013-12-31"))) #selecting period of available gauged data

Reservoir_flow[2:49] <- (Reservoir_flow[2:49]/35.314666)*86400/1000000 #as reservoir data is in million cubic meters (mcm), let's convert cfs to mcm
Reservoir_flow
```

```{r}
#Water balance of reservoirs to determine inflows
#Subbasin 38 - Florence
Florence_outflow <- Reservoir_flow[c(1,4,5)] %>% #gauges 11229500 and USGS 11230000, outflow of Florence
  filter(between(Date, as.Date("1950-01-01"), as.Date("1980-09-29")))

Sub38 <- Reservoir_stor[c(1,8)] %>% #selecting gauge with storage data
filter(between(Date, as.Date("1950-01-01"), as.Date("1980-09-29"))) %>%
   mutate(ChangeStorage = `USGS 11229600 FLORENCE LK NR BIG C CA` - lag(`USGS 11229600 FLORENCE LK NR BIG C CA`, default = first(`USGS 11229600 FLORENCE LK NR BIG C CA`)), #getting daily change in storage
          Outflow = rowSums(Florence_outflow[2:3], na.rm=T), #sum of USGS 11229500 and USGS 11230000, which correspond to flow below the reservoir and Ward Tunnel intake, respectively
          flow = ifelse((ChangeStorage + Outflow) < 0, 0, (ChangeStorage + Outflow))) %>% #getting inflow, in which inflow = change in storage (storage of the day minus the storage of the previous day + outflow of the day
           dplyr::select(Date, flow)#remove the already used column
Sub38

#Subbasins 31 - Thomas Edison
Thomas_outflow <- Reservoir_flow[c(1,14)] %>%
  filter(between(Date, as.Date("1954-10-12"), as.Date("2013-12-31")))

Sub31 <- Reservoir_stor[c(1,6)] %>% #selecting gauge with storage data
  filter(between(Date, as.Date("1954-10-12"), as.Date("2013-12-31"))) %>%
   mutate(ChangeStorage = `USGS 11231000 LK THOMAS A. EDISON NR BIG C CA` - lag(`USGS 11231000 LK THOMAS A. EDISON NR BIG C CA`, default = first(`USGS 11231000 LK THOMAS A. EDISON NR BIG C CA`)), #getting daily change in storage
          Outflow = Thomas_outflow$`USGS 11231500 MONO C BL LK THOMAS A EDISON CA`, #getting outflow
          flow = ifelse((ChangeStorage + Outflow) < 0, 0, (ChangeStorage + Outflow))) %>%
           dplyr::select(Date, flow)#remove the already used column
Sub31

```


```{r}
#Subbasin 22
Sub22 <- Reservoir_flow[c(1,3)] %>% #Granite Creek
   filter(between(Date, as.Date("1965-10-01"), as.Date("1986-09-29"))) %>%
  rename(flow = 2)
Sub22

#Subbasin 23
Sub23 <- Reservoir_flow[c(1,2)]%>% #San Joaquin River
    filter(between(Date, as.Date("1951-10-01"), as.Date("1991-09-29")))%>%
  rename(flow = 2)
Sub23

#Subbasin 36
Sub36 <- Reservoir_flow[c(1,8)]%>% #Bear Creek
   filter(between(Date, as.Date("1950-01-01"), as.Date("2013-12-31")))%>%
  rename(flow = 2)
Sub36

#Subbasin 37
#If we consider these gauges, around 24% of the days will have negative values

#Sub37 <- Reservoir_flow[c(1,5,7)] %>% #SF San Joaquin, with inflow from Florence Lake
#   filter(between(Date, as.Date("1975-10-01"), as.Date("1980-09-29"))) %>%
#  mutate(Inflow = (`USGS 11230215 SF SAN JOAQUIN R BL HOOPER C NR FLORENCE LAKE CA` - `USGS 11230000 SF SAN JOAQUIN R NR FLORENCE LK CA`))
#Sub37

#Considering this slightly shorter record, we only have 4 days with negative flow
Sub37 <- Reservoir_flow[c(1,6,7)] %>% #SF San Joaquin, with inflow from Florence Lake
   filter(between(Date, as.Date("2009-10-01"), as.Date("2013-12-31"))) %>% #but this is actually a very short period of time for bias correction... should we really do it? Probably not
  mutate(flow = ifelse((`USGS 11230215 SF SAN JOAQUIN R BL HOOPER C NR FLORENCE LAKE CA` - `USGS 11230070 SF SAN JOAQUIN R AB HOOPER C NR FLORENCE LAKE CA`) < 0,0, (`USGS 11230215 SF SAN JOAQUIN R BL HOOPER C NR FLORENCE LAKE CA` - `USGS 11230070 SF SAN JOAQUIN R AB HOOPER C NR FLORENCE LAKE CA`)))%>%
           dplyr::select(Date, flow)#remove the already used column
Sub37

```