%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%
%     Compute the frequency response of the
%     Dolph-Chebyshev Filter.   
%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

Dt = 300;
M = 18;
N = 2*M+1;
TauS= 3*3600;
omegaS = (2*pi)/TauS;
thetaS = omegaS*Dt;
Tspan = (N-1)*Dt/3600;

fprintf(' Delta t  %g \n', Dt) 
fprintf(' M  N  %i   %i  \n', M,N) 
fprintf(' TauS     %g \n', TauS) 
fprintf(' omegaS   %g \n', omegaS) 
fprintf(' thetaS  %g  \n', thetaS) 
fprintf(' Tspan   %g  \n', Tspan) 

Dolph(1:19)   = [   0.03380; 
                    0.03370;
                    0.03342; 
                    0.03295; 
                    0.03230; 
                    0.03149;  
                    0.03049;  
                    0.02936;  
                    0.02809;  
                    0.02671;  
                    0.02522;  
                    0.02365;  
                    0.02201;  
                    0.02032;  
                    0.01860;  
                    0.01688;  
                    0.01517;  
                    0.01348;  
                    0.04928  ];

%%%   Compute the frequency response.
h = Dolph;
ntheta = 501;
theta = linspace(0,pi,ntheta);
H = zeros(1,ntheta);
H = H + h(1);
for m=1:M
  H = H + 2*h(m+1)*cos(m*theta);
end

r = max(H(100:ntheta));
fprintf(' Ripple   %g \n', r) 

%%%   Plot the frequency response.
plot(theta,H,'LineWidth',2)
pause
thnum = 160;
plot(theta(1:thnum),abs(H(1:thnum)),'LineWidth',2);
grid on;
axis([0,1,0,1]);
xlabel('Frequency   \theta','FontSize',12);
ylabel('Response   H','FontSize',12);
title('Frequency Response of Dolph-Chebyshev Filter','FontSize',16);

%  print -deps 'dolph.eps'

return
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
