From aba3997b20e7cab43560a70848e20d052a4157d2 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Ludger=20K=C3=BCperkoch?= Date: Mon, 22 Jun 2015 12:39:29 +0200 Subject: [PATCH] Modified earllatepicker: Mean is removed from trace calculated from noise + signal window. --- pylot/core/pick/utils.py | 81 ++++++++++++++++++---------------------- 1 file changed, 36 insertions(+), 45 deletions(-) diff --git a/pylot/core/pick/utils.py b/pylot/core/pick/utils.py index 25eef111..47d06b69 100644 --- a/pylot/core/pick/utils.py +++ b/pylot/core/pick/utils.py @@ -53,6 +53,9 @@ def earllatepicker(X, nfac, TSNR, Pick1, iplot=None): inoise = getnoisewin(t, Pick1, TSNR[0], TSNR[1]) # get signal window isignal = getsignalwin(t, Pick1, TSNR[2]) + # remove mean + meanwin = np.hstack((inoise, isignal)) + x = x - np.mean(x[meanwin]) # calculate noise level nlevel = np.sqrt(np.mean(np.square(x[inoise]))) * nfac # get time where signal exceeds nlevel @@ -412,34 +415,28 @@ def wadaticheck(pickdic, dttolerance, iplot): SPtimes = [] for key in pickdic: if pickdic[key]['P']['weight'] < 4 and pickdic[key]['S']['weight'] < 4: - # calculate S-P time - spt = pickdic[key]['S']['mpp'] - pickdic[key]['P']['mpp'] - # add S-P time to dictionary - pickdic[key]['SPt'] = spt - # add P onsets and corresponding S-P times to list - UTCPpick = UTCDateTime(pickdic[key]['P']['mpp']) - UTCDateTime(1970, - 1, 1, - 0, 0, - 0) - UTCSpick = UTCDateTime(pickdic[key]['S']['mpp']) - UTCDateTime(1970, - 1, 1, - 0, 0, - 0) - Ppicks.append(UTCPpick) - Spicks.append(UTCSpick) - SPtimes.append(spt) - vpvs.append(UTCPpick / UTCSpick) + # calculate S-P time + spt = pickdic[key]['S']['mpp'] - pickdic[key]['P']['mpp'] + # add S-P time to dictionary + pickdic[key]['SPt'] = spt + # add P onsets and corresponding S-P times to list + UTCPpick = UTCDateTime(pickdic[key]['P']['mpp']) + UTCSpick = UTCDateTime(pickdic[key]['S']['mpp']) + Ppicks.append(UTCPpick.timestamp) + Spicks.append(UTCSpick.timestamp) + SPtimes.append(spt) + if len(SPtimes) >= 3: - # calculate slope - p1 = np.polyfit(Ppicks, SPtimes, 1) - wdfit = np.polyval(p1, Ppicks) + # calculate slope + p1 = np.polyfit(Ppicks, SPtimes, 1) + wdfit = np.polyval(p1, Ppicks) wfitflag = 0 # calculate vp/vs ratio before check vpvsr = p1[0] + 1 print 'wadaticheck: Average Vp/Vs ratio before check:', vpvsr - + checkedPpicks = [] checkedSpicks = [] checkedSPtimes = [] @@ -457,23 +454,19 @@ def wadaticheck(pickdic, dttolerance, iplot): pickdic[key]['S']['weight'] = 9 else: marker = 'goodWadatiCheck' - checkedPpick = UTCDateTime(pickdic[key]['P']['mpp']) - \ - UTCDateTime(1970, 1, 1, 0, 0, 0) - checkedPpicks.append(checkedPpick) - checkedSpick = UTCDateTime(pickdic[key]['S']['mpp']) - \ - UTCDateTime(1970, 1, 1, 0, 0, 0) - checkedSpicks.append(checkedSpick) - checkedSPtime = pickdic[key]['S']['mpp'] - \ - pickdic[key]['P']['mpp'] + checkedPpick = UTCDateTime(pickdic[key]['P']['mpp']) + checkedPpicks.append(checkedPpick.timestamp) + checkedSpick = UTCDateTime(pickdic[key]['S']['mpp']) + checkedSpicks.append(checkedSpick.timestamp) + checkedSPtime = pickdic[key]['S']['mpp'] - pickdic[key]['P']['mpp'] checkedSPtimes.append(checkedSPtime) - checkedvpvs.append(checkedPpick / checkedSpick) pickdic[key]['S']['marked'] = marker - # calculate new slope - p2 = np.polyfit(checkedPpicks, checkedSPtimes, 1) - wdfit2 = np.polyval(p2, checkedPpicks) + # calculate new slope + p2 = np.polyfit(checkedPpicks, checkedSPtimes, 1) + wdfit2 = np.polyval(p2, checkedPpicks) # calculate vp/vs ratio after check cvpvsr = p2[0] + 1 @@ -482,26 +475,24 @@ def wadaticheck(pickdic, dttolerance, iplot): checkedonsets = pickdic else: - print 'wadaticheck: Not enough S-P times available for reliable regression!' + print 'wadaticheck: Not enough S-P times available for reliable regression!' print 'Skip wadati check!' wfitflag = 1 # plot results if iplot > 1: - plt.figure(iplot) - f1, = plt.plot(Ppicks, SPtimes, 'ro') + plt.figure(iplot) + f1, = plt.plot(Ppicks, SPtimes, 'ro') if wfitflag == 0: - f2, = plt.plot(Ppicks, wdfit, 'k') - f3, = plt.plot(checkedPpicks, checkedSPtimes, 'ko') - f4, = plt.plot(checkedPpicks, wdfit2, 'g') + f2, = plt.plot(Ppicks, wdfit, 'k') + f3, = plt.plot(checkedPpicks, checkedSPtimes, 'ko') + f4, = plt.plot(checkedPpicks, wdfit2, 'g') plt.ylabel('S-P Times [s]') plt.xlabel('P Times [s]') - plt.title( - 'Wadati-Diagram, %d S-P Times, Vp/Vs(old)=%5.2f, Vp/Vs(checked)=%5.2f' \ - % (len(SPtimes), vpvsr, cvpvsr)) - plt.legend([f1, f2, f3, f4], - ['Skipped S-Picks', 'Wadati 1', 'Reliable S-Picks', \ - 'Wadati 2'], loc='best') + plt.title('Wadati-Diagram, %d S-P Times, Vp/Vs(raw)=%5.2f, Vp/Vs(checked)=%5.2f' \ + % (len(SPtimes), vpvsr, cvpvsr)) + plt.legend([f1, f2, f3, f4], ['Skipped S-Picks', 'Wadati 1', 'Reliable S-Picks', \ + 'Wadati 2'], loc='best') plt.show() raw_input() plt.close(iplot)