forked from HongjianFang/global_tomo_taup
-
Notifications
You must be signed in to change notification settings - Fork 0
/
velocity_layer.py
148 lines (127 loc) · 5.25 KB
/
velocity_layer.py
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
139
140
141
142
143
144
145
146
147
148
#!/usr/bin/env python
# -*- coding: utf-8 -*-
"""
Functionality for dealing with a single velocity layer.
"""
from __future__ import (absolute_import, division, print_function,
unicode_literals)
from future.builtins import * # NOQA
from future.utils import native_str
import numpy as np
#: The VelocityLayer dtype stores a single layer. An entire velocity model is
#: implemented as an array of layers. The elements are:
#:
#: * ``top_depth``: The top depth of the layer.
#: * ``bot_depth``: The bottom depth of the layer.
#: * ``top_p_velocity``: The compressional (P) wave velocity at the top.
#: * ``bot_p_velocity``: The compressional (P) wave velocity at the bottom.
#: * ``top_s_velocity``: The shear (S) wave velocity at the top.
#: * ``bot_s_velocity``: The shear (S) wave velocity at the bottom.
#: * ``top_density``: The density at the top.
#: * ``bot_density``: The density at the bottom.
#: * ``top_qp``: The P wave attenuation at the top.
#: * ``bot_qp``: The P wave attenuation at the bottom.
#: * ``top_qs``: The S wave attenuation at the top.
#: * ``bot_qs``: The S wave attenuation at the bottom.
VelocityLayer = np.dtype([
(native_str('top_depth'), np.float_),
(native_str('bot_depth'), np.float_),
(native_str('top_p_velocity'), np.float_),
(native_str('bot_p_velocity'), np.float_),
(native_str('top_s_velocity'), np.float_),
(native_str('bot_s_velocity'), np.float_),
(native_str('top_density'), np.float_),
(native_str('bot_density'), np.float_),
(native_str('top_qp'), np.float_),
(native_str('bot_qp'), np.float_),
(native_str('top_qs'), np.float_),
(native_str('bot_qs'), np.float_),
])
def evaluate_velocity_at_bottom(layer, prop):
"""
Evaluate material properties at bottom of a velocity layer.
.. seealso:: :func:`evaluate_velocity_at_top`, :func:`evaluate_velocity_at`
:param layer: The velocity layer to use for evaluation.
:type layer: :class:`~numpy.ndarray`, dtype = :py:const:`.VelocityLayer`
:param prop: The material property to evaluate. One of:
* ``p``
Compressional (P) velocity (km/s)
* ``s``
Shear (S) velocity (km/s)
* ``r`` or ``d``
Density (g/cm^3)
:type prop: str
:returns: The value of the material property requested.
:rtype: :class:`~numpy.ndarray` (dtype = :class:`float`, shape equivalent
to ``layer``)
"""
prop = prop.lower()
if prop == "p":
return layer['bot_p_velocity']
elif prop == "s":
return layer['bot_s_velocity']
elif prop in "rd":
return layer['bot_density']
raise ValueError("Unknown material property, use p, s, or d.")
def evaluate_velocity_at_top(layer, prop):
"""
Evaluate material properties at top of a velocity layer.
.. seealso:: :func:`evaluate_velocity_at_bottom`,
:func:`evaluate_velocity_at`
:param layer: The velocity layer to use for evaluation.
:type layer: :class:`~numpy.ndarray`, dtype = :py:const:`.VelocityLayer`
:param prop: The material property to evaluate. One of:
* ``p``
Compressional (P) velocity (km/s)
* ``s``
Shear (S) velocity (km/s)
* ``r`` or ``d``
Density (g/cm^3)
:type prop: str
:returns: The value of the material property requested.
:rtype: :class:`~numpy.ndarray` (dtype = :class:`float`, shape equivalent
to ``layer``)
"""
prop = prop.lower()
if prop == "p":
return layer['top_p_velocity']
elif prop == "s":
return layer['top_s_velocity']
elif prop in "rd":
return layer['top_density']
raise ValueError("Unknown material property, use p, s, or d.")
def evaluate_velocity_at(layer, depth, prop):
"""
Evaluate material properties at some depth in a velocity layer.
.. seealso:: :func:`evaluate_velocity_at_top`,
:func:`evaluate_velocity_at_bottom`
:param layer: The velocity layer to use for evaluation.
:type layer: :class:`~numpy.ndarray`, dtype = :py:const:`.VelocityLayer`
:param depth: The depth at which the material property should be
evaluated. Must be within the bounds of the layer or results will be
undefined.
:type depth: float
:param prop: The material property to evaluate. One of:
* ``p``
Compressional (P) velocity (km/s)
* ``s``
Shear (S) velocity (km/s)
* ``r`` or ``d``
Density (g/cm^3)
:type prop: str
:returns: The value of the material property requested.
:rtype: :class:`~numpy.ndarray` (dtype = :class:`float`, shape equivalent
to ``layer``)
"""
thick = layer['bot_depth'] - layer['top_depth']
prop = prop.lower()
if prop == "p":
slope = (layer['bot_p_velocity'] - layer['top_p_velocity']) / thick
return slope * (depth - layer['top_depth']) + layer['top_p_velocity']
elif prop == "s":
slope = (layer['bot_s_velocity'] - layer['top_s_velocity']) / thick
return slope * (depth - layer['top_depth']) + layer['top_s_velocity']
elif prop in "rd":
slope = (layer['bot_density'] - layer['top_density']) / thick
return slope * (depth - layer['top_depth']) + layer['top_density']
raise ValueError("Unknown material property, use p, s, or d.")