Merge branch 'develop' of ariadne.geophysik.ruhr-uni-bochum.de:/data/git/pylot into develop

autoPyLoT.py

@@ -111,7 +111,7 @@ def autoPyLoT(inputfile):
                     # write phases to NLLoc-phase file
                     writephases(picks, 'NLLoc', phasefile)
 
-                    # For locating the events we have to modify the NLLoc-control file!
+                    # For locating the event the NLLoc-control file has to be modified!
                     # create comment line for NLLoc-control file
                     # NLLoc-output file
                     nllocout = '%s/loc/%s_%s' % (nllocroot, event, nllocoutpatter)
@@ -131,6 +131,10 @@ def autoPyLoT(inputfile):
                     # locate the event
                     subprocess.call([nlloccall, locfile])
                 ##########################################################
+                # write phase files for various location routines
+                # HYPO71
+                hypo71file = '%s/%s/autoPyLoT_HYPO71.pha' % (datapath, eventID)
+                writephases(picks, 'HYPO71', hypo71file)
 
                 print '------------------------------------------'
                 print '-----Finished event %s!-----' % event
@@ -153,7 +157,7 @@ def autoPyLoT(inputfile):
                 # write phases to NLLoc-phase file
                 writephases(picks, 'NLLoc', phasefile)
 
-                # For locating the event we have to modify the NLLoc-control file!
+                # For locating the event the NLLoc-control file has to be modified!
                 # create comment line for NLLoc-control file NLLoc-output file
                 nllocout = '%s/loc/%s_%s' % (nllocroot, parameter.getParam('eventID'), nllocoutpatter)
                 locfiles = 'LOCFILES %s NLLOC_OBS %s %s 0' % (phasefile, ttpatter, nllocout)
@@ -172,7 +176,10 @@ def autoPyLoT(inputfile):
                 # locate the event
                 subprocess.call([nlloccall, locfile])
             ##########################################################
-
+            # write phase files for various location routines
+            # HYPO71
+            hypo71file = '%s/%s/autoPyLoT_HYPO71.pha' % (datapath, parameter.getParam('eventID'))
+            writephases(picks, 'HYPO71', hypo71file)
            
                    
             print '------------------------------------------'

autoPyLoT_local.in

@@ -10,14 +10,14 @@ EVENT_DATA/LOCAL		    #datapath#     %data path
 e0019.048.13                        #eventID#      %event ID for single event processing
 /DATA/Insheim/STAT_INFO             #invdir#       %full path to inventory or dataless-seed file
 PILOT				    #datastructure#%choose data structure
-0                                   #iplot#        %flag for plotting: 0 none, 1, partly, >1 everything
+0                                   #iplot#        %flag for plotting: 0 none, 1 partly, >1 everything
 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
-#NLLoc settings
+#NLLoc settings#
 /home/ludger/NLLOC                  #nllocbin#     %path to NLLoc executable
 /home/ludger/NLLOC/Insheim          #nllocroot#    %root of NLLoc-processing directory
 AUTOPHASES.obs                      #phasefile#    %name of autoPyLoT-output phase file for NLLoc
                                                    %(in nllocroot/obs)
-Insheim_min1d2015.in                #locfile#      %name of autoPyLoT-output control file for NLLoc
+Insheim_min1d2015_auto.in           #locfile#      %name of autoPyLoT-output control file for NLLoc
                                                    %(in nllocroot/run)
 ttime                               #ttpatter#     %pattern of NLLoc ttimes from grid
                                                    %(in nllocroot/times)

autoPyLoT_regional.in

@@ -88,7 +88,7 @@ ARH                                 #algoS#        %choose algorithm for S-onset
 2.5                                 #noisefactor#  %noiselevel*noisefactor=threshold
 60                                  #minpercent#   %required percentage of samples higher than threshold
 #check for spuriously picked S-onsets#
-1.0                                 #zfac#         %P-amplitude must exceed at least zfac times RMS-S amplitude 
+0.5                                 #zfac#         %P-amplitude must exceed at least zfac times RMS-S amplitude 
 #check statistics of P onsets#
 45                                  #mdttolerance#  %maximum allowed deviation of P picks from median [s]
 #wadati check#

pylot/RELEASE-VERSION

@@ -1 +1 @@
-ac7d-dirty
+a31e-dirty

pylot/core/pick/autopick.py

@@ -335,7 +335,8 @@ def autopickstation(wfstream, pickparam):
                                "no zero crossings derived!")
                         print ("Cannot calculate source spectrum!")
                     else:
-                        calcwin = (zc[3] - zc[0]) * z_copy[0].stats.delta
+                        index = min([3, len(zc) - 1])
+                        calcwin = (zc[index] - zc[0]) * z_copy[0].stats.delta
                     # calculate source spectrum and get w0 and fc
                     specpara = DCfc(z_copy, mpickP, calcwin, iplot)
                     w0 = specpara.getw0()

pylot/core/pick/utils.py

@@ -8,7 +8,7 @@
 
    :author: Ludger Kueperkoch / MAGS2 EP3 working group
 """
-
+import pdb
 import numpy as np
 import matplotlib.pyplot as plt
 from obspy.core import Stream, UTCDateTime
@@ -937,7 +937,7 @@ def writephases(arrivals, fformat, filename):
     Function of methods to write phases to the following standard file
     formats used for locating earthquakes:
 
-    HYPO71, NLLoc, VELEST, HYPOSAT, HYPOINVERSE and hypoDD
+    HYPO71, NLLoc, VELEST, HYPOSAT, and hypoDD
 
     :param: arrivals
     :type: dictionary containing all phase information including
@@ -999,11 +999,80 @@ def writephases(arrivals, fformat, filename):
                  mm,
                  ss_ms))
 
-
         fid.close()
 
+    elif fformat == 'HYPO71':
+    	print ("Writing phases to %s for HYPO71" % filename)
+        fid = open("%s" % filename, 'w')
+        # write header
+        fid.write('                                                              EQ001\n')
+        for key in arrivals:
+            if arrivals[key]['P']['weight'] < 4:
+                Ponset = arrivals[key]['P']['mpp']
+                Sonset = arrivals[key]['S']['mpp']
+                pweight = arrivals[key]['P']['weight'] 
+                sweight = arrivals[key]['S']['weight'] 
+                fm = arrivals[key]['P']['fm']
+                if fm is None:
+                    fm = '-'
+                Ao = arrivals[key]['S']['Ao']
+                if Ao is None:
+                    Ao = ''
+                else:
+                    Ao = str('%7.2f' % Ao)
+                year = Ponset.year
+                if year >= 2000:
+                   year = year -2000
+                else:
+                   year = year - 1900
+                month = Ponset.month
+                day = Ponset.day
+                hh = Ponset.hour
+                mm = Ponset.minute
+                ss = Ponset.second
+                ms = Ponset.microsecond
+                ss_ms = ss + ms / 1000000.0
+                if pweight < 2:
+                     pstr = 'I'
+                elif pweight >= 2:
+                     pstr = 'E'
+                if arrivals[key]['S']['weight'] < 4:
+                    Sss = Sonset.second
+                    Sms = Sonset.microsecond
+                    Sss_ms = Sss + Sms / 1000000.0
+                    Sss_ms = str('%5.02f' % Sss_ms)
+                    if sweight < 2:
+                        sstr = 'I'
+                    elif sweight >= 2:
+                        sstr = 'E'
+                    fid.write('%s%sP%s%d %02d%02d%02d%02d%02d%5.2f       %s%sS %d   %s\n' % (key,
+                     pstr,
+                     fm,
+                     pweight,
+                     year,
+                     month,
+                     day,
+                     hh,
+                     mm,
+                     ss_ms,
+                     Sss_ms,
+                     sstr,
+                     sweight,
+                     Ao))
+                else:
+                    fid.write('%s%sP%s%d %02d%02d%02d%02d%02d%5.2f                  %s\n' % (key,
+                     pstr,
+                     fm,
+                     pweight,
+                     year,
+                     month,
+                     day,
+                     hh,
+                     mm,
+                     ss_ms,
+                     Ao))
 
-
+        fid.close()
         
 
 if __name__ == '__main__':