clc;
clearvars; close all;

%% Compute

lambda = 30;
alpha = 0.25;
sum = 0;
confidence = 1 - alpha;
k = 0; 

while (sum < confidence)
    sum = sum + ((lambda^k)/(factorial(k))) * exp(-lambda);
    k = k+1;
end;

fprintf('n = %d, q_x(%f) = %f\n', k-1, confidence, sum);

%% Show Poisson distribution pdf

ssum = 0;
kk = 0;
px = [];

while (ssum < 0.9999)
    px(end+1) = ((lambda^kk)/(factorial(kk))) * exp(-lambda);
    ssum = ssum + ((lambda^kk)/(factorial(kk))) * exp(-lambda);
    kk = kk+1;
end;

figure;
stem(px); hold on;
set(gca, 'XTickLabel', -1:size(px, 2)+1);
set(gca, 'XTick', 0:size(px, 2)+1);
xlim([1 length(px)]);
stem(px(1:k), 'color', 'r');
%vline(k, 'r', sprintf('%d', k-1));
