New function invoked by autoPyLoT for automated picking of onset times. Main tool for automatic picking!
This commit is contained in:
Ludger Küperkoch
parent
5be662524f
commit
74682952e7
401
pylot/core/pick/run_autopicking.py
Executable file
401
pylot/core/pick/run_autopicking.py
Executable file
@ -0,0 +1,401 @@
|
|||||||
|
#!/usr/bin/env python
|
||||||
|
# -*- coding: utf-8 -*-
|
||||||
|
|
||||||
|
"""
|
||||||
|
Function to run automated picking algorithms using AIC,
|
||||||
|
HOS and AR prediction. Uses object CharFuns and Picker and
|
||||||
|
function conglomerate utils.
|
||||||
|
|
||||||
|
:author: MAGS2 EP3 working group / Ludger Kueperkoch
|
||||||
|
"""
|
||||||
|
|
||||||
|
from obspy.core import read
|
||||||
|
import matplotlib.pyplot as plt
|
||||||
|
import numpy as np
|
||||||
|
from pylot.core.pick.CharFuns import *
|
||||||
|
from pylot.core.pick.Picker import *
|
||||||
|
from pylot.core.pick.CharFuns import *
|
||||||
|
from pylot.core.pick import utils
|
||||||
|
|
||||||
|
import pdb
|
||||||
|
|
||||||
|
def run_autopicking(wfstream, pickparam):
|
||||||
|
|
||||||
|
'''
|
||||||
|
param: wfstream
|
||||||
|
:type: `~obspy.core.stream.Stream`
|
||||||
|
|
||||||
|
param: pickparam
|
||||||
|
:type: container of picking parameters from input file,
|
||||||
|
usually autoPyLoT.in
|
||||||
|
'''
|
||||||
|
|
||||||
|
# declaring pickparam variables (only for convenience)
|
||||||
|
# read your autoPyLoT.in for details!
|
||||||
|
|
||||||
|
#special parameters for P picking
|
||||||
|
algoP = pickparam.getParam('algoP')
|
||||||
|
iplot = pickparam.getParam('iplot')
|
||||||
|
pstart = pickparam.getParam('pstart')
|
||||||
|
pstop = pickparam.getParam('pstop')
|
||||||
|
thosmw = pickparam.getParam('tlta')
|
||||||
|
hosorder = pickparam.getParam('hosorder')
|
||||||
|
tsnrz = pickparam.getParam('tsnrz')
|
||||||
|
hosorder = pickparam.getParam('hosorder')
|
||||||
|
bpz1 = pickparam.getParam('bpz1')
|
||||||
|
bpz2 = pickparam.getParam('bpz2')
|
||||||
|
pickwinP = pickparam.getParam('pickwinP')
|
||||||
|
tsmoothP = pickparam.getParam('tsmoothP')
|
||||||
|
ausP = pickparam.getParam('ausP')
|
||||||
|
nfacP = pickparam.getParam('nfacP')
|
||||||
|
tpred1z = pickparam.getParam('tpred1z')
|
||||||
|
tdet1z = pickparam.getParam('tdet1z')
|
||||||
|
Parorder = pickparam.getParam('Parorder')
|
||||||
|
addnoise = pickparam.getParam('addnoise')
|
||||||
|
Precalcwin = pickparam.getParam('Precalcwin')
|
||||||
|
minAICPslope = pickparam.getParam('minAICPslope')
|
||||||
|
minAICPSNR = pickparam.getParam('minAICPSNR')
|
||||||
|
timeerrorsP = pickparam.getParam('timeerrorsP')
|
||||||
|
#special parameters for S picking
|
||||||
|
algoS = pickparam.getParam('algoS')
|
||||||
|
sstart = pickparam.getParam('sstart')
|
||||||
|
sstop = pickparam.getParam('sstop')
|
||||||
|
bph1 = pickparam.getParam('bph1')
|
||||||
|
bph2 = pickparam.getParam('bph2')
|
||||||
|
tsnrh = pickparam.getParam('tsnrh')
|
||||||
|
pickwinS = pickparam.getParam('pickwinS')
|
||||||
|
tpred1h = pickparam.getParam('tpred1h')
|
||||||
|
tdet1h = pickparam.getParam('tdet1h')
|
||||||
|
tpred2h = pickparam.getParam('tpred2h')
|
||||||
|
tdet2h = pickparam.getParam('tdet2h')
|
||||||
|
Sarorder = pickparam.getParam('Sarorder')
|
||||||
|
aictsmoothS = pickparam.getParam('aictsmoothS')
|
||||||
|
tsmoothS = pickparam.getParam('tsmoothS')
|
||||||
|
ausS = pickparam.getParam('ausS')
|
||||||
|
minAICSslope = pickparam.getParam('minAICSslope')
|
||||||
|
minAICSSNR = pickparam.getParam('minAICSSNR')
|
||||||
|
Srecalcwin = pickparam.getParam('Srecalcwin')
|
||||||
|
nfacS = pickparam.getParam('nfacS')
|
||||||
|
timeerrorsS = pickparam.getParam('timeerrorsS')
|
||||||
|
|
||||||
|
|
||||||
|
# split components
|
||||||
|
zdat = wfstream.select(component="Z")
|
||||||
|
edat = wfstream.select(component="E")
|
||||||
|
ndat = wfstream.select(component="N")
|
||||||
|
|
||||||
|
if algoP == 'HOS' or algoP == 'ARZ' and zdat is not None:
|
||||||
|
print '##########################################'
|
||||||
|
print 'run_autopicking: Working on P onset of station %s' % zdat[0].stats.station
|
||||||
|
print 'Filtering vertical trace ...'
|
||||||
|
print zdat
|
||||||
|
z_copy = zdat.copy()
|
||||||
|
#filter and taper data
|
||||||
|
tr_filt = zdat[0].copy()
|
||||||
|
tr_filt.filter('bandpass', freqmin=bpz1[0], freqmax=bpz1[1], zerophase=False)
|
||||||
|
tr_filt.taper(max_percentage=0.05, type='hann')
|
||||||
|
z_copy[0].data = tr_filt.data
|
||||||
|
##############################################################
|
||||||
|
#check length of waveform and compare with cut times
|
||||||
|
Lc = pstop - pstart
|
||||||
|
Lwf = zdat[0].stats.endtime - zdat[0].stats.starttime
|
||||||
|
Ldiff = Lwf - Lc
|
||||||
|
if Ldiff < 0:
|
||||||
|
print 'run_autopicking: Cutting times are too large for actual waveform!'
|
||||||
|
print 'Use entire waveform instead!'
|
||||||
|
pstart = 0
|
||||||
|
pstop = len(zdat[0].data) * zdat[0].stats.delta
|
||||||
|
cuttimes = [pstart, pstop]
|
||||||
|
if algoP == 'HOS':
|
||||||
|
#calculate HOS-CF using subclass HOScf of class CharacteristicFunction
|
||||||
|
cf1 = HOScf(z_copy, cuttimes, thosmw, hosorder) #instance of HOScf
|
||||||
|
elif algoP == 'ARZ':
|
||||||
|
#calculate ARZ-CF using subclass ARZcf of class CharcteristicFunction
|
||||||
|
cf1 = ARZcf(z_copy, cuttimes, tpred1z, Parorder, tdet1z, addnoise) #instance of ARZcf
|
||||||
|
##############################################################
|
||||||
|
#calculate AIC-HOS-CF using subclass AICcf of class CharacteristicFunction
|
||||||
|
#class needs stream object => build it
|
||||||
|
tr_aic = tr_filt.copy()
|
||||||
|
tr_aic.data =cf1.getCF()
|
||||||
|
z_copy[0].data = tr_aic.data
|
||||||
|
aiccf = AICcf(z_copy, cuttimes) #instance of AICcf
|
||||||
|
##############################################################
|
||||||
|
#get prelimenary onset time from AIC-HOS-CF using subclass AICPicker of class AutoPicking
|
||||||
|
aicpick = AICPicker(aiccf, tsnrz, pickwinP, iplot, None, tsmoothP)
|
||||||
|
##############################################################
|
||||||
|
#go on with processing if AIC onset passes quality control
|
||||||
|
if aicpick.getSlope() >= minAICPslope and aicpick.getSNR() >= minAICPSNR:
|
||||||
|
print 'AIC P-pick passes quality control: Slope: %f, SNR: %f' % \
|
||||||
|
(aicpick.getSlope(), aicpick.getSNR())
|
||||||
|
print 'Go on with refined picking ...'
|
||||||
|
#re-filter waveform with larger bandpass
|
||||||
|
print 'run_autopicking: re-filtering vertical trace ...'
|
||||||
|
z_copy = zdat.copy()
|
||||||
|
tr_filt = zdat[0].copy()
|
||||||
|
tr_filt.filter('bandpass', freqmin=bpz2[0], freqmax=bpz2[1], zerophase=False)
|
||||||
|
tr_filt.taper(max_percentage=0.05, type='hann')
|
||||||
|
z_copy[0].data = tr_filt.data
|
||||||
|
#############################################################
|
||||||
|
#re-calculate CF from re-filtered trace in vicinity of initial onset
|
||||||
|
cuttimes2 = [round(max([aicpick.getpick() - Precalcwin, 0])), \
|
||||||
|
round(min([len(zdat[0].data) * zdat[0].stats.delta, \
|
||||||
|
aicpick.getpick() + Precalcwin]))]
|
||||||
|
if algoP == 'HOS':
|
||||||
|
#calculate HOS-CF using subclass HOScf of class CharacteristicFunction
|
||||||
|
cf2 = HOScf(z_copy, cuttimes2, thosmw, hosorder) #instance of HOScf
|
||||||
|
elif algoP == 'ARZ':
|
||||||
|
#calculate ARZ-CF using subclass ARZcf of class CharcteristicFunction
|
||||||
|
cf2 = ARZcf(z_copy, cuttimes2, tpred1z, Parorder, tdet1z, addnoise) #instance of ARZcf
|
||||||
|
##############################################################
|
||||||
|
#get refined onset time from CF2 using class Picker
|
||||||
|
refPpick = PragPicker(cf2, tsnrz, pickwinP, iplot, ausP, tsmoothP, aicpick.getpick())
|
||||||
|
#############################################################
|
||||||
|
#quality assessment
|
||||||
|
#get earliest and latest possible pick and symmetrized uncertainty
|
||||||
|
[lpickP, epickP, Perror] = earllatepicker(z_copy, nfacP, tsnrz, refPpick.getpick(), iplot)
|
||||||
|
|
||||||
|
#get SNR
|
||||||
|
[SNRP, SNRPdB, Pnoiselevel] = getSNR(z_copy, tsnrz, refPpick.getpick())
|
||||||
|
|
||||||
|
#weight P-onset using symmetric error
|
||||||
|
if Perror <= timeerrorsP[0]:
|
||||||
|
Pweight = 0
|
||||||
|
elif Perror > timeerrorsP[0] and Perror <= timeerrorsP[1]:
|
||||||
|
Pweight = 1
|
||||||
|
elif Perror > timeerrorsP[1] and Perror <= timeerrorsP[2]:
|
||||||
|
Pweight = 2
|
||||||
|
elif Perror > timeerrorsP[2] and Perror <= timeerrorsP[3]:
|
||||||
|
Pweight = 3
|
||||||
|
elif Perror > timeerrorsP[3]:
|
||||||
|
Pweight = 4
|
||||||
|
|
||||||
|
print 'run_autopicking: P-weight: %d, SNR: %f, SNR[dB]: %f' % (Pweight, SNRP, SNRPdB)
|
||||||
|
|
||||||
|
else:
|
||||||
|
print 'Bad initial (AIC) P-pick, skip this onset!'
|
||||||
|
Pweight = 4
|
||||||
|
else:
|
||||||
|
print 'run_autopicking: No vertical component data available, skipping station!'
|
||||||
|
return
|
||||||
|
|
||||||
|
if edat is not None and ndat is not None and Pweight < 4:
|
||||||
|
print 'Go on picking S onset ...'
|
||||||
|
print '##################################################'
|
||||||
|
print 'Working on S onset of station %s' % edat[0].stats.station
|
||||||
|
print 'Filtering horizontal traces ...'
|
||||||
|
|
||||||
|
#determine time window for calculating CF after P onset
|
||||||
|
#cuttimesh = [round(refPpick.getpick() + sstart), round(refPpick.getpick() + sstop)]
|
||||||
|
cuttimesh = [round(max([refPpick.getpick() + sstart, 0])), \
|
||||||
|
round(min([refPpick.getpick() + sstop, Lwf]))]
|
||||||
|
|
||||||
|
if algoS == 'ARH':
|
||||||
|
print edat, ndat
|
||||||
|
#re-create stream object including both horizontal components
|
||||||
|
hdat = edat.copy()
|
||||||
|
hdat += ndat
|
||||||
|
h_copy = hdat.copy()
|
||||||
|
#filter and taper data
|
||||||
|
trH1_filt = hdat[0].copy()
|
||||||
|
trH2_filt = hdat[1].copy()
|
||||||
|
trH1_filt.filter('bandpass', freqmin=bph1[0], freqmax=bph1[1], zerophase=False)
|
||||||
|
trH2_filt.filter('bandpass', freqmin=bph1[0], freqmax=bph1[1], zerophase=False)
|
||||||
|
trH1_filt.taper(max_percentage=0.05, type='hann')
|
||||||
|
trH2_filt.taper(max_percentage=0.05, type='hann')
|
||||||
|
h_copy[0].data = trH1_filt.data
|
||||||
|
h_copy[1].data = trH2_filt.data
|
||||||
|
elif algoS == 'AR3':
|
||||||
|
print zdat, edat, ndat
|
||||||
|
#re-create stream object including both horizontal components
|
||||||
|
hdat = zdat.copy()
|
||||||
|
hdat += edat
|
||||||
|
hdat += ndat
|
||||||
|
h_copy = hdat.copy()
|
||||||
|
#filter and taper data
|
||||||
|
trH1_filt = hdat[0].copy()
|
||||||
|
trH2_filt = hdat[1].copy()
|
||||||
|
trH3_filt = hdat[2].copy()
|
||||||
|
trH1_filt.filter('bandpass', freqmin=bph1[0], freqmax=bph1[1], zerophase=False)
|
||||||
|
trH2_filt.filter('bandpass', freqmin=bph1[0], freqmax=bph1[1], zerophase=False)
|
||||||
|
trH3_filt.filter('bandpass', freqmin=bph1[0], freqmax=bph1[1], zerophase=False)
|
||||||
|
trH1_filt.taper(max_percentage=0.05, type='hann')
|
||||||
|
trH2_filt.taper(max_percentage=0.05, type='hann')
|
||||||
|
trH3_filt.taper(max_percentage=0.05, type='hann')
|
||||||
|
h_copy[0].data = trH1_filt.data
|
||||||
|
h_copy[1].data = trH2_filt.data
|
||||||
|
h_copy[2].data = trH3_filt.data
|
||||||
|
##############################################################
|
||||||
|
if algoS == 'ARH':
|
||||||
|
#calculate ARH-CF using subclass ARHcf of class CharcteristicFunction
|
||||||
|
arhcf1 = ARHcf(h_copy, cuttimesh, tpred1h, Sarorder, tdet1h, addnoise) #instance of ARHcf
|
||||||
|
elif algoS == 'AR3':
|
||||||
|
#calculate ARH-CF using subclass AR3cf of class CharcteristicFunction
|
||||||
|
arhcf1 = AR3Ccf(h_copy, cuttimesh, tpred1h, Sarorder, tdet1h, addnoise) #instance of ARHcf
|
||||||
|
##############################################################
|
||||||
|
#calculate AIC-ARH-CF using subclass AICcf of class CharacteristicFunction
|
||||||
|
#class needs stream object => build it
|
||||||
|
tr_arhaic = trH1_filt.copy()
|
||||||
|
tr_arhaic.data = arhcf1.getCF()
|
||||||
|
h_copy[0].data = tr_arhaic.data
|
||||||
|
#calculate ARH-AIC-CF
|
||||||
|
haiccf = AICcf(h_copy, cuttimesh) #instance of AICcf
|
||||||
|
##############################################################
|
||||||
|
#get prelimenary onset time from AIC-HOS-CF using subclass AICPicker of class AutoPicking
|
||||||
|
aicarhpick = AICPicker(haiccf, tsnrh, pickwinS, iplot, None, aictsmoothS)
|
||||||
|
###############################################################
|
||||||
|
#go on with processing if AIC onset passes quality control
|
||||||
|
if aicarhpick.getSlope() >= minAICSslope and aicarhpick.getSNR() >= minAICSSNR:
|
||||||
|
print 'AIC S-pick passes quality control: Slope: %f, SNR: %f' \
|
||||||
|
% (aicarhpick.getSlope(), aicarhpick.getSNR())
|
||||||
|
print 'Go on with refined picking ...'
|
||||||
|
#re-calculate CF from re-filtered trace in vicinity of initial onset
|
||||||
|
cuttimesh2 = [round(aicarhpick.getpick() - Srecalcwin), \
|
||||||
|
round(aicarhpick.getpick() + Srecalcwin)]
|
||||||
|
#re-filter waveform with larger bandpass
|
||||||
|
print 'run_autopicking: re-filtering horizontal traces...'
|
||||||
|
h_copy = hdat.copy()
|
||||||
|
#filter and taper data
|
||||||
|
if algoS == 'ARH':
|
||||||
|
trH1_filt = hdat[0].copy()
|
||||||
|
trH2_filt = hdat[1].copy()
|
||||||
|
trH1_filt.filter('bandpass', freqmin=bph2[0], freqmax=bph2[1], zerophase=False)
|
||||||
|
trH2_filt.filter('bandpass', freqmin=bph2[0], freqmax=bph2[1], zerophase=False)
|
||||||
|
trH1_filt.taper(max_percentage=0.05, type='hann')
|
||||||
|
trH2_filt.taper(max_percentage=0.05, type='hann')
|
||||||
|
h_copy[0].data = trH1_filt.data
|
||||||
|
h_copy[1].data = trH2_filt.data
|
||||||
|
#############################################################
|
||||||
|
arhcf2 = ARHcf(h_copy, cuttimesh2, tpred2h, Sarorder, tdet2h, addnoise) #instance of ARHcf
|
||||||
|
elif algoS == 'AR3':
|
||||||
|
trH1_filt = hdat[0].copy()
|
||||||
|
trH2_filt = hdat[1].copy()
|
||||||
|
trH3_filt = hdat[2].copy()
|
||||||
|
trH1_filt.filter('bandpass', freqmin=bph2[0], freqmax=bph2[1], zerophase=False)
|
||||||
|
trH2_filt.filter('bandpass', freqmin=bph2[0], freqmax=bph2[1], zerophase=False)
|
||||||
|
trH3_filt.filter('bandpass', freqmin=bph2[0], freqmax=bph2[1], zerophase=False)
|
||||||
|
trH1_filt.taper(max_percentage=0.05, type='hann')
|
||||||
|
trH2_filt.taper(max_percentage=0.05, type='hann')
|
||||||
|
trH3_filt.taper(max_percentage=0.05, type='hann')
|
||||||
|
h_copy[0].data = trH1_filt.data
|
||||||
|
h_copy[1].data = trH2_filt.data
|
||||||
|
h_copy[2].data = trH3_filt.data
|
||||||
|
#############################################################
|
||||||
|
arhcf2 = AR3Ccf(h_copy, cuttimesh2, tpred2h, Sarorder, tdet2h, addnoise) #instance of ARHcf
|
||||||
|
|
||||||
|
#get refined onset time from CF2 using class Picker
|
||||||
|
refSpick = PragPicker(arhcf2, tsnrh, pickwinS, iplot, ausS, tsmoothS, aicarhpick.getpick())
|
||||||
|
#############################################################
|
||||||
|
#quality assessment
|
||||||
|
#get earliest and latest possible pick and symmetrized uncertainty
|
||||||
|
h_copy[0].data = trH1_filt.data
|
||||||
|
[lpickS1, epickS1, Serror1] = earllatepicker(h_copy, nfacS, tsnrh, refSpick.getpick(), iplot)
|
||||||
|
h_copy[0].data = trH2_filt.data
|
||||||
|
[lpickS2, epickS2, Serror2] = earllatepicker(h_copy, nfacS, tsnrh, refSpick.getpick(), iplot)
|
||||||
|
if algoS == 'ARH':
|
||||||
|
#get earliest pick of both earliest possible picks
|
||||||
|
epick = [epickS1, epickS2]
|
||||||
|
lpick = [lpickS1, lpickS2]
|
||||||
|
pickerr = [Serror1, Serror2]
|
||||||
|
ipick =np.argmin([epickS1, epickS2])
|
||||||
|
elif algoS == 'AR3':
|
||||||
|
[lpickS3, epickS3, Serror3] = earllatepicker(h_copy, nfacS, tsnrh, refSpick.getpick(), iplot)
|
||||||
|
#get earliest pick of all three picks
|
||||||
|
epick = [epickS1, epickS2, epickS3]
|
||||||
|
lpick = [lpickS1, lpickS2, lpickS3]
|
||||||
|
pickerr = [Serror1, Serror2, Serror3]
|
||||||
|
ipick =np.argmin([epickS1, epickS2, epickS3])
|
||||||
|
epickS = epick[ipick]
|
||||||
|
lpickS = lpick[ipick]
|
||||||
|
Serror = pickerr[ipick]
|
||||||
|
|
||||||
|
#get SNR
|
||||||
|
[SNRS, SNRSdB, Snoiselevel] = getSNR(h_copy, tsnrh, refSpick.getpick())
|
||||||
|
|
||||||
|
#weight S-onset using symmetric error
|
||||||
|
if Serror <= timeerrorsS[0]:
|
||||||
|
Sweight = 0
|
||||||
|
elif Serror > timeerrorsS[0] and Serror <= timeerrorsS[1]:
|
||||||
|
Sweight = 1
|
||||||
|
elif Perror > timeerrorsS[1] and Serror <= timeerrorsS[2]:
|
||||||
|
Sweight = 2
|
||||||
|
elif Serror > timeerrorsS[2] and Serror <= timeerrorsS[3]:
|
||||||
|
Sweight = 3
|
||||||
|
elif Serror > timeerrorsS[3]:
|
||||||
|
Sweight = 4
|
||||||
|
|
||||||
|
print 'run_autopicking: S-weight: %d, SNR: %f, SNR[dB]: %f' % (Sweight, SNRS, SNRSdB)
|
||||||
|
|
||||||
|
else:
|
||||||
|
print 'Bad initial (AIC) S-pick, skip this onset!'
|
||||||
|
Sweight = 4
|
||||||
|
|
||||||
|
else:
|
||||||
|
print 'run_autopicking: No horizontal component data available, skipping S picking!'
|
||||||
|
return
|
||||||
|
|
||||||
|
##############################################################
|
||||||
|
if iplot > 0:
|
||||||
|
#plot vertical trace
|
||||||
|
plt.figure()
|
||||||
|
plt.subplot(3,1,1)
|
||||||
|
tdata = np.arange(0, tr_filt.stats.npts / tr_filt.stats.sampling_rate, tr_filt.stats.delta)
|
||||||
|
p1, = plt.plot(tdata, tr_filt.data/max(tr_filt.data), 'k')
|
||||||
|
p2, = plt.plot(cf1.getTimeArray(), cf1.getCF() / max(cf1.getCF()), 'b')
|
||||||
|
p3, = plt.plot(cf2.getTimeArray(), cf2.getCF() / max(cf2.getCF()), 'm')
|
||||||
|
p4, = plt.plot([aicpick.getpick(), aicpick.getpick()], [-1, 1], 'r')
|
||||||
|
plt.plot([aicpick.getpick()-0.5, aicpick.getpick()+0.5], [1, 1], 'r')
|
||||||
|
plt.plot([aicpick.getpick()-0.5, aicpick.getpick()+0.5], [-1, -1], 'r')
|
||||||
|
p5, = plt.plot([refPpick.getpick(), refPpick.getpick()], [-1.3, 1.3], 'r', linewidth=2)
|
||||||
|
plt.plot([refPpick.getpick()-0.5, refPpick.getpick()+0.5], [1.3, 1.3], 'r', linewidth=2)
|
||||||
|
plt.plot([refPpick.getpick()-0.5, refPpick.getpick()+0.5], [-1.3, -1.3], 'r', linewidth=2)
|
||||||
|
plt.plot([lpickP, lpickP], [-1.1, 1.1], 'r--')
|
||||||
|
plt.plot([epickP, epickP], [-1.1, 1.1], 'r--')
|
||||||
|
plt.yticks([])
|
||||||
|
plt.ylim([-1.5, 1.5])
|
||||||
|
plt.ylabel('Normalized Counts')
|
||||||
|
plt.title('%s, %s, P Weight=%d, SNR=%7.2f, SNR[dB]=%7.2f' % (tr_filt.stats.station, \
|
||||||
|
tr_filt.stats.channel, Pweight, SNRP, SNRPdB))
|
||||||
|
plt.suptitle(tr_filt.stats.starttime)
|
||||||
|
plt.legend([p1, p2, p3, p4, p5], ['Data', 'CF1', 'CF2', 'Initial P Onset', 'Final P Pick'])
|
||||||
|
#plot horizontal traces
|
||||||
|
plt.subplot(3,1,2)
|
||||||
|
th1data = np.arange(0, trH1_filt.stats.npts / trH1_filt.stats.sampling_rate, trH1_filt.stats.delta)
|
||||||
|
p21, = plt.plot(th1data, trH1_filt.data/max(trH1_filt.data), 'k')
|
||||||
|
p22, = plt.plot(arhcf1.getTimeArray(), arhcf1.getCF()/max(arhcf1.getCF()), 'b')
|
||||||
|
p23, = plt.plot(arhcf2.getTimeArray(), arhcf2.getCF()/max(arhcf2.getCF()), 'm')
|
||||||
|
p24, = plt.plot([aicarhpick.getpick(), aicarhpick.getpick()], [-1, 1], 'g')
|
||||||
|
plt.plot([aicarhpick.getpick() - 0.5, aicarhpick.getpick() + 0.5], [1, 1], 'g')
|
||||||
|
plt.plot([aicarhpick.getpick() - 0.5, aicarhpick.getpick() + 0.5], [-1, -1], 'g')
|
||||||
|
p25, = plt.plot([refSpick.getpick(), refSpick.getpick()], [-1.3, 1.3], 'g', linewidth=2)
|
||||||
|
plt.plot([refSpick.getpick() - 0.5, refSpick.getpick() + 0.5], [1.3, 1.3], 'g', linewidth=2)
|
||||||
|
plt.plot([refSpick.getpick() - 0.5, refSpick.getpick() + 0.5], [-1.3, -1.3], 'g', linewidth=2)
|
||||||
|
plt.plot([lpickS, lpickS], [-1.1, 1.1], 'g--')
|
||||||
|
plt.plot([epickS, epickS], [-1.1, 1.1], 'g--')
|
||||||
|
plt.yticks([])
|
||||||
|
plt.ylim([-1.5, 1.5])
|
||||||
|
plt.ylabel('Normalized Counts')
|
||||||
|
plt.title('%s, S Weight=%d, SNR=%7.2f, SNR[dB]=%7.2f' % (trH1_filt.stats.channel, \
|
||||||
|
Sweight, SNRS, SNRSdB))
|
||||||
|
plt.suptitle(trH1_filt.stats.starttime)
|
||||||
|
plt.legend([p21, p22, p23, p24, p25], ['Data', 'CF1', 'CF2', 'Initial S Onset', 'Final S Pick'])
|
||||||
|
plt.subplot(3,1,3)
|
||||||
|
th2data = np.arange(0, trH2_filt.stats.npts / trH2_filt.stats.sampling_rate, trH2_filt.stats.delta)
|
||||||
|
plt.plot(th2data, trH2_filt.data/max(trH2_filt.data), 'k')
|
||||||
|
plt.plot(arhcf1.getTimeArray(), arhcf1.getCF()/max(arhcf1.getCF()), 'b')
|
||||||
|
plt.plot(arhcf2.getTimeArray(), arhcf2.getCF()/max(arhcf2.getCF()), 'm')
|
||||||
|
plt.plot([aicarhpick.getpick(), aicarhpick.getpick()], [-1, 1], 'g')
|
||||||
|
plt.plot([aicarhpick.getpick() - 0.5, aicarhpick.getpick() + 0.5], [1, 1], 'g')
|
||||||
|
plt.plot([aicarhpick.getpick() - 0.5, aicarhpick.getpick() + 0.5], [-1, -1], 'g')
|
||||||
|
plt.plot([refSpick.getpick(), refSpick.getpick()], [-1.3, 1.3], 'g', linewidth=2)
|
||||||
|
plt.plot([refSpick.getpick() - 0.5, refSpick.getpick() + 0.5], [1.3, 1.3], 'g', linewidth=2)
|
||||||
|
plt.plot([refSpick.getpick() - 0.5, refSpick.getpick() + 0.5], [-1.3, -1.3], 'g', linewidth=2)
|
||||||
|
plt.plot([lpickS, lpickS], [-1.1, 1.1], 'g--')
|
||||||
|
plt.plot([epickS, epickS], [-1.1, 1.1], 'g--')
|
||||||
|
plt.yticks([])
|
||||||
|
plt.ylim([-1.5, 1.5])
|
||||||
|
plt.xlabel('Time [s] after %s' % tr_filt.stats.starttime)
|
||||||
|
plt.ylabel('Normalized Counts')
|
||||||
|
plt.title(trH2_filt.stats.channel)
|
||||||
|
plt.show()
|
||||||
|
raw_input()
|
||||||
|
plt.close()
|
||||||
Loading…
Reference in New Issue
Block a user