function [y,index] = markovsim(z,Pi,T)


t=1;
ns = length(z);
ssz = ns/2 + .5;

zyes = ssz;

y(t) = z(zyes);
index(t) = ssz;
   while (t+1<=T)
       u = rand(1);
       ztoday = 0;
       sumz=0;
       while (u>sumz && ztoday<=ns)
           sumz=sumz+Pi(zyes,ztoday+1);
           ztoday=ztoday+1;
       end;
       
       y(t+1) = z(ztoday);
       index(t+1) = ztoday;
       zyes = ztoday;
       t=t+1;
   end;