%Demo done during Session 6 illustrating
%the utility of a simple second-order difference
%equation to notch out an undesired tone while
%preserving the integrity of the desired signal
%as best as possible.  This demo also illustrates
%the effects of zeroes and poles on the frequency
%response of an LTI System described by a difference 
%eqn. FIR vs IIR notch filtered is explored as well.
%Although this code is not commented, the input lines
%serve as comments that guide one through the code.
clear
clf
set(0,'defaultaxesfontsize',20);
data=getspeech('erf1s1t0.wav');
%data=getspeech('rpf1s1t0.wav');
Fs=6250;
omega0=pi/4;
NL=100;
plot(data)
end1=input('end1?');
end2=input('end2?');
datar=data(end1:end2);
[d,dsize]=size(datar);
datar2=datar(1:2:dsize);
dsize=dsize/2
x=datar2+500*cos(omega0*(1:length(datar2)));
input('Utterance played back at orignal sampling rate, Fs = 6.25 KHz');
soundsc(datar2,Fs)
plot(x)
input('Utterance with 1.5625 KHz tone interference played at Fs = 6.25 KHz');
soundsc(x,Fs)
r=0.95; c0=2*cos(omega0); rs=r^2; c0r=r*c0;
h2poles=r.^(1:NL).*sin(omega0*(1:NL));
hiir=conv(h2poles,[1 -c0 1]);
poles=roots([1 -c0r rs]);
zeroes=roots([1 -c0 1]);
polar(1,0,'.');
hold on
plot(real(poles),imag(poles),'x','Linewidth',4,'MarkerSize',18)
plot(real(zeroes),imag(zeroes),'o','Linewidth',4,'MarkerSize',18)
hold
input('Pole-zero plot -- NEXT: magnitude of freq. response for notch filter')
domega=2*pi/256;
omega=-pi:domega:pi-domega;
hf1=abs(fftshift(fft(hiir,256)));
plot(omega,hf1,'LineWidth',4)
axis([-pi pi 0 max(hf1)]);
title('Magnitude of DTFT of IIR notch filter');
xlabel('omega (radians/s)');
%
r=0.95; c0=2*cos(omega0); rs=r^2; c0r=r*c0;
input('run tone corrupted speech through simple second');
input('order difference equation');
y(1)=0; y(2)=0;
for n=3:length(x),
y(n)=c0r*y(n-1)-rs*y(n-2)+x(n)-c0*x(n-1)+x(n-2);
end
%
%y=conv(x,hiir);
input('IIR Notch-filtered utterance plus 1.5625 KHz tone interference');
soundsc(y,Fs)
z=conv(x,[1 -c0 1]);
hf2=abs(fftshift(fft([1 -c0 1],256)));
plot(omega,hf2,'LineWidth',4)
axis([-pi pi 0 max(hf2)]);
title('Magnitude of DTFT of FIR notch filter');
xlabel('omega (radians/s)');
input('FIR Notch-filtered utterance plus 1.5625 KHz tone interference');
soundsc(z,Fs)
input('Pole-zero plot -- NEXT: magnitude of freq. response for notch filter')
domega=2*pi/8192;
omega=-pi:domega:pi-domega;
daf1=abs(fftshift(fft(datar,8192)));
plot(omega,daf1)
axis([-pi pi 0 max(daf1)]);
title('Magnitude of DTFT of original utterance');
xlabel('omega (radians/s)');

