annotate extra codes/stpt2.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=stpt2(x,w,N,H,t)
katkost@3 2 %Analysis/synthesis of a sound using the peaks
katkost@3 3 %of the short-time Fourier transform
katkost@3 4 %x: input sound, w: analysis window (odd size), N: FFT size, H: hop
katkost@3 5 %size, t: threshold in negative dB, y: output sound
katkost@3 6
katkost@3 7 %Internal variables
katkost@3 8 M=length(w);
katkost@3 9 N2=N/2+1;
katkost@3 10 soundlength=length(x);
katkost@3 11 hM=(M-1)/2;
katkost@3 12
katkost@3 13 %Pointers
katkost@3 14 pin=1+hM;
katkost@3 15 pend=soundlength-hM;
katkost@3 16
katkost@3 17 %Initialize buffers
katkost@3 18 fftbuffer=zeros(N,1);
katkost@3 19 yw=zeros(M,1);
katkost@3 20 y=zeros(soundlength,1);
katkost@3 21
katkost@3 22 %Analysis/Synthesis windows
katkost@3 23 sw=hanning(M);
katkost@3 24 sw=sw./sum(sw);
katkost@3 25 w=w./sum(w);
katkost@3 26
katkost@3 27 while pin<pend
katkost@3 28 %Analysis
katkost@3 29 xw=x(pin-hM:pin+hM).*w(1:M); %Frame selector
katkost@3 30 fftbuffer=zeros(N,1); %Zero-phase windowing
katkost@3 31 fftbuffer(1:(M+1)/2)=xw((M+1)/2:M);
katkost@3 32 fftbuffer(N-(M-1)/2+1:N)=xw(1:(M-1)/2);
katkost@3 33
katkost@3 34 X=fft(fftbuffer); %FFT performance
katkost@3 35 mX=20*log10(abs(X(1:N2))); %Compute the FFT Magnitude (dB)
katkost@3 36 pX=angle(X(1:N2)); %Compute the FFT phase
katkost@3 37
katkost@3 38 ploc=1+find((mX(2:N2-1)>t).*(mX(2:N2-1)>mX(3:N2)).*(mX(2:N2-1)>mX(1:N2-2))); %Peak finder above threshold (t)
katkost@3 39 pmag=mX(ploc); %Peak magnitude recovery
katkost@3 40 pphase=pX(ploc); %Peak phase recovery
katkost@3 41
katkost@3 42 %Synthesis
katkost@3 43 Y=zeros(N,1); %Memory reserve
katkost@3 44 Y(ploc)=10.^(pmag/2).*exp(1i.*pphase); %Spectrum reconstruction from peak info
katkost@3 45 Y(N+2-ploc)=10.^(pmag/20).*exp(-1i.*pphase); %Spectrum reconstruction negative freq.
katkost@3 46 fftbuffer=real(ifft(Y)); %Compute inverse FFT
katkost@3 47 yw((M+1)/2:M)=fftbuffer(1:(M+1)/2); %Zero-phase windowing undo
katkost@3 48 yw(1:(M-1)/2)=fftbuffer(N-(M-1)/2+1:N);
katkost@3 49 y(pin-hM:pin+hM)=y(pin-hM:pin+hM)+H*N*sw.*yw(1:M); %overlap-add
katkost@3 50 pin=pin+H; %Pointer advance
katkost@3 51 end
katkost@3 52
katkost@3 53