Scipy and Matlab Spectrogram not matching

Question

I am trying to replicate a spectrogram from MATLAB in Python. I've read other posts but they either don't use complex data or the data doesn't match between languages.

I've defined my FFT length (nff), sample rate (fs), overlap (nov) and window (hamming) to all be the same between MATLAB and Python.

Here is the spectrogram from scipy and you can see the issue is the scaling is incorrect. I see that Matlab handles windowing differently than scipy. In scipy your window length and segment value (nperseg) must be the same. This is not the case in Matlab. Also, in Matlab your FFT length does not have to be larger than your window length however that is a requirement in scipy.

I've tried to keep the parameters the same between Matlab and Python but I still am getting the wrong matrix values for plotting. I've been banging my head against the wall on this for longer than I'd like to admit, any help would be great.

MATLAB:

fs = 50000;
t = (0:1/fs:10);
x = sin(2*pi*10000*t)+sin(2*pi*20000*t);
nsc = floor(50e-3*fs);           % Section Length 
nov = floor(nsc*0.90);           % Percent Overlap 
nff = 2^12;                      % FFT Length 

%Plot Spectrogram 
[Sxx, spec_f,t] = spectrogram(x,hamming(nsc),nov,nff,fs);

Sxx_dB = 10*log10(abs(Sxx));

figure;
imagesc(t,spec_f,Sxx_dB)
ylabel('Frequency (kHz)','Fontsize',12,'Interpreter','latex')
xlabel('Time (sec)','Fontsize',12,'Interpreter','latex')
colorbar

colormap('jet')
set(gca,'Fontsize',12);
set(gca, 'YDir','normal')

Python:

import numpy as np 
from scipy import signal
import math
from matplotlib import pyplot as plt

fs = 50000
t = np.arange(0, 10, 1/fs)
x = np.sin(2*np.pi*10000*t) + np.sin(2*np.pi*20000*t)
nsc = math.floor(50e-3*fs);           
nov = math.floor(nsc*0.90);         
nff = 2**12;     

f,t,Sxx = signal.spectrogram(x, fs=fs, window=np.hamming(nsc), 
                             noverlap=nov, nperseg=nsc, nfft=nff, mode='complex')

fig1 = plt.figure()
im1= plt.pcolormesh(t,f,10*np.log10(np.abs(Sxx)), cmap='jet')
fig1.colorbar(im1).set_label('Intensity (dB)')
plt.xlabel('Time [sec]')
plt.ylabel('Frequency [Hz]')
plt.show()

Answer 1

It’s a scaling problem, for which there are two options:

Option 1: Use PSD scaling for both Python and Matlab

1. Use mode='psd' in the scipy.spectrogram call.

2. Scale the matlab version:

   win = hamming(nsc); 
   scaling_factor = 2/(fs*sum(win.^2)); % PSD scaling factor
   
   % Spectrogram not scaled
   [Sxx, spec_f, t] = spectrogram(x,win,nov,nff,fs);
   
   % Correct PSD scaling
   Sxx_psd = scaling_factor .* abs(Sxx).^2;
   Sxx_psd_dB = 10*log10(Sxx_psd);
   
   % Plot
   figure;
   imagesc(t,spec_f,Sxx_psd_dB)
   

   Alternatively, Matlab can take care of the scaling:

   % Spectrogram with PSD scaling
   [~, spec_f, t, Pxx] = spectrogram(x,win,nov,nff,fs);  
   Pxx_dB = 10*log10(Pxx);
   
   % Plot
   figure;
   imagesc(t,spec_f, Pxx_dB);
   

Option 2: Scale Matlab Only

If you want to keep mode = 'complex' in scipy.spectrogram, then modify the Matlab scaling_factor from option 1 to:

scaling_factor = 1/sqrt(fs*sum(win.^2));

That's the scaling factor used by scipy.

Additionally, do not square abs(Sxx):

Sxx_complex = scaling_factor .* abs(Sxx);
Sxx_complex_dB = 10*log10(Sxx_complex);

Answer 2

I was able to resolve the issue. It seems that the scaling is incorrect when generating the complex matrix so instead compute the matrix manually and stack the FFT results. A great resource is this SDR website.

Here is a small Python code snippet showing how to calculate your spectrogram manually where nff is the number of FFT's, Y is the signal you want to plot:

#Calculate Spectrogram Matrix
num_rows = int(np.floor(len(Y)/nff))
spectrogram = np.zeros((num_rows, nff))
for i in range(num_rows):
    spectrogram[i,:] = 10*np.log10(np.abs(np.fft.fftshift(np.fft.fft(Y[i*nff:(i+1)*nff]))))         #create spectrogram, fill in rows with FFT, shift to center, calculate dB magnitude
spectrogram = spectrogram[:,nff//2:] # get rid of negative freqs because we simulated a real signal