# Build 2 experiments with observed data, and calculate
# probability for the given limit.  The user can  iterate towards 
# whatever limit type they choose, e.g., for a 90% upper 
# limit they adjust the limit guess until the probability is shown as 0.10; 
# for a 90% lower limit they would aim for 0.90.

from jhpro.stat.limit import *

sensitivity=1;
sigma=0;
limitguess=90;
MCevents=10000;
s=StatConfidence(sensitivity,sigma,limitguess,MCevents);

# build an experimental data
# - 200 observed events
# - 140 expected background
# - 20  error on expected background 
exp1=ExpData(200,140,20)
s.addData(exp1)

# second experiment with 150 observed events
exp2=ExpData(190,190,20)
s.addData(exp2)

print "-> Get probabilities for upper limits:"
s.runUpper()
print "Prob (Cousins+Highland)=",s.getProbabilityCH()   
print "Prob (BaBar SWG)=",s.getProbabilitySWG()
print "Prob (Jeffreys)=",s.getProbabilityJ()

# run for lower limit
print "\n->  Get probabilities for lower limits:"
s.runLower()
print "Prob (Cousins+Highland)=",s.getProbabilityCH()              
print "Prob (BaBar SWG)=",s.getProbabilitySWG()
print "Prob (Jeffreys)=",s.getProbabilityJ()