IRISSeismic (version 1.0.5)

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)

Arguments

PSDs

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

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

References

http://pubs.usgs.gov/of/2005/1438/pdf/OFR-1438.pdf{Seismic Noise Analysis System Using Power Spectral Density Probability Density Functions} (McNamara and Boaz 2005)

See Also

McNamaraPSD, psdList, psdPlot

Examples

# 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(-250,-50), las=1,
     xlab="log10(period)", yalb="Power (dB)",
     main="Model 'quiet' 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)