Spectrogram output versus figure
17 views (last 30 days)
Show older comments
When using spectrogram to calculate power spectral density for a large dataset, the plot produced by directly running the command with no output arguments is slightly different from what I get when I run it with an output argument and plot that output myself.
For example, the following plots spectrograms for a subset of my data and compares the direct result of spectrogram with the output raster:
load('exampledata.mat')
window = 10000;
noverlap = round(0.1*window);
fs = 10000;
figure(1);clf
t = tiledlayout(3,1, 'TileIndexing', 'columnmajor');
% Plot spectrogram directly from spectrogram function
s(1) = nexttile;
spectrogram(exampledata,window,noverlap,window,fs);
title('Original spectrogram plot (real data)')
PSDdB_orig = get(get(gca,'Children'),'CData');
caxis([min(PSDdB_orig,[],'all') max(PSDdB_orig,[],'all')])
% Calculate power spectral density and plot as decibels
[~,F,T,PSD] = spectrogram(exampledata,window,noverlap,window,fs);
PSDdB = (10.*log10(PSD))';
s(2) = nexttile;
imagesc(F./1000,T./60,PSDdB)
title('Spectrogram output raster (real data)')
c = colorbar;
caxis([min(PSDdB,[],'all') max(PSDdB,[],'all')])
set(get(c,'label'),'string','Power/frequency (dB/Hz)','Rotation',90);
% Calculate ratio of the two approaches above and plot that
PSDdiff = PSDdB_orig./PSDdB;
s(3) = nexttile;
imagesc(F./1000,T./60,PSDdiff)
title('Original plot data divided by output raster data')
c = colorbar;
caxis([min(PSDdiff,[],'all') max(PSDdiff,[],'all')])
set(get(c,'label'),'string','original/output','Rotation',90);
set(s(2:3),'YDir','normal')
xlabel(s(2:3), s(1).XLabel.String)
ylabel(s(2:3), s(1).YLabel.String)
From the figures above, it looks like they produce nearly the same result, except at high frequencies (the middle figure looks more different than the first mostly due to the different color scales).
Comparing the ranges and plotting a cross-section (example spectra at a single time) from each one shows that the difference seems to be that spectrogram's direct plotting method is cutting off local minima:
disp(['min/max of original spectrogram plot: ',num2str(min(PSDdB_orig,[],'all')),', ',num2str(max(PSDdB_orig,[],'all'))])
disp(['min/max of spectrogram output raster: ',num2str(min(PSDdB,[],'all')),', ',num2str(max(PSDdB,[],'all'))])
% Plot example spectra from each approach
figure(2);clf
subplot(2,1,1);
plot(F./1000,PSDdB(10,:),F./1000,PSDdB_orig(10,:),'--');
legend('Spectrogram output raster','Original spectrogram')
title('Example spectra')
xlabel('Frequency (kHz)')
ylabel('Power/frequency (dB/Hz)')
% zoom in on high frequency range where values differ
subplot(2,1,2,copyobj(gca,gcf))
xlim([4.75 4.85])
title('high frequency close-up')
But weirdly, I can't reproduce this with randomly generated synthetic data:
x = randn(length(exampledata),1);
figure(3);clf
t = tiledlayout(3,1, 'TileIndexing', 'columnmajor');
s(1) = nexttile;
spectrogram(x,window,noverlap,window,fs);
title('Original spectrogram plot (synthetic data)')
PSDdB_orig = get(get(gca,'Children'),'CData');
caxis([min(PSDdB_orig,[],'all') max(PSDdB_orig,[],'all')])
[~,F,T,PSD] = spectrogram(x,window,noverlap,window,fs);
PSDdB = (10.*log10(PSD))';
PSDdiff = PSDdB_orig./PSDdB;
s(2) = nexttile;
imagesc(F./1000,T./60,PSDdB)
title('Spectrogram output raster (synthetic data)')
c = colorbar;
caxis([min(PSDdB,[],'all') max(PSDdB,[],'all')])
set(get(c,'label'),'string','Power/frequency (dB/Hz)','Rotation',90);
s(3) = nexttile;
imagesc(F./1000,T./60,PSDdiff)
title('Original plot data divided by output raster data')
c = colorbar;
caxis([min(PSDdiff,[],'all') max(PSDdiff,[],'all')])
set(get(c,'label'),'string','original/output','Rotation',90);
set(s(2:3),'YDir','normal')
xlabel(s(2:3), s(1).XLabel.String)
ylabel(s(2:3), s(1).YLabel.String)
disp(['min/max of original spectrogram plot: ',num2str(min(PSDdB_orig,[],'all')),', ',num2str(max(PSDdB_orig,[],'all'))])
disp(['min/max of spectrogram output raster: ',num2str(min(PSDdB,[],'all')),', ',num2str(max(PSDdB,[],'all'))])
I do notice that the minimum PSD value for my data (-196 dB) is significantly lower than the minimum of the synthetic data. Is there a lower limit where spectrogram starts to censor spectra (say, around -156.5351 dB)? Or is something else going on here? Any help would be much appreciated!! Thanks!
0 Comments
Answers (1)
Paul
on 5 Sep 2024
Hi DLR,
Stepping through spectrogram shows that, when output argruments aren't specified, the PSD data that's plotted is offset by eps("double") when creating the plot. See modified code below.
load('exampledata.mat')
window = 10000;
noverlap = round(0.1*window);
fs = 10000;
figure(1);clf
t = tiledlayout(3,1, 'TileIndexing', 'columnmajor');
% Plot spectrogram directly from spectrogram function
s(1) = nexttile;
spectrogram(exampledata,window,noverlap,window,fs);
title('Original spectrogram plot (real data)')
PSDdB_orig = get(get(gca,'Children'),'CData');
caxis([min(PSDdB_orig,[],'all') max(PSDdB_orig,[],'all')])
% Calculate power spectral density and plot as decibels
[~,F,T,PSD] = spectrogram(exampledata,window,noverlap,window,fs);
%PSDdB = (10.*log10(PSD))';
PSDdB = 10*log10(PSD + eps("double"))'; % what's actually plotted when output arguments not specified.
s(2) = nexttile;
imagesc(F./1000,T./60,PSDdB)
title('Spectrogram output raster (real data)')
c = colorbar;
caxis([min(PSDdB,[],'all') max(PSDdB,[],'all')])
set(get(c,'label'),'string','Power/frequency (dB/Hz)','Rotation',90);
% Calculate ratio of the two approaches above and plot that
PSDdiff = PSDdB_orig./PSDdB;
s(3) = nexttile;
imagesc(F./1000,T./60,PSDdiff)
title('Original plot data divided by output raster data')
c = colorbar;
% caxis([min(PSDdiff,[],'all') max(PSDdiff,[],'all')])
set(get(c,'label'),'string','original/output','Rotation',90);
set(s(2:3),'YDir','normal')
xlabel(s(2:3), s(1).XLabel.String)
ylabel(s(2:3), s(1).YLabel.String)
2 Comments
Paul
on 13 Sep 2024 at 16:42
eps is not the double precision limit, at least that's not how I would state it.
As best I can tell, the purpose of adding eps to PSD for purposes of plotting is just to limit the range of the colors on the plot. On the other hand, there may have been better ways to do that besides adding a small value to every single element of PSD, so maybe there really is something more to it. But I doubt it.
"Do you know whether this implies that working with values below 10^-16 should generally be avoided?"
No I don't. Double precision floating point can represent much small numbers:
realmin("double")
Whether or not it's sensible to use such small numbers is beyond my knowledge.
Maybe you can change the units on your data or use the outputs from spectrogram and roll your own plots if needed.
See Also
Categories
Find more on AI for Signals in Help Center and File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!