24 views (last 30 days)

I hope this will be fun for you to investigate and play with. Weird stuff happens. If weird stuff doesn't happen for you, then maybe my version of MATLAB (R2013a) is to blame.

I want to make a spectrogram of a noisy signal. Due to the signal's complexity, I've been messing around with a matlab example to understand how the inputs affect the output and what I expect to happen based on the help info vs. what I see actually plotted. This is the example (I added the 'yaxis' bit):

t=0:0.001:2; % 2 secs @ 1kHz sample rate

y=chirp(t,100,1,200,'q'); % Start @ 100Hz, cross 200Hz at t=1sec

spectrogram(y,128,120,128,1E3,'yaxis'); % Display the spectrogram

title('Quadratic Chirp: start at 100Hz and cross 200Hz at t=1sec');

I can tell, then, that they're using the form "spectrogram(X,WINDOW,NOVERLAP,NFFT,Fs,FREQLOCATION)."

So my understanding is that the 2nd argument is the number of time-axis segments to split the data into, 3rd is not super clear actually but the default is 50%. 50% of, I initially assumed, the adjacent two segments. The later statement of "NOVERLAP is the number of samples each segment of X overlaps." Ok, well in the default case of 8 segments and 50% of the samples, how can you make 8 segments when each segment is supposed to have 50% of the samples in it?? And why does it show 7 segments? This will show you what I mean:

spectrogram(y,[],[],128,1E3,'yaxis');

Then lets decide we want 128 segments like they did, and we'll just be ok with the default NOVERLAP. Well, plot that and you'll see 29 segments, which isn't 50% of anything I can see.

Now put NOVERLAP to 20, and we get 17 segments. Change it to 5 and you'll get 15! Plot the original again, then lets look at what happens when we double the number of segments and bump the NOVERLAP a bit:

spectrogram(y,256,200,128,1E3,'yaxis');

That doesn't look right! Maybe what I thought about NFFT was unclear. Lets double that...nope...at least that seems to behave properly. That's the number of segments the sample frequency spectrum is broken in to. Just for fun you can defy logic one last time with this:

spectrogram(y,256,200,50,1E3,'yaxis');

How am I supposed to select the appropriate values for this function? Is there something wrong with it? For my real case I'll have a sampling frequency of 72 kHz and will have applied a zero-phase Butterworth low pass filter with cutoff at about 7 kHz.

Shane L
on 26 Jun 2017

The 2nd argument (if it is a single integer), WINDOW, is not the number of time segments into which to split the data, but the length of each time segment. To see this, try, for example, a modified version of the example code:

t = 0:0.001:2-0.001; % 2 secs @ 1kHz sample rate

y = chirp(t,100,1,200,'q'); % start @ 100Hz, cross 200Hz at t = 1 sec

spectrogram(y,200,0,200,1e3,'yaxis'); % display the spectrogram

The signal has a length of 2000 samples. We set the second argument, WINDOW, to 200, which specifies that each segment contains 200 samples. We also set the 3rd input, NOVERLAP, to zero for this example. Since the signal with 2000 samples is split into segments of 200 samples, there are 10 time segments in the resulting spectrogram, as shown in the figure below.

The 3rd argument, NOVERLAP, specifies the number of samples that overlap between two adjacent segments in the calculation of the short-time Fourier transform. In the example above, we set NOVERLAP to zero, so no samples were overlapped. By default, 50% of the samples are overlapped. To see this, try, for example:

t = 0:0.001:2-0.001; % 2 secs @ 1kHz sample rate

y = chirp(t,100,1,200,'q'); % start @ 100Hz, cross 200Hz at t = 1 sec

spectrogram(y,200,[],200,1e3,'yaxis'); % display the spectrogram

Now, you will see 19 time segments in the resulting spectrogram, as shown in the figure below. The reason is because of the overlap. The first segment includes samples 1-200; the second segment includes samples 101-300 (50% overlap of 200 samples is 100 samples); the third segment includes samples 201-400; up to the 19th segment, which includes samples 1801-2000.

I don't know if this was part of your question, but the 4th argument, NFFT, is the number of points used in the calculation of the FFT for each segment. If NFFT is less than WINDOW, then not all of the points in each segment are used in the FFT. If NFFT is greater than WINDOW, then the segment is zero-padded before computing the FFT. Try changing the value of NFFT (for example, 10 vs. 1000) and plotting the spectrogram, and you will see the frequency resolution change, as shown in the two figures below.

For more information, you can consult the MATLAB help page for spectrogram ; several examples demonstrate how to compute the number of resulting time and frequency bins based on WINDOW, NOVERLAP, and NFFT.

Shane L
on 27 Jun 2017

Hi Mike, I'm glad that the first part was helpful for you. As to why you don't get a color bar with units, this functionality was introduced in R2015a (I was using R2015b to generate the plots). However, I do not know why your version is only plotting N-1 segments. Could you try running spectrogram with output arguments, as follows:

t = 0:0.001:2-0.001; % 2 secs @ 1kHz sample rate

y = chirp(t,100,1,200,'q'); % start @ 100Hz, cross 200Hz at t = 1 sec

[s,f,t] = spectrogram(y,200,0,200,1e3,'yaxis'); % display the spectrogram

What is the size of s (or alternatively, the length of t)?

Shane L
on 27 Jun 2017

Hi Mike, I got a chance to access R2014b. MATLAB introduced additional code into the R2015a version of spectrogram that shifts the time series by a half-interval and adds an additional point at the end, which explains the differences in the plots. See below for the block of code in R2015a:

% Surf encodes data points at the edges and takes the color from

% the last edge so we need to add an additional point so that surf

% does the right thing. This is important especially when

% spectrogram has only one estimate (e.g. window length = signal

% length). Although this issue also exists on the frequency

% direction we will not add an extra row since we will never

% encounter a single frequency point estimate. For the plot we set

% time values to be at: nwin/2-a/2, nwin/2+a/2, nwin/2+3a/2,

% nwin/2+5a/2 ... where a is the number of new samples for each

% segment (i.e. nwind-noverlap). For the case of zero overlap this

% corresponds to 0, nwind, 2*nwind, ...

a = nwind - noverlap;

t = [(nwind/2-a/2)/Fs t+((a/2)/Fs)];

Pxx = [Pxx Pxx(:,end)];

If you would like your plots to be in style of R2015a and later, then edit spectrogram and insert that block of code in between lines 195 and 196 (if the code in your version is the same as in R2014b, otherwise, find the code as shown below:

[Pxx,W] = compute_PSD(win,y,nfft,f,options,Fs,esttype); %#ok<NASGU>

%%%INSERT LINES HERE

displayspectrogram(t,f,double(Pxx),isFsnormalized,faxisloc);

I tested my first example code using R2014b and it now matches R2015a. I hope this helps!

Sign in to comment.

Sign in to answer this question.

Opportunities for recent engineering grads.

Apply Today
## 0 Comments

Sign in to comment.