Model the seismic spectrum due to bedload transport in rivers
The function calculates a seismic spectrum as predicted by the model of Tsai et al. (2012) for river bedload transport. The code was written to R by Sophie Lagarde and integrated to the R package 'eseis' by Michael Dietze.
model_bedload( gsd, d_s, s_s, r_s, q_s, h_w, w_w, a_w, f = c(1, 100), r_0, f_0, q_0, e_0, v_0, p_0, n_0, n_c, res = 100, adjust = TRUE, eseis = FALSE, ... )
gsd: data frame grain-size distribution function. Must be provided as data frame with two variables: grain-size class (in m, first column) and wgt/vol percentage per class (second column). See examples for details.d_s: Numeric value, mean sediment grain diameter (m). Alternative to gsd.s_s: Numeric value, standard deviation of sediment grain diameter (m). Alternative to gsd.r_s: Numeric value, specific sediment density (kg / m^3)q_s: Numeric value, unit sediment flux (m^2 / s)h_w: Numeric value, fluid flow depth (m)w_w: Numeric value, fluid flow width (m)a_w: Numeric value, fluid flow inclination angle (radians)f: Numeric vector, frequency range to be modelled. If of length two the argument is interpreted as representing the lower and upper limit and the final length of the frequency vector is set by the argument res. If f contains more than two values it is interpreted as the actual frequency vector and the value of res is ignored.r_0: Numeric value, distance of seismic station to sourcef_0: Numeric value, reference frequency (Hz)q_0: Numeric value, ground quality factor at f_0. "Reasonable value may be 20" (Tsai et al. 2012).e_0: Numeric value, exponent characterizing quality factor increase with frequency (dimensionless). "Reasonable value may be 0" (Tsai et al. 2012).v_0: Numeric value, phase speed of the Rayleigh wave at f_0 (m/s). Assuming a shear wave velocity of about 2200 m/s, Tsai et al. (2012) yield a value of 1295 m/s for this parameter.p_0: Numeric value, exponent of the power law variation of Rayleigh wave velocities with frequency. Attention, previous package versions used the argument name x_0, which is still silently supported for compatibility reasons, only giving a warning.n_0: Numeric vector of length two, Greens function displacement amplitude coefficients. Cf. N_ij in eq. 36 in Gimbert et al. (2014)n_c: Numeric value, option to include single particle hops coherent in time, causing spectrum modulation due to secondary effects. Omitted is no value is specified, here. Usual values may be between 2 and 4.res: Numeric value, output resolution, i.e. length of the spectrum vector. Default is 1000.adjust: Logical value, option to adjust PSD for wide grain-size distributions, according to implementation by Tsai et al. (2012).eseis: Character value, option to return an eseis object instead of a data frame. Default is FALSE....: Further arguments passed to the function.eseis object containing the modelled spectrum.
The model uses a set of predefined constants. These can also be changed by the user, using the ... argument:
g = 9.81, gravitational acceleration (m/s^2)r_w = 1000, fluid specific density (kg/m^3)k_s = 3 * d_50, roughness length (m)log_lim = c(0.0001, 100), limits of grain-size distribution functiontemplate (m)log_length = 10000, length of grain-size distribution function templatenu = 10^(-6), specific density of the fluid (kg/m^3)power_d = 3, grain-size power exponentgamma = 0.9, gamma parameter, after Parker (1990)s_c = 0.8, drag coefficient parameters_p = 3.5, drag coefficient parameterc_1 = 2 / 3, inter-impact time scaling, after Sklar & Dietrich (2004)When no user defined grain-size distribution function is provided,the function calculates the raised cosine distribution function as defined in Tsai et al. (2012) using the default range and resolution as specified by log_lim and log_length (see additional arguments list above). These default values are appropriate for mean sediment sizes between 0.001 and 10 m and log standard deivations between 0.05 and 1. When more extreme distributions are to be used, it is necessary to either adjust the arguments log_lim and log_length or use a user defined distribution function.
The adjustment option (implemented with package version 0.6.0) is only relevant for wide grain-size distributions, i.e., s_s > 0.2. In such cases, the unadjusted version tends to underestimate seismic power.
## calculate spectrum (i.e., fig. 1b in Tsai et al., 2012) p_bedload <- model_bedload(d_s = 0.7, s_s = 0.1, r_s = 2650, q_s = 0.001, h_w = 4, w_w = 50, a_w = 0.005, f = c(0.1, 20), r_0 = 600, f_0 = 1, q_0 = 20, e_0 = 0, v_0 = 1295, x_0 = 0.374, n_0 = 1, res = 100, eseis = TRUE) ## plot spectrum plot_spectrum(data = p_bedload, ylim = c(-170, -110)) ## define empiric grain-size distribution gsd_empiric <- data.frame(d = c(0.70, 0.82, 0.94, 1.06, 1.18, 1.30), p = c(0.02, 0.25, 0.45, 0.23, 0.04, 0.00)) ## calculate spectrum p_bedload <- model_bedload(gsd = gsd_empiric, r_s = 2650, q_s = 0.001, h_w = 4, w_w = 50, a_w = 0.005, f = c(0.1, 20), r_0 = 600, f_0 = 1, q_0 = 20, e_0 = 0, v_0 = 1295, x_0 = 0.374, n_0 = 1, res = 100, eseis = TRUE) ## plot spectrum plot_spectrum(data = p_bedload, ylim = c(-170, -110)) ## define mean and sigma for parametric distribution function d_50 <- 1 sigma <- 0.1 ## define raised cosine distribution function following Tsai et al. (2012) d_1 <- 10^seq(log10(d_50 - 5 * sigma), log10(d_50 + 5 * sigma), length.out = 20) sigma_star <- sigma / sqrt(1 / 3 - 2 / pi^2) p_1 <- (1 / (2 * sigma_star) * (1 + cos(pi * (log(d_1) - log(d_50)) / sigma_star))) / d_1 p_1[log(d_1) - log(d_50) > sigma_star] <- 0 p_1[log(d_1) - log(d_50) < -sigma_star] <- 0 p_1 <- p_1 / sum(p_1) gsd_raised_cos <- data.frame(d = d_1, p = p_1)
Tsai, V. C., B. Minchew, M. P. Lamb, and J.-P. Ampuero (2012), A physical model for seismic noise generation from sediment transport in rivers, Geophys. Res. Lett., 39, L02404, doi:10.1029/2011GL050255.
Sophie Lagarde, Michael Dietze
Useful links