% dft signalu - pocitani a zobrazeni. 

% chci 1 znely recovy ramec: 
s = wavread('test.wav')'; 
sfr = frame (s,160,80); 
x = sfr(:,13); 
plot (x); print -depsc FIG2/sig.eps
Fs = 8000; f = (0:159) / 160 * Fs
X = fft(x); 
subplot (211); plot(f,abs(X)); 
subplot (212); plot(f,angle(X)); 
print -depsc FIG2/do_fs.eps
%
Fs = 8000; f = (0:79) / 160 * Fs
X = fft(x); X = X(1:80); 
subplot (211); plot(f,abs(X)); 
subplot (212); plot(f,angle(X)); 
print -depsc FIG2/do_fs2.eps
% vice bodu ve spek. 
Fs = 8000; f = (0:511) / 1024 * Fs
X = fft([x' zeros(1,1024-160)]); X = X(1:512); 
subplot (211); plot(f,abs(X)); 
subplot (212); plot(f,angle(X)); 
print -depsc FIG2/padded.eps

%%%  PSD - stejne triky jako pro DFT
Fs = 8000; f = (0:511) / 1024 * Fs
X = fft([x' zeros(1,1024-160)]); X = X(1:512); Gdft= 1/160*abs(X) .^ 2;
plot(f,Gdft); 
print -depsc FIG2/Gdft.eps

plot(f,10 * log10 (Gdft)); 
print -depsc FIG2/Gdftlog.eps



% case study - potrebuji omezeni audio signalu na tlf. kanal. 
% signal z radia ... 

% - linkuji sem kus sinfonietty. 
x = wavread ('sinf.wav'); 
% jen 1 kanal a radeji radkovy vektor ... 
x = (x(:,1) + x(:,2))'; 
wavwrite(x,44100,16,'janacek.wav'); 

%%% navrh filtru pro omezeni na tlf. kanal 300-3400 Hz. 
Fs = 44100; Fs2 = Fs/2; 
Wp = [300/Fs2 3400/Fs2]; 
Ws = [200/Fs2 3500/Fs2];
Rp = 3;
Rs = 30;
[N, Wn] = ellipord(Wp, Ws, Rp, Rs)
[B,A] = ellip(N,Rp,Rs,Wn)
% overeni 
freqz (B,A,512,Fs); 
slidify(1,3,1,15); print -depsc FIG2/tel_filter.eps
zplane (B,A); 
slidify(1,3,1,15); print -depsc FIG2/tel_np.eps

y=filter(B,A,x); 
wavwrite(y,44100,16,'out.wav');