Commit 2701c51c authored by Semih Demir's avatar Semih Demir
Browse files

delete double events and plotting solution for electronic interference

parent f39d44f2
......@@ -23,42 +23,106 @@ def dug_trigger(sta_total, tparam, event_nr, event_nr_s):
t = 1 / sta_total[0].stats.sampling_rate
trig = coincidence_trigger(tparam['algorithm'], tparam['sta_lta']['threshold_on'],
tparam['sta_lta']['threshold_off'], sta_total, tparam['coincidence'],
sta=tparam['sta_lta']['st_window'] * t, lta=tparam['sta_lta']['lt_window'] * t,
sta=tparam['sta_lta']['st_window'] * t, lta=tparam['sta_lta']['lt_window'] * t,trigger_off_extension=0.01,
details=True)
# print(trig)
time_s = []
classification_s = []
# coins_sum_s=[]
time_min_s = []
for f in range(len(trig)):
event_nr += 1 # update cummulative event_id
event_nr_s = event_nr_s + [event_nr] # create a list of event_ids
event_nr += 1 # update cummulative event_id
event_nr_s = event_nr_s + [event_nr] # create a list of event_ids
coins_sum = int(trig[f]['coincidence_sum'])
time_min = trig[f]['time'][0] # trigger time of earliest arrival of one event
time_max = max(trig[f]['time']) # trigger time of latest arrival of one event
time_min = trig[f]['time'][0] # trigger time of earliest arrival of one event
time_min_s = time_min_s + [trig[f]['time'][0]] # list of trigger time of earliest arrival of one event
time_max = max(trig[f]['time']) # trigger time of latest arrival of one event
diff = time_max-time_min # difference between earliest and latest arrival of one event
trace_id = trig[f]['trace_ids'] # find out trace ids of one event
diff = time_max - time_min # difference between earliest and latest arrival of one event
trace_id = trig[f]['trace_ids'] # find out trace ids of one event
id_int = [int(i[3:6]) for i in trace_id]
id_int.sort() # sort the trace ids of one event
id_int.sort() # sort the trace ids of one event
# if the difference between earliest and latest arrival is too small, classification=electronic interference
if diff < float(tparam['classification']['spread_int']):
cla_s = 'electronic'
elif id_int[0] == 1: # if the event has an arrival on the first trace, classification=active
elif id_int[0] == 1: # if the event has an arrival on the first trace, classification=active
cla_s = 'active'
else:
cla_s = 'passive' # if none of the above are true, classification=passive
cla_s = 'passive' # if none of the above are true, classification=passive
classification_s.append(cla_s) # make a list of classifications for each event
time_s.append(UTCDateTime(time_min.isoformat())) # make a list of earliest arrival times for each event
diff = []
delete = []
for i in range(1, len(trig)):
diff = time_min_s[i] - time_min_s[i - 1]
if diff < tparam['time_between_events']:
delete = delete + [i]
if len(delete):
n_del = 0
delete.reverse()
for i in delete:
del trig[i]
del event_nr_s[i]
event_nr_s[i:] = [x - 1 for x in event_nr_s[i:]]
del time_min_s[i]
del time_s[i]
del classification_s[i]
n_del = n_del + 1
coins_sum = coins_sum - n_del
event_nr = event_nr - n_del
for f in range(len(trig)):
# write a log file
data = [event_nr] + [coins_sum] + [time_min.isoformat()] + [cla_s]
data = [event_nr_s[f]] + [coins_sum] + [time_min_s[f].isoformat()] + [classification_s[f]]
cols = pd.Index(['Event_id', 'Coincidence_sum', 'Time', 'Classification'], name='cols')
df = pd.DataFrame(data=[data], columns=cols)
log_file = 'trigger.csv'
df.to_csv(log_file, mode='a', header=False, index=False)
classification_s.append(cla_s) # make a list of classifications for each event
time_s.append(UTCDateTime(time_min.isoformat())) # make a list of earliest arrival times for each event
# set up data frame containing list of event ids, trigger times and classifications
trigger_out = pd.DataFrame({'Event_id': event_nr_s, 'Time': time_s, 'Classification': classification_s})
return trigger_out, event_nr
return trigger_out, event_nr
\ No newline at end of file
# # print(trig)
# time_s = []
# classification_s = []
# for f in range(len(trig)):
# event_nr += 1 # update cummulative event_id
# event_nr_s = event_nr_s + [event_nr] # create a list of event_ids
# coins_sum = int(trig[f]['coincidence_sum'])
#
# time_min = trig[f]['time'][0] # trigger time of earliest arrival of one event
# time_max = max(trig[f]['time']) # trigger time of latest arrival of one event
#
# diff = time_max-time_min # difference between earliest and latest arrival of one event
# trace_id = trig[f]['trace_ids'] # find out trace ids of one event
# id_int = [int(i[3:6]) for i in trace_id]
# id_int.sort() # sort the trace ids of one event
# # if the difference between earliest and latest arrival is too small, classification=electronic interference
# if diff < float(tparam['classification']['spread_int']):
# cla_s = 'electronic'
# elif id_int[0] == 1: # if the event has an arrival on the first trace, classification=active
# cla_s = 'active'
# else:
# cla_s = 'passive' # if none of the above are true, classification=passive
#
# # write a log file
# data = [event_nr] + [coins_sum] + [time_min.isoformat()] + [cla_s]
# cols = pd.Index(['Event_id', 'Coincidence_sum', 'Time', 'Classification'], name='cols')
# df = pd.DataFrame(data=[data], columns=cols)
# log_file = 'trigger.csv'
# df.to_csv(log_file, mode='a', header=False, index=False)
#
# classification_s.append(cla_s) # make a list of classifications for each event
# time_s.append(UTCDateTime(time_min.isoformat())) # make a list of earliest arrival times for each event
#
# # set up data frame containing list of event ids, trigger times and classifications
# trigger_out = pd.DataFrame({'Event_id': event_nr_s, 'Time': time_s, 'Classification': classification_s})
#
# return trigger_out, event_nr
\ No newline at end of file
......@@ -258,29 +258,55 @@ class Event(ObsPyEvent):
if i < len(self.wf_stream) - 1:
axs[i].set_xticklabels([])
if self.classification != 'noise':
if self.classification != 'noise' and self.classification != 'electronic':
trig_ch = [int(i.waveform_id['station_code']) for i in self.picks]
ch_in = [int(i.stats['station']) for i in self.wf_stream.traces]
# finding the correct axes(stations)
# print(self.loc_ind)
ind_trig = np.where(np.isin(ch_in, trig_ch)) # Finds correct axis
trig_time = [(i.time - self.wf_stream.traces[0].stats["starttime"]) * 1000 for i in self.picks]
loc_ind_trig = ind_trig[0] #why was this? if this line gives error switch between ind_trig[0] and self.loc_ind[-1][:]
# plotting pick time
for m in range(len(trig_time)):
axs[ind_trig[0][m]].plot((trig_time[m], trig_time[m]),
(axs[ind_trig[0][m]].get_ylim()[0], axs[ind_trig[0][m]].get_ylim()[1]),
'r-', linewidth=1.0)
for m in range(len(loc_ind_trig)):
if len(self.origins):
axs[loc_ind_trig[m]].plot((self.tcalc[len(self.origins)-1][m], self.tcalc[len(self.origins)-1][m]),
(axs[loc_ind_trig[m]].get_ylim()[0], axs[loc_ind_trig[m]].get_ylim()[1]),
'g-', linewidth=1.0)
if len(self.magnitudes):
axs[loc_ind_trig[m]].plot((self.ts_approx[m], self.ts_approx[m]),
(axs[loc_ind_trig[m]].get_ylim()[0], axs[loc_ind_trig[m]].get_ylim()[1]),
'y-', linewidth=1.0)
if len(self.origins):
loc_ind_trig = self.loc_ind[-1][
:] # why was this? if this line gives error switch between ind_trig[0] and self.loc_ind[-1][:]
for m in range(len(loc_ind_trig)):
axs[loc_ind_trig[m]].plot(
(self.tcalc[len(self.origins) - 1][m], self.tcalc[len(self.origins) - 1][m]),
(axs[loc_ind_trig[m]].get_ylim()[0], axs[loc_ind_trig[m]].get_ylim()[1]),
'g-', linewidth=1.0)
if len(self.magnitudes):
axs[loc_ind_trig[m]].plot((self.ts_approx[m], self.ts_approx[m]),
(axs[loc_ind_trig[m]].get_ylim()[0],
axs[loc_ind_trig[m]].get_ylim()[1]),
'y-', linewidth=1.0)
# if self.classification != 'noise':
# trig_ch = [int(i.waveform_id['station_code']) for i in self.picks]
# ch_in = [int(i.stats['station']) for i in self.wf_stream.traces]
# # finding the correct axes(stations)
# ind_trig = np.where(np.isin(ch_in, trig_ch)) # Finds correct axis
# trig_time = [(i.time - self.wf_stream.traces[0].stats["starttime"]) * 1000 for i in self.picks]
# loc_ind_trig = ind_trig[0] #why was this? if this line gives error switch between ind_trig[0] and self.loc_ind[-1][:]
#
# # plotting pick time
# for m in range(len(trig_time)):
# axs[ind_trig[0][m]].plot((trig_time[m], trig_time[m]),
# (axs[ind_trig[0][m]].get_ylim()[0], axs[ind_trig[0][m]].get_ylim()[1]),
# 'r-', linewidth=1.0)
# for m in range(len(loc_ind_trig)):
# if len(self.origins):
# axs[loc_ind_trig[m]].plot((self.tcalc[len(self.origins)-1][m], self.tcalc[len(self.origins)-1][m]),
# (axs[loc_ind_trig[m]].get_ylim()[0], axs[loc_ind_trig[m]].get_ylim()[1]),
# 'g-', linewidth=1.0)
# if len(self.magnitudes):
# axs[loc_ind_trig[m]].plot((self.ts_approx[m], self.ts_approx[m]),
# (axs[loc_ind_trig[m]].get_ylim()[0], axs[loc_ind_trig[m]].get_ylim()[1]),
# 'y-', linewidth=1.0)
plt.suptitle(t_title, fontsize=17)
fig.text(0.49, 0.035, 'time [ms]', ha='center', fontsize=14)
......
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment