In mkoohafkan/rremat: Russian River Estuary data collected by the Largier Lab

source(purl("rremat/vignettes/habitat-analysis.Rmd", tempfile()))

add some identifiers for scenarios

isfall = function(x){
  xx = as.numeric(format(x,"%m"))
  ifelse(xx > 8 & xx < 11, TRUE, FALSE)
}
data(inflows)
data(tides)
data(wll)

data(closures)
closures["duration"] = closures$end - closures$start
closures["isFall"] = isfall(closures$start)
closures["mean.flow"] = NA
closures["min.flow"] = NA
closures["max.flow"] = NA
closures["mean.tide"]= NA
closures["min.tide"] = NA
closures["max.tide"] = NA
closures["range.tide"] = NA
closures["mean.wll"]= NA
closures["min.wll"] = NA
closures["max.wll"] = NA
closures["range.wll"] = NA

for(i in seq(nrow(closures))){
  sd = closures[i, "start"]
  ed = closures[i, "end"]
  thisflow = filter(inflows, as.Date(datetime) >= sd & 
    as.Date(datetime) <= ed, gauge == "russian river")
  closures[i, "mean.flow"]= mean(thisflow$flow)
  closures[i, "min.flow"] = min(thisflow$flow)
  closures[i, "max.flow"] = max(thisflow$flow) 
  thistide = filter(tides, as.Date(datetime) >= sd & 
    as.Date(datetime) <= ed)  
  closures[i, "mean.tide"] = mean(thistide$height)
  closures[i, "min.tide"] = min(thistide$height)
  closures[i, "max.tide"] = max(thistide$height)
  closures[i, "range.tide"] = max(thistide$height) - min(thistide$height)
  thiswll = filter(wll, as.Date(mtime) >= sd & as.Date(mtime) <= ed, 
    site == "river mouth")  
  closures[i, "mean.wll"] = mean(thiswll$depth)
  closures[i, "min.wll"] = min(thiswll$depth)
  closures[i, "max.wll"] = max(thiswll$depth)
  closures[i, "range.wll"] = max(thiswll$depth) - min(thiswll$depth)
}

data(ctdmeta)
ctdmeta["isFall"] = isfall(ctdmeta$start)
ctdmeta["mean.flow"] = NA
ctdmeta["min.flow"] = NA
ctdmeta["max.flow"] = NA
ctdmeta["mean.tide"]= NA
ctdmeta["min.tide"] = NA
ctdmeta["max.tide"] = NA
ctdmeta["mean.wll"]= NA
ctdmeta["min.wll"] = NA
ctdmeta["max.wll"] = NA
for(i in seq(nrow(ctdmeta))){
  sd = ctdmeta[i, "start"]
  thisflow = filter(inflows, as.Date(datetime) == as.Date(sd), 
    gauge == "russian river")
  ctdmeta[i, "mean.flow"]= mean(thisflow$flow)
  ctdmeta[i, "min.flow"] = min(thisflow$flow)
  ctdmeta[i, "max.flow"] = max(thisflow$flow) 
  thistide = filter(tides, as.Date(datetime) == as.Date(sd))
  ctdmeta[i, "mean.tide"]= mean(thistide$height)
  ctdmeta[i, "min.tide"] = min(thistide$height)
  ctdmeta[i, "max.tide"] = max(thistide$height)   
}

# analysis of single closure (case 2)
thisclosure = closures[29,]

thesegrids = filter(habgrids, 
  date >= thisclosure$start & date <= thisclosure$end)

theseflow = filter(inflows, as.Date(datetime) >=thisclosure$start & 
  as.Date(datetime) <= thisclosure$end, gauge == "russian river")
theseflow["date"] = as.Date(theseflow$datetime)
thesetransects = filter(ctdmeta, as.Date(start) >= thisclosure$start & as.Date(start) <= thisclosure$end, numcasts == 12)["start"]
thesetransects["date"] = as.Date(thesetransects$start)
theseflow = left_join(theseflow, thesetransects, "date")

ggplot(theseflow, aes(x = datetime, y = flow)) + 
  geom_line(size = 1, color = "darkblue") + scale_x_datetime("") +
  geom_vline(aes(xintercept = as.numeric(start)), linetype = "dashed", size = 1) + 
  ylab("inflow (m3/s)") 
ggsave(paste0("C:/Users/Michael/Desktop/john_plots/", "closureflows", 
  ".png"), width = 10, height = 8)

closureplot = TRUE
#load(purl("rremat/vignettes/analyze-scenario.Rmd", tempfile())  

analysis of Open conditions

#opengeom = geom_bar()
thisopen = filter(ctdmeta, code == "O", isFall, mean.flow < 101, 
  numcasts == 12)
opendates = as.Date(thisopen$start)

thesegrids = filter(habgrids, date %in% opendates)
thesetides = filter(tides, as.Date(datetime) %in% opendates)  
thesetides["date"] = as.Date(thesetides$datetime)
thesetransects = thisopen["start"]
thesetransects["date"] = as.Date(thesetransects$start)
thesetides = left_join(thesetides, thesetransects, "date")
ggplot(thesetides, aes(x = datetime, y = height, group = date)) + 
  geom_line(size = 1, color = "darkblue") + scale_x_datetime("") +
  geom_vline(aes(xintercept = as.numeric(start)), linetype = "dashed", size = 1) + 
  facet_wrap(~date, scales = "free_x", ncol = 1) + ylab("Tide height (MLLW, m)")
ggsave(paste0("C:/Users/Michael/Desktop/john_plots/", "opentides", 
  ".png"), width = 10, height = 8)

closureplot = TRUE
#load(purl("rremat/vignettes/analyze-scenario.Rmd", tempfile())  

mkoohafkan/rremat documentation built on July 3, 2021, 12:06 p.m.