function  pz_plot(p,z,figno,dt)
% pz_plot(p,z,figno,dt)
% plot the poles and zeros of a MIMO system in the Laplace domain
% in figure number figno
% if dt is provided the poles and zeros are interpreted as 
% discrete-time poles and zeros and are plotted with respect to the unit circle


    if nargin == 4, dt=1; else dt=0; end

    p = p(:);
    z = z(:);

    ff = 1.20;
    if dt
	ax = [ -1 , 1 , -1 , 1 ]*1.1;
    else
	xr = (max([real(p);real(z)]) - min([real(p);real(z)]))/2; 
	yr = (max([imag(p);imag(z)]) - min([imag(p);imag(z)]))/2; 
	xa = (max([real(p);real(z)]) + min([real(p);real(z)]))/2; 
	ya = (max([imag(p);imag(z)]) + min([imag(p);imag(z)]))/2; 
	ax = [ xa-ff*xr , max(xa+ff*xr,0.1*xr) , ya-ff*yr , ya+ff*yr ];
    end

   figure(figno);
     clf
     hold on
     plot(real(p), imag(p), 'xr', 'MarkerSize',5);
     plot(real(z), imag(z), 'ob', 'MarkerSize',5);
     axis(ax);
     axis('square')
     plot([ax(1), ax(2)],[0,0],'-k', [0,0],[ax(3),ax(4)],'-k')
     if dt,
	plot(cosd(0:360),sind(0:360),'-k'); % unit circle
     end
     set(gca,'box','off');
     set(gca,'XAxisLocation','origin');
     set(gca,'YAxisLocation','origin');
     set(gca,'LineWidth',2);

     if dt
	text(0.7,-0.1,'Re(\lambda_i)')
	text(0.1,0.7,'Im(\lambda_i)','Rotation',90)
     else
	text(0.8*ax(1),0.1*ax(4),'Re(\lambda_i) = \zeta_i \omega_{ni}')
	text(0.5*ax(2),0.5*ax(4),'Im(\lambda_i) = \omega_{di}','Rotation',90)
     end
     hold off

% --------------------------------------------------------------------------
% H.P. Gavin, Duke Univ.,  2017-10-07
