[x,fs]=audioread('flute.wav');

n = [0:length(x)-1];
t = n/fs;

ax1 = subplot(311);
plot(t,x);
Nfft = 256;
Nstep = 1024;
win = rectwin(Nfft);
%win = hamming(Nfft);
freq = [0:Nfft-1]/Nfft-0.5;
for i=2e4:Nstep:length(x)
    xlim(ax1,[i i+Nfft-1]/fs);
    y = x(i:i+Nfft-1).*win;
    Y = fft(y);
    [Ys,I] = sort(abs(Y),'descend');
    for k=2:2:Nfft-1
        subplot(312);
        plot(freq,abs(Y));
        hold on;
        plot(freq(I(1:k)),abs(Y(I(1:k))),'r.');
        hold off;
        subplot(313);
        Z = Y;
        Z(I(k+1:end)) = 0;
        z = real(ifft(Z));
        plot([0:Nfft-1],[y,z]);
        drawnow;
        pause;
    end
    
end

