r1=1/76.8*10^6;
r2=1/230.4*10^6;
T=3.84*10^-6;
lan=1/T;
N=5000;
tin=zeros(1,N);
mu1=1/2*(lan+r1+r2+sqrt((lan+r1+r2)^2-4*lan*r2));
mu2=1/2*(lan+r1+r2-sqrt((lan+r1+r2)^2-4*lan*r2));
p=(lan-mu2)/(mu1-mu2);

for i=1:N
    if rand<=p
        tin(i)=genindex(mu1,1);
    else
        tin(i)=genindex(mu2,1);
    end
end

M=700;
P1=zeros(1,M);
P2=zeros(1,M);

[P1,x]=hist(tin,M);
dt=x(2)-x(1);
P1=P1/length(tin);

for i=1:M
    P2(i)=hyperexpcdf(x(i)+dt/2,mu1,mu2,p)-hyperexpcdf(x(i)-dt/2,mu1,mu2,p);
end

plot(x,P1,'r');
hold;
plot(x,P2,'b');
hold;

K=100;
Tin=cumsum(tin);
figure(2)
h=plot(0,0,'.');
axis([0,Tin(K),0,K]);
grid off;
set(h,'EraseMode','none','MarkerSize',5);
figure(3)
h2=plot(0,1,'.');
axis([0,Tin(K),0,2]);
grid off;
set(h2,'EraseMode','none','MarkerSize',5);

k=1;
while(k<=K)
    t=timer('StartDelay',tin(k),'TimerFcn','disp(Tin(k))');
    start(t);
    wait(t);
    drawnow;
    set(h,'XData',Tin(k),'YData',k);
    set(h2,'XData',Tin(k),'YData',1);
    k=k+1;
end
delete(t);