function [k,t] = erlbes(A,h,N,I,K0)
% Call
%Where A = offered traffic in Erlang (for M/x queue lambda = A/h)
%      h = holding time is seconds   (For M/x queye mu =1/h)
%      N = Number of resources
%      I = number of attempts/arrivals to process
%      K0 = warm up/cool down period
%      K = sequence of states visited

if (nargin == 4)
  K0 = 0; 
endif
I = I + K0
g = h/A; % Average Call Gap

% Initialise start times 
ts =  cumsum([0;g*(-log(rand(I-1,1)))]);
te =  ts + h*(-log(rand(I,1)));

[t,it] = sort([ts;te]); %it = original location

ic = ones(2*I,1);
ic(find(it>I)) =-1;


% To deal with blocked calls
% Find location of any time when state k > N
% Find corresponding call end  and put ic =0 

k=0;
for i=1:2*I
  k=k+ic(i);
  if (k>N)
    k=N;
    ic(i) =0;
    ic(find(it==(it(i)+I)))=0;
  endif
[i,ic(i),k];
endfor
k=cumsum(ic);
I=I-K0
k=k(K0+1:2*I+K0);
t=t(K0+1:2*I+K0);