why does spectrogram of stft shows different magnitude while changing the window size?

40 Ansichten (letzte 30 Tage)
fs=5000; % Sampling rate
N=1000; % Number of data points
y=sin(500*pi*[0:1:N-1]/fs)...
+2*sin(500*1.5*pi*[0:1:N-1]/fs)...
+3*sin(500*3*pi*[0:1:N-1]/fs)...
+4*sin(500*4*pi*[0:1:N-1]/fs);t=[0:1:N-1]/fs;
figure(1);
windowsize = 200;
window = boxcar(windowsize);
nfft = windowsize;
noverlap = windowsize-1;
[S,F,T] = spectrogram(y,window,noverlap,nfft,fs);
imagesc(T,F,(abs(S)))
set(gca,'YDir','Normal')
xlabel('Time (secs)')
ylabel('Freq (Hz)')
title('Short-time Fourier Transform spectrum')
For the above signal, the expected amplitudes of frequencies are in the range 1-4. once we i do the stft with window size 200, i am getting amplitudes of frequency in the range of 0-400 Hz. if i reduce the window size to 16, the amplitudes moving to the range of 0-45 Hz. Did i do anything wrong in my code or am i conceptually wrong somewhere?

Akzeptierte Antwort

Mathieu NOE
Mathieu NOE am 17 Mär. 2021
hello
I was tired of not being able to understand why some matlab functions (mostly based on fft) would not scale the output as I was expecting.
Finally decided (long time ago) to rewritte the spectrogram function myself - now the amplitude of the spectrogram are ok with the input data. I only decided to simplify the function, so it's only hanning window for the time being. If you change the type of window , you have to correct for the new correction coefficient.
also the overlap is a percentage of nfft (most of the time I leave it at 75%)
hope it helps
fs=5000; % Sampling rate
N=1000; % Number of data points
t=[0:1:N-1]/fs;
y=sin(500*pi*t)...
+2*sin(500*1.5*pi*t)...
+3*sin(500*3*pi*t)...
+4*sin(500*4*pi*t);
figure(1);
windowsize = 200;
window = boxcar(windowsize);
nfft = windowsize;
noverlap = windowsize-1;
% [S,F,T] = spectrogram(y,window,noverlap,nfft,fs); % S = SPECTROGRAM(X,WINDOW,NOVERLAP,NFFT,Fs)
[S,F,T] = myspecgram(y, fs, nfft, 0.75); % overlap is 75% of nfft here
imagesc(T,F,(abs(S)))
colorbar('vert');
set(gca,'YDir','Normal')
xlabel('Time (secs)')
ylabel('Freq (Hz)')
title('Short-time Fourier Transform spectrum')
colormap('jet');
function [fft_specgram,freq_vector,time] = myspecgram(signal, Fs, nfft, Overlap)
% FFT peak spectrogram of signal (example sinus amplitude 1 = 0 dB after fft).
% signal - input signal,
% Fs - Sampling frequency (Hz).
% nfft - FFT window size
% Overlap - buffer overlap % (between 0 and 0.95)
signal = signal(:);
samples = length(signal);
% fill signal with zeros if its length is lower than nfft
if samples<nfft
s_tmp = zeros(nfft,1);
s_tmp((1:samples)) = signal;
signal = s_tmp;
samples = nfft;
end
% window : hanning
window = hanning(nfft);
window = window(:);
% compute fft with overlap
offset = fix((1-Overlap)*nfft);
spectnum = 1+ fix((samples-nfft)/offset); % Number of windows
% % for info is equivalent to :
% noverlap = Overlap*nfft;
% spectnum = fix((samples-noverlap)/(nfft-noverlap)); % Number of windows
% main loop
fft_specgram = [];
for ci=1:spectnum
start = (ci-1)*offset;
sw = signal((1+start):(start+nfft)).*window;
fft_specgram = [fft_specgram abs(fft(sw))*4/nfft]; % X=fft(x.*hanning(N))*4/N; % hanning only
end
% one sidded fft spectrum % Select first half
if rem(nfft,2) % nfft odd
select = (1:(nfft+1)/2)';
else
select = (1:nfft/2+1)';
end
fft_specgram = fft_specgram(select,:);
freq_vector = (select - 1)*Fs/nfft;
% time vector
% time stamps are defined in the middle of the buffer
time = ((0:spectnum-1)*offset + round(offset/2))/Fs;
end

Weitere Antworten (0)

Kategorien

Mehr zu Time-Frequency Analysis finden Sie in Help Center und File Exchange

Tags

Community Treasure Hunt

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

Start Hunting!

Translated by