greenfn
- 简介:
使用广义反射透射系数矩阵法计算动态全波解格林函数
语法
grt greenfn -Mmodel -Ddepsrc/deprcv -Nnt/dt[+wzeta][+nfac][+a][+f] -Rr1/r2/dr|r1,r2,…|file -Ooutdir [ -Bf|F|r|R|h|H ] [ -Hf1/f2 ] [ -Llength[+lFlength][+aFtol][+ooffset] ] [ -Cd|p|n ] [ -K[+kk0][+sampk][+ekeps][+vvmin] ] [ -E[p]t0[/v0] ] [ -Pnthreads ] [ -Ge|v|h|s ] [ -S[i1,i2,…] ] [ -e ] [ -s ] [ -h ]
描述
greenfn 模块计算每个震中距 r 的格林函数波形保存路径为
{outdir}/{model}_{depsrc}_{deprcv}_{r}/{stype}.sac ，
其中支持的格林函数类型 stype 详见 格林函数分类 。
执行 greenfn 模块的原始命令会附加式地保存到 {outdir}/command 文件。
不同震源的格林函数单位为：
爆炸源： \(10^{-20} \, \frac{\text{cm}}{\text{dyne} \cdot \text{cm}}\)
单力源： \(10^{-15} \, \frac{\text{cm}}{\text{dyne}}\)
剪切源： \(10^{-20} \, \frac{\text{cm}}{\text{dyne} \cdot \text{cm}}\)
必选选项
- -Mmodel
半空间层状模型文件的路径。模型格式如下:
Thickness(km) Vp(km/s) Vs(km/s) Rho(g/cm^3) [Qp] [Qs]
后两列的 Qp, Qs 可省略。
如果首列首行的值为正数，则首列表示 每层厚度 (km) ，最后一层为半空间，其厚度设置为 0 即可（不会被使用）。
如果首列首行的值为 0，则首列表示 每层层顶深度 (km) 。
如果某层设置 Vs = 0，则该层为液体层。
- -Ddepsrc/deprcv
震源深度 depsrc (km) 和台站深度 deprcv (km)。 如果是在使用波数积分法求解格林函数，则当二者深度差小于 1 km 时，自动使用快速收敛算法。
- -Nnt/dt[+wzeta][+nfac][+a][+f]
采样点数 nt 和采样时间间隔 dt (secs) ，这将决定计算的最高频。还可设置：
+wzeta - 虚频率系数 \(\zeta\) [0.8]。 Bouchon, 1981 提出在频率上添加一个虚频率偏移， \(\omega \leftarrow \omega - i \zeta \dfrac{\pi}{T}\) ， 其中 \(T\) 为时窗长度 nt*dt，以使波数积分适当偏移实轴上的极点。
+nfac - 频率域插值倍数[1]。即在做逆傅里叶变换时，在频域上最高频后补零， 相当于 nt ← nt * fac , dt ← dt / fac ，使计算的波形更平滑。
+a - 计算所有频点，不论频率多低。除非做数值实验，否则不建议使用该选项。 默认情况下，程序会跳过非常低频的几个点以避免引入误差， 详见 解决“波形漂移”问题 的介绍。
+f - 不进行由虚频率引起的振幅补偿。 默认情况下，程序会对时间域结果乘上 \(\exp(\frac{\zeta \pi}{T} t)\) 进行振幅补偿， 但当时窗长度远超有效信号的长度时，这种补偿会放大尾部噪声。
注： 当时窗长度 nt*dt 太小“包不住”有效信号，或时窗长度足够但时延(-E)不合适，输出的波形会发生混叠， 此时需调整 -N 和 -E 。
- -Rr1/r2/dr|r1,r2,…|file
多个震中距 (km)。有三种设置方式：
-Rr1/r2/dr - 在 [r1, r2] 范围内设置间隔为 dr 的一系列等距点。
-Rr1,r2,… - 逗号分隔的多个震中距
-Rfile - 仅含一列震中距数值的文件 file 。
- -Ooutdir
输出目录名，不存在会自动新建。若存在同名保存路径则直接覆盖已有结果。
可选选项
- -Bf|F|r|R|h|H
设置边界条件，小写和大写字母分别控制顶界面和底界面的边界条件， 分别为自由边界(f|F)，刚性边界(r|R)和半空间(h|H)， 默认顶界面自由边界，底界面半空间。 用法类似 -BrH , -BhF 等。
- -Hf1/f2
待计算的频率范围（闭区间）(Hz) [默认计算由 -D 指定的全部范围]。 f1 和 f2 分别为最小最大频率，若取 -1 则分别对应低通和高通。
- -Llength[+lFlength][+aFtol][+ooffset]
控制波数积分间隔。默认使用传统的离散波数积分，自动确定 length ， 准则详见 Bouchon (1981) 或代码源文件， 对应波数间隔 \(\Delta k =\frac{2\pi}{\textit{<length>}\cdot R_{\text{max}}}\) ， 其中 \(R_{\text{max}}\) 为多个震中距的最大值。
支持以下子选项：
+lFlength - 使用 固定间隔的Filon积分法 ，间隔由 Flength 按相同公式给定。
+aFtol - 使用 自适应Filon积分法 ， 收敛精度 Ftol 一般取 1e-2 即可。
+ooffset - 将积分区间 \([0, k_{\text{max}}]\) 划分为两段， \([0, k^*]\) 和 \([k^*, k_{\text{max}}]\) ， 其中 \(k^*=\frac{\text{<offset>}}{R_{\text{max}}}\) 前段仍使用离散波数法求解积分，间隔由 length 控制， 后段则使用以上用户设置的方法。默认 offset 为0。
例如：
-L20 - 手动设置离散波数积分间隔。
-L+a1e-2 - 在全区间 \([0, k_{\text{max}}]\) 使用自适应Filon积分。
-L20+a1e-2+o10 - 分两段计算积分。
备注
固定间隔和自适应 Filon 积分不计算近场项，即 计算动态格林函数 中介绍的 \(p=1\) 对应的积分项。
- -K[+kk0][+sampk][+ekeps][+vvmin]
控制波数积分上限
\[k_{\text{max}} = \sqrt{ k_0 \cdot \dfrac{\pi}{\Delta h} + \textit{<ampk>} \cdot \left(\dfrac{\omega}{v_{\text{min}}}\right)^2}\]+kk0 - 控制零频的积分上限 [5.0]，其中深度差 \(\Delta h = \max(|z_s - z_r|, 1.0)\) 。
+sampk - 放大倍数 [1.15] 。
+ekeps - 用于判断提前结束波数积分的收敛精度[0.0, 默认不使用]，详见 Yao and Harkrider (1983) 和 控制波数积分上限 。
+vvmin - 参考最小速度，默认 \(\max{(\min\limits_{i} (\alpha_i \cup \beta_i), 0.1)}\) 。
- -Cd|p|n
设置波数积分收敛方法。默认当 \(|z_s - z_r| <= 1.0\) km 时，自动使用直接收敛法。 支持以下子选项：
n - 强制禁用任何收敛方法。
- -E[p]t0[/v0]
增加时延以调整波形起始时刻。输出的波形默认起始时刻为 0 时刻，即发震时刻。 对于大震中距的记录，P 波初至走时较大，此时可以适当调整波形起始时刻， 使 -N 定义的时窗长度能框住整段波形， 而不必使用 -N 定义从发震时刻开始的超宽时窗，增加不必要的计算量。支持以下用法：
-Et0 - 起始时刻调整为 \(t_0\) (secs) ， e.g. -E20 。
-Et0/v0 - 每个震中距 \(r\) 的记录起始时刻不同， 调整为 \(t_0 + \frac{r}{v_0}\) ， v0 为参考速度 (km/s) ， e.g. -E-10/9。
-Ept0 - 每个震中距 \(r\) 的记录起始时刻不同， 调整为 \(t_0 + t_p\) ，其中 \(t_p\) 为每个震中距的初至 P 到时， e.g. -Ep-20 。
- -Pnthreads
多线程数，默认使用全部线程。
- -Ge|v|h|s
控制输出的震源类型 [evhs]，不会减少计算量。
e - 爆炸源
v - 垂直力源
h - 水平力源
s - 剪切源
- -S[i1,i2,…]
指定频率索引值，输出对应频率下波数积分过程中的核函数文件，多个索引值用逗号分隔，文件保存目录为
{outdir}_grtstats。若使用 -S 则保存所有频点。 关于文件格式及其读取详见 积分收敛性 。注意：
如果使用了波数积分收敛方法(-C)，输出的核函数会有变动。 使用 -Cd 会保存改动后的核函数；使用 -Cp 会额外输出应用 PTAM 过程中的核函数和积分峰谷值。
如果使用 -e 计算格林函数空间导数，则输出核函数对 z 的偏导。
- -e
在计算格林函数的同时，也计算其空间导数 \(\dfrac{\partial (u_z, u_r, u_\theta)}{\partial (z, r)}\) 。 偏导对应在文件名/变量名开头添加了
z和r。 关于 \(\theta\) 的偏导与方向因子有关，这将由 syn 和 static_syn 模块计算。不同震源的格林函数空间导数单位为:
爆炸源： \(10^{-25} \, /\text{dyne} \cdot \text{cm}\)
单力源： \(10^{-20} \, /\text{dyne}\)
剪切源： \(10^{-25} \, /\text{dyne} \cdot \text{cm}\)
- -s
静默输出，不在终端打印参数和进度条。
- -h
打印帮助文档。
示例
详见教程：