IRISSeismic (version 1.6.6)

psdStatistics: Return statistics for a set of PSDs

Description

The psdStatistics function calculates a variety of information associated with the incoming set of PSDs.

Usage

psdStatistics(PSDs, evalresp=NULL)

Value

A list of elements:

• noiseMatrix -- matrix of corrected power levels; rows=PSDs, columns=frequencies

• pdfMatrix -- matrix of probability density values; rows=dB level, columns=frequencies

• freq -- vector of frequencies associated statistics vectors and with matrix columns

• pdfBins -- vector of power values (dB) associated with pdfMatrix rows

• max -- maximum power level at each frequency

• min -- minimum power level at each frequency

• mean -- mean power level at each frequency

• median -- median power level at each frequency

• mode -- mode of power level at each frequency (obtained from pdfMatrix)

• nlnm -- low noise model power level at each frequency

• nhnm -- high noise model power level at each frequency

• pct_above -- percent of PSDs above the high noise model at each frequency

• pct_below -- percent of PSDS below the low noise model at each frequency

A variety of plots can be generated form the information in this list.

Arguments

PSDs

either a list as returned by psdList or a dataframe of PSD values obtained from the BSS

evalresp

dataframe of freq, amp, phase information matching output of getEvalresp, optional

Author

Jonathan Callahan

References

Seismic Noise Analysis System Using Power Spectral Density Probability Density Functions (McNamara and Boaz 2005)

See Also

McNamaraPSD, psdList, psdPlot

Examples

if (FALSE) {
# Create a new IrisClient
iris <- new("IrisClient", debug=TRUE)

# Get seismic data
starttime <- as.POSIXct("2011-05-05", tz="GMT") # 2011.125
endtime <- starttime + 1*24*3600
st <- getDataselect(iris,"IU","GRFO","--","BHE",starttime,endtime)

# Generate power spectral density for each hour long segment
psdList <- psdList(st)

# Generate Statistics
stats <- psdStatistics(psdList)

# Just for fun plot
logPeriod <- log10(1/stats$freq)
plot(logPeriod,stats$max,ylim=c(-200,-50), las=1,
     xlab="log10(period)", ylab="Power (dB)",
     main="Model 'normal background noise' area and area of seismic signal.")
points(logPeriod,stats$min)

# Overlay a polygon showing the range between the noise models
x <- c(logPeriod,rev(logPeriod),logPeriod[1])
y <- c(stats$nhnm,rev(stats$nlnm),stats$nhnm[1])
transparentBlack <- adjustcolor('black',0.4)
polygon(x,y,col=transparentBlack)

# Overlay a polygon showing the range of measured values
y <- c(stats$max,rev(stats$min),stats$max[1])
transparentBlue <- adjustcolor('blue',0.6)
polygon(x,y,col=transparentBlue)  
}