``` 1 # Statistics. Computation of 95 % C.L. limits using stat. and systematical errors. 2 # This program demonstrates the computation of 95 % C.L. limits. 3 # Statistical errors are included in the treatment. 4 # 5 # Signal hypothesis is excluded at the 95% CL if CLs = 0.05 6 # and at more than the 95% CL if CLs < 0.05, assuming that signal is present 7 # 8 # Authors: Christophe.Delaere@cern.ch on 21.08.02 and Sergei Chekanov 9 10 from java.awt import Color,Font 11 from java.util import Random 12 from jhplot import * 13 from jhpro.stat import * 14 15 16 # just convinient print 17 def printCL(mess,confidence ): 18 print "\n------- "+mess+" -----------" 19 print "CLs : " ,confidence.getCLs() 20 print "CLb : " ,confidence.getCLb() 21 print "CLsb : " ,confidence.getCLsb() 22 print "expected : " ,confidence.getExpectedCLs_b() 23 print "expected : " ,confidence.getExpectedCLb_b() 24 print "expected : " ,confidence.getExpectedCLb_b() 25 print "Prob for background fluctuation", 1-confidence.getCLb() 26 27 28 c1 = HPlot("Canvas") 29 c1.setGTitle("Statistical significance of a signal") 30 c1.visible() 31 c1.setRange(-4,4,0.0,300) 32 c1.setNameX("Variable") 33 c1.setNameY("Events") 34 35 # set 36 background = H1D("Background",30,-4.0,4.0) 37 background.setColor(Color.green) 38 background.setFill(1) 39 background.setFillColor(Color.green) 40 background.setErrAll(0) 41 42 signal = H1D("Signal",30,-4.0,4.0) 43 signal.setFill(1) 44 signal.setFillColor(Color.red) 45 signal.setColor(Color.red) 46 47 data = H1D("Data",30,-4.0,4.0) 48 data.setColor(Color.black) 49 data.setStyle("p") 50 51 r=Random() 52 for i in range(1000): 53 background.fill(r.nextGaussian(),1.0) 54 data.fill(r.nextGaussian(),1.0) 55 56 for i in range(200): 57 signal.fill(1+0.2*r.nextGaussian(),1.0) 58 59 for i in range(70): 60 data.fill(1+0.2*r.nextGaussian(),1.0) 61 62 63 sigback=background.oper(signal,"Signal+Background","+") 64 sigback.setErrAll(0) 65 66 c1.cd(1,1) 67 c1.draw(signal) 68 c1.draw(background) 69 c1.draw(sigback) 70 c1.draw(data) 71 72 # add systematic errors from 2 sources 73 # err_names=["source1","source2"] 74 # err_b=[0.05, 0.0] # error on background is 5% (source1) and 0% (source2) 75 # err_s=[0.0, 0.01] # error on signal is 0% (source1) and 1% (source2) 76 77 print "Wait.. Calculating.." 78 # assume first that we do not have any systematics, only statistical uncertanties 79 d=DataSource() 80 d.addChannel(signal,background,data) 81 climit=CLimits(d, 100000) 82 printCL("No systematics", climit.getLimit()) 83 84 print "Wait.. Calculating.." 85 # assume data source has +/-10% systematics 86 d=DataSource() 87 d.addChannel(signal,background,data,0,0,0.1,"source") 88 d.addChannel(signal,background,data,0,0,-0.1,"source") 89 climit=CLimits(d, 100000) 90 printCL("Data with +/-5% systematics", climit.getLimit()) 91 92 print "Wait.. Calculating.." 93 # assume data source has +/-2% systematics and background has 10% systematics 94 d=DataSource() 95 d.addChannel(signal,background,data,0,0.1,0.02,"source") 96 d.addChannel(signal,background,data,0,-0.1,-0.02,"source") 97 climit=CLimits(d, 100000) 98 printCL("Data with +/-5% and background with +/-4% systematics", climit.getLimit()) 99 100 101 102 103 # export to some image (png,eps,pdf,jpeg...) 104 # c1.export(Editor.DocMasterName()+".png"); 105 # edit the image 106 # IEditor(Editor.DocMasterName()+".png"); ### © jHepWork. S.Chekanov ### ```