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

autoPyLoT.py

@@ -32,19 +32,17 @@ def autoPyLoT(inputfile):
     .. rubric:: Example
 
     """
-    print '************************************'
+    splash = '''************************************\n
-    print '*********autoPyLoT starting*********'
+                *********autoPyLoT starting*********\n
-    print 'The Python picking and Location Tool'
+                The Python picking and Location Tool\n
-    print '      Version ', _getVersionString(), '2015'
+                Version {version} 2015\n
-    print ' '
+                \n
-    print 'Authors:'
+                Authors:\n
-    print 'S. Wehling-Benatelli'
+                S. Wehling-Benatelli (Ruhr-Universität Bochum)\n
-    print '   Ruhr-Universität Bochum'
+                L. Küperkoch (BESTEC GmbH, Landau i. d. Pfalz)\n
-    print 'L. Küperkoch'
+                K. Olbert (Christian-Albrechts Universität zu Kiel)\n
-    print '   BESTEC GmbH, Landau (Pfalz)'
+                ***********************************'''.format(version=_getVersionString())
-    print 'K. Olbert'
+    print(splash)
-    print '   Christian-Albrechts Universität Kiel'
-    print '************************************'
 
     # reading parameter file
 
@@ -77,7 +75,7 @@ def autoPyLoT(inputfile):
             # get path to NLLoc executable
             nllocbin = parameter.getParam('nllocbin')
             nlloccall = '%s/NLLoc' % nllocbin
-            # get name of phase file  
+            # get name of phase file
             phasef = parameter.getParam('phasefile')
             phasefile = '%s/obs/%s' % (nllocroot, phasef)
             # get name of NLLoc-control file
@@ -99,8 +97,8 @@ def autoPyLoT(inputfile):
         if not parameter.hasParam('eventID'):
             for event in [events for events in glob.glob(os.path.join(datapath, '*')) if os.path.isdir(events)]:
                 data.setWFData(glob.glob(os.path.join(datapath, event, '*')))
-                print 'Working on event %s' % event
+                print('Working on event %s' % event)
-                print data
+                print(data)
 
                 wfdat = data.getWFData()  # all available streams
                 ##########################################################
@@ -113,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)
@@ -133,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
@@ -148,14 +150,14 @@ def autoPyLoT(inputfile):
             ##########################################################
             # !automated picking starts here!
             picks = autopickevent(wfdat, parameter)
-            
+
             ##########################################################
             # locating
             if locflag == 1:
                 # 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)
@@ -166,17 +168,20 @@ def autoPyLoT(inputfile):
                 filedata = nllfile.read()
                 if filedata.find(locfiles) < 0:
                     # replace old command
-            	    filedata = filedata.replace('LOCFILES', locfiles)
+                    filedata = filedata.replace('LOCFILES', locfiles)
-            	    nllfile = open(locfile, 'w')
+                    nllfile = open(locfile, 'w')
-            	    nllfile.write(filedata)
+                    nllfile.write(filedata)
-            	    nllfile.close()
+                    nllfile.close()
 
                 # 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 '------------------------------------------'
             print '-------Finished event %s!-------' % parameter.getParam('eventID')
             print '------------------------------------------'

autoPyLoT_local.in

@@ -9,15 +9,15 @@ EVENT_DATA/LOCAL		    #datapath#     %data path
 2013.02_Insheim                     #database#     %name of data base
 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
+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/loc/__init__.py

@@ -0,0 +1,2 @@
+#!/usr/bin/env python
+# -*- coding: utf-8 -*-

pylot/core/loc/hsat.py

@@ -0,0 +1,2 @@
+#!/usr/bin/env python
+# -*- coding: utf-8 -*-

pylot/core/loc/nll.py

@@ -0,0 +1,60 @@
+#!/usr/bin/env python
+# -*- coding: utf-8 -*-
+
+import subprocess
+from obspy.core.event import readEvents
+from pylot.core.pick.utils import writephases
+
+def picksExport(picks, phasefile):
+    '''
+    Take <picks> dictionary and exports picking data to a NLLOC-obs <phasefile> without creating an ObsPy event object.
+    :param picks: picking data dictionary
+    :type picks: dict
+    :param phasefile: complete path to the exporting obs file
+    :type phasefile: str
+    '''
+    # write phases to NLLoc-phase file
+    writephases(picks, 'NLLoc', phasefile)
+
+def modfiyInputFile(fn, root, outpath, phasefile, tttable):
+    '''
+
+    :param fn:
+    :param root:
+    :param outpath:
+    :param phasefile:
+    :param tttable:
+    :return:
+    '''
+    # For locating the event we have to modify the NLLoc-control file!
+    # create comment line for NLLoc-control file NLLoc-output file
+    print ("Modifying  NLLoc-control file %s ..." % fn)
+    nllocout = '%s/loc/%s_%s' % (root, outpath)
+    locfiles = 'LOCFILES %s NLLOC_OBS %s %s 0' % (phasefile, tttable, nllocout)
+
+    # modification of NLLoc-control file
+    nllfile = open(fn, 'r')
+    filedata = nllfile.read()
+    if filedata.find(locfiles) < 0:
+        # replace old command
+        filedata = filedata.replace('LOCFILES', locfiles)
+        nllfile = open(fn, 'w')
+        nllfile.write(filedata)
+        nllfile.close()
+
+def locate(call, fnin):
+    '''
+    Takes paths to NLLoc executable <call> and input parameter file <fnin> and starts the location calculation.
+    :param call: full path to NLLoc executable
+    :type call: str
+    :param fnin: full path to input parameter file
+    :type fnin: str
+    '''
+    # locate the event
+    subprocess.call([call, fnin])
+
+def readLocation(fn):
+    pass
+
+if __name__=='__main__':
+    pass

pylot/core/loc/velest.py

@@ -0,0 +1,2 @@
+#!/usr/bin/env python
+# -*- coding: utf-8 -*-

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__':