Skip to content

Commit

Permalink
Changed the SNR detector of impulse detections, added the detection b…
Browse files Browse the repository at this point in the history 
…and and alaysis band selection
  • Loading branch information
Clea Parcerisas committed Jun 25, 2021
1 parent 32650ba commit ee759f3
Show file tree
Hide file tree
Showing 3 changed files with 65 additions and 33 deletions.
15 changes: 7 additions & 8 deletions pypam/acoustic_file.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -812,7 +812,7 @@ def correlation(self, signal, fs_signal, binsize=1.0):
return fs

def detect_piling_events(self, min_separation, max_duration, threshold, dt, binsize=None,
verbose=False, **kwargs):
verbose=False, save_path=None, band=None, method=None):
"""
Detect piling events
Expand All @@ -839,20 +839,19 @@ def detect_piling_events(self, min_separation, max_duration, threshold, dt, bins

detector = impulse_detector.PilingDetector(min_separation=min_separation,
max_duration=max_duration,
threshold=threshold, dt=dt)
threshold=threshold, dt=dt, detection_band=band,
analysis_band=self.band)
total_events = pd.DataFrame()
if 'save_path' in kwargs.keys():
save_path = kwargs['save_path']
else:
save_path = None
for i, block in enumerate(sf.blocks(self.file_path, blocksize=blocksize, start=self._start_frame)):
time_bin = self.time_bin(blocksize, i)
print('bin %s' % time_bin)
signal_upa = self.wav2upa(wav=block)
signal = Signal(signal=signal_upa, fs=self.fs, channel=self.channel)
signal.set_band(band=self.band)
events_df = detector.detect_events(signal, method='snr', verbose=verbose,
save_path=save_path)
events_df = detector.detect_events(signal, method=method, verbose=verbose,
save_path=save_path.joinpath('%s.png' %
datetime.datetime.strftime(time_bin,
"%y%m%d_%H%M%S")))
events_df['datetime'] = pd.to_timedelta(events_df[('temporal', 'start_seconds')],
unit='seconds') + time_bin
events_df = events_df.set_index('datetime')
Expand Down
4 changes: 2 additions & 2 deletions pypam/acoustic_survey.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -368,7 +368,7 @@ def cut_and_place_files_period(self, period, folder_name, extensions=None):
pass
return 0

def detect_piling_events(self, min_separation, max_duration, threshold, dt=None, verbose=False):
def detect_piling_events(self, min_separation, max_duration, threshold, dt=None, verbose=False, **kwargs):
"""
Return a DataFrame with all the piling events and their rms, sel and peak values
Expand Down Expand Up @@ -398,7 +398,7 @@ def detect_piling_events(self, min_separation, max_duration, threshold, dt=None,
threshold=threshold,
max_duration=max_duration,
dt=dt, binsize=self.binsize,
verbose=verbose)
verbose=verbose, **kwargs)
df = df.append(df_output)
return df

Expand Down
79 changes: 56 additions & 23 deletions pypam/impulse_detector.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -26,7 +26,7 @@


class ImpulseDetector:
def __init__(self, min_separation, max_duration, band, threshold=150, dt=None):
def __init__(self, min_separation, max_duration, detection_band, analysis_band, threshold=150, dt=None):
"""
Impulse events detector
Expand All @@ -36,6 +36,10 @@ def __init__(self, min_separation, max_duration, band, threshold=150, dt=None):
Minimum separation of the event, in seconds
max_duration : float
Maxium duration of the event, in seconds
detection_band : list
List or tuple of (min_freq, max_freq) of the frequency band used to detect the signals
analysis_band : list
List or tuple of (min_freq, max_freq) of the frequency band used to characterize the impulses
threshold : float
Threshold above ref value which one it is considered piling, in db
dt : float
Expand All @@ -46,7 +50,8 @@ def __init__(self, min_separation, max_duration, band, threshold=150, dt=None):
self.max_duration = max_duration
self.threshold = threshold
self.dt = dt
self.band = band
self.detection_band = detection_band
self.analysis_band = analysis_band

def detect_events(self, signal, method='dt', verbose=False, save_path=None):
"""
Expand Down Expand Up @@ -89,7 +94,7 @@ def detect_events_dt(self, signal, verbose=False, save_path=None):
save_path : string or Path
Where to save the image. Set to None if it should not be saved
"""
signal.set_band(self.band)
signal.set_band(self.detection_band)
levels = []
blocksize = int(self.dt * signal.fs)
for block in signal.Blocks(blocksize=blocksize):
Expand All @@ -116,10 +121,8 @@ def detect_events_envelope(self, signal, verbose=False, save_path=None):
save_path : string or Path
Where to save the image. Set to None if it should not be saved
"""
signal.set_band(band=self.band)
signal.set_band(band=self.detection_band)
envelope = signal.envelope()
envelope = utils.to_db(envelope, ref=1.0, square=True)

times_events = events_times_diff(signal=envelope, fs=signal.fs, threshold=self.threshold,
max_duration=self.max_duration,
min_separation=self.min_separation)
Expand All @@ -144,12 +147,12 @@ def detect_events_snr(self, signal, verbose=False, save_path=None):
Where to save the image. Set to None if it should not be saved
"""
blocksize = int(self.dt * signal.fs)
signal.set_band(band=self.band)
signal.set_band(band=self.detection_band)
envelope = signal.envelope()
envelope = utils.to_db(envelope, ref=1.0, square=True)
times_events = events_times_snr(signal=envelope, blocksize=blocksize, fs=signal.fs,
threshold=self.threshold, max_duration=self.max_duration,
min_separation=self.min_separation)
min_separation=self.min_separation, original_sig=signal.signal)

events_df = self.load_all_times_events(times_events, signal, verbose=verbose)

if verbose:
Expand Down Expand Up @@ -183,10 +186,11 @@ def load_event(self, s, t, duration, removenoise=True, verbose=False):
n2 = min(int((t + duration + self.min_separation) * s.fs), s.signal.shape[0])
event = Event(s.signal[n1:n2], s.fs)
noise_clip = np.concatenate((s.signal[n1:start_n], s.signal[end_n:n2]))
event.reduce_noise(noise_clip=noise_clip, nfft=512, verbose=verbose)
event.reduce_noise(noise_clip=noise_clip, nfft=512, verbose=False)
event.signal = event.signal[start_n - n1:end_n - n1]
else:
event = Event(s.signal[start_n:end_n], s.fs)

return event

def load_all_times_events(self, times_events, signal, verbose=False):
Expand All @@ -202,10 +206,10 @@ def load_all_times_events(self, times_events, signal, verbose=False):
verbose : bool
Set to True to plot all the events of the signal
"""
signal.set_band([10, 20000])
signal.set_band(self.analysis_band)
columns_temp = ['start_seconds', 'end_seconds', 'duration', 'rms', 'sel', 'peak']
columns_df = pd.DataFrame({'variable': 'temporal', 'value': columns_temp})
freq = sci.fft.rfftfreq(128) * 40000
freq = sci.fft.rfftfreq(128) * signal.fs
columns_df = pd.concat([columns_df, pd.DataFrame({'variable': 'psd', 'value': freq})])
columns = pd.MultiIndex.from_frame(columns_df)
events_df = pd.DataFrame(columns=columns)
Expand All @@ -231,7 +235,7 @@ def plot_all_events(self, signal, events_df, save_path=None):
save_path : string or Path
Where to save the image. Set to None if it should not be saved
"""
signal.set_band(band=self.band)
signal.set_band(band=self.detection_band)
fbands, t, sxx = signal.spectrogram(nfft=512, scaling='spectrum', db=True, mode='fast')
fig, ax = plt.subplots(3, 1, sharex='col')
im = ax[0].pcolormesh(t, fbands, sxx, shading='auto')
Expand Down Expand Up @@ -275,9 +279,13 @@ def plot_all_events(self, signal, events_df, save_path=None):


class PilingDetector(ImpulseDetector):
def __init__(self, min_separation, max_duration, threshold, dt):
def __init__(self, min_separation, max_duration, threshold, dt, detection_band=None, analysis_band=None):
if detection_band is None:
detection_band = [500, 1000]
if analysis_band is None:
analysis_band = [100, 20000]
super().__init__(min_separation=min_separation, max_duration=max_duration,
band=[5000, 10000], threshold=threshold, dt=dt)
detection_band=detection_band, analysis_band=analysis_band, threshold=threshold, dt=dt)


@nb.jit
Expand All @@ -304,21 +312,30 @@ def events_times_diff(signal, fs, threshold, max_duration, min_separation):
event_max_val = 0
last_xi = 0
i = 0
threshold_upa = 10 ** (threshold/20.0)
while i < len(signal):
xi = signal[i]

if event_on:
duration = (i - event_start) / fs
if duration >= max_duration or xi < (event_max_val - threshold):
# Event finished, too long! Or event detected!
if duration >= max_duration or (event_max_val - xi) >= threshold_upa:
# Event finished, too long! Or event end detected
event_on = False
event_end = i
times_events.append([event_start / fs, duration, event_end / fs])
i += min_separation_samples
event_max_val = 0
# plt.Figure()
# plt.plot(signal[event_start - min_separation_samples:event_end + min_separation_samples])
# plt.axvline(min_separation_samples, color='green', label='Start')
# plt.axvline(min_separation_samples + event_end - event_start, color='red', label='End')
# plt.show()
# plt.close()
elif xi > event_max_val:
event_max_val = xi
else:
if i != 0 and (xi - last_xi) >= threshold:
if i != 0 and (xi - last_xi) >= threshold_upa:
# There is a big jump, start event!
event_on = True
event_start = i
event_max_val = xi
Expand All @@ -328,31 +345,47 @@ def events_times_diff(signal, fs, threshold, max_duration, min_separation):


@nb.jit
def events_times_snr(signal, fs, blocksize, threshold, max_duration, min_separation):
def events_times_snr(signal, fs, blocksize, threshold, max_duration, min_separation, original_sig):
times_events = []
min_separation_samples = int(min_separation * fs)
event_on = False
event_start = 0
event_end = 0
j = 0
threshold_upa = 10 ** (threshold/20.0)
while j < len(signal):
if j + blocksize > len(signal):
blocksize = len(signal) - j
noise = np.mean(signal[j:j + blocksize])
noise = np.sqrt(np.mean(signal[j:j + blocksize]**2))
max_value = noise
for i in np.arange(blocksize - 1) + j:
xi = signal[i]
snr = xi - noise
if event_on:
duration = (i - event_start) / fs
if duration >= max_duration or snr < threshold:
if duration >= max_duration or (xi - noise) < 10**(6.0/20.0):
# Event finished, too long! Or event detected!
event_on = False
event_end = i
times_events.append([event_start / fs, duration, event_end / fs])
# plt.Figure()
# plt.plot(original_sig[event_start-min_separation_samples:event_end+min_separation_samples],
# label='Signal')
# plt.plot(signal[event_start-min_separation_samples:event_end+min_separation_samples],
# label='Envelope')
# plt.axhline(noise + threshold_upa, label='Threshold')
# plt.axhline(noise + threshold_upa/2.0, label='Threshold')
# plt.axvline(min_separation_samples, color='green', label='Start')
# plt.axvline(min_separation_samples+event_end-event_start, color='red', label='End')
# plt.show()
# plt.close()
if xi > max_value:
max_value = xi
else:
if snr >= threshold:
if (xi - noise) >= threshold_upa:
if len(times_events) == 0 or (i - event_end) >= min_separation_samples:
event_on = True
event_start = i
max_value = xi
j += blocksize

return times_events

0 comments on commit ee759f3

Please sign in to comment.