0
$\begingroup$

I want to make a plot of the power spectrum for a .wav sound file, over the frequency range from -2000 to 2000 Hz.

Attempt:

I used the following code, but the resulting plot is one sided. It doesn't plot the negative frequencies. But the graph I need to produce has to contain negative frequencies. What do I need to do?

Any help would be appreciated.

Here is my code:

[s, Fs] = wavread('chord.wav');
nfft = 2^nextpow2(length(s));
Pxx = abs(fft(s,nfft)).^2/length(s)/Fs;
Hpsd = dspdata.psd(Pxx(1:length(Pxx)/2),'Fs',Fs);  
plot(Hpsd)
axis([-2 2 -150 -35])  %change axis range

Produces:

enter image description here

$\endgroup$
6
  • $\begingroup$ Check this link from Mathworks for dspdata.psd. $\endgroup$
    – Gilles
    Apr 26, 2016 at 9:03
  • $\begingroup$ I have already checked the link, and I used the syntax they give for a double-sided spectrum and I still get the same one sided plot. I used Hpsd=dspdata.psd(Pxx(1:length(Pxx)/2),'Fs',Fs,'SpectrumType','twosided'); . What is wrong here? $\endgroup$
    – Merin
    Apr 26, 2016 at 9:24
  • $\begingroup$ You're plotting your Pxx from 1 to length(Pxx)/2. You should plot the full of it. i.e. from 1 to length(Pxx) or simply Pxx with no indexes and use the twosided argument. $\endgroup$
    – Gilles
    Apr 26, 2016 at 9:28
  • $\begingroup$ Okay. So I tried dspdata.psd(Pxx,'Fs',Fs,'SpectrumType','twosided');. But it still returns the same one sided graph. Why is that? So, if we can't do this using psd, can we get the power spectrum by directly plotting $Pxx^2$ versus frequency? $\endgroup$
    – Merin
    Apr 26, 2016 at 10:19
  • $\begingroup$ Your negative frequencies are now from $2-4 \textrm{ kHz}$. Remove the line axis([-2 2 -150 -35]). $\endgroup$
    – Gilles
    Apr 26, 2016 at 10:23

1 Answer 1

1
$\begingroup$

First of all, wavread has been removed in the 2015 & 2016 versions of MATLAB. So, if you plan on using your code in a newer version of MATLAB, I recommend that you use the audioread function instead.

Now, to answer your question. MATLAB already have a built-in function for plotting the power spectral density (Pxx) of a given input signal - periodogram. I've taken the liberty to copy a piece of sample code from MATLAB's documentation that should help you plot the Pxx you want (it will contain negative frequencies as well).

periodogram(x,[],NFFT,Fs,'centered')

All you need to keep from your code snippet are the input values seen here. The rest, like the axis function, aren't necessary to get the plot you want if you use this periodogram function. The above function will plot the power spectral density of a given input signal x with no specified window ([] which just means the function defaults to a simple rectangular window). The final argument is the one that I presume you are most interested in: freqrange. It can have several different inputs, including twosided and centered. You should try both options and see which output you find more useful.

I hope this helps you!

$\endgroup$
4
  • $\begingroup$ Thank you. I am actually using Matlab 2012 version so wavread is usable. Here's my code now: [s, Fs] = wavread('chord.wav'); nfft = 2^nextpow2(length(s)); periodogram(s,[-2000:2000],nfft,Fs,'centered'). But I get the error: "Ambiguous or invalid string option specified". Did I use an incorrect syntax? $\endgroup$
    – Merin
    Apr 26, 2016 at 14:38
  • $\begingroup$ My apologies, I never checked if the commands I gave above are valid in earlier MATLAB versions. So, the problem you just ran into is not your fault: what happened is that the 'centered' option does not exist in the MATLAB version you are using, so as expected an error was thrown. Here is how you fix that error: plot(psd(spectrum.periodogram, x, 'Fs', Fs, 'NFFT', nfft, 'spectrumtype', 'twosided', 'centerdc', true)) ... you can find the full list of syntax replacements MATLAB made here: mathworks.com/help/signal/ug/periodogram-psd-.html $\endgroup$ Apr 26, 2016 at 14:46
  • $\begingroup$ Thank you so much, it worked perfectly. I have one last question: In plot(psd(spectrum.periodogram, x, 'Fs', Fs, 'NFFT', nfft, 'spectrumtype', 'twosided', 'centerdc', true)), where can you specify the x and y range? I want it to be plotted from -2k to 2k, rather than me zooming in manually. $\endgroup$
    – Merin
    Apr 26, 2016 at 20:22
  • $\begingroup$ I know it may be tempting to set the window parameter within the psd function to the x-axis and y-axis ranges you want to limit your plot to. However, that parameter is there to specify a window function (read up about their application to spectral analysis: en.wikipedia.org/wiki/Window_function#Spectral_analysis). You don't want that, all you need to do is specify the axis limits. You already had that set up fine in the code snippet you provide in the original question! Add this line back into your code: axis([-2 2 -150 -30]);. Zooming manually is a bother :P $\endgroup$ Apr 26, 2016 at 20:55

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.