In [254]:
matplotlib inline
In [255]:
from obspy.core import read
from obspy import UTCDateTime
#from obspy.xseed import Parser
#from obspy.signal.rotate import rotate2ZNE
from obspy.core import Stream, Trace, read, AttribDict
from obspy.signal.invsim import paz_2_amplitude_value_of_freq_resp
In [256]:
import matplotlib.pyplot as plt
import pandas as pd
In [257]:
from __future__ import print_function
%matplotlib inline
import matplotlib.pyplot as plt
#plt.style.use('ggplot')
plt.rcParams['figure.figsize'] = 18, 6
In [258]:
import numpy as np
import pandas as pd
from obspy import UTCDateTime
import matplotlib.pyplot as plt
import math
from scipy.optimize import curve_fit
import numpy.random as npr
import os.path
In [259]:
#print(np.__version__)
In [260]:
import obspy as op
#print(op.__version__)
In [261]:
import sys
#print(sys.version)
In [262]:
import warnings
warnings.filterwarnings('ignore')
In [263]:
from IPython.display import HTML
HTML('''<script>
code_show=true; 
function code_toggle() {
 if (code_show){
 $('div.input').hide();
 } else {
 $('div.input').show();
 }
 code_show = !code_show
} 
$( document ).ready(code_toggle);
</script>
<form action="javascript:code_toggle()"><input type="submit" value="Click here to toggle on/off the raw code."></form>''')
Out[263]:
In [264]:
hnm_pdf =pd.read_csv("http://ncedc.org/ftp/outgoing/taira/WQC/ByerlyVault_PSD/HNM.dat", sep=" ",header=None)
In [265]:
lnm_pdf =pd.read_csv("http://ncedc.org/ftp/outgoing/taira/WQC/ByerlyVault_PSD/LNM.dat", sep=" ",header=None)
In [266]:
data_set_name = "run_071919to091719"
In [267]:
bb_100s =pd.read_csv("http://ncedc.org/ftp/outgoing/taira/WQC/ByerlyVault_PSD/"+data_set_name+"/bktest_BB_0.009413_.out", 
                          sep=" ",names=["sta", "zval_med", "nval_med", "eval_med", "nz_val_med", "ez_val_med", "en_val_med"],
                          header=None, skiprows=1)
In [268]:
#bb_200s =pd.read_csv("http://ncedc.org/ftp/outgoing/taira/WQC/ByerlyVault_PSD/bktest_BB_0.005132_.out", 
#                          sep=" ",names=["sta", "zval_med", "nval_med", "eval_med", "nz_val_med", "ez_val_med", "en_val_med"],
#                          header=None, skiprows=1)
In [269]:
bb_50s =pd.read_csv("http://ncedc.org/ftp/outgoing/taira/WQC/ByerlyVault_PSD/"+data_set_name+"/bktest_BB_0.020530_.out", 
                          sep=" ",names=["sta", "zval_med", "nval_med", "eval_med", "nz_val_med", "ez_val_med", "en_val_med"],
                          header=None, skiprows=1)
In [270]:
bb_1s =pd.read_csv("http://ncedc.org/ftp/outgoing/taira/WQC/ByerlyVault_PSD/"+data_set_name+"/bktest_BB_1.013156_.out", 
                          sep=" ",names=["sta", "zval_med", "nval_med", "eval_med", "nz_val_med", "ez_val_med", "en_val_med"],
                          header=None, skiprows=1)
In [271]:
bks_hhz_psd =pd.read_csv("http://ncedc.org/ftp/outgoing/taira/WQC/ByerlyVault_PSD/"+data_set_name+"/freq.stack_BKS.BK.HHZ.00_.out", 
                          sep=" ",names=["freq", "med", "mean", "mode", "l95", "h95", "min", "max", "sd","var","sta","net","com","loc"],
                          header=None, skiprows=1)
bks_hhz_psd_nona = bks_hhz_psd.dropna()
bks_hhn_psd =pd.read_csv("http://ncedc.org/ftp/outgoing/taira/WQC/ByerlyVault_PSD/"+data_set_name+"/freq.stack_BKS.BK.HHN.00_.out", 
                          sep=" ",names=["freq", "med", "mean", "mode", "l95", "h95", "min", "max", "sd","var","sta","net","com","loc"],
                          header=None, skiprows=1)
bks_hhn_psd_nona = bks_hhn_psd.dropna()
bks_hhe_psd =pd.read_csv("http://ncedc.org/ftp/outgoing/taira/WQC/ByerlyVault_PSD/"+data_set_name+"/freq.stack_BKS.BK.HHE.00_.out", 
                          sep=" ",names=["freq", "med", "mean", "mode", "l95", "h95", "min", "max", "sd","var","sta","net","com","loc"],
                          header=None, skiprows=1)
bks_hhe_psd_nona = bks_hhe_psd.dropna()
In [272]:
bk82_hhz_psd =pd.read_csv("http://ncedc.org/ftp/outgoing/taira/WQC/ByerlyVault_PSD/"+data_set_name+"/freq.stack_BK82.BK.HHZ.00_.out", 
                          sep=" ",names=["freq", "med", "mean", "mode", "l95", "h95", "min", "max", "sd","var","sta","net","com","loc"],
                          header=None, skiprows=1)
bk82_hhz_psd_nona = bk82_hhz_psd.dropna()
bk82_hhn_psd =pd.read_csv("http://ncedc.org/ftp/outgoing/taira/WQC/ByerlyVault_PSD/"+data_set_name+"/freq.stack_BK82.BK.HHN.00_.out", 
                          sep=" ",names=["freq", "med", "mean", "mode", "l95", "h95", "min", "max", "sd","var","sta","net","com","loc"],
                          header=None, skiprows=1)
bk82_hhn_psd_nona = bk82_hhn_psd.dropna()
bk82_hhe_psd =pd.read_csv("http://ncedc.org/ftp/outgoing/taira/WQC/ByerlyVault_PSD/"+data_set_name+"/freq.stack_BK82.BK.HHE.00_.out", 
                          sep=" ",names=["freq", "med", "mean", "mode", "l95", "h95", "min", "max", "sd","var","sta","net","com","loc"],
                          header=None, skiprows=1)
bk82_hhe_psd_nona = bk82_hhe_psd.dropna()
In [273]:
bk80_hhz_psd =pd.read_csv("http://ncedc.org/ftp/outgoing/taira/WQC/ByerlyVault_PSD/"+data_set_name+"/freq.stack_BK80.BK.HHZ.00_.out", 
                          sep=" ",names=["freq", "med", "mean", "mode", "l95", "h95", "min", "max", "sd","var","sta","net","com","loc"],
                          header=None, skiprows=1)
bk80_hhz_psd_nona = bk80_hhz_psd.dropna()
bk80_hhn_psd =pd.read_csv("http://ncedc.org/ftp/outgoing/taira/WQC/ByerlyVault_PSD/"+data_set_name+"/freq.stack_BK80.BK.HHN.00_.out", 
                          sep=" ",names=["freq", "med", "mean", "mode", "l95", "h95", "min", "max", "sd","var","sta","net","com","loc"],
                          header=None, skiprows=1)
bk80_hhn_psd_nona = bk80_hhn_psd.dropna()
bk80_hhe_psd =pd.read_csv("http://ncedc.org/ftp/outgoing/taira/WQC/ByerlyVault_PSD/"+data_set_name+"/freq.stack_BK80.BK.HHE.00_.out", 
                          sep=" ",names=["freq", "med", "mean", "mode", "l95", "h95", "min", "max", "sd","var","sta","net","com","loc"],
                          header=None, skiprows=1)
bk80_hhe_psd_nona = bk80_hhe_psd.dropna()
In [274]:
bk69_hhz_psd =pd.read_csv("http://ncedc.org/ftp/outgoing/taira/WQC/ByerlyVault_PSD/"+data_set_name+"/freq.stack_BK69.BK.HHZ.00_.out", 
                          sep=" ",names=["freq", "med", "mean", "mode", "l95", "h95", "min", "max", "sd","var","sta","net","com","loc"],
                          header=None, skiprows=1)
bk69_hhz_psd_nona = bk69_hhz_psd.dropna()
bk69_hhn_psd =pd.read_csv("http://ncedc.org/ftp/outgoing/taira/WQC/ByerlyVault_PSD/"+data_set_name+"/freq.stack_BK69.BK.HHN.00_.out", 
                          sep=" ",names=["freq", "med", "mean", "mode", "l95", "h95", "min", "max", "sd","var","sta","net","com","loc"],
                          header=None, skiprows=1)
bk69_hhn_psd_nona = bk69_hhn_psd.dropna()
bk69_hhe_psd =pd.read_csv("http://ncedc.org/ftp/outgoing/taira/WQC/ByerlyVault_PSD/"+data_set_name+"/freq.stack_BK69.BK.HHE.00_.out", 
                          sep=" ",names=["freq", "med", "mean", "mode", "l95", "h95", "min", "max", "sd","var","sta","net","com","loc"],
                          header=None, skiprows=1)
bk69_hhe_psd_nona = bk69_hhe_psd.dropna()
In [275]:
bk70_hhz_psd =pd.read_csv("http://ncedc.org/ftp/outgoing/taira/WQC/ByerlyVault_PSD/"+data_set_name+"/freq.stack_BK70.BK.HHZ.00_.out", 
                          sep=" ",names=["freq", "med", "mean", "mode", "l95", "h95", "min", "max", "sd","var","sta","net","com","loc"],
                          header=None, skiprows=1)
bk70_hhz_psd_nona = bk70_hhz_psd.dropna()
bk70_hhn_psd =pd.read_csv("http://ncedc.org/ftp/outgoing/taira/WQC/ByerlyVault_PSD/"+data_set_name+"/freq.stack_BK70.BK.HHN.00_.out", 
                          sep=" ",names=["freq", "med", "mean", "mode", "l95", "h95", "min", "max", "sd","var","sta","net","com","loc"],
                          header=None, skiprows=1)
bk70_hhn_psd_nona = bk70_hhn_psd.dropna()
bk70_hhe_psd =pd.read_csv("http://ncedc.org/ftp/outgoing/taira/WQC/ByerlyVault_PSD/"+data_set_name+"/freq.stack_BK70.BK.HHE.00_.out", 
                          sep=" ",names=["freq", "med", "mean", "mode", "l95", "h95", "min", "max", "sd","var","sta","net","com","loc"],
                          header=None, skiprows=1)
bk70_hhe_psd_nona = bk70_hhe_psd.dropna()
In [276]:
bk71_hhz_psd =pd.read_csv("http://ncedc.org/ftp/outgoing/taira/WQC/ByerlyVault_PSD/"+data_set_name+"/freq.stack_BK71.BK.HHZ.00_.out", 
                          sep=" ",names=["freq", "med", "mean", "mode", "l95", "h95", "min", "max", "sd","var","sta","net","com","loc"],
                          header=None, skiprows=1)
bk71_hhz_psd_nona = bk71_hhz_psd.dropna()
bk71_hhn_psd =pd.read_csv("http://ncedc.org/ftp/outgoing/taira/WQC/ByerlyVault_PSD/"+data_set_name+"/freq.stack_BK71.BK.HHN.00_.out", 
                          sep=" ",names=["freq", "med", "mean", "mode", "l95", "h95", "min", "max", "sd","var","sta","net","com","loc"],
                          header=None, skiprows=1)
bk71_hhn_psd_nona = bk71_hhn_psd.dropna()
bk71_hhe_psd =pd.read_csv("http://ncedc.org/ftp/outgoing/taira/WQC/ByerlyVault_PSD/"+data_set_name+"/freq.stack_BK71.BK.HHE.00_.out", 
                          sep=" ",names=["freq", "med", "mean", "mode", "l95", "h95", "min", "max", "sd","var","sta","net","com","loc"],
                          header=None, skiprows=1)
bk71_hhe_psd_nona = bk71_hhe_psd.dropna()
In [277]:
bk72_hhz_psd =pd.read_csv("http://ncedc.org/ftp/outgoing/taira/WQC/ByerlyVault_PSD/"+data_set_name+"/freq.stack_BK72.BK.HHZ.00_.out", 
                          sep=" ",names=["freq", "med", "mean", "mode", "l95", "h95", "min", "max", "sd","var","sta","net","com","loc"],
                          header=None, skiprows=1)
bk72_hhz_psd_nona = bk72_hhz_psd.dropna()
bk72_hhn_psd =pd.read_csv("http://ncedc.org/ftp/outgoing/taira/WQC/ByerlyVault_PSD/"+data_set_name+"/freq.stack_BK72.BK.HHN.00_.out", 
                          sep=" ",names=["freq", "med", "mean", "mode", "l95", "h95", "min", "max", "sd","var","sta","net","com","loc"],
                          header=None, skiprows=1)
bk72_hhn_psd_nona = bk72_hhn_psd.dropna()
bk72_hhe_psd =pd.read_csv("http://ncedc.org/ftp/outgoing/taira/WQC/ByerlyVault_PSD/"+data_set_name+"/freq.stack_BK72.BK.HHE.00_.out", 
                          sep=" ",names=["freq", "med", "mean", "mode", "l95", "h95", "min", "max", "sd","var","sta","net","com","loc"],
                          header=None, skiprows=1)
bk72_hhe_psd_nona = bk72_hhe_psd.dropna()
In [278]:
bk63_hhz_psd =pd.read_csv("http://ncedc.org/ftp/outgoing/taira/WQC/ByerlyVault_PSD/"+data_set_name+"/freq.stack_BK63.BK.HHZ.--_.out", 
                          sep=" ",names=["freq", "med", "mean", "mode", "l95", "h95", "min", "max", "sd","var","sta","net","com","loc"],
                          header=None, skiprows=1)
bk63_hhz_psd_nona = bk63_hhz_psd.dropna()
bk63_hhn_psd =pd.read_csv("http://ncedc.org/ftp/outgoing/taira/WQC/ByerlyVault_PSD/"+data_set_name+"/freq.stack_BK63.BK.HHN.--_.out", 
                          sep=" ",names=["freq", "med", "mean", "mode", "l95", "h95", "min", "max", "sd","var","sta","net","com","loc"],
                          header=None, skiprows=1)
bk63_hhn_psd_nona = bk63_hhn_psd.dropna()
bk63_hhe_psd =pd.read_csv("http://ncedc.org/ftp/outgoing/taira/WQC/ByerlyVault_PSD/"+data_set_name+"/freq.stack_BK63.BK.HHE.--_.out", 
                          sep=" ",names=["freq", "med", "mean", "mode", "l95", "h95", "min", "max", "sd","var","sta","net","com","loc"],
                          header=None, skiprows=1)
bk63_hhe_psd_nona = bk63_hhe_psd.dropna()
In [279]:
bk64_hhz_psd =pd.read_csv("http://ncedc.org/ftp/outgoing/taira/WQC/ByerlyVault_PSD/"+data_set_name+"/freq.stack_BK64.BK.HHZ.--_.out", 
                          sep=" ",names=["freq", "med", "mean", "mode", "l95", "h95", "min", "max", "sd","var","sta","net","com","loc"],
                          header=None, skiprows=1)
bk64_hhz_psd_nona = bk64_hhz_psd.dropna()
bk64_hhn_psd =pd.read_csv("http://ncedc.org/ftp/outgoing/taira/WQC/ByerlyVault_PSD/"+data_set_name+"/freq.stack_BK64.BK.HHN.--_.out", 
                          sep=" ",names=["freq", "med", "mean", "mode", "l95", "h95", "min", "max", "sd","var","sta","net","com","loc"],
                          header=None, skiprows=1)
bk64_hhn_psd_nona = bk64_hhn_psd.dropna()
bk64_hhe_psd =pd.read_csv("http://ncedc.org/ftp/outgoing/taira/WQC/ByerlyVault_PSD/"+data_set_name+"/freq.stack_BK64.BK.HHE.--_.out", 
                          sep=" ",names=["freq", "med", "mean", "mode", "l95", "h95", "min", "max", "sd","var","sta","net","com","loc"],
                          header=None, skiprows=1)
bk64_hhe_psd_nona = bk64_hhe_psd.dropna()
In [280]:
bk65_hhz_psd =pd.read_csv("http://ncedc.org/ftp/outgoing/taira/WQC/ByerlyVault_PSD/"+data_set_name+"/freq.stack_BK65.BK.HHZ.--_.out", 
                          sep=" ",names=["freq", "med", "mean", "mode", "l95", "h95", "min", "max", "sd","var","sta","net","com","loc"],
                          header=None, skiprows=1)
bk65_hhz_psd_nona = bk65_hhz_psd.dropna()
bk65_hhn_psd =pd.read_csv("http://ncedc.org/ftp/outgoing/taira/WQC/ByerlyVault_PSD/"+data_set_name+"/freq.stack_BK65.BK.HHN.--_.out", 
                          sep=" ",names=["freq", "med", "mean", "mode", "l95", "h95", "min", "max", "sd","var","sta","net","com","loc"],
                          header=None, skiprows=1)
bk65_hhn_psd_nona = bk65_hhn_psd.dropna()
bk65_hhe_psd =pd.read_csv("http://ncedc.org/ftp/outgoing/taira/WQC/ByerlyVault_PSD/"+data_set_name+"/freq.stack_BK65.BK.HHE.--_.out", 
                          sep=" ",names=["freq", "med", "mean", "mode", "l95", "h95", "min", "max", "sd","var","sta","net","com","loc"],
                          header=None, skiprows=1)
bk65_hhe_psd_nona = bk65_hhe_psd.dropna()
In [ ]:
In [281]:
ref_hhz_pasd_nona = bks_hhz_psd_nona.copy()
ref_hhn_pasd_nona = bks_hhn_psd_nona.copy()
ref_hhe_pasd_nona = bks_hhe_psd_nona.copy()
In [282]:
bk63_hhz_psd_nona['rmode'] = bk63_hhz_psd_nona['mode'] - ref_hhz_pasd_nona['mode']
bk63_hhn_psd_nona['rmode'] = bk63_hhn_psd_nona['mode'] - ref_hhn_pasd_nona['mode']
bk63_hhe_psd_nona['rmode'] = bk63_hhe_psd_nona['mode'] - ref_hhe_pasd_nona['mode']
In [283]:
bk64_hhz_psd_nona['rmode'] = bk64_hhz_psd_nona['mode'] - ref_hhz_pasd_nona['mode']
bk64_hhn_psd_nona['rmode'] = bk64_hhn_psd_nona['mode'] - ref_hhn_pasd_nona['mode']
bk64_hhe_psd_nona['rmode'] = bk64_hhe_psd_nona['mode'] - ref_hhe_pasd_nona['mode']
In [284]:
bk65_hhz_psd_nona['rmode'] = bk65_hhz_psd_nona['mode'] - ref_hhz_pasd_nona['mode']
bk65_hhn_psd_nona['rmode'] = bk65_hhn_psd_nona['mode'] - ref_hhn_pasd_nona['mode']
bk65_hhe_psd_nona['rmode'] = bk65_hhe_psd_nona['mode'] - ref_hhe_pasd_nona['mode']
In [285]:
bk69_hhz_psd_nona['rmode'] = bk69_hhz_psd_nona['mode'] - ref_hhz_pasd_nona['mode']
bk69_hhn_psd_nona['rmode'] = bk69_hhn_psd_nona['mode'] - ref_hhn_pasd_nona['mode']
bk69_hhe_psd_nona['rmode'] = bk69_hhe_psd_nona['mode'] - ref_hhe_pasd_nona['mode']
In [286]:
bk70_hhz_psd_nona['rmode'] = bk70_hhz_psd_nona['mode'] - ref_hhz_pasd_nona['mode']
bk70_hhn_psd_nona['rmode'] = bk70_hhn_psd_nona['mode'] - ref_hhn_pasd_nona['mode']
bk70_hhe_psd_nona['rmode'] = bk70_hhe_psd_nona['mode'] - ref_hhe_pasd_nona['mode']
In [287]:
bk71_hhz_psd_nona['rmode'] = bk71_hhz_psd_nona['mode'] - ref_hhz_pasd_nona['mode']
bk71_hhn_psd_nona['rmode'] = bk71_hhn_psd_nona['mode'] - ref_hhn_pasd_nona['mode']
bk71_hhe_psd_nona['rmode'] = bk71_hhe_psd_nona['mode'] - ref_hhe_pasd_nona['mode']
In [288]:
bk72_hhz_psd_nona['rmode'] = bk72_hhz_psd_nona['mode'] - ref_hhz_pasd_nona['mode']
bk72_hhn_psd_nona['rmode'] = bk72_hhn_psd_nona['mode'] - ref_hhn_pasd_nona['mode']
bk72_hhe_psd_nona['rmode'] = bk72_hhe_psd_nona['mode'] - ref_hhe_pasd_nona['mode']
In [289]:
bk80_hhz_psd_nona['rmode'] = bk80_hhz_psd_nona['mode'] - ref_hhz_pasd_nona['mode']
bk80_hhn_psd_nona['rmode'] = bk80_hhn_psd_nona['mode'] - ref_hhn_pasd_nona['mode']
bk80_hhe_psd_nona['rmode'] = bk80_hhe_psd_nona['mode'] - ref_hhe_pasd_nona['mode']
In [290]:
bk82_hhz_psd_nona['rmode'] = bk82_hhz_psd_nona['mode'] - ref_hhz_pasd_nona['mode']
bk82_hhn_psd_nona['rmode'] = bk82_hhn_psd_nona['mode'] - ref_hhn_pasd_nona['mode']
bk82_hhe_psd_nona['rmode'] = bk82_hhe_psd_nona['mode'] - ref_hhe_pasd_nona['mode']
In [ ]:
In [291]:
plt.rcParams['figure.figsize'] = 18, 8
plt.rcParams['figure.figsize'] = 18, 10

import matplotlib.pyplot as plt
import matplotlib.dates as mdates
fig, ax = plt.subplots()
ticks = [10,20,50,100,200,500,1000]
ax.set_xticks(ticks, minor=False)
cmap = plt.get_cmap("tab20c") 


plt.plot(1/bks_hhz_psd_nona.freq, bks_hhz_psd_nona['mode'], color=cmap(0), linestyle='solid',linewidth = 3.0, label='BKS STS-1')

plt.plot(1/bk63_hhz_psd_nona.freq, bk63_hhz_psd_nona['mode'], color=cmap(1), linestyle='solid',linewidth = 3.0, label='BK63 Trillium Compact')
plt.plot(1/bk64_hhz_psd_nona.freq, bk64_hhz_psd_nona['mode'], color=cmap(2), linestyle='solid',linewidth = 3.0, label='BK64 Trillium 120 PH')
plt.plot(1/bk65_hhz_psd_nona.freq, bk65_hhz_psd_nona['mode'], color=cmap(3), linestyle='solid',linewidth = 3.0, label='BK65 Trillium Horizon 120')
#plt.plot(1/bk66_hhz_psd_nona.freq, bk66_hhz_psd_nona['mode'], color=cmap(3), linestyle='solid',linewidth = 3.0, label='BK66')
#plt.plot(1/bk67_hhz_psd_nona.freq, bk67_hhz_psd_nona['mode'], color=cmap(4), linestyle='solid',linewidth = 3.0, label='BK67')
#plt.plot(1/bk68_hhz_psd_nona.freq, bk68_hhz_psd_nona['mode'], color=cmap(5), linestyle='solid',linewidth = 3.0, label='BK68')

plt.plot(1/bk69_hhz_psd_nona.freq, bk69_hhz_psd_nona['mode'], color=cmap(4), linestyle='solid',linewidth = 3.0, label='BK69 MBB-2')
#plt.plot(1/bk70_hhz_psd_nona.freq, bk70_hhz_psd_nona['mode'], color=cmap(5), linestyle='solid',linewidth = 3.0, label='BK70 MBB-2')
plt.plot(1/bk71_hhz_psd_nona.freq, bk71_hhz_psd_nona['mode'], color=cmap(6), linestyle='solid',linewidth = 3.0, label='BK71 MBB-2')
plt.plot(1/bk72_hhz_psd_nona.freq, bk72_hhz_psd_nona['mode'], color=cmap(7), linestyle='solid',linewidth = 3.0, label='BK72 MBB-2')

plt.plot(1/bk80_hhz_psd_nona.freq, bk80_hhz_psd_nona['mode'], color=cmap(8), linestyle='solid',linewidth = 3.0, label='BK80 CMG-3T (sand small beads)')
#plt.plot(1/bk81_hhz_psd_nona.freq, bk81_hhz_psd_nona['mode'], color=cmap(8), linestyle='solid',linewidth = 3.0, label='BK81')
plt.plot(1/bk82_hhz_psd_nona.freq, bk82_hhz_psd_nona['mode'], color=cmap(9), linestyle='solid',linewidth = 3.0, label='BK82 CMG-3T (sand 3/5mm mix)')
plt.plot(hnm_pdf[0], hnm_pdf[1], color='black',linestyle='dashed',linewidth = 2.0,label='HNM')
plt.plot(lnm_pdf[0], lnm_pdf[1], color='black',linestyle='solid',linewidth = 2.0,label='LNM')
plt.xscale("log")
plt.xlim(0.01,10000)
plt.xlim(0.01,100000)
plt.xlim(0.01,2000)

plt.ylim(-190,-80)
plt.ylim(-190,-40)

plt.grid(which='major',color='black',linestyle='-',linewidth = 0.5)
plt.grid(which='minor',color='black',linestyle='-',linewidth = 0.25)

plt.xlabel("Periods (s)", fontsize=16)
plt.ylabel("PSD Power [10log(m**2/sec**4/Hz)] (dB)", fontsize=16)
plt.title("Byerly Vault PSD Com=HHZ",fontsize=16)

plt.tick_params(labelsize=14)
#plt.xticks([0.01,0.02,0.05,0.1,0.2,0.5,1,2,5,10,20, 50,100,200,500,1000,2000,5000,10000], 
#           ["0.01","0.02","0.05","0.1","0.2","0.5","1","2","5","10","20","50","100","200","500","1000","2000","5000","10000"])
plt.xticks([0.01,0.02,0.05,0.1,0.2,0.5,1,2,5,10,20, 50,100,200,500,1000,2000], 
           ["0.01","0.02","0.05","0.1","0.2","0.5","1","2","5","10","20","50","100","200","500","1000","2000"])

#plt.legend([p3, p2, p1, p4, p5, p6, p7], ["BK.JEPS.LHN", "BK.JEPS.LHE", "High Noise Model", "TA-Mean (LHN,LHE)", "BDSN-RefTek-Mean (LHN,LHE)", "USGS-Crest-Mean (HHN,HHE)","BDSN-STS1-Mean (LHN,LHE)"], loc="upper right", fontsize=14)
#plt.legend([p33, p22, p3, p2, p1, p4, p5], ["BK.OAKV.LHN", "BK.OAKV.LHE", "BK.JEPS.LHN", "BK.JEPS.LHE", "High Noise Model", "TA-Mean (LHN,LHE)", "BDSN-RefTek-Mean (LHN,LHE)"], loc="upper right", fontsize=14)

plt.legend(loc="upper left", fontsize=14) # 凡例を表示

print("Figure 1: PSD of broadband data in the vertical component")

Figure 1: PSD of broadband data in the vertical component
In [ ]:
In [292]:
plt.rcParams['figure.figsize'] = 18, 10

import matplotlib.pyplot as plt
import matplotlib.dates as mdates
fig, ax = plt.subplots()
ticks = [10,20,50,100,200,500,1000]
ax.set_xticks(ticks, minor=False)
cmap = plt.get_cmap("tab20c") 


plt.plot(1/bks_hhn_psd_nona.freq, bks_hhn_psd_nona['mode'], color=cmap(0), linestyle='solid',linewidth = 3.0, label='BKS STS-1')

plt.plot(1/bk63_hhn_psd_nona.freq, bk63_hhn_psd_nona['mode'], color=cmap(1), linestyle='solid',linewidth = 3.0, label='BK63 Trillium Compact')
plt.plot(1/bk64_hhn_psd_nona.freq, bk64_hhn_psd_nona['mode'], color=cmap(2), linestyle='solid',linewidth = 3.0, label='BK64 Trillium 120 PH')
plt.plot(1/bk65_hhn_psd_nona.freq, bk65_hhn_psd_nona['mode'], color=cmap(3), linestyle='solid',linewidth = 3.0, label='BK65 Trillium Horizon 120')

plt.plot(1/bk69_hhn_psd_nona.freq, bk69_hhn_psd_nona['mode'], color=cmap(4), linestyle='solid',linewidth = 3.0, label='BK69 MBB-2')
plt.plot(1/bk71_hhn_psd_nona.freq, bk71_hhn_psd_nona['mode'], color=cmap(6), linestyle='solid',linewidth = 3.0, label='BK71 MBB-2')
plt.plot(1/bk72_hhn_psd_nona.freq, bk72_hhn_psd_nona['mode'], color=cmap(7), linestyle='solid',linewidth = 3.0, label='BK72 MBB-2')

plt.plot(1/bk80_hhn_psd_nona.freq, bk80_hhn_psd_nona['mode'], color=cmap(8), linestyle='solid',linewidth = 3.0, label='BK80 CMG-3T (sand small beads)')
plt.plot(1/bk82_hhn_psd_nona.freq, bk82_hhn_psd_nona['mode'], color=cmap(9), linestyle='solid',linewidth = 3.0, label='BK82 CMG-3T (sand 3/5mm mix)')



plt.plot(hnm_pdf[0], hnm_pdf[1], color='black',linestyle='dashed',linewidth = 2.0,label='HNM')
plt.plot(lnm_pdf[0], lnm_pdf[1], color='black',linestyle='solid',linewidth = 2.0,label='LNM')
plt.xscale("log")
plt.xlim(0.01,10000)
plt.xlim(0.01,100000)
plt.xlim(0.01,2000)

plt.ylim(-190,-80)
plt.ylim(-190,-40)

plt.grid(which='major',color='black',linestyle='-',linewidth = 0.5)
plt.grid(which='minor',color='black',linestyle='-',linewidth = 0.25)

plt.xlabel("Periods (s)", fontsize=16)
plt.ylabel("PSD Power [10log(m**2/sec**4/Hz)] (dB)", fontsize=16)
plt.title("Byerly Vault PSD Com=HHN",fontsize=16)

plt.tick_params(labelsize=14)
#plt.xticks([0.01,0.02,0.05,0.1,0.2,0.5,1,2,5,10,20, 50,100,200,500,1000,2000,5000,10000], 
#           ["0.01","0.02","0.05","0.1","0.2","0.5","1","2","5","10","20","50","100","200","500","1000","2000","5000","10000"])
plt.xticks([0.01,0.02,0.05,0.1,0.2,0.5,1,2,5,10,20, 50,100,200,500,1000,2000], 
           ["0.01","0.02","0.05","0.1","0.2","0.5","1","2","5","10","20","50","100","200","500","1000","2000"])

#plt.legend([p3, p2, p1, p4, p5, p6, p7], ["BK.JEPS.LHN", "BK.JEPS.LHE", "High Noise Model", "TA-Mean (LHN,LHE)", "BDSN-RefTek-Mean (LHN,LHE)", "USGS-Crest-Mean (HHN,HHE)","BDSN-STS1-Mean (LHN,LHE)"], loc="upper right", fontsize=14)
#plt.legend([p33, p22, p3, p2, p1, p4, p5], ["BK.OAKV.LHN", "BK.OAKV.LHE", "BK.JEPS.LHN", "BK.JEPS.LHE", "High Noise Model", "TA-Mean (LHN,LHE)", "BDSN-RefTek-Mean (LHN,LHE)"], loc="upper right", fontsize=14)

plt.legend(loc="upper left", fontsize=14) # 凡例を表示

print("Figure 2: PSD of broadband data in the N-S component")

Figure 2: PSD of broadband data in the N-S component
In [ ]:
In [293]:
plt.rcParams['figure.figsize'] = 18, 10
#plt.rcParams['figure.figsize'] = 14, 6

import matplotlib.pyplot as plt
import matplotlib.dates as mdates
fig, ax = plt.subplots()
ticks = [10,20,50,100,200,500,1000]
ax.set_xticks(ticks, minor=False)
cmap = plt.get_cmap("tab20c") 



plt.plot(1/bks_hhe_psd_nona.freq, bks_hhe_psd_nona['mode'], color=cmap(0), linestyle='solid',linewidth = 3.0, label='BKS STS-1')

plt.plot(1/bk63_hhe_psd_nona.freq, bk63_hhe_psd_nona['mode'], color=cmap(1), linestyle='solid',linewidth = 3.0, label='BK63 Trillium Compact')
plt.plot(1/bk64_hhe_psd_nona.freq, bk64_hhe_psd_nona['mode'], color=cmap(2), linestyle='solid',linewidth = 3.0, label='BK64 Trillium 120 PH')
plt.plot(1/bk65_hhe_psd_nona.freq, bk65_hhe_psd_nona['mode'], color=cmap(3), linestyle='solid',linewidth = 3.0, label='BK65 Trillium Horizon 120')

plt.plot(1/bk69_hhe_psd_nona.freq, bk69_hhe_psd_nona['mode'], color=cmap(4), linestyle='solid',linewidth = 3.0, label='BK69 MBB-2')
plt.plot(1/bk71_hhe_psd_nona.freq, bk71_hhe_psd_nona['mode'], color=cmap(6), linestyle='solid',linewidth = 3.0, label='BK71 MBB-2')
plt.plot(1/bk72_hhe_psd_nona.freq, bk72_hhe_psd_nona['mode'], color=cmap(7), linestyle='solid',linewidth = 3.0, label='BK72 MBB-2')

plt.plot(1/bk80_hhe_psd_nona.freq, bk80_hhe_psd_nona['mode'], color=cmap(8), linestyle='solid',linewidth = 3.0, label='BK80 CMG-3T (sand small beads)')
plt.plot(1/bk82_hhe_psd_nona.freq, bk82_hhe_psd_nona['mode'], color=cmap(9), linestyle='solid',linewidth = 3.0, label='BK82 CMG-3T (sand 3/5mm mix)')




plt.plot(hnm_pdf[0], hnm_pdf[1], color='black',linestyle='dashed',linewidth = 2.0,label='HNM')
plt.plot(lnm_pdf[0], lnm_pdf[1], color='black',linestyle='solid',linewidth = 2.0,label='LNM')
plt.xscale("log")

#plt.xlim(0.01,10000)
#plt.xlim(0.01,100000)
plt.xlim(0.01,2000)

plt.ylim(-190,-80)
plt.ylim(-190,-40)

plt.grid(which='major',color='black',linestyle='-',linewidth = 0.5)
plt.grid(which='minor',color='black',linestyle='-',linewidth = 0.25)

plt.xlabel("Periods (s)", fontsize=16)
plt.ylabel("PSD Power [10log(m**2/sec**4/Hz)] (dB)", fontsize=16)
plt.title("Byerly Vault PSD Com=HHE",fontsize=16)

plt.tick_params(labelsize=14)
#plt.xticks([0.01,0.02,0.05,0.1,0.2,0.5,1,2,5,10,20, 50,100,200,500,1000,2000,5000,10000], 
#           ["0.01","0.02","0.05","0.1","0.2","0.5","1","2","5","10","20","50","100","200","500","1000","2000","5000","10000"])
plt.xticks([0.01,0.02,0.05,0.1,0.2,0.5,1,2,5,10,20, 50,100,200,500,1000,2000], 
           ["0.01","0.02","0.05","0.1","0.2","0.5","1","2","5","10","20","50","100","200","500","1000","2000"])

#plt.legend([p3, p2, p1, p4, p5, p6, p7], ["BK.JEPS.LHN", "BK.JEPS.LHE", "High Noise Model", "TA-Mean (LHN,LHE)", "BDSN-RefTek-Mean (LHN,LHE)", "USGS-Crest-Mean (hhe,HHE)","BDSN-STS1-Mean (LHN,LHE)"], loc="upper right", fontsize=14)
#plt.legend([p33, p22, p3, p2, p1, p4, p5], ["BK.OAKV.LHN", "BK.OAKV.LHE", "BK.JEPS.LHN", "BK.JEPS.LHE", "High Noise Model", "TA-Mean (LHN,LHE)", "BDSN-RefTek-Mean (LHN,LHE)"], loc="upper right", fontsize=14)

plt.legend(loc="upper left", fontsize=14) # 凡例を表示

print("Figure 3: PSD of broadband data in the E-W component")

Figure 3: PSD of broadband data in the E-W component
In [ ]:
In [294]:
plt.rcParams['figure.figsize'] = 14, 6
plt.rcParams['figure.figsize'] = 18, 10

import matplotlib.pyplot as plt
import matplotlib.dates as mdates
fig, ax = plt.subplots()
ticks = [10,20,50,100,200,500,1000]
ax.set_xticks(ticks, minor=False)
cmap = plt.get_cmap("tab20c") 


#plt.plot(1/bks_hhz_psd_nona.freq, bks_hhz_psd_nona.rmed, color=cmap(0), linestyle='solid',linewidth = 3.0, label='BKS')
plt.plot(1/bk63_hhz_psd_nona.freq, bk63_hhz_psd_nona['rmode'], color=cmap(1), linestyle='solid',linewidth = 3.0, label='BK63 Trillium Compact')
plt.plot(1/bk64_hhz_psd_nona.freq, bk64_hhz_psd_nona['rmode'], color=cmap(2), linestyle='solid',linewidth = 3.0, label='BK64 Trillium 120 PH')
plt.plot(1/bk65_hhz_psd_nona.freq, bk65_hhz_psd_nona['rmode'], color=cmap(3), linestyle='solid',linewidth = 3.0, label='BK65 Trillium Horizon 120')



plt.plot(1/bk69_hhz_psd_nona.freq, bk69_hhz_psd_nona['rmode'], color=cmap(4), linestyle='solid',linewidth = 3.0, label='BK69 MBB-2')
plt.plot(1/bk71_hhz_psd_nona.freq, bk71_hhz_psd_nona['rmode'], color=cmap(6), linestyle='solid',linewidth = 3.0, label='BK71 MBB-2')
plt.plot(1/bk72_hhz_psd_nona.freq, bk72_hhz_psd_nona['rmode'], color=cmap(7), linestyle='solid',linewidth = 3.0, label='BK72 MBB-2')

plt.plot(1/bk80_hhz_psd_nona.freq, bk80_hhz_psd_nona['rmode'], color=cmap(8), linestyle='solid',linewidth = 3.0, label='BK80 CMG-3T (sand small beads)')
plt.plot(1/bk82_hhz_psd_nona.freq, bk82_hhz_psd_nona['rmode'], color=cmap(9), linestyle='solid',linewidth = 3.0, label='BK82 CMG-3T (sand 3/5mm mix)')


#plt.plot(hnm_pdf[0], hnm_pdf[1], color='black',linestyle='dashed',linewidth = 2.0,label='HNM')
#plt.plot(lnm_pdf[0], lnm_pdf[1], color='black',linestyle='solid',linewidth = 2.0,label='LNM')
plt.xscale("log")
plt.xlim(0.01,10000)
plt.xlim(0.01,2000)
plt.ylim(-20,80)

plt.grid(which='major',color='black',linestyle='-',linewidth = 0.5)
plt.grid(which='minor',color='black',linestyle='-',linewidth = 0.25)

plt.xlabel("Periods (s)", fontsize=16)
plt.ylabel("Relative PSD Power [10log(m**2/sec**4/Hz)] (dB)", fontsize=16)
plt.title("Byerly Vault rPSD (ref:BKS) Com=HHZ",fontsize=16)

plt.tick_params(labelsize=14)
#plt.xticks([0.01,0.02,0.05,0.1,0.2,0.5,1,2,5,10,20, 50,100,200,500,1000,2000,5000,10000], 
#           ["0.01","0.02","0.05","0.1","0.2","0.5","1","2","5","10","20","50","100","200","500","1000","2000","5000","10000"])
plt.xticks([0.01,0.02,0.05,0.1,0.2,0.5,1,2,5,10,20, 50,100,200,500,1000,2000], 
           ["0.01","0.02","0.05","0.1","0.2","0.5","1","2","5","10","20","50","100","200","500","1000","2000"])

#plt.legend([p3, p2, p1, p4, p5, p6, p7], ["BK.JEPS.LHN", "BK.JEPS.LHE", "High Noise Model", "TA-Mean (LHN,LHE)", "BDSN-RefTek-Mean (LHN,LHE)", "USGS-Crest-Mean (HHN,HHE)","BDSN-STS1-Mean (LHN,LHE)"], loc="upper right", fontsize=14)
#plt.legend([p33, p22, p3, p2, p1, p4, p5], ["BK.OAKV.LHN", "BK.OAKV.LHE", "BK.JEPS.LHN", "BK.JEPS.LHE", "High Noise Model", "TA-Mean (LHN,LHE)", "BDSN-RefTek-Mean (LHN,LHE)"], loc="upper right", fontsize=14)

plt.legend(loc="upper left", fontsize=14) # 凡例を表示

print("Figure 4: Relative PSD (BKS as reference) of broadband data in the vertical component")

Figure 4: Relative PSD (BKS as reference) of broadband data in the vertical component
In [ ]:
In [ ]:
In [295]:
plt.rcParams['figure.figsize'] = 14, 6
plt.rcParams['figure.figsize'] = 18, 10

import matplotlib.pyplot as plt
import matplotlib.dates as mdates
fig, ax = plt.subplots()
ticks = [10,20,50,100,200,500,1000]
ax.set_xticks(ticks, minor=False)
cmap = plt.get_cmap("tab20c") 


#plt.plot(1/bks_hhn_psd_nona.freq, bks_hhn_psd_nona.rmed, color=cmap(0), linestyle='solid',linewidth = 3.0, label='BKS')
plt.plot(1/bk63_hhn_psd_nona.freq, bk63_hhn_psd_nona['rmode'], color=cmap(1), linestyle='solid',linewidth = 3.0, label='BK63 Trillium Compact')
plt.plot(1/bk64_hhn_psd_nona.freq, bk64_hhn_psd_nona['rmode'], color=cmap(2), linestyle='solid',linewidth = 3.0, label='BK64 Trillium 120 PH')
plt.plot(1/bk65_hhn_psd_nona.freq, bk65_hhn_psd_nona['rmode'], color=cmap(3), linestyle='solid',linewidth = 3.0, label='BK65 Trillium Horizon 120')



plt.plot(1/bk69_hhn_psd_nona.freq, bk69_hhn_psd_nona['rmode'], color=cmap(4), linestyle='solid',linewidth = 3.0, label='BK69 MBB-2')
plt.plot(1/bk71_hhn_psd_nona.freq, bk71_hhn_psd_nona['rmode'], color=cmap(6), linestyle='solid',linewidth = 3.0, label='BK71 MBB-2')
plt.plot(1/bk72_hhn_psd_nona.freq, bk72_hhn_psd_nona['rmode'], color=cmap(7), linestyle='solid',linewidth = 3.0, label='BK72 MBB-2')

plt.plot(1/bk80_hhn_psd_nona.freq, bk80_hhn_psd_nona['rmode'], color=cmap(8), linestyle='solid',linewidth = 3.0, label='BK80 CMG-3T (sand small beads)')
plt.plot(1/bk82_hhn_psd_nona.freq, bk82_hhn_psd_nona['rmode'], color=cmap(9), linestyle='solid',linewidth = 3.0, label='BK82 CMG-3T (sand 3/5mm mix)')


#plt.plot(hnm_pdf[0], hnm_pdf[1], color='black',linestyle='dashed',linewidth = 2.0,label='HNM')
#plt.plot(lnm_pdf[0], lnm_pdf[1], color='black',linestyle='solid',linewidth = 2.0,label='LNM')
plt.xscale("log")
plt.xlim(0.01,10000)
plt.xlim(0.01,2000)
plt.ylim(-10,50)
plt.ylim(-20,80)

plt.grid(which='major',color='black',linestyle='-',linewidth = 0.5)
plt.grid(which='minor',color='black',linestyle='-',linewidth = 0.25)

plt.xlabel("Periods (s)", fontsize=16)
plt.ylabel("Relative PSD Power [10log(m**2/sec**4/Hz)] (dB)", fontsize=16)
plt.title("Byerly Vault rPSD (ref:BKS) Com=HHN",fontsize=16)

plt.tick_params(labelsize=14)
#plt.xticks([0.01,0.02,0.05,0.1,0.2,0.5,1,2,5,10,20, 50,100,200,500,1000,2000,5000,10000], 
#           ["0.01","0.02","0.05","0.1","0.2","0.5","1","2","5","10","20","50","100","200","500","1000","2000","5000","10000"])
plt.xticks([0.01,0.02,0.05,0.1,0.2,0.5,1,2,5,10,20, 50,100,200,500,1000,2000], 
           ["0.01","0.02","0.05","0.1","0.2","0.5","1","2","5","10","20","50","100","200","500","1000","2000"])

#plt.legend([p3, p2, p1, p4, p5, p6, p7], ["BK.JEPS.LHN", "BK.JEPS.LHE", "High Noise Model", "TA-Mean (LHN,LHE)", "BDSN-RefTek-Mean (LHN,LHE)", "USGS-Crest-Mean (HHN,HHE)","BDSN-STS1-Mean (LHN,LHE)"], loc="upper right", fontsize=14)
#plt.legend([p33, p22, p3, p2, p1, p4, p5], ["BK.OAKV.LHN", "BK.OAKV.LHE", "BK.JEPS.LHN", "BK.JEPS.LHE", "High Noise Model", "TA-Mean (LHN,LHE)", "BDSN-RefTek-Mean (LHN,LHE)"], loc="upper right", fontsize=14)

plt.legend(loc="upper left", fontsize=14) # 凡例を表示

print("Figure 5: Relative PSD (BKS as reference) of broadband data in the N-S component")

Figure 5: Relative PSD (BKS as reference) of broadband data in the N-S component
In [ ]:
In [296]:
plt.rcParams['figure.figsize'] = 14, 6
plt.rcParams['figure.figsize'] = 18, 10

import matplotlib.pyplot as plt
import matplotlib.dates as mdates
fig, ax = plt.subplots()
ticks = [10,20,50,100,200,500,1000]
ax.set_xticks(ticks, minor=False)
cmap = plt.get_cmap("tab20c") 


#plt.plot(1/bks_hhe_psd_nona.freq, bks_hhe_psd_nona.rmed, color=cmap(0), linestyle='solid',linewidth = 3.0, label='BKS')
plt.plot(1/bk63_hhe_psd_nona.freq, bk63_hhe_psd_nona['rmode'], color=cmap(1), linestyle='solid',linewidth = 3.0, label='BK63 Trillium Compact')
plt.plot(1/bk64_hhe_psd_nona.freq, bk64_hhe_psd_nona['rmode'], color=cmap(2), linestyle='solid',linewidth = 3.0, label='BK64 Trillium 120 PH')
plt.plot(1/bk65_hhe_psd_nona.freq, bk65_hhe_psd_nona['rmode'], color=cmap(3), linestyle='solid',linewidth = 3.0, label='BK65 Trillium Horizon 120')



plt.plot(1/bk69_hhe_psd_nona.freq, bk69_hhe_psd_nona['rmode'], color=cmap(4), linestyle='solid',linewidth = 3.0, label='BK69 MBB-2')
plt.plot(1/bk71_hhe_psd_nona.freq, bk71_hhe_psd_nona['rmode'], color=cmap(6), linestyle='solid',linewidth = 3.0, label='BK71 MBB-2')
plt.plot(1/bk72_hhe_psd_nona.freq, bk72_hhe_psd_nona['rmode'], color=cmap(7), linestyle='solid',linewidth = 3.0, label='BK72 MBB-2')

plt.plot(1/bk80_hhe_psd_nona.freq, bk80_hhe_psd_nona['rmode'], color=cmap(8), linestyle='solid',linewidth = 3.0, label='BK80 CMG-3T (sand small beads)')
plt.plot(1/bk82_hhe_psd_nona.freq, bk82_hhe_psd_nona['rmode'], color=cmap(9), linestyle='solid',linewidth = 3.0, label='BK82 CMG-3T (sand 3/5mm mix)')


#plt.plot(hnm_pdf[0], hnm_pdf[1], color='black',linestyle='dashed',linewidth = 2.0,label='HNM')
#plt.plot(lnm_pdf[0], lnm_pdf[1], color='black',linestyle='solid',linewidth = 2.0,label='LNM')
plt.xscale("log")
plt.xlim(0.01,10000)
plt.xlim(0.01,2000)
plt.ylim(-10,50)
plt.ylim(-20,80)
plt.grid(which='major',color='black',linestyle='-',linewidth = 0.5)
plt.grid(which='minor',color='black',linestyle='-',linewidth = 0.25)

plt.xlabel("Periods (s)", fontsize=16)
plt.ylabel("Relative PSD Power [10log(m**2/sec**4/Hz)] (dB)", fontsize=16)
plt.title("Byerly Vault rPSD (ref:BKS) Com=HHE",fontsize=16)

plt.tick_params(labelsize=14)
#plt.xticks([0.01,0.02,0.05,0.1,0.2,0.5,1,2,5,10,20, 50,100,200,500,1000,2000,5000,10000], 
#           ["0.01","0.02","0.05","0.1","0.2","0.5","1","2","5","10","20","50","100","200","500","1000","2000","5000","10000"])
plt.xticks([0.01,0.02,0.05,0.1,0.2,0.5,1,2,5,10,20, 50,100,200,500,1000,2000], 
           ["0.01","0.02","0.05","0.1","0.2","0.5","1","2","5","10","20","50","100","200","500","1000","2000"])

#plt.legend([p3, p2, p1, p4, p5, p6, p7], ["BK.JEPS.LHN", "BK.JEPS.LHE", "High Noise Model", "TA-Mean (LHN,LHE)", "BDSN-RefTek-Mean (LHN,LHE)", "USGS-Crest-Mean (hhe,HHE)","BDSN-STS1-Mean (LHN,LHE)"], loc="upper right", fontsize=14)
#plt.legend([p33, p22, p3, p2, p1, p4, p5], ["BK.OAKV.LHN", "BK.OAKV.LHE", "BK.JEPS.LHN", "BK.JEPS.LHE", "High Noise Model", "TA-Mean (LHN,LHE)", "BDSN-RefTek-Mean (LHN,LHE)"], loc="upper right", fontsize=14)

plt.legend(loc="upper left", fontsize=14) # 凡例を表示

print("Figure 6: Relative PSD (BKS as reference) of broadband data in the E-W component")

Figure 6: Relative PSD (BKS as reference) of broadband data in the E-W component
In [ ]:
In [297]:
#plt.rcParams['figure.figsize'] = 11, 6
plt.rcParams['figure.figsize'] = 15, 5
plt.xlabel('Station',fontsize=16)

sta_list = bb_100s.sta.as_matrix()
zval_med = bb_100s.zval_med.as_matrix()
nval_med = bb_100s.nval_med.as_matrix()
eval_med = bb_100s.eval_med.as_matrix()

x_position = np.arange(len(sta_list))
plt.bar(x_position - 0.3, zval_med, width=0.3, label='HHZ')
plt.bar(x_position, nval_med, width=0.3, label='HHN')                    
plt.bar(x_position + 0.3, eval_med, width=0.3, label='HHE')

plt.legend(fontsize=16)
#plt.xticks(x_position + 0.2, x)
plt.xticks(x_position + 0.1, sta_list)

plt.ylabel('PSD at 100-s period', fontsize=16)
#plt.ylabel('PSD at 200-s period', fontsize=16)
plt.title("Byerly Vault PSDs at 100 s period",fontsize=16)

# HN
#plt.ylim(-85, -120)
# HH
plt.ylim(-120, -180)
plt.ylim(-120, -190)
plt.legend(loc="upper right", fontsize=14) # 凡例を表示
plt.tick_params(labelsize=14)


print("Figure 7: PSD of broadband data at 100-s period")

plt.show()

Figure 7: PSD of broadband data at 100-s period
In [ ]:
In [ ]:
In [298]:
#plt.rcParams['figure.figsize'] = 11, 6
plt.rcParams['figure.figsize'] = 15, 5
plt.xlabel('Station',fontsize=16)

sta_list = bb_50s.sta.as_matrix()
zval_med = bb_50s.zval_med.as_matrix()
nval_med = bb_50s.nval_med.as_matrix()
eval_med = bb_50s.eval_med.as_matrix()

x_position = np.arange(len(sta_list))
plt.bar(x_position - 0.3, zval_med, width=0.3, label='HHZ')
plt.bar(x_position, nval_med, width=0.3, label='HHN')                    
plt.bar(x_position + 0.3, eval_med, width=0.3, label='HHE')

plt.legend(fontsize=16)
#plt.xticks(x_position + 0.2, x)
plt.xticks(x_position + 0.1, sta_list)

plt.ylabel('PSD at 50-s period', fontsize=16)
#plt.ylabel('PSD at 200-s period', fontsize=16)
plt.title("Byerly Vault PSDs at 50 s period",fontsize=16)

# HN
#plt.ylim(-85, -120)
# HH
plt.ylim(-120, -190)
plt.legend(loc="upper right", fontsize=14) # 凡例を表示
plt.tick_params(labelsize=14)


print("Figure 8: PSD of broadband data at 50-s period")


plt.show()

Figure 8: PSD of broadband data at 50-s period
In [299]:
#plt.rcParams['figure.figsize'] = 11, 6
plt.rcParams['figure.figsize'] = 15, 5
plt.xlabel('Station',fontsize=16)

sta_list = bb_1s.sta.as_matrix()
zval_med = bb_1s.zval_med.as_matrix()
nval_med = bb_1s.nval_med.as_matrix()
eval_med = bb_1s.eval_med.as_matrix()

x_position = np.arange(len(sta_list))
plt.bar(x_position - 0.3, zval_med, width=0.3, label='HHZ')
plt.bar(x_position, nval_med, width=0.3, label='HHN')                    
plt.bar(x_position + 0.3, eval_med, width=0.3, label='HHE')

plt.legend(fontsize=16)
#plt.xticks(x_position + 0.2, x)
plt.xticks(x_position + 0.1, sta_list)

plt.ylabel('PSD at 1-s period', fontsize=16)
#plt.ylabel('PSD at 200-s period', fontsize=16)
plt.title("Byerly Vault PSDs at 1 s period",fontsize=16)

# HN
#plt.ylim(-85, -120)
# HH
plt.ylim(-120, -140)
plt.legend(loc="upper right", fontsize=14) # 凡例を表示
plt.tick_params(labelsize=14)


print("Figure 9: PSD of broadband data at 1-s period")


plt.show()

Figure 9: PSD of broadband data at 1-s period
In [ ]:
In [ ]: