[update] adding possibility to display other waveform data (e.g. denoised/synthetic) together with genuine data for comparison

closes #24

PyLoT.py

@@ -60,7 +60,7 @@ except ImportError:
     from matplotlib.backends.backend_qt5agg import FigureCanvasQTAgg as FigureCanvas
 from matplotlib.figure import Figure
 
-from pylot.core.analysis.magnitude import LocalMagnitude, MomentMagnitude
+from pylot.core.analysis.magnitude import LocalMagnitude, MomentMagnitude, calcsourcespec
 from pylot.core.io.data import Data
 from pylot.core.io.inputs import FilterOptions, PylotParameter
 from autoPyLoT import autoPyLoT
@@ -76,20 +76,21 @@ from pylot.core.util.utils import fnConstructor, getLogin, \
     full_range, readFilterInformation, pick_color_plt, \
     pick_linestyle_plt, identifyPhaseID, excludeQualityClasses, \
     transform_colors_mpl, transform_colors_mpl_str, getAutoFilteroptions, check_all_obspy, \
-    check_all_pylot, get_Bool, get_None
+    check_all_pylot, get_bool, get_None
 from pylot.core.util.gui import make_pen
 from pylot.core.util.event import Event
 from pylot.core.io.location import create_creation_info, create_event
 from pylot.core.util.widgets import FilterOptionsDialog, NewEventDlg, \
     PylotCanvas, WaveformWidgetPG, PropertiesDlg, HelpForm, createAction, PickDlg, \
     ComparisonWidget, TuneAutopicker, PylotParaBox, AutoPickDlg, CanvasWidget, AutoPickWidget, \
-    CompareEventsWidget, ProgressBarWidget, AddMetadataWidget, SingleTextLineDialog, LogWidget
+    CompareEventsWidget, ProgressBarWidget, AddMetadataWidget, SingleTextLineDialog, LogWidget, PickQualitiesFromXml, \
+    SourceSpecWindow, ChooseWaveFormWindow, SpectrogramTab
 from pylot.core.util.array_map import Array_map
 from pylot.core.util.structure import DATASTRUCTURE
 from pylot.core.util.thread import Thread, Worker
 from pylot.core.util.version import get_git_version as _getVersionString
 from pylot.core.io.getEventListFromXML import geteventlistfromxml
-from pylot.core.io.getQualitiesfromxml import getQualitiesfromxml
+from pylot.core.io.phases import getQualitiesfromxml
 
 from pylot.styles import style_settings
 
@@ -116,17 +117,19 @@ class MainWindow(QMainWindow):
         if not infile:
             infile = os.path.join(os.path.expanduser('~'), '.pylot', 'pylot.in')
             print('Using default input file {}'.format(infile))
-        if os.path.isfile(infile) == False:
+        if os.path.isfile(infile) is False:
-            infile = QFileDialog().getOpenFileName(caption='Choose PyLoT-input file')
+            infile = QFileDialog().getOpenFileName(caption='Choose PyLoT-input file')[0]
 
-            if not os.path.exists(infile[0]):
+            if not os.path.exists(infile):
                 QMessageBox.warning(self, "PyLoT Warning",
-                                    "No PyLoT-input file declared!")
+                                    "No PyLoT-input file declared! Using default parameters!")
-                sys.exit(0)
+                infile = None
-            self.infile = infile[0]
+
-        else:
-            self.infile = infile
         self._inputs = PylotParameter(infile)
+        if not infile:
+            self._inputs.reset_defaults()
+
+        self.infile = infile
         self._props = None
 
         self.gain = 1.
@@ -178,6 +181,7 @@ class MainWindow(QMainWindow):
         self.autodata = Data(self)
 
         self.fnames = None
+        self.fnames_comp = None
         self._stime = None
 
         # track deleted picks for logging
@@ -254,7 +258,7 @@ class MainWindow(QMainWindow):
             self._inputs.export2File(infile)
         self.infile = infile
 
-    def setupUi(self, use_logwidget=True):
+    def setupUi(self, use_logwidget=False):
         try:
             self.startTime = min(
                 [tr.stats.starttime for tr in self.data.wfdata])
@@ -480,7 +484,7 @@ class MainWindow(QMainWindow):
                                                               "automatic pick "
                                                               "data.", False)
         self.compare_action.setEnabled(False)
-        self.qualities_action = self.createAction(parent=self, text='Show pick qualitites...',
+        self.qualities_action = self.createAction(parent=self, text='Show pick qualities...',
                                                   slot=self.pickQualities, shortcut='Alt+Q',
                                                   icon=qualities_icon, tip='Histogram of pick qualities')
         self.qualities_action.setEnabled(False)
@@ -493,7 +497,6 @@ class MainWindow(QMainWindow):
                                                       icon=eventlist_xml_icon,
                                                       tip='Create an Eventlist from a XML File')
         self.eventlist_xml_action.setEnabled(False)
-
         printAction = self.createAction(self, "&Print event ...",
                                         self.show_event_information, QKeySequence.Print,
                                         print_icon,
@@ -506,6 +509,8 @@ class MainWindow(QMainWindow):
         logAction = self.createAction(self, "&Show Log", self.showLogWidget,
                                       tip="""Display Log""")
 
+        logAction.setEnabled(use_logwidget)
+
         # create button group for component selection
 
         componentGroup = QActionGroup(self)
@@ -695,14 +700,17 @@ class MainWindow(QMainWindow):
         wf_tab = QtWidgets.QWidget(self)
         array_tab = QtWidgets.QWidget(self)
         events_tab = QtWidgets.QWidget(self)
+        spectro_tab = QtWidgets.QWidget(self)
 
         # init main widgets layouts
         self.wf_layout = QtWidgets.QVBoxLayout()
         self.array_layout = QtWidgets.QVBoxLayout()
         self.events_layout = QtWidgets.QVBoxLayout()
+        self.spectro_layout = QtWidgets.QVBoxLayout()
         wf_tab.setLayout(self.wf_layout)
         array_tab.setLayout(self.array_layout)
         events_tab.setLayout(self.events_layout)
+        spectro_tab.setLayout(self.spectro_layout)
 
         # tighten up layouts inside tabs
         for layout in [self.wf_layout, self.array_layout, self.events_layout]:
@@ -713,12 +721,14 @@ class MainWindow(QMainWindow):
         self.tabs.addTab(wf_tab, 'Waveform Plot')
         self.tabs.addTab(array_tab, 'Array Map')
         self.tabs.addTab(events_tab, 'Eventlist')
+        self.tabs.addTab(spectro_tab, 'Spectro')
 
         self.wf_layout.addWidget(self.no_data_label)
         self.wf_layout.addWidget(self.wf_scroll_area)
         self.wf_scroll_area.setWidgetResizable(True)
         self.init_array_tab()
         self.init_event_table()
+        self.init_spectro_tab()
         self.tabs.setCurrentIndex(0)
 
         self.eventLabel = QLabel()
@@ -731,17 +741,19 @@ class MainWindow(QMainWindow):
         _widget.setLayout(self._main_layout)
         _widget.showFullScreen()
 
+
         if use_logwidget:
             self.logwidget = LogWidget(parent=None)
             self.logwidget.show()
             self.stdout = self.logwidget.stdout
             self.stderr = self.logwidget.stderr
 
-            sys.stdout = self.stdout
-            sys.stderr = self.stderr
-
         self.setCentralWidget(_widget)
 
+        # Need to store PickQualities Window somewhere so it doesnt disappear
+        self.pickQualitiesWindow = None
+
+
     def init_wfWidget(self):
         xlab = self.startTime.strftime('seconds since %Y/%m/%d %H:%M:%S (%Z)')
         plottitle = None  # "Overview: {0} components ".format(self.getComponent())
@@ -1119,16 +1131,19 @@ class MainWindow(QMainWindow):
             else:
                 return
 
-    def getWFFnames_from_eventbox(self, eventbox=None):
+    def getWFFnames_from_eventbox(self, eventbox: str = None, subpath: str = None) -> list:
         '''
         Return waveform filenames from event in eventbox.
         '''
         # TODO: add dataStructure class for obspyDMT here, this is just a workaround!
         eventpath = self.get_current_event_path(eventbox)
-        basepath = eventpath.split(os.path.basename(eventpath))[0]
+        if subpath:
+            eventpath = os.path.join(eventpath, subpath)
+        if not os.path.isdir(eventpath):
+            return []
         if self.dataStructure:
             if not eventpath:
-                return
+                return []
             fnames = [os.path.join(eventpath, f) for f in os.listdir(eventpath)]
         else:
             raise DatastructureError('not specified')
@@ -1172,7 +1187,7 @@ class MainWindow(QMainWindow):
         '''
         if not self.project:
             self.createNewProject()
-        ed = getExistingDirectories(self, 'Select event directories...')
+        ed = GetExistingDirectories(self, 'Select event directories...')
         if ed.exec_():
             eventlist = [event for event in ed.selectedFiles() if not event.endswith('EVENTS-INFO')]
             basepath = eventlist[0].split(os.path.basename(eventlist[0]))[0]
@@ -1373,7 +1388,7 @@ class MainWindow(QMainWindow):
         index = eventBox.currentIndex()
         tv = QtWidgets.QTableView()
         header = tv.horizontalHeader()
-        header.setResizeMode(QtWidgets.QHeaderView.ResizeToContents)
+        header.setSectionResizeMode(QtWidgets.QHeaderView.ResizeToContents)
         header.setStretchLastSection(True)
         header.hide()
         tv.verticalHeader().hide()
@@ -1669,16 +1684,55 @@ class MainWindow(QMainWindow):
         self.cmpw.show()
 
     def pickQualities(self):
-        path = self._inputs['rootpath'] + '/' + self._inputs['datapath'] + '/' + self._inputs['database']
+        path = self.get_current_event_path()
-        getQualitiesfromxml(path)
+        (_, plot) = getQualitiesfromxml(path, self._inputs.get('timeerrorsP'), self._inputs.get('timeerrorsS'),plotflag=1)
+        self.pickQualitiesWindow = PickQualitiesFromXml(figure=plot, path=self.get_current_event_path(),inputVar=self._inputs)
+        self.pickQualitiesWindow.showUI()
         return 
 
+    # WIP JG
     def eventlistXml(self):
         path = self._inputs['rootpath'] + '/' + self._inputs['datapath'] + '/' + self._inputs['database']
         outpath = self.project.location[:self.project.location.rfind('/')]
         geteventlistfromxml(path, outpath)
         return
 
+        # WIP JG
+    def spectogramView(self):
+        global test
+        stations = []
+        names = []
+        traces = {}
+        for tr in self.get_data().wfdata.traces:
+            if not tr.stats.station in stations:
+                stations.append(tr.stats.station)
+                names.append(tr.stats.network + '.' + tr.stats.station)
+        for station in stations:
+            traces[station] = {}
+        for ch in ['Z', 'N', 'E']:
+            for tr in self.get_data().wfdata.select(component=ch).traces:
+                traces[tr.stats.station][ch] = tr
+
+        names.sort()
+        a = self.get_current_event()
+
+        print (self.get_data().wfdata.traces[0])
+
+        test = SpectrogramTab(traces, self.get_data().wfdata)
+        height = self.tabs.widget(0).height()
+        width = self.tabs.widget(0).width()
+        self.tabs.setCurrentIndex(3)
+        figCanvas = test.makeSpecFig(direction=self.dispComponent, height = height, width = width, parent = self.tabs.widget)
+        return figCanvas
+        #self.spectro_layout.addWidget()
+        # self.get_data().wfdata.spectrogram()
+        # self.tabs.addTab(figCanvas, 'Spectrogram')
+        # self.tabs[3] = figCanvas
+        # self.refreshTabs()
+        # test.show()
+
+
+
     def compareMulti(self):
         if not self.compareoptions:
             return
@@ -1864,6 +1918,13 @@ class MainWindow(QMainWindow):
                     self.newWF(plot=False)
                 self.update_obspy_dmt()
                 self.refresh_array_map()
+        if self.tabs.currentIndex() == 3:
+            if self.spectroWidget != None:
+                self.spectro_layout.removeWidget(self.spectroWidget)
+            newSpectroWidget = self.spectogramView()
+            self.spectro_layout.addWidget(newSpectroWidget)
+            self.spectroWidget = newSpectroWidget
+
 
     def newWF(self, event=None, plot=True):
         '''
@@ -1894,13 +1955,20 @@ class MainWindow(QMainWindow):
 
     def prepareLoadWaveformData(self):
         self.fnames = self.getWFFnames_from_eventbox()
-        self.fnames_syn = []
+        self.fnames_comp = []
+        fnames_comp = self.getWFFnames_from_eventbox(subpath='compare')
+        self.dataPlot.activateCompareOptions(bool(fnames_comp))
+        if fnames_comp:
+            if self.dataPlot.comp_checkbox.isChecked():
+                self.fnames_comp = fnames_comp
+
         eventpath = self.get_current_event_path()
         basepath = eventpath.split(os.path.basename(eventpath))[0]
         self.obspy_dmt = check_obspydmt_structure(basepath)
         self.dataPlot.activateObspyDMToptions(self.obspy_dmt)
         if self.obspy_dmt:
             self.prepareObspyDMT_data(eventpath)
+            self.dataPlot.activateCompareOptions(True)
 
     def loadWaveformData(self):
         '''
@@ -1931,7 +1999,7 @@ class MainWindow(QMainWindow):
             tstop = None
 
         self.data.setWFData(self.fnames,
-                            self.fnames_syn,
+                            self.fnames_comp,
                             checkRotated=True,
                             metadata=self.metadata,
                             tstart=tstart,
@@ -1939,7 +2007,7 @@ class MainWindow(QMainWindow):
 
     def prepareObspyDMT_data(self, eventpath):
         qcbox_processed = self.dataPlot.qcombo_processed
-        qcheckb_syn = self.dataPlot.syn_checkbox
+        qcheckb_syn = self.dataPlot.comp_checkbox
         qcbox_processed.setEnabled(False)
         qcheckb_syn.setEnabled(False)
         for fpath in os.listdir(eventpath):
@@ -1947,8 +2015,8 @@ class MainWindow(QMainWindow):
             if 'syngine' in fpath:
                 eventpath_syn = os.path.join(eventpath, fpath)
                 qcheckb_syn.setEnabled(True)
-                if self.dataPlot.syn_checkbox.isChecked():
+                if self.dataPlot.comp_checkbox.isChecked():
-                    self.fnames_syn = [os.path.join(eventpath_syn, filename) for filename in os.listdir(eventpath_syn)]
+                    self.fnames_comp = [os.path.join(eventpath_syn, filename) for filename in os.listdir(eventpath_syn)]
             if 'processed' in fpath:
                 qcbox_processed.setEnabled(True)
         if qcbox_processed.isEnabled():
@@ -2121,6 +2189,7 @@ class MainWindow(QMainWindow):
             self.locateEventAction.setEnabled(True)
             self.qualities_action.setEnabled(True)
             self.eventlist_xml_action.setEnabled(True)
+
         if True in self.comparable.values():
             self.compare_action.setEnabled(True)
         self.draw()
@@ -2228,7 +2297,7 @@ class MainWindow(QMainWindow):
         comp = self.getComponent()
         title = 'section: {0} components'.format(zne_text[comp])
         wfst = self.get_data().getWFData()
-        wfsyn = self.get_data().getSynWFData()
+        wfsyn = self.get_data().getAltWFdata()
         if self.filterActionP.isChecked() and filter:
             self.filterWaveformData(plot=False, phase='P')
         elif self.filterActionS.isChecked() and filter:
@@ -2237,7 +2306,7 @@ class MainWindow(QMainWindow):
         # wfst += self.get_data().getWFData().select(component=alter_comp)
         plotWidget = self.getPlotWidget()
         self.adjustPlotHeight()
-        if get_Bool(settings.value('large_dataset')) == True:
+        if get_bool(settings.value('large_dataset')):
             self.plot_method = 'fast'
         else:
             self.plot_method = 'normal'
@@ -2540,17 +2609,21 @@ class MainWindow(QMainWindow):
             print("Warning! No network, station, and location info available!")
             return
         self.update_status('picking on station {0}'.format(station))
-        data = self.get_data().getOriginalWFData().copy()
+        wfdata = self.get_data().getOriginalWFData().copy()
+        wfdata_comp = self.get_data().getAltWFdata().copy()
         event = self.get_current_event()
         wftype = self.dataPlot.qcombo_processed.currentText() if self.obspy_dmt else None
         pickDlg = PickDlg(self, parameter=self._inputs,
-                          data=data.select(station=station),
+                          data=wfdata.select(station=station),
+                          data_compare=wfdata_comp.select(station=station),
                           station=station, network=network,
                           location=location,
                           picks=self.getPicksOnStation(station, 'manual'),
                           autopicks=self.getPicksOnStation(station, 'auto'),
                           metadata=self.metadata, event=event,
-                          filteroptions=self.filteroptions, wftype=wftype)
+                          model=self.inputs.get('taup_model'),
+                          filteroptions=self.filteroptions, wftype=wftype,
+                          show_comp_data=self.dataPlot.comp_checkbox.isChecked())
         if self.filterActionP.isChecked():
             pickDlg.currentPhase = "P"
             pickDlg.filterWFData()
@@ -2663,7 +2736,7 @@ class MainWindow(QMainWindow):
         self.init_fig_dict()
         # if not self.tap:
         # init TuneAutopicker object
-        self.tap = TuneAutopicker(self)
+        self.tap = TuneAutopicker(self, self.obspy_dmt)
         # first call of update to init tabs with empty canvas
         self.update_autopicker()
         # connect update signal of TuneAutopicker with update function
@@ -2862,7 +2935,9 @@ class MainWindow(QMainWindow):
         self.log_deleted_picks([deleted_pick])
 
     def log_deleted_picks(self, deleted_picks, event_path=None):
-        ''' Log deleted picks to list self.deleted_picks '''
+        '''
+        Log deleted picks to list self.deleted_picks
+        '''
         if not event_path:
             event_path = self.get_current_event_path()
         for deleted_pick in deleted_picks:
@@ -2876,7 +2951,9 @@ class MainWindow(QMainWindow):
             self.deleted_picks[event_path].append(deleted_pick)
 
     def dump_deleted_picks(self, event_path):
-        ''' Save deleted picks to json file for event in event_path. Load old file before and merge'''
+        '''
+        Save deleted picks to json file for event in event_path. Load old file before and merge
+        '''
         try:
             deleted_picks_from_file = self.load_deleted_picks(event_path)
         except Exception as e:
@@ -3072,11 +3149,13 @@ class MainWindow(QMainWindow):
         lt = locateTool[loctool]
         # get working directory
         locroot = parameter['nllocroot']
+        #locroot = 'E:/NLL/src/Insheim'
         if locroot is None:
             self.PyLoTprefs()
             self.locate_event()
 
         ctrfile = os.path.join(locroot, 'run', parameter['ctrfile'])
+        #ctrfile = 'E:/NLL/src/Insheim/run/Insheim_min1d032016.in'
         ttt = parameter['ttpatter']
         outfile = parameter['outpatter']
         eventname = self.get_current_event_name()
@@ -3087,7 +3166,7 @@ class MainWindow(QMainWindow):
         phasefile = os.path.join(obsdir, filename + '.obs')
         lt.modify_inputs(ctrfile, locroot, filename, phasefile, ttt)
         try:
-            lt.locate(ctrfile)
+            lt.locate(ctrfile, self._inputs)
         except RuntimeError as e:
             print(e.message)
         # finally:
@@ -3141,7 +3220,7 @@ class MainWindow(QMainWindow):
         '''
         if checked: pass  # dummy argument for QAction trigger signal
         self.tabs.setCurrentIndex(1)
-        # if there is no metadata (invetories is an empty list), just initialize the default empty tab
+        # if there is no metadata (inventories is an empty list), just initialize the default empty tab
         if not self.metadata.inventories:
             self.init_array_tab()
             return
@@ -3163,6 +3242,15 @@ class MainWindow(QMainWindow):
         self.tabs.setCurrentIndex(index)
         self.refresh_array_map()
 
+    def init_spectro_tab(self):
+        '''
+        Init spectrogram tab with currently selected event.
+        '''
+        self.spectroWidget = None
+        #self.spectro_layout.addWidget( self.spectogramView() )
+        pass
+
+
     def array_map_thread(self):
         '''
         Start modal thread to init the array_map object.
@@ -3385,7 +3473,7 @@ class MainWindow(QMainWindow):
                     self.event_table.setCellWidget(r_index, c_index, item)
 
         header = self.event_table.horizontalHeader()
-        header.setResizeMode(QtWidgets.QHeaderView.ResizeToContents)
+        header.setSectionResizeMode(QtWidgets.QHeaderView.ResizeToContents)
         header.setStretchLastSection(True)
         self.event_table.cellChanged[int, int].connect(cell_changed)
         self.event_table.cellClicked[int, int].connect(cell_clicked)
@@ -3492,7 +3580,7 @@ class MainWindow(QMainWindow):
 
         wf_select = Stream()
         # restitute only picked traces
-        for station in np.unique(self.getPicks('manual').keys() + self.getPicks('auto').keys()):
+        for station in np.unique(list(self.getPicks('manual').keys()) + list(self.getPicks('auto').keys())):
             wf_select += wf_copy.select(station=station)
 
         corr_wf = restitute_data(wf_select, self.metadata)
@@ -3671,7 +3759,7 @@ class MainWindow(QMainWindow):
             filename = fnm[0] + '.plp'
         self.project.parameter = self._inputs
         settings = QSettings()
-        autosaveXML = get_Bool(settings.value('autosaveXML', True))
+        autosaveXML = get_bool(settings.value('autosaveXML', True))
         if autosaveXML:
             self.exportEvents()
         if not self.project.save(filename): return False
@@ -3695,7 +3783,7 @@ class MainWindow(QMainWindow):
                 self.metadata.clear_inventory()
                 self.project.parameter = self._inputs
                 settings = QSettings()
-                autosaveXML = get_Bool(settings.value('autosaveXML', True))
+                autosaveXML = get_bool(settings.value('autosaveXML', True))
                 if autosaveXML:
                     self.exportEvents()
                 if not self.project.save(): return False
@@ -3784,7 +3872,7 @@ class MainWindow(QMainWindow):
     def helpHelp(self):
         if checkurl():
             form = HelpForm(self,
-                            'https://ariadne.geophysik.ruhr-uni-bochum.de/trac/PyLoT/wiki')
+                            'https://github.com/seismology-RUB/PyLoT')
         else:
             form = HelpForm(self, ':/help.html')
         form.show()
@@ -3964,13 +4052,13 @@ class Project(object):
         return project
 
 
-class getExistingDirectories(QFileDialog):
+class GetExistingDirectories(QFileDialog):
     '''
     File dialog with possibility to select multiple folders.
     '''
 
     def __init__(self, *args):
-        super(getExistingDirectories, self).__init__(*args)
+        super(GetExistingDirectories, self).__init__(*args)
         self.setOption(self.DontUseNativeDialog, True)
         self.setOption(self.ReadOnly, True)
         self.setFileMode(self.Directory)
@@ -3996,16 +4084,7 @@ def create_window():
     return app, app_created
 
 
-def main(args=None):
+def main(project_filename=None, pylot_infile=None, reset_qsettings=False):
-    project_filename = None
-    # args.project_filename = 'C:/Shared/AlpArray/alparray_data/project_alparray_test.plp'
-    pylot_infile = None
-    if args:
-        if args.project_filename:
-            project_filename = args.project_filename
-        if args.input_filename:
-            pylot_infile = args.input_filename
-        reset_qsettings = args.reset_qsettings
 
     # create the Qt application
     pylot_app, app_created = create_window()
@@ -4054,4 +4133,5 @@ if __name__ == "__main__":
     parser.add_argument('--reset_qsettings', default=False, action='store_true',
                         help='reset qsettings (debug option)')
     args = parser.parse_args()
-    sys.exit(main(args))
+    sys.exit(main(project_filename=args.project_filename, pylot_infile=args.input_filename,
+                  reset_qsettings=args.reset_qsettings))

README.md

@@ -24,10 +24,10 @@ It is highly recommended to use Anaconda for a simple creation of a Python insta
 
 Afterwards run (from the PyLoT main directory where the files *requirements.txt* and *pylot.yml* are located)
 
-    conda create --name pylot_38 --file requirements.txt
+    conda env create -f pylot.yml
 or
     
-    conda env create -f pylot.yml
+    conda create --name pylot_38 --file requirements.txt
 
 to create a new Anaconda environment called "pylot_38".
 
@@ -99,13 +99,13 @@ We hope to solve these with the next release.
 
 ## Staff
 
-Original author(s): L. Kueperkoch, S. Wehling-Benatelli, M. Bischoff (PILOT)
+Original author(s): M. Rische, S. Wehling-Benatelli, L. Kueperkoch, M. Bischoff (PILOT)
 
-Developer(s): S. Wehling-Benatelli, L. Kueperkoch, K. Olbert, M. Bischoff, C. Wollin, M. Rische, M. Paffrath
+Developer(s): S. Wehling-Benatelli, M. Paffrath, L. Kueperkoch, K. Olbert, M. Bischoff, C. Wollin, M. Rische, D. Arnold, K. Cökerim, S. Zimmermann
 
 Others: A. Bruestle, T. Meier, W. Friederich
 
 
 [ObsPy]: http://github.com/obspy/obspy/wiki
 
-September 2017
+April 2022

autopylot.sh

@@ -7,4 +7,6 @@
 #$ -l h_vmem=2G
 #$ -l os=*stretch
 
-python ./autoPyLoT.py -i /home/marcel/.pylot/pylot_alparray_mantle_corr_stack_0.03-0.5.in -dmt processed -c $NSLOTS
+conda activate pylot_38
+
+python ./autoPyLoT.py -i /home/marcel/.pylot/pylot_janis_noisy.in -c $NSLOTS

pylot.yml

@@ -8,7 +8,7 @@ dependencies:
   - numpy=1.22.3
   - obspy=1.3.0
   - pyqtgraph=0.12.4
-  - pyside2=5.13.2
+  - pyside2>=5.13.2
   - python=3.8.12
-  - qt=5.12.9
+  - qt>=5.12.9
   - scipy=1.8.0

pylot/core/analysis/magnitude.py

@@ -144,6 +144,10 @@ class Magnitude(object):
                     azimuthal_gap=self.origin_id.get_referred_object().quality.azimuthal_gap)
             else:
                 # no scaling necessary
+                # Temporary fix needs rework
+                if (len(self.magnitudes.keys()) == 0):
+                    print("Error in local magnitude calculation ")
+                    return None
                 mag = ope.Magnitude(
                     mag=np.median([M.mag for M in self.magnitudes.values()]),
                     magnitude_type=self.type,
@@ -414,6 +418,10 @@ class MomentMagnitude(Magnitude):
             distance = degrees2kilometers(a.distance)
             azimuth = a.azimuth
             incidence = a.takeoff_angle
+            if not 0. <= incidence <= 360.:
+                if self.verbose:
+                    print(f'WARNING: Incidence angle outside bounds - {incidence}')
+                return
             w0, fc = calcsourcespec(scopy, onset, self.p_velocity, distance,
                                     azimuth, incidence, self.p_attenuation,
                                     self.plot_flag, self.verbose)
@@ -432,6 +440,40 @@ class MomentMagnitude(Magnitude):
             self.event.station_magnitudes.append(magnitude)
             self.magnitudes = (station, magnitude)
 
+    # WIP JG
+    def getSourceSpec(self):
+        for a in self.arrivals:
+            if a.phase not in 'pP':
+                continue
+            # make sure calculating Mo only from reliable onsets
+            # NLLoc: time_weight = 0 => do not use onset!
+            if a.time_weight == 0:
+                continue
+            pick = a.pick_id.get_referred_object()
+            station = pick.waveform_id.station_code
+            if len(self.stream) <= 2:
+                print("Station:" '{0}'.format(station))
+                print("WARNING: No instrument corrected data available,"
+                      " no magnitude calculation possible! Go on.")
+                continue
+            wf = self.stream.select(station=station)
+            if not wf:
+                continue
+            try:
+                scopy = wf.copy()
+            except AssertionError:
+                print("WARNING: Something's wrong with the data,"
+                      "station {},"
+                      "no calculation of moment magnitude possible! Go on.".format(station))
+                continue
+            onset = pick.time
+            distance = degrees2kilometers(a.distance)
+            azimuth = a.azimuth
+            incidence = a.takeoff_angle
+            w0, fc, plt = calcsourcespec(scopy, onset, self.p_velocity, distance,
+                                    azimuth, incidence, self.p_attenuation,
+                                    3, self.verbose)
+            return w0, fc, plt
 
 def calcMoMw(wfstream, w0, rho, vp, delta, verbosity=False):
     '''
@@ -602,7 +644,7 @@ def calcsourcespec(wfstream, onset, vp, delta, azimuth, incidence,
         N = min(int(np.power(m, 2)), 16384)
         # N = int(np.power(m, 2))
         y = dt * np.fft.fft(xdat, N)
-        Y = abs(y[: N / 2])
+        Y = abs(y[: int(N / 2)])
         L = (N - 1) / freq
         f = np.arange(0, fny, 1 / L)
 
@@ -679,6 +721,8 @@ def calcsourcespec(wfstream, onset, vp, delta, azimuth, incidence,
             plt.xlabel('Frequency [Hz]')
             plt.ylabel('Amplitude [m/Hz]')
             plt.grid()
+            if iplot == 3:
+                return w0, Fc, plt
             plt.show()
             try:
                 input()

pylot/core/io/data.py

@@ -260,7 +260,6 @@ class Data(object):
         can be a str or a list of strings of ['manual', 'auto', 'origin', 'magnitude']
         """
         from pylot.core.util.defaults import OUTPUTFORMATS
-
         if not type(fcheck) == list:
             fcheck = [fcheck]
 
@@ -321,35 +320,61 @@ class Data(object):
             if lendiff != 0:
                 print("Manual as well as automatic picks available. Prefered the {} manual ones!".format(lendiff))
 
+
+            no_uncertainties_p = []
+            no_uncertainties_s = []
             if upperErrors:
                 # check for pick uncertainties exceeding adjusted upper errors
                 # Picks with larger uncertainties will not be saved in output file!
                 for j in range(len(picks)):
                     for i in range(len(picks_copy)):
                         if picks_copy[i].phase_hint[0] == 'P':
-                            if (picks_copy[i].time_errors['upper_uncertainty'] >= upperErrors[0]) or \
+                            # Skipping pick if no upper_uncertainty is found and warning user
-                                    (picks_copy[i].time_errors['uncertainty'] is None):
+                            if picks_copy[i].time_errors['upper_uncertainty'] is None:
+                                #print("{1} P-Pick of station {0} does not have upper_uncertainty and cant be checked".format(
+                                 #   picks_copy[i].waveform_id.station_code,
+                                  #  picks_copy[i].method_id))
+                                if not picks_copy[i].waveform_id.station_code in no_uncertainties_p:
+                                    no_uncertainties_p.append(picks_copy[i].waveform_id.station_code)
+                                continue
+
+                            #print ("checking for upper_uncertainty")
+                            if  (picks_copy[i].time_errors['uncertainty'] is None) or \
+                                    (picks_copy[i].time_errors['upper_uncertainty'] >= upperErrors[0]):
                                 print("Uncertainty exceeds or equal adjusted upper time error!")
                                 print("Adjusted uncertainty: {}".format(upperErrors[0]))
                                 print("Pick uncertainty: {}".format(picks_copy[i].time_errors['uncertainty']))
                                 print("{1} P-Pick of station {0} will not be saved in outputfile".format(
                                     picks_copy[i].waveform_id.station_code,
                                     picks_copy[i].method_id))
-                                print("#")
                                 del picks_copy[i]
                                 break
                         if picks_copy[i].phase_hint[0] == 'S':
-                            if (picks_copy[i].time_errors['upper_uncertainty'] >= upperErrors[1]) or \
+
-                                    (picks_copy[i].time_errors['uncertainty'] is None):
+                            # Skipping pick if no upper_uncertainty is found and warning user
+                            if picks_copy[i].time_errors['upper_uncertainty'] is None:
+                                #print("{1} S-Pick of station {0} does not have upper_uncertainty and cant be checked".format(
+                                    #picks_copy[i].waveform_id.station_code,
+                                    #picks_copy[i].method_id))
+                                if not picks_copy[i].waveform_id.station_code in no_uncertainties_s:
+                                    no_uncertainties_s.append(picks_copy[i].waveform_id.station_code)
+                                continue
+
+
+                            if (picks_copy[i].time_errors['uncertainty'] is None) or \
+                                (picks_copy[i].time_errors['upper_uncertainty'] >= upperErrors[1]):
                                 print("Uncertainty exceeds or equal adjusted upper time error!")
                                 print("Adjusted uncertainty: {}".format(upperErrors[1]))
                                 print("Pick uncertainty: {}".format(picks_copy[i].time_errors['uncertainty']))
                                 print("{1} S-Pick of station {0} will not be saved in outputfile".format(
                                     picks_copy[i].waveform_id.station_code,
                                     picks_copy[i].method_id))
-                                print("#")
                                 del picks_copy[i]
                                 break
+            for s in no_uncertainties_p:
+                print("P-Pick of station {0} does not have upper_uncertainty and cant be checked".format(s))
+            for s in no_uncertainties_s:
+                print("S-Pick of station {0} does not have upper_uncertainty and cant be checked".format(s))
 
             if fnext == '.obs':
                 try:
@@ -426,20 +451,25 @@ class Data(object):
         data.filter(**kwargs)
         self.dirty = True
 
-    def setWFData(self, fnames, fnames_syn=None, checkRotated=False, metadata=None, tstart=0, tstop=0):
+    def setWFData(self, fnames, fnames_alt=None, checkRotated=False, metadata=None, tstart=0, tstop=0):
         """
         Clear current waveform data and set given waveform data
         :param fnames: waveform data names to append
+        :param fnames_alt: alternative data to show (e.g. synthetic/processed)
         :type fnames: list
         """
         self.wfdata = Stream()
         self.wforiginal = None
-        self.wfsyn = Stream()
+        self.wf_alt = Stream()
         if tstart == tstop:
             tstart = tstop = None
         self.tstart = tstart
         self.tstop = tstop
 
+        # remove directories
+        fnames = [fname for fname in fnames if not os.path.isdir(fname)]
+        fnames_alt = [fname for fname in fnames_alt if not os.path.isdir(fname)]
+
         # if obspy_dmt:
         #     wfdir = 'raw'
         #     self.processed = False
@@ -457,8 +487,8 @@ class Data(object):
         #     wffnames = fnames
         if fnames is not None:
             self.appendWFData(fnames)
-            if fnames_syn is not None:
+            if fnames_alt is not None:
-                self.appendWFData(fnames_syn, synthetic=True)
+                self.appendWFData(fnames_alt, alternative=True)
         else:
             return False
 
@@ -478,7 +508,7 @@ class Data(object):
         self.dirty = False
         return True
 
-    def appendWFData(self, fnames, synthetic=False):
+    def appendWFData(self, fnames, alternative=False):
         """
         Read waveform data from fnames and append it to current wf data
         :param fnames: waveform data to append
@@ -491,20 +521,20 @@ class Data(object):
         if self.dirty:
             self.resetWFData()
 
-        real_or_syn_data = {True: self.wfsyn,
+        orig_or_alternative_data = {True: self.wf_alt,
-                            False: self.wfdata}
+                                    False: self.wfdata}
 
         warnmsg = ''
         for fname in set(fnames):
             try:
-                real_or_syn_data[synthetic] += read(fname, starttime=self.tstart, endtime=self.tstop)
+                orig_or_alternative_data[alternative] += read(fname, starttime=self.tstart, endtime=self.tstop)
             except TypeError:
                 try:
-                    real_or_syn_data[synthetic] += read(fname, format='GSE2', starttime=self.tstart, endtime=self.tstop)
+                    orig_or_alternative_data[alternative] += read(fname, format='GSE2', starttime=self.tstart, endtime=self.tstop)
                 except Exception as e:
                     try:
-                        real_or_syn_data[synthetic] += read(fname, format='SEGY', starttime=self.tstart,
+                        orig_or_alternative_data[alternative] += read(fname, format='SEGY', starttime=self.tstart,
-                                                            endtime=self.tstop)
+                                                                      endtime=self.tstop)
                     except Exception as e:
                         warnmsg += '{0}\n{1}\n'.format(fname, e)
             except SacIOError as se:
@@ -519,8 +549,8 @@ class Data(object):
     def getOriginalWFData(self):
         return self.wforiginal
 
-    def getSynWFData(self):
+    def getAltWFdata(self):
-        return self.wfsyn
+        return self.wf_alt
 
     def resetWFData(self):
         """

pylot/core/io/getEventListFromXML.py

@@ -8,6 +8,11 @@
    Edited for use in PyLoT
    JG, igem, 01/2022
 """
+
+import os
+import argparse
+import numpy as np
+import matplotlib.pyplot as plt
 import glob
 
 from obspy.core.event import read_events

pylot/core/io/getQualitiesfromxml.py

@@ -1,138 +0,0 @@
-#!/usr/bin/python
-# -*- coding: utf-8 -*-
-
-"""
-   Script to get onset uncertainties from Quakeml.xml files created by PyLoT.
-   Uncertainties are tranformed into quality classes and visualized via histogram if desired.
-   Ludger Küperkoch, BESTEC GmbH, 07/2017
-   rev.: Ludger Küperkoch, igem, 10/2020
-   Edited for usage in PyLoT: Jeldrik Gaal, igem, 01/2022
-"""
-
-import glob
-
-import matplotlib.pyplot as plt
-import numpy as np
-from obspy.core.event import read_events
-
-
-def getQualitiesfromxml(path):
-    # uncertainties
-    ErrorsP = [0.02, 0.04, 0.08, 0.16]
-    ErrorsS = [0.04, 0.08, 0.16, 0.32]
-
-    Pw0 = []
-    Pw1 = []
-    Pw2 = []
-    Pw3 = []
-    Pw4 = []
-    Sw0 = []
-    Sw1 = []
-    Sw2 = []
-    Sw3 = []
-    Sw4 = []
-
-    # data path
-    dp = path + '/e*/*.xml'
-    # list of all available xml-files
-    xmlnames = glob.glob(dp)
-
-    # read all onset weights
-    for names in xmlnames:
-        print("Getting onset weights from {}".format(names))
-        cat = read_events(names)
-        arrivals = cat.events[0].picks
-        for Pick in arrivals:
-            if Pick.phase_hint[0] == 'P':
-                if Pick.time_errors.uncertainty <= ErrorsP[0]:
-                    Pw0.append(Pick.time_errors.uncertainty)
-                elif Pick.time_errors.uncertainty > ErrorsP[0] and \
-                        Pick.time_errors.uncertainty <= ErrorsP[1]:
-                    Pw1.append(Pick.time_errors.uncertainty)
-                elif Pick.time_errors.uncertainty > ErrorsP[1] and \
-                        Pick.time_errors.uncertainty <= ErrorsP[2]:
-                    Pw2.append(Pick.time_errors.uncertainty)
-                elif Pick.time_errors.uncertainty > ErrorsP[2] and \
-                        Pick.time_errors.uncertainty <= ErrorsP[3]:
-                    Pw3.append(Pick.time_errors.uncertainty)
-                elif Pick.time_errors.uncertainty > ErrorsP[3]:
-                    Pw4.append(Pick.time_errors.uncertainty)
-                else:
-                    pass
-            elif Pick.phase_hint[0] == 'S':
-                if Pick.time_errors.uncertainty <= ErrorsS[0]:
-                    Sw0.append(Pick.time_errors.uncertainty)
-                elif Pick.time_errors.uncertainty > ErrorsS[0] and \
-                        Pick.time_errors.uncertainty <= ErrorsS[1]:
-                    Sw1.append(Pick.time_errors.uncertainty)
-                elif Pick.time_errors.uncertainty > ErrorsS[1] and \
-                        Pick.time_errors.uncertainty <= ErrorsS[2]:
-                    Sw2.append(Pick.time_errors.uncertainty)
-                elif Pick.time_errors.uncertainty > ErrorsS[2] and \
-                        Pick.time_errors.uncertainty <= ErrorsS[3]:
-                    Sw3.append(Pick.time_errors.uncertainty)
-                elif Pick.time_errors.uncertainty > ErrorsS[3]:
-                    Sw4.append(Pick.time_errors.uncertainty)
-                else:
-                    pass
-            else:
-                print("Phase hint not defined for picking!")
-                pass
-    # get percentage of weights
-    numPweights = np.sum([len(Pw0), len(Pw1), len(Pw2), len(Pw3), len(Pw4)])
-    numSweights = np.sum([len(Sw0), len(Sw1), len(Sw2), len(Sw3), len(Sw4)])
-    try:
-        P0perc = 100.0 / numPweights * len(Pw0)
-    except:
-        P0perc = 0
-    try:
-        P1perc = 100.0 / numPweights * len(Pw1)
-    except:
-        P1perc = 0
-    try:
-        P2perc = 100.0 / numPweights * len(Pw2)
-    except:
-        P2perc = 0
-    try:
-        P3perc = 100.0 / numPweights * len(Pw3)
-    except:
-        P3perc = 0
-    try:
-        P4perc = 100.0 / numPweights * len(Pw4)
-    except:
-        P4perc = 0
-    try:
-        S0perc = 100.0 / numSweights * len(Sw0)
-    except:
-        Soperc = 0
-    try:
-        S1perc = 100.0 / numSweights * len(Sw1)
-    except:
-        S1perc = 0
-    try:
-        S2perc = 100.0 / numSweights * len(Sw2)
-    except:
-        S2perc = 0
-    try:
-        S3perc = 100.0 / numSweights * len(Sw3)
-    except:
-        S3perc = 0
-    try:
-        S4perc = 100.0 / numSweights * len(Sw4)
-    except:
-        S4perc = 0
-
-    weights = ('0', '1', '2', '3', '4')
-    y_pos = np.arange(len(weights))
-    width = 0.34
-    p1, = plt.bar(0 - width, P0perc, width, color='black')
-    p2, = plt.bar(0, S0perc, width, color='red')
-    plt.bar(y_pos - width, [P0perc, P1perc, P2perc, P3perc, P4perc], width, color='black')
-    plt.bar(y_pos, [S0perc, S1perc, S2perc, S3perc, S4perc], width, color='red')
-    plt.ylabel('%')
-    plt.xticks(y_pos, weights)
-    plt.xlim([-0.5, 4.5])
-    plt.xlabel('Qualities')
-    plt.title('{0} P-Qualities, {1} S-Qualities'.format(numPweights, numSweights))
-    plt.legend([p1, p2], ['P-Weights', 'S-Weights'])
-    plt.show()

pylot/core/io/phases.py

@@ -17,7 +17,7 @@ from pylot.core.io.location import create_event, \
     create_magnitude
 from pylot.core.pick.utils import select_for_phase, get_quality_class
 from pylot.core.util.utils import getOwner, full_range, four_digits, transformFilterString4Export, \
-    backtransformFilterString
+    backtransformFilterString, loopIdentifyPhase, identifyPhase
 
 
 def add_amplitudes(event, amplitudes):
@@ -280,6 +280,7 @@ def picksdict_from_picks(evt):
             infile = os.path.join(os.path.expanduser('~'), '.pylot', 'pylot.in')
             print('Using default input file {}'.format(infile))
             parameter = PylotParameter(infile)
+            pick.phase_hint = identifyPhase(pick.phase_hint)
             if pick.phase_hint == 'P':
                 errors = parameter['timeerrorsP']
             elif pick.phase_hint == 'S':
@@ -375,7 +376,6 @@ def picks_from_picksdict(picks, creation_info=None):
 
 
 def reassess_pilot_db(root_dir, db_dir, out_dir=None, fn_param=None, verbosity=0):
-    import glob
     # TODO: change root to datapath
     db_root = os.path.join(root_dir, db_dir)
     evt_list = glob.glob1(db_root, 'e????.???.??')
@@ -518,7 +518,7 @@ def writephases(arrivals, fformat, filename, parameter=None, eventinfo=None):
         arrivals = chooseArrivals(arrivals)  # MP MP what is chooseArrivals? It is not defined anywhere
         for key in arrivals:
             # P onsets
-            if arrivals[key].has_key('P'):
+            if 'P' in arrivals[key]:
                 try:
                     fm = arrivals[key]['P']['fm']
                 except KeyError as e:
@@ -552,7 +552,7 @@ def writephases(arrivals, fformat, filename, parameter=None, eventinfo=None):
                                                                                            ss_ms,
                                                                                            pweight))
             # S onsets
-            if arrivals[key].has_key('S') and arrivals[key]['S']['mpp'] is not None:
+            if 'S' in arrivals[key] and arrivals[key]['S']['mpp'] is not None:
                 fm = '?'
                 onset = arrivals[key]['S']['mpp']
                 year = onset.year
@@ -671,7 +671,7 @@ def writephases(arrivals, fformat, filename, parameter=None, eventinfo=None):
         arrivals = chooseArrivals(arrivals)  # MP MP what is chooseArrivals? It is not defined anywhere
         for key in arrivals:
             # P onsets
-            if arrivals[key].has_key('P') and arrivals[key]['P']['mpp'] is not None:
+            if 'P' in  arrivals[key] and arrivals[key]['P']['mpp'] is not None:
                 if arrivals[key]['P']['weight'] < 4:
                     Ponset = arrivals[key]['P']['mpp']
                     pyear = Ponset.year
@@ -700,7 +700,7 @@ def writephases(arrivals, fformat, filename, parameter=None, eventinfo=None):
                     fid.write('%-5s P1       %4.0f %02d %02d %02d %02d %05.02f   %5.3f -999.   0.00 -999.  0.00\n'
                               % (key, pyear, pmonth, pday, phh, pmm, Pss, pstd))
             # S onsets
-            if arrivals[key].has_key('S') and arrivals[key]['S']['mpp'] is not None:
+            if 'S' in arrivals[key] and arrivals[key]['S']['mpp'] is not None:
                 if arrivals[key]['S']['weight'] < 4:
                     Sonset = arrivals[key]['S']['mpp']
                     syear = Sonset.year
@@ -769,7 +769,7 @@ def writephases(arrivals, fformat, filename, parameter=None, eventinfo=None):
         usedarrivals = chooseArrivals(arrivals)
         for key in usedarrivals:
             # P onsets
-            if usedarrivals[key].has_key('P'):
+            if 'P' in usedarrivals[key]:
                 if usedarrivals[key]['P']['weight'] < 4:
                     n += 1
                     stat = key
@@ -783,7 +783,7 @@ def writephases(arrivals, fformat, filename, parameter=None, eventinfo=None):
                     else:
                         fid.write('%-4sP%d%6.2f\n' % (stat, Pweight, Prt))
                         # S onsets
-            if usedarrivals[key].has_key('S'):
+            if 'S' in usedarrivals[key]:
                 if usedarrivals[key]['S']['weight'] < 4:
                     n += 1
                     stat = key
@@ -829,13 +829,13 @@ def writephases(arrivals, fformat, filename, parameter=None, eventinfo=None):
         # prefer manual picks
         usedarrivals = chooseArrivals(arrivals)
         for key in usedarrivals:
-            if usedarrivals[key].has_key('P'):
+            if 'P' in usedarrivals[key]:
                 # P onsets
                 if usedarrivals[key]['P']['weight'] < 4:
                     Ponset = usedarrivals[key]['P']['mpp']
                     Prt = Ponset - stime  # onset time relative to source time
                     fid.write('%s    %6.3f  1  P\n' % (key, Prt))
-            if usedarrivals[key].has_key('S'):
+            if 'S' in usedarrivals[key]:
                 # S onsets
                 if usedarrivals[key]['S']['weight'] < 4:
                     Sonset = usedarrivals[key]['S']['mpp']
@@ -880,7 +880,7 @@ def writephases(arrivals, fformat, filename, parameter=None, eventinfo=None):
         # prefer manual picks
         usedarrivals = chooseArrivals(arrivals)
         for key in usedarrivals:
-            if usedarrivals[key].has_key('P'):
+            if 'P' in usedarrivals[key]:
                 if usedarrivals[key]['P']['weight'] < 4 and usedarrivals[key]['P']['fm'] is not None:
                     stat = key
                     for i in range(len(picks)):
@@ -962,7 +962,7 @@ def writephases(arrivals, fformat, filename, parameter=None, eventinfo=None):
         arrivals = chooseArrivals(arrivals)  # MP MP what is chooseArrivals? It is not defined anywhere
         # write phase lines
         for key in arrivals:
-            if arrivals[key].has_key('P'):
+            if 'P' in arrivals[key]:
                 if arrivals[key]['P']['weight'] < 4 and arrivals[key]['P']['fm'] is not None:
                     stat = key
                     ccode = arrivals[key]['P']['channel']
@@ -1015,7 +1015,6 @@ def chooseArrivals(arrivals):
     :return: arrivals but with the manual picks prefered if possible
     """
     # If len of arrivals is greater than 2 it comes from autopicking so only autopicks are available
-    print("=== CHOOSE ===")
     if len(arrivals) > 2:
         return arrivals
     if arrivals['auto'] and arrivals['manual']:
@@ -1056,37 +1055,63 @@ def merge_picks(event, picks):
     return event
 
 
-def getQualitiesfromxml(xmlnames, ErrorsP, ErrorsS, plotflag=1):
+def getQualitiesfromxml(path, errorsP, errorsS, plotflag=1, figure=None, verbosity=0):
     """
     Script to get onset uncertainties from Quakeml.xml files created by PyLoT.
     Uncertainties are tranformed into quality classes and visualized via histogram if desired.
     Ludger Küperkoch, BESTEC GmbH, 07/2017
-    :param xmlnames: list of xml obspy event files containing picks
+    :param path: path containing xml files
-    :type xmlnames: list
+    :type path: str
-    :param ErrorsP: time errors of P waves for the four discrete quality classes
+    :param errorsP: time errors of P waves for the four discrete quality classes
-    :type ErrorsP:
+    :type errorsP:
-    :param ErrorsS: time errors of S waves for the four discrete quality classes
+    :param errorsS: time errors of S waves for the four discrete quality classes
-    :type ErrorsS:
+    :type errorsS:
     :param plotflag:
     :type plotflag:
     :return:
     :rtype:
     """
 
-    from pylot.core.pick.utils import get_quality_class
+    def calc_perc(uncertainties, ntotal):
-    from pylot.core.util.utils import loopIdentifyPhase, identifyPhase
+        ''' simple function that calculates percentage of number of uncertainties (list length)'''
+        if len(uncertainties) == 0:
+            return 0
+        else:
+            return 100. / ntotal * len(uncertainties)
+
+    def calc_weight_perc(psweights, weight_ids):
+        ''' calculate percentages of different weights (pick classes!?) of total number of uncertainties of a phase'''
+        # count total number of list items for this phase
+        numWeights = np.sum([len(weight) for weight in psweights.values()])
+
+        # iterate over all available weights to return a list with percentages for plotting
+        plot_list = []
+        for weight_id in weight_ids:
+            plot_list.append(calc_perc(psweights[weight_id], numWeights))
+
+        return plot_list, numWeights
+
+    # get all xmlfiles in path (maybe this should be changed to one xml file for this function, selectable via GUI?)
+    xmlnames = glob.glob(os.path.join(path, '*.xml'))
+    if len(xmlnames) == 0:
+        print(f'No files found in path {path}.')
+        return False
+
+    # first define possible phases here
+    phases = ['P', 'S']
+
+    # define possible weights (0-4)
+    weight_ids = list(range(5))
+
+    # put both error lists in a dictionary with P/S key so that amount of code can be halfed by simply using P/S as key
+    errors = dict(P=errorsP, S=errorsS)
+
+    # create dictionaries for each phase (P/S) with a dictionary of empty list for each weight defined in weights
+    # tuple above
+    weights = {}
+    for phase in phases:
+        weights[phase] = {weight_id: [] for weight_id in weight_ids}
 
-    # read all onset weights
-    Pw0 = []
-    Pw1 = []
-    Pw2 = []
-    Pw3 = []
-    Pw4 = []
-    Sw0 = []
-    Sw1 = []
-    Sw2 = []
-    Sw3 = []
-    Sw4 = []
     for names in xmlnames:
         print("Getting onset weights from {}".format(names))
         cat = read_events(names)
@@ -1094,119 +1119,60 @@ def getQualitiesfromxml(xmlnames, ErrorsP, ErrorsS, plotflag=1):
         arrivals = cat.events[0].picks
         arrivals_copy = cat_copy.events[0].picks
         # Prefere manual picks if qualities are sufficient!
-        for Pick in arrivals:
+        for pick in arrivals:
-            if Pick.method_id.id.split('/')[1] == 'manual':
+            if pick.method_id.id.split('/')[1] == 'manual':
-                mstation = Pick.waveform_id.station_code
+                mstation = pick.waveform_id.station_code
                 mstation_ext = mstation + '_'
                 for mpick in arrivals_copy:
-                    phase = identifyPhase(loopIdentifyPhase(Pick.phase_hint))
+                    phase = identifyPhase(loopIdentifyPhase(pick.phase_hint)) # MP MP catch if this fails?
-                    if phase == 'P':
+                    if ((mpick.waveform_id.station_code == mstation) or
-                        if ((mpick.waveform_id.station_code == mstation) or
+                        (mpick.waveform_id.station_code == mstation_ext)) and \
-                            (mpick.waveform_id.station_code == mstation_ext)) and \
+                            (mpick.method_id.id.split('/')[1] == 'auto') and \
-                                (mpick.method_id.id.split('/')[1] == 'auto') and \
+                            (mpick.time_errors['uncertainty'] <= errors[phase][3]):
-                                (mpick.time_errors['uncertainty'] <= ErrorsP[3]):
+                        del mpick
-                            del mpick
+                        break
-                            break
-                    elif phase == 'S':
-                        if ((mpick.waveform_id.station_code == mstation) or
-                            (mpick.waveform_id.station_code == mstation_ext)) and \
-                                (mpick.method_id.id.split('/')[1] == 'auto') and \
-                                (mpick.time_errors['uncertainty'] <= ErrorsS[3]):
-                            del mpick
-                            break
         lendiff = len(arrivals) - len(arrivals_copy)
         if lendiff != 0:
             print("Found manual as well as automatic picks, prefered the {} manual ones!".format(lendiff))
 
-        for Pick in arrivals_copy:
+        for pick in arrivals_copy:
-            phase = identifyPhase(loopIdentifyPhase(Pick.phase_hint))
+            phase = identifyPhase(loopIdentifyPhase(pick.phase_hint))
-            if phase == 'P':
+            uncertainty = pick.time_errors.uncertainty
-                Pqual = get_quality_class(Pick.time_errors.uncertainty, ErrorsP)
+            if not uncertainty:
-                if Pqual == 0:
+                if verbosity > 0:
-                    Pw0.append(Pick.time_errors.uncertainty)
+                    print('No uncertainty, pick {} invalid!'.format(pick.method_id.id))
-                elif Pqual == 1:
+                continue
-                    Pw1.append(Pick.time_errors.uncertainty)
+            # check P/S phase
-                elif Pqual == 2:
+            if phase not in phases:
-                    Pw2.append(Pick.time_errors.uncertainty)
-                elif Pqual == 3:
-                    Pw3.append(Pick.time_errors.uncertainty)
-                elif Pqual == 4:
-                    Pw4.append(Pick.time_errors.uncertainty)
-            elif phase == 'S':
-                Squal = get_quality_class(Pick.time_errors.uncertainty, ErrorsS)
-                if Squal == 0:
-                    Sw0.append(Pick.time_errors.uncertainty)
-                elif Squal == 1:
-                    Sw1.append(Pick.time_errors.uncertainty)
-                elif Squal == 2:
-                    Sw2.append(Pick.time_errors.uncertainty)
-                elif Squal == 3:
-                    Sw3.append(Pick.time_errors.uncertainty)
-                elif Squal == 4:
-                    Sw4.append(Pick.time_errors.uncertainty)
-            else:
                 print("Phase hint not defined for picking!")
-                pass
+                continue
+
+            qual = get_quality_class(uncertainty, errors[phase])
+            weights[phase][qual].append(uncertainty)
 
     if plotflag == 0:
-        Punc = [Pw0, Pw1, Pw2, Pw3, Pw4]
+        p_unc = [weights['P'][weight_id] for weight_id in weight_ids]
-        Sunc = [Sw0, Sw1, Sw2, Sw3, Sw4]
+        s_unc = [weights['S'][weight_id] for weight_id in weight_ids]
-        return Punc, Sunc
+        return p_unc, s_unc
     else:
+        if not figure:
+            fig = plt.figure()
+        ax = fig.add_subplot(111)
         # get percentage of weights
-        numPweights = np.sum([len(Pw0), len(Pw1), len(Pw2), len(Pw3), len(Pw4)])
+        listP, numPweights = calc_weight_perc(weights['P'], weight_ids)
-        numSweights = np.sum([len(Sw0), len(Sw1), len(Sw2), len(Sw3), len(Sw4)])
+        listS, numSweights = calc_weight_perc(weights['S'], weight_ids)
-        if len(Pw0) > 0:
-            P0perc = 100 / numPweights * len(Pw0)
-        else:
-            P0perc = 0
-        if len(Pw1) > 0:
-            P1perc = 100 / numPweights * len(Pw1)
-        else:
-            P1perc = 0
-        if len(Pw2) > 0:
-            P2perc = 100 / numPweights * len(Pw2)
-        else:
-            P2perc = 0
-        if len(Pw3) > 0:
-            P3perc = 100 / numPweights * len(Pw3)
-        else:
-            P3perc = 0
-        if len(Pw4) > 0:
-            P4perc = 100 / numPweights * len(Pw4)
-        else:
-            P4perc = 0
-        if len(Sw0) > 0:
-            S0perc = 100 / numSweights * len(Sw0)
-        else:
-            S0perc = 0
-        if len(Sw1) > 0:
-            S1perc = 100 / numSweights * len(Sw1)
-        else:
-            S1perc = 0
-        if len(Sw2) > 0:
-            S2perc = 100 / numSweights * len(Sw2)
-        else:
-            S2perc = 0
-        if len(Sw3) > 0:
-            S3perc = 100 / numSweights * len(Sw3)
-        else:
-            S3perc = 0
-        if len(Sw4) > 0:
-            S4perc = 100 / numSweights * len(Sw4)
-        else:
-            S4perc = 0
 
-        weights = ('0', '1', '2', '3', '4')
+        y_pos = np.arange(len(weight_ids))
-        y_pos = np.arange(len(weights))
         width = 0.34
-        plt.bar(y_pos - width, [P0perc, P1perc, P2perc, P3perc, P4perc], width, color='black')
+        ax.bar(y_pos - width, listP, width, color='black')
-        plt.bar(y_pos, [S0perc, S1perc, S2perc, S3perc, S4perc], width, color='red')
+        ax.bar(y_pos, listS, width, color='red')
-        plt.ylabel('%')
+        ax.set_ylabel('%')
-        plt.xticks(y_pos, weights)
+        ax.set_xticks(y_pos, weight_ids)
-        plt.xlim([-0.5, 4.5])
+        ax.set_xlim([-0.5, 4.5])
-        plt.xlabel('Qualities')
+        ax.set_xlabel('Qualities')
-        plt.title('{0} P-Qualities, {1} S-Qualities'.format(numPweights, numSweights))
+        ax.set_title('{0} P-Qualities, {1} S-Qualities'.format(numPweights, numSweights))
-        plt.show()
 
-        return [P0perc, P1perc, P2perc, P3perc, P4perc], [S0perc, S1perc, S2perc, S3perc, S4perc]
+        if not figure:
+            fig.show()
+
+        return listP, listS

pylot/core/loc/nll.py

@@ -82,8 +82,8 @@ def locate(fnin, parameter=None):
     :param fnin: external program name
     :return: None
     """
-    exe_path = which('NLLoc', parameter)
+    exe_path = os.path.join(parameter['nllocbin'], 'NLLoc')
-    if exe_path is None:
+    if not os.path.isfile(exe_path):
         raise NLLocError('NonLinLoc executable not found; check your '
                          'environment variables')
 

pylot/core/pick/autopick.py

@@ -22,7 +22,7 @@ from pylot.core.pick.picker import AICPicker, PragPicker
 from pylot.core.pick.utils import checksignallength, checkZ4S, earllatepicker, \
     getSNR, fmpicker, checkPonsets, wadaticheck, get_quality_class
 from pylot.core.util.utils import getPatternLine, gen_Pool, \
-    get_Bool, identifyPhaseID, get_None, correct_iplot
+    get_bool, identifyPhaseID, get_None, correct_iplot
 
 
 def autopickevent(data, param, iplot=0, fig_dict=None, fig_dict_wadatijack=None, ncores=0, metadata=None, origin=None):
@@ -429,7 +429,7 @@ class AutopickStation(object):
             if station_coords is None:
                 exit_taupy()
                 raise AttributeError('Warning: Could not find station in metadata')
-            # TODO raise when metadata.get_coordinates returns None
+            # TODO: raise when metadata.get_coordinates returns None
             source_origin = origin[0]
             model = TauPyModel(taup_model)
             taup_phases = self.pickparams['taup_phases']
@@ -471,11 +471,13 @@ class AutopickStation(object):
             """If taupy failed to calculate theoretical starttimes, picking continues.
             For this a clean exit is required, since the P starttime is no longer relative to the theoretic onset but
             to the vertical trace starttime, eg. it can't be < 0."""
+            # TODO here the pickparams is modified, instead of a copy
             if self.pickparams["pstart"] < 0:
-                # TODO here the pickparams is modified, instead of a copy
                 self.pickparams["pstart"] = 0
+            if self.pickparams["sstart"] < 0:
+                self.pickparams["sstart"] = 0
 
-        if self.pickparams["use_taup"] is False:
+        if get_bool(self.pickparams["use_taup"]) is False:
             # correct user mistake where a relative cuttime is selected (pstart < 0) but use of taupy is disabled/ has
             # not the required parameters
             exit_taupy()
@@ -499,6 +501,20 @@ class AutopickStation(object):
         self.pickparams["pstart"] = max(self.pickparams["pstart"], 0)
         self.pickparams["pstop"] = min(self.pickparams["pstop"], len(self.ztrace) * self.ztrace.stats.delta)
 
+        if self.horizontal_traces_exist():
+            # for the two horizontal components take earliest and latest time to make sure that the s onset is not clipped
+            # if start and endtime of horizontal traces differ, the s windowsize will automatically increase
+            trace_s_start = min([self.etrace.stats.starttime, self.ntrace.stats.starttime])
+            # modifiy sstart and sstop relative to estimated first S arrival (relative to station time axis)
+            self.pickparams["sstart"] += (self.origin[0].time + estFirstS) - trace_s_start
+            self.pickparams["sstop"] += (self.origin[0].time + estFirstS) - trace_s_start
+            print('autopick: CF calculation times respectively:'
+                  ' sstart: {} s, sstop: {} s'.format(self.pickparams["sstart"], self.pickparams["sstop"]))
+            # make sure pstart and pstop are inside the starttime/endtime of horizontal traces
+            self.pickparams["sstart"] = max(self.pickparams["sstart"], 0)
+            self.pickparams["sstop"] = min(self.pickparams["sstop"], len(self.ntrace) * self.ntrace.stats.delta,
+                                           len(self.etrace) * self.etrace.stats.delta)
+
     def autopickstation(self):
         """
         Main function of autopickstation, which calculates P and S picks and returns them in a dictionary.
@@ -508,6 +524,17 @@ class AutopickStation(object):
         station's value is the station name on which the picks were calculated.
         :rtype: dict
         """
+
+        if get_bool(self.pickparams['use_taup']) is True and self.origin is not None:
+            try:
+                # modify pstart, pstop, sstart, sstop to be around theoretical onset if taupy should be used,
+                # else do nothing
+                self.modify_starttimes_taupy()
+            except AttributeError as ae:
+                print(ae)
+            except MissingTraceException as mte:
+                print(mte)
+
         try:
             self.pick_p_phase()
         except MissingTraceException as mte:
@@ -515,13 +542,15 @@ class AutopickStation(object):
         except PickingFailedException as pfe:
             print(pfe)
 
-        if self.horizontal_traces_exist() and self.p_results.weight is not None and self.p_results.weight < 4:
+        if self.horizontal_traces_exist():
-            try:
+            if (self.p_results.weight is not None and self.p_results.weight < 4) or \
-                self.pick_s_phase()
+                    get_bool(self.pickparams.get('use_taup')):
-            except MissingTraceException as mte:
+                try:
-                print(mte)
+                    self.pick_s_phase()
-            except PickingFailedException as pfe:
+                except MissingTraceException as mte:
-                print(pfe)
+                    print(mte)
+                except PickingFailedException as pfe:
+                    print(pfe)
 
         self.plot_pick_results()
         self.finish_picking()
@@ -661,13 +690,14 @@ class AutopickStation(object):
                         ax2.plot([self.s_results.lpp, self.s_results.lpp], [-1.1, 1.1], 'g--', label='lpp')
                         ax2.plot([self.s_results.epp, self.s_results.epp], [-1.1, 1.1], 'g--', label='epp')
                         title = '{channel}, S weight={sweight}, SNR={snr:7.2}, SNR[dB]={snrdb:7.2}'
-                        ax2.set_title(title.format(channel=self.etrace.stats.channel,
+                        ax2.set_title(title.format(channel=str(self.etrace.stats.channel),
-                                                   sweight=self.s_results.weight,
+                                                   sweight=str(self.s_results.weight),
-                                                   snr=self.s_results.snr,
+                                                   snr=str(self.s_results.snr),
-                                                   snrdb=self.s_results.snrdb))
+                                                   snrdb=str(self.s_results.snrdb)))
                     else:
                         title = '{channel}, S weight={sweight}, SNR=None, SNR[dB]=None'
-                        ax2.set_title(title.format(channel=self.etrace.stats.channel, sweight=self.s_results.weight))
+                        ax2.set_title(title.format(channel=str(self.etrace.stats.channel),
+                                                   sweight=str(self.s_results.weight)))
                 ax2.legend(loc=1)
                 ax2.set_yticks([])
                 ax2.set_ylim([-1.5, 1.5])
@@ -795,15 +825,6 @@ class AutopickStation(object):
         # save filtered trace in instance for later plotting
         self.tr_filt_z_bpz2 = tr_filt
 
-        if get_Bool(self.pickparams['use_taup']) is True and self.origin is not None:
-            try:
-                # modify pstart, pstop to be around theoretical onset if taupy should be used, else does nothing
-                self.modify_starttimes_taupy()
-            except AttributeError as ae:
-                print(ae)
-            except MissingTraceException as mte:
-                print(mte)
-
         Lc = self.pickparams['pstop'] - self.pickparams['pstart']
 
         Lwf = self.ztrace.stats.endtime - self.ztrace.stats.starttime
@@ -915,7 +936,7 @@ class AutopickStation(object):
             "minFMSNR"]:
             # if SNR is high enough, try to determine first motion of onset
             self.set_current_figure('fm_picker')
-            self.p_results.fm = fmpicker(self.zstream, z_copy, self.pickparams["fmpickwin"], self.p_results.mpp,
+            self.p_results.fm = fmpicker(self.zstream.copy(), z_copy, self.pickparams["fmpickwin"], self.p_results.mpp,
                                          self.iplot, self.current_figure, self.current_linecolor)
         msg = "autopickstation: P-weight: {}, SNR: {}, SNR[dB]: {}, Polarity: {}"
         msg = msg.format(self.p_results.weight, self.p_results.snr, self.p_results.snrdb, self.p_results.fm)
@@ -1127,9 +1148,11 @@ class AutopickStation(object):
                   ''.format(self.s_results.weight, self.s_results.snr, self.s_results.snrdb))
 
     def pick_s_phase(self):
-
+        if get_bool(self.pickparams.get('use_taup')) is True:
-        # determine time window for calculating CF after P onset
+            cuttimesh = (self.pickparams.get('sstart'), self.pickparams.get('sstop'))
-        cuttimesh = self._calculate_cuttimes(type='S', iteration=1)
+        else:
+            # determine time window for calculating CF after P onset
+            cuttimesh = self._calculate_cuttimes(type='S', iteration=1)
 
         # calculate autoregressive CF
         self.arhcf1 = self._calculate_autoregressive_cf_s_pick(cuttimesh)

pylot/core/pick/charfuns.py

@@ -16,10 +16,9 @@ autoregressive prediction: application ot local and regional distances, Geophys.
 
 :author: MAGS2 EP3 working group
 """
-
 import numpy as np
-from obspy.core import Stream
 from scipy import signal
+from obspy.core import Stream
 
 
 class CharacteristicFunction(object):
@@ -259,7 +258,7 @@ class HOScf(CharacteristicFunction):
         """
         Function to calculate skewness (statistics of order 3) or kurtosis
         (statistics of order 4), using one long moving window, as published
-        in Kueperkoch et al. (2010).
+        in Kueperkoch et al. (2010), or order 2, i.e. STA/LTA.
         :param data: data, time series (whether seismogram or CF)
         :type data: tuple
         :return: HOS cf
@@ -276,28 +275,47 @@ class HOScf(CharacteristicFunction):
         elif self.getOrder() == 4:  # this is kurtosis
             y = np.power(xnp, 4)
             y1 = np.power(xnp, 2)
+        elif self.getOrder() == 2:  # this is variance, used for STA/LTA processing
+            y = np.power(xnp, 2)
+            y1 = np.power(xnp, 2)
 
         # Initialisation
         # t2: long term moving window
         ilta = int(round(self.getTime2() / self.getIncrement()))
+        ista = int(round((self.getTime2() / 10) / self.getIncrement()))  # TODO: still hard coded!!
         lta = y[0]
         lta1 = y1[0]
+        sta = y[0]
         # moving windows
         LTA = np.zeros(len(xnp))
+        STA = np.zeros(len(xnp))
         for j in range(0, len(xnp)):
             if j < 4:
                 LTA[j] = 0
+                STA[j] = 0
+            elif j <= ista:
+                lta = (y[j] + lta * (j - 1)) / j
+                if self.getOrder() == 2:
+                    sta = (y[j] + sta * (j - 1)) / j
+            # elif j < 4:
             elif j <= ilta:
                 lta = (y[j] + lta * (j - 1)) / j
                 lta1 = (y1[j] + lta1 * (j - 1)) / j
+                if self.getOrder() == 2:
+                    sta = (y[j] - y[j - ista]) / ista + sta
             else:
                 lta = (y[j] - y[j - ilta]) / ilta + lta
                 lta1 = (y1[j] - y1[j - ilta]) / ilta + lta1
+                if self.getOrder() == 2:
+                    sta = (y[j] - y[j - ista]) / ista + sta
             # define LTA
             if self.getOrder() == 3:
                 LTA[j] = lta / np.power(lta1, 1.5)
             elif self.getOrder() == 4:
                 LTA[j] = lta / np.power(lta1, 2)
+            else:
+                LTA[j] = lta
+                STA[j] = sta
 
         # remove NaN's with first not-NaN-value,
         # so autopicker doesnt pick discontinuity at start of the trace
@@ -306,7 +324,10 @@ class HOScf(CharacteristicFunction):
             first = ind[0]
             LTA[:first] = LTA[first]
 
-        self.cf = LTA
+        if self.getOrder() > 2:
+            self.cf = LTA
+        else:  # order 2 means STA/LTA!
+            self.cf = STA / LTA
         self.xcf = x
 
 

pylot/core/pick/picker.py

@@ -336,7 +336,7 @@ class AICPicker(AutoPicker):
                     self.slope = 1 / (len(dataslope) * self.Data[0].stats.delta) * (datafit[-1] - datafit[0])
                     # normalize slope to maximum of cf to make it unit independent
                     self.slope /= aicsmooth[iaicmax]
-            except ValueError as e:
+            except Exception as e:
                 print("AICPicker: Problems with data fitting! {}".format(e))
 
         else:

pylot/core/pick/utils.py

@@ -15,7 +15,7 @@ import numpy as np
 from obspy.core import Stream, UTCDateTime
 from scipy.signal import argrelmax
 
-from pylot.core.util.utils import get_Bool, get_None, SetChannelComponents
+from pylot.core.util.utils import get_bool, get_None, SetChannelComponents
 
 
 def earllatepicker(X, nfac, TSNR, Pick1, iplot=0, verbosity=1, fig=None, linecolor='k'):
@@ -62,8 +62,8 @@ def earllatepicker(X, nfac, TSNR, Pick1, iplot=0, verbosity=1, fig=None, linecol
     plt_flag = 0
     try:
         iplot = int(iplot)
-    except:
+    except ValueError:
-        if get_Bool(iplot):
+        if get_bool(iplot):
             iplot = 2
         else:
             iplot = 0
@@ -275,7 +275,7 @@ def fmpicker(Xraw, Xfilt, pickwin, Pick, iplot=0, fig=None, linecolor='k'):
             try:
                 P1 = np.polyfit(xslope1, xraw[islope1], 1)
                 datafit1 = np.polyval(P1, xslope1)
-            except ValueError as e:
+            except Exception as e:
                 print("fmpicker: Problems with data fit! {}".format(e))
                 print("Skip first motion determination!")
                 return FM
@@ -321,7 +321,7 @@ def fmpicker(Xraw, Xfilt, pickwin, Pick, iplot=0, fig=None, linecolor='k'):
             try:
                 P2 = np.polyfit(xslope2, xfilt[islope2], 1)
                 datafit2 = np.polyval(P2, xslope2)
-            except ValueError as e:
+            except Exception as e:
                 emsg = 'fmpicker: polyfit failed: {}'.format(e)
                 print(emsg)
                 return FM
@@ -816,7 +816,7 @@ def checksignallength(X, pick, minsiglength, pickparams, iplot=0, fig=None, line
     try:
         iplot = int(iplot)
     except:
-        if get_Bool(iplot):
+        if get_bool(iplot):
             iplot = 2
         else:
             iplot = 0
@@ -1130,7 +1130,7 @@ def checkZ4S(X, pick, pickparams, iplot, fig=None, linecolor='k'):
     try:
         iplot = int(iplot)
     except:
-        if get_Bool(iplot):
+        if get_bool(iplot):
             iplot = 2
         else:
             iplot = 0
@@ -1320,7 +1320,7 @@ def get_quality_class(uncertainty, weight_classes):
     :return: quality of pick (0-4)
     :rtype: int
     """
-    if not uncertainty: return max(weight_classes)
+    if not uncertainty: return len(weight_classes)
     try:
         # create generator expression containing all indices of values in weight classes that are >= than uncertainty.
         # call next on it once to receive first value
@@ -1499,7 +1499,7 @@ def get_pickparams(pickparam):
     first_motion_params = dict(zip(first_motion_names, fm_parameter_values))
     signal_length_params = dict(zip(signal_length_names, sl_parameter_values))
 
-    p_params['use_taup'] = get_Bool(p_params['use_taup'])
+    p_params['use_taup'] = get_bool(p_params['use_taup'])
 
     return p_params, s_params, first_motion_params, signal_length_params
 

pylot/core/util/array_map.py

@@ -75,7 +75,6 @@ class Array_map(QtWidgets.QWidget):
 
         self._style = None if not hasattr(parent, '_style') else parent._style
 
-        self.show()
 
     def init_map(self):
         self.init_colormap()
@@ -464,20 +463,23 @@ class Array_map(QtWidgets.QWidget):
                                                                   transform=ccrs.PlateCarree(), label='deleted'))
 
     def openPickDlg(self, ind):
-        data = self._parent.get_data().getWFData()
+        wfdata = self._parent.get_data().getWFData()
+        wfdata_comp = self._parent.get_data().getWFDataComp()
         for index in ind:
             network, station = self._station_onpick_ids[index].split('.')[:2]
             pyl_mw = self._parent
             try:
-                data = data.select(station=station)
+                wfdata = wfdata.select(station=station)
-                if not data:
+                wfdata_comp = wfdata_comp.select(station=station)
+                if not wfdata:
                     self._warn('No data for station {}'.format(station))
                     return
                 pickDlg = PickDlg(self._parent, parameter=self.parameter,
-                                  data=data, network=network, station=station,
+                                  data=wfdata.copy(), data_compare=wfdata_comp.copy(), network=network, station=station,
                                   picks=self._parent.get_current_event().getPick(station),
                                   autopicks=self._parent.get_current_event().getAutopick(station),
                                   filteroptions=self._parent.filteroptions, metadata=self.metadata,
+                                  model=self.parameter.get('taup_model'),
                                   event=pyl_mw.get_current_event())
             except Exception as e:
                 message = 'Could not generate Plot for station {st}.\n {er}'.format(st=station, er=e)

pylot/core/util/dataprocessing.py

@@ -271,7 +271,7 @@ class Metadata(object):
                     continue
             invtype, robj = self._read_metadata_file(os.path.join(path_to_inventory, fname))
             try:
-                robj.get_coordinates(station_seed_id)
+                # robj.get_coordinates(station_seed_id)  # TODO: Commented out, failed with Parser, is this needed?
                 self.inventory_files[fname] = {'invtype': invtype,
                                                'data': robj}
                 if station_seed_id in self.seed_ids.keys():

pylot/core/util/defaults.py

@@ -26,9 +26,7 @@ elif system_name == "Windows":
 # suffix for phase name if not phase identified by last letter (P, p, etc.)
 ALTSUFFIX = ['diff', 'n', 'g', '1', '2', '3']
 
-FILTERDEFAULTS = readDefaultFilterInformation(os.path.join(os.path.expanduser('~'),
+FILTERDEFAULTS = readDefaultFilterInformation()
-                                                           '.pylot',
-                                                           'pylot.in'))
 
 TIMEERROR_DEFAULTS = os.path.join(os.path.expanduser('~'),
                                   '.pylot',

pylot/core/util/gui.py

@@ -56,7 +56,7 @@ def which(program, parameter):
         nllocpath = ":" + parameter.get('nllocbin')
         os.environ['PATH'] += nllocpath
     except Exception as e:
-        print(e.message)
+        print(e)
 
     def is_exe(fpath):
         return os.path.exists(fpath) and os.access(fpath, os.X_OK)

pylot/core/util/thread.py

@@ -41,7 +41,8 @@ class Thread(QThread):
             exctype, value = sys.exc_info()[:2]
             self._executedErrorInfo = '{} {} {}'. \
                 format(exctype, value, traceback.format_exc())
-        sys.stdout = self.old_stdout
+        if self.redirect_stdout:
+            sys.stdout = self.old_stdout
 
     def showProgressbar(self):
         if self.progressText:

pylot/core/util/utils.py

@@ -37,15 +37,14 @@ def getAutoFilteroptions(phase, parameter):
     return filteroptions
 
 
-def readDefaultFilterInformation(fname):
+def readDefaultFilterInformation():
     """
     Read default filter information from pylot.in file
-    :param fname: path to pylot.in file
-    :type fname: str
     :return: dictionary containing the defailt filter information
     :rtype: dict
     """
-    pparam = PylotParameter(fname)
+    pparam = PylotParameter()
+    pparam.reset_defaults()
     return readFilterInformation(pparam)
 
 
@@ -328,7 +327,7 @@ def get_None(value):
         return value
 
 
-def get_Bool(value):
+def get_bool(value):
     """
     Convert string representations of bools to their true boolean value
     :param value:
@@ -336,13 +335,14 @@ def get_Bool(value):
     :return: true boolean value
     :rtype: bool
     """
-    if value in ['True', 'true']:
+    if type(value) is bool:
+        return value
+    elif value in ['True', 'true']:
         return True
     elif value in ['False', 'false']:
         return False
     else:
-        return value
+        return bool(value)
-
 
 def four_digits(year):
     """
@@ -1169,7 +1169,7 @@ def correct_iplot(iplot):
     try:
         iplot = int(iplot)
     except ValueError:
-        if get_Bool(iplot):
+        if get_bool(iplot):
             iplot = 2
         else:
             iplot = 0

pylot/core/util/widgets.py

@@ -1,7 +1,6 @@
 # -*- coding: utf-8 -*-
 """
 Created on Wed Mar 19 11:27:35 2014
-
 @author: sebastianw
 """
 
@@ -54,7 +53,7 @@ from pylot.core.util.utils import prepTimeAxis, full_range, demeanTrace, isSorte
     pick_linestyle_plt, pick_color_plt, \
     check4rotated, check4doubled, merge_stream, identifyPhase, \
     loopIdentifyPhase, trim_station_components, transformFilteroptions2String, \
-    identifyPhaseID, get_Bool, get_None, pick_color, getAutoFilteroptions, SetChannelComponents, \
+    identifyPhaseID, get_bool, get_None, pick_color, getAutoFilteroptions, SetChannelComponents, \
     station_id_remove_channel
 from autoPyLoT import autoPyLoT
 from pylot.core.util.thread import Thread
@@ -90,6 +89,7 @@ class TextLogWidget(QtWidgets.QTextEdit):
     def __init__(self, parent, highlight_input=False):
         super(TextLogWidget, self).__init__(parent)
         self.highlight_input = highlight_input
+        self.append('DUMMY TEXT\n')
 
     def write(self, text):
         self.append(text)
@@ -107,7 +107,9 @@ class LogWidget(QtWidgets.QWidget):
         self.setMinimumHeight(600)
 
         self.stdout = TextLogWidget(self)
+        self.stdout.write('DUMMY TEXT 2\n')
         self.stderr = TextLogWidget(self, highlight_input=True)
+        self.stderr.write('DUMMY TEXT 2\n')
         self.stderr.highlight.connect(self.active_error)
 
         self.tabs = QTabWidget()
@@ -120,7 +122,7 @@ class LogWidget(QtWidgets.QWidget):
         self.layout.addWidget(self.tabs)
 
     def active_error(self):
-        if self.current_active_error == False:
+        if not self.current_active_error:
             self.current_active_error = True
             self.show()
             self.activateWindow()
@@ -221,6 +223,7 @@ class AddMetadataWidget(QWidget):
         self.show()
         # self.__test__()
 
+    # TODO: what is this for? Remove?
     def __test__(self):
         self.add_item(r'/rscratch/minos14/marcel/git/pylot/tests')
         self.add_item(r'/rscratch/minos14/marcel/git/pylot/inputs')
@@ -790,7 +793,7 @@ class WaveformWidgetPG(QtWidgets.QWidget):
 
     def connect_signals(self):
         self.qcombo_processed.activated.connect(self.parent().newWF)
-        self.syn_checkbox.clicked.connect(self.parent().newWF)
+        self.comp_checkbox.clicked.connect(self.parent().newWF)
 
     def init_labels(self):
         self.label_layout.addWidget(self.status_label)
@@ -801,13 +804,13 @@ class WaveformWidgetPG(QtWidgets.QWidget):
         # use widgets as placeholder, so that child widgets keep position when others are hidden
         mid_layout = QHBoxLayout()
         right_layout = QHBoxLayout()
-        mid_layout.addWidget(self.syn_checkbox)
+        mid_layout.addWidget(self.comp_checkbox)
         right_layout.addWidget(self.qcombo_processed)
         mid_widget.setLayout(mid_layout)
         right_widget.setLayout(right_layout)
         self.label_layout.addWidget(mid_widget)
         self.label_layout.addWidget(right_widget)
-        self.syn_checkbox.setLayoutDirection(Qt.RightToLeft)
+        self.comp_checkbox.setLayoutDirection(Qt.RightToLeft)
         self.label_layout.setStretch(0, 4)
         self.label_layout.setStretch(1, 0)
         self.label_layout.setStretch(2, 0)
@@ -822,7 +825,7 @@ class WaveformWidgetPG(QtWidgets.QWidget):
             label = QtWidgets.QLabel()
             self.perm_labels.append(label)
         self.qcombo_processed = QtWidgets.QComboBox()
-        self.syn_checkbox = QtWidgets.QCheckBox('synthetics')
+        self.comp_checkbox = QtWidgets.QCheckBox('Load comparison data')
         self.addQCboxItem('processed', 'green')
         self.addQCboxItem('raw', 'black')
         # self.perm_qcbox_right.setAlignment(2)
@@ -831,9 +834,11 @@ class WaveformWidgetPG(QtWidgets.QWidget):
     def getPlotDict(self):
         return self.plotdict
 
-    def activateObspyDMToptions(self, activate):
+    def activateObspyDMToptions(self, activate: bool) -> None:
-        self.syn_checkbox.setVisible(activate)
+        self.qcombo_processed.setEnabled(activate)
-        self.qcombo_processed.setVisible(activate)
+
+    def activateCompareOptions(self, activate: bool) -> None:
+        self.comp_checkbox.setEnabled(activate)
 
     def setPermText(self, number, text=None, color='black'):
         if not 0 <= number < len(self.perm_labels):
@@ -857,10 +862,8 @@ class WaveformWidgetPG(QtWidgets.QWidget):
     def clearPlotDict(self):
         self.plotdict = dict()
 
-    def plotWFData(self, wfdata, wfsyn=None, title=None, zoomx=None, zoomy=None,
+    def plotWFData(self, wfdata, wfsyn=None, title=None, scaleddata=False, mapping=True,
-                   noiselevel=None, scaleddata=False, mapping=True,
+                   component='*', nth_sample=1, verbosity=0, method='normal', gain=1., shift_syn=0.2):
-                   component='*', nth_sample=1, iniPick=None, verbosity=0,
-                   method='normal', gain=1.):
         if not wfdata:
             print('Nothing to plot.')
             return
@@ -943,7 +946,7 @@ class WaveformWidgetPG(QtWidgets.QWidget):
                     [time for index, time in enumerate(time_ax_syn) if not index % nth_sample] if st_syn else [])
                 trace.data = np.array(
                     [datum * gain + n for index, datum in enumerate(trace.data) if not index % nth_sample])
-                trace_syn.data = np.array([datum + n for index, datum in enumerate(trace_syn.data)
+                trace_syn.data = np.array([datum + n + shift_syn for index, datum in enumerate(trace_syn.data)
                                            if not index % nth_sample] if st_syn else [])
                 plots.append((times, trace.data,
                               times_syn, trace_syn.data))
@@ -1132,12 +1135,12 @@ class PylotCanvas(FigureCanvas):
         ax.set_xlim(self.cur_xlim)
         ax.set_ylim(self.cur_ylim)
         self.refreshPickDlgText()
-        ax.figure.canvas.draw()
+        ax.figure.canvas.draw_idle()
 
     def panRelease(self, gui_event):
         self.press = None
         self.press_rel = None
-        self.figure.canvas.draw()
+        self.figure.canvas.draw_idle()
 
     def panZoom(self, gui_event, threshold=2., factor=1.1):
         if not gui_event.x and not gui_event.y:
@@ -1355,11 +1358,15 @@ class PylotCanvas(FigureCanvas):
             plot_positions[channel] = plot_pos
         return plot_positions
 
-    def plotWFData(self, wfdata, title=None, zoomx=None, zoomy=None,
+    def plotWFData(self, wfdata, wfdata_compare=None, title=None, zoomx=None, zoomy=None,
                    noiselevel=None, scaleddata=False, mapping=True,
                    component='*', nth_sample=1, iniPick=None, verbosity=0,
                    plot_additional=False, additional_channel=None, scaleToChannel=None,
                    snr=None):
+        def get_wf_dict(data: Stream = Stream(), linecolor = 'k', offset: float = 0., **plot_kwargs):
+            return dict(data=data, linecolor=linecolor, offset=offset, plot_kwargs=plot_kwargs)
+
+
         ax = self.axes[0]
         ax.cla()
 
@@ -1370,21 +1377,33 @@ class PylotCanvas(FigureCanvas):
         settings = QSettings()
         compclass = SetChannelComponents.from_qsettings(settings)
 
+        linecolor = (0., 0., 0., 1.) if not self.style else self.style['linecolor']['rgba_mpl']
+
+        plot_streams = dict(wfdata=get_wf_dict(linecolor=linecolor, linewidth=0.7),
+                            wfdata_comp=get_wf_dict(offset=0.1, linecolor='b', alpha=0.7, linewidth=0.5))
+
         if not component == '*':
             alter_comp = compclass.getCompPosition(component)
             # alter_comp = str(alter_comp[0])
 
-            st_select = wfdata.select(component=component)
+            plot_streams['wfdata']['data'] = wfdata.select(component=component)
-            st_select += wfdata.select(component=alter_comp)
+            plot_streams['wfdata']['data'] += wfdata.select(component=alter_comp)
+            if wfdata_compare:
+                plot_streams['wfdata_comp']['data'] = wfdata_compare.select(component=component)
+                plot_streams['wfdata_comp']['data'] += wfdata_compare.select(component=alter_comp)
         else:
-            st_select = wfdata
+            plot_streams['wfdata']['data'] = wfdata
+            if wfdata_compare:
+                plot_streams['wfdata_comp']['data'] = wfdata_compare
+            
+        st_main = plot_streams['wfdata']['data']
 
         if mapping:
-            plot_positions = self.calcPlotPositions(st_select, compclass)
+            plot_positions = self.calcPlotPositions(st_main, compclass)
 
         # list containing tuples of network, station, channel and plot position (for sorting)
         nslc = []
-        for plot_pos, trace in enumerate(st_select):
+        for plot_pos, trace in enumerate(st_main):
             if not trace.stats.channel[-1] in ['Z', 'N', 'E', '1', '2', '3']:
                 print('Warning: Unrecognized channel {}'.format(trace.stats.channel))
                 continue
@@ -1392,44 +1411,48 @@ class PylotCanvas(FigureCanvas):
         nslc.sort()
         nslc.reverse()
 
-        linecolor = (0., 0., 0., 1.) if not self.style else self.style['linecolor']['rgba_mpl']
-
         for n, seed_id in enumerate(nslc):
             network, station, location, channel = seed_id.split('.')
-            st = st_select.select(id=seed_id)
+            for wf_name, wf_dict in plot_streams.items():
-            trace = st[0].copy()
+                st_select = wf_dict.get('data')
-            if mapping:
+                if not st_select:
-                n = plot_positions[trace.stats.channel]
+                    continue
-            if n > nmax:
+                st = st_select.select(id=seed_id)
-                nmax = n
+                trace = st[0].copy()
-            if verbosity:
+                if mapping:
-                msg = 'plotting %s channel of station %s' % (channel, station)
+                    n = plot_positions[trace.stats.channel]
-                print(msg)
+                if n > nmax:
-            stime = trace.stats.starttime - wfstart
+                    nmax = n
-            time_ax = prepTimeAxis(stime, trace)
+                if verbosity:
-            if time_ax is not None:
+                    msg = 'plotting %s channel of station %s' % (channel, station)
-                if scaleToChannel:
+                    print(msg)
-                    st_scale = wfdata.select(channel=scaleToChannel)
+                stime = trace.stats.starttime - wfstart
-                    if st_scale:
+                time_ax = prepTimeAxis(stime, trace)
-                        tr = st_scale[0]
+                if time_ax is not None:
+                    if scaleToChannel:
+                        st_scale = wfdata.select(channel=scaleToChannel)
+                        if st_scale:
+                            tr = st_scale[0]
+                            trace.detrend('constant')
+                            trace.normalize(np.max(np.abs(tr.data)) * 2)
+                            scaleddata = True
+                    if not scaleddata:
                         trace.detrend('constant')
-                        trace.normalize(np.max(np.abs(tr.data)) * 2)
+                        trace.normalize(np.max(np.abs(trace.data)) * 2)
-                        scaleddata = True
-                if not scaleddata:
-                    trace.detrend('constant')
-                    trace.normalize(np.max(np.abs(trace.data)) * 2)
 
-                times = [time for index, time in enumerate(time_ax) if not index % nth_sample]
+                    offset = wf_dict.get('offset')
-                data = [datum + n for index, datum in enumerate(trace.data) if not index % nth_sample]
+
-                ax.axhline(n, color="0.5", lw=0.5)
+                    times = [time for index, time in enumerate(time_ax) if not index % nth_sample]
-                ax.plot(times, data, color=linecolor, linewidth=0.7)
+                    data = [datum + n + offset for index, datum in enumerate(trace.data) if not index % nth_sample]
-                if noiselevel is not None:
+                    ax.axhline(n, color="0.5", lw=0.5)
-                    for level in [-noiselevel[channel], noiselevel[channel]]:
+                    ax.plot(times, data, color=wf_dict.get('linecolor'), **wf_dict.get('plot_kwargs'))
-                        ax.plot([time_ax[0], time_ax[-1]],
+                    if noiselevel is not None:
-                                [n + level, n + level],
+                        for level in [-noiselevel[channel], noiselevel[channel]]:
-                                color=linecolor,
+                            ax.plot([time_ax[0], time_ax[-1]],
-                                linestyle='dashed')
+                                    [n + level, n + level],
-                self.setPlotDict(n, seed_id)
+                                    color=wf_dict.get('linecolor'),
+                                    linestyle='dashed')
+                    self.setPlotDict(n, seed_id)
         if plot_additional and additional_channel:
             compare_stream = wfdata.select(channel=additional_channel)
             if compare_stream:
@@ -1656,8 +1679,8 @@ class PhaseDefaults(QtWidgets.QDialog):
 class PickDlg(QDialog):
     update_picks = QtCore.Signal(dict)
 
-    def __init__(self, parent=None, data=None, station=None, network=None, location=None, picks=None,
+    def __init__(self, parent=None, data=None, data_compare=None, station=None, network=None, location=None, picks=None,
-                 autopicks=None, rotate=False, parameter=None, embedded=False, metadata=None,
+                 autopicks=None, rotate=False, parameter=None, embedded=False, metadata=None, show_comp_data=False,
                  event=None, filteroptions=None, model=None, wftype=None):
         super(PickDlg, self).__init__(parent, Qt.Window)
         self.orig_parent = parent
@@ -1666,6 +1689,7 @@ class PickDlg(QDialog):
         # initialize attributes
         self.parameter = parameter
         self._embedded = embedded
+        self.showCompData = show_comp_data
         self.station = station
         self.network = network
         self.location = location
@@ -1704,22 +1728,6 @@ class PickDlg(QDialog):
         else:
             self.filteroptions = FILTERDEFAULTS
         self.pick_block = False
-        self.nextStation = QtWidgets.QCheckBox('Continue with next station ')
-
-        # comparison channel
-        self.compareChannel = QtWidgets.QComboBox()
-        self.compareChannel.activated.connect(self.resetPlot)
-
-        # scale channel
-        self.scaleChannel = QtWidgets.QComboBox()
-        self.scaleChannel.activated.connect(self.resetPlot)
-
-        # initialize panning attributes
-        self.press = None
-        self.xpress = None
-        self.ypress = None
-        self.cur_xlim = None
-        self.cur_ylim = None
 
         # set attribute holding data
         if data is None or not data:
@@ -1732,6 +1740,31 @@ class PickDlg(QDialog):
                 raise Exception(errmsg)
         else:
             self.data = data
+        self.data_compare = data_compare
+
+        self.nextStation = QtWidgets.QCheckBox('Continue with next station ')
+
+        # comparison channel
+        self.referenceChannel = QtWidgets.QComboBox()
+        self.referenceChannel.activated.connect(self.resetPlot)
+
+        # comparison channel
+        self.compareCB = QtWidgets.QCheckBox()
+        self.compareCB.setChecked(self.showCompData)
+        self.compareCB.clicked.connect(self.switchCompData)
+        self.compareCB.clicked.connect(self.resetPlot)
+        self.compareCB.setVisible(bool(self.data_compare))
+
+        # scale channel
+        self.scaleChannel = QtWidgets.QComboBox()
+        self.scaleChannel.activated.connect(self.resetPlot)
+
+        # initialize panning attributes
+        self.press = None
+        self.xpress = None
+        self.ypress = None
+        self.cur_xlim = None
+        self.cur_ylim = None
 
         self.stime, self.etime = full_range(self.getWFData())
 
@@ -1744,12 +1777,12 @@ class PickDlg(QDialog):
         self.setupUi()
 
         # fill compare and scale channels
-        self.compareChannel.addItem('-', None)
+        self.referenceChannel.addItem('-', None)
         self.scaleChannel.addItem('individual', None)
 
         for trace in self.getWFData():
             channel = trace.stats.channel
-            self.compareChannel.addItem(channel, trace)
+            self.referenceChannel.addItem(channel, trace)
             if not channel[-1] in ['Z', 'N', 'E', '1', '2', '3']:
                 print('Skipping unknown channel for scaling: {}'.format(channel))
                 continue
@@ -1766,7 +1799,7 @@ class PickDlg(QDialog):
         if self.wftype is not None:
             title += ' | ({})'.format(self.wftype)
 
-        self.multicompfig.plotWFData(wfdata=self.getWFData(),
+        self.multicompfig.plotWFData(wfdata=self.getWFData(), wfdata_compare=self.getWFDataComp(),
                                      title=title)
 
         self.multicompfig.setZoomBorders2content()
@@ -1783,7 +1816,7 @@ class PickDlg(QDialog):
 
         # init expected picks using obspy Taup
         try:
-            if self.metadata and model is not None:
+            if self.metadata and model != "None":
                 self.model = TauPyModel(model)
                 self.get_arrivals()
                 self.drawArrivals()
@@ -1873,7 +1906,7 @@ class PickDlg(QDialog):
         self.sChannels.triggered.connect(self.updateChannelSettingsS)
 
         settings = QSettings()
-        self.autoFilterAction.setChecked(get_Bool(settings.value('autoFilter')))
+        self.autoFilterAction.setChecked(get_bool(settings.value('autoFilter')))
 
         # create other widget elements
         phaseitems = [None] + list(FILTERDEFAULTS.keys())
@@ -1942,8 +1975,11 @@ class PickDlg(QDialog):
         _dialtoolbar.addWidget(est_label)
         _dialtoolbar.addWidget(self.plot_arrivals_button)
         _dialtoolbar.addSeparator()
-        _dialtoolbar.addWidget(QtWidgets.QLabel('Compare to channel: '))
+        _dialtoolbar.addWidget(QtWidgets.QLabel('Plot reference channel: '))
-        _dialtoolbar.addWidget(self.compareChannel)
+        _dialtoolbar.addWidget(self.referenceChannel)
+        _dialtoolbar.addSeparator()
+        _dialtoolbar.addWidget(QtWidgets.QLabel('Compare: '))
+        _dialtoolbar.addWidget(self.compareCB)
         _dialtoolbar.addSeparator()
         _dialtoolbar.addWidget(QtWidgets.QLabel('Scaling: '))
         _dialtoolbar.addWidget(self.scaleChannel)
@@ -2278,7 +2314,7 @@ class PickDlg(QDialog):
     def activatePicking(self):
         self.leave_rename_phase()
         self.renamePhaseAction.setEnabled(False)
-        self.compareChannel.setEnabled(False)
+        self.referenceChannel.setEnabled(False)
         self.scaleChannel.setEnabled(False)
         phase = self.currentPhase
         phaseID = self.getPhaseID(phase)
@@ -2310,7 +2346,7 @@ class PickDlg(QDialog):
         self.disconnectPressEvent()
         self.multicompfig.connectEvents()
         self.renamePhaseAction.setEnabled(True)
-        self.compareChannel.setEnabled(True)
+        self.referenceChannel.setEnabled(True)
         self.scaleChannel.setEnabled(True)
         self.connect_pick_delete()
         self.draw()
@@ -2352,7 +2388,7 @@ class PickDlg(QDialog):
         settings = QSettings()
         phaseID = self.getPhaseID(phase)
 
-        if get_Bool(settings.value('useGuiFilter')) or gui_filter:
+        if get_bool(settings.value('useGuiFilter')) or gui_filter:
             filteroptions = self.filteroptions[phaseID]
         else:
             filteroptions = getAutoFilteroptions(phaseID, self.parameter)
@@ -2383,6 +2419,12 @@ class PickDlg(QDialog):
     def getWFData(self):
         return self.data
 
+    def getWFDataComp(self):
+        if self.showCompData:
+            return self.data_compare
+        else:
+            return Stream()
+
     def selectWFData(self, channel):
         component = channel[-1].upper()
         wfdata = Stream()
@@ -2504,24 +2546,33 @@ class PickDlg(QDialog):
 
         stime = self.getStartTime()
 
-        # copy data for plotting
+        # copy wfdata for plotting
-        data = self.getWFData().copy()
+        wfdata = self.getWFData().copy()
-        data = self.getPickPhases(data, phase)
+        wfdata_comp = self.getWFDataComp().copy()
-        data.normalize()
+        wfdata = self.getPickPhases(wfdata, phase)
-        if not data:
+        wfdata_comp = self.getPickPhases(wfdata_comp, phase)
+        for wfd in [wfdata, wfdata_comp]:
+            if wfd:
+                wfd.normalize()
+
+        if not wfdata:
             QtWidgets.QMessageBox.warning(self, 'No channel to plot',
-                                          'No channel to plot for phase: {}.'.format(phase))
+                                          'No channel to plot for phase: {}. '
+                                          'Make sure to select the correct channels for P and S '
+                                          'in the menu in the top panel.'.format(phase))
             self.leave_picking_mode()
             return
 
-        # filter data and trace on which is picked prior to determination of SNR
+        # filter wfdata and trace on which is picked prior to determination of SNR
         filterphase = self.currentFilterPhase()
         if filterphase:
             filteroptions = self.getFilterOptions(filterphase).parseFilterOptions()
             try:
-                data.detrend('linear')
+                for wfd in [wfdata, wfdata_comp]:
-                data.filter(**filteroptions)
+                    if wfd:
-                # wfdata.filter(**filteroptions)# MP MP removed filtering of original data
+                        wfd.detrend('linear')
+                        wfd.filter(**filteroptions)
+                # wfdata.filter(**filteroptions)# MP MP removed filtering of original wfdata
             except ValueError as e:
                 self.qmb = QtWidgets.QMessageBox(QtWidgets.QMessageBox.Icon.Information,
                                                  'Denied',
@@ -2531,8 +2582,8 @@ class PickDlg(QDialog):
         snr = []
         noiselevels = {}
         # determine SNR and noiselevel
-        for trace in data.traces:
+        for trace in wfdata.traces:
-            st = data.select(channel=trace.stats.channel)
+            st = wfdata.select(channel=trace.stats.channel)
             stime_diff = trace.stats.starttime - stime
             result = getSNR(st, (noise_win, gap_win, signal_win), ini_pick - stime_diff)
             snr.append(result[0])
@@ -2543,23 +2594,25 @@ class PickDlg(QDialog):
                 noiselevel = nfac
             noiselevels[trace.stats.channel] = noiselevel
 
-        # prepare plotting of data
+        # prepare plotting of wfdata
-        for trace in data:
+        for wfd in [wfdata, wfdata_comp]:
-            t = prepTimeAxis(trace.stats.starttime - stime, trace)
+            if wfd:
-            inoise = getnoisewin(t, ini_pick, noise_win, gap_win)
+                for trace in wfd:
-            trace = demeanTrace(trace, inoise)
+                    t = prepTimeAxis(trace.stats.starttime - stime, trace)
-            # upscale trace data in a way that each trace is vertically zoomed to noiselevel*factor
+                    inoise = getnoisewin(t, ini_pick, noise_win, gap_win)
-            channel = trace.stats.channel
+                    trace = demeanTrace(trace, inoise)
-            noiselevel = noiselevels[channel]
+                    # upscale trace wfdata in a way that each trace is vertically zoomed to noiselevel*factor
-            noiseScaleFactor = self.calcNoiseScaleFactor(noiselevel, zoomfactor=5.)
+                    channel = trace.stats.channel
-            trace.data *= noiseScaleFactor
+                    noiselevel = noiselevels[channel]
-            noiselevels[channel] *= noiseScaleFactor
+                    noiseScaleFactor = self.calcNoiseScaleFactor(noiselevel, zoomfactor=5.)
+                    trace.data *= noiseScaleFactor
+                    noiselevels[channel] *= noiseScaleFactor
 
         mean_snr = np.mean(snr)
         x_res = getResolutionWindow(mean_snr, parameter.get('extent'))
 
         xlims = [ini_pick - x_res, ini_pick + x_res]
-        ylims = list(np.array([-.5, .5]) + [0, len(data) - 1])
+        ylims = list(np.array([-.5, .5]) + [0, len(wfdata) - 1])
 
         title = self.getStation() + ' picking mode'
         title += ' | SNR: {}'.format(mean_snr)
@@ -2567,9 +2620,10 @@ class PickDlg(QDialog):
             filtops_str = transformFilteroptions2String(filteroptions)
             title += ' | Filteroptions: {}'.format(filtops_str)
 
-        plot_additional = bool(self.compareChannel.currentText())
+        plot_additional = bool(self.referenceChannel.currentText())
-        additional_channel = self.compareChannel.currentText()
+        additional_channel = self.referenceChannel.currentText()
-        self.multicompfig.plotWFData(wfdata=data,
+        self.multicompfig.plotWFData(wfdata=wfdata,
+                                     wfdata_compare=wfdata_comp,
                                      title=title,
                                      zoomx=xlims,
                                      zoomy=ylims,
@@ -2650,7 +2704,7 @@ class PickDlg(QDialog):
             minFMSNR = parameter.get('minFMSNR')
             quality = get_quality_class(spe, parameter.get('timeerrorsP'))
             if quality <= minFMweight and snr >= minFMSNR:
-                FM = fmpicker(self.getWFData().select(channel=channel), wfdata, parameter.get('fmpickwin'),
+                FM = fmpicker(self.getWFData().select(channel=channel).copy(), wfdata.copy(), parameter.get('fmpickwin'),
                               pick - stime_diff)
 
         # save pick times for actual phase
@@ -2960,7 +3014,8 @@ class PickDlg(QDialog):
         self.cur_xlim = self.multicompfig.axes[0].get_xlim()
         self.cur_ylim = self.multicompfig.axes[0].get_ylim()
         # self.multicompfig.updateCurrentLimits()
-        data = self.getWFData().copy()
+        wfdata = self.getWFData().copy()
+        wfdata_comp = self.getWFDataComp().copy()
         title = self.getStation()
         if filter:
             filtoptions = None
@@ -2968,19 +3023,22 @@ class PickDlg(QDialog):
                 filtoptions = self.getFilterOptions(self.getPhaseID(phase), gui_filter=True).parseFilterOptions()
 
             if filtoptions is not None:
-                data.detrend('linear')
+                for wfd in [wfdata, wfdata_comp]:
-                data.taper(0.02, type='cosine')
+                    if wfd:
-                data.filter(**filtoptions)
+                        wfd.detrend('linear')
+                        wfd.taper(0.02, type='cosine')
+                        wfd.filter(**filtoptions)
                 filtops_str = transformFilteroptions2String(filtoptions)
                 title += ' | Filteroptions: {}'.format(filtops_str)
 
         if self.wftype is not None:
             title += ' | ({})'.format(self.wftype)
 
-        plot_additional = bool(self.compareChannel.currentText())
+        plot_additional = bool(self.referenceChannel.currentText())
-        additional_channel = self.compareChannel.currentText()
+        additional_channel = self.referenceChannel.currentText()
         scale_channel = self.scaleChannel.currentText()
-        self.multicompfig.plotWFData(wfdata=data, title=title,
+        self.multicompfig.plotWFData(wfdata=wfdata, wfdata_compare=wfdata_comp,
+                                     title=title,
                                      zoomx=self.getXLims(),
                                      zoomy=self.getYLims(),
                                      plot_additional=plot_additional,
@@ -3022,7 +3080,7 @@ class PickDlg(QDialog):
     @staticmethod
     def getChannelSettingsP(channel):
         settings = QSettings()
-        rval = get_Bool(settings.value('p_channel_{}'.format(channel)))
+        rval = get_bool(settings.value('p_channel_{}'.format(channel)))
         compclass = SetChannelComponents.from_qsettings(settings)
         components = ['Z']
         for component in components[:]:
@@ -3037,7 +3095,7 @@ class PickDlg(QDialog):
     @staticmethod
     def getChannelSettingsS(channel):
         settings = QSettings()
-        rval = get_Bool(settings.value('s_channel_{}'.format(channel)))
+        rval = get_bool(settings.value('s_channel_{}'.format(channel)))
         compclass = SetChannelComponents.from_qsettings(settings)
         components = ['N', 'E']
         for component in components[:]:
@@ -3053,6 +3111,9 @@ class PickDlg(QDialog):
         self.resetZoom()
         self.refreshPlot()
 
+    def switchCompData(self):
+        self.showCompData = self.compareCB.isChecked()
+
     def refreshPlot(self):
         if self.autoFilterAction.isChecked():
             self.filterActionP.setChecked(False)
@@ -3157,7 +3218,6 @@ class CanvasWidget(QWidget):
 class MultiEventWidget(QWidget):
     start = Signal()
     '''
-
     '''
 
     def __init__(self, options=None, parent=None, windowflag=Qt.Window):
@@ -3434,7 +3494,6 @@ class TuneAutopicker(QWidget):
     update = QtCore.Signal(str)
     '''
     QWidget used to modifiy and test picking parameters for autopicking algorithm.
-
     :param: parent
     :type: PyLoT Mainwindow
     '''
@@ -3610,14 +3669,14 @@ class TuneAutopicker(QWidget):
         self.listWidget.scrollToBottom()
 
     def get_current_event(self):
-        path = self.eventBox.currentText()
+        path = self.get_current_event_fp()
         return self.parent().project.getEventFromPath(path)
 
     def get_current_event_name(self):
-        return self.eventBox.currentText().split('/')[-1].split('*')[0]
+        return self.eventBox.currentText().split('/')[-1].rstrip('*')
 
     def get_current_event_fp(self):
-        return self.eventBox.currentText().split('*')[0]
+        return self.eventBox.currentText().rstrip('*')
 
     def get_current_event_picks(self, station):
         event = self.get_current_event()
@@ -3665,16 +3724,19 @@ class TuneAutopicker(QWidget):
             location = None
 
         wfdata = self.data.getWFData()
+        wfdata_comp = self.data.getWFDataComp()
         metadata = self.parent().metadata
         event = self.get_current_event()
         filteroptions = self.parent().filteroptions
         wftype = self.wftype if self.obspy_dmt else ''
         self.pickDlg = PickDlg(self.parent(), data=wfdata.select(station=station).copy(),
+                               data_comp=wfdata_comp.select(station=station).copy(),
                                station=station, network=network,
                                location=location, parameter=self.parameter,
                                picks=self.get_current_event_picks(station),
                                autopicks=self.get_current_event_autopicks(station),
                                metadata=metadata, event=event, filteroptions=filteroptions,
+                               model=self.parameter.get('taup_model'),
                                embedded=True, wftype=wftype)
         self.pickDlg.update_picks.connect(self.picks_from_pickdlg)
         self.pickDlg.update_picks.connect(self.fill_eventbox)
@@ -3752,7 +3814,7 @@ class TuneAutopicker(QWidget):
             self.plot_manual_pick_to_ax(ax=ax, picks=picks, phase='S',
                                         starttime=starttime, quality=qualitySpick)
         for canvas in self.parent().canvas_dict.values():
-            canvas.draw()
+            canvas.draw_idle()
 
     def plot_manual_pick_to_ax(self, ax, picks, phase, starttime, quality):
         mpp = picks[phase]['mpp'] - starttime
@@ -3954,10 +4016,8 @@ class PylotParaBox(QtWidgets.QWidget):
     def __init__(self, parameter, parent=None, windowflag=Qt.Window):
         '''
         Generate Widget containing parameters for PyLoT.
-
         :param: parameter
         :type: PylotParameter (object)
-
         '''
         QtWidgets.QWidget.__init__(self, parent, windowflag)
         self.parameter = parameter
@@ -5194,7 +5254,7 @@ class FilterOptionsDialog(QDialog):
                            'S': QtWidgets.QGroupBox('S Filter')}
 
         settings = QSettings()
-        overwriteFilter = get_Bool(settings.value('useGuiFilter'))
+        overwriteFilter = get_bool(settings.value('useGuiFilter'))
 
         self.overwriteFilterCheckbox = QCheckBox('Overwrite filteroptions')
         self.overwriteFilterCheckbox.setToolTip('Overwrite filter settings for refined pick with GUI settings')
@@ -5556,6 +5616,251 @@ class HelpForm(QDialog):
         self.pageLabel.setText(self.webBrowser.title())
 
 
+class PickQualitiesFromXml(QWidget):
+    """
+            PyLoT widget PickQualitiesFromXml is a QWidget object. It is an UI that enables the user
+            to create a plot showing the pick qualities in the event selected inside the QComboBox created
+            by this Widget. The user can also choose to select all Events.
+            The plot is being shown by a FigureCanvas from matlplotlib.
+    """
+
+    def __init__(self, parent=None, figure=Figure(), path="", inputVar=None):
+        super(PickQualitiesFromXml, self).__init__(parent)
+
+        self.fig = figure
+        self.chooseBox = QComboBox()
+        self.path = path
+        self.inputs = inputVar
+        self.setupUi()
+
+    def setupUi(self):
+        self.setWindowTitle("Get pick qualities from xml files")
+        self.main_layout = QtWidgets.QVBoxLayout()
+        self.figureC = FigureCanvas(self.fig)
+        self.chooseBox = self.createComboBox()
+        if self.chooseBox:
+            self.main_layout.addWidget(self.chooseBox)
+        self.main_layout.addWidget(self.figureC)
+        self.setLayout(self.main_layout)
+
+    def showUI(self):
+        self.show()
+
+    # Creates a QComboBox and adds all events in the current folder as options. Also gives the option to choose all events
+    def createComboBox(self):
+        self.chooseBox.addItems(glob.glob(os.path.join(os.path.dirname(self.path) + "/*/", '*.xml')))
+        self.chooseBox.addItem("All")
+        self.chooseBox.currentIndexChanged.connect(self.selectionChanged)
+        self.chooseBox.setCurrentIndex(self.chooseBox.count() - 1)
+        return self.chooseBox
+
+    # Function that gets called when the user changes the current selection in the QComboBox. Redraws the plot with the new data
+    def selectionChanged(self):
+        self.figureC.setParent(None)
+        if self.chooseBox.currentIndex() == self.chooseBox.count() - 1:
+            (_, _, plot) = getQualitiesfromxml(self.path, self.inputs.get('timeerrorsP'),
+                                               self.inputs.get('timeerrorsS'), plotflag=1)
+            self.figureC = FigureCanvas(plot)
+        else:
+            (_, _, plot) = getQualitiesfromxml(self.path, self.inputs.get('timeerrorsP'),
+                                               self.inputs.get('timeerrorsS'), plotflag=1,
+                                               xmlnames=[self.chooseBox.currentText()])
+            self.figureC = FigureCanvas(plot)
+        self.figureC.draw()
+        self.main_layout.addWidget(self.figureC)
+        self.setLayout(self.main_layout)
+
+
+class SourceSpecWindow(QWidget):
+    def __init__(self, parent=None, figure=Figure()):
+        super(SourceSpecWindow, self).__init__(parent)
+        self.main_layout = QVBoxLayout()
+        self.setWindowTitle("Display source spectrum from selected trace")
+        self.fig = figure
+
+    def setupUi(self):
+        self.figureC = FigureCanvas(self.fig)
+        self.main_layout.addWidget(self.figureC)
+        self.setLayout(self.main_layout)
+
+
+class ChooseWaveFormWindow(QWidget):
+    def __init__(self, parent=None, WaveForms=[], traces=[], stream=None, chooseB=False):
+        super(ChooseWaveFormWindow, self).__init__(parent)
+        self.main_layout = QVBoxLayout()
+        self.setWindowTitle("Choose trace to display source spectrum")
+        self.wFs = WaveForms
+        self.chooseBoxTraces = QComboBox()
+        self.chooseBoxComponent = QComboBox()
+        self.submitButton = QPushButton(text='test')
+
+        self.chooseB = chooseB
+
+        self.verticalButton = QPushButton(text='Z')
+        self.northButton = QPushButton(text='N')
+        self.eastButton = QPushButton(text='E')
+        self.component = ''
+
+        self.traces = traces
+        self.stream = stream
+        self.setupUI()
+        self.currentSpectro = Figure()
+
+    def setupUI(self):
+        self.submitButton.clicked.connect(self.submit)
+        self.createComboBoxTraces()
+        self.main_layout.addWidget(self.chooseBoxTraces)
+
+        if self.chooseB:
+            self.createComboBoxComponent()
+            self.main_layout.addWidget(self.submitButton)
+            self.main_layout.addWidget(self.chooseBoxComponent)
+        else:
+            self.createButtonsComponent()
+
+        self.setLayout(self.main_layout)
+
+    def submit(self):
+        matplotlib.pyplot.close(self.currentSpectro)
+        t = self.chooseBoxTraces.currentText() + " " + self.chooseBoxComponent.currentText()
+        #self.currentSpectro = self.traces[
+        #    self.chooseBoxTraces.currentText()[3:]][self.chooseBoxComponent.currentText()].spectrogram(show=False, title=t)
+        #self.currentSpectro.show()
+        applyFFT()
+
+    def applyFFT(self, trace):
+        tra = self.traces[self.chooseBoxTraces.currentText()[3:]]['Z']
+        transformed = abs(np.fft.rfft(tra.data))
+        print ( transformed )
+        matplotlib.pyplot.plot ( transformed )
+        matplotlib.pyplot.show()
+
+    def applyFFTs(self, tra):
+        transformed = abs(np.fft.rfft(tra.data))
+        print ( transformed )
+        matplotlib.pyplot.plot ( transformed )
+        matplotlib.pyplot.show()
+
+    def submitN(self):
+        matplotlib.pyplot.close(self.currentSpectro)
+        t = self.chooseBoxTraces.currentText() + " " + self.chooseBoxComponent.currentText()
+        self.currentSpectro = self.traces[
+            self.chooseBoxTraces.currentText()[3:]]['N'].spectrogram(show=False, title=t)
+        self.currentSpectro.show()
+
+    def submitE(self):
+        matplotlib.pyplot.close(self.currentSpectro)
+        t = self.chooseBoxTraces.currentText() + " " + self.chooseBoxComponent.currentText()
+        self.currentSpectro = self.traces[
+            self.chooseBoxTraces.currentText()[3:]]['E'].spectrogram(show=False, title=t)
+        self.currentSpectro.show()
+
+    # Creates a QComboBox and adds all traces provided
+    def createComboBoxTraces(self):
+        if len(self.wFs) <= 0:
+            raise 'No traces provided'
+        self.chooseBoxTraces.addItems(self.wFs)
+        self.chooseBoxTraces.currentIndexChanged.connect(self.selectionChanged)
+        self.chooseBoxTraces.setCurrentIndex(0)
+        return self.chooseBoxTraces
+
+    def createButtonsComponent(self):
+        self.northButton.clicked.connect(self.submitN)
+        self.eastButton.clicked.connect(self.submitE)
+        self.verticalButton.clicked.connect(self.submitZ)
+
+        self.main_layout.addWidget(self.verticalButton)
+        self.main_layout.addWidget(self.northButton)
+        self.main_layout.addWidget(self.eastButton)
+
+    def createComboBoxComponent(self):
+        self.chooseBoxComponent.addItems(['Z', 'N', 'E'])
+
+    # Function that gets called when the user changes the current selection in the QComboBox.
+    def selectionChanged(self):
+        pass
+
+
+class SpectrogramTab(QWidget):
+    def __init__(self, traces, wfdata, parent=None):
+        super(SpectrogramTab, self).__init__(parent)
+        self.setupUi()
+        self.traces = traces
+        self.wfdata = wfdata
+
+    def setupUi(self):
+        pass
+    def makeSpecFig(self, direction = 'Z', height = 0, width = 0, parent = None):
+
+        i = 0
+        grams = []
+        figure, axis = matplotlib.pyplot.subplots(len(self.traces), sharex=True)
+
+        start, end = full_range(self.wfdata)
+
+        if height != 0 and width != 0:
+            figure.figsize = (width, height)
+            figure.set_figwidth = width
+            figure.set_figheight = height
+
+        #figure.tight_layout()
+
+        for t in self.traces:
+            tra = self.traces[t][direction]
+            #print(start, end)
+
+            # Set Title
+            if i == 0:
+                if direction == 'Z':
+                    figure.suptitle("section: vertical components")
+                elif direction == 'E':
+                    figure.suptitle("section: east-west components")
+                elif direction == 'N':
+                    figure.suptitle("section: north-south components")
+                axis[i].vlines(0, axis[i].get_ylim()[0], axis[i].get_ylim()[1],
+                          colors='m', linestyles='dashed',
+                          linewidth=2)
+
+            # Different axis settings for visual improvements
+            # axis[i].set_xlim(left=0, right=end - start)
+            # axis[i].spines['top'].set_visible(False)
+            # axis[i].spines['right'].set_visible(False)
+            # # axis[i].spines['left'].set_visible(False)
+
+            # axis[i].tick_params(axis='x', which='both', bottom=False, top=False, labelbottom=False)
+            # if not (len(self.traces) == i - 1):
+                # axis[i].spines['bottom'].set_visible(False)
+            # axis[i].set_yticks([])
+            # axis[i].set_ylabel(t, loc='center', rotation='horizontal')
+
+            #ax.axhline(n, color="0.5", lw=0.5)
+
+
+            grams.append(tra.spectrogram(show=False, axes=axis[i]))
+            i+=1
+
+        #figure.setXLims([0, end - start])
+        figure.set_tight_layout(True)
+        fC = FigureCanvas(figure)
+        return fC
+
+    #for t in self.traces:
+       # tra = self.traces[t]['Z']
+       # transformed = abs(np.fft.rfft(tra.data))
+       # axis[i].plot(transformed, label=t)
+       # # axis[i].tick_params(labelbottom=False)
+       # axis[i].spines['top'].set_visible(False)
+       # axis[i].spines['right'].set_visible(False)
+       # axis[i].spines['left'].set_visible(False)
+       # if not (len(self.traces) == i - 1):
+       #     axis[i].spines['bottom'].set_visible(False)
+       # axis[i].set_yticks([])
+       # axis[i].set_ylabel(t, loc='center', rotation='horizontal')
+       # # axis[i].axis('off')
+       # i += 1
+       # # self.applyFFTs(t)
+
+
 if __name__ == '__main__':
     import doctest
 

tests/test_get_qualities_from_xml.py

@@ -5,7 +5,7 @@ from pylot.core.io.phases import getQualitiesfromxml
 
 class TestQualityFromXML(unittest.TestCase):
     def setUp(self):
-        self.xmlpaths = ['PyLoT_e0019.048.13.xml']
+        self.path = '.'
         self.ErrorsP = [0.02, 0.04, 0.08, 0.16]
         self.ErrorsS = [0.04, 0.08, 0.16, 0.32]
         self.test0_result = [[0.0136956521739, 0.0126, 0.0101612903226, 0.00734848484849, 0.0135069444444,
@@ -23,10 +23,10 @@ class TestQualityFromXML(unittest.TestCase):
                             [92.0, 4.0, 4.0, 0, 0]
 
     def test_result_plotflag0(self):
-        self.assertEqual(getQualitiesfromxml(self.xmlpaths, self.ErrorsP, self.ErrorsS, 0), self.test0_result)
+        self.assertEqual(getQualitiesfromxml(self.path, self.ErrorsP, self.ErrorsS, 0), self.test0_result)
 
     def test_result_plotflag1(self):
-        self.assertEqual(getQualitiesfromxml(self.xmlpaths, self.ErrorsP, self.ErrorsS, 1), self.test1_result)
+        self.assertEqual(getQualitiesfromxml(self.path, self.ErrorsP, self.ErrorsS, 1), self.test1_result)
 
 
 if __name__ == '__main__':