plot_spectrogram() doesn't show dB re. 2e-5 Pa #211

Activity

1. topherbuckley commented on Feb 28, 2018

   

   topherbuckley

   on Feb 28, 2018

   More actions

   @duracell013 can you provide steps (or more preferably code) as to how you are creating your plot? It looks like the plot_spectrogram() method within _signal.pyx doesn't modify the data array at all, in other words, the input data to this method should already represent the SPL in dB, if that is what you are seeking to plot.

   @FRidh Can you elaborate as to where the "data" variable from line 679 of _signal.pyx is generated? np.squeeze(self) is not very descriptive, and without understanding the whole Signal class or how the input attribute "numpy.ndarry" on line 16 is provided, I don't think I'd be able to help troubleshoot this very far. In short, how can we modify the "data" ndarray or where is it being modified prior to being sent to the Signal class?
2. duracell013 commented on Feb 28, 2018

   

   duracell013

   on Feb 28, 2018

   Author

   More actions

   @topherbuckley Here is a sample code which plots the spectrogram with the acoustics package and with scipy.signal.spectrogram().

   import numpy as np
   import matplotlib.pyplot as plt
   import scipy.signal as ss
   import acoustics
   
   N_POINTS = 200000
   MAIN_FREQUENCY = 500
   SAMPLING_FREQUENCY = 5000
   
   FFT_POINTS = 4096
   FFT_OVERLAP = 128
   PLOT_CLIM = (40, 100)
   
   P_REF = 2e-5
   
   plt.ion()
   
   T = np.arange(N_POINTS)/SAMPLING_FREQUENCY
   Y_main = np.sin(2*np.pi*MAIN_FREQUENCY*T)
   Y = Y_main + np.random.rand(N_POINTS)
   Y = ss.detrend(Y)
   
   rms = np.sqrt(np.mean(Y_main**2))
   rms_dB = 20*np.log10(rms/P_REF)
   print('Amplitude of main frequency = {:.2f} dB'.format(rms_dB))
   
   # Plot spectrogram using acoustics
   signal = acoustics.Signal(Y, SAMPLING_FREQUENCY)
   signal.plot_spectrogram(title='acoustics.Signal.plot_spectrogram()')
   
   # Spectrogram with scipy
   f, t, Sxx = ss.spectrogram(Y, fs=SAMPLING_FREQUENCY,
                              window='hanning',
                              nperseg=FFT_POINTS,
                              noverlap=FFT_OVERLAP,
                              mode='psd',
                              scaling='spectrum')
   Sxx_dB = 10*np.log10(Sxx/P_REF**2)
   print('Max amplitude using Scipy = {:.2f} dB'.format(Sxx_dB.max()))
   
   fig = plt.figure()
   pc = plt.pcolormesh(t, f, Sxx_dB)
   pc.set_clim(PLOT_CLIM)
   cb = fig.colorbar(pc)
   plt.ylabel('Frequency [Hz]')
   plt.xlabel('Time [sec]')
   plt.title('scipy.signal.spectrogram()')
   cb.set_label('SPL [dB re. 20 $\mathrm{\mu{}Pa}$]')
3. FRidh commented on Feb 28, 2018

   

   FRidh

   on Feb 28, 2018

   Member

   More actions

   The Signal class is a subclass of the numpy array and thus self is also an array. Squeezing is done because plot_spectrogram can handle only one-channel signals.
4. FRidh commented on Feb 28, 2018

   

   FRidh

   on Feb 28, 2018

   Member

   More actions

   I believe the plot_spectrogram() function doesn't output dB relative to the reference acoustic pressure (2e-5 Pa). This can be misleading as the colorbar label is "SPL in dB".

   Correct, the values are re. 1, so there is a shift in decibel.

   In my opinion the way forward is that the spectrogram method gets a reference kwarg (see how it is done with the other methods) which defaults to 1 (so we don't break anything, and remain neutral as to what type of signal we're talking about). The plot_spectrogram method ideally reuses the spectrogram method but it should be investigated how this impacts performance.
5. duracell013 commented on Feb 28, 2018

   

   duracell013

   on Feb 28, 2018

   Last edited by duracell013

   Author

   More actions

   @FRidh I dont see the conversion to dB in the plot_spectrogram function. I understand that the input data of specgram is in Pa, so I would expect the output to be in Pa. Am I missing something here?
6. FRidh commented on Feb 28, 2018

   

   FRidh

   on Feb 28, 2018

   Member

   More actions

   The specgram function in matplotlib does the conversion to decibel.
7. duracell013 commented on Feb 28, 2018

   

   duracell013

   on Feb 28, 2018

   Author

   More actions

   I see, thanks!
8. topherbuckley commented on Mar 1, 2018

   

   topherbuckley

   on Mar 1, 2018

   More actions

   @duracell013 , maybe @FRidh has already answered your question, in which case no need to further debug, but your code doesn't seem to show any plot when run on my machine. It flashes a plot window and immediately closes. I tried tinkering a bit (adding plt.show()), but nothing showing yet. One thing I did notice was that your scipy spectrogram is in 'psd' mode, while this library is using 'magnitude'. In addition to the reference pressure being different, this may also add to any differences you are seeing.

   @FRidh , yes I'm sorry, I definitely mistook the ndarray inheritance as an input parameter. Thanks for pointing that out kindly :) . I'm not a native python programmer (or O.O. programmer for that matter) so these things trick me once and a while.

   I think I understand your method for moving this forward. I'll send a PR to see if it agrees with what you were thinking.