plot fails with -Inf

9 visualizzazioni (ultimi 30 giorni)
Jonathan
Jonathan il 23 Mar 2020
Commentato: Jonathan il 23 Mar 2020
I need to plot a fft of an audio file. It appears that if any of the fft bins are zero, the conversion to dB produces an -Inf. Then the plot function does draw the figure but does not plot any data.
Is there an efficent way to trap the 0 value bins in the conversion to dB? For my work -Inf PSD would be equiv to zero.
The test audio file was a code generated audio file, 100Hz Sine for 10 seconds. Since only one fft bin will have a value, all the other bins are zero.
Thanks.
%Code Fragment:
[yS, Fs] = audioread(PathFN,'double');
L = length(yS);
Y = fft(yS,L);
% Version One
Px = Y.*conj(Y); % PSD Power of each freq components
Px = Px(1:L/2);
Px = 10*log10(Px); % Convert to dB
fVals = Fs * (0:(L/2)-1)/L;
plot(fVals, Px); % Plot
  3 Commenti
Jonathan
Jonathan il 23 Mar 2020
Hi, thank you for the response. The problem is that the plot function seems to failing because some of the y data is -inf
Peng Li
Peng Li il 23 Mar 2020
did you try my solution? If you add a small number like eps to Px, there won't be any -inf.

Accedi per commentare.

Risposte (1)

Stijn Haenen
Stijn Haenen il 23 Mar 2020
you can use Px(isinf(abs(Px)))=0; to set all -Inf values to 0
  3 Commenti
Stijn Haenen
Stijn Haenen il 23 Mar 2020
Modificato: Stijn Haenen il 23 Mar 2020
did you place it after the convertion to dB?
like this:
Px = Y.*conj(Y); % PSD Power of each freq components
Px = Px(1:L/2);
Px = 10*log10(Px); % Convert to dB
Px(isinf(abs(Px)))=0;
fVals = Fs * (0:(L/2)-1)/L;
plot(fVals, Px); % Plot
Does your Px also contains non real values?
Jonathan
Jonathan il 23 Mar 2020
Yes, it was placed there and no just real values. aside from a typo, it does plot now. The single positive is at 100Hz as expected. However, it plots also at every 40Hz a spectra that is at/about -250db for the full range of Fs/2 except where the -Inf were. I didn't think it possible to have a negative PSD. Round off/precision errors from the program that generated the audio file? Thanks.

Accedi per commentare.

Prodotti


Release

R2019b

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by