katkost@3: function y=stpt2(x,w,N,H,t) 
katkost@3: %Analysis/synthesis of a sound using the peaks 
katkost@3: %of the short-time Fourier transform 
katkost@3: %x: input sound, w: analysis window (odd size), N: FFT size, H: hop 
katkost@3: %size, t: threshold in negative dB, y: output sound 
katkost@3:   
katkost@3: %Internal variables 
katkost@3: M=length(w); 
katkost@3: N2=N/2+1; 
katkost@3: soundlength=length(x); 
katkost@3: hM=(M-1)/2; 
katkost@3:   
katkost@3: %Pointers 
katkost@3: pin=1+hM; 
katkost@3: pend=soundlength-hM; 
katkost@3:   
katkost@3: %Initialize buffers 
katkost@3: fftbuffer=zeros(N,1); 
katkost@3: yw=zeros(M,1); 
katkost@3: y=zeros(soundlength,1); 
katkost@3:   
katkost@3: %Analysis/Synthesis windows 
katkost@3: sw=hanning(M); 
katkost@3: sw=sw./sum(sw); 
katkost@3: w=w./sum(w); 
katkost@3:   
katkost@3: while pin<pend 
katkost@3:     %Analysis 
katkost@3:     xw=x(pin-hM:pin+hM).*w(1:M); %Frame selector 
katkost@3:     fftbuffer=zeros(N,1);        %Zero-phase windowing 
katkost@3:     fftbuffer(1:(M+1)/2)=xw((M+1)/2:M); 
katkost@3:     fftbuffer(N-(M-1)/2+1:N)=xw(1:(M-1)/2); 
katkost@3:      
katkost@3:     X=fft(fftbuffer);   %FFT performance 
katkost@3:     mX=20*log10(abs(X(1:N2))); %Compute the FFT Magnitude (dB) 
katkost@3:     pX=angle(X(1:N2));  %Compute the FFT phase 
katkost@3:      
katkost@3:     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:     pmag=mX(ploc);     %Peak magnitude recovery 
katkost@3:     pphase=pX(ploc);   %Peak phase recovery 
katkost@3:  
katkost@3:     %Synthesis 
katkost@3:     Y=zeros(N,1);    %Memory reserve 
katkost@3:     Y(ploc)=10.^(pmag/2).*exp(1i.*pphase); %Spectrum reconstruction from peak info 
katkost@3:     Y(N+2-ploc)=10.^(pmag/20).*exp(-1i.*pphase); %Spectrum reconstruction negative freq. 
katkost@3:     fftbuffer=real(ifft(Y));     %Compute inverse FFT 
katkost@3:     yw((M+1)/2:M)=fftbuffer(1:(M+1)/2); %Zero-phase windowing undo 
katkost@3:     yw(1:(M-1)/2)=fftbuffer(N-(M-1)/2+1:N); 
katkost@3:     y(pin-hM:pin+hM)=y(pin-hM:pin+hM)+H*N*sw.*yw(1:M); %overlap-add 
katkost@3:     pin=pin+H;                   %Pointer advance 
katkost@3: end 
katkost@3:      
katkost@3: