annotate extra codes/stft_primer_ejemplo.m @ 3:1c0f36c348d4

extra code for matlab
author Katerina <katkost@gmail.com>
date Sat, 20 Apr 2013 13:03:01 +0100
parents
children
rev   line source
katkost@3 1 function y = stft(x, w, N, H)
katkost@3 2 % Analysis/synthesis of a sound using the short-time fourier transform
katkost@3 3 % x: input sound, w: analysis window (odd size), N: FFT size, H: hop size
katkost@3 4 % y: output sound
katkost@3 5 M = length(w); % analysis window size
katkost@3 6 N2 = N/2+1; % size of positive spectrum
katkost@3 7 soundlength = length(x); % length of input sound array
katkost@3 8 hM = (M-1)/2; % half analysis window size
katkost@3 9 pin = 1+hM; % initialize sound pointer in middle of analysis window
katkost@3 10 pend = soundlength-hM; % last sample to start a frame
katkost@3 11 fftbuffer = zeros(N,1); % initialize buffer for FFT
katkost@3 12 yw = zeros(M,1); % initialize output sound frame
katkost@3 13 y = zeros(soundlength,1); % initialize output array
katkost@3 14 w = w/sum(w); % normalize analysis window
katkost@3 15 while pin<pend
katkost@3 16 %-----analysis-----%
katkost@3 17 xw = x(pin-hM:pin+hM).*w(1:M); % window the input sound
katkost@3 18 fftbuffer(:) = 0; % reset buffer
katkost@3 19 fftbuffer(1:(M+1)/2) = xw((M+1)/2:M); % zero-phase window in fftbuffer
katkost@3 20 fftbuffer(N-(M-1)/2+1:N) = xw(1:(M-1)/2);
katkost@3 21 X = fft(fftbuffer); % compute FFT
katkost@3 22 mX = 20*log10(abs(X(1:N2))); % magnitude spectrum of positive frequencies
katkost@3 23 pX = unwrap(angle(X(1:N2))); % unwrapped phase spect. of positive freq.
katkost@3 24 %-----synthesis-----%
katkost@3 25 Y = zeros(N,1); % initialize output spectrum
katkost@3 26 Y(1:N2) = 10.^(mX/20).*exp(i.*pX); % generate positive freq.
katkost@3 27 Y(N2+1:N) = 10.^(mX(N2-1:-1:2)/20).*exp(-i.*pX(N2-1:-1:2));
katkost@3 28 % generate neg.freq.
katkost@3 29 fftbuffer = real(ifft(Y)); % inverse FFT
katkost@3 30 yw(1:(M-1)/2) = fftbuffer(N-(M-1)/2+1:N); % undo zero-phase window
katkost@3 31 yw((M+1)/2:M) = fftbuffer(1:(M+1)/2);
katkost@3 32 y(pin-hM:pin+hM) = y(pin-hM:pin+hM) + H*yw(1:M); % overlap-add
katkost@3 33 pin = pin+H; % advance sound pointer
katkost@3 34 end