Manual Implementation of STFT of an audio signal.
18 次查看(过去 30 天)
显示 更早的评论
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 个评论
采纳的回答
Wayne King
2011-9-26
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 个评论
更多回答(5 个)
Wayne King
2011-9-22
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 个评论
Daniel Shub
2011-9-26
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
2011-9-23
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 个评论
Wayne King
2011-9-26
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 个评论
UJJWAL
2011-9-26
2 个评论
Daniel Shub
2011-9-26
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
2011-9-26
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
2011-11-6
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 个评论
另请参阅
类别
在 Help Center 和 File Exchange 中查找有关 Time-Frequency Analysis 的更多信息
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!