r params$region
"library(tidyverse) library(tmap) library(bcmaps) library(tidyhydat) library(lubridate) library(sf) library(rmapshaper) library(kableExtra) library(knitr) library(RcppRoll)
## Knitr options dependent on type of document if(params$table_format == "latex"){ opts_chunk$set(echo = FALSE, warning = FALSE, message = FALSE, fig.width = 8, fig.height = 8, fig.path = file.path("report/regional_streamflow/")) } if(params$table_format == "html"){ opts_chunk$set(echo = FALSE, warning = FALSE, message = FALSE, fig.width = 8, fig.height = 8) } options(knitr.table.format = params$table_format) #params <- tibble(region = "Skeena Natural Resource Region")
## Calculate percentiles from tidyhydat extracted data ## Limited to current day calculation calc_percentiles <- function(historical_flow, realtime_data) { if(class(realtime_data$Date) != "Date") stop("Date column is not in date class") historical_flow %>% filter(yday(Date) == yday(Sys.Date())) %>% group_by(STATION_NUMBER) %>% nest() %>% left_join(realtime_data, by = c("STATION_NUMBER")) %>% mutate(prctile = map2_dbl(data, Value, ~ecdf(.x$Value)(.y))) %>% left_join(allstations, by = c("STATION_NUMBER")) %>% mutate(pct_bin = case_when( is.na(prctile) ~ "Not ranked", prctile >= 0 & prctile <= 0.01 ~ "Low", prctile > 0.01 & prctile <= 0.10 ~ "Much below normal (<10)", prctile > 0.10 & prctile <= 0.24 ~ "Below Normal (10-24)", prctile > 0.24 & prctile <= 0.75 ~ "Normal (25-75)", prctile > 0.75 & prctile <= 0.90 ~ "Above normal (76-90)", prctile > 0.90 & prctile < 1 ~ "Much above normal (>90)", prctile == 1 ~ "High" )) %>% mutate(pct_bin = factor(pct_bin, levels = expected)) %>% mutate(prctile = prctile*100) %>% st_as_sf(coords = c("LONGITUDE", "LATITUDE"), crs= 4326, agr= "constant") %>% transform_bc_albers() }
## For testing #params <- tibble(date = Sys.Date(),region = "West Coast Natural Resource Region") ## Load realtime stations and convert to spatial object stns_sf <- st_as_sf( realtime_stations(prov_terr_state_loc = "BC"), coords = c("LONGITUDE", "LATITUDE"), crs = 4326, agr = "constant" ) %>% transform_bc_albers()
### Spatial Cleaning ## Choose the correct region nr_regions_sub <- filter(nr_regions(), ORG_UNIT_NAME == params$region) ## Find the inters region_borders <- st_intersection(nr_regions_sub, bc_bound_hres()) region_rivers <- st_intersection(nr_regions_sub, watercourses_5M()) wsc_drainages <- ms_simplify(wsc_drainages()) region_watershed <- st_intersection(select(nr_regions_sub, SHAPE), wsc_drainages()) %>% group_by(SUB_DRAINAGE_AREA_NAME) %>% summarise() %>% ms_simplify() stns_region <- st_join(stns_sf, nr_regions_sub, left = FALSE) city <- bc_cities() %>% st_join(region_borders, left = FALSE) %>% filter(POP_2000 == max(POP_2000)) %>% rename(geometry = SHAPE) %>% mutate( CENTROID = map(geometry, st_centroid), COORDS = map(CENTROID, st_coordinates), COORDS_X = map_dbl(COORDS, 1), COORDS_Y = map_dbl(COORDS, 2) ) %>% as_tibble() %>% st_sf()
## All stations from the specified region, that are flow and ## have a record longer than 20 years q_stns <- unique(stns_region$STATION_NUMBER) %>% hy_stn_data_range() %>% filter(DATA_TYPE == "Q") %>% filter(RECORD_LENGTH >=20) %>% pull(STATION_NUMBER) ## Query realtime data rl_data <- realtime_dd(q_stns) ## Find most recent instantaneous discharge value rl_data_instant <- rl_data %>% filter(Parameter == "Flow") %>% group_by(STATION_NUMBER) %>% filter(Date == max(Date)) %>% select(STATION_NUMBER, Date, Value) %>% mutate(Date = as.Date(Date)) %>% filter(Date == Sys.Date()) %>% ## drop max values that aren't today ungroup() ## Queery historical data ## NOTE: Should this be done with rl_data_recent$STATION_NUMBER? hist_flow <- hy_daily_flows(q_stns)
This report was generated on
cat(format(Sys.Date(), '%d %B, %Y'))
The information is taken from real time data recorded at WSC stations and is based on percentiles of long term historical flows for that day. A percentile is a value on a scale of one hundred that indicates the percent of a distribution that is equal to or below it. In general:
Below are the maps with percentile data for WSC stations:
## Find the average of the last 24 hours rl_data_last24 <- rl_data %>% filter(Parameter == "Flow") %>% group_by(STATION_NUMBER) %>% filter(Date >= Sys.time() - 60*60*24) %>% ## all data from last 24 hours select(STATION_NUMBER, Date, Value) %>% mutate(Date = Sys.Date()) %>% ## label last twenty four hours as from today group_by(Date, STATION_NUMBER) %>% summarise(Value = mean(Value)) %>% ungroup()
## Realtime 7 day average ## I took the realtime daily averages rl_data_7day_mean <- rl_data %>% filter(Parameter == "Flow") %>% mutate(Date = as.Date(Date)) %>% group_by(Date, STATION_NUMBER) %>% summarise(Value = mean(Value, na.rm = TRUE)) %>% mutate(Value = roll_mean(Value, n = 7,na.rm = TRUE, fill = "right")) %>% #filter(Date == Sys.Date()) %>% ungroup() hist_flow_7day_mean <- hist_flow %>% group_by(STATION_NUMBER) %>% mutate(Value = roll_mean(Value, n = 7, na.rm = TRUE, fill = "right")) %>% ungroup()
## Expected NAWW percentile bins expected <- c("Not ranked", "Low","Much below normal (<10)", "Below Normal (10-24)", "Normal (25-75)", "Above normal (76-90)","Much above normal (>90)", "High") ## Calculate instantaneous percentiles pct_flow_instant <- calc_percentiles(hist_flow, rl_data_instant) ## Calculate 24 hours percentiles pct_flow_last24 <- calc_percentiles(hist_flow, rl_data_last24) pct_flow_7day_mean <- calc_percentiles(hist_flow_7day_mean, rl_data_7day_mean)
naww_pal <- c("#9c0ad1","#F90010","#B31F23","#FEA116","#08FA09","#46DED2","#0003F6","#000000") ## Create a named vector named_naww_pal <- setNames(naww_pal, expected) tm_shape(region_watershed) + tm_polygons(col = "grey90") + tm_shape(region_rivers) + tm_lines(col = "blue") + tm_shape(pct_flow_last24) + tm_dots(col = "pct_bin", size = 0.4, title = "Streamflow conditions \n based percentile quantities", palette = naww_pal) + tm_shape(city) + tm_dots(col = "black", size = 0.4, shape = 17) + tm_text("NAME", auto.placement = TRUE) + tm_style_white(main.title = paste0(params$region))
The following are comparisons of two percentile summaries:
These are evaluated separately in this document.
num_year_data <- hy_stn_data_range(unique(pct_flow_instant$STATION_NUMBER)) %>% filter(DATA_TYPE == "Q") %>% select(STATION_NUMBER, RECORD_LENGTH) %>% rename(`Record Length` = RECORD_LENGTH) ## Grab only the latest flow and merge into one data frame pct_flow_instant_tbl_data <- pct_flow_instant %>% st_set_geometry(NULL) %>% select(STATION_NAME, STATION_NUMBER, Value, prctile) %>% rename(`%tile-instant` = prctile, `Latest Q` = Value) %>% select(-`%tile-instant`) %>% ##remove instant %tile left_join(pct_flow_last24 %>% st_set_geometry(NULL) %>% select(STATION_NUMBER, prctile) %>% rename(`%tile-last24` = prctile)) %>% left_join(pct_flow_7day_mean %>% st_set_geometry(NULL) %>% filter(Date == Sys.Date()) %>% select(STATION_NUMBER, prctile) %>% rename(`%tile-7day_mean` = prctile)) %>% left_join(num_year_data, by = c("STATION_NUMBER")) %>% arrange(STATION_NUMBER) if(params$table_format == "latex"){ pct_flow_instant_tbl_data %>% kable(format = "latex", booktabs = T, digits = 1) %>% kable_styling(font_size = 5,latex_options = c("HOLD_position")) } if(params$table_format == "html"){ pct_flow_instant_tbl_data %>% kable(format = "html", digits = 2) %>% kable_styling(bootstrap_options = c("striped", "hover", "condensed", "responsive")) }
If a location is not listed here there are no immediate low streamflow concerns.
much_below_normal_tbl <- pct_flow_instant %>% filter(pct_bin %in% c("Low","Much below normal (<10)")) %>% st_set_geometry(NULL) %>% select(STATION_NAME, STATION_NUMBER, Value, prctile) %>% rename(`Latest Q` = Value, `%tile` = prctile) %>% arrange(STATION_NUMBER) if(params$table_format == "latex"){ much_below_normal_tbl %>% kable(format = "latex", booktabs = T, digits = 2) %>% kable_styling(font_size = 7,latex_options = c("HOLD_position")) } if(params$table_format == "html"){ much_below_normal_tbl %>% kable(format = "html", digits = 2) %>% kable_styling(bootstrap_options = c("striped", "hover", "condensed", "responsive")) }
Flow needs during the summer will range from juvenile rearing flows near 10-20 % (preferred) MAD to adult salmon/char passage flows (>20%LT MAD) required prior to spawning. Sub-standard flows (<10%LT MAD) can affect fish populations by reducing the area and quality of riffle habitats that generate fish food and aeration. Flows nearing 5%LT MAD are considered suboptimal for fish rearing and migration, indicate a degradation of instream habitat and may require restrictions on water use.
## Calculate %MAD long_term_mad <- hy_annual_stats(q_stns) %>% filter(Sum_stat == "MEAN") %>% spread(Sum_stat, Value) %>% group_by(STATION_NUMBER) %>% summarise(`MEAN MAD (m^3/s)` = mean(MEAN, na.rm = TRUE)) %>% left_join(rl_data_instant, by = c("STATION_NUMBER")) %>% mutate(`% MAD` = (Value/`MEAN MAD (m^3/s)`)*100) %>% left_join(allstations[,1:2], by = c("STATION_NUMBER")) %>% select(STATION_NAME, STATION_NUMBER, Value, `MEAN MAD (m^3/s)`, `% MAD`) %>% rename(`Latest Q` = Value) %>% arrange(STATION_NUMBER) if(params$table_format == "latex"){ long_term_mad %>% kable(format = "latex", booktabs = T, digits = 2) %>% kable_styling(font_size = 7,latex_options = c("HOLD_position")) } if(params$table_format == "html"){ long_term_mad %>% kable(format = "html", digits = 2) %>% kable_styling(bootstrap_options = c("striped", "hover", "condensed", "responsive")) }
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.