Debugging: If no P-Pick was determined, no plot of of horizontal data comes up.

This commit is contained in:
Ludger Küperkoch 2015-06-22 14:59:57 +02:00
parent eb153679ba
commit 99a5a4499a

View File

@ -14,13 +14,12 @@ import numpy as np
from pylot.core.pick.Picker import * from pylot.core.pick.Picker import *
from pylot.core.pick.CharFuns import * from pylot.core.pick.CharFuns import *
def run_autopicking(wfstream, pickparam): def run_autopicking(wfstream, pickparam):
""" """
param: wfstream :param: wfstream
:type: `~obspy.core.stream.Stream` :type: `~obspy.core.stream.Stream`
param: pickparam :param: pickparam
:type: container of picking parameters from input file, :type: container of picking parameters from input file,
usually autoPyLoT.in usually autoPyLoT.in
""" """
@ -95,6 +94,7 @@ def run_autopicking(wfstream, pickparam):
aicSflag = 0 aicSflag = 0
aicPflag = 0 aicPflag = 0
Sflag = 0
# split components # split components
zdat = wfstream.select(component="Z") zdat = wfstream.select(component="Z")
@ -218,10 +218,12 @@ def run_autopicking(wfstream, pickparam):
print 'run_autopicking: P-weight: %d, SNR: %f, SNR[dB]: %f, ' \ print 'run_autopicking: P-weight: %d, SNR: %f, SNR[dB]: %f, ' \
'Polarity: %s' % (Pweight, SNRP, SNRPdB, FM) 'Polarity: %s' % (Pweight, SNRP, SNRPdB, FM)
Sflag = 1
else: else:
print 'Bad initial (AIC) P-pick, skip this onset!' print 'Bad initial (AIC) P-pick, skip this onset!'
print 'AIC-SNR=', aicpick.getSNR(), 'AIC-Slope=', aicpick.getSlope() print 'AIC-SNR=', aicpick.getSNR(), 'AIC-Slope=', aicpick.getSlope()
Sflag = 0
else: else:
print 'run_autopicking: No vertical component data available, ' \ print 'run_autopicking: No vertical component data available, ' \
@ -295,7 +297,7 @@ def run_autopicking(wfstream, pickparam):
# class needs stream object => build it # class needs stream object => build it
tr_arhaic = trH1_filt.copy() tr_arhaic = trH1_filt.copy()
tr_arhaic.data = arhcf1.getCF() tr_arhaic.data = arhcf1.getCF()
h_copy[0].data = tr_arhaic.data h_copy[0].data = tr_arhaic.data
# calculate ARH-AIC-CF # calculate ARH-AIC-CF
haiccf = AICcf(h_copy, cuttimesh) # instance of AICcf haiccf = AICcf(h_copy, cuttimesh) # instance of AICcf
############################################################## ##############################################################
@ -463,108 +465,118 @@ def run_autopicking(wfstream, pickparam):
plt.ylabel('Normalized Counts') plt.ylabel('Normalized Counts')
plt.suptitle(tr_filt.stats.starttime) plt.suptitle(tr_filt.stats.starttime)
# plot horizontal traces if len(edat) > 1 and len(ndat) > 1 and Sflag == 1:
plt.subplot(3, 1, 2) # plot horizontal traces
th1data = np.arange(0, plt.subplot(3, 1, 2)
trH1_filt.stats.npts / th1data = np.arange(0,
trH1_filt.stats.sampling_rate, trH1_filt.stats.npts /
trH1_filt.stats.delta) trH1_filt.stats.sampling_rate,
# check equal length of arrays, sometimes they are different!? trH1_filt.stats.delta)
wfldiff = len(trH1_filt.data) - len(th1data) # check equal length of arrays, sometimes they are different!?
if wfldiff < 0: wfldiff = len(trH1_filt.data) - len(th1data)
th1data = th1data[0:len(th1data) - abs(wfldiff)] if wfldiff < 0:
p21, = plt.plot(th1data, trH1_filt.data / max(trH1_filt.data), 'k') th1data = th1data[0:len(th1data) - abs(wfldiff)]
if Pweight < 4: p21, = plt.plot(th1data, trH1_filt.data / max(trH1_filt.data), 'k')
p22, = plt.plot(arhcf1.getTimeArray(), if Pweight < 4:
arhcf1.getCF() / max(arhcf1.getCF()), 'b') p22, = plt.plot(arhcf1.getTimeArray(),
if aicSflag == 1: arhcf1.getCF() / max(arhcf1.getCF()), 'b')
p23, = plt.plot(arhcf2.getTimeArray(), if aicSflag == 1:
arhcf2.getCF() / max(arhcf2.getCF()), 'm') p23, = plt.plot(arhcf2.getTimeArray(),
p24, = plt.plot([aicarhpick.getpick(), aicarhpick.getpick()], arhcf2.getCF() / max(arhcf2.getCF()), 'm')
[-1, 1], 'g') p24, = plt.plot([aicarhpick.getpick(), aicarhpick.getpick()],
plt.plot( [-1, 1], 'g')
[aicarhpick.getpick() - 0.5, aicarhpick.getpick() + 0.5], plt.plot(
[1, 1], 'g') [aicarhpick.getpick() - 0.5, aicarhpick.getpick() + 0.5],
plt.plot( [1, 1], 'g')
[aicarhpick.getpick() - 0.5, aicarhpick.getpick() + 0.5], plt.plot(
[-1, -1], 'g') [aicarhpick.getpick() - 0.5, aicarhpick.getpick() + 0.5],
p25, = plt.plot([refSpick.getpick(), refSpick.getpick()], [-1, -1], 'g')
[-1.3, 1.3], 'g', linewidth=2) p25, = plt.plot([refSpick.getpick(), refSpick.getpick()],
plt.plot([refSpick.getpick() - 0.5, refSpick.getpick() + 0.5], [-1.3, 1.3], 'g', linewidth=2)
[1.3, 1.3], 'g', linewidth=2) plt.plot([refSpick.getpick() - 0.5, refSpick.getpick() + 0.5],
plt.plot([refSpick.getpick() - 0.5, refSpick.getpick() + 0.5], [1.3, 1.3], 'g', linewidth=2)
[-1.3, -1.3], 'g', linewidth=2) plt.plot([refSpick.getpick() - 0.5, refSpick.getpick() + 0.5],
plt.plot([lpickS, lpickS], [-1.1, 1.1], 'g--') [-1.3, -1.3], 'g', linewidth=2)
plt.plot([epickS, epickS], [-1.1, 1.1], 'g--') plt.plot([lpickS, lpickS], [-1.1, 1.1], 'g--')
plt.legend([p21, p22, p23, p24, p25], plt.plot([epickS, epickS], [-1.1, 1.1], 'g--')
['Data', 'CF1', 'CF2', 'Initial S Onset', plt.legend([p21, p22, p23, p24, p25],
'Final S Pick']) ['Data', 'CF1', 'CF2', 'Initial S Onset',
plt.title('%s, S Weight=%d, SNR=%7.2f, SNR[dB]=%7.2f' % ( 'Final S Pick'])
trH1_filt.stats.channel, plt.title('%s, S Weight=%d, SNR=%7.2f, SNR[dB]=%7.2f' % (
Sweight, SNRS, SNRSdB)) trH1_filt.stats.channel,
else: Sweight, SNRS, SNRSdB))
plt.legend([p21, p22], ['Data', 'CF1']) else:
plt.title('%s, S Weight=%d, SNR=None, SNRdB=None' % ( plt.legend([p21, p22], ['Data', 'CF1'])
trH1_filt.stats.channel, Sweight)) plt.title('%s, S Weight=%d, SNR=None, SNRdB=None' % (
plt.yticks([]) trH1_filt.stats.channel, Sweight))
plt.ylim([-1.5, 1.5]) plt.yticks([])
plt.ylabel('Normalized Counts') plt.ylim([-1.5, 1.5])
plt.suptitle(trH1_filt.stats.starttime) plt.ylabel('Normalized Counts')
plt.suptitle(trH1_filt.stats.starttime)
plt.subplot(3, 1, 3) plt.subplot(3, 1, 3)
th2data = np.arange(0, th2data = np.arange(0,
trH2_filt.stats.npts / trH2_filt.stats.npts /
trH2_filt.stats.sampling_rate, trH2_filt.stats.sampling_rate,
trH2_filt.stats.delta) trH2_filt.stats.delta)
# check equal length of arrays, sometimes they are different!? # check equal length of arrays, sometimes they are different!?
wfldiff = len(trH2_filt.data) - len(th2data) wfldiff = len(trH2_filt.data) - len(th2data)
if wfldiff < 0: if wfldiff < 0:
th2data = th2data[0:len(th2data) - abs(wfldiff)] th2data = th2data[0:len(th2data) - abs(wfldiff)]
plt.plot(th2data, trH2_filt.data / max(trH2_filt.data), 'k') plt.plot(th2data, trH2_filt.data / max(trH2_filt.data), 'k')
if Pweight < 4: if Pweight < 4:
p22, = plt.plot(arhcf1.getTimeArray(), p22, = plt.plot(arhcf1.getTimeArray(),
arhcf1.getCF() / max(arhcf1.getCF()), 'b') arhcf1.getCF() / max(arhcf1.getCF()), 'b')
if aicSflag == 1: if aicSflag == 1:
p23, = plt.plot(arhcf2.getTimeArray(), p23, = plt.plot(arhcf2.getTimeArray(),
arhcf2.getCF() / max(arhcf2.getCF()), 'm') arhcf2.getCF() / max(arhcf2.getCF()), 'm')
p24, = plt.plot([aicarhpick.getpick(), aicarhpick.getpick()], p24, = plt.plot([aicarhpick.getpick(), aicarhpick.getpick()],
[-1, 1], 'g') [-1, 1], 'g')
plt.plot( plt.plot(
[aicarhpick.getpick() - 0.5, aicarhpick.getpick() + 0.5], [aicarhpick.getpick() - 0.5, aicarhpick.getpick() + 0.5],
[1, 1], 'g') [1, 1], 'g')
plt.plot( plt.plot(
[aicarhpick.getpick() - 0.5, aicarhpick.getpick() + 0.5], [aicarhpick.getpick() - 0.5, aicarhpick.getpick() + 0.5],
[-1, -1], 'g') [-1, -1], 'g')
p25, = plt.plot([refSpick.getpick(), refSpick.getpick()], p25, = plt.plot([refSpick.getpick(), refSpick.getpick()],
[-1.3, 1.3], 'g', linewidth=2) [-1.3, 1.3], 'g', linewidth=2)
plt.plot([refSpick.getpick() - 0.5, refSpick.getpick() + 0.5], plt.plot([refSpick.getpick() - 0.5, refSpick.getpick() + 0.5],
[1.3, 1.3], 'g', linewidth=2) [1.3, 1.3], 'g', linewidth=2)
plt.plot([refSpick.getpick() - 0.5, refSpick.getpick() + 0.5], plt.plot([refSpick.getpick() - 0.5, refSpick.getpick() + 0.5],
[-1.3, -1.3], 'g', linewidth=2) [-1.3, -1.3], 'g', linewidth=2)
plt.plot([lpickS, lpickS], [-1.1, 1.1], 'g--') plt.plot([lpickS, lpickS], [-1.1, 1.1], 'g--')
plt.plot([epickS, epickS], [-1.1, 1.1], 'g--') plt.plot([epickS, epickS], [-1.1, 1.1], 'g--')
plt.legend([p21, p22, p23, p24, p25], plt.legend([p21, p22, p23, p24, p25],
['Data', 'CF1', 'CF2', 'Initial S Onset', ['Data', 'CF1', 'CF2', 'Initial S Onset',
'Final S Pick']) 'Final S Pick'])
else: else:
plt.legend([p21, p22], ['Data', 'CF1']) plt.legend([p21, p22], ['Data', 'CF1'])
plt.yticks([]) plt.yticks([])
plt.ylim([-1.5, 1.5]) plt.ylim([-1.5, 1.5])
plt.xlabel('Time [s] after %s' % tr_filt.stats.starttime) plt.xlabel('Time [s] after %s' % tr_filt.stats.starttime)
plt.ylabel('Normalized Counts') plt.ylabel('Normalized Counts')
plt.title(trH2_filt.stats.channel) plt.title(trH2_filt.stats.channel)
plt.show() plt.show()
raw_input() raw_input()
plt.close() plt.close()
########################################################################## ##########################################################################
# calculate "real" onset times
if mpickP is not None:
lpickP = zdat[0].stats.starttime + lpickP
epickP = zdat[0].stats.starttime + epickP
mpickP = zdat[0].stats.starttime + mpickP
if mpickS is not None:
lpickS = edat[0].stats.starttime + lpickS
epickS = edat[0].stats.starttime + epickS
mpickS = edat[0].stats.starttime + mpickS
# create dictionary # create dictionary
# for P phase # for P phase
phase = 'P' phase = 'P'
phasepick = {'lpp': lpickP, 'epp': epickP, 'mpp': mpickP, 'spe': Perror, \ phasepick = {'lpp': lpickP, 'epp': epickP, 'mpp': mpickP, 'spe': Perror, \
'snr': SNRP, 'snrdb': SNRPdB, 'weight': Pweight, 'fm': FM} 'snr': SNRP, 'snrdb': SNRPdB, 'weight': Pweight, 'fm': FM}
picks = {phase: phasepick} picks = {phase: phasepick}
#stationID = tr_filt.stats.station
#onsets = {stationID: picks}
# add S phase # add S phase
phase = 'S' phase = 'S'
phasepick = {'lpp': lpickS, 'epp': epickS, 'mpp': mpickS, 'spe': Serror, \ phasepick = {'lpp': lpickS, 'epp': epickS, 'mpp': mpickS, 'spe': Serror, \