###############################################################################
# Pathname to location of main source codes:
source("P:/Research_Tools/EMR-MRC/codes-EMR/EMR_calc_v39.1.r")
source("P:/current_development-EMR-MRC/codes-EMR-calc-plot/EMR_plot_v36.3b.r")
# Folders for input data and output:
inputlocation = "P:/current_hydrograph-sites/SH_CZO/inputs/"
outputlocation = "P:/current_hydrograph-sites/SH_CZO/outputs/"
###############################################################################
# File names for input data and for output (use "" for no output file):
inputfile = "SH-Q.2009.7mo.cd.csv" # units are d, m^3/d, mm
outputfile = "SH-Q.2009.7mo.out.80QC.csv"
# Type of hydrograph:
# "WT" for water table, "QS" for streamflow with episode end set by the given
# stabilization time, "QC" for streamflow with episode end set at the
# intersection of the constant-slope hydrograph separation with the hydrograph,
# "QM" for the start set by a moving average criterion and end as for "QC"
hytype = "QC"
# Site and input information--input units, watershed area (for QS), labeling:
inpunits = c("d", "m^3/d", "mm")
wsarea = .079 # km^2
probeloc = "Shale Hills"
timeper = "1/20/09-8/17/09"
# MRC coefficient matrix: Each row represents a different part (portion of
# range) of the MRC, with the intercept in the first column and the slope in
# the second.
mrc_coef = matrix(nrow=2, ncol=2)
mrc_coef[1,] = c(0, -0.207971474833795) # 2018-07-17 09:37:46
mrc_coef[2,] = c(45.6112458589481, -0.810849050006656) # 2018-07-17 09:39:58
# Specific yield (when data are with water-table levels) or other appropriate
# measure of response sensitivity (unimportant for streamflow):
capacity = 1
# Lag time between input and response:
lag_time = .7
# Fluctuation tolerance (the maximum dH/dt value allowable as system noise
# rather than response to incoming water flux):
fluc_tol = 70
# Minimum amount of precip within precip window to allow inclusion as an
# episode: Use a negative number to allow no exclusions on this basis.
minprecip = .5
# Parameters for moving minimum criterion: For use with QM type only. For
# WT, QS, and QC, these parameters will be ignored and may be omitted.
# The first parameter is the width of the window (units of time) preceding
# an episode that will be examined for its minimum H. The second is the
# threshold (units of H) by which H must exceed the moving minimum to indicate
# that an episode has begun.
epst_par = c(.2, .5)
# Parameter used in determining the end of an episode:
# For QS, time allowed for stabilization of return to MRC behavior
# For WT, as contingency, time allowed for stabilization
# For QC or QM, slope [L^3 T^-2] of hydrograph separation line
epend_par = 1.5
# For smoothing of hydrograph-slope values, number of data points to left and
# to right of the selected point that will be averaged (weighted in proportion
# to proximity to that point): Use 0 for no smoothing.
Nsmooth = 8
# Axis limits (set = for autogenerated limits):
limt = c(20, 230)
limresp = c(0, 500)
limprecip = c(0, 0)
limsl = c(-100, 200)
# End of frequently-changed parameters ++++++++++++++++++++++++++++++++++++++++
# Units for output, which can be different from input units:
outunits = vector("character")
outunits["precip"] = "mm"
outunits["prectime"] = "h"
outunits["vol"] = "m^3"
# Conversion factors, input to output units:
unicon = vector("numeric")
unicon["precip"] = 1
unicon["prectime"] = 1/24
unicon["epis_product"] = 1 # Product is typically recharge or stormflow
unicon["volarea"] = (1/wsarea) * 1.e-3
###############################################################################
# Read data and convert labeling:
rawdata = read.csv(paste(inputlocation, inputfile, sep=""), header=TRUE)
axislabels = read.csv(paste(inputlocation, inputfile, sep=""),
nrows = 1, header=FALSE)
# Calculate:
outfile_csv = paste0(outputlocation, outputfile)
if (outputfile == "") outfile_csv = NULL
toplabel = paste(probeloc, timeper, hytype, fluc_tol, lag_time, epend_par,
sep=" ")
EMR_results = applyEMR(inputdata=rawdata, hytype=hytype, fluc_tol=fluc_tol,
lag_time=lag_time, epend_par=epend_par, epst_par=epst_par,
capacity=capacity, mrc_coef=mrc_coef, Nsmooth=Nsmooth,
precip_min=minprecip, outfile_csv=outfile_csv, toplabel=toplabel,
inputfile=inputfile, unicon=unicon, inpunits=inpunits, outunits=outunits)
# Plot results:
plot_emr(EMR_results, plot_type="mrc", limt=limt, limresp=limresp,
limprecip=limprecip, limsl=limsl)
plot_emr(EMR_results, plot_type="hydro", limt=limt, limresp=limresp,
limprecip=limprecip, limsl=limsl)
###############################################################################