function [IMP,prumer_IMP,prumer_DO,prumer_CZ,prumer_CD,prumer_NLO,IP,p,N,Nq,Ns]=simulace_1(CPS, lambda,mi)
%funkce simulace_1 simuluje cinnost systemu M/M/1/Inf/FIFO
%syntaxe:
%[IMP,prumer_IMP,prumer_DO,prumer_CZ,prumer_CD,prumer_NLO,IP,p,N,Nq,Ns]=simulace_1(CPS, lambda,mi)
%vstupni parametry:
 %CPS ... celkovy pocet simulaci
 %lambda ... intenzita vstupniho proudu zakazniku
 %tj. pocet zakazniku za jednotku casu
 %mi ... intenzita obsluhy
 %tj. pocet obslouzenych zakazniku za jednotku casu
 %lambda musi byt mensi nez mi, jinak se system nestabilizuje
%vystupni parametry:
 %IMP ... delky intervalu mezi prichody zakazniku
 %prumer_IMP ... prumerna doba mezi prichody zakazniku
 %prumer_DO ... prumerna doba obsluhy
 %prumer_CZ ... prumerna doba cekani ve fronte
 %prumer_CD ... prumerna doba stravena v systemu
 %prumer_NLO ... prumerna doba nevyuziti linky obsluhy
 %IP ... intenzita provozu
 %p ... odhady pravdepodobnosti poctu zakazniku v systemu
 %N ... prumerny pocet zakazniku v systemu
%Nq ... prumerny pocet zakazniku ve fronte
%Ns ... prumerny pocet obsluhovanych zakazniku
if lambda/mi>=1 error('System se nemuze stabilizovat. Intenzita provozu je vetsi nez 1.') 
end
IMP=exprnd(1/lambda,CPS,1);
%IMP ... intervaly mezi prichody zakazniku
DO=exprnd(1/mi,CPS,1);
%DO ... doba obsluhy zakaznika
OPZ=cumsum(IMP);
%OPZ ... okamziky prichodu zakazniku
OOZ(1)=OPZ(1)+DO(1);
%OOZ ... okamziky odchodu zakazniku
%ZO ... zacatek obsluhy zakaznika
for i=2:CPS
if OPZ(i)<OOZ(i-1)
    ZO(i)=OOZ(i-1);
    OOZ(i)=ZO(i)+DO(i);
else
    ZO(i)=OPZ(i);
    OOZ(i)=ZO(i)+DO(i);
end
end
OOZ=OOZ';
CZ(1)=0;
%CZ ... cekani zakaznika na obsluhu
for i=2:CPS
    if ZO(i)>OPZ(i)
        CZ(i)=ZO(i)-OPZ(i);
    else
        CZ(i)=0;
    end
end
ZO=ZO';
CZ=CZ';
%PCZ ... celkovy pocet cekajicich zakazniku
PCZ=0;
for i=1:CPS
    if CZ(i)>0 PCZ=PCZ+1; end
end
%NLO ... nevyuziti linky obsluhy
NLO(1)=OPZ(1);
for i=2:CPS
    if OPZ(i)>OOZ(i-1)
        NLO(i)=ZO(i)-OOZ(i-1);
    else
        NLO(i)=0;
    end
end
NLO=NLO';
prumer_IMP=mean(IMP);
prumer_DO=mean(DO);
prumer_CZ=mean(CZ);
prumer_CD=prumer_CZ+prumer_DO;
prumer_NLO=mean(NLO);
%DATA=[IMP OPZ DO ZO OOZ CZ NLO];
cas=[OPZ;OOZ];
cas=sort(cas);
%cas ... serazene casy prichodu a odchodu zakazniku
%DCI ... delky casovych intervalu
DCI=diff(cas);
PZ(1)=1;
%PZ ... pocet zakazniku
for i=2:2*CPS 
    citac=0;
    for j=1:CPS
    if (OPZ(j)==cas(i)) citac=citac+1;end
    if (OOZ(j)==cas(i)) citac=citac-1;end 
    end
    PZ(i)=PZ(i-1)+citac;
end
PZ=PZ(1:2*CPS-1);
PZ=PZ';
[PZ,index]=sort(PZ);
DCI=DCI(index);
pocet=length(PZ);
%v ... varianty vektoru poctu zakazniku
%abscet ... abs. cetnosti jednotlivych variant
j=1;
v(j)=PZ(j);
pv(j)=1;
for i=1:(pocet-1)
  if (PZ(i)==PZ(i+1)) pv(j)=pv(j)+1;
  else j=j+1;
       pv(j)=1;
       v(j)=PZ(i+1);
  end
end
r=length(v);
v=v';
k=cumsum(pv);
%CDCI ... celkove delky casovych intervalu prislusne urcitemu poctu
%zakazniku 
CDCI(1)=cas(1)+sum(DCI(1:k(1)));
for i=2:r
    CDCI(i)=sum(DCI(k(i-1)+1:k(i)));
end
%p ... odhady pravdepodobnosti poctu zakazniku v systemu
p=CDCI/cas(length(cas));
%N ... prumerny pocet zakazniku v systemu
%Nq ... prumerny pocet zakazniku ve fronte
%Ns ... prumerny pocet obsluhovanych zakazniku
N=p*v;
Nq=p(2:r)*v(1:r-1);
Ns=N-Nq;
%IP ... intenzita provozu
IP=PCZ/CPS;
