Taup初探

Posted on | In

描述

What and WHere

README.md 开始信息
LICENSE GNU LGPL licence
StdModels 在/lib/TauP-2.6.1.jar/edu/sc/seis/TauP/StdModels中存在的模型

创建模型

.nd模型

模型中，每一列为：
depth（km）；P波速度（km/s）；S波速度(km/s)；密度（g/cm3）；Qp；Qs。

但是只有前三列为有效值，其余三列无效，只是为了兼容之后的某些情况。

注意！模型最后一行留一空白行。

使用英文单词表示分界面：

  • mantle 莫霍界面
  • outer-core CMB
  • inner-core ICB
    示例：
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
/* below is a simple named discontinuities model. */
0.0 5.0 3.0 2.7
20 5.0 3.0 2.7
20 6.5 3.7 2.9
33 6.5 3.7 2.9
mantle # the word "mantle" designates that this is the moho
33 7.8 4.4 3.3
410 8.9 4.7 3.5
410 9.1 4.9 3.7
670 10.2 5.5 4.0
670 10.7 5.9 4.4
2891 13.7 7.2 5.6
outer-core # "outer-core" designates that this is the core mantle boundary
2891 8.0 0.0 9.9
5149.5 10.3 0.0 12.2
inner-core # "inner-core" makes this the inner-outer core boundary
5149.5 11 3.5 12.7
6371 11.3 3.7 13

.vel模型

前两行为模型描述，两行之后为模型。
间断面使用两行相同的深度表示。每一行为：
depth；P波速度；S波速度；密度 （单位同上）

示例：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
ak135 - P
ak135 - S
     0.000      5.8000      3.4600     2.7200
    20.000      5.8000      3.4600      2.7200
    20.000      6.5000      3.8500      2.9200
    35.000      6.5000      3.8500      2.9200
    35.000      8.0400      4.4800      3.3198
    77.500      8.0450      4.4900      3.3455
   120.000      8.0500      4.5000      3.3713
   165.000      8.1750      4.5090      3.3985
   210.000      8.3000      4.5180      3.4258
   210.000      8.3000      4.5230      3.4258
   260.000      8.4825      4.6090      3.4561
   310.000      8.6650      4.6960      3.4864
   360.000      8.8475      4.7830      3.5167
   410.000      9.0300      4.8700      3.5470
   410.000      9.3600      5.0800      3.7557
   460.000      9.5280      5.1860      3.8175
   510.000      9.6960      5.2920      3.8793
   560.000      9.8640      5.3980      3.9410
   610.000     10.0320      5.5040      4.0028
   660.000     10.2000      5.6100      4.0646
   660.000     10.7900      5.9600      4.3714
   710.000     10.9229      6.0897      4.4010
   760.000     11.0558      6.2095      4.4305
   809.500     11.1353      6.2426      4.4596
   859.000     11.2221      6.2798      4.4885
   908.500     11.3068      6.3160      4.5173
   958.000     11.3896      6.3512      4.5459
  1007.500     11.4705      6.3854      4.5744
  1057.000     11.5495      6.4187      4.6028
  1106.500     11.6269      6.4510      4.6310
  1156.000     11.7026      6.4828      4.6591
  1205.500     11.7766      6.5138      4.6870
  1255.000     11.8491      6.5439      4.7148
  1304.500     11.9200      6.5727      4.7424
  1354.000     11.9895      6.6008      4.7699
  1403.500     12.0577      6.6285      4.7973
  1453.000     12.1245      6.6555      4.8245
  1502.500     12.1912      6.6815      4.8515
  1552.000     12.2550      6.7073      4.8785
  1601.500     12.3185      6.7326      4.9052
  1651.000     12.3819      6.7573      4.9319
  1700.500     12.4426      6.7815      4.9584
  1750.000     12.5031      6.8052      4.9847
  1799.500     12.5631      6.8286      5.0109
  1849.000     12.6221      6.8515      5.0370
  1898.500     12.6804      6.8742      5.0629
  1948.000     12.7382      6.8972      5.0887
  1997.500     12.7956      6.9194      5.1143
  2047.000     12.8526      6.9418      5.1398
  2096.500     12.9096      6.9627      5.1652
  2146.000     12.9668      6.9855      5.1904
  2195.500     13.0222      7.0063      5.2154
  2245.000     13.0783      7.0281      5.2403
  2294.500     13.1336      7.0500      5.2651
  2344.000     13.1894      7.0720      5.2898
  2393.500     13.2465      7.0931      5.3142
  2443.000     13.3018      7.1144      5.3386
  2492.500     13.3585      7.1369      5.3628
  2542.000     13.4156      7.1586      5.3869
  2591.500     13.4741      7.1807      5.4108
  2640.000     13.5312      7.2031      5.4345
  2690.000     13.5900      7.2258      5.4582
  2740.000     13.6494      7.2490      5.4817
  2740.000     13.6494      7.2490      5.4817
  2789.670     13.6530      7.2597      5.5051
  2839.330     13.6566      7.2704      5.5284
  2891.500     13.6602      7.2811      5.5515
  2891.500      8.0000      0.0000      9.9145
  2939.330      8.0382      0.0000      9.9942
  2989.660      8.1283      0.0000     10.0722
  3039.990      8.2213      0.0000     10.1485
  3090.320      8.3122      0.0000     10.2233
  3140.660      8.4001      0.0000     10.2964
  3190.990      8.4861      0.0000     10.3679
  3241.320      8.5692      0.0000     10.4378
  3291.650      8.6496      0.0000     10.5062
  3341.980      8.7283      0.0000     10.5731
  3392.310      8.8036      0.0000     10.6385
  3442.640      8.8761      0.0000     10.7023
  3492.970      8.9461      0.0000     10.7647
  3543.300      9.0138      0.0000     10.8257
  3593.640      9.0792      0.0000     10.8852
  3643.970      9.1426      0.0000     10.9434
  3694.300      9.2042      0.0000     11.0001
  3744.630      9.2634      0.0000     11.0555
  3794.960      9.3205      0.0000     11.1095
  3845.290      9.3760      0.0000     11.1623
  3895.620      9.4297      0.0000     11.2137
  3945.950      9.4814      0.0000     11.2639
  3996.280      9.5306      0.0000     11.3127
  4046.620      9.5777      0.0000     11.3604
  4096.950      9.6232      0.0000     11.40690
  4147.280      9.6673      0.0000     11.4521
  4197.610      9.7100      0.0000     11.4962
  4247.940      9.7513      0.0000     11.5391
  4298.270      9.7914      0.0000     11.5809
  4348.600      9.8304      0.0000     11.6216
  4398.930      9.8682      0.0000     11.6612
  4449.260      9.9051      0.0000     11.6998
  4499.600      9.9410      0.0000     11.7373
  4549.930      9.9761      0.0000     11.7737
  4600.260     10.0103      0.0000     11.8092
  4650.590     10.0439      0.0000     11.8437
  4700.920     10.0768      0.0000     11.8772
  4801.580     10.1415      0.0000     11.94140
  4851.910     10.1739      0.0000     11.9722
  4902.240     10.2049      0.0000     12.0001
  4952.580     10.2329      0.0000     12.0311
  5002.910     10.2565      0.0000     12.0593
  5053.240     10.2745      0.0000     12.0867
  5103.570     10.2854      0.0000     12.1133
  5153.500     10.2890      0.0000     12.1391
  5153.500     11.0427      3.5043     12.7037
  5204.610     11.0585      3.5187     12.7289
  5255.320     11.0718      3.5314     12.7530
  5306.040     11.0850      3.5435     12.7760
  5356.750     11.0983      3.5551     12.7980
  5407.460     11.1166      3.5661     12.8188
  5458.170     11.1316      3.5765     12.8387
  5508.890     11.1457      3.5864     12.8574
  5559.600     11.1590      3.5957     12.8751
  5610.310     11.1715      3.6044     12.8917
  5661.020     11.1832      3.6126     12.9072
  5711.740     11.1941      3.6202     12.9217
  5813.160     11.2134      3.6337     12.9474
  5863.870     11.2219      3.6396     12.9586
  5914.590     11.2295      3.6450     12.9688
  5965.300     11.2364      3.6498     12.9779
  6016.010     11.2424      3.6540     12.9859
  6066.720     11.2477      3.6577     12.9929
  6117.440     11.2521      3.6608     12.9988
  6168.150     11.2557      3.6633     13.0036
  6218.860     11.2586      3.6653     13.0074
  6269.570     11.2606      3.6667     13.0100
  6320.290     11.2618      3.6675     13.0117
  6371.000     11.2622      3.6678     13.0122

创建taup模型

创建完成模型后，可以通过语句完成taup模型创建

taup create -nd prem_normal.nd --verbose

生成的模型名称为prem_normal.taup

错误信息

创建模型时，注意地表处S波速度不能为0,否则报错。

1
Error starting tool: create edu.sc.seis.TauP.SlownessModelException: TauP Exception: SlownessModel Exception: Unable to handle zero S velocity layers at surface. This should be fixed at some point, but is a limitation of TauP at this point

S波速度与P波速度关系

S波速度小于P波速度。因此大于时会报错。

1
Caused by: edu.sc.seis.TauP.VelocityModelException: TauP Exception: VelocityModel Exception: S velocity, 3.19101 at depth 0.0 is greater than the P velocity, 1.44996

计算震相走时

交互输入

taup time --mod prem_normal.nd

示例：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
Enter:
h for new depth
r to recalculate
p to append phases, 
c to clear phases
l to list phases
s for new station lat lon
e for new event lat lon
a for new azimuth
b for new back azimuth
t for takeoff angle
m for new model or 
q to quit.

Enter Distance or Option [hrpclseabmqt]: h
Enter Depth: 140.0
Enter Distance or Option [hrpclseabmqt]: l
18 phases.
p,s,P,S,pP,sS,Pn,Sn,PcP,ScS,Pdiff,Sdiff,PKP,SKS,PKiKP,SKiKS,PKIKP,SKIKS
Enter Distance or Option [hrpclseabmqt]: 100

Model: prem_normal
Distance   Depth   Phase   Travel    Ray Param  Takeoff  Incident  Purist    Purist
  (deg)     (km)   Name    Time (s)  p (s/deg)   (deg)    (deg)   Distance   Name 
-----------------------------------------------------------------------------------
  100.00   140.0   Pdiff    812.04     4.437     18.61     3.32   100.00   = Pdiff
  100.00   140.0   PKiKP   1078.11     1.797      7.43     1.34   100.00   = PKiKP
  100.00   140.0   SKS     1443.09     4.955     11.38     2.55   100.00   = SKS  
  100.00   140.0   SKiKS   1502.74     1.888      4.31     0.97   100.00   = SKiKS
  100.00   140.0   S       1503.79     8.537     19.87     4.40   100.00   = S    
  100.00   140.0   ScS     1503.93     8.405     19.55     4.33   100.00   = ScS  
  100.00   140.0   sS      1571.60     8.584    160.02     4.43   100.00   = sS   

Enter Distance or Option [hrpclseabmqt]: q

也可以使用命令行命令直接进行计算：

-mod 输入模型名称
-h 震源深度 km
–stadepth 地震仪深度 km

-deg 震中距 度
-km 震中距 千米
或者给出经纬度（假设地球为球体）
-sta lat lon 台站经纬度
-evt lat lon 事件经纬度

–rayp 只输出射线参数
–time 只输出走时
–rel phasename 同时输出相对走时？
–json json格式输出走时

-ph 震相名称
-pf 震相文件，即包含震相名称的文件

示例：

taup time -mod prem_normal.nd -h 200 -deg 170 --rel P --json

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
{
  "model": "prem_normal",
  "sourcedepth": 200.0,
  "receiverdepth": 0.0,
  "phases": [ "p", "s", "P", "S", "pP", "sS", "Pn", "Sn", "PcP", "ScS", "Pdiff", "Sdiff", "PKP", "SKS", "PKiKP", "SKiKS", "PKIKP", "SKIKS" ],
  "arrivals": [
    {
      "distdeg": 170.0,
      "phase": "PKIKP",
      "time": 1186.1327,
      "rayparam": 0.5805277,
      "takeoff": 2.3948026,
      "incident": 0.43373135,
      "puristdist": 170.0,
      "puristname": "PKIKP"
    },
    {
      "distdeg": 170.0,
      "phase": "PKP",
      "time": 1261.7003,
      "rayparam": 4.42584,
      "takeoff": 18.575956,
      "incident": 3.308497,
      "puristdist": 170.0,
      "puristname": "PKP"
    },
    {
      "distdeg": 170.0,
      "phase": "SKIKS",
      "time": 1596.8145,
      "rayparam": 0.45553252,
      "takeoff": 1.0524957,
      "incident": 0.23472442,
      "puristdist": 170.0,
      "puristname": "SKIKS"
    }
  ]
}

截至此时，其只能在命令行中输出信息。

绘制速度剖面图

taup velplot -nd prem_normal.nd
sh prem_normal_vel.gmt

即可使用gmt绘制出漂亮的速度剖面图。

  • 模型类型
  • -mod 直接输入Taup中自带的模型名称。不需要调用其他文件。例如：iasp91 prem等。可以在安装目录StdModel文件夹中查看。
  • -nd 输入文件名称。
  • -tvel
  • –version 在命令行输出Taup版本。
  • –debug 可能是出现bug会在命令行有输出？目前还没有创建bug测试。
  • –verbose 未知

注意：必须已安装gmt。安装见《GMT安装》。

绘制模型慢度剖面图

绘制走时曲线

1
2
3
4
export model='premoonaver'
export phases='p,s,PKP,ScS,PcP,P,S'
taup_curve -mod $model.nd -h 867.0 -ph ${phases} -o time_curve -gmt
sh time_curve.gmt

注意：震相s与S不同，需要单独绘制。