How to set STFT-parameters to visualize woodpecker / Sawtooth Signal? ​

Question

I always wanted to work with FFTs and Spectrograms to characterise sounds and their frequencies using python.

Lately I recorded a woodpecker at his work.

Woodpecker Wav File

To know where to look at I compared the sound with a generated sawtooth from Onlinegenerator.com. I think 17 Hz fits well.

17 Hz Sawtooth Wav File

Checking the code with the generated file worked good.

But sadly I cannot find the frequency in the DFT/FFT or spectrogram generated with the code below. No matter which settings I tried, I cannot visualise the frequencies below 20 Hz properly. No major signals cuts on or off. Only a signal at 23 Hz, which seems to be continuous.

​from scipy import signal
import scipy.io.wavfile as wav
import matplotlib.pyplot as plt
import numpy as np

​# Read WAV File
the_file = 'IMG_0710_short.wav'
samplerate, samples = wav.read(the_file)

#samples_left = samples[:,0]
samples = samples[:,1] # reduce to right channel

timevec = range(len(samples)) # Time vector for plot
timevec = [x / samplerate for x in timevec]
t_max =timevec[-1]

​# Lineplot of Signal
dpi = 200
plt.rcParams['figure.dpi']= dpi

plt.plot(timevec,samples)
plt.title(the_file)
plt.ylabel('Amplitude [-]')
plt.xlabel('Time [sec]')

​# STFT Settings
nperseg  = 0.5 * samplerate # Window Size 1
noverlap = nperseg*0.95 # Overlap
nfft     = 1 * samplerate # STFT 2
window   = 'hann' # Window Type

timespan = [0, 2] # Calculation Window
fromm = int(len(samples)/t_max*timespan[0])
too   = int(len(samples)/t_max*timespan[1])

f, t, Zxx = signal.stft(samples[fromm:too], samplerate, nperseg=nperseg, window=window, noverlap=noverlap, nfft=nfft)
t = t + timespan[0]
cmap=plt.cm.nipy_spectral
vmin = 10
vmax = 18
fig = plt.figure(figsize=(7, 5))
pcm = plt.pcolormesh(t, f, np.log(np.abs(Zxx)), cmap=cmap, vmin=vmin, vmax=vmax)

plt.ylim(0,50)
fig.colorbar(pcm)
plt.title(the_file)
plt.ylabel('Frequency [Hz]')
plt.xlabel('Time [sec]')

# FFT
​n = len(samples) # length of the signal
k = np.arange(n)
T = n/samplerate
frq = k/T # two sides frequency range
frq = frq[range(int((n/2)))] # one side frequency range

Y = np.fft.fft(samples, norm='ortho')#/n # fft computing and normalization
Y = Y[range(int(n/2))]

fig, ax = plt.subplots(2, 1)
ax[0].plot(timevec,samples) # plotting the signal
ax[0].set_xlabel('Time')
ax[0].set_ylabel('Amplitude')
ax[1].plot(frq,abs(Y),'r') # plotting the spectrum
ax[1].set_xlabel('Frequency (Hz)')
ax[1].set_ylabel('|Y(freq)|')
ax[1].set_xlim(0,50)

How do I need to set the stft parameters in order to visualise the woodpecker spectrum correctly? Or is there a fundamental problem with stft here?

Answer 1

Your first problem is that you have your axes mixed up.

Your "17 Hz" picture shows this quite clearly. The red bars are the impulses from the signal. The bars should indicate an impulse (wide range of frequencies) but your axes indicate they are across time.

This is what it should look like:

Notice that the bars are along the frequency axis rather than time as in yours.

This is your woodpecker recording:

It contains bars, too. But, they are not as wide band and are much closer (in intensity) to the background noise.

The bars repeat at just above 20Hz - you can see that the frequency is around 23Hz in your last plot (the red line plot below the blue oscilloscope style plot.)

Your overlap of 0.95 amounts to a time resolution of about just about 40Hz. You need a higher overlap of maybe 0.99 to separate the bars.

You could downsample the recording down to like 11025 Hz sampling rate. That would get you better frequency resolution and make calculating the spectrogram faster. You can do that because there's not really much of interest in your signal above about 4500Hz.