-
Notifications
You must be signed in to change notification settings - Fork 0
/
TestSPlot.C
154 lines (138 loc) · 5.89 KB
/
TestSPlot.C
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
#include "TSPlot.h"
#include "TTree.h"
#include "TH1.h"
#include "TCanvas.h"
#include "TFile.h"
#include "TPaveLabel.h"
#include "TPad.h"
#include "TPaveText.h"
#include "Riostream.h"
void TestSPlot()
{
//This tutorial illustrates the use of class TSPlot and of the sPlots method
//
//It is an example of analysis of charmless B decays, performed for the BABAR experiment.
//One is dealing with a data sample in which two species are present:
//the first is termed signal and the second background.
//A maximum Likelihood fit is performed to obtain the two yields N1 and N2
//The fit relies on two discriminating variables collectively denoted y,
//which are chosen within three possible variables denoted Mes, dE and F.
//The variable which is not incorporated in y, is used as the control variable x.
//The distributions of discriminating variables and more details about the method
//can be found in the TSPlot class description
//
// NOTE: This script requires a data file "TestSPlot_toyMC.dat".
// This data file is not distributed in the standard ROOT binary tar file.
// You can download it from ftp://root.cern.ch/root/TestSPlot_toyMC.dat
ifstream toymc;
toymc.open("TestSPlot_toyMC.dat");
if (!toymc.good()) {
printf("Before executing this script you must download the data file from\n");
printf(" ftp://root.cern.ch/root/TestSPlot_toyMC.dat\n");
return;
}
toymc.close();
//Read the data and initialize a TSPlot object
TTree *datatree = new TTree("datatree", "datatree");
datatree->ReadFile("TestSPlot_toyMC.dat", "Mes/D:dE/D:F/D:MesSignal/D:MesBackground/D:dESignal/D:dEBackground/D:FSignal/D:FBackground/D");
TSPlot *splot = new TSPlot(0, 3, 5420, 2, datatree);
//Set the selection for data tree
//Note the order of the variables: first the control variables (not presented in this example),
//then the 3 discriminating variables, then their probability distribution functions for
//the first species(signal) and then their pdfs for the second species(background)
splot->SetTreeSelection("Mes:dE:F:MesSignal:dESignal:FSignal:MesBackground:dEBackground:FBackground");
//Set the initial estimates of the number of events in each species - used as initial
//parameter values for the Minuit likelihood fit
Int_t ne[2];
ne[0]=500; ne[1]=5000;
splot->SetInitialNumbersOfSpecies(ne);
//Compute the weights
splot->MakeSPlot();
//Fill the sPlots
splot->FillSWeightsHists(25);
//Now let's look at the sPlots
//The first two histograms are sPlots for the Mes variable signal and background.
//dE and F were chosen as discriminating variables to determine N1 and N2, through a
//maximum Likelihood fit, and thus the sPlots for the control variable Mes, unknown
//to the fit, was contructed.
//One can see that the sPlot for signal reproduces the PDF correctly,
//even when the latter vanishes.
//
//The lower two histograms are sPlots for the F variables signal and background.
//dE and Mes were chosen as discriminating variables to determine N1 and N2, through a
//maximum Likelihood fit, and thus the sPlots for the control variable F, unknown
//to the fit, was contructed.
TCanvas *myc = new TCanvas("myc", "sPlots of Mes and F signal and background", 800, 600);
myc->SetFillColor(40);
TPaveText *pt = new TPaveText(0.02,0.85,0.98,0.98);
pt->SetFillColor(18);
pt->SetTextFont(20);
pt->SetTextColor(4);
pt->AddText("sPlots of Mes and F signal and background,");
pt->AddText("obtained by running the tutorial TestSPlot.C on BABAR Monte Carlo data (sPlot_toyMC.fit)");
TText *t3=pt->AddText("M. Pivk and F. R. Le Diberder, Nucl.Inst.Meth.A (in press), physics/0402083");
t3->SetTextColor(1);
t3->SetTextFont(30);
pt->Draw();
TPad* pad1 = new TPad("pad1","Mes signal",0.02,0.43,0.48,0.83,33);
TPad* pad2 = new TPad("pad2","Mes background",0.5,0.43,0.98,0.83,33);
TPad* pad3 = new TPad("pad3", "F signal", 0.02, 0.02, 0.48, 0.41,33);
TPad* pad4 = new TPad("pad4", "F background", 0.5, 0.02, 0.98, 0.41,33);
pad1->Draw();
pad2->Draw();
pad3->Draw();
pad4->Draw();
pad1->cd();
pad1->SetGrid();
TH1D *sweight00 = splot->GetSWeightsHist(-1, 0, 0);
sweight00->SetTitle("Mes signal");
sweight00->SetStats(kFALSE);
sweight00->Draw("e");
sweight00->SetMarkerStyle(21);
sweight00->SetMarkerSize(0.7);
sweight00->SetMarkerColor(2);
sweight00->SetLineColor(2);
sweight00->GetXaxis()->SetLabelSize(0.05);
sweight00->GetYaxis()->SetLabelSize(0.06);
sweight00->GetXaxis()->SetLabelOffset(0.02);
pad2->cd();
pad2->SetGrid();
TH1D *sweight10 = splot->GetSWeightsHist(-1, 1, 0);
sweight10->SetTitle("Mes background");
sweight10->SetStats(kFALSE);
sweight10->Draw("e");
sweight10->SetMarkerStyle(21);
sweight10->SetMarkerSize(0.7);
sweight10->SetMarkerColor(2);
sweight10->SetLineColor(2);
sweight10->GetXaxis()->SetLabelSize(0.05);
sweight10->GetYaxis()->SetLabelSize(0.06);
sweight10->GetXaxis()->SetLabelOffset(0.02);
pad3->cd();
pad3->SetGrid();
TH1D *sweight02 = splot->GetSWeightsHist(-1, 0, 2);
sweight02->SetTitle("F signal");
sweight02->SetStats(kFALSE);
sweight02->Draw("e");
sweight02->SetMarkerStyle(21);
sweight02->SetMarkerSize(0.7);
sweight02->SetMarkerColor(2);
sweight02->SetLineColor(2);
sweight02->GetXaxis()->SetLabelSize(0.06);
sweight02->GetYaxis()->SetLabelSize(0.06);
sweight02->GetXaxis()->SetLabelOffset(0.01);
pad4->cd();
pad4->SetGrid();
TH1D *sweight12 = splot->GetSWeightsHist(-1, 1, 2);
sweight12->SetTitle("F background");
sweight12->SetStats(kFALSE);
sweight12->Draw("e");
sweight12->SetMarkerStyle(21);
sweight12->SetMarkerSize(0.7);
sweight12->SetMarkerColor(2);
sweight12->SetLineColor(2);
sweight12->GetXaxis()->SetLabelSize(0.06);
sweight12->GetYaxis()->SetLabelSize(0.06);
sweight02->GetXaxis()->SetLabelOffset(0.01);
myc->cd();
}