function [uplim,pmu,ompb] = CLsUpperLimit(x,CL)
%CLsUpperLimit Compute CLs upper limit.
%   [uplim,pmu,ompb] = CLsUpperLimit(x,CL) finds the upper limit
%   at the CL confidence level using CLs=pmu/(1-pb) for the unit
%   gaussian measurement x with the null (background-only) hypothesis
%   at 0 and requiring the signal to be positive. The p-value for the
%   signal plus background hypothesis, pmu, and ompb, the p-value for the
%   null hypothesis, often called 1-p_b, are also returned. 
%
%   If x is a list of observations, lists of upper limits and
%   p-values are returned.

ompb = normcdf(x,0);
for i=1:length(x)
    uplim(i) = fzero(@(mu) normcdf(x(i),mu)/ompb(i)-(1-CL),x(i));
    pmu(i) = normcdf(x(i),uplim(i));
end

end