NAME
spectrum1d - compute auto- [and cross- ] spectra from one
[or two] timeseries.
SYNOPSIS
spectrum1d [ x[y]file ] -Ssegment_size] [ -C[xycnpago] ] [
-Ddt ] [ -Nname_stem ] [ -V ] [ -W ] [ -bi[s][n] ] [ -bo[s]
]
DESCRIPTION
spectrum1d reads X [and Y] values from the first [and
second] columns on standard input [or x[y]file]. These
values are treated as timeseries X(t) [Y(t)] sampled at
equal intervals spaced dt units apart. There may be any
number of lines of input. spectrum1d will create file[s]
containing auto- [and cross- ] spectral density estimates by
Welch's method of ensemble averaging of multiple overlapped
windows, using standard error estimates from Bendat and
Piersol.
The output files have 3 columns: f or w, p, and e. f or w
is the frequency or wavelength, p is the spectral density
estimate, and e is the one standard deviation error bar
size. These files are named based on name_stem. If the -C
option is used, up to eight files are created; otherwise
only one (xpower) is written. The files (which are ASCII
unless -bo is set) are as follows:
name_stem.xpower
Power spectral density of X(t). Units of X * X * dt.
name_stem.ypower
Power spectral density of Y(t). Units of Y * Y * dt.
name_stem.cpower
Power spectral density of the coherent output. Units
same as ypower.
name_stem.npower
Power spectral density of the noise output. Units same
as ypower.
name_stem.gain
Gain spectrum, or modulus of the transfer function.
Units of (Y / X).
name_stem.phase
Phase spectrum, or phase of the transfer function.
Units are radians.
name_stem.admit
Admittance spectrum, or real part of the transfer
function. Units of (Y / X).
name_stem.coh
(Squared) coherency spectrum, or linear correlation
coefficient as a function of frequency. Dimensionless
number in [0, 1]. The Signal-to-Noise-Ratio (SNR) is
coh / (1 - coh). SNR = 1 when coh = 0.5.
REQUIRED ARGUMENTS
x[y]file
ASCII (or binary, see -bi) file holding X(t) [Y(t)]
samples in the first 1 [or 2] columns. If no file is
specified, spectrum1d will read from standard input.
-S segment_size is a radix-2 number of samples per window
for ensemble averaging. The smallest frequency
estimated is 1.0/(segment_size * dt), while the largest
is 1.0/(2 * dt). One standard error in power spectral
density is approximately 1.0 / sqrt(n_data /
segment_size), so if segment_size = 256, you need
25,600 data to get a one standard error bar of 10%.
Cross-spectral error bars are larger and more compli-
cated, being a function also of the coherency.
OPTIONS
-C Read the first two columns of input as samples of two
timeseries, X(t) and Y(t).
Consider Y(t) to be the output and X(t) the input in
a linear system with noise. Estimate the optimum f
requency response function by least squares, such that
the noise output is minimized and the coherent outpu t
and the noise output are uncorrelated. Optionally
specify up to 8 letters from the set { x y c n p a g o
} in any order to create only those output files
instead of the default [all]. x = xpower, y = ypower,
c = cpower, n = npower, p = phase, a = admit, g = gain,
o = coh.
-D dt Set the spacing between samples in the timeseries
[Default = 1].
-N name_stem Supply the name stem to be used for output
files [Default = "spectrum"].
-V Selects verbose mode, which will send progress reports
to stderr [Default runs "silently"].
-W Write Wavelength rather than frequency in column 1 of
the output file[s] [Default = frequency, (cycles /
dt)].
-bi Selects binary input. Append s for single precision
[Default is double]. Append n for the number of
columns in the binary file(s). [Default is 2 input
columns].
-bo Selects binary output. Append s for single precision
[Default is double].
EXAMPLES
Suppose data.g is gravity data in mGal, sampled every 1.5
km. To write its power spectrum, in mGal**2-km, to the file
data.xpower, try
spectrum1d data.g -S256 -D1.5 -Ndata
Suppose in addition to data.g you have data.t, which is
topography in meters sampled at the same points as data.g.
To estimate various features of the transfer function, con-
sidering data.t as input and data.g as output, try
paste data.t data.g | spectrum1d -S256 -D1.5 -Ndata -C
SEE ALSO
gmt(l), grdfft(l)
REFERENCES
Bendat, J. S., and A. G. Piersol, 1986, Random Data, 2nd
revised ed., John Wiley & Sons.
Welch, P. D., 1967, "The use of Fast Fourier Transform for
the estimation of power spectra: a method based on time
averaging over short, modified periodograms", IEEE Transac-
tions on Audio and Electroacoustics, Vol AU-15, No 2.