Loading [MathJax]/jax/output/HTML-CSS/config.js
Logo ROOT  
Reference Guide
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Modules Pages
limit.C File Reference

Detailed Description

View in nbviewer Open in SWAN This program demonstrates the computation of 95 % C.L.

limits. It uses a set of randomly created histograms.

Computing limits...
CLs : 0.0179535
CLsb : 0.00930313
CLb : 0.51818
< CLs > : 0.0165344
< CLsb > : 0.00826754
< CLb > : 0.50002
Computing limits with stat systematics...
CLs : 0.0189723
CLsb : 0.00989482
CLb : 0.52154
< CLs > : 0.0172591
< CLsb > : 0.00862989
< CLb > : 0.50002
Computing limits with systematics...
CLs : 0.0352037
CLsb : 0.0184327
CLb : 0.5236
< CLs > : 0.0328702
< CLsb > : 0.0164358
< CLb > : 0.50002
#include <iostream>
#include "TH1.h"
#include "THStack.h"
#include "TCanvas.h"
#include "TFrame.h"
#include "TRandom2.h"
#include "TSystem.h"
#include "TVector.h"
#include "TObjArray.h"
#include "TLimit.h"
using std::cout;
using std::endl;
void limit() {
// Create a new canvas.
TCanvas *c1 = new TCanvas("c1","Dynamic Filling Example",200,10,700,500);
c1->SetFillColor(42);
// Create some histograms
TH1D* backgroundHist = new TH1D("background","The expected background",30,-4,4);
TH1D* signalHist = new TH1D("signal","the expected signal",30,-4,4);
TH1D* dataHist = new TH1D("data","some fake data points",30,-4,4);
backgroundHist->SetFillColor(48);
signalHist->SetFillColor(41);
dataHist->SetMarkerStyle(21);
dataHist->SetMarkerColor(kBlue);
backgroundHist->Sumw2(); // needed for stat uncertainty
signalHist->Sumw2(); // needed for stat uncertainty
// Fill histograms randomly
Float_t bg,sig,dt;
for (Int_t i = 0; i < 25000; i++) {
bg = r.Gaus(0,1);
sig = r.Gaus(1,.2);
backgroundHist->Fill(bg,0.02);
signalHist->Fill(sig,0.001);
}
for (Int_t i = 0; i < 500; i++) {
dt = r.Gaus(0,1);
dataHist->Fill(dt);
}
THStack *hs = new THStack("hs","Signal and background compared to data...");
hs->Add(backgroundHist);
hs->Add(signalHist);
hs->Draw("hist");
dataHist->Draw("PE1,Same");
c1->Modified();
c1->Update();
c1->GetFrame()->SetFillColor(21);
c1->GetFrame()->SetBorderSize(6);
c1->GetFrame()->SetBorderMode(-1);
c1->Modified();
c1->Update();
// Compute the limits
cout << "Computing limits... " << endl;
TLimitDataSource* mydatasource = new TLimitDataSource(signalHist,backgroundHist,dataHist);
TConfidenceLevel *myconfidence = TLimit::ComputeLimit(mydatasource,50000);
cout << "CLs : " << myconfidence->CLs() << endl;
cout << "CLsb : " << myconfidence->CLsb() << endl;
cout << "CLb : " << myconfidence->CLb() << endl;
cout << "< CLs > : " << myconfidence->GetExpectedCLs_b() << endl;
cout << "< CLsb > : " << myconfidence->GetExpectedCLsb_b() << endl;
cout << "< CLb > : " << myconfidence->GetExpectedCLb_b() << endl;
// Add stat uncertainty
cout << endl << "Computing limits with stat systematics... " << endl;
TConfidenceLevel *mystatconfidence = TLimit::ComputeLimit(mydatasource,50000,true);
cout << "CLs : " << mystatconfidence->CLs() << endl;
cout << "CLsb : " << mystatconfidence->CLsb() << endl;
cout << "CLb : " << mystatconfidence->CLb() << endl;
cout << "< CLs > : " << mystatconfidence->GetExpectedCLs_b() << endl;
cout << "< CLsb > : " << mystatconfidence->GetExpectedCLsb_b() << endl;
cout << "< CLb > : " << mystatconfidence->GetExpectedCLb_b() << endl;
// Add some systematics
cout << endl << "Computing limits with systematics... " << endl;
TVectorD errorb(2);
TVectorD errors(2);
TObjArray* names = new TObjArray();
TObjString name1("bg uncertainty");
TObjString name2("sig uncertainty");
names->AddLast(&name1);
names->AddLast(&name2);
errorb[0]=0.05; // error source 1: 5%
errorb[1]=0; // error source 2: 0%
errors[0]=0; // error source 1: 0%
errors[1]=0.01; // error source 2: 1%
TLimitDataSource* mynewdatasource = new TLimitDataSource();
mynewdatasource->AddChannel(signalHist,backgroundHist,dataHist,&errors,&errorb,names);
TConfidenceLevel *mynewconfidence = TLimit::ComputeLimit(mynewdatasource,50000,true);
cout << "CLs : " << mynewconfidence->CLs() << endl;
cout << "CLsb : " << mynewconfidence->CLsb() << endl;
cout << "CLb : " << mynewconfidence->CLb() << endl;
cout << "< CLs > : " << mynewconfidence->GetExpectedCLs_b() << endl;
cout << "< CLsb > : " << mynewconfidence->GetExpectedCLsb_b() << endl;
cout << "< CLb > : " << mynewconfidence->GetExpectedCLb_b() << endl;
// show canonical -2lnQ plots in a new canvas
// - The histogram of -2lnQ for background hypothesis (full)
// - The histogram of -2lnQ for signal and background hypothesis (dashed)
TCanvas *c2 = new TCanvas("c2");
myconfidence->Draw();
// clean up (except histograms and canvas)
delete myconfidence;
delete mydatasource;
delete mystatconfidence;
delete mynewconfidence;
delete mynewdatasource;
}
ROOT::R::TRInterface & r
Definition: Object.C:4
int Int_t
Definition: RtypesCore.h:43
float Float_t
Definition: RtypesCore.h:55
@ kBlue
Definition: Rtypes.h:64
R__EXTERN TSystem * gSystem
Definition: TSystem.h:556
virtual void SetFillColor(Color_t fcolor)
Set the fill area color.
Definition: TAttFill.h:37
virtual void SetMarkerColor(Color_t mcolor=1)
Set the marker color.
Definition: TAttMarker.h:38
virtual void SetMarkerStyle(Style_t mstyle=1)
Set the marker style.
Definition: TAttMarker.h:40
The Canvas class.
Definition: TCanvas.h:27
Class to compute 95% CL limits.
Double_t GetExpectedCLs_b(Int_t sigma=0) const
Double_t GetExpectedCLb_b(Int_t sigma=0) const
Get the expected Confidence Level for the background only if there is only background.
Double_t CLsb(bool use_sMC=kFALSE) const
Get the Confidence Level for the signal plus background hypothesis.
Double_t CLb(bool use_sMC=kFALSE) const
Get the Confidence Level for the background only.
Double_t GetExpectedCLsb_b(Int_t sigma=0) const
Get the expected Confidence Level for the signal plus background hypothesis if there is only backgrou...
Double_t CLs(bool use_sMC=kFALSE) const
Get the Confidence Level defined by CLs = CLsb/CLb.
void Draw(const Option_t *option="")
Display sort of a "canonical" -2lnQ plot.
1-D histogram with a double per channel (see TH1 documentation)}
Definition: TH1.h:614
virtual Int_t Fill(Double_t x)
Increment bin with abscissa X by 1.
Definition: TH1.cxx:3275
virtual void Draw(Option_t *option="")
Draw this histogram with options.
Definition: TH1.cxx:2998
virtual void Sumw2(Bool_t flag=kTRUE)
Create structure to store sum of squares of weights.
Definition: TH1.cxx:8476
The Histogram stack class.
Definition: THStack.h:31
virtual void Draw(Option_t *chopt="")
Draw this multihist with its current attributes.
Definition: THStack.cxx:445
virtual void Add(TH1 *h, Option_t *option="")
add a new histogram to the list Only 1-d and 2-d histograms currently supported.
Definition: THStack.cxx:359
This class serves as interface to feed data into the TLimit routines.
virtual void AddChannel(TH1 *, TH1 *, TH1 *)
static TConfidenceLevel * ComputeLimit(TLimitDataSource *data, Int_t nmc=50000, bool stat=false, TRandom *generator=0)
Definition: TLimit.cxx:103
An array of TObjects.
Definition: TObjArray.h:37
virtual void AddLast(TObject *obj)
Add object in the next empty slot in the array.
Definition: TObjArray.cxx:178
Collectable string class.
Definition: TObjString.h:28
Random number generator class based on the maximally quidistributed combined Tausworthe generator by ...
Definition: TRandom2.h:27
virtual Bool_t ProcessEvents()
Process pending events (GUI, timers, sockets).
Definition: TSystem.cxx:414
return c1
Definition: legend1.C:41
return c2
Definition: legend2.C:14
Author
Christophe Delaere

Definition in file limit.C.