fg_06_05.m

英文书《Digital Signal Processing with Examples in MATLAB》附带的MATLAB实例

% fg_06_05  Analog Chebyshev poles.
L=[2 4 6];
clr=['r' 'g' 'b'];
sp_fig(1,8.5,3);
for plt=1:3,
   [z,p]=chebyshev_roots(L(plt),1,2);
   b=1;
   a=poly([p,-p]);
   subplot(1,3,plt);
%
% Altered version of function pz_plot(b,a) follows:
% pz_plot(b,a)
% Plots zeros and poles on the z plane in figure 1.
% b and a are row-vector coefficients in
%        b(1)+b(2)*z^(-1)+...+b(N)*z^(-(N-1))
% H(z) = ------------------------------------
%        a(1)+a(2)*z^(-1)+...+a(M)*z^(-(M-1))
% b and a can have multiple rows if H(z) is in cascade form.
% zeros ('o') and poles ('x') are plotted on the z plane.
[Nsect,N]=size(b);
[Nsecta,M]=size(a);
if(Nsect~=Nsecta)
   error('pz_plot error: b and a must have same no. rows.')
end
% Plot the coordinates.
circle=exp(j*2*pi*linspace(0,1,1000));
cla; plot(circle,'k'); hold on;
set(gca,'fontname','times','fontsize',16);
if(plt==2),
   xlabel('Re({s})');
elseif(plt==1),
   ylabel('Im({s})');
end
n=linspace(-1,1,1000);
Ncircles=5;
for c=2:Ncircles
   plot(circle*(c-1)/Ncircles, ':k')
end
Nlines=4;
for q=1:Nlines
   plot(n*cos(pi*q/Nlines),n*sin(pi*q/Nlines),':k')
end
% Plot the zeroes. lim= outer limit of plot axes.
lim=1;
if N>1
   for s=1:Nsect
      z=roots(b(s,:))';
      lim=max([lim,abs(real(z)),abs(imag(z))]);
      for n=1:N-1
         hz(n)=plot(real(z(n)),imag(z(n)),'o');
         set(hz(n),'MarkerSize',12,'color',clr(plt));
         set(hz(n),'LineWidth',1.5);   
      end
   end
end
% Plot the poles.
if M>1
   for s=1:Nsect
      p=roots(a(s,:))';
      lim=max([lim,abs(real(p)),abs(imag(p))]);
      for n=1:M-1
         hp=plot(real(p(n)),imag(p(n)),'x');
         set(hp,'MarkerSize',12,'color',clr(plt));
         set(hp,'LineWidth',1.5);   
      end
   end
end
% Adjust the axes so |z|=1 is a circle.
axis([-lim,lim,-lim,lim]);
axis square
hold off;
% End of altered function
title(['L=',int2str(L(plt))])
end