| Title: | Fourier Analysis of Discharge Data |
|---|---|
| Description: | Computes discrete fast Fourier transform of river discharge data and the derived metrics. The methods are described in J. L. Sabo, D. M. Post (2008) <doi:10.1890/06-1340.1> and J. L. Sabo, A. Ruhi, G. W. Holtgrieve, V. Elliott, M. E. Arias, P. B. Ngor, T. A. Räsänsen, S. Nam (2017) <doi:10.1126/science.aao1053>. |
| Authors: | Samarth Shah [aut, cre], Albert Ruhi [aut], Future H2O [cph] |
| Maintainer: | Samarth Shah <[email protected]> |
| License: | GPL-3 |
| Version: | 1.0.0 |
| Built: | 2024-10-27 04:34:32 UTC |
| Source: | https://github.com/cran/discharge |
Calculates all parameters at once for a single data file. The output gives numerical results from the functions annualnoise, fourierAnalysis, lp3Events,sigmaHighFlows, and sigmaLowFlows.
allstats(file.name, river.name, file.type="txt", date.col=3, discharge.col=4, skipped.rows=28)allstats(file.name, river.name, file.type="txt", date.col=3, discharge.col=4, skipped.rows=28)
file.name |
Character string of the form "file.txt" or "file.csv". |
river.name |
Character string specifying river name. |
file.type |
Character string, "txt" or "csv". Defaults to "txt". |
date.col |
Numeric specifying column containing date in "MM-DD-YYYY" format. Defaults to 3. |
discharge.col |
Numeric specifying column containing discharge data. Defaults to 4. |
skipped.rows |
Numeric indicating number of rows to skip at beginning of file. |
A data frame with columns
a.rms |
Root mean squared amplitude. |
n.rms |
Root mean squared noise. |
snr |
Signal-to-noise ratio. |
theta.d |
Daily noise color. |
theta.a |
Annual noise color. |
sigma.lf |
Sigma for low flow events. |
sigma.hf |
Sigma for high flow events. |
q2 |
2-year return level (flood). |
q10 |
10-year return level (flood). |
l2 |
2-year return level (drought). |
l10 |
10-year return level (drought). |
The arguments "date.col", "discharge.col", and "skipped.rows" are designed to give some flexibility in file input; however, tab-delimited text without extra columns will likely work best.
# allstats function works on files # read R data into temporary file handle data(sycamore) f = tempfile(fileext="txt") write.table(sycamore, file=f, sep="\t") # print all statistics for this river allstats(f,river.name="sycamore", date.col = 2, discharge.col = 3, skipped.rows = 1)# allstats function works on files # read R data into temporary file handle data(sycamore) f = tempfile(fileext="txt") write.table(sycamore, file=f, sep="\t") # print all statistics for this river allstats(f,river.name="sycamore", date.col = 2, discharge.col = 3, skipped.rows = 1)
Calculates annual extreme events for all years in the series. By default, the function finds annual extreme discharge in streamflow object, but any matrix or data.frame may be used.
annualExtremes(x,data.col=NULL, year.col=NULL, moving.avg=FALSE)annualExtremes(x,data.col=NULL, year.col=NULL, moving.avg=FALSE)
x |
Object from which to extract extremes. Should be of class |
data.col |
Optional. If input is a matrix or |
year.col |
Optional. If input is a matrix or |
moving.avg |
Logical; defaults to FALSE. Can be specified TRUE to use 7-day moving average discharge when input is of class " |
A list with items
annual.max |
Matrix giving maximum flow for each year in series. Each row contains the maximum values and all corresponding variables from that observation. |
annual.min |
Matrix giving minimum flow for each year in series. Each row contains the minimum values and all corresponding variables from that observation. |
data(sycamore) sycamore.flows<-asStreamflow(sycamore,river.name="Sycamore Creek") syc.extremes<-annualExtremes(sycamore.flows) names(syc.extremes) syc.extremes$annual.max[1:3,]data(sycamore) sycamore.flows<-asStreamflow(sycamore,river.name="Sycamore Creek") syc.extremes<-annualExtremes(sycamore.flows) names(syc.extremes) syc.extremes$annual.max[1:3,]
The autocovariance function is estimated for the annual maxima in the series. An autoregressive model of the order of the highest significant lag is fit, using the Yule-Walker method to estimate the parameters. The function is transformed into the frequency domain, yielding an estimate theta.a of the annual noise color.
annualnoise(x)annualnoise(x)
x |
A numeric vector of annual extremes. A |
To determine the order of the AR model, the ACF is calculated at all lags less than or equal to the highest power of 2 less than the length of the series. The order of the AR model is the lag with the highest significantly non-zero autocorrelation.
An object of S3 class annualnoise with the following attributes:
auto.corr |
Sample autocorrelation. |
lm.fit |
|
interval |
Upper and lower bounds of a 95% acceptance region when |
log.log |
Matrix with log frequency and log power spectrum. |
reg.stats |
Slope and intercept of regression of log power spectrum on log frequency, where slope is the annual noise color ( |
order |
Indicates order of fitted AR model. |
fit.ar |
Object of class |
data(sycamore) sycamore.flows<-asStreamflow(sycamore,river.name="Sycamore Creek") syc.ar<-annualnoise(sycamore.flows)data(sycamore) sycamore.flows<-asStreamflow(sycamore,river.name="Sycamore Creek") syc.ar<-annualnoise(sycamore.flows)
Ensure that the given arguments are equal length vectors
assert.equal.length(arg0, ...)assert.equal.length(arg0, ...)
arg0 |
At least one input argument to be checked needs to be given |
... |
Other inputs whose length equality needs to be checked |
Stops execution in case of failed check
for.year argumentEnsure that the for.year argument is a numeric scalar value
assert.for.year(for.year)assert.for.year(for.year)
for.year |
input argument to be checked |
Stops execution in case of failed check
Ensure that the given argument is a numeric vector
assert.numeric.vector(arg)assert.numeric.vector(arg)
arg |
input argument to be checked |
Stops execution in case of failed check
This function converts a dataset to S3 object of class "streamflow". Data in this format can be provided as arguments to other functions to call default procedures.
asStreamflow(x,river.name=NULL, start.date=NULL, end.date=NULL,max.na=10) ## S3 method for class 'streamflow' print(x, ...) ## S3 method for class 'streamflow' summary(object, ...)asStreamflow(x,river.name=NULL, start.date=NULL, end.date=NULL,max.na=10) ## S3 method for class 'streamflow' print(x, ...) ## S3 method for class 'streamflow' summary(object, ...)
x, object
|
A matrix with dates in the first column and discharge values in the second column. Dates should be of the format "YYYY-MM-DD". |
river.name |
A character vector listing the river name. |
start.date |
Optional character string giving date to start analysis, of the format "YYYY-MM-DD" |
end.date |
Optional date to start analysis, of the format "YYYY-MM-DD" |
max.na |
Optional number specifying maximum NA values to allow. |
... |
Other parameters. |
An object of class streamflow containing the following items:
data |
Data frame with 8 columns |
n |
Number of observations in series. |
n.nas |
Number of NA observations in series. |
start |
Start date. |
end |
End date. |
name |
Name of river. |
na.info |
Matrix containing the index and start date of all blocks of more than one NA observation. |
data(sycamore) sycamore.flows<-asStreamflow(sycamore,start.date="1965-01-01", end.date="2010-03-15",river.name="Sycamore Creek")data(sycamore) sycamore.flows<-asStreamflow(sycamore,start.date="1965-01-01", end.date="2010-03-15",river.name="Sycamore Creek")
Uses the package "CircStats" to find the mean direction, rho, and kappa of the von Mises distribution summarizing the ordinal day of high- and low-flow events.
circ.s(x)circ.s(x)
x |
Output from "compare.periods" function (of class " |
circstats |
|
# load data data("sycamore") # compare for periods from 1960 to 1979 and 1980 to 1999 comp = compare.periods(c("1960-01-01", "1979-12-31"), c("1980-01-01", "1999-12-31"), sycamore, plot=FALSE) circ.s(comp)# load data data("sycamore") # compare for periods from 1960 to 1979 and 1980 to 1999 comp = compare.periods(c("1960-01-01", "1979-12-31"), c("1980-01-01", "1999-12-31"), sycamore, plot=FALSE) circ.s(comp)
The function finds the seasonal signal in the first of two time periods within the data series. This signal is used to calculate residual flows in the first time period and "historical residuals", the residuals when the second period is compared to the signal in the first period.
The output of the function gives event information for all residual events greater than as calculated from the functions sigmaHighFlows and sigmaLowFlows.
compare.periods(p1, p2, x, plot=T) ## S3 method for class 'compflows' plot(x, ...)compare.periods(p1, p2, x, plot=T) ## S3 method for class 'compflows' plot(x, ...)
p1 |
Character vector specifying start and end date of the first period. "YYYY-MM-DD" format. |
p2 |
Character vector specifying start and end date of the second period. "YYYY-MM-DD" format. |
x |
Matrix with first column specifying dates and second column specifying raw discharge data. |
plot |
Logical; defaults to TRUE. If TRUE, the seasonal signal for both periods will be plotted. |
... |
Other parameters. |
Object of S3 ckass compflows with the following items:
sigma.low |
Estimate of |
sigma.high |
Estimate of |
p1.levents |
Matrix of low flow events for period 1 |
p1.hevents |
Matrix of high flow events for period 1 |
p2.levents |
Matrix of low flow events for period 2 |
p2.hevents |
Matrix of high flow events for period 2 |
# load data data("sycamore") # compare for periods from 1960 to 1979 and 1980 to 1999 compare.periods(c("1960-01-01", "1979-12-31"), c("1980-01-01", "1999-12-31"), sycamore)# load data data("sycamore") # compare for periods from 1960 to 1979 and 1980 to 1999 compare.periods(c("1960-01-01", "1979-12-31"), c("1980-01-01", "1999-12-31"), sycamore)
A data series from USGS station 02392000, Etowah River (Canton,GA).
data(etowah)data(etowah)
The data are supplied as a matrix with years in the first column and average daily discharge in the second column. The series contains 18282 observations starting Jan 1, 1960 and ending Dec 31, 2009.
USGS National Water Information Sysem: http://waterdata.usgs.gov/usa/nwis/uv?site_no=02392000
This is a wrapper function to calculate all the DFFT metrics for the given input signal
fft_metrics(data, candmin, candmax, river.name = "", baseline.signal = NULL)fft_metrics(data, candmin, candmax, river.name = "", baseline.signal = NULL)
data |
A matrix with dates in the first column and discharge values in the second column. Dates should be of the format "YYYY-MM-DD" |
candmin |
numeric vector of possible ordinal days in which the predicted signal is lowest. This range need not be narrow, but a string of consecutive days should not include more than only local minimum. Used for calculating the high- and low-flow windows |
candmax |
numeric vector of possible ordinal days in which the predicted signal is highest. This range need not be narrow, but a string of consecutive days should not include more than only local maximum. |
river.name |
A character vector listing the river name. |
baseline.signal |
If |
A list containing 2 data frames:
high.level.metrics |
Data frame containing NAA and FPExt values for each year in the given series |
naa.shape.components |
Data frame containing HSAM, LSAM, Transition time, HSAF, LSAF, timing of HSAM, timing of LSAM, IFI, IDI |
# load sample data data("sycamore") x = sycamore # get streamflow object for the sample data x.streamflow = asStreamflow(x) # prepare baseline signal x.bl = prepareBaseline(x.streamflow) # fetch the DFFT metrics for this sample data # "candmax" chosen because preliminary analysis (e.g. with fourierAnalysis # output) shows the signal is highest sometime between # day 190 and day 330 # "candmin" can be estimated analogously. x.fftmetrics = fft_metrics(x, river.name = "Sycamore", candmin = c(40:125), candmax = c(190:330), baseline.signal = x.bl)# load sample data data("sycamore") x = sycamore # get streamflow object for the sample data x.streamflow = asStreamflow(x) # prepare baseline signal x.bl = prepareBaseline(x.streamflow) # fetch the DFFT metrics for this sample data # "candmax" chosen because preliminary analysis (e.g. with fourierAnalysis # output) shows the signal is highest sometime between # day 190 and day 330 # "candmin" can be estimated analogously. x.fftmetrics = fft_metrics(x, river.name = "Sycamore", candmin = c(40:125), candmax = c(190:330), baseline.signal = x.bl)
This function takes as input a streamFlow object and a year, and outputs the timing and magnitude of noteworthy spectral anomolies.
fftmetrics(x,year,candmin,candmax)fftmetrics(x,year,candmin,candmax)
x |
streamflow object, as output from the asStreamflow function. |
year |
integer; data from this year will be analyzed. The plotted hydrograph will include data from August of the previous year to April of the following year. |
candmin |
numeric vector of possible ordinal days in which the predicted signal is lowest. This range need not be narrow, but a string of consecutive days should not include more than only local minimum. Used for calculating the high- and low-flow windows. |
candmax |
numeric vector of possible ordinal days in which the predicted signal is highest. This range need not be narrow, but a string of consecutive days should not include more than only local maximum. |
A list object with the following components;
sam |
dataframe where row 1 corresponds to the largest low-residual event in the year, and row 2 corresponds to the largest high-residual. The dataframe contains each event's date, ordinal day, magnitude of residual and signal, index in the original data, and timing (in days) relative to the reference point. |
ref.point |
date used as a reference point for the timing of max and min events. If not given by the user, it is the first local maximum of the signal, within the year specified. |
events |
dataframe with rows corresponding to start dates of the longest low-flow and high-flow events. The data frame contains signal and residual data for the start date and the lenght of the run. |
auc |
dataframe with three values: 'net.auc', 'rel.auc.low', and rel.auc.high'. 'net.auc' is the sum of positive residuals in the high flow window divided by the sum of negative residuals in the low flow window. 'rel.auc.low' is postive residuals divided by the sum of negative residuals |
noise.color |
numeric, "theta.d" as calculated as in the fourierAnalysis function. |
data(etowah) etowah.flows=asStreamflow(etowah, river.name="Etowah") # "candmax" chosen because preliminary analysis (e.g., with fourierAnalysis output) # shows the signal is highest sometime between day 40 and day 125. # "candmin" can be estimated analogously. fftmetrics(etowah.flows,2002,candmin=c(190:330),candmax=c(40:125))data(etowah) etowah.flows=asStreamflow(etowah, river.name="Etowah") # "candmax" chosen because preliminary analysis (e.g., with fourierAnalysis output) # shows the signal is highest sometime between day 40 and day 125. # "candmin" can be estimated analogously. fftmetrics(etowah.flows,2002,candmin=c(190:330),candmax=c(40:125))
Filter the baseline signal for a given time window
filterBaseline(bl, filter.date.start, filter.date.end, date.format = "%Y-%m-%d")filterBaseline(bl, filter.date.start, filter.date.end, date.format = "%Y-%m-%d")
bl |
baseline signal as returned from the function |
filter.date.start |
start date of the filtering window |
filter.date.end |
end date of the filtering window |
date.format |
format of date specified in |
baseline signal filtered for the given date window
# load sample data data("sycamore") x = sycamore # get streamflow object for the sample data x.streamflow = asStreamflow(x) # baseline for single run for all the years in input signal x.bl = prepareBaseline(x.streamflow) # filter the baseline signal between years 1993 and 2000 x.bl.filtered = filterBaseline(x.bl, filter.date.start = "1993-01-01", filter.date.end = "2000-12-31")# load sample data data("sycamore") x = sycamore # get streamflow object for the sample data x.streamflow = asStreamflow(x) # baseline for single run for all the years in input signal x.bl = prepareBaseline(x.streamflow) # filter the baseline signal between years 1993 and 2000 x.bl.filtered = filterBaseline(x.bl, filter.date.start = "1993-01-01", filter.date.end = "2000-12-31")
Find median value from the filtered part of the input vector. This function is used as a helper to the bootstrapping routine in prepare baseline
findMed(data, indices)findMed(data, indices)
data |
complete input vector |
indices |
logical vector to filter indices from |
median value of the filtered vector
The fast Fourier transform is used to extract the seasonal signal of a time series. The significant frequencies are found from among periods of length 2-, 3-, 4-, 6-, 12-, and 18-months.
The signal may be specified as stationary or non-stationary. If a non-stationary fit is allowed, simple linear regression estimates the long term linear trend. The seasonal signal is calculcated from the residuals.
Predicted flow (and corresponding residual) at each time point is calculated from seasonal signal and, if non-stationary, long term trend coefficient.
fourierAnalysis(x, stationary=F) ## S3 method for class 'ssignal' plot(x, plot.type="hydrograph", ...)fourierAnalysis(x, stationary=F) ## S3 method for class 'ssignal' plot(x, plot.type="hydrograph", ...)
x |
An object of class |
stationary |
Logical; defaults to FALSE. |
plot.type |
Indicates the type of plot to create. The default "hydrograph" produces a plot of ordinary day and log normalized discharge, with the seasonal signal overlaid. "auto.corr" produces a plot of daily autocorrelation as calculated from the residual flows. |
... |
Other parameters. |
An object of class ssignal with items
signal |
Data matrix augmented to included predicted and residual values. |
terms |
Matrix containing amplitude, phase, and frequency of seasonal signal. |
detrend.fit |
An |
logps.regression |
An |
rms |
|
data(sycamore) sycamore.flows<-asStreamflow(sycamore,river.name="Sycamore Creek") syc.seas<-fourierAnalysis(sycamore.flows) summary(syc.seas)data(sycamore) sycamore.flows<-asStreamflow(sycamore,river.name="Sycamore Creek") syc.seas<-fourierAnalysis(sycamore.flows) summary(syc.seas)
Calculate the Flood Pulse Extent (FPExt) from the given residual values.
getFPExt(resid, years, for.year = NULL)getFPExt(resid, years, for.year = NULL)
resid |
A vector of residual values generated with respect to the baseline signal |
years |
A vector of years corrosponding to the residual values |
for.year |
(optional) Calculate FPExt values only for the given year in this argument. If argument is omitted, NAA values for all years are calculated. |
Data frame containing two columns:
year |
First column, represents year |
FPExt |
Second column, represents FPExt values |
# load sample data data("sycamore") x = sycamore # get streamflow object for the sample data x.streamflow = asStreamflow(x) # prepare baseline signal x.bl = prepareBaseline(x.streamflow) # FPExt fpext = getFPExt(x.bl$resid.sig, x.streamflow$data$year)# load sample data data("sycamore") x = sycamore # get streamflow object for the sample data x.streamflow = asStreamflow(x) # prepare baseline signal x.bl = prepareBaseline(x.streamflow) # FPExt fpext = getFPExt(x.bl$resid.sig, x.streamflow$data$year)
Compute High Spectral Anomaly Frequency (HSAF) from the given residual values.
getHSAF(resid, years, for.year = NULL)getHSAF(resid, years, for.year = NULL)
resid |
A vector of residual values generated with respect to the baseline signal |
years |
A vector of years corrosponding to the residual values |
for.year |
(optional) Calculate HSAF values only for the given year in this argument. If argument is omitted, HSAF values for all years are calculated. |
Data frame containing two Columns:
year |
First column, represents year |
HSAF |
Second column, represents HSAF values |
# load sample data data("sycamore") x = sycamore # get streamflow object for the sample data x.streamflow = asStreamflow(x) # prepare baseline signal x.bl = prepareBaseline(x.streamflow) # HSAF hsaf = getHSAF(x.bl$resid.sig, x.streamflow$data$year)# load sample data data("sycamore") x = sycamore # get streamflow object for the sample data x.streamflow = asStreamflow(x) # prepare baseline signal x.bl = prepareBaseline(x.streamflow) # HSAF hsaf = getHSAF(x.bl$resid.sig, x.streamflow$data$year)
Compute High Spectral Anomaly Magnitude (HSAM) from the given residual values each year
getHSAM(resid, years, for.year = NULL)getHSAM(resid, years, for.year = NULL)
resid |
A vector of residual values generated with respect to the baseline signal |
years |
A vector of years corrosponding to the residual values |
for.year |
(optional) Calculate HSAM values only for the given year in this argument. If argument is omitted, HSAM values for all years are calculated. |
Data frame containing four columns:
year |
First column, represents year |
HSAM |
Second column, represents HSAM values |
index.year |
Third column, representing index of HSAM value in that year |
index.all |
Fourth column, representing index of HSAM value in the input resid |
# load sample data data("sycamore") x = sycamore # get streamflow object for the sample data x.streamflow = asStreamflow(x) # prepare baseline signal x.bl = prepareBaseline(x.streamflow) # HSAM hsam = getHSAM(x.bl$resid.sig, x.streamflow$data$year)# load sample data data("sycamore") x = sycamore # get streamflow object for the sample data x.streamflow = asStreamflow(x) # prepare baseline signal x.bl = prepareBaseline(x.streamflow) # HSAM hsam = getHSAM(x.bl$resid.sig, x.streamflow$data$year)
Compute Inter-Draught Interval (IDI) from the given residual values.
getIDI(resid, years, highflow.start, highflow.end, unique.years, for.year = NULL)getIDI(resid, years, highflow.start, highflow.end, unique.years, for.year = NULL)
resid |
A vector of residual values generated with respect to the baseline signal |
years |
A vector of years corrosponding to the residual values |
highflow.start |
A vector giving start index of high-flow window in each year |
highflow.end |
A vector giving end index of high-flow window in each year |
unique.years |
A vector or year values corresponding to the |
for.year |
(optional) Calculate IDI values only for the given year in this argument. If argument is omitted, IDI values for all years are calculated. |
Data frame containing two columns:
year |
First column, represents year |
IDI |
Second column, represents IDI values |
# load sample data data("sycamore") x = sycamore # get streamflow object for the sample data x.streamflow = asStreamflow(x) # prepare baseline signal x.bl = prepareBaseline(x.streamflow) # signal parts x.sp = getSignalParts(x.bl$pred2, candmin = c(40:125), candmax = c(190:330), years = x.streamflow$data$year, months = x.streamflow$data$month, jdays = x.streamflow$data$jday) # IDI idi = getIDI(x.bl$resid.sig, x.streamflow$data$year, x.sp$HF.window.start, x.sp$HF.window.end, x.sp$year)# load sample data data("sycamore") x = sycamore # get streamflow object for the sample data x.streamflow = asStreamflow(x) # prepare baseline signal x.bl = prepareBaseline(x.streamflow) # signal parts x.sp = getSignalParts(x.bl$pred2, candmin = c(40:125), candmax = c(190:330), years = x.streamflow$data$year, months = x.streamflow$data$month, jdays = x.streamflow$data$jday) # IDI idi = getIDI(x.bl$resid.sig, x.streamflow$data$year, x.sp$HF.window.start, x.sp$HF.window.end, x.sp$year)
Compute Inter-Flood Interval (IFI) from the given residual values.
getIFI(resid, years, lowflow.start, lowflow.end, unique.years, for.year = NULL)getIFI(resid, years, lowflow.start, lowflow.end, unique.years, for.year = NULL)
resid |
A vector of residual values generated with respect to the baseline signal |
years |
A vector of years corrosponding to the residual values |
lowflow.start |
A vector giving start index of low-flow window in each year |
lowflow.end |
A vector giving end index of low-flow window in each year |
unique.years |
A vector or year values corresponding to the |
for.year |
(optional) Calculate IFI values only for the given year in this argument. If argument is omitted, IFI values for all years are calculated. |
Data frame containing two columns:
year |
First column, represents year |
IFI |
Second column, represents IFI values |
# load sample data data("sycamore") x = sycamore # get streamflow object for the sample data x.streamflow = asStreamflow(x) # prepare baseline signal x.bl = prepareBaseline(x.streamflow) # signal parts x.sp = getSignalParts(x.bl$pred2, candmin = c(40:125), candmax = c(190:330), years = x.streamflow$data$year, months = x.streamflow$data$month, jdays = x.streamflow$data$jday) # IFI ifi = getIFI(x.bl$resid.sig, x.streamflow$data$year, x.sp$LF.window.start, x.sp$LF.window.end, x.sp$year)# load sample data data("sycamore") x = sycamore # get streamflow object for the sample data x.streamflow = asStreamflow(x) # prepare baseline signal x.bl = prepareBaseline(x.streamflow) # signal parts x.sp = getSignalParts(x.bl$pred2, candmin = c(40:125), candmax = c(190:330), years = x.streamflow$data$year, months = x.streamflow$data$month, jdays = x.streamflow$data$jday) # IFI ifi = getIFI(x.bl$resid.sig, x.streamflow$data$year, x.sp$LF.window.start, x.sp$LF.window.end, x.sp$year)
Compute Low Spectral Anomaly Frequency (LSAF) from the given residual values.
getLSAF(resid, years, for.year = NULL)getLSAF(resid, years, for.year = NULL)
resid |
A vector of residual values generated with respect to the baseline signal |
years |
A vector of years corrosponding to the residual values |
for.year |
(optional) Calculate LSAF values only for the given year in this argument. If argument is omitted, LSAF values for all years are calculated. |
Data frame containing two Columns:
year |
First column, represents year |
LSAF |
Second column, represents LSAF values |
# load sample data data("sycamore") x = sycamore # get streamflow object for the sample data x.streamflow = asStreamflow(x) # prepare baseline signal x.bl = prepareBaseline(x.streamflow) # LSAF lsaf = getLSAF(x.bl$resid.sig, x.streamflow$data$year)# load sample data data("sycamore") x = sycamore # get streamflow object for the sample data x.streamflow = asStreamflow(x) # prepare baseline signal x.bl = prepareBaseline(x.streamflow) # LSAF lsaf = getLSAF(x.bl$resid.sig, x.streamflow$data$year)
Compute Low Spectral Anomaly Magnitude (LSAM) from the given residual values each year
getLSAM(resid, years, for.year = NULL)getLSAM(resid, years, for.year = NULL)
resid |
A vector of residual values generated with respect to the baseline signal |
years |
A vector of years corrosponding to the residual values |
for.year |
(optional) Calculate LSAM values only for the given year in this argument. If argument is omitted, LSAM values for all years are calculated. |
Data frame containing four columns:
year |
First column, represents year |
LSAM |
Second column, represents LSAM values |
index.year |
Third column, representing index of LSAM value in that year |
index.all |
Fourth column, representing index of LSAM value in the input resid |
# load sample data data("sycamore") x = sycamore # get streamflow object for the sample data x.streamflow = asStreamflow(x) # prepare baseline signal x.bl = prepareBaseline(x.streamflow) # LSAM lsam = getLSAM(x.bl$resid.sig, x.streamflow$data$year)# load sample data data("sycamore") x = sycamore # get streamflow object for the sample data x.streamflow = asStreamflow(x) # prepare baseline signal x.bl = prepareBaseline(x.streamflow) # LSAM lsam = getLSAM(x.bl$resid.sig, x.streamflow$data$year)
Calculate Net Annual Anomaly (NAA) from the given residual values.
getNAA(resid, years, for.year = NULL)getNAA(resid, years, for.year = NULL)
resid |
A vector of residual values generated with respect to the baseline signal |
years |
A vector of years corrosponding to the residual values |
for.year |
(optional) Calculate NAA values only for the given year in this argument. If argument is omitted, NAA values for all years are calculated. |
Data frame containing two columns:
year |
First column, represents year |
NAA |
Second column, represents NAA values |
# load sample data data("sycamore") x = sycamore # get streamflow object for the sample data x.streamflow = asStreamflow(x) # prepare baseline signal x.bl = prepareBaseline(x.streamflow) # NAA naa = getNAA(x.bl$resid.sig, x.streamflow$data$year)# load sample data data("sycamore") x = sycamore # get streamflow object for the sample data x.streamflow = asStreamflow(x) # prepare baseline signal x.bl = prepareBaseline(x.streamflow) # NAA naa = getNAA(x.bl$resid.sig, x.streamflow$data$year)
This function computes high flow and low flow window of seasonal signal, and the peak max and peak min values.
getSignalParts(seas.sig, candmin, candmax, years, months, jdays, for.year = NULL)getSignalParts(seas.sig, candmin, candmax, years, months, jdays, for.year = NULL)
seas.sig |
Seasonal signal as generated from DFFT methods |
candmin |
numeric vector of possible ordinal days in which the predicted signal is lowest. This range need not be narrow, but a string of consecutive days should not include more than one local minimum. Used for calculating the high- and low-flow windows. |
candmax |
numeric vector of possible ordinal days in which the predicted signal is highest. This range need not be narrow, but a string of consecutive days should not include more than one local maximum. |
years |
A vector of years corrosponding to the seasonal signal values |
months |
A vector of months corrosponding to the seasonal signal values |
jdays |
A vector of julian days corrosponding to the seasonal signal values |
for.year |
(optional) Calculate signal parts only for the given year in this argument. If argument is omitted, all years are considered. |
Data frame containing following columns.
year |
represents year |
max.peak.index.all |
represents index value within the entire vector |
max.peak.value |
represents value of max peak |
highwind.start.index.all |
start index of high flow window within the entire vector |
highwind.end.index.all |
end index of high flow window within the entire vector |
lowwind.start.index.all |
start index of low flow window within the entire vector |
lowwind.end.index.all |
end index of low flow window within the entire vector |
# load sample data data("sycamore") x = sycamore # get streamflow object for the sample data x.streamflow = asStreamflow(x) # prepare baseline signal x.bl = prepareBaseline(x.streamflow) # signal parts x.sp = getSignalParts(x.bl$pred2, candmin = c(40:125), candmax = c(190:330), years = x.streamflow$data$year, months = x.streamflow$data$month, jdays = x.streamflow$data$jday)# load sample data data("sycamore") x = sycamore # get streamflow object for the sample data x.streamflow = asStreamflow(x) # prepare baseline signal x.bl = prepareBaseline(x.streamflow) # signal parts x.sp = getSignalParts(x.bl$pred2, candmin = c(40:125), candmax = c(190:330), years = x.streamflow$data$year, months = x.streamflow$data$month, jdays = x.streamflow$data$jday)
Compute the number of days separating HSAM and reference point for each year.
getTimingHSAM(index.hsam, index.ref, years, for.year = NULL)getTimingHSAM(index.hsam, index.ref, years, for.year = NULL)
index.hsam |
A scalar/vector of index of HSAM values in given year/years |
index.ref |
A scalar/vector of index of reference point in given year/years |
years |
A vector of years corresponding to HSAM and ref values. This argument can be NULL if the HSAM and ref values are scalars. |
for.year |
(optional) Calculate timing (HSAM) only for the given year in this argument. If argument is omitted, timing (HSAM) values for all years are calculated. |
Scalar timing HSAM value if the inputs are scalars, or a Data frame containing two Columns:
year |
First column, represents year |
timing.hsam |
Second column, represents hsam timing values |
# load sample data data("sycamore") x = sycamore # get streamflow object for the sample data x.streamflow = asStreamflow(x) # prepare baseline signal x.bl = prepareBaseline(x.streamflow) # get signal parts x.sp = getSignalParts(x.bl$pred2, candmin = c(40:125), candmax = c(190:330), years = x.streamflow$data$year, months = x.streamflow$data$month, jdays = x.streamflow$data$jday) # get HSAM values hsam = getHSAM(x.bl$resid.sig, x.streamflow$data$year) # timing HSAM thsam = getTimingHSAM(hsam$Index.all, x.sp$peak.index, x.sp$year)# load sample data data("sycamore") x = sycamore # get streamflow object for the sample data x.streamflow = asStreamflow(x) # prepare baseline signal x.bl = prepareBaseline(x.streamflow) # get signal parts x.sp = getSignalParts(x.bl$pred2, candmin = c(40:125), candmax = c(190:330), years = x.streamflow$data$year, months = x.streamflow$data$month, jdays = x.streamflow$data$jday) # get HSAM values hsam = getHSAM(x.bl$resid.sig, x.streamflow$data$year) # timing HSAM thsam = getTimingHSAM(hsam$Index.all, x.sp$peak.index, x.sp$year)
Compute the number of days separating LSAM and reference point for each year.
getTimingLSAM(index.lsam, index.ref, years = NULL, for.year = NULL)getTimingLSAM(index.lsam, index.ref, years = NULL, for.year = NULL)
index.lsam |
A scalar/vector of index of LSAM values in given year/years |
index.ref |
A scalar/vector of index of reference point in given year/years |
years |
(optional) A vector of years corresponding to LSAM and ref values. This argument can be NULL if the LSAM and ref values are scalars. |
for.year |
(optional) Calculate timing (LSAM) only for the given year in this argument. If argument is omitted, timing (LSAM) values for all years are calculated. |
Scalar timing LSAM value if the inputs are scalars, or a Data frame containing two Columns:
year |
First column, represents year |
timing.lsam |
Second column, represents lsam timing values |
# load sample data data("sycamore") x = sycamore # get streamflow object for the sample data x.streamflow = asStreamflow(x) # prepare baseline signal x.bl = prepareBaseline(x.streamflow) # get signal parts x.sp = getSignalParts(x.bl$pred2, candmin = c(40:125), candmax = c(190:330), years = x.streamflow$data$year, months = x.streamflow$data$month, jdays = x.streamflow$data$jday) # get LSAM values lsam = getLSAM(x.bl$resid.sig, x.streamflow$data$year) # timing LSAM tlsam = getTimingLSAM(lsam$Index.all, x.sp$peak.index, x.sp$year)# load sample data data("sycamore") x = sycamore # get streamflow object for the sample data x.streamflow = asStreamflow(x) # prepare baseline signal x.bl = prepareBaseline(x.streamflow) # get signal parts x.sp = getSignalParts(x.bl$pred2, candmin = c(40:125), candmax = c(190:330), years = x.streamflow$data$year, months = x.streamflow$data$month, jdays = x.streamflow$data$jday) # get LSAM values lsam = getLSAM(x.bl$resid.sig, x.streamflow$data$year) # timing LSAM tlsam = getTimingLSAM(lsam$Index.all, x.sp$peak.index, x.sp$year)
Compute the number of days separating HSAM and LSAM for the given year/years.
getTransitionTime(index.hsam, index.lsam, years, for.year = NULL)getTransitionTime(index.hsam, index.lsam, years, for.year = NULL)
index.hsam |
A scalar/vector of index of HSAM values in given year/years |
index.lsam |
A scalar/vector of index of LSAM values in given year/years |
years |
A vector of years corresponding to HSAM and LSAM values. This argument can be NULL if the HSAM and LSAM values are scalars. |
for.year |
(optional) Calculate transition time only for the given year in this argument. If argument is omitted, transition times for all years are calculated. |
Scalar transition time if the inputs are scalars, or a Data frame containing two Columns:
year |
First column, represents year |
transition.time |
Second column, represents transition times |
# load sample data data("sycamore") x = sycamore # get streamflow object for the sample data x.streamflow = asStreamflow(x) # prepare baseline signal x.bl = prepareBaseline(x.streamflow) # get HSAM and LSAM values hsam = getHSAM(x.bl$resid.sig, x.streamflow$data$year) lsam = getLSAM(x.bl$resid.sig, x.streamflow$data$year) # transition time tt = getTransitionTime(hsam$Index.all, lsam$Index.all, hsam$year)# load sample data data("sycamore") x = sycamore # get streamflow object for the sample data x.streamflow = asStreamflow(x) # prepare baseline signal x.bl = prepareBaseline(x.streamflow) # get HSAM and LSAM values hsam = getHSAM(x.bl$resid.sig, x.streamflow$data$year) lsam = getLSAM(x.bl$resid.sig, x.streamflow$data$year) # transition time tt = getTransitionTime(hsam$Index.all, lsam$Index.all, hsam$year)
Finds independent events greater than or less than a specified criterion. High (or low) -flow days occurring on consecutive days are considered part of one event. This function can be used to find events exceeding 2- or 10- year return levels (as calculated in lp3Events function, for example), or to find residual flows of a certain magnitude.
independentEvents(cutoff.val, data, data.column, below.cutoff=FALSE)independentEvents(cutoff.val, data, data.column, below.cutoff=FALSE)
cutoff.val |
Numeric specifying event criterion. |
data |
Data matrix or data frame with one column of streamflow data. |
data.column |
Numeric; specifies column in which to look for events. |
below.cutoff |
Logical. TRUE to find events less than the |
A data.frame with columns
events.starts |
Index of event start. |
events.ends |
Index of event end. |
events.duration |
Length (days) of event. |
extreme.this.events |
Maximum or minimum flow for this event. |
ind.extreme |
Index of maximum or minimum flow for this event. If extreme is not unique, the chronologically first index is given. |
... |
All columns of original data, corresponding to max or min flow. These columns will have the same column names as the original data. |
duplicates |
0 if the extreme is unique, 1 if it is not unique. |
data(sycamore) syc.sf<-asStreamflow(sycamore) #find 10-year flood q10<-lp3Events(syc.sf)$Q10 #find all events greater than 10-year flood independentEvents(q10,syc.sf$data, data.col=8 , below.cutoff=FALSE)data(sycamore) syc.sf<-asStreamflow(sycamore) #find 10-year flood q10<-lp3Events(syc.sf)$Q10 #find all events greater than 10-year flood independentEvents(q10,syc.sf$data, data.col=8 , below.cutoff=FALSE)
Uses the method of moments to find the 2- and 10-year droughts and floods under the Log-III Pearson distribution.
lp3Events(x)lp3Events(x)
x |
Object of class " |
Return levels are calculated using the method of moments through the package "lmom". High return levels (floods) are calculated using annual maxima of the raw (log normalized) data, while the low return levels (droughts) are calculated using annual mimima of the 7-day moving averages.
A list with items
Q2 |
2-year high return level. |
Q10 |
10-year high return level. |
L2 |
2-year low return level. |
L10 |
10-year low return level. |
data(sycamore) syc.sf<-asStreamflow(sycamore) lp3Events(syc.sf)data(sycamore) syc.sf<-asStreamflow(sycamore) lp3Events(syc.sf)
Function takes a vector of file names and returns seasonal signal-to-noise ratio, daily and annual noise color, 2- and 10-year return levels, and for high- and low-flow events.
All files in the vector should be of the same format.
parameters.list(x, names=NULL, file.type="txt", date.col=3, dis.col=4, skipped.rows=28)parameters.list(x, names=NULL, file.type="txt", date.col=3, dis.col=4, skipped.rows=28)
x |
Character vector containing file names. |
names |
Optional character vector with names of sites. |
file.type |
Character string, "txt" or "csv". Defaults to "txt". |
date.col |
Numeric specifying column containing date in "MM-DD-YYYY" format. Defaults to 3. |
dis.col |
Numeric specifying column containing discharge data. Defaults to 4. |
skipped.rows |
Numeric indicating number of rows to skip at beginning of file. |
A data frame with one row for each file and the following columns:
a.rms |
Root mean squared amplitude. |
n.rms |
Root mean squared noise. |
snr |
Signal-to-noise ratio. |
theta.d |
Daily noise color. |
theta.a |
Annual noise color. |
sigma.lf |
Sigma for low flow events. |
sigma.hf |
Sigma for high flow events. |
q2 |
2-year return level (flood). |
q10 |
10-year return level (flood). |
l2 |
2-year return level (drought). |
l10 |
10-year return level (drought). |
The arguments "date.col", "discharge.col", and "skipped.rows" are designed to give some flexibility in file input; however, tab-delimited text without extra columns will work best.
# this function works on list of files # read R data into temporary file handle data(sycamore) f = tempfile(fileext="txt") write.table(sycamore, file=f, sep="\t") # print all statistics for the list of rivers parameters.list(c(f), names=c("sycamore"), date.col=2,dis.col=3,skipped.rows = 1)# this function works on list of files # read R data into temporary file handle data(sycamore) f = tempfile(fileext="txt") write.table(sycamore, file=f, sep="\t") # print all statistics for the list of rivers parameters.list(c(f), names=c("sycamore"), date.col=2,dis.col=3,skipped.rows = 1)
Runs fourier analysis on the input signal, to build baseline signal.
prepareBaseline(x, year.start = NULL, year.end = NULL, window.20 = FALSE)prepareBaseline(x, year.start = NULL, year.end = NULL, window.20 = FALSE)
x |
streamflow object, as output from the |
year.start |
Start of the year for estimating baseline, or |
year.end |
End of the year for estimating baseline, or |
window.20 |
If |
ssignal object containing the baseline signal
# load sample data data("sycamore") x = sycamore # get streamflow object for the sample data x.streamflow = asStreamflow(x) # baseline for single run for all the years in input signal bl.singlerun.all = prepareBaseline(x.streamflow) # baseline for singlerun between the given start and end years bl.singlerun.filtered = prepareBaseline(x.streamflow, year.start = 1961, year.end = 2000) # baseline with windowinng and bootstrapping on all years in the input signal bl.windowed.all = prepareBaseline(x.streamflow, window.20 = TRUE) # baseline with windowing and bootstrapping on given start year # with end year inferred from singal bl.windowed.filtered = prepareBaseline(x.streamflow, year.start = 1961, window.20 = TRUE)# load sample data data("sycamore") x = sycamore # get streamflow object for the sample data x.streamflow = asStreamflow(x) # baseline for single run for all the years in input signal bl.singlerun.all = prepareBaseline(x.streamflow) # baseline for singlerun between the given start and end years bl.singlerun.filtered = prepareBaseline(x.streamflow, year.start = 1961, year.end = 2000) # baseline with windowinng and bootstrapping on all years in the input signal bl.windowed.all = prepareBaseline(x.streamflow, window.20 = TRUE) # baseline with windowing and bootstrapping on given start year # with end year inferred from singal bl.windowed.filtered = prepareBaseline(x.streamflow, year.start = 1961, window.20 = TRUE)
Creates a plot with the maximum annual low- and high residuals for each year in the series.
residplot.extreme(x, text=FALSE, data=FALSE)residplot.extreme(x, text=FALSE, data=FALSE)
x |
Object of class |
text |
Logical. If true, points corresponding to flows greater than |
data |
Logical. If true, the extreme residuals are returned in the output. |
Plot with year on the x-axis and the maximum residual magnitude for that year on the y-axis.
If data=TRUE, output includes a list with the following components:
annual.max |
Matrix with data corresponding to the maximum residual flow for each year in series. |
annual.min |
Matrix with data corresponding to the minimum residual flow for each year in series. |
# load data data(sycamore) # plot residplot.extreme(asStreamflow(sycamore))# load data data(sycamore) # plot residplot.extreme(asStreamflow(sycamore))
Creates a . This is similar to the plots produced in the "compare.periods" function, but only displays the data for a single time period.
sigmaeventsplot(x)sigmaeventsplot(x)
x |
Object of class "streamflow". |
A "ggplot2" plot depicting frequency of events greater than sigma, organized circularly by ordinal day of the year.
sigmaHighFlows
sigmaLowFlows
compare.periods
Calculates catastrophic variability for high flow events. Positive residuals from the seasonal signal are used to calculate , the standard deviation of high-flow events.
sigmaHighFlows(x, resid.column)sigmaHighFlows(x, resid.column)
x |
An object of class |
resid.column |
Optional numeric specifiying which column contains residuals. Required if |
An object of class list with items
n.floods |
Number of independent events with positive residuals. |
sigma.hfa |
Estimated sigma using the y-intercept. |
sigma.hfb |
Estimated sigma using the slope ( |
flood.line |
Matrix containing fitted, observed, and residual values from regression of log counts on bin midpoints. |
onesigma.events |
matrix containing information for all events below |
twosigma.events |
matrix containing information for all events below |
independentEvents
sigmaLowFlows
# load data data(sycamore) # get streamflow object sf = asStreamflow(sycamore) # estimate catastrophic high flow variability sigmaHighFlows(sf)# load data data(sycamore) # get streamflow object sf = asStreamflow(sycamore) # estimate catastrophic high flow variability sigmaHighFlows(sf)
Calculates catastrophic variability for low flow events. Negative residuals from the seasonal signal are used to calculate , the standard deviation of low-flow events.
sigmaLowFlows(x, resid.column)sigmaLowFlows(x, resid.column)
x |
An object of class |
resid.column |
Optional numeric specifiying which column contains residuals. Required if |
An object of class list with items
n.droughts |
Number of independent events with negative residuals. |
sigma.lfa |
Estimated sigma using the y-intercept. |
sigma.lfb |
Estimated sigma using the slope ( |
drought.line |
Matrix containing fitted, observed, and residual values from regression of log counts on bin midpoints. |
onesigma.events |
matrix containing information for all events below |
twosigma.events |
matrix containing information for all events below |
independentEvents
sigmaHighFlows
# load data data(sycamore) # get streamflow object sf = asStreamflow(sycamore) # estimate catastrophic low flow variability sigmaLowFlows(sf)# load data data(sycamore) # get streamflow object sf = asStreamflow(sycamore) # estimate catastrophic low flow variability sigmaLowFlows(sf)
A data series from USGS station 09510200, Sycamore Creek (Fort McDowell, AZ).
data(sycamore)data(sycamore)
The data are supplied as a matrix with years in the first column and average daily discharge in the second column. The series contains 18263 observations starting Jan 1, 1961 and ending Jan 1, 2011.
USGS National Water Information Sysem: http://waterdata.usgs.gov/usa/nwis/uv?09510200