GF3D/make_station_section.py at main · lsawade/GF3D · GitHub

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
#!/bin/env python
"""
Station Section
===============

The tutorial will go over the reading of a subset file by loading one that is
included in the directory. At the end we plot a station section of waveforms
using built-in plotting tools.

Loading all modules
-------------------


"""
# sphinx_gallery_thumbnail_number = 1
# sphinx_gallery_dummy_images = 1

# External
import matplotlib.pyplot as plt
from obspy import read, read_inventory, Stream

# Internal
from gf3d.source import CMTSOLUTION
from gf3d.seismograms import GFManager
from gf3d.process import process_stream, select_pairs
from gf3d.plot.section import plotsection

# %%

# CMTSOLUTION
cmt = CMTSOLUTION.read('examples/DATA/single_element_read/CMTSOLUTION')

# Load subset
gfsub = GFManager("temp_subset.h5")
gfsub.load()

# Load Observed Data
raw = read("examples/DATA/single_element_read/traces/*.sac")
inv = read_inventory("examples/DATA/single_element_read/station.xml")

# %% Get seismograms from the database

rp = gfsub.get_seismograms(cmt)


# %% Process

obs = process_stream(raw, inv=inv, cmt=cmt, duration=10800)
syn = process_stream(rp, cmt=cmt, duration=10800)

# %%

pobs, psyn = select_pairs(obs, syn)

# %% Plot section

starttime = psyn[0].stats.starttime + 0
endtime = starttime + 10800
limits = (starttime, endtime)

# Plots a section of observed and synthetic
plotsection(pobs, psyn, cmt, comp='Z', lw=0.75, limits=limits)

plt.show()