From ca809c467376f214b5fc981bdb5c7801558c5cfc Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Ludger=20K=C3=BCperkoch?= Date: Fri, 5 Dec 2014 16:14:03 +0100 Subject: [PATCH] Initial version of new class of methods for automatic picking, AICPicker is running but without quality attributes --- pylot/core/pick/Picker.py | 183 ++++++++++++++++++++++++++++++++++++++ 1 file changed, 183 insertions(+) create mode 100644 pylot/core/pick/Picker.py diff --git a/pylot/core/pick/Picker.py b/pylot/core/pick/Picker.py new file mode 100644 index 00000000..5e2f6f6d --- /dev/null +++ b/pylot/core/pick/Picker.py @@ -0,0 +1,183 @@ +# -*- coding: utf-8 -*- +""" +Created Dec 2014 +Implementation of the picking algorithms published and described in: + +Küperkoch, L., Meier, T., Lee, J., Friederich, W., & Egelados Working Group, 2010: +Automated determination of P-phase arrival times at regional and local distances +using higher order statistics, Geophys. J. Int., 181, 1159-1170 + +Küperkoch, L., Meier, T., Brüstle, A., Lee, J., Friederich, W., & Egelados +Working Group, 2012: Automated determination of S-phase arrival times using +autoregressive prediction: application ot local and regional distances, Geophys. J. Int., +188, 687-702. + +:author: MAGS2 EP3 working group / Ludger Küperkoch +""" +import numpy as np +import matplotlib.pyplot as plt +import pdb + +class AutoPicking(object): + ''' + Superclass of different, automated picking algorithms applied on a CF determined + using AIC, HOS, or AR prediction. + ''' + def __init__(self, cf, Tcf, dt, Tslope, aerr, TSNR, PickWindow, peps, Tsmooth): + ''' + :param: cf, characteristic function, on which the picking algorithm is applied + :type: array + + :param: Tcf, corresponding time array + :type: array + + :param: dt, sampling interval [s] + :type: float + + :param: Tslope, length of time window after pick used to determine slope + for quality estimation [s] + :type: float + + :param: aerr (adjusted error), percentage of maximum of CF to determine slope for quality estimation + :type: int + + :param: TSNR, length of time windows around pick used to determine SNR [s] + :type: tuple (T_noise, T_gap, T_signal) + + :param: PickWindow, length of pick window [s] + :type: float + + :param: peps, find local minimum at i if aic(i-1)*(1+peps) >= aic(i) + :type: float + + :param: Tsmooth, length of moving smoothing window to calculate smoothed CF [s] + :type: float + ''' + + self.cf = cf + self.Tcf = Tcf + self.dt = dt + self.setTslope(Tslope) + self.setaerr(aerr) + self.setTSNR(TSNR) + self.setPickWindow(PickWindow) + self.setpeps(peps) + self.setTsmooth(Tsmooth) + self.calcPick() + + def __str__(self): + return '''\n\t{name} object:\n + TSlope:\t{Tslope}\n + aerr:\t{aerr}\n + TSNR:\t\t\t{TSNR}\n + PickWindow:\t{PickWindow}\n + peps:\t{peps}\n + Tsmooth:\t{Tsmooth}\n + '''.format(name=type(self).__name__, + Tslope=self.getTslope(), + aerr=self.getaerr(), + TSNR=self.getTSNR(), + PickWindow=self.getPickWindow(), + peps=self.getpeps(), + Tsmooth=self.getTsmooth()) + + def getTslope(self): + return self.Tslope + + def setTslope(self, Tslope): + self.Tslope = Tslope + + def getaerr(self): + return self.aerr + + def setaerr(self, aerr): + self.aerr = aerr + + def getTSNR(self): + return self.TSNR + + def setTSNR(self, TSNR): + self.TSNR = TSNR + + def getPickWindow(self): + return self.PickWindow + + def setPickWindow(self, PickWindow): + self.PickWindow = PickWindow + + def getpeps(self): + return self.peps + + def setpeps(self, peps): + self.peps = peps + + def setTsmooth(self, Tsmooth): + self.Tsmooth = Tsmooth + + def getTsmooth(self): + return self.Tsmooth + + def getpick(self): + return self.Pick + + def calcPick(self): + self.Pick = Pick + +class AICPicker(AutoPicking): + ''' + Method to derive onset time of arriving phase based on CF + derived from AIC. + ''' + + def calcPick(self): + + print 'Get onset (pick) from AIC-CF ...' + + #taper AIC-CF to get rid off side maxima + tap = np.hanning(len(self.cf)) + aic = tap * self.cf + max(abs(self.cf)) + #get maximum of CF as starting point + icfmax = np.argmax(aic) + + #smooth CF + aicsmooth = np.zeros(len(aic)) + ismooth = round(self.Tsmooth / self.dt) + if len(aic) < ismooth: + print 'AICPicker: Tsmooth larger than AIC function!' + self.Pick = -1 + return self.Pick + else: + self.Pick = -1 + for i in range(1, len(aic)): + if i > ismooth: + ii1 = i - ismooth + aicsmooth[i] = aicsmooth[i - 1] + (aic[i] - aic[ii1]) / ismooth + else: + aicsmooth[i] = np.mean(aic[0:i]) + + #find common, local minimum in front of maximum + #of smoothed and unsmoothed AIC-CF + lpickwindow = int(round(self.PickWindow / self.dt)) + for i in range(icfmax - 1, max([icfmax - lpickwindow, 2]), -1): + if aic[i - 1] * (1 + self.peps) >= aic[i]: + if aicsmooth[i - 1] * (1 + self.peps) >= aicsmooth[i]: + self.Pick = self.Tcf[i] + break + + #try again with larger peps if picking failed + if self.Pick < 0: + peps2 = self.peps + 0.01 + for i in range(icfmax - 1, max([icfmax - lpickwindow, 2]), -1): + if aic[i - 1] * (1 + peps2) >= aic[i]: + if aicsmooth[i - 1] * (1 + peps2) >= aicsmooth[i]: + self.Pick = self.Tcf[i] + break + +class PragPicker(AutoPicking): + ''' + Method of pragmatic picking exploiting information given by CF. + ''' + + def calcPick(self): + + print 'Get onset (pick) from HOS- or AR-CF using pragmatic picking algorithm ...'