Description of the data and problem:
I have a signal sampled at 1000 Hz. I'm low-pass filtering it at 120 Hz, and want to make spectrograms of the frequencies below this threshold. I'm using the scipy functions
fftpack.fftfreq to get the Fourier coefficients, then
fftpack.fft to do the actual transform.
The signal is quite long, and I want the spectrograms to be about 5 seconds in length using 50 millisecond windows. I also filter out negative frequencies for plotting. Given the number of samples in this window (50), I get 50 Fourier coefficients. However, half of these are negative, and another portion goes between 120 Hz and 500 Hz (naturally, because the sampling rate is 1000 Hz).
This leaves me with pretty low-resolution spectrograms. I only have about 6 blocks in the frequency area of interest (0, 20, 40, 60, 100, and 120 Hz), then 18 blocks showing low activity in the frequencies up to 500 Hz, since they were filtered out.
How could I, for instance, do the FFT specifying something like 50 frequency bins between only 0 and 120 Hz?
What I tried:
I tried doing something like this by using
np.linspace to specify the frequencies (instead of
fft.fftfreq), but this introduced some other bug, namely the spectrograms always looked like a mirror, with high power at the highest and lowest frequencies and low power in the middle of the graph, regardless of the range. I'm honestly not sure why this happens with linspace and not with fftfreq. They both return arrays of floats.
Code sample included below. Thanks for any help in advance! Cheers
X = samples_filt # long 1-D vector of low-pass filtered samples fs = 1000 time = .05 # window length in sec N = int(fs * time) # num samples tot_len = 5 * fs # 5 sec of the whole signal X = X[:tot_len] f = fftpack.fftfreq(N, 1.0/fs) # f = np.ceil(np.linspace(0, max_freq, 52)[1:51]) # introduced error described above mask = (f > 0) # mask for positive freqs n_max = int(np.ceil(X.shape / N)) # the number of segments of length N in the sample array data f_values = np.sum(1 * mask) # how many values meet mask reqs spectogram_data = np.zeros((f_values, n_max)) window = sp.signal.blackman(N) # taper used to improve contrast of spectrogram for n in range(0, n_max): subdata = X[(N * n):(N * (n + 1))] F = fftpack.fft(subdata * window) spectogram_data[:, n] = np.log(abs(F[mask])) fig, ax = plt.subplots(1, 1, figsize=(8, 6)) p = ax.imshow(spectogram_data, origin='lower', extent=(0, X.shape / fs, 0, max(f)), aspect='auto', cmap=mpl.cm.RdBu_r) cb = fig.colorbar(p, ax=ax)
For instance, it would be great to have this same graph with 24 frequency bins going up to 120 Hz.