[bugfix] fixes and modifications for (inhomogeneous) global datasets, WIP

This commit is contained in:
2017-05-30 16:09:48 +02:00
parent ceed165307
commit f493359c26
4 changed files with 66 additions and 30 deletions

View File

@@ -1 +1 @@
a745-dirty
ceed1-dirty

View File

@@ -343,6 +343,10 @@ def getSNR(X, TSNR, t1, tracenum=0):
assert isinstance(X, Stream), "%s is not a stream object" % str(X)
SNR = None
SNRdb = None
noiselevel = None
x = X[tracenum].data
npts = X[tracenum].stats.npts
sr = X[tracenum].stats.sampling_rate
@@ -356,10 +360,7 @@ def getSNR(X, TSNR, t1, tracenum=0):
isignal = getsignalwin(t, t1, TSNR[2])
if np.size(inoise) < 1:
print ("getSNR: Empty array inoise, check noise window!")
return
elif np.size(isignal) < 1:
print ("getSNR: Empty array isignal, check signal window!")
return
return SNR, SNRdB, noiselevel
# demean over entire waveform
x = x - np.mean(x[inoise])
@@ -368,6 +369,10 @@ def getSNR(X, TSNR, t1, tracenum=0):
noiselevel = np.sqrt(np.mean(np.square(x[inoise])))
#signallevel = np.sqrt(np.mean(np.square(x[isignal])))
if np.size(isignal) < 1:
print ("getSNR: Empty array isignal, check signal window!")
return SNR, SNRdB, noiselevel
#noiselevel = np.abs(x[inoise]).max()
signallevel = np.abs(x[isignal]).max()
@@ -431,7 +436,7 @@ def getsignalwin(t, t1, tsignal):
return isignal
def getResolutionWindow(snr):
def getResolutionWindow(snr, aperture):
"""
Number -> Float
produce the half of the time resolution window width from given SNR
@@ -453,18 +458,24 @@ def getResolutionWindow(snr):
>>> getResolutionWindow(2)
2.5
"""
res_wins = {'HRW': 2., 'MRW': 5., 'LRW': 10., 'VLRW': 15.}
res_wins = {
'regional': {'HRW': 2., 'MRW': 5., 'LRW': 10., 'VLRW': 15.},
'local': {'HRW': 2., 'MRW': 5., 'LRW': 10., 'VLRW': 15.},
'global': {'HRW': 40., 'MRW': 100., 'LRW': 200., 'VLRW': 300.}
}
if snr < 1.5:
time_resolution = res_wins['VLRW']
time_resolution = res_wins[aperture]['VLRW']
elif snr < 2.:
time_resolution = res_wins['LRW']
time_resolution = res_wins[aperture]['LRW']
elif snr < 3.:
time_resolution = res_wins['MRW']
time_resolution = res_wins[aperture]['MRW']
elif snr >3.:
time_resolution = res_wins[aperture]['HRW']
else:
time_resolution = res_wins['HRW']
time_resolution = res_wins[aperture]['VLRW']
return time_resolution / 2

View File

@@ -431,13 +431,19 @@ class WaveformWidget(FigureCanvas):
def plotWFData(self, wfdata, title=None, zoomx=None, zoomy=None,
noiselevel=None, scaleddata=False, mapping=True,
nth_sample=1):
component='*', nth_sample=1):
self.getAxes().cla()
self.clearPlotDict()
wfstart, wfend = full_range(wfdata)
nmax = 0
compclass = SetChannelComponents()
alter_comp = compclass.getCompPosition(component)
alter_comp = str(alter_comp[0])
wfdata = wfdata.select(component=component)
wfdata += wfdata.select(component=alter_comp)
# list containing tuples of network, station, channel (for sorting)
nsc = []
for trace in wfdata:
@@ -970,9 +976,13 @@ class PickDlg(QDialog):
result = getSNR(wfdata, (noise_win, gap_win, signal_win), ini_pick, itrace)
snr = result[0]
noiselevel = result[2] * nfac
noiselevel = result[2]
if noiselevel:
noiselevel *= nfac
else:
noiselevel = nfac
x_res = getResolutionWindow(snr)
x_res = getResolutionWindow(snr, 'regional')
# remove mean noise level from waveforms
for trace in data:
@@ -1015,7 +1025,12 @@ class PickDlg(QDialog):
# determine SNR and noiselevel
result = getSNR(wfdata, (noise_win, gap_win, signal_win), ini_pick)
snr = result[0]
noiselevel = result[2] * nfac
noiselevel = result[2]
if noiselevel:
noiselevel *= nfac
else:
noiselevel = nfac
# prepare plotting of data
for trace in data:
@@ -1027,7 +1042,7 @@ class PickDlg(QDialog):
horiz_comp = find_horizontals(data)
data = scaleWFData(data, noiselevel * 2.5, horiz_comp)
x_res = getResolutionWindow(snr)
x_res = getResolutionWindow(snr, 'regional')
self.setXLims(tuple([ini_pick - x_res, ini_pick + x_res]))
traces = self.getTraceID(horiz_comp)