Manual Implementation of STFT of an audio signal.
18 Ansichten (letzte 30 Tage)
Ältere Kommentare anzeigen
Hi ,
I have an Audio File whose spectrogram I want to find manually by calculating its STFT without using the inbuilt Spectrogram function. My steps can be mentioned as follows :-
a) I take samples out of the file in chunks of 1000 samples each and process them one by one. b) For each chunk , I find its STFT . c) I plot the magnitute squared of the STFT and hold that plot. d) I repeat the whole process again and again till all the samples of the audio file have been used.
However I face the following big problem :-
The process is never complete as everytime the computer runs out of memory. I try to use as few variables as possible and i clear the variables as soon as they are used. I also use the pack command before executing the program. However each time the computer almost hangs after some iterations. I think this is because the hold on command causes MATLAB to keep storing its previous plot data.
Is there any way around this issue ? Please enlighten. I am using an audio file of length 3.5 seconds at a sampling frequency of 44100 Hz.
0 Kommentare
Akzeptierte Antwort
Wayne King
am 26 Sep. 2011
Here's what you have to do, something like this:
dt = 1/44100;
t=0:dt:(117204*dt)-dt;
x = chirp(t,1500,1,8000);
S = zeros(501,117);
windowlength = 1e3;
k = 1;
for jj = 1:117
y = x(k:k+windowlength-1)';
ydft = fft(y).*gausswin(1e3);
S(:,jj) = ydft(1:501);
k = k+windowlength;
end
F = 0:44100/1000:44100/2;
T = 0:(1e3*dt):(117000*dt)-(1e3*dt);
surf(T,F,20*log10(abs(S)),'EdgeColor','none')
axis xy; axis tight; colormap(jet); view(0,90);
xlabel('Time');
ylabel('Frequency (Hz)');
Like I said, the segments do not overlap and I haven't been too careful with somethings, like whether I skip a sample between adjacent segments, but this is the general idea.
0 Kommentare
Weitere Antworten (5)
Wayne King
am 22 Sep. 2011
Why are you plotting each segment and holding it? I think you should store the DFT of the individual segments as column vectors (or row vectors) in a matrix and then make a surf plot when you are done.
You will also want to overlap your segments and not block them.
6 Kommentare
Daniel Shub
am 26 Sep. 2011
Have you looked at my answer? It seems to do what you want minus the windowing. Adding the window in is relatively easy (you need to multiple each column of my z with your window vector) or replace ones(1000,1) in the spectrogram method with a window vector.
Daniel Shub
am 23 Sep. 2011
I am confused ...
Is this roughly what you are trying to do?
x = randn(117424,1); % Make some fake data of the correct size
y = [x; zeros((1000-mod(length(x), 1000)), 1)];
z = reshape(y, 1000, length(y)/1000);
Z = fft(z, 2^16);
mesh(abs(Z))
Then again, I think that is roughly
spectrogram(x, ones(1000, 1), 0, 2^16, 44.1e3)
0 Kommentare
Wayne King
am 26 Sep. 2011
Hi, I think you are misunderstanding some things about the STFT. Why do you think the sampling frequency has anything to do with the length of the short-time Fourier transforms? The length of the short-time Fourier transforms depends on the length of the segments you choose to extract from your data.
Also, the number of columns is NOT equal to the number of points in your data vector.
Below I use a Gaussian window of length 1000 and overlap the segments by 900 points. Look at the size of S, the spectrogram matrix.
win = gausswin(1e3);
dt = 1/44100;
t=0:dt:(117204*dt)-dt;
x = chirp(t,1500,1,8000);
[S,F,T,P] = spectrogram(x,win,900,length(win),44100);
surf(T,F,10*log10(P),'EdgeColor','none')
axis xy; axis tight; colormap(jet); view(0,90);
xlabel('Time');
ylabel('Frequency (Hz)');
0 Kommentare
UJJWAL
am 26 Sep. 2011
2 Kommentare
Daniel Shub
am 26 Sep. 2011
Again, have you looked at my answer? I am not sure why you don't want to use spectrogram, but I gave you a solution that doesn't require it.
Wayne King
am 26 Sep. 2011
Irrespective of whether you want to write your own STFT or not, you are not correctly understanding what the size of the output matrix is.
Let's say you don't overlap your time segments (which will result in a blocky-looking spectogram, but ok), if you want to take the data 1,000 samples points at a time, then you will only have about 117 columns in your matrix. If your data is real-valued, then each column will only have 501 points, the positive frequencies of the fft() of each of your 1,000 sample segments.
raj
am 6 Nov. 2011
hello friends i too have the same question, i have an audio file and i want to implement stft manually but i have a problem with the overlap i implemented an overlap using if elseif branch but the window size doesnt change ...please help me with it i am pasting my code here
clear all
close all
clc
fs = 50e3;
%t = 1/fs;
w1 = 28000;
feed = w1/2;
load srtm
x1 = (srtm(8,1:10000));
nx = length(x1);% size of signal
w = hamming(w1)';% hamming window
nw = length(w);
T = (1:nx/fs); %size of window
pos=1;
n = 1;
%Null Matrices
W_total1 =[];%nullmatrix
while (pos+nw <= nx)
% while enough signal left
if n == 1
w1start = 0;
w1(end) = w1;
elseif n==2
w1start = (n-1)*w1 - feed;
w1(end) = (w1 + feed);
else
w1start = (w1start + w1(end))/2;
w1(end) = w1start + w1;
end
z(n,:) = x1(pos:pos+nw-1).*w(:,28000);
n = n+1;
pos = pos + nw/2;
W = fft(z(n-1,:));
W = W(1:end/2);
W_total1 = [W_total1 W.'];
end
figure;
imagesc(T,1:fs/2,10*log10(abs(W_total1/2e-5)))
axis xy; axis tight; colormap(jet);
ylim([0 1000])
xlabel('Time');
ylabel('Frequency (Hz)');
colorbar
0 Kommentare
Siehe auch
Kategorien
Mehr zu Time-Frequency Analysis finden Sie in Help Center und File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!