% http://en.wikipedia.org/wiki/Exponential_distribution
%
% -ln(u)/lambda
%
function p = poissonMaxTime(lambda, maxTime)
p = cumsum(-log(rand(ceil(lambda*maxTime),1))/lambda);
t = p(end) - log(rand(1))/lambda;
while(t < maxTime)
p = [p; t];
t = t - log(rand(1))/lambda;
end
p = p(find(p < maxTime),1);