Why are the results of my spectrogram sinusoidal?

6 次查看（过去 30 天）
Harvey
Harvey 2016-4-21
评论： Harvey 2016-4-22
I am trying to determine the amplitude envelope of specific frequencies over time, from a sample of an instrument (a trumpet). I use the spectrogram function to find the amplitude of each frequency bin over time. Knowing that a harmonic lies within a certain bin, I can plot that particular bin against time to produce the following figure:
Although I can make out the general envelope, as marked in red, this method is much less accurate than what I was expecting.
I can't understand why the results are sinusoidal. Is this a result of aliasing, or something else entirely?
If we compare this plot, to a plot of all frequencies, we notice that there is no sinusoidal fluctuation in the power of the particular harmonic (The harmonic plotted in the first figure has been boxed):
So where is difference in the values used to plot the consistent power we see here, and the values that produce the sinusoidal plot?
For reference I have attached the audio file I am using. Bellow are the functions I used to get these results. Row 13 refers to the bin containing the specific harmonic;
[y,fs]=audioread('Trumpet_C4.wav');
[s,f,t]=spectrogram(y,hann(2048),1024,2048,fs,'yaxis');
s_real=real(s);
plot(t,s_real(13,:))
figure(2)
spectrogram(y,hann(2048),1024,2048,fs,'taxi's')
Any help would be extremely gratefully received!

请先登录，再回答此问题。

采纳的回答

J. Webster
J. Webster 2016-4-21
编辑：J. Webster 2016-4-21
You don't want just the real part of the fft, you want the magnitude of it.
change
s_real = real(s)
to
s_mag = abs(s)
and use that in your plot.
  3 个评论
显示 1更早的评论
J. Webster
J. Webster 2016-4-21
编辑：J. Webster 2016-4-21
I believe that what's happening is that when you call spectrogram with no output arguments, it plots the PSD of the signal, but when you call it with output arguments, it returns the complex sffts of the signal.
Note that the abs(s) is not equal to the PSD, it just has the same behavior. If you want the PSD instead of just the magnitude, I think you can just use
PSD = 1/(fs*2048) .* s_mag.^2 %2048 is the number of points in the fft

请先登录，再进行评论。

更多回答（1 个）

Star Strider
Star Strider 2016-4-21
I am not certain what you want to do. Just filtering the fundamental frequency gives it as a function of time, and the envelope is relatively easy to calculate:
[y,fs]=audioread('Trumpet_C4.wav');
t = linspace(0, length(y), length(y))/fs; % Time Vector
Fn = fs/2; % Nyquist Frequency
Wp = [0.2 0.4]*1E+3/Fn; % Filter Passband
Ws = [0.5 2.0].*Wp; % Filter Stopband
Rp = 10; % Passband Ripple
Rs = 45; % Stopband Ripple
[n,Wn] = buttord(Wp, Ws, Rp, Rs); % Filter Order
[b,a] = butter(n,Wn); % Transfer Function Coefficient Vectors
[sos,g] = tf2sos(b,a); % Second-Order-Section (For Stability)
figure(1)
freqz(sos, 4028, fs); % Filter Bode Plot
yf = filtfilt(sos,g, y); % Filter Signal (Phase-Neutral)
env = abs(hilbert(yf)); % Calculate Envelope
figure(2)
subplot(2,1,1)
plot(t, y)
title('Unfiltered Signal')
grid
axis([0 0.2 ylim]) % Axis Resolution (Comment-Out To See Entire Signal)
subplot(2,1,2)
plot(t, yf)
hold on
plot(t, env, '-r', t, -env, '-r')
hold off
title('Bandpass-Filtered Signal')
grid
axis([0 0.2 ylim]) % Axis Resolution (Comment-Out To See Entire Signal)
xlabel('Time (s)')
It’s a signal with a fundamental frequency and (apparently 13) harmonics, so you would have to filter each one, and calculate the envelope for each one. They seem to be evenly-spaced, so you could determine the spacing and then design a bank of filters (see How to shape the spectrum of an audio waveform? to filter each frequency.
  3 个评论
显示 1更早的评论
Star Strider
Star Strider 2016-4-21
My pleasure.
There are several ways to get the individual harmonics. I would use the filter-bank method to get them all as individual records. I’ve had success with it in other applications, and it probably produces the most accurate results for what you want to do.
I went through the spectrogram documentation and couldn’t figure out the reason it gave you the results it did. I suspect it has to do with windowing, overlap and other necessary operations it has to do to calculate its results.
Harvey
Harvey 2016-4-21
I will try the filter bank method. Strange that the spectrogram appears steady in figure 2, but could be the windowing.
Thanks

请先登录，再进行评论。

类别

在 Help Center 和 File Exchange 中查找有关 Spectral Measurements 的更多信息

标签

产品

Community Treasure Hunt

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

Start Hunting!

Translated by