Streamflow Analysis
Pablo E. GutiƩrrez-Fonseca
04 February, 2023
Step 1
Load libraries.
Cargar las librerias que necesitas.
libraries <- c('hydroTSM','dataRetrieval','data.table','dplyr','ggplot2','lubridate','raster','rgdal','waterData','zoo')
lapply(libraries, require, character.only = TRUE)## [[1]]
## [1] TRUE
##
## [[2]]
## [1] TRUE
##
## [[3]]
## [1] TRUE
##
## [[4]]
## [1] TRUE
##
## [[5]]
## [1] TRUE
##
## [[6]]
## [1] TRUE
##
## [[7]]
## [1] TRUE
##
## [[8]]
## [1] TRUE
##
## [[9]]
## [1] TRUE
##
## [[10]]
## [1] TRUElibrary(dataRetrieval): USGS package that gets streamflow data direct from the USGS website
Step 2
Data minning.
Try the code below with the site.code here, then
use the site code for your watershed.
site.code = "04290500" # The USGS streamgage code for Winnoski River at EssesStep 3
Download data from the USGS.
parameter.code = "00060" # this is the code for stream discharge.
start.date = "1929-01-01" # Blanks get all of the data
end.date = "2022-12-31"
# store data in "winnoski dataframe
winooski = readNWISdv("04290500", parameter.code, start.date, end.date)Step 4
Explore date.
head(winooski)## agency_cd site_no Date X_00060_00003 X_00060_00003_cd
## 1 USGS 04290500 1929-01-01 630 A
## 2 USGS 04290500 1929-01-02 540 A
## 3 USGS 04290500 1929-01-03 550 A
## 4 USGS 04290500 1929-01-04 570 A
## 5 USGS 04290500 1929-01-05 540 A
## 6 USGS 04290500 1929-01-06 480 Atail(winooski)## agency_cd site_no Date X_00060_00003 X_00060_00003_cd
## 34328 USGS 04290500 2022-12-26 3750 P
## 34329 USGS 04290500 2022-12-27 3280 P
## 34330 USGS 04290500 2022-12-28 2840 P
## 34331 USGS 04290500 2022-12-29 2610 P
## 34332 USGS 04290500 2022-12-30 2380 P
## 34333 USGS 04290500 2022-12-31 4270 PStep 5
The names for the discharge and QA columns arenāt very nice, so
rename them:
names(winooski)[c(4,5)] = c("Q.ft.s","QA")
head(winooski)## agency_cd site_no Date Q.ft.s QA
## 1 USGS 04290500 1929-01-01 630 A
## 2 USGS 04290500 1929-01-02 540 A
## 3 USGS 04290500 1929-01-03 550 A
## 4 USGS 04290500 1929-01-04 570 A
## 5 USGS 04290500 1929-01-05 540 A
## 6 USGS 04290500 1929-01-06 480 AStep 6
Explore data:
plot(winooski$Date,winooski$Q.ft.s, type='l')
Step 7
Summarize data and find extreme values.
summary(winooski)## agency_cd site_no Date Q.ft.s
## Length:34333 Length:34333 Min. :1929-01-01 Min. : 24
## Class :character Class :character 1st Qu.:1952-07-02 1st Qu.: 608
## Mode :character Mode :character Median :1976-01-01 Median : 1060
## Mean :1976-01-01 Mean : 1819
## 3rd Qu.:1999-07-02 3rd Qu.: 2050
## Max. :2022-12-31 Max. :41600
## QA
## Length:34333
## Class :character
## Mode :character
##
##
## winooski[which.max(winooski$Q.ft.s), ] # Identify the maximum flood events in the record## agency_cd site_no Date Q.ft.s QA
## 2635 USGS 04290500 1936-03-19 41600 Awinooski[which.min(winooski$Q.ft.s), ] # Identify the minimum discharge events in the record## agency_cd site_no Date Q.ft.s QA
## 14495 USGS 04290500 1968-09-07 24 AStep 8
Data visualization.
winooski.df <- winooski %>%
mutate(year = factor(year(Date)), # use year to define separate curves
Date = update(Date, year = 1)) # use a constant year for the x-axis
p <- ggplot(winooski.df, aes(Date, Q.ft.s, color = year)) +
scale_x_date(date_breaks = "1 month", date_labels = "%b") +
theme_bw()
p + geom_line() + theme(legend.position = "none") # Add the lines. 
Step 9
Data visualization.
p + geom_line(aes(group = year), color = "gray20", alpha = 0.1) +
geom_line(data = function(x) filter(x, year == 1936), color='blue' , size = 1) +
theme_bw()
Step 10
Contrasting 1936 with the rest of the years.
year.1936 <- setDT(winooski)[Date %between% c('1936-01-01', '1936-12-31')]
no.2023 <- setDT(winooski)[(Date < '1936-01-01' | Date > '1936-12-31'),]
summary(year.1936)## agency_cd site_no Date Q.ft.s
## Length:366 Length:366 Min. :1936-01-01 Min. : 116.0
## Class :character Class :character 1st Qu.:1936-04-01 1st Qu.: 534.8
## Mode :character Mode :character Median :1936-07-01 Median : 891.5
## Mean :1936-07-01 Mean : 2223.6
## 3rd Qu.:1936-09-30 3rd Qu.: 2192.5
## Max. :1936-12-31 Max. :41600.0
## QA
## Length:366
## Class :character
## Mode :character
##
##
## summary(no.2023)## agency_cd site_no Date Q.ft.s
## Length:33967 Length:33967 Min. :1929-01-01 Min. : 24
## Class :character Class :character 1st Qu.:1953-04-02 1st Qu.: 610
## Mode :character Mode :character Median :1976-07-02 Median : 1070
## Mean :1976-06-04 Mean : 1815
## 3rd Qu.:1999-10-01 3rd Qu.: 2050
## Max. :2022-12-31 Max. :29800
## QA
## Length:33967
## Class :character
## Mode :character
##
##
## Step 11.1
Calculate cumulative discharge for each year by first grouping by
water year.
cum.data <- addWaterYear(winooski) # Add a new column with Year
cumulative_dat <- group_by(cum.data, waterYear) %>% # Group by Year
mutate(cumulative_dis = cumsum(Q.ft.s), # Sum daily discharge (function cumsum())
wy_doy = seq(1:n())) # Add consecutive dayStep 11.1
First graph.
q <- ggplot(cumulative_dat, aes(x = wy_doy, y = cumulative_dis, group = waterYear)) +
geom_line(lwd = 0.6) +
xlab("Julian Day") + ylab("Cumulative dischage (ft.s)") +
ylim(c(0, 1300000)) +
xlim(0,366) +
theme_bw()
q
Step 11.2
Highlight the rainiest year
q + geom_line(data=subset(cumulative_dat, waterYear == "2011"), colour="blue", size=0.9)
Step 11.3
highlight the driest year
q + geom_line(data=subset(cumulative_dat, waterYear == "1965"), colour="red", size=0.9) 
References