marcel/pylot
3
0
Fork 0
Code Issues 7 Pull Requests 1 Releases Wiki Activity

Merge branch 'develop' into feature/parameter_limits

Browse Source
This commit is contained in:
Darius Arnold 2018-03-07 11:16:37 +01:00
commit ef8c19b747
33 changed files with 123003 additions and 105727 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

212
PyLoT.py Normal file → Executable file
View File

@ -75,7 +75,7 @@ from pylot.core.util.dataprocessing import read_metadata, restitute_data
from pylot.core.util.utils import fnConstructor, getLogin, \
full_range, readFilterInformation, trim_station_components, check4gaps, make_pen, pick_color_plt, \
pick_linestyle_plt, remove_underscores, check4doubled, identifyPhaseID, excludeQualityClasses, has_spe, \
check4rotated, transform_colors_mpl, transform_colors_mpl_str
check4rotated, transform_colors_mpl, transform_colors_mpl_str, getAutoFilteroptions
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, \
@ -182,25 +182,22 @@ class MainWindow(QMainWindow):
settings.setValue("data/dataRoot", dirname)
if settings.value('compclass', None) is None:
settings.setValue('compclass', SetChannelComponents())
if settings.value('useGuiFilter') is None:
settings.setValue('useGuiFilter', False)
if settings.value('output/Format', None) is None:
outformat = QInputDialog.getText(self,
"Enter output format (*.xml, *.cnv, *.obs):",
"Format")
settings.setValue("output/Format", outformat)
if settings.value('autoFilter', None) is None:
settings.setValue('autoFilter', True)
settings.sync()
# setup UI
self.setupUi()
filter_info = readFilterInformation(self._inputs)
p_filter = filter_info['P']
s_filter = filter_info['S']
self.filteroptions = {'P': FilterOptions(p_filter['filtertype'],
p_filter['freq'],
p_filter['order']),
'S': FilterOptions(s_filter['filtertype'],
s_filter['freq'],
s_filter['order'])}
self.updateFilteroptions()
self.loc = False
def setupUi(self):
@ -270,8 +267,12 @@ class MainWindow(QMainWindow):
s_icon.addPixmap(QPixmap(':/icons/key_S.png'))
print_icon = QIcon()
print_icon.addPixmap(QPixmap(':/icons/printer.png'))
filter_icon = QIcon()
filter_icon.addPixmap(QPixmap(':/icons/filter.png'))
self.filter_icon = QIcon()
self.filter_icon.addPixmap(QPixmap(':/icons/filter.png'))
self.filter_icon_p = QIcon()
self.filter_icon_p.addPixmap(QPixmap(':/icons/filter_p.png'))
self.filter_icon_s = QIcon()
self.filter_icon_s.addPixmap(QPixmap(':/icons/filter_s.png'))
z_icon = QIcon()
z_icon.addPixmap(QPixmap(':/icons/key_Z.png'))
n_icon = QIcon()
@ -286,6 +287,8 @@ class MainWindow(QMainWindow):
locate_icon.addPixmap(QPixmap(':/icons/locate_button.png'))
compare_icon = QIcon()
compare_icon.addPixmap(QPixmap(':/icons/compare_button.png'))
qualities_icon = QIcon()
qualities_icon.addPixmap(QPixmap(':/icons/pick_qualities_button.png'))
self.newProjectAction = self.createAction(self, "&New project ...",
self.createNewProject,
QKeySequence.New, newIcon,
@ -366,15 +369,23 @@ class MainWindow(QMainWindow):
self.setParameter,
None, paraIcon,
"Modify Parameter")
self.filterAction = self.createAction(self, "&Filter ...",
self.filterWaveformData,
"Ctrl+F", filter_icon,
"""Toggle un-/filtered waveforms
to be displayed, according to the
desired seismic phase.""", True)
self.filterActionP = createAction(parent=self, text='Apply P Filter',
slot=self.filterP,
icon=self.filter_icon_p,
tip='Toggle filtered/original'
' waveforms',
checkable=True,
shortcut='P')
self.filterActionS = createAction(parent=self, text='Apply S Filter',
slot=self.filterS,
icon=self.filter_icon_s,
tip='Toggle filtered/original'
' waveforms',
checkable=True,
shortcut='S')
filterEditAction = self.createAction(self, "&Filter parameter ...",
self.adjustFilterOptions,
"Alt+F", filter_icon,
"Ctrl+F", self.filter_icon,
"""Adjust filter parameters.""")
self.inventoryAction = self.createAction(self, "Select &Inventory ...",
self.get_new_metadata,
@ -402,6 +413,10 @@ class MainWindow(QMainWindow):
"automatic pick "
"data.", False)
self.compare_action.setEnabled(False)
self.qualities_action = self.createAction(parent=self, text='Show pick qualitites...',
slot=self.pickQualities, shortcut='Alt+Q',
icon=qualities_icon, tip='Histogram of pick qualities')
self.qualities_action.setEnabled(False)
printAction = self.createAction(self, "&Print event ...",
self.show_event_information, QKeySequence.Print,
@ -472,7 +487,7 @@ class MainWindow(QMainWindow):
' the complete project on grid engine.')
self.auto_pick_sge.setEnabled(False)
pickActions = (self.auto_tune, self.auto_pick, self.compare_action)
pickActions = (self.auto_tune, self.auto_pick, self.compare_action, self.qualities_action)
# pickToolBar = self.addToolBar("PickTools")
# pickToolActions = (selectStation, )
@ -500,8 +515,8 @@ class MainWindow(QMainWindow):
self.updateFileMenu()
self.editMenu = self.menuBar().addMenu('&Edit')
editActions = (self.filterAction, filterEditAction, None,
self.selectPAction, self.selectSAction, None,
editActions = (self.filterActionP, self.filterActionS, filterEditAction, None,
#self.selectPAction, self.selectSAction, None,
self.inventoryAction, self.initMapAction, None,
prefsEventAction)
#printAction) #TODO: print event?
@ -561,20 +576,29 @@ class MainWindow(QMainWindow):
style = settings.value('style')
self.set_style(style)
# add event combo box and ref/test buttons
# add event combo box, forward, backward and ref/test buttons
self.eventBox = self.createEventBox()
self.eventBox.setMaxVisibleItems(30)
self.eventBox.setEnabled(False)
self.previous_button = QPushButton('<')
self.next_button = QPushButton('>')
self.init_ref_test_buttons()
self._event_layout = QHBoxLayout()
self._event_layout.addWidget(QLabel('Event: '))
self._event_layout.addWidget(self.eventBox)
self._event_layout.addWidget(self.previous_button)
self._event_layout.addWidget(self.next_button)
self._event_layout.addWidget(self.ref_event_button)
self._event_layout.addWidget(self.test_event_button)
self._event_layout.setStretch(1, 1) # set stretch of item 1 to 1
self._main_layout.addLayout(self._event_layout)
self.eventBox.activated.connect(self.refreshEvents)
self.previous_button.clicked.connect(self.previous_event)
self.next_button.clicked.connect(self.next_event)
self.previous_button.setEnabled(False)
self.next_button.setEnabled(False)
# add main tab widget
self.tabs = QTabWidget(self)
self._main_layout.addWidget(self.tabs)
@ -650,6 +674,8 @@ class MainWindow(QMainWindow):
'''
self.ref_event_button = QtGui.QPushButton('Tune')
self.test_event_button = QtGui.QPushButton('Test')
self.ref_event_button.setMinimumWidth(100)
self.test_event_button.setMinimumWidth(100)
self.ref_event_button.setToolTip('Set manual picks of current ' +
'event as reference picks for autopicker tuning.')
self.test_event_button.setToolTip('Set manual picks of current ' +
@ -763,6 +789,17 @@ class MainWindow(QMainWindow):
def metadata(self, value):
self._metadata = value
def updateFilteroptions(self):
filter_info = readFilterInformation(self._inputs)
p_filter = filter_info['P']
s_filter = filter_info['S']
self.filteroptions = {'P': FilterOptions(p_filter['filtertype'],
p_filter['freq'],
p_filter['order']),
'S': FilterOptions(s_filter['filtertype'],
s_filter['freq'],
s_filter['order'])}
def updateFileMenu(self):
self.fileMenu.clear()
@ -1004,7 +1041,8 @@ class MainWindow(QMainWindow):
and len(item.split('/')[-1].split('.')) == 3
and len(item.split('/')[-1]) == 12]
if not eventlist:
print('No events found! Expected structure for event folders: [evID.DOY.YR]')
print('No events found! Expected structure for event folders: [eEVID.DOY.YR],\n'
' e.g. eventID=1, doy=2, yr=2016: e0001.002.16')
return
else:
return
@ -1353,6 +1391,9 @@ class MainWindow(QMainWindow):
self.cmpw.refresh_tooltips()
self.cmpw.show()
def pickQualities(self):
return
def compareMulti(self):
if not self.compareoptions:
return
@ -1401,7 +1442,8 @@ class MainWindow(QMainWindow):
return None
def getStime(self):
return self._stime
if self.get_data():
return full_range(self.get_data().getWFData())[0]
def addActions(self, target, actions):
for action in actions:
@ -1470,6 +1512,26 @@ class MainWindow(QMainWindow):
if self.tap:
self.tap.fill_eventbox()
def checkEventButtons(self):
if self.eventBox.currentIndex() == 0:
prev_state = False
else:
prev_state = True
if self.eventBox.currentIndex() == len(self.project.eventlist) - 1:
next_state = False
else:
next_state = True
self.previous_button.setEnabled(prev_state)
self.next_button.setEnabled(next_state)
def previous_event(self):
self.eventBox.setCurrentIndex(self.eventBox.currentIndex() - 1)
self.eventBox.activated.emit(-1)
def next_event(self):
self.eventBox.setCurrentIndex(self.eventBox.currentIndex() + 1)
self.eventBox.activated.emit(+1)
def refreshEvents(self):
'''
Refresh GUI when events get changed.
@ -1480,6 +1542,7 @@ class MainWindow(QMainWindow):
# array_map refresh is not necessary when changing event in waveform plot tab,
# but gets necessary when switching from one to another after changing an event.
self._eventChanged = [True, True]
self.checkEventButtons()
self.refreshTabs()
def refreshTabs(self):
@ -1554,18 +1617,9 @@ class MainWindow(QMainWindow):
# else:
# ans = False
self.fnames = self.getWFFnames_from_eventbox()
self.data.setWFData(self.fnames)
wfdat = self.data.getWFData() # all available streams
# remove possible underscores in station names
wfdat = remove_underscores(wfdat)
# check for gaps and doubled channels
check4gaps(wfdat)
check4doubled(wfdat)
# check for stations with rotated components
wfdat = check4rotated(wfdat, self.metadata, verbosity=0)
# trim station components to same start value
trim_station_components(wfdat, trim_start=True, trim_end=False)
self._stime = full_range(self.get_data().getWFData())[0]
self.data.setWFData(self.fnames,
checkRotated=True,
metadata=self.metadata)
def connectWFplotEvents(self):
'''
@ -1633,7 +1687,7 @@ class MainWindow(QMainWindow):
if self.pg:
self.finish_pg_plot()
else:
self._max_xlims = self.dataPlot.getXLims()
self._max_xlims = self.dataPlot.getXLims(self.dataPlot.axes[0])
plotWidget = self.getPlotWidget()
plotDict = plotWidget.getPlotDict()
pos = plotDict.keys()
@ -1661,8 +1715,10 @@ class MainWindow(QMainWindow):
if event.pylot_picks:
self.drawPicks(picktype='manual')
self.locateEvent.setEnabled(True)
self.qualities_action.setEnabled(True)
if event.pylot_autopicks:
self.drawPicks(picktype='auto')
self.qualities_action.setEnabled(True)
if True in self.comparable.values():
self.compare_action.setEnabled(True)
self.draw()
@ -1707,7 +1763,8 @@ class MainWindow(QMainWindow):
self.dataPlot.plotWidget.hideAxis('bottom')
self.dataPlot.plotWidget.hideAxis('left')
else:
self.dataPlot.getAxes().cla()
for ax in self.dataPlot.axes:
ax.cla()
self.loadlocationaction.setEnabled(False)
self.auto_tune.setEnabled(False)
self.auto_pick.setEnabled(False)
@ -1723,18 +1780,19 @@ class MainWindow(QMainWindow):
self.disableSaveEventAction()
self.draw()
def plotWaveformDataThread(self):
def plotWaveformDataThread(self, filter=True):
'''
Open a modal thread to plot current waveform data.
'''
self.clearWaveformDataPlot()
self.wfp_thread = Thread(self, self.plotWaveformData,
arg=filter,
progressText='Plotting waveform data...',
pb_widget=self.mainProgressBarWidget)
self.wfp_thread.finished.connect(self.finishWaveformDataPlot)
self.wfp_thread.start()
def plotWaveformData(self):
def plotWaveformData(self, filter=True):
'''
Plot waveform data to current plotWidget.
'''
@ -1746,6 +1804,10 @@ class MainWindow(QMainWindow):
comp = self.getComponent()
title = 'section: {0} components'.format(zne_text[comp])
wfst = self.get_data().getWFData()
if self.filterActionP.isChecked() and filter:
self.filterWaveformData(plot=False, phase='P')
elif self.filterActionS.isChecked() and filter:
self.filterWaveformData(plot=False, phase='S')
# wfst = self.get_data().getWFData().select(component=comp)
# wfst += self.get_data().getWFData().select(component=alter_comp)
plotWidget = self.getPlotWidget()
@ -1785,18 +1847,51 @@ class MainWindow(QMainWindow):
def pushFilterWF(self, param_args):
self.get_data().filterWFData(param_args)
def filterWaveformData(self):
def filterP(self):
self.filterActionS.setChecked(False)
if self.filterActionP.isChecked():
self.filterWaveformData(phase='P')
else:
self.resetWFData()
def filterS(self):
self.filterActionP.setChecked(False)
if self.filterActionS.isChecked():
self.filterWaveformData(phase='S')
else:
self.resetWFData()
def resetWFData(self):
self.get_data().resetWFData()
self.plotWaveformDataThread()
def filterWaveformData(self, plot=True, phase=None):
if not self.get_current_event():
return
if self.get_data():
if self.getFilterOptions() and self.filterAction.isChecked():
kwargs = self.getFilterOptions()[self.getSeismicPhase()].parseFilterOptions()
self.pushFilterWF(kwargs)
elif self.filterAction.isChecked():
if not phase:
if self.filterActionP.isChecked():
phase = 'P'
elif self.filterActionS.isChecked():
phase = 'S'
if self.getFilterOptions():
if (phase == 'P' and self.filterActionP.isChecked()) or (phase == 'S' and self.filterActionS.isChecked()):
kwargs = self.getFilterOptions()[phase].parseFilterOptions()
self.pushFilterWF(kwargs)
else:
self.get_data().resetWFData()
elif self.filterActionP.isChecked() or self.filterActionS.isChecked():
self.adjustFilterOptions()
else:
self.get_data().resetWFData()
self.plotWaveformDataThread()
self.drawPicks()
self.draw()
if plot:
self.plotWaveformDataThread(filter=False)
#self.drawPicks()
#self.draw()
def getAutoFilteroptions(self, phase):
return getAutoFilteroptions(phase, self._inputs)
def adjustFilterOptions(self):
fstring = "Filter Options"
@ -1805,10 +1900,11 @@ class MainWindow(QMainWindow):
if self.filterDlg.exec_():
filteroptions = self.filterDlg.getFilterOptions()
self.setFilterOptions(filteroptions)
if self.filterAction.isChecked():
if self.filterActionP.isChecked() or self.filterActionS.isChecked():
kwargs = self.getFilterOptions()[self.getSeismicPhase()].parseFilterOptions()
self.pushFilterWF(kwargs)
self.plotWaveformDataThread()
return True
def checkFilterOptions(self):
fstring = "Filter Options"
@ -1885,7 +1981,7 @@ class MainWindow(QMainWindow):
# '[{0}: {1} Hz]'.format(
# self.getFilterOptions().getFilterType(),
# self.getFilterOptions().getFreq()))
# if self.filterAction.isChecked():
# if self.filterActionP.isChecked() or self.filterActionS.isChecked():
# self.filterWaveformData()
def getSeismicPhase(self):
@ -1933,7 +2029,7 @@ class MainWindow(QMainWindow):
if self._shift:
factor = {'up': 1. / 2.,
'down': 2.}
xlims = self.dataPlot.getXLims()
xlims = self.dataPlot.getXLims(self.dataPlot.axes[0])
xdiff = xlims[1] - xlims[0]
xdiff *= factor[button]
xl = x - 0.5 * xdiff
@ -1942,7 +2038,7 @@ class MainWindow(QMainWindow):
xl = self._max_xlims[0]
if xr > self._max_xlims[1]:
xr = self._max_xlims[1]
self.dataPlot.setXLims((xl, xr))
self.dataPlot.setXLims(self.dataPlot.axes[0], (xl, xr))
self.dataPlot.draw()
def pickOnStation(self, gui_event):
@ -1969,7 +2065,7 @@ class MainWindow(QMainWindow):
if not station:
return
self.update_status('picking on station {0}'.format(station))
data = self.get_data().getWFData()
data = self.get_data().getOriginalWFData().copy()
event = self.get_current_event()
pickDlg = PickDlg(self, parameter=self._inputs,
data=data.select(station=station),
@ -1978,6 +2074,9 @@ class MainWindow(QMainWindow):
autopicks=self.getPicksOnStation(station, 'auto'),
metadata=self.metadata, event=event,
filteroptions=self.filteroptions)
if self.filterActionP.isChecked() or self.filterActionS.isChecked():
pickDlg.currentPhase = self.getSeismicPhase()
pickDlg.filterWFData()
pickDlg.nextStation.setChecked(nextStation)
if pickDlg.exec_():
if pickDlg._dirty:
@ -2314,7 +2413,7 @@ class MainWindow(QMainWindow):
if self.pg:
pw = self.getPlotWidget().plotWidget
else:
ax = self.getPlotWidget().axes
ax = self.getPlotWidget().axes[0]
ylims = np.array([-.5, +.5]) + plotID
stat_picks = self.getPicks(type=picktype)[station]
@ -2384,7 +2483,7 @@ class MainWindow(QMainWindow):
else:
if picktype == 'manual':
linestyle_mpp, width_mpp = pick_linestyle_plt(picktype, 'mpp')
color = pick_color_plt(picktype, phase, quality)
color = pick_color_plt(picktype, self.getPhaseID(phase), quality)
if picks['epp'] and picks['lpp']:
ax.fill_between([epp, lpp], ylims[0], ylims[1],
alpha=.25, color=color, label='EPP, LPP')
@ -2926,6 +3025,7 @@ class MainWindow(QMainWindow):
if hasattr(self.project, 'parameter'):
if self.project.parameter:
self._inputs = self.project.parameter
self.updateFilteroptions()
self.tabs.setCurrentIndex(0) # implemented to prevent double-loading of waveform data
self.init_events(new=True)
self.setDirty(False)

4
autoPyLoT.py
View File

@ -338,6 +338,7 @@ def autoPyLoT(input_dict=None, parameter=None, inputfile=None, fnames=None, even
local_mag = LocalMagnitude(corr_dat, evt,
parameter.get('sstop'),
WAscaling, True, iplot)
# update pick with local magnitude property values
for stats, amplitude in local_mag.amplitudes.items():
picks[stats]['S']['Ao'] = amplitude.generic_amplitude
print("Local station magnitudes scaled with:")
@ -412,6 +413,7 @@ def autoPyLoT(input_dict=None, parameter=None, inputfile=None, fnames=None, even
local_mag = LocalMagnitude(corr_dat, evt,
parameter.get('sstop'),
WAscaling, True, iplot)
# update pick with local magnitude property values
for stats, amplitude in local_mag.amplitudes.items():
if stats in picks:
picks[stats]['S']['Ao'] = amplitude.generic_amplitude
@ -496,7 +498,7 @@ if __name__ == "__main__":
help='''full path to the file containing the input
parameters for autoPyLoT''')
parser.add_argument('-p', '-P', '--iplot', type=int,
action='store',
action='store', default=0,
help='''optional, logical variable for plotting: 0=none, 1=partial, 2=all''')
parser.add_argument('-f', '-F', '--fnames', type=str,
action='store',

4
icons.qrc
View File

@ -29,10 +29,12 @@
<file>icons/map.png</file>
<file>icons/openloc.png</file>
<file>icons/compare_button.png</file>
<file>icons/pick_qualities_button.png</file>
<file>icons/locate_button.png</file>
<file>icons/Matlab_PILOT_icon.png</file>
<file>icons/printer.png</file>
<file>icons/delete.png</file>
<file>icons/key_A.png</file>
<file>icons/key_E.png</file>
<file>icons/key_N.png</file>
<file>icons/key_P.png</file>
@ -45,6 +47,8 @@
<file>icons/key_W.png</file>
<file>icons/key_Z.png</file>
<file>icons/filter.png</file>
<file>icons/filter_p.png</file>
<file>icons/filter_s.png</file>
<file>icons/sync.png</file>
<file>icons/zoom_0.png</file>
<file>icons/zoom_in.png</file>

BIN
icons/filter_p.png Normal file
View File

Binary file not shown.
Side by Side

After

Width:  |  Height:  |  Size: 5.4 KiB

BIN
icons/filter_s.png Normal file
View File

Binary file not shown.
Side by Side

After

Width:  |  Height:  |  Size: 7.8 KiB

BIN
icons/key_A.png Normal file
View File

Binary file not shown.
Side by Side

After

Width:  |  Height:  |  Size: 48 KiB

11574
icons_rc_2.py
View File

File diff suppressed because it is too large Load Diff

4580
icons_rc_3.py
View File

File diff suppressed because it is too large Load Diff

1
pylot/RELEASE-VERSION Normal file
View File

@ -0,0 +1 @@
9e74-dirty

39
pylot/core/analysis/magnitude.py
View File

@ -16,7 +16,6 @@ from pylot.core.util.utils import common_range, fit_curve
from scipy import integrate, signal
from scipy.optimize import curve_fit
def richter_magnitude_scaling(delta):
distance = np.array([0, 10, 20, 25, 30, 35, 40, 45, 50, 60, 70, 75, 85, 90, 100, 110,
120, 130, 140, 150, 160, 170, 180, 190, 200, 210, 230, 240, 250,
@ -133,7 +132,7 @@ class Magnitude(object):
station_count=len(self.magnitudes),
azimuthal_gap=self.origin_id.get_referred_object().quality.azimuthal_gap)
else:
# no saling necessary
# no scaling necessary
mag = ope.Magnitude(
mag=np.median([M.mag for M in self.magnitudes.values()]),
magnitude_type=self.type,
@ -233,17 +232,33 @@ class LocalMagnitude(Magnitude):
# check for plot flag (for debugging only)
fig = None
if iplot > 1:
st.plot()
fig = plt.figure()
ax = fig.add_subplot(111)
ax = fig.add_subplot(211)
ax.plot(th, st[0].data, 'k')
ax.plot(th, sqH)
ax.plot(th[iwin], sqH[iwin], 'g')
ax.plot([t0, t0], [0, max(sqH)], 'r', linewidth=2)
ax.title(
'Station %s, RMS Horizontal Traces, WA-peak-to-peak=%4.1f mm' \
% (st[0].stats.station, wapp))
ax.plot([t0 - stime, t0 - stime], [0, max(sqH)], 'r', linewidth=2)
ax.set_title('Station %s, Channel %s, RMS Horizontal Trace, '
'WA-peak-to-peak=%6.3f mm' % (st[0].stats.station,
st[0].stats.channel,
wapp))
ax.set_xlabel('Time [s]')
ax.set_ylabel('Displacement [mm]')
ax = fig.add_subplot(212)
ax.plot(th, st[1].data, 'k')
ax.plot(th, sqH)
ax.plot(th[iwin], sqH[iwin], 'g')
ax.plot([t0 - stime, t0 - stime], [0, max(sqH)], 'r', linewidth=2)
ax.set_title('Channel %s, RMS Horizontal Trace, '
'WA-peak-to-peak=%6.3f mm' % (st[1].stats.channel,
wapp))
ax.set_xlabel('Time [s]')
ax.set_ylabel('Displacement [mm]')
fig.show()
try: input()
except SyntaxError: pass
plt.close(fig)
return wapp, fig
@ -251,6 +266,10 @@ class LocalMagnitude(Magnitude):
for a in self.arrivals:
if a.phase not in 'sS':
continue
# make sure calculating Ml 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
wf = select_for_phase(self.stream.select(
@ -349,6 +368,10 @@ class MomentMagnitude(Magnitude):
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
scopy = self.stream.copy()

168
pylot/core/io/data.py
View File

@ -12,7 +12,8 @@ from pylot.core.io.phases import readPILOTEvent, picks_from_picksdict, \
picksdict_from_pilot, merge_picks
from pylot.core.util.errors import FormatError, OverwriteError
from pylot.core.util.event import Event
from pylot.core.util.utils import fnConstructor, full_range
from pylot.core.util.utils import fnConstructor, full_range, remove_underscores, check4gaps, check4doubled, \
check4rotated, trim_station_components
import pylot.core.loc.velest as velest
@ -99,6 +100,11 @@ class Data(object):
return self
def getPicksStr(self):
"""
Return picks in event data
:return: picks seperated by newlines
:rtype: str
"""
picks_str = ''
for pick in self.get_evt_data().picks:
picks_str += str(pick) + '\n'
@ -106,18 +112,11 @@ class Data(object):
def getParent(self):
"""
:return:
Get PySide.QtGui.QWidget parent object
"""
return self._parent
def isNew(self):
"""
:return:
"""
return self._new
def setNew(self):
@ -125,9 +124,9 @@ class Data(object):
def getCutTimes(self):
"""
:return:
Returns earliest start and latest end of all waveform data
:return: minimum start time and maximum end time as a tuple
:rtype: (UTCDateTime, UTCDateTime)
"""
if self.cuttimes is None:
self.updateCutTimes()
@ -135,22 +134,34 @@ class Data(object):
def updateCutTimes(self):
"""
Update cuttimes to contain earliest start and latest end time
of all waveform data
:rtype: None
"""
self.cuttimes = full_range(self.getWFData())
def getEventFileName(self):
"""
:return:
"""
ID = self.getID()
# handle forbidden filenames especially on windows systems
return fnConstructor(str(ID))
def checkEvent(self, event, fcheck, forceOverwrite=False):
"""
Check information in supplied event and own event and replace own
information with supplied information if own information not exiisting
or forced by forceOverwrite
:param event: Event that supplies information for comparison
:type event: pylot.core.util.event.Event
:param fcheck: check and delete existing information
can be a str or a list of strings of ['manual', 'auto', 'origin', 'magnitude']
:type fcheck: str, [str]
:param forceOverwrite: Set to true to force overwrite own information. If false,
supplied information from event is only used if there is no own information in that
category (given in fcheck: manual, auto, origin, magnitude)
:type forceOverwrite: bool
:return:
:rtype: None
"""
if 'origin' in fcheck:
self.replaceOrigin(event, forceOverwrite)
if 'magnitude' in fcheck:
@ -161,18 +172,47 @@ class Data(object):
self.replacePicks(event, 'manual')
def replaceOrigin(self, event, forceOverwrite=False):
"""
Replace own origin with the one supplied in event if own origin is not
existing or forced by forceOverwrite = True
:param event: Event that supplies information for comparison
:type event: pylot.core.util.event.Event
:param forceOverwrite: always replace own information with supplied one if true
:type forceOverwrite: bool
:return:
:rtype: None
"""
if self.get_evt_data().origins or forceOverwrite:
if event.origins:
print("Found origin, replace it by new origin.")
event.origins = self.get_evt_data().origins
def replaceMagnitude(self, event, forceOverwrite=False):
"""
Replace own magnitude with the one supplied in event if own magnitude is not
existing or forced by forceOverwrite = True
:param event: Event that supplies information for comparison
:type event: pylot.core.util.event.Event
:param forceOverwrite: always replace own information with supplied one if true
:type forceOverwrite: bool
:return:
:rtype: None
"""
if self.get_evt_data().magnitudes or forceOverwrite:
if event.magnitudes:
print("Found magnitude, replace it by new magnitude")
event.magnitudes = self.get_evt_data().magnitudes
def replacePicks(self, event, picktype):
"""
Replace own picks with the one in event
:param event: Event that supplies information for comparison
:type event: pylot.core.util.event.Event
:param picktype: 'auto' or 'manual' picks
:type picktype: str
:return:
:rtype: None
"""
checkflag = 0
picks = event.picks
# remove existing picks
@ -189,10 +229,10 @@ class Data(object):
picks.append(pick)
def exportEvent(self, fnout, fnext='.xml', fcheck='auto', upperErrors=None):
"""
Export event to file
:param fnout: basename of file
:param fnext: file extension
:param fnext: file extension, xml, cnv, obs
:param fcheck: check and delete existing information
can be a str or a list of strings of ['manual', 'auto', 'origin', 'magnitude']
"""
@ -304,17 +344,13 @@ class Data(object):
def getComp(self):
"""
:return:
Get component (ZNE)
"""
return self.comp
def getID(self):
"""
:return:
Get unique resource id
"""
try:
return self.evtdata.get('resource_id').id
@ -323,16 +359,20 @@ class Data(object):
def filterWFData(self, kwargs):
"""
:param kwargs:
Filter waveform data
:param kwargs: arguments to pass through to filter function
"""
self.getWFData().filter(**kwargs)
data = self.getWFData()
data.detrend('linear')
data.taper(0.02, type='cosine')
data.filter(**kwargs)
self.dirty = True
def setWFData(self, fnames):
def setWFData(self, fnames, checkRotated=False, metadata=None):
"""
:param fnames:
Clear current waveform data and set given waveform data
:param fnames: waveform data names to append
:type fnames: list
"""
self.wfdata = Stream()
self.wforiginal = None
@ -340,14 +380,31 @@ class Data(object):
self.appendWFData(fnames)
else:
return False
# various pre-processing steps:
# remove possible underscores in station names
self.wfdata = remove_underscores(self.wfdata)
# check for gaps and doubled channels
check4gaps(self.wfdata)
check4doubled(self.wfdata)
# check for stations with rotated components
if checkRotated and metadata is not None:
self.wfdata = check4rotated(self.wfdata, metadata, verbosity=0)
# trim station components to same start value
trim_station_components(self.wfdata, trim_start=True, trim_end=False)
# make a copy of original data
self.wforiginal = self.getWFData().copy()
self.dirty = False
return True
def appendWFData(self, fnames):
"""
:param fnames:
Read waveform data from fnames and append it to current wf data
:param fnames: waveform data to append
:type fnames: list
"""
assert isinstance(fnames, list), "input parameter 'fnames' is " \
"supposed to be of type 'list' " \
@ -372,54 +429,45 @@ class Data(object):
print(warnmsg)
def getWFData(self):
"""
:return:
"""
return self.wfdata
def getOriginalWFData(self):
"""
:return:
"""
return self.wforiginal
def resetWFData(self):
"""
Set waveform data to original waveform data
"""
self.wfdata = self.getOriginalWFData().copy()
if self.getOriginalWFData():
self.wfdata = self.getOriginalWFData().copy()
else:
self.wfdata = Stream()
self.dirty = False
def resetPicks(self):
"""
Clear all picks from event
"""
self.get_evt_data().picks = []
def get_evt_data(self):
"""
:return:
"""
return self.evtdata
def setEvtData(self, event):
self.evtdata = event
def applyEVTData(self, data, typ='pick', authority_id='rub'):
"""
:param data:
:param typ:
:param authority_id:
Either takes an `obspy.core.event.Event` object and applies all new
information on the event to the actual data if typ is 'event or
creates ObsPy pick objects and append it to the picks list from the
PyLoT dictionary contain all picks if type is pick
:param data: data to apply, either picks or complete event
:type data:
:param typ: which event data to apply, 'pick' or 'event'
:type typ: str
:param authority_id: (currently unused)
:type: str
:raise OverwriteError:
"""

10
pylot/core/io/default_parameters.py
View File

@ -2,6 +2,10 @@
# -*- coding: utf-8 -*-
import numpy as np
"""
Default parameters used for picking
"""
defaults = {'rootpath': {'type': str,
'tooltip': 'project path',
'value': '',
@ -239,7 +243,7 @@ defaults = {'rootpath': {'type': str,
'namestring': 'AIC smooth P'},
'tsmoothP': {'type': float,
'tooltip': 'for HOS/AR, take average of samples for smoothing CF [s]',
'tooltip': 'for HOS/AR, take average of samples in this time window for smoothing CF [s]',
'value': 0.1,
'min:': 0.,
'max': np.inf,
@ -321,7 +325,7 @@ defaults = {'rootpath': {'type': str,
'namestring': ('SNR windows S', 'Noise', 'Safety', 'Signal', 'Slope')},
'aictsmoothS': {'type': float,
'tooltip': 'for AIC-picker, take average of samples for smoothing of AIC-function [s]',
'tooltip': 'for AIC-picker, take average of samples in this time window for smoothing of AIC-function [s]',
'value': 0.5,
'min:': 0.,
'max': np.inf,
@ -361,7 +365,7 @@ defaults = {'rootpath': {'type': str,
'namestring': 'Min SNR'},
'fmpickwin': {'type': float,
'tooltip': 'pick window around P onset for calculating zero crossings',
'tooltip': 'pick window [s] around P onset for calculating zero crossings',
'value': 0.2,
'min:': 0.,
'max': np.inf,

102
pylot/core/io/inputs.py
View File

@ -70,6 +70,7 @@ class PylotParameter(object):
# Set default values of parameter names
def __init_default_paras(self):
"""set default values of parameter names"""
parameters = default_parameters.defaults
self.__defaults = parameters
@ -104,15 +105,34 @@ class PylotParameter(object):
return len(self.__parameter.keys())
def iteritems(self):
"""
Iterate over parameters
:return: key, value tupel
:rtype:
"""
for key, value in self.__parameter.items():
yield key, value
def hasParam(self, parameter):
"""
Check if parameter is in keys
:param parameter: parameter to look for in keys
:type parameter:
:return:
:rtype: bool
"""
if parameter in self.__parameter.keys():
return True
return False
def get(self, *args):
"""
Get first available parameter in args
:param args:
:type args:
:return:
:rtype:
"""
try:
for param in args:
try:
@ -128,15 +148,35 @@ class PylotParameter(object):
raise ParameterError(e)
def get_defaults(self):
"""
get default parameters
:return:
:rtype: dict
"""
return self.__defaults
def get_main_para_names(self):
"""
Get main parameter names
:return: list of keys available in parameters
:rtype:
"""
return self._settings_main
def get_special_para_names(self):
"""
Get pick parameter names
:return: list of keys available in parameters
:rtype:
"""
return self._settings_special_pick
def get_all_para_names(self):
"""
Get all parameter names
:return:
:rtype: list
"""
all_names = []
all_names += self.get_main_para_names()['dirs']
all_names += self.get_main_para_names()['nlloc']
@ -151,6 +191,16 @@ class PylotParameter(object):
return all_names
def checkValue(self, param, value):
"""
Check type of value against expected type of param.
Print warning message if type check fails
:param param:
:type param:
:param value:
:type value:
:return:
:rtype:
"""
is_type = type(value)
expect_type = self.get_defaults()[param]['type']
if not is_type == expect_type and not is_type == tuple:
@ -159,9 +209,25 @@ class PylotParameter(object):
print(Warning(message))
def setParamKV(self, param, value):
"""
set parameter param to value
:param param:
:type param:
:param value:
:type value:
:return:
:rtype: None
"""
self.__setitem__(param, value)
def setParam(self, **kwargs):
"""
Set multiple parameters
:param kwargs:
:type kwargs:
:return:
:rtype: None
"""
for key in kwargs:
self.__setitem__(key, kwargs[key])
@ -170,11 +236,23 @@ class PylotParameter(object):
print('ParameterError:\n non-existent parameter %s' % errmsg)
def reset_defaults(self):
"""
Reset current parameters to default parameters
:return:
:rtype: None
"""
defaults = self.get_defaults()
for param in defaults:
self.setParamKV(param, defaults[param]['value'])
def from_file(self, fnin=None):
"""
read parameters from file and set values to read values
:param fnin: filename
:type fnin:
:return:
:rtype: None
"""
if not fnin:
if self.__filename is not None:
fnin = self.__filename
@ -221,6 +299,13 @@ class PylotParameter(object):
self.__parameter = self._parFileCont
def export2File(self, fnout):
"""
Export parameters to file
:param fnout: Filename of export file
:type fnout: str
:return:
:rtype:
"""
fid_out = open(fnout, 'w')
lines = []
# for key, value in self.iteritems():
@ -257,6 +342,19 @@ class PylotParameter(object):
'quality assessment', None)
def write_section(self, fid, names, title, separator):
"""
write a section of parameters to file
:param fid: File object to write to
:type fid:
:param names: which parameter names to write to file
:type names:
:param title: title of section
:type title: str
:param separator: section separator, written at start of section
:type separator: str
:return:
:rtype:
"""
if separator:
fid.write(separator)
fid.write('#{}#\n'.format(title))
@ -341,7 +439,9 @@ class FilterOptions(object):
def parseFilterOptions(self):
if self:
robject = {'type': self.getFilterType(), 'corners': self.getOrder()}
robject = {'type': self.getFilterType(),
'corners': self.getOrder(),
'zerophase': False}
if not self.getFilterType() in ['highpass', 'lowpass']:
robject['freqmin'] = self.getFreq()[0]
robject['freqmax'] = self.getFreq()[1]

8
pylot/core/io/location.py
View File

@ -54,7 +54,7 @@ def create_arrival(pickresID, cinfo, phase, azimuth=None, dist=None):
def create_creation_info(agency_id=None, creation_time=None, author=None):
'''
get creation info of obspy event
:param agency_id:
:param creation_time:
:param author:
@ -197,9 +197,9 @@ def create_pick(origintime, picknum, picktime, eventnum, cinfo, phase, station,
def create_resourceID(timetohash, restype, authority_id=None, hrstr=None):
'''
:param timetohash:
:type timetohash
create unique resource id
:param timetohash: event origin time to hash
:type timetohash: class: `~obspy.core.utcdatetime.UTCDateTime` object
:param restype: type of the resource, e.g. 'orig', 'earthquake' ...
:type restype: str
:param authority_id: name of the institution carrying out the processing

106
pylot/core/io/phases.py
View File

@ -16,7 +16,8 @@ from pylot.core.io.inputs import PylotParameter
from pylot.core.io.location import create_event, \
create_magnitude
from pylot.core.pick.utils import select_for_phase
from pylot.core.util.utils import getOwner, full_range, four_digits
from pylot.core.util.utils import getOwner, full_range, four_digits, transformFilteroptions2String, \
transformFilterString4Export, backtransformFilterString
def add_amplitudes(event, amplitudes):
@ -118,6 +119,13 @@ def readPILOTEvent(phasfn=None, locfn=None, authority_id='RUB', **kwargs):
def picksdict_from_pilot(fn):
"""
Create pick dictionary from matlab file
:param fn: matlab file
:type fn:
:return: pick dictionary
:rtype: dict
"""
from pylot.core.util.defaults import TIMEERROR_DEFAULTS
picks = dict()
phases_pilot = sio.loadmat(fn)
@ -147,6 +155,13 @@ def picksdict_from_pilot(fn):
def stations_from_pilot(stat_array):
"""
Create stations list from pilot station array
:param stat_array:
:type stat_array:
:return:
:rtype: list
"""
stations = list()
cur_stat = None
for stat in stat_array:
@ -164,6 +179,13 @@ def stations_from_pilot(stat_array):
def convert_pilot_times(time_array):
"""
Convert pilot times to UTCDateTimes
:param time_array: pilot times
:type time_array:
:return:
:rtype:
"""
times = [int(time) for time in time_array]
microseconds = int((time_array[-1] - times[-1]) * 1e6)
times.append(microseconds)
@ -171,6 +193,13 @@ def convert_pilot_times(time_array):
def picksdict_from_obs(fn):
"""
create pick dictionary from obs file
:param fn: filename
:type fn:
:return:
:rtype:
"""
picks = dict()
station_name = str()
for line in open(fn, 'r'):
@ -207,6 +236,10 @@ def picksdict_from_picks(evt):
network = pick.waveform_id.network_code
mpp = pick.time
spe = pick.time_errors.uncertainty
if pick.filter_id:
filter_id = backtransformFilterString(str(pick.filter_id.id))
else:
filter_id = None
try:
picker = str(pick.method_id)
if picker.startswith('smi:local/'):
@ -222,10 +255,15 @@ def picksdict_from_picks(evt):
lpp = mpp + pick.time_errors.upper_uncertainty
epp = mpp - pick.time_errors.lower_uncertainty
except TypeError as e:
msg = e + ',\n falling back to symmetric uncertainties'
if not spe:
msg = 'No uncertainties found for pick: {}. Uncertainty set to 0'.format(pick)
lpp = mpp
epp = mpp
else:
msg = str(e) + ',\n falling back to symmetric uncertainties'
lpp = mpp + spe
epp = mpp - spe
warnings.warn(msg)
lpp = mpp + spe
epp = mpp - spe
phase['mpp'] = mpp
phase['epp'] = epp
phase['lpp'] = lpp
@ -233,6 +271,7 @@ def picksdict_from_picks(evt):
phase['channel'] = channel
phase['network'] = network
phase['picker'] = picker
phase['filter_id'] = filter_id if filter_id is not None else ''
onsets[pick.phase_hint] = phase.copy()
picksdict[picker][station] = onsets.copy()
@ -240,6 +279,16 @@ def picksdict_from_picks(evt):
def picks_from_picksdict(picks, creation_info=None):
"""
Create a list of picks out of a pick dictionary
:param picks: pick dictionary
:type picks: dict
:param creation_info: obspy creation information to apply to picks
:type creation_info:
:param creation_info: obspy creation information to apply to picks
:return: list of picks
:rtype: list
"""
picks_list = list()
for station, onsets in picks.items():
for label, phase in onsets.items():
@ -274,6 +323,13 @@ def picks_from_picksdict(picks, creation_info=None):
pick.waveform_id = ope.WaveformStreamID(station_code=station,
channel_code=ccode,
network_code=ncode)
try:
filter_id = phase['filteroptions']
filter_id = transformFilterString4Export(filter_id)
except KeyError as e:
warnings.warn(e.message, RuntimeWarning)
filter_id = ''
pick.filter_id = filter_id
try:
polarity = phase['fm']
if polarity == 'U' or '+':
@ -290,7 +346,6 @@ def picks_from_picksdict(picks, creation_info=None):
picks_list.append(pick)
return picks_list
def reassess_pilot_db(root_dir, db_dir, out_dir=None, fn_param=None, verbosity=0):
import glob
@ -410,25 +465,24 @@ def writephases(arrivals, fformat, filename, parameter=None, eventinfo=None):
HYPO71, NLLoc, VELEST, HYPOSAT, and hypoDD
:param: arrivals
:type: dictionary containing all phase information including
station ID, phase, first motion, weight (uncertainty),
....
:param arrivals:dictionary containing all phase information including
station ID, phase, first motion, weight (uncertainty), ...
:type arrivals: dict
:param: fformat
:type: string, chosen file format (location routine),
choose between NLLoc, HYPO71, HYPOSAT, VELEST,
HYPOINVERSE, and hypoDD
:param fformat: chosen file format (location routine),
choose between NLLoc, HYPO71, HYPOSAT, VELEST,
HYPOINVERSE, and hypoDD
:type fformat: str
:param: filename, full path and name of phase file
:type: string
:param filename: full path and name of phase file
:type filename: string
:param: parameter, all input information
:type: object
:param parameter: all input information
:type parameter: object
:param: eventinfo, optional, needed for VELEST-cnv file
:param eventinfo: optional, needed for VELEST-cnv file
and FOCMEC- and HASH-input files
:type: `obspy.core.event.Event` object
:type eventinfo: `obspy.core.event.Event` object
"""
if fformat == 'NLLoc':
@ -874,10 +928,20 @@ def merge_picks(event, picks):
def getQualitiesfromxml(xmlnames, ErrorsP, ErrorsS, plotflag=1):
"""
Script to get onset uncertainties from Quakeml.xml files created by PyLoT.
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
:type xmlnames: list
:param ErrorsP: time errors of P waves for the four discrete quality classes
:type ErrorsP:
:param ErrorsS: time errors of S waves for the four discrete quality classes
:type ErrorsS:
:param plotflag:
:type plotflag:
:return:
:rtype:
"""
from pylot.core.pick.utils import getQualityFromUncertainty
from pylot.core.util.utils import loopIdentifyPhase, identifyPhase

8
pylot/core/loc/focmec.py
View File

@ -18,11 +18,11 @@ def export(picks, fnout, parameter, eventinfo):
:param fnout: complete path to the exporting obs file
:type fnout: str
:param: parameter, all input information
:type: object
:param parameter: all input information
:type parameter: object
:param: eventinfo, source information needed for focmec format
:type: list object
:param eventinfo: source information needed for focmec format
:type eventinfo: list object
'''
# write phases to FOCMEC-phase file
writephases(picks, 'FOCMEC', fnout, parameter, eventinfo)

8
pylot/core/loc/hash.py
View File

@ -18,11 +18,11 @@ def export(picks, fnout, parameter, eventinfo):
:param fnout: complete path to the exporting obs file
:type fnout: str
:param: parameter, all input information
:type: object
:param parameter: all input information
:type parameter: object
:param: eventinfo, source information needed for HASH format
:type: list object
:param eventinfo: source information needed for HASH format
:type eventinfo: list object
'''
# write phases to HASH-phase file
writephases(picks, 'HASH', fnout, parameter, eventinfo)

4
pylot/core/loc/hypo71.py
View File

@ -18,8 +18,8 @@ def export(picks, fnout, parameter):
:param fnout: complete path to the exporting obs file
:type fnout: str
:param: parameter, all input information
:type: object
:param parameter: all input information
:type parameter: object
'''
# write phases to HYPO71-phase file
writephases(picks, 'HYPO71', fnout, parameter)

8
pylot/core/loc/hypodd.py
View File

@ -18,11 +18,11 @@ def export(picks, fnout, parameter, eventinfo):
:param fnout: complete path to the exporting obs file
:type fnout: str
:param: parameter, all input information
:type: object
:param parameter: all input information
:type parameter: object
:param: eventinfo, source information needed for hypoDD format
:type: list object
:param eventinfo: source information needed for hypoDD format
:type eventinfo: list object
'''
# write phases to hypoDD-phase file
writephases(picks, 'hypoDD', fnout, parameter, eventinfo)

4
pylot/core/loc/hyposat.py
View File

@ -18,8 +18,8 @@ def export(picks, fnout, parameter):
:param fnout: complete path to the exporting obs file
:type fnout: str
:param: parameter, all input information
:type: object
:param parameter: all input information
:type parameter: object
'''
# write phases to HYPOSAT-phase file
writephases(picks, 'HYPOSAT', fnout, parameter)

20
pylot/core/loc/nll.py
View File

@ -28,8 +28,8 @@ def export(picks, fnout, parameter):
:param fnout: complete path to the exporting obs file
:type fnout: str
:param: parameter, all input information
:type: object
:param parameter: all input information
:type parameter: object
'''
# write phases to NLLoc-phase file
writephases(picks, 'NLLoc', fnout, parameter)
@ -38,19 +38,19 @@ def export(picks, fnout, parameter):
def modify_inputs(ctrfn, root, nllocoutn, phasefn, tttn):
'''
:param ctrfn: name of NLLoc-control file
:type: str
:type ctrfn: str
:param root: root path to NLLoc working directory
:type: str
:type root: str
:param nllocoutn: name of NLLoc-location output file
:type: str
:type nllocoutn: str
:param phasefn: name of NLLoc-input phase file
:type: str
:type phasefn: str
:param tttn: pattern of precalculated NLLoc traveltime tables
:type: str
:type tttn: str
'''
# For locating the event the NLLoc-control file has to be modified!
# create comment line for NLLoc-control file NLLoc-output file
@ -75,9 +75,9 @@ def modify_inputs(ctrfn, root, nllocoutn, phasefn, tttn):
def locate(fnin, infile=None):
"""
takes an external program name
:param fnin:
:return:
takes an external program name and tries to run it
:param fnin: external program name
:return: None
"""
if infile is None:

8
pylot/core/loc/velest.py
View File

@ -18,11 +18,11 @@ def export(picks, fnout, eventinfo, parameter=None):
:param fnout: complete path to the exporting obs file
:type fnout: str
:param: eventinfo, source time needed for VELEST-cnv format
:type: list object
:param eventinfo: source time needed for VELEST-cnv format
:type eventinfo: list object
:param: parameter, all input information
:type: object
:param parameter: all input information
:type parameter: object
'''
# write phases to VELEST-phase file
writephases(picks, 'VELEST', fnout, parameter, eventinfo)

2
pylot/core/pick/autopick.py
View File

@ -155,7 +155,7 @@ def autopickstation(wfstream, pickparam, verbose=False,
:type metadata: tuple (str, ~obspy.io.xseed.parser.Parser)
:param origin: list containing origin objects representing origins for all events
:type origin: list(~obspy.core.event.origin)
:return: :dictionary containing P pick, S pick and station name
:return: dictionary containing P pick, S pick and station name
:rtype: dict
"""

40
pylot/core/pick/compare.py
View File

@ -109,10 +109,9 @@ class Comparison(object):
Comparison is carried out with the help of pdf representation of the picks
and a probabilistic approach to the time difference of two onset
measurements.
:param a: filename for pickset A
:type a: str
:param b: filename for pickset B
:type b: str
:param type: type of the returned `~pylot.core.util.pdf.ProbabilityDensityFunction` object.
Possible values: 'exp' and 'gauss', representing the type of branches of the PDF
:type type: str
:return: dictionary containing the resulting comparison pdfs for all picks
:rtype: dict
"""
@ -142,8 +141,7 @@ class Comparison(object):
istations = range(nstations)
fig, axarr = plt.subplots(nstations, 2, sharex='col', sharey='row')
for n in istations:
station = stations[n]
for n, station in enumerate(stations):
if station not in self.comparison.keys():
continue
compare_pdf = self.comparison[station]
@ -190,6 +188,20 @@ class Comparison(object):
return self.get_array(phase, 'standard_deviation')
def hist_expectation(self, phases='all', bins=20, normed=False):
"""
Plot a histogram of the expectation values of the PDFs.
Expectation represents the time difference between two most likely arrival times
:param phases: type of phases to compare
:type phases: str
:param bins: number of bins in histogram
:type bins: int
:param normed: Normalize histogram
:type normed: bool
:return: None
:rtype: None
"""
phases.strip()
if phases.find('all') is 0:
phases = 'ps'
@ -210,6 +222,20 @@ class Comparison(object):
plt.show()
def hist_standard_deviation(self, phases='all', bins=20, normed=False):
"""
Plot a histogram of the compared standard deviation values of two arrivals.
Standard deviation of two compared picks represents the combined uncertainties/pick errors
(earliest possible pick, latest possible pick)
:param phases: type of phases to compare
:type phases: str
:param bins: number of bins in histogram
:type bins: int
:param normed: Normalize histogram
:type normed: bool
:return: None
:rtype: None
"""
phases.strip()
if phases.find('all') == 0:
phases = 'ps'
@ -370,7 +396,7 @@ class PDFDictionary(object):
class PDFstatistics(object):
"""
This object can be used to get various statistic values from probabillity density functions.
This object can be used to get various statistic values from probability density functions.
Takes a path as argument.
"""

70
pylot/core/pick/picker.py
View File

@ -29,38 +29,36 @@ from pylot.core.pick.utils import getnoisewin, getsignalwin
class AutoPicker(object):
'''
"""
Superclass of different, automated picking algorithms applied on a CF determined
using AIC, HOS, or AR prediction.
'''
"""
warnings.simplefilter('ignore')
def __init__(self, cf, TSNR, PickWindow, iplot=0, aus=None, Tsmooth=None, Pick1=None, fig=None, linecolor='k'):
'''
:param: cf, characteristic function, on which the picking algorithm is applied
:type: `~pylot.core.pick.CharFuns.CharacteristicFunction` object
:param: TSNR, length of time windows around pick used to determine SNR [s]
:type: tuple (T_noise, T_gap, T_signal)
:param: PickWindow, length of pick window [s]
:type: float
:param: iplot, no. of figure window for plotting interims results
:type: integer
:param: aus ("artificial uplift of samples"), find local minimum at i if aic(i-1)*(1+aus) >= aic(i)
:type: float
:param: Tsmooth, length of moving smoothing window to calculate smoothed CF [s]
:type: float
:param: Pick1, initial (prelimenary) onset time, starting point for PragPicker and
EarlLatePicker
:type: float
'''
"""
Create AutoPicker object
:param cf: characteristic function, on which the picking algorithm is applied
:type cf: `~pylot.core.pick.CharFuns.CharacteristicFunction`
:param TSNR: length of time windows around pick used to determine SNR [s], tuple (T_noise, T_gap, T_signal)
:type TSNR: (float, float, float)
:param PickWindow: length of pick window [s]
:type PickWindow: float
:param iplot: flag used for plotting, if > 1, results will be plotted. Use iplot = 0 to disable plotting
:type iplot: int
:param aus: ("artificial uplift of samples"), find local minimum at i if aic(i-1)*(1+aus) >= aic(i)
:type aus: float
:param Tsmooth: length of moving smoothing window to calculate smoothed CF [s]
:type Tsmooth: float
:param Pick1: initial (preliminary) onset time, starting point for PragPicker and EarlLatePicker
:type Pick1: float
:param fig: matplotlib figure used for plotting. If not given and plotting is enabled, a new figure will
be created
:type fig: `~matplotlib.figure.Figure`
:param linecolor: matplotlib line color string
:type linecolor: str
"""
assert isinstance(cf, CharacteristicFunction), "%s is not a CharacteristicFunction object" % str(cf)
self._linecolor = linecolor
@ -79,6 +77,11 @@ class AutoPicker(object):
self.calcPick()
def __str__(self):
"""
String representation of AutoPicker object
:return:
:rtype: str
"""
return '''\n\t{name} object:\n
TSNR:\t\t\t{TSNR}\n
PickWindow:\t{PickWindow}\n
@ -142,12 +145,12 @@ class AutoPicker(object):
class AICPicker(AutoPicker):
'''
"""
Method to derive the onset time of an arriving phase based on CF
derived from AIC. In order to get an impression of the quality of this inital pick,
derived from AIC. In order to get an impression of the quality of this initial pick,
a quality assessment is applied based on SNR and slope determination derived from the CF,
from which the AIC has been calculated.
'''
"""
def calcPick(self):
@ -214,7 +217,6 @@ class AICPicker(AutoPicker):
self.Pick = self.Tcf[i]
break
# quality assessment using SNR and slope from CF
if self.Pick is not None:
# get noise window
inoise = getnoisewin(self.Tcf, self.Pick, self.TSNR[0], self.TSNR[1])
@ -251,8 +253,8 @@ class AICPicker(AutoPicker):
except IndexError:
print("Slope Calculation: empty array islope, check signal window")
return
if len(dataslope) < 1:
print('No data in slope window found!')
if len(dataslope) <= 1:
print('No or not enough data in slope window found!')
return
imaxs, = argrelmax(dataslope)
if imaxs.size:
@ -364,9 +366,9 @@ class AICPicker(AutoPicker):
class PragPicker(AutoPicker):
'''
"""
Method of pragmatic picking exploiting information given by CF.
'''
"""
def calcPick(self):

444
pylot/core/pick/utils.py
View File

@ -13,31 +13,39 @@ import warnings
import matplotlib.pyplot as plt
import numpy as np
from obspy.core import Stream, UTCDateTime
from pylot.core.util.utils import real_Bool, real_None
def earllatepicker(X, nfac, TSNR, Pick1, iplot=0, verbosity=1, fig=None, linecolor='k'):
'''
"""
Function to derive earliest and latest possible pick after Diehl & Kissling (2009)
as reasonable uncertainties. Latest possible pick is based on noise level,
earliest possible pick is half a signal wavelength in front of most likely
pick given by PragPicker or manually set by analyst. Most likely pick
(initial pick Pick1) must be given.
:param: X, time series (seismogram)
:type: `~obspy.core.stream.Stream`
:param: nfac (noise factor), nfac times noise level to calculate latest possible pick
:type: int
:param: TSNR, length of time windows around pick used to determine SNR [s]
:type: tuple (T_noise, T_gap, T_signal)
:param: Pick1, initial (most likely) onset time, starting point for earllatepicker
:type: float
:param: iplot, if given, results are plotted in figure(iplot)
:type: int
'''
:param X: time series (seismogram)
:type X: `~obspy.core.stream.Stream`
:param nfac: (noise factor), nfac times noise level to calculate latest possible pick
:type nfac: int
:param TSNR: length of time windows around pick used to determine SNR [s]
:type TSNR: tuple (T_noise, T_gap, T_signal)
:param Pick1: initial (most likely) onset time, starting point for earllatepicker
:type Pick1: float
:param iplot: if given, results are plotted in figure(iplot)
:type iplot: int
:param verbosity: amount of displayed information about the process:
2 = all
1 = default
0 = none
:type verbosity: int
:param fig: Matplotlib figure ised for plotting
:type fig: `~matplotlib.figure.Figure`
:param linecolor: color for plotting results
:type linecolor: str
:return: tuple containing earliest possible pick, latest possible pick and pick error
:rtype: (float, float, float)
"""
assert isinstance(X, Stream), "%s is not a stream object" % str(X)
@ -54,10 +62,10 @@ def earllatepicker(X, nfac, TSNR, Pick1, iplot=0, verbosity=1, fig=None, linecol
try:
iplot = int(iplot)
except:
if iplot == True or iplot == 'True':
iplot = 2
if real_Bool(iplot):
iplot = 2
else:
iplot = 0
iplot = 0
if verbosity:
print('earllatepicker: Get earliest and latest possible pick'
@ -127,7 +135,7 @@ def earllatepicker(X, nfac, TSNR, Pick1, iplot=0, verbosity=1, fig=None, linecol
PickError = symmetrize_error(diffti_te, diffti_tl)
if iplot > 1:
if fig == None or fig == 'None':
if real_None(fig) is None:
fig = plt.figure() # iplot)
plt_flag = 1
fig._tight = True
@ -162,26 +170,27 @@ def earllatepicker(X, nfac, TSNR, Pick1, iplot=0, verbosity=1, fig=None, linecol
def fmpicker(Xraw, Xfilt, pickwin, Pick, iplot=0, fig=None, linecolor='k'):
'''
"""
Function to derive first motion (polarity) of given phase onset Pick.
Calculation is based on zero crossings determined within time window pickwin
after given onset time.
:param: Xraw, unfiltered time series (seismogram)
:type: `~obspy.core.stream.Stream`
:param: Xfilt, filtered time series (seismogram)
:type: `~obspy.core.stream.Stream`
:param: pickwin, time window after onset Pick within zero crossings are calculated
:type: float
:param: Pick, initial (most likely) onset time, starting point for fmpicker
:type: float
:param: iplot, if given, results are plotted in figure(iplot)
:type: int
'''
:param Xraw: unfiltered time series (seismogram)
:type Xraw: `~obspy.core.stream.Stream`
:param Xfilt: filtered time series (seismogram)
:type Xfilt: `~obspy.core.stream.Stream`
:param pickwin: time window after onset Pick within zero crossings are calculated
:type pickwin: float
:param Pick: initial (most likely) onset time, starting point for fmpicker
:type Pick: float
:param iplot: if given, results are plotted in figure(iplot)
:type iplot: int
:param fig: Matplotlib figure ised for plotting
:type fig: `~matplotlib.figure.Figure`
:param linecolor: color for plotting results
:type linecolor: str
:return: None if first motion detection was skipped, otherwise string indicating the polarity
:rtype: None, str
"""
plt_flag = 0
try:
@ -268,7 +277,7 @@ def fmpicker(Xraw, Xfilt, pickwin, Pick, iplot=0, fig=None, linecolor='k'):
index2 = []
i = 0
for j in range(ipick[0][1], ipick[0][len(t[ipick]) - 1]):
i = i + 1
i += 1
if xfilt[j - 1] <= 0 <= xfilt[j]:
zc2.append(t[ipick][i])
index2.append(i)
@ -321,7 +330,7 @@ def fmpicker(Xraw, Xfilt, pickwin, Pick, iplot=0, fig=None, linecolor='k'):
print("fmpicker: Found polarity %s" % FM)
if iplot > 1:
if fig == None or fig == 'None':
if real_None(fig) is None:
fig = plt.figure() # iplot)
plt_flag = 1
fig._tight = True
@ -360,8 +369,15 @@ def fmpicker(Xraw, Xfilt, pickwin, Pick, iplot=0, fig=None, linecolor='k'):
def crossings_nonzero_all(data):
"""
Returns the indices of zero crossings the data array
:param data: data array (seismic trace)
:type data: `~numpy.ndarray`
:return: array containing indices of zero crossings in the data array.
:rtype: `~numpy.ndarray`
"""
pos = data > 0
npos = ~pos
npos = ~pos # get positions of negative values
return ((pos[:-1] & npos[1:]) | (npos[:-1] & pos[1:])).nonzero()[0]
@ -372,26 +388,28 @@ def symmetrize_error(dte, dtl):
:param dte: relative lower uncertainty
:param dtl: relative upper uncertainty
:return: symmetrized error
:rtype: float
"""
return (dte + 2 * dtl) / 3
def getSNR(X, TSNR, t1, tracenum=0):
'''
"""
Function to calculate SNR of certain part of seismogram relative to
given time (onset) out of given noise and signal windows. A safety gap
between noise and signal part can be set. Returns SNR and SNR [dB] and
noiselevel.
:param: X, time series (seismogram)
:type: `~obspy.core.stream.Stream`
:param: TSNR, length of time windows [s] around t1 (onset) used to determine SNR
:type: tuple (T_noise, T_gap, T_signal)
:param: t1, initial time (onset) from which noise and signal windows are calculated
:type: float
'''
:param X: time series (seismogram)
:type X: `~obspy.core.stream.Stream`
:param TSNR: length of time windows [s] around t1 (onset) used to determine SNR
:type TSNR: (T_noise, T_gap, T_signal)
:param t1: initial time (onset) from which noise and signal windows are calculated
:type t1: float
:param tracenum: used to select the trace in stream X
:type tracenum: int
:return: tuple containing SNR, SNRdB and noise level
:rtype: (float, float, float)
"""
assert isinstance(X, Stream), "%s is not a stream object" % str(X)
@ -435,26 +453,24 @@ def getSNR(X, TSNR, t1, tracenum=0):
def getnoisewin(t, t1, tnoise, tgap):
'''
Function to extract indeces of data out of time series for noise calculation.
Returns an array of indeces.
:param: t, array of time stamps
:type: numpy array
:param: t1, time from which relativ to it noise window is extracted
:type: float
:param: tnoise, length of time window [s] for noise part extraction
:type: float
:param: tgap, safety gap between t1 (onset) and noise window to
ensure, that noise window contains no signal
:type: float
'''
"""
Function to extract indices of data out of time series for noise calculation.
Returns an array of indices.
:param t: array of time stamps
:type t: `numpy.ndarray`
:param t1: time from which relative to it noise window is extracted
:type t1: float
:param tnoise: length of time window [s] for noise part extraction
:type tnoise: float
:param tgap: safety gap between t1 (onset) and noise window to ensure, that
the noise window contains no signal
:type tgap: float
:return: indices of noise window in t
:rtype: `~numpy.ndarray`
"""
# get noise window
inoise, = np.where((t <= max([t1 - tgap, 0])) \
inoise, = np.where((t <= max([t1 - tgap, 0]))
& (t >= max([t1 - tnoise - tgap, 0])))
if np.size(inoise) < 1:
inoise, = np.where((t >= t[0]) & (t <= t1))
@ -465,22 +481,21 @@ def getnoisewin(t, t1, tnoise, tgap):
def getsignalwin(t, t1, tsignal):
'''
"""
Function to extract data out of time series for signal level calculation.
Returns an array of indeces.
:param: t, array of time stamps
:type: numpy array
:param: t1, time from which relativ to it signal window is extracted
:type: float
:param: tsignal, length of time window [s] for signal level calculation
:type: float
'''
Returns an array of indices.
:param t: array of time stamps
:type t: `numpy.ndarray`
:param t1: time from which relative to it signal window is extracted
:type t1: float
:param tsignal: length of time window [s] for signal level calculation
:type tsignal: float
:return: indices of signal window i t
:rtype: `numpy.ndarray`
"""
# get signal window
isignal, = np.where((t <= min([t1 + tsignal, t[-1]])) \
isignal, = np.where((t <= min([t1 + tsignal, t[-1]]))
& (t >= t1))
if np.size(isignal) < 1:
print("getsignalwin: Empty array isignal, check signal window!")
@ -490,27 +505,18 @@ def getsignalwin(t, t1, tsignal):
def getResolutionWindow(snr, extent):
"""
Number -> Float
produce the half of the time resolution window width from given SNR
value
Produce the half of the time resolution window width from given SNR value
SNR >= 3 -> 2 sec HRW
3 > SNR >= 2 -> 5 sec MRW
2 > SNR >= 1.5 -> 10 sec LRW
1.5 > SNR -> 15 sec VLRW
see also Diehl et al. 2009
:parameter: extent, can be 'local', 'regional', 'global'
>>> getResolutionWindow(0.5)
7.5
>>> getResolutionWindow(1.8)
5.0
>>> getResolutionWindow(2.3)
2.5
>>> getResolutionWindow(4)
1.0
>>> getResolutionWindow(2)
2.5
:param snr: Signal to noise ration which decides the witdth of the resolution window
:type snr: float
:param extent: can be 'local', 'regional', 'global'
:type extent: str
:return: half width of the resolution window
:rtype: float
"""
res_wins = {
@ -535,17 +541,17 @@ def getResolutionWindow(snr, extent):
def select_for_phase(st, phase):
'''
takes a STream object and a phase name and returns that particular component
"""
Takes a Stream object and a phase name and returns that particular component
which presumably shows the chosen PHASE best
:param st: stream object containing one or more component[s]
:type st: `~obspy.core.stream.Stream`
:param phase: label of the phase for which the stream selection is carried
out; 'P' or 'S'
:param phase: label of the phase for which the stream selection is carried out; 'P' or 'S'
:type phase: str
:return:
'''
:return: stream object containing the selected phase
:rtype: `~obspy.core.stream.Stream`
"""
from pylot.core.util.defaults import SetChannelComponents
sel_st = Stream()
@ -569,21 +575,22 @@ def select_for_phase(st, phase):
def wadaticheck(pickdic, dttolerance, iplot=0, fig_dict=None):
'''
"""
Function to calculate Wadati-diagram from given P and S onsets in order
to detect S pick outliers. If a certain S-P time deviates by dttolerance
from regression of S-P time the S pick is marked and down graded.
: param: pickdic, dictionary containing picks and quality parameters
: type: dictionary
: param: dttolerance, maximum adjusted deviation of S-P time from
S-P time regression
: type: float
: param: iplot, if iplot > 1, Wadati diagram is shown
: type: int
'''
:param pickdic: dictionary containing picks and quality parameters
:type pickdic: dict
:param dttolerance: dttolerance, maximum adjusted deviation of S-P time from
S-P time regression
:type dttolerance: float
:param iplot: iplot, if iplot > 1, Wadati diagram is shown
:type iplot: int
:param fig_dict: Matplotlib figure used for plotting
:type fig_dict: `~matplotlib.figure.Figure`
:return: dictionary containing all onsets that passed the wadati check
:rtype: dict
"""
checkedonsets = pickdic
@ -699,8 +706,8 @@ def wadaticheck(pickdic, dttolerance, iplot=0, fig_dict=None):
for Ppick, SPtime, station in zip(Ppicks, SPtimes, stations):
ax.text(Ppick, SPtime + 0.01, '{0}'.format(station), color='0.25')
ax.set_title('Wadati-Diagram, %d S-P Times, Vp/Vs(raw)=%5.2f,' \
'Vp/Vs(checked)=%5.2f' % (len(SPtimes), vpvsr, cvpvsr))
ax.set_title('Wadati-Diagram, %d S-P Times, Vp/Vs(raw)=%5.2f,'
'Vp/Vs(checked)=%5.2f' % (len(SPtimes), vpvsr, cvpvsr))
ax.legend(loc=1, numpoints=1)
else:
ax.set_title('Wadati-Diagram, %d S-P Times' % len(SPtimes))
@ -714,48 +721,51 @@ def wadaticheck(pickdic, dttolerance, iplot=0, fig_dict=None):
def RMS(X):
'''
Function returns root mean square of a given array X
'''
"""
Returns root mean square of a given array X
:param X: Array
:type X: `~numpy.ndarray`
:return: root mean square value of given array
:rtype: float
"""
return np.sqrt(np.sum(np.power(X, 2)) / len(X))
def checksignallength(X, pick, TSNR, minsiglength, nfac, minpercent, iplot=0, fig=None, linecolor='k'):
'''
"""
Function to detect spuriously picked noise peaks.
Uses RMS trace of all 3 components (if available) to determine,
how many samples [per cent] after P onset are below certain
threshold, calculated from noise level times noise factor.
: param: X, time series (seismogram)
: type: `~obspy.core.stream.Stream`
: param: pick, initial (AIC) P onset time
: type: float
: param: TSNR, length of time windows around initial pick [s]
: type: tuple (T_noise, T_gap, T_signal)
: param: minsiglength, minium required signal length [s] to
declare pick as P onset
: type: float
: param: nfac, noise factor (nfac * noise level = threshold)
: type: float
: param: minpercent, minimum required percentage of samples
above calculated threshold
: type: float
: param: iplot, if iplot > 1, results are shown in figure
: type: int
'''
:param X: time series (seismogram)
:type X: `~obspy.core.stream.Stream`
:param pick: initial (AIC) P onset time
:type pick: float
:param TSNR: length of time windows around initial pick [s]
:type TSNR: (T_noise, T_gap, T_signal)
:param minsiglength: minium required signal length [s] to declare pick as P onset
:type minsiglength: float
:param nfac: noise factor (nfac * noise level = threshold)
:type nfac: float
:param minpercent: minimum required percentage of samples above calculated threshold
:type minpercent: float
:param iplot: iplot, if iplot > 1, results are shown in figure
:type iplot: int
:param fig: Matplotlib figure to plot results in
:type fig: `~matplotlib.figure.Figure`
:param linecolor: color of seismic traces
:type linecolor: str
:return: flag, value of 1 if signal reached required length, 0 if signal is shorter than
required length
:rtype: int
"""
plt_flag = 0
try:
iplot = int(iplot)
except:
if iplot == True or iplot == 'True':
if real_Bool(iplot):
iplot = 2
else:
iplot = 0
@ -802,7 +812,7 @@ def checksignallength(X, pick, TSNR, minsiglength, nfac, minpercent, iplot=0, fi
returnflag = 0
if iplot > 1:
if fig == None or fig == 'None':
if real_None(fig) is None:
fig = plt.figure() # iplot)
plt_flag = 1
fig._tight = True
@ -828,21 +838,26 @@ def checksignallength(X, pick, TSNR, minsiglength, nfac, minpercent, iplot=0, fi
def checkPonsets(pickdic, dttolerance, jackfactor=5, iplot=0, fig_dict=None):
'''
"""
Function to check statistics of P-onset times: Control deviation from
median (maximum adjusted deviation = dttolerance) and apply pseudo-
bootstrapping jackknife.
: param: pickdic, dictionary containing picks and quality parameters
: type: dictionary
: param: dttolerance, maximum adjusted deviation of P-onset time from
median of all P onsets
: type: float
: param: iplot, if iplot > 1, Wadati diagram is shown
: type: int
'''
:param pickdic: dictionary containing picks and quality parameters
:type pickdic: dict
:param dttolerance: maximum adjusted deviation of P onset time from the
median of all P onsets
:type dttolerance: float
:param jackfactor: if pseudo value is larger than jackfactor * Jackknife estimator,
the distorts the estimator too much and will be removed
:type jackfactor: int
:param iplot: if iplot > 1, Wadati diagram is shown
:type iplot: int
:param fig_dict: Matplotlib figure used for plotting
:type fig_dict: `~matplotlib.figure.Figure`
:return: dictionary containing all onsets that passed the jackknife and the
median test
:rtype: dict
"""
checkedonsets = pickdic
@ -883,7 +898,7 @@ def checkPonsets(pickdic, dttolerance, jackfactor=5, iplot=0, fig_dict=None):
print("checkPonsets: %d pick(s) deviate too much from median!" % len(ibad))
print(badstations)
print("checkPonsets: Skipped %d P pick(s) out of %d" % (len(badstations) \
print("checkPonsets: Skipped %d P pick(s) out of %d" % (len(badstations)
+ len(badjkstations), len(stations)))
goodmarker = 'goodPonsetcheck'
@ -921,7 +936,7 @@ def checkPonsets(pickdic, dttolerance, jackfactor=5, iplot=0, fig_dict=None):
ax = fig.add_subplot(111)
if len(badstations) > 0:
ax.plot(ibad, np.array(Ppicks)[ibad], marker ='o', markerfacecolor='orange', markersize=14,
ax.plot(ibad, np.array(Ppicks)[ibad], marker='o', markerfacecolor='orange', markersize=14,
linestyle='None', label='Median Skipped P Picks')
if len(badjkstations) > 0:
ax.plot(badjk[0], np.array(Ppicks)[badjk], 'ro', markersize=14, label='Jackknife Skipped P Picks')
@ -934,7 +949,7 @@ def checkPonsets(pickdic, dttolerance, jackfactor=5, iplot=0, fig_dict=None):
for index, pick in enumerate(Ppicks):
ax.text(index, pick + 0.01, '{0}'.format(stations[index]), color='0.25')
ax.set_xlabel('Number of P Picks')
ax.set_ylabel('Onset Time [s] from 1.1.1970') # MP MP Improve this?
ax.set_ylabel('Onset Time [s] from 1.1.1970') # MP MP Improve this?
ax.legend(loc=1, numpoints=1)
ax.set_title('Jackknifing and Median Tests on P Onsets')
if plt_flag:
@ -943,24 +958,27 @@ def checkPonsets(pickdic, dttolerance, jackfactor=5, iplot=0, fig_dict=None):
return checkedonsets
def jackknife(X, phi, h):
'''
def jackknife(X, phi, h=1):
"""
Function to calculate the Jackknife Estimator for a given quantity,
special type of boot strapping. Returns the jackknife estimator PHI_jack
the pseudo values PHI_pseudo and the subgroup parameters PHI_sub.
special type of boot strapping.
: param: X, given quantity
: type: list
: param: phi, chosen estimator, choose between:
Returns the jackknife estimator PHI_jack the pseudo values PHI_pseudo
and the subgroup parameters PHI_sub.
:param X: list containing UTCDateTime objcects representing P onsets
:type X: list (`~obspy.core.utcdatetime.UTCDateTime`)
:param phi:phi, chosen estimator, choose between:
"MED" for median
"MEA" for arithmetic mean
"VAR" for variance
: type: string
: param: h, size of subgroups, optinal, default = 1
: type: integer
'''
:type phi: str
:param h: size of subgroups, optional (default = 1)
:type h: int
:return: Tuple containing the Jackknife estimator PHI_jack, a list of jackknife pseudo
values (PHI_pseudo) and the Jackknife estimators (PHI_sub) of the subgroups.
Will return (None, None, None) if X cannot be divided in h subgroups of equals size
:rtype: (float, list, list)
"""
PHI_jack = None
PHI_pseudo = None
@ -1008,40 +1026,42 @@ def jackknife(X, phi, h):
def checkZ4S(X, pick, zfac, checkwin, iplot, fig=None, linecolor='k'):
'''
"""
Function to compare energy content of vertical trace with
energy content of horizontal traces to detect spuriously
picked S onsets instead of P onsets. Usually, P coda shows
larger longitudal energy on vertical trace than on horizontal
traces, where the transversal energy is larger within S coda.
Be careful: there are special circumstances, where this is not
the case!
picked S onsets instead of P onsets.
: param: X, fitered(!) time series, three traces
: type: `~obspy.core.stream.Stream`
Usually, P coda shows larger longitudal energy on vertical trace
than on horizontal traces, where the transversal energy is larger
within S coda. Be careful: there are special circumstances, where
this is not the case!
To pass the test, vertical P-coda level must exceed horizontal P-coda level
zfac times EN-coda level
: param: pick, initial (AIC) P onset time
: type: float
: param: zfac, factor for threshold determination,
vertical energy must exceed coda level times zfac
to declare a pick as P onset
: type: float
: param: checkwin, window length [s] for calculating P-coda
energy content
: type: float
: param: iplot, if iplot > 1, energy content and threshold
are shown
: type: int
'''
:param X: fitered(!) time series, three traces
:type X: `~obspy.core.stream.Stream`
:param pick: initial (AIC) P onset time
:type pick: float
:param zfac: factor for threshold determination, vertical energy must
exceed coda level times zfac to declare a pick as P onset
:type zfac: float
:param checkwin: window length [s] for calculating P-coda engergy content
:type checkwin: float
:param iplot: if iplot > 1, energy content and threshold are shown
:type iplot: int
:param fig: Matplotlib figure to plot results in
:type fig: `~matplotlib.figure.Figure`
:param linecolor: color of seismic traces
:type linecolor: str
:return: returnflag; 0 if onset failed test, 1 if onset passed test
:rtype: int
"""
plt_flag = 0
try:
iplot = int(iplot)
except:
if iplot == True or iplot == 'True':
if real_Bool(iplot):
iplot = 2
else:
iplot = 0
@ -1123,14 +1143,14 @@ def checkZ4S(X, pick, zfac, checkwin, iplot, fig=None, linecolor='k'):
t = np.arange(diff_dict[key], trace.stats.npts / trace.stats.sampling_rate + diff_dict[key],
trace.stats.delta)
if i == 0:
if fig == None or fig == 'None':
if real_None(fig) is None:
fig = plt.figure() # self.iplot) ### WHY? MP MP
plt_flag = 1
ax1 = fig.add_subplot(3, 1, i + 1)
ax = ax1
ax.set_title('CheckZ4S, Station %s' % zdat[0].stats.station)
else:
if fig == None or fig == 'None':
if real_None(fig) is None:
fig = plt.figure() # self.iplot) ### WHY? MP MP
plt_flag = 1
ax = fig.add_subplot(3, 1, i + 1, sharex=ax1)
@ -1154,14 +1174,20 @@ def checkZ4S(X, pick, zfac, checkwin, iplot, fig=None, linecolor='k'):
def getQualityFromUncertainty(uncertainty, Errors):
'''Script to transform uncertainty into quality classes 0-4
regarding adjusted time errors Errors.
'''
"""
Script to transform uncertainty into quality classes 0-4 regarding adjusted time errors
:param uncertainty: symmetric picking error of picks
:type uncertainty: float
:param Errors: Width of uncertainty classes 0-4 in seconds
:type Errors: list
:return: quality of pick (0-4)
:rtype: int
"""
# set initial quality to 4 (worst) and change only if one condition is hit
quality = 4
if uncertainty == None or uncertainty == 'None':
if real_None(uncertainty) is None:
return quality
if uncertainty <= Errors[0]:

7
pylot/core/util/connection.py
View File

@ -8,6 +8,13 @@ except:
def checkurl(url='https://ariadne.geophysik.ruhr-uni-bochum.de/trac/PyLoT/'):
"""
check if URL is available
:param url: url
:type url: str
:return: available: True/False
:rtype: bool
"""
try:
urlopen(url, timeout=1)
return True

154
pylot/core/util/event.py
View File

@ -15,6 +15,11 @@ class Event(ObsPyEvent):
'''
def __init__(self, path):
"""
Initialize event by event directory
:param path: path to event directory
:type path: str
"""
self.pylot_id = path.split('/')[-1]
# initialize super class
super(Event, self).__init__(resource_id=ResourceIdentifier('smi:local/' + self.pylot_id))
@ -30,10 +35,20 @@ class Event(ObsPyEvent):
self.get_notes()
def get_notes_path(self):
"""
Notes files is freely editable by the user and can contain notes regarding the event
:return: path to notes file
:rtype: str
"""
notesfile = os.path.join(self.path, 'notes.txt')
return notesfile
def get_notes(self):
"""
set self.note attribute to content of notes file
:return:
:rtype: None
"""
notesfile = self.get_notes_path()
if os.path.isfile(notesfile):
with open(notesfile) as infile:
@ -48,34 +63,81 @@ class Event(ObsPyEvent):
pass
def addNotes(self, notes):
"""
Set new notes string
:param notes: notes to save in Event object
:type notes: str
:return:
:rtype: None
"""
self.notes = str(notes)
def clearNotes(self):
"""
Clear event notes
:return:
:rtype: None
"""
self.notes = None
def isRefEvent(self):
"""
Return reference event flag
:return: True if event is refence event
:rtype: bool
"""
return self._refEvent
def isTestEvent(self):
"""
Return test event flag
:return: True if event is test event
:rtype: bool
"""
return self._testEvent
def setRefEvent(self, bool):
"""
Set reference event flag
:param bool: new reference event flag
:type bool: bool
:return:
:rtype: None
"""
self._refEvent = bool
if bool: self._testEvent = False
def setTestEvent(self, bool):
"""
Set test event flag
:param bool: new test event flag
:type bool: bool
:return:
:rtype: None
"""
self._testEvent = bool
if bool: self._refEvent = False
def clearObsPyPicks(self, picktype):
"""
Remove picks of a certain type from event
:param picktype: type of picks to remove, 'auto' or 'manual'
:type picktype: str
:return:
:rtype: None
"""
for index, pick in reversed(list(enumerate(self.picks))):
if picktype in str(pick.method_id):
self.picks.pop(index)
def addPicks(self, picks):
'''
"""
add pylot picks and overwrite existing ones
'''
:param picks: picks to add to event in pick dictionary
:type picks: dict
:return:
:rtype: None
"""
for station in picks:
self.pylot_picks[station] = picks[station]
# add ObsPy picks (clear old manual and copy all new manual from pylot)
@ -83,6 +145,13 @@ class Event(ObsPyEvent):
self.picks += picks_from_picksdict(self.pylot_picks)
def addAutopicks(self, autopicks):
"""
Add automatic picks to event
:param autopicks: automatic picks to add to event
:type autopicks dict:
:return:
:rtype: None
"""
for station in autopicks:
self.pylot_autopicks[station] = autopicks[station]
# add ObsPy picks (clear old auto and copy all new auto from pylot)
@ -90,6 +159,15 @@ class Event(ObsPyEvent):
self.picks += picks_from_picksdict(self.pylot_autopicks)
def setPick(self, station, pick):
"""
Set pick for a station
:param station: station name
:type station: str
:param pick:
:type pick: dict
:return:
:rtype:
"""
if pick:
self.pylot_picks[station] = pick
else:
@ -101,21 +179,46 @@ class Event(ObsPyEvent):
self.picks += picks_from_picksdict(self.pylot_picks)
def setPicks(self, picks):
'''
set pylot picks and delete and overwrite all existing
'''
"""
Set pylot picks and delete and overwrite all existing
:param picks: new picks
:type picks: dict
:return:
:rtype: None
"""
self.pylot_picks = picks
self.clearObsPyPicks('manual')
self.picks += picks_from_picksdict(self.pylot_picks)
def getPick(self, station):
"""
Get pick at station
:param station: station name
:type station: str
:return: pick dictionary of station
:rtype: dict
"""
if station in self.pylot_picks.keys():
return self.pylot_picks[station]
def getPicks(self):
"""
Return pylot picks
:return:
:rtype: dict
"""
return self.pylot_picks
def setAutopick(self, station, pick):
"""
Set pick at station
:param station: station name
:type station: str
:param pick:
:type pick: dict
:return:
:rtype: None
"""
if pick:
self.pylot_autopicks[station] = pick
else:
@ -127,25 +230,46 @@ class Event(ObsPyEvent):
self.picks += picks_from_picksdict(self.pylot_autopicks)
def setAutopicks(self, picks):
'''
set pylot picks and delete and overwrite all existing
'''
"""
Set pylot picks and delete and overwrite all existing
:param picks: new picks
:type picks: dict
:return:
:rtype: None
"""
self.pylot_autopicks = picks
self.clearObsPyPicks('auto')
self.picks += picks_from_picksdict(self.pylot_autopicks)
def getAutopick(self, station):
"""
Return autopick at station
:param station: station name
:type station: str
:return: pick dictionary
:rtype: dict
"""
if station in self.pylot_autopicks.keys():
return self.pylot_autopicks[station]
def getAutopicks(self):
"""
Get autopicks of event
:return: dict containing automatic picks
:rtype: dict
"""
return self.pylot_autopicks
def save(self, filename):
'''
"""
Save PyLoT Event to a file.
Can be loaded by using event.load(filename).
'''
Uses pickling to save event object to file
:param filename: filename to save project under
:type filename: str
:return:
:rtype: None
"""
try:
import cPickle
except ImportError:
@ -159,9 +283,13 @@ class Event(ObsPyEvent):
@staticmethod
def load(filename):
'''
Load project from filename.
'''
"""
Load project from filename
:param filename: to load event file
:type filename: str
:return: event loaded from file
:rtype: Event
"""
try:
import cPickle
except ImportError:

2
pylot/core/util/thread.py
View File

@ -24,7 +24,7 @@ class Thread(QThread):
if self.redirect_stdout:
sys.stdout = self
try:
if self.arg:
if self.arg is not None:
self.data = self.func(self.arg)
else:
self.data = self.func()

450
pylot/core/util/utils.py
View File

@ -13,7 +13,7 @@ from obspy.core import AttribDict
from obspy.signal.rotate import rotate2zne
from obspy.io.xseed.utils import SEEDParserException
from pylot.core.io.inputs import PylotParameter
from pylot.core.io.inputs import PylotParameter, FilterOptions
from pylot.styles import style_settings
from scipy.interpolate import splrep, splev
@ -31,13 +31,36 @@ def _pickle_method(m):
else:
return getattr, (m.im_self, m.im_func.func_name)
def getAutoFilteroptions(phase, parameter):
filtername = {'P': 'bpz2',
'S': 'bph2'}
if not phase in filtername.keys():
print('autoPickParameter: No filter options for phase {}.'.format(phase))
return
freqmin, freqmax = parameter.get(filtername[phase])
filteroptions = FilterOptions(type='bandpass', freq=[freqmin, freqmax], order=4) # order=4 default from obspy
return filteroptions
def readDefaultFilterInformation(fname):
"""
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)
return readFilterInformation(pparam)
def readFilterInformation(pylot_parameter):
"""
Read filter information from PylotParameter object into a dictionary
:param pylot_parameter: PylotParameter object
:type pylot_parameter: `~pylot.pylot.core.io.inputs.PylotParameter`
:return: dictionary containing the filter information
:rtype: dict
"""
p_filter = {'filtertype': pylot_parameter['filter_type'][0],
'freq': [pylot_parameter['minfreq'][0], pylot_parameter['maxfreq'][0]],
'order': int(pylot_parameter['filter_order'][0])}
@ -50,23 +73,46 @@ def readFilterInformation(pylot_parameter):
def fit_curve(x, y):
"""
:param x: data points defining a curve y = f(x)
:type x: array_like
:param y: data points defining a curve y = f(x)
:type y: array_like
:return: tuple containing a function to evaluate a B-spline and
(the vector of knots, B-spline coefficients, degree of the spline)
:rtype: (func, (t, c, k))
"""
return splev, splrep(x, y)
def getindexbounds(f, eta):
mi = f.argmax()
m = max(f)
b = m * eta
l = find_nearest(f[:mi], b)
u = find_nearest(f[mi:], b) + mi
"""
Get indices of values closest below and above maximum value in an array
:param f: array
:type f: `~numpy.ndarray`
:param eta: look for value in array that is closes to max_value * eta
:type eta: float
:return: tuple containing index of max value, index of value closest below max value,
index of value closest above max value
:rtype: (int, int, int)
"""
mi = f.argmax() # get indices of max values
m = max(f) # get maximum value
b = m * eta #
l = find_nearest(f[:mi], b) # find closest value below max value
u = find_nearest(f[mi:], b) + mi # find closest value above max value
return mi, l, u
def gen_Pool(ncores=0):
'''
"""
Generate mulitprocessing pool object utilizing ncores amount of cores
:param ncores: number of CPU cores for multiprocessing.Pool, if ncores == 0 use all available
:type ncores: int
:return: multiprocessing.Pool object
'''
:rtype: `~multiprocessing.Pool`
"""
import multiprocessing
if ncores == 0:
@ -79,14 +125,19 @@ def gen_Pool(ncores=0):
def excludeQualityClasses(picks, qClasses, timeerrorsP, timeerrorsS):
'''
"""
takes PyLoT picks dictionary and returns a new dictionary with certain classes excluded.
:param picks: PyLoT picks dictionary
:type picks: dict
:param qClasses: list (or int) of quality classes (0-4) to exclude
:param timeerrorsP: time errors for classes (0-4) for P
:param timeerrorsS: time errors for classes (0-4) for S
:return: new picks dictionary
'''
:type qClasses: [int]
:param timeerrorsP: width of quality classes for P onsets in seconds
:type timeerrorsP: (float, float, float, float)
:param timeerrorsS: width of quality classes for S onsets in seconds
:type timeerrorsS: (float, float, float, float])
:return: dictionary containing only picks above the excluded quality class(es)
:rtype: dict
"""
from pylot.core.pick.utils import getQualityFromUncertainty
if type(qClasses) in [int, float]:
@ -114,9 +165,12 @@ def excludeQualityClasses(picks, qClasses, timeerrorsP, timeerrorsS):
def clims(lim1, lim2):
"""
takes two pairs of limits and returns one pair of common limts
:param lim1:
:param lim2:
:return:
:param lim1: limit 1
:type lim1: int
:param lim2: limit 2
:type lim2: int
:return: new upper and lower limit common to both given limits
:rtype: [int, int]
>>> clims([0, 4], [1, 3])
[0, 4]
@ -155,7 +209,7 @@ def demeanTrace(trace, window):
demeaned within a certain time window
:param trace: waveform trace object
:type trace: `~obspy.core.stream.Trace`
:param window:
:param window: time window whitin which data is demeaned
:type window: tuple
:return: trace
:rtype: `~obspy.core.stream.Trace`
@ -182,8 +236,11 @@ def find_in_list(list, str):
takes a list of strings and a string and returns the first list item
matching the string pattern
:param list: list to search in
:type list: list
:param str: pattern to search for
:type str: str
:return: first list item containing pattern
:rtype: str
.. example::
@ -208,13 +265,15 @@ def find_in_list(list, str):
def find_nearest(array, value):
'''
"""
function find_nearest takes an array and a value and returns the
index of the nearest value found in the array
:param array: array containing values
:type array: `~numpy.ndarray`
:param value: number searched for
:type value: float
:return: index of the array item being nearest to the value
:rtype: int
>>> a = np.array([ 1.80339578, -0.72546654, 0.95769195, -0.98320759, 0.85922623])
>>> find_nearest(a, 1.3)
@ -228,17 +287,18 @@ def find_nearest(array, value):
>>> a = np.array([ 1.1, -0.7, 0.9, -0.9, 0.8])
>>> find_nearest(a, 0.849)
4
'''
"""
return (np.abs(array - value)).argmin()
def fnConstructor(s):
'''
"""
takes a string and returns a valid filename (especially on windows machines)
:param s: desired filename
:type s: str
:return: valid filename
'''
:rtype: str
"""
if type(s) is str:
s = s.split(':')[-1]
else:
@ -255,6 +315,13 @@ def fnConstructor(s):
def real_None(value):
"""
Convert "None" to None
:param value:
:type value: str, bool
:return:
:rtype: bool
"""
if value == 'None':
return None
else:
@ -262,9 +329,16 @@ def real_None(value):
def real_Bool(value):
if value == 'True':
"""
Convert string representations of bools to their true boolean value
:param value:
:type value: str, bool
:return: true boolean value
:rtype: bool
"""
if value in ['True', 'true']:
return True
elif value == 'False':
elif value in ['False', 'false']:
return False
else:
return value
@ -276,7 +350,8 @@ def four_digits(year):
from the last 100 years
:param year: two digit year
:type year: int
:return: four digit year correspondant
:return: four digit year correspondent
:rtype: int
>>> four_digits(20)
1920
@ -293,13 +368,13 @@ def four_digits(year):
def common_range(stream):
'''
takes a stream object and returns the earliest end and the latest start
time of all contained trace objects
:param stream: seismological data stream
"""
takes a stream object and returns the earliest end and the latest start time of all contained trace objects
:param stream: seismological data stream
:type stream: `~obspy.core.stream.Stream`
:return: maximum start time and minimum end time
'''
:rtype: (`~maximum start time and minimum end time`, maximum start time and minimum end time)
"""
max_start = None
min_end = None
for trace in stream:
@ -311,51 +386,82 @@ def common_range(stream):
def full_range(stream):
'''
"""
takes a stream object and returns the latest end and the earliest start
time of all contained trace objects
:param stream: seismological data stream
:type stream: `~obspy.core.stream.Stream`
:return: minimum start time and maximum end time
'''
min_start = UTCDateTime()
max_end = None
for trace in stream:
if trace.stats.starttime < min_start:
min_start = trace.stats.starttime
if max_end is None or trace.stats.endtime > max_end:
max_end = trace.stats.endtime
:rtype: (`~maximum start time and minimum end time`, maximum start time and minimum end time)
"""
min_start = min([trace.stats.starttime for trace in stream])
max_end = max([trace.stats.endtime for trace in stream])
return min_start, max_end
def transformFilteroptions2String(filtopts):
st = ''
if not filtopts:
return st
if 'type' in filtopts.keys():
st += '{}'.format(filtopts['type'])
if 'freq' in filtopts.keys():
st += ' | freq: {}'.format(filtopts['freq'])
elif 'freqmin' in filtopts.keys() and 'freqmax' in filtopts.keys():
st += ' | freqmin: {} | freqmax: {}'.format(filtopts['freqmin'], filtopts['freqmax'])
for key, value in filtopts.items():
if key in ['type', 'freq', 'freqmin', 'freqmax']:
continue
st += ' | {}: {}'.format(key, value)
return st
def transformFilterString4Export(st):
st = st.replace('|', '//')
st = st.replace(':', '/')
st = st.replace(' ', '')
return st
def backtransformFilterString(st):
st = st.split('smi:local/')
st = st[1] if len(st) > 1 else st[0]
st = st.replace('//', ' | ')
st = st.replace('/', ': ')
return st
def getHash(time):
'''
takes a time object and returns the corresponding SHA1 hash of the
formatted date string
"""
takes a time object and returns the corresponding SHA1 hash of the formatted date string
:param time: time object for which a hash should be calculated
:type time: :class: `~obspy.core.utcdatetime.UTCDateTime` object
:return: str
'''
:type time: `~obspy.core.utcdatetime.UTCDateTime`
:return: SHA1 hash
:rtype: str
"""
hg = hashlib.sha1()
hg.update(time.strftime('%Y-%m-%d %H:%M:%S.%f'))
return hg.hexdigest()
def getLogin():
'''
"""
returns the actual user's login ID
:return: login ID
'''
:rtype: str
"""
return os.getlogin()
def getOwner(fn):
'''
"""
takes a filename and return the login ID of the actual owner of the file
:param fn: filename of the file tested
:type fn: str
:return: login ID of the file's owner
'''
:rtype: str
"""
system_name = platform.system()
if system_name in ["Linux", "Darwin"]:
import pwd
@ -376,6 +482,8 @@ def getPatternLine(fn, pattern):
:param pattern: pattern string to search for
:type pattern: str
:return: the complete line containing the pattern string or None
:rtype: int, None
>>> getPatternLine('utils.py', 'python')
'#!/usr/bin/env python\\n'
@ -397,18 +505,20 @@ def is_executable(fn):
and False otherwise
:param fn: path to the file to be tested
:return: True or False
:rtype: bool
"""
return os.path.isfile(fn) and os.access(fn, os.X_OK)
def isSorted(iterable):
'''
takes an iterable and returns 'True' if the items are in order otherwise
'False'
"""
takes an iterable and returns True if the items are in order otherwise False
:param iterable: an iterable object
:type iterable:
:return: Boolean
:rtype: bool
..example::
>>> isSorted(1)
Traceback (most recent call last):
...
@ -421,7 +531,7 @@ def isSorted(iterable):
False
>>> isSorted([2,3,1,4])
False
'''
"""
assert isIterable(iterable), 'object is not iterable; object: {' \
'0}'.format(iterable)
if type(iterable) is str:
@ -431,10 +541,12 @@ def isSorted(iterable):
def isIterable(obj):
"""
takes a python object and returns 'True' is the object is iterable and
'False' otherwise
takes a python object and returns True is the object is iterable and
False otherwise
:param obj: a python object
:type obj: object
:return: True of False
:rtype: bool
"""
try:
iterator = iter(obj)
@ -451,6 +563,7 @@ def key_for_set_value(d):
:type d: dict
:return: key to the first non-False value found; None if no value's
boolean equals True
:rtype:
"""
r = None
for k, v in d.items():
@ -460,13 +573,18 @@ def key_for_set_value(d):
def prepTimeAxis(stime, trace, verbosity=0):
'''
"""
takes a starttime and a trace object and returns a valid time axis for
plotting
:param stime: start time of the actual seismogram as UTCDateTime
:type stime: `~obspy.core.utcdatetime.UTCDateTime`
:param trace: seismic trace object
:type trace: `~obspy.core.trace.Trace`
:param verbosity: if != 0, debug output will be written to console
:type verbosity: int
:return: valid numpy array with time stamps for plotting
'''
:rtype: `~numpy.ndarray`
"""
nsamp = trace.stats.npts
srate = trace.stats.sampling_rate
tincr = trace.stats.delta
@ -513,6 +631,20 @@ def find_horizontals(data):
def make_pen(picktype, phase, key, quality):
"""
Make PyQtGraph.QPen
:param picktype: 'manual' or 'automatic'
:type picktype: str
:param phase: 'P' or 'S'
:type phase: str
:param key: 'mpp', 'epp', 'lpp' or 'spe', (earliest/latest possible pick, symmetric picking error or
most probable pick)
:type key: str
:param quality: quality class of pick, decides color modifier
:type quality: int
:return: PyQtGraph QPen
:rtype: `~QPen`
"""
if pg:
rgba = pick_color(picktype, phase, quality)
linestyle, width = pick_linestyle_pg(picktype, key)
@ -521,8 +653,22 @@ def make_pen(picktype, phase, key, quality):
def pick_color(picktype, phase, quality=0):
"""
Create pick color by modifying the base color by the quality.
Returns rgba values in a range of [0, 255]. picktype, phase decide the base color,
quality decides the applied modifier
:param picktype: 'manual' or 'automatic'
:type picktype: str
:param phase: 'P' or 'S'
:type phase: str
:param quality: quality of pick. Decides the new intensity of the modifier color
:type quality: int
:return: tuple containing modified rgba color values
:rtype: (int, int, int, int)
"""
min_quality = 3
bpc = base_phase_colors(picktype, phase)
bpc = base_phase_colors(picktype, phase) # returns dict like {'modifier': 'g', 'rgba': (0, 0, 255, 255)}
rgba = bpc['rgba']
modifier = bpc['modifier']
intensity = 255.*quality/min_quality
@ -531,6 +677,20 @@ def pick_color(picktype, phase, quality=0):
def pick_color_plt(picktype, phase, quality=0):
"""
Create pick color by modifying the base color by the quality.
Returns rgba values in a range of [0, 1]. picktype, phase decide the base color,
quality decides the applied modifier
:param picktype: 'manual' or 'automatic'
:type picktype: str
:param phase: 'P' or 'S'
:type phase: str
:param quality: quality of pick. Decides the new intensity of the modifier color
:type quality: int
:return: tuple containing rgba values matplotlib style, ranging from [0, 1]
:rtype: (float, float, float, float)
"""
rgba = list(pick_color(picktype, phase, quality))
for index, val in enumerate(rgba):
rgba[index] /= 255.
@ -538,6 +698,16 @@ def pick_color_plt(picktype, phase, quality=0):
def pick_linestyle_plt(picktype, key):
"""
Get matplotlib line style for plotting by pick type and pick parameter (earliest/latest possible pick,
symmetric picking error or most probable pick).
:param picktype: 'manual' or 'automatic'
:type picktype: str
:param key: which pick parameter should be plotted, 'mpp', 'epp', 'lpp' or 'spe'
:type key: str
:return: tuple containing matplotlib line style string and line thicknes
:rtype: (str, float)
"""
linestyles_manu = {'mpp': ('solid', 2.),
'epp': ('dashed', 1.),
'lpp': ('dashed', 1.),
@ -552,6 +722,16 @@ def pick_linestyle_plt(picktype, key):
def pick_linestyle_pg(picktype, key):
"""
Get Qt line style by picktype and pick parameter (earliest/latest possible pick, symmetric picking error or
most probable pick)
:param picktype: 'manual' or 'automatic'
:type picktype: str
:param key: which pick parameter should be plotted, 'mpp', 'epp', 'lpp' or 'spe'
:type key: str
:return: Qt line style parameters
:rtype:
"""
linestyles_manu = {'mpp': (QtCore.Qt.SolidLine, 2.),
'epp': (QtCore.Qt.DashLine, 1.),
'lpp': (QtCore.Qt.DashLine, 1.),
@ -566,6 +746,17 @@ def pick_linestyle_pg(picktype, key):
def modify_rgba(rgba, modifier, intensity):
"""
Modify rgba color by adding the given intensity to the modifier color
:param rgba: tuple containing rgba values
:type rgba: (int, int, int, int)
:param modifier: which color should be modified, eg. 'r', 'g', 'b'
:type modifier: str
:param intensity: intensity to be added to selected color
:type intensity: float
:return: tuple containing rgba values
:rtype: (int, int, int, int)
"""
rgba = list(rgba)
index = {'r': 0,
'g': 1,
@ -580,10 +771,28 @@ def modify_rgba(rgba, modifier, intensity):
def base_phase_colors(picktype, phase):
"""
Get base color for a phase from style settings
:param picktype: 'manual' or 'automatic' picks
:type picktype: str
:param phase: Phase to select color for, 'P' or 'S'
:type phase: str
:return: dictionary {'modifier': 'g', 'rgba': (0, 0, 255, 255)}
:rtype: dict
"""
phasecolors = style_settings.phasecolors
return phasecolors[picktype][phase]
def transform_colors_mpl_str(colors, no_alpha=False):
"""
Transforms rgba color values to a matplotlib string of color values with a range of [0, 1]
:param colors: tuple of rgba color values ranging from [0, 255]
:type colors: (float, float, float, float)
:param no_alpha: Wether to return a alpha value in the matplotlib color string
:type no_alpha: bool
:return: String containing r, g, b values and alpha value if no_alpha is False (default)
:rtype: str
"""
colors = list(colors)
colors_mpl = tuple([color / 255. for color in colors])
if no_alpha:
@ -593,6 +802,13 @@ def transform_colors_mpl_str(colors, no_alpha=False):
return colors_mpl
def transform_colors_mpl(colors):
"""
Transform rgba colors from [0, 255] to [0, 1]
:param colors: tuple of rgba color values ranging from [0, 255]
:type colors: (float, float, float, float)
:return: tuple of rgba color values ranging from [0, 1]
:rtype: (float, float, float, float)
"""
colors = list(colors)
colors_mpl = tuple([color / 255. for color in colors])
return colors_mpl
@ -600,10 +816,11 @@ def transform_colors_mpl(colors):
def remove_underscores(data):
"""
takes a `obspy.core.stream.Stream` object and removes all underscores
from stationnames
from station names
:param data: stream of seismic data
:type data: `obspy.core.stream.Stream`
:type data: `~obspy.core.stream.Stream`
:return: data stream
:rtype: `~obspy.core.stream.Stream`
"""
for tr in data:
# remove underscores
@ -612,16 +829,17 @@ def remove_underscores(data):
def trim_station_components(data, trim_start=True, trim_end=True):
'''
"""
cut a stream so only the part common to all three traces is kept to avoid dealing with offsets
:param data: stream of seismic data
:type data: `obspy.core.stream.Stream`
:type data: `~obspy.core.stream.Stream`
:param trim_start: trim start of stream
:type trim_start: bool
:param trim_end: trim end of stream
:type trim_end: bool
:return: data stream
'''
:rtype: `~obspy.core.stream.Stream`
"""
starttime = {False: None}
endtime = {False: None}
@ -639,11 +857,13 @@ def trim_station_components(data, trim_start=True, trim_end=True):
def check4gaps(data):
'''
"""
check for gaps in Stream and remove them
:param data: stream of seismic data
:type data: `~obspy.core.stream.Stream`
:return: data stream
'''
:rtype: `~obspy.core.stream.Stream`
"""
stations = get_stations(data)
for station in stations:
@ -657,11 +877,13 @@ def check4gaps(data):
def check4doubled(data):
'''
"""
check for doubled stations for same channel in Stream and take only the first one
:param data: stream of seismic data
:type data: `~obspy.core.stream.Stream`
:return: data stream
'''
:rtype: `~obspy.core.stream.Stream`
"""
stations = get_stations(data)
for station in stations:
@ -682,6 +904,13 @@ def check4doubled(data):
def get_stations(data):
"""
Get list of all station names in data stream
:param data: stream containing seismic traces
:type data: `~obspy.core.stream.Stream`
:return: list of all station names in data, no duplicates
:rtype: list of str
"""
stations = []
for tr in data:
station = tr.stats.station
@ -692,14 +921,35 @@ def get_stations(data):
def check4rotated(data, metadata=None, verbosity=1):
"""
:param data: stream object containing seismic traces
:type data: `~obspy.core.stream.Stream`
:param metadata: tuple containing metadata type string and metadata parser object
:type metadata: (str, `~obspy.io.xseed.parser.Parser`)
:param verbosity: if 0 print no information at runtime
:type verbosity: int
:return: stream object with traditionally oriented traces (ZNE) for stations that had misaligned traces (123) before
:rtype: `~obspy.core.stream.Stream`
"""
def rotate_components(wfstream, metadata=None):
"""rotates components if orientation code is numeric.
azimut and dip are fetched from metadata"""
"""
Rotate components if orientation code is numeric (= non traditional orientation).
Azimut and dip are fetched from metadata. To be rotated, traces of a station have to be cut to the same length.
Returns unrotated traces of no metadata is provided
:param wfstream: stream containing seismic traces
:type wfstream: `~obspy.core.stream.Stream`
:param metadata: tuple containing metadata type string and metadata parser object
:type metadata: (str, `~obspy.io.xseed.parser.Parser`)
:return: stream object with traditionally oriented traces (ZNE)
:rtype: `~obspy.core.stream.Stream`
"""
try:
# indexing fails if metadata is None
metadata[0]
except:
except TypeError:
if verbosity:
msg = 'Warning: could not rotate traces since no metadata was given\nset Inventory file!'
print(msg)
@ -714,7 +964,15 @@ def check4rotated(data, metadata=None, verbosity=1):
parser = metadata[1]
def get_dip_azimut(parser, trace_id):
"""gets azimut and dip for a trace out of the metadata parser"""
"""
Gets azimuth and dip by trace id out of the metadata parser
:param parser: metadata parser object
:type parser: `~obspy.io.xseed.parser.Parser`
:param trace_id: eg. 'BW.RJOB..EHZ',
:type trace_id: str
:return: tuple containing dip and azimuth of the trace corresponding to trace_id
:rtype: (float, float)
"""
dip = None
azimut = None
try:
@ -728,14 +986,14 @@ def check4rotated(data, metadata=None, verbosity=1):
dip = blockette_.dip
azimut = blockette_.azimuth
break
if dip is None or azimut is None:
if (dip is None or azimut is None) or (dip == 0 and azimut == 0):
error_msg = 'Dip and azimuth not available for trace_id {}'.format(trace_id)
raise ValueError(error_msg)
return dip, azimut
trace_ids = [trace.id for trace in wfstream]
for trace_id in trace_ids:
orientation = trace_id[-1]
orientation = trace_id[-1] # last letter if trace id is orientation code, ZNE or 123
if orientation.isnumeric():
# misaligned channels have a number as orientation
azimuts = []
@ -752,10 +1010,10 @@ def check4rotated(data, metadata=None, verbosity=1):
# to rotate all traces must have same length
wfstream = trim_station_components(wfstream, trim_start=True, trim_end=True)
z, n, e = rotate2zne(wfstream[0], azimuts[0], dips[0],
wfstream[1], azimuts[1], dips[1],
wfstream[2], azimuts[2], dips[2])
wfstream[1], azimuts[1], dips[1],
wfstream[2], azimuts[2], dips[2])
print('check4rotated: rotated station {} to ZNE'.format(trace_id))
z_index = dips.index(min(dips)) # get z-trace index (dip is measured from 0 to -90
z_index = dips.index(min(dips)) # get z-trace index (dip is measured from 0 to -90)
wfstream[z_index].data = z
wfstream[z_index].stats.channel = wfstream[z_index].stats.channel[0:-1] + 'Z'
del trace_ids[z_index]
@ -775,9 +1033,9 @@ def check4rotated(data, metadata=None, verbosity=1):
stations = get_stations(data)
for station in stations:
for station in stations: # loop through all stations and rotate data if neccessary
wf_station = data.select(station=station)
wf_station = rotate_components(wf_station, metadata)
rotate_components(wf_station, metadata)
return data
@ -837,7 +1095,9 @@ def which(program, infile=None):
takes a program name and returns the full path to the executable or None
modified after: http://stackoverflow.com/questions/377017/test-if-executable-exists-in-python
:param program: name of the desired external program
:type program: str
:return: full path of the executable file
:rtype: str
"""
try:
from PySide.QtCore import QSettings
@ -880,14 +1140,16 @@ def which(program, infile=None):
def loopIdentifyPhase(phase):
'''
"""
Loop through phase string and try to recognize its type (P or S wave).
Global variable ALTSUFFIX gives alternative suffix for phases if they do not end with P, p or S, s.
If ALTSUFFIX is not given, the function will cut the last letter of the phase string until string ends
with P or S.
:param phase: phase name (str)
:return:
'''
:param phase: phase name
:type phase: str
:return: str of phase ending with identified type, None if phase could not be identified
:rtype: str or None
"""
from pylot.core.util.defaults import ALTSUFFIX
phase_copy = phase
@ -906,14 +1168,18 @@ def loopIdentifyPhase(phase):
def identifyPhase(phase):
'''
"""
Returns capital P or S if phase string is identified by last letter. Else returns False.
:param phase: phase name (str)
:param phase: phase name
:type phase: str
:return: 'P', 'S' or False
'''
:rtype: str or bool
"""
# common phase suffix for P and S
common_P = ['P', 'p']
common_P = ['P', 'p', 'R']
common_S = ['S', 's']
if phase is None:
return False
if phase[-1] in common_P:
return 'P'
if phase[-1] in common_S:
@ -923,10 +1189,24 @@ def identifyPhase(phase):
def identifyPhaseID(phase):
"""
Returns phase id (capital P or S)
:param phase: phase name
:type phase: str
:return: phase type string
:rtype: str
"""
return identifyPhase(loopIdentifyPhase(phase))
def has_spe(pick):
"""
Check for 'spe' key (symmetric picking error) in dict and return its value if found, else return None
:param pick: pick dictionary
:type pick: dict
:return: value of 'spe' key
:rtype: float or None
"""
if not 'spe' in pick.keys():
return None
else:

1119
pylot/core/util/widgets.py
View File

File diff suppressed because it is too large Load Diff

1
pylot/styles/bright.qss
View File

@ -124,7 +124,6 @@ color:rgba(0, 0, 0, 255);
border-style: outset;
border-width: 1px;
border-color: rgba(100, 100, 120, 255);
min-width: 6em;
padding: 4px;
padding-left:5px;
padding-right:5px;

1
pylot/styles/dark.qss
View File

@ -123,7 +123,6 @@ color:rgba(255, 255, 255, 255);
border-style: outset;
border-width: 2px;
border-color: rgba(50, 50, 60, 255);
min-width: 6em;
padding: 4px;
padding-left:5px;
padding-right:5px;