clear all; close all;

i=sqrt(-1);

theta=0:0.1:2*pi+0.1;

Z = rho(exp(i*theta))./sigma(exp(i*theta));
plot(Z)
axis([-3.5 0.5 -2 2])

tp = -1; 
tp = 0.5;

tpol1=[ 1 -1 0 0 ]; tpol2=[9 19 -5 1]/24;
cpol = tpol1 - tp*tpol2;
tz = roots(cpol);
max(abs(tz))