Merge remote-tracking branch 'origin/develop' into develop

pylot/core/pick/picker.py

@@ -173,21 +173,14 @@ class AICPicker(AutoPicker):
                     aicsmooth[i] = aicsmooth[i - 1] + (aic[i] - aic[ii1]) / ismooth
                 else:
                     aicsmooth[i] = np.mean(aic[1: i])
-        # remove offset
+        # remove offset in AIC function
         offset = abs(min(aic) - min(aicsmooth))
         aicsmooth = aicsmooth - offset
-        # get maximum of 1st derivative of AIC-CF (more stable!) as starting point
-        diffcf = np.diff(aicsmooth)
-        # find NaN's
-        nn = np.isnan(diffcf)
-        if len(nn) > 1:
-            diffcf[nn] = 0
-        # taper CF to get rid off side maxima
-        tap = np.hanning(len(diffcf))
-        diffcf = tap * diffcf * max(abs(aicsmooth))
-        icfmax = np.argmax(diffcf)
+        # get maximum of HOS/AR-CF as startimg point for searching
+        # minimum in AIC function 
+        icfmax = np.argmax(self.Data[0].data)
 
-        # find minimum in AIC-CF front of maximum
+        # find minimum in AIC-CF front of maximum of HOS/AR-CF
         lpickwindow = int(round(self.PickWindow / self.dt))
         for i in range(icfmax - 1, max([icfmax - lpickwindow, 2]), -1):
             if aicsmooth[i - 1] >= aicsmooth[i]:
@@ -196,6 +189,14 @@ class AICPicker(AutoPicker):
         # if no minimum could be found:
         # search in 1st derivative of AIC-CF
         if self.Pick is None:
+            diffcf = np.diff(aicsmooth)
+            # find NaN's
+            nn = np.isnan(diffcf)
+            if len(nn) > 1:
+                diffcf[nn] = 0
+            # taper CF to get rid off side maxima
+            tap = np.hanning(len(diffcf))
+            diffcf = tap * diffcf * max(abs(aicsmooth))
             for i in range(icfmax - 1, max([icfmax - lpickwindow, 2]), -1):
                 if diffcf[i - 1] >= diffcf[i]:
                     self.Pick = self.Tcf[i]

pylot/core/pick/utils.py

@@ -9,7 +9,6 @@
 """
 
 import warnings
-
 import matplotlib.pyplot as plt
 import numpy as np
 from obspy.core import Stream, UTCDateTime
@@ -405,7 +404,9 @@ def getnoisewin(t, t1, tnoise, tgap):
     inoise, = np.where((t <= max([t1 - tgap, 0])) \
                        & (t >= max([t1 - tnoise - tgap, 0])))
     if np.size(inoise) < 1:
-        print ("getnoisewin: Empty array inoise, check noise window!")
+        inoise, = np.where((t>=t[0]) & (t<=t1))
+        if np.size(inoise) < 1:
+            print ("getnoisewin: Empty array inoise, check noise window!")
 
     return inoise
 
@@ -599,7 +600,7 @@ def wadaticheck(pickdic, dttolerance, iplot):
         wfitflag = 1
 
     # plot results
-    if iplot > 1:
+    if iplot > 0:
         plt.figure()#iplot)
         f1, = plt.plot(Ppicks, SPtimes, 'ro')
         if wfitflag == 0:
@@ -743,11 +744,12 @@ def checkPonsets(pickdic, dttolerance, iplot):
     [xjack, PHI_pseudo, PHI_sub] = jackknife(Ppicks, 'VAR', 1)
     # get pseudo variances smaller than average variances
     # (times safety factor), these picks passed jackknife test
-    ij = np.where(PHI_pseudo <= 2 * xjack)
+    ij = np.where(PHI_pseudo <= 5 * xjack)
     # these picks did not pass jackknife test
-    badjk = np.where(PHI_pseudo > 2 * xjack)
+    badjk = np.where(PHI_pseudo > 5 * xjack)
     badjkstations = np.array(stations)[badjk]
     print ("checkPonsets: %d pick(s) did not pass jackknife test!" % len(badjkstations))
+    print(badjkstations)
 
     # calculate median from these picks
     pmedian = np.median(np.array(Ppicks)[ij])
@@ -782,19 +784,22 @@ def checkPonsets(pickdic, dttolerance, iplot):
 
     checkedonsets = pickdic
 
-    if iplot > 1:
-        p1, = plt.plot(np.arange(0, len(Ppicks)), Ppicks, 'r+', markersize=14)
-        p2, = plt.plot(igood, np.array(Ppicks)[igood], 'g*', markersize=14)
+    if iplot > 0:
+        p1, = plt.plot(np.arange(0, len(Ppicks)), Ppicks, 'ro', markersize=14)
+        if len(badstations) < 1 and len(badjkstations) < 1:
+            p2, = plt.plot(np.arange(0, len(Ppicks)), Ppicks, 'go', markersize=14)
+        else:
+            p2, = plt.plot(igood, np.array(Ppicks)[igood], 'go', markersize=14)
         p3, = plt.plot([0, len(Ppicks) - 1], [pmedian, pmedian], 'g',
                        linewidth=2)
         for i in range(0, len(Ppicks)):
-            plt.text(i, Ppicks[i] + 0.2, stations[i])
+            plt.text(i, Ppicks[i] + 0.01, '{0}'.format(stations[i]))
 
         plt.xlabel('Number of P Picks')
         plt.ylabel('Onset Time [s] from 1.1.1970')
         plt.legend([p1, p2, p3], ['Skipped P Picks', 'Good P Picks', 'Median'],
                    loc='best')
-        plt.title('Check P Onsets')
+        plt.title('Jackknifing and Median Tests on P Onsets')
 
     return checkedonsets
 
@@ -928,8 +933,18 @@ def checkZ4S(X, pick, zfac, checkwin, iplot):
     isignal = getsignalwin(tz, pick, checkwin)
 
     # calculate energy levels
-    zcodalevel = max(absz[isignal])
-    encodalevel = max(absen[isignal])
+    try:
+       zcodalevel = max(absz[isignal])
+    except:
+       ii = np.where(isignal <= len(absz))
+       isignal = isignal[ii]
+       zcodalevel = max(absz[isignal - 1])
+    try:
+       encodalevel = max(absen[isignal])
+    except:
+       ii = np.where(isignal <= len(absen))
+       isignal = isignal[ii]
+       encodalevel = max(absen[isignal - 1])
 
     # calculate threshold
     minsiglevel = encodalevel * zfac