spectrum1d
- 官方文档
- 简介
计算一个时间序列的自功率谱,或两个时间序列的互功率谱
spectrum1d
从标准输入流或数据文件中读取一列或两列数据。这些数据被当作是采样间隔为 <dt>
的等间隔的时间序列。spectrum1d 采用Welch方法,即加窗多段平均周期图法,计算输出自功率或互功率谱密度。其输出的功率谱的标准差,是用
Bendat和Piersol提供的算法。
spectrum1d
的输出文件有三列:
f|w p e
其中, f代表频率,w代表波长,p代表计算的功率谱密度,e代表一个标准差的值。
spectrum1d
的输出文件的文件名是使用统一的前缀 name_stem
。如果使用了
-C
选项,那么将会有8个文件输出,否则只生成一个功率谱文件( .xpower
)。这些文件默认是以ASCII码格式,除非用 -bo
选项指定为二进制格式输出。这8个文件介绍如下:
name_stem.xpower
: X(t)的功率谱。单位是X*X*dt
。name_stem.ypower
: Y(t)的功率谱。单位是Y*Y*dt
。name_stem.cpower
: 一致性(coherent)的功率谱。单位和 ypower 一样。name_stem.npower
: 噪声的功率谱。单位和 ypower 一样。name_stem.gain
: 增益谱,或传输函数的模。单位是Y/X
name_stem.phase
: 相位谱,或传输函数的相位。单位是弧度。name_stem.admit
: 导纳(Admittance)谱,或传输函数的实部。单位是Y/X
。name_stem.coh
: (平方)相干谱,或者线性相关系数(它是频率的函数)。i 无单位,取值范围为[0,1]
。信噪比SNR=coh/(1-coh)
。当 coh=0.5 时,SNR=1。
除非使用 -T
选项,否则以上文件会以单个文件单列的形式输出。
选项
-S<segment_size>
<segment_size>
是一个2的指数数值,用于控制Welch方法中分段平均时的窗口长度。它也决定了功率谱密度的最小频率分辨率和最大频率分辨率,即1.0/(segment_size*dt)
和1.0/(2*dt)
(即Nyquist频率)。在功率谱密度中的一个标准误差大约为1.0/(n_data/segment_size)
,比如 segment_size=256,那么就需要25600个数据点去计算一个误差棒的10%。互功率谱误差棒的计算则需要更多数据点,而且是相干性的函数,比较复杂。table
输入文件名。它是ASCII类型的一列数据或两列数据。如果是一列数据文件,就计算自功率谱;如果是两列,就计算互功率谱。若未指定文件名,
spectrum1d
会从标准输入流中读取数据。-C[xycnpago]
默认会输出全部8个文件。使用该选项可以指定输出8个文件中的某些文件。x=xpower、y=ypower、c=cpower、n=npower、p=phase、a=admit、g=gain、o=coh。
-D<dt>
设置读入的时间序列的时间采样间隔,默认值是1。
-L[m|h]
不去除信号中的线性趋势。默认情况下,在对信号进行变换处理前会先去掉其中的线性趋势。
m
表示去掉数据的均值,h
表示去掉数据的中值。-N[name_stem]
输出文件名的前缀,默认为
spectrum
。若不使用此选项,则输出的8个文件会合到一个文件里。-T
不让单个分量的结果输出到标准输出流。
-W
输出文件中第一列是波长而不是频率。默认输出时第一列是频率。
示例
假设 data.g
是重力数据,单位为 mGal,空间采样间隔为1.5 km。如下命令会输出数据的功率谱,单位为 mGal^2 km 表示:
gmt spectrum1d data.g -S256 -D1.5 -Ndata
假设你除了有重力数据 data.g
之外,还有在相同地点测得的地形数据 data.t
,单位为 m。计算二者之间的传输函数,即 data.t
是输入, data.g
是输出:
paste data.t data.g | gmt spectrum1d -S256 -D1.5 -Ndata -C > results.txt