MATLAB and C Programs - Digital Signal Processing: A Practitioner's Approach - page 149

subplot(211)
plot(f,mn);
xlabel ('Normalised Frequency')
ylabel('Magnitude')
grid;
pause;
%print -depsc f31_7
%gtext('r = 0.98');

A.3.12 Program f 31_90

%Digital Signal Processing:A Practitioner's Approach
%Dr.Kaluri Venkata Ranga Rao, [e-mail]
%Generates figures 3.23,3.24

clear
close
range=2;
N=200;
delp=range/N;
p=-1;
r=0.95;
r1=0.5*(1-r*r);
k=1:150;
f=0.2;
power=0.5;
u=sin(2*pi*f*k);
noise=power*randn(size(u));
u=u+0*noise;
for i=1:N
b1=[r1 0 -r1];
b2=[0 2*r];
a=[1 -2*r*p r*r];
p=p+delp;
fhat=acos(p)/(2*pi);
x=filter(b1,a,u);
sk=filter(b2,a,x);
ek=x-u;
amp(i) = mean(ek.*ek);

% F(i)=p;

F(i)=fhat;
slope(i)=2*mean(ek.*sk);

end
subplot(211)