clear
%24 settembre 2019
%dati
durata=.02; %durata pacchetto in secondi
pacchetti=100000; %pacchetti totali simulati
step=10; %secondi; intervallo temporale verifica Poisson 

G=.05:.05:1; %carichi simulati
for j=1:length(G)
    lambda=G(j)/durata; %pacchetti/secondo
    interarrivo=-1/lambda*log(1-rand(1,pacchetti)); %esponenziale
    arrivo=cumsum(interarrivo); %temi di arrivo dei pacchetti
    intervals=step*(1:ceil(max(arrivo)/step));  %delimitatori intervalli
    %arrivi all'interno di cisacun intervallo
    arrivi=[];
    arrivi(1)=sum(arrivo<=intervals(1));
    for i=2:length(intervals)
        arrivi(i)=sum(arrivo<=intervals(i))-sum(arrivi(1:i-1));
    end
%verifica proprietà arrivi e interarrivi
%numero medio e varianza di arrivi per intervallo teorico: lambda*step
    risultati=sprintf('media e varianza teorica arrivi=%.2f, media=%.2f varianza=%.2f',lambda*step,mean(arrivi),var(arrivi));
    disp(risultati);
%tempo medio teorico di interarrivo 1/lambda, varianza 1/lambda^2
    risultati=sprintf('interarrivo media t.=%.1e, misura=%.1e, var. t.=%.1e, misura=%.1e',1/lambda,mean(interarrivo),1/lambda^2, var(interarrivo));
    disp(risultati);
%stima probabilità d'errore
    successi=0;
    if arrivo(1)+durata<arrivo(2)
        successi=1;
    end
    for i=2:pacchetti-1;
        if arrivo(i)>arrivo(i-1)+durata  && arrivo(i)+durata < arrivo(i+1)
            successi=successi+1;
        end
    end
    if arrivo(pacchetti)>arrivo(pacchetti-1)+durata
        successi=successi+1;
    end
    Ps(j)=successi/pacchetti;
%visualizza risultati
    risultati=sprintf('G=%.2f, tasso successo teorico = %.3f, stimato = %.3f S=%.3f',G(j),exp(-2*G(j)),Ps(j),G(j)*Ps(j));
    disp(risultati);
    input('premere un tasto');
end
%grafico
plot(G,G.*exp(-2*G),'k'); %teorico
hold on
plot(G,Ps.*G,'k*'); %simulato
xlabel('G');
ylabel('S');