1 |
//-----------------------------------------------
|
2 |
// Takes the text files and dumps the contents
|
3 |
// in an ntuple.
|
4 |
// It's OK if a text file does not exist.
|
5 |
// However, the text files need to have the
|
6 |
// same order of gluino-mass LSP-mass etc...
|
7 |
// Meaning that the lines must correspond
|
8 |
// If they do not, the program will barf
|
9 |
//
|
10 |
// It expects that there are files for 7 signal regions.
|
11 |
//
|
12 |
// The txt files are supposed to be "sorted".
|
13 |
// If they are not, you may want to sort them
|
14 |
// ahead of time, eg:
|
15 |
//
|
16 |
// sort t1tttt_sr1.txt > SR1_sorted.txt
|
17 |
// sort t1tttt_sr2.txt > SR2_sorted.txt
|
18 |
// sort t1tttt_sr3.txt > SR3_sorted.txt
|
19 |
// sort t1tttt_sr4.txt > SR4_sorted.txt
|
20 |
// sort t1tttt_sr5.txt > SR5_sorted.txt
|
21 |
// sort t1tttt_sr6.txt > SR6_sorted.txt
|
22 |
// sort t1tttt_sr7.txt > SR7_sorted.txt
|
23 |
//
|
24 |
// Make sure that the file names, which are
|
25 |
// hardwired (sorry) are correct.
|
26 |
//
|
27 |
// Also reads in the gluino xsection txt file.
|
28 |
// Assumes that the gluino mass increases monotonically
|
29 |
// from one line to the next
|
30 |
//
|
31 |
//
|
32 |
// Usage:
|
33 |
// root> .L makeNtuple()
|
34 |
// root> makeNtuple()
|
35 |
//-------------------------------------------
|
36 |
|
37 |
//
|
38 |
//
|
39 |
|
40 |
void makeNtuple(){
|
41 |
|
42 |
//-----------------------------------------------
|
43 |
// The name of the file with the output ntuple
|
44 |
//-----------------------------------------------
|
45 |
char* ntFile = "ntuple.root";
|
46 |
|
47 |
|
48 |
//-----------------------------------------------
|
49 |
// Raed in the gluino pair production xsection
|
50 |
//-----------------------------------------------
|
51 |
float xsec_mass[2000];
|
52 |
float xsec[2000];
|
53 |
float xsecup[2000];
|
54 |
float xsecdwn[2000];
|
55 |
int nxsec=0;
|
56 |
ifstream fromfile("gluino_xsec_7Tev.txt");
|
57 |
if ( !fromfile.good() ) {
|
58 |
cout << "gluino xsec file does not exist" << endl ;
|
59 |
return;
|
60 |
}
|
61 |
cout << "Reading gluino xsection" << endl;
|
62 |
int xs=0;
|
63 |
while (!fromfile.eof()) {
|
64 |
TString d1, d2, d3, d4;
|
65 |
fromfile >> d1 >> xsec_mass[xs] >> d2 >> xsec[xs] >> d3
|
66 |
>> xsecup[xs] >> d4 >> xsecdwn[xs];
|
67 |
xs++;
|
68 |
nxsec++;
|
69 |
}
|
70 |
fromfile.close();
|
71 |
nxsec = nxsec-1; //last line is junk
|
72 |
|
73 |
for (int i=0; i<nxsec; i++) {
|
74 |
cout<<xsec_mass[i]<<" "<<xsec[i]<<" "<<xsecup[i]<<" "<< xsecdwn[i]<< endl;
|
75 |
}
|
76 |
|
77 |
//--------------------------------------------------------
|
78 |
// Check that the gluino mass in the xsection file increases monotonically
|
79 |
//---------------------------------------------------------
|
80 |
float previous = xsec_mass[0];
|
81 |
for (int i=1; i<nxsec; i++) {
|
82 |
if (xsec_mass[i] < previous) {
|
83 |
cout << "Mass in xsection file does not increase monotonically..Abort"<<endl;
|
84 |
return;
|
85 |
}
|
86 |
previous = xsec_mass[i];
|
87 |
}
|
88 |
|
89 |
//-----------------------------------------------
|
90 |
// Define the variables in the txt files
|
91 |
//-----------------------------------------------
|
92 |
float glmass[7][1000];
|
93 |
float lspmass[7][1000];
|
94 |
float eff[7][1000];
|
95 |
float err[7][1000];
|
96 |
float obslim[7][1000];
|
97 |
float explim[7][1000];
|
98 |
bool usedflag[7][1000]; // a flag to be used later
|
99 |
int nentries[7]; // number of entries in the file
|
100 |
|
101 |
//-----------------------------------------------
|
102 |
// The file names, ordered by signal region...
|
103 |
//-----------------------------------------------
|
104 |
vector<TString> filenames;
|
105 |
filenames.push_back("t1tttt_sr1.txt");
|
106 |
filenames.push_back("t1tttt_sr2.txt");
|
107 |
filenames.push_back("t1tttt_sr3.txt");
|
108 |
filenames.push_back("t1tttt_sr4.txt");
|
109 |
filenames.push_back("t1tttt_sr5.txt");
|
110 |
filenames.push_back("t1tttt_sr6.txt");
|
111 |
filenames.push_back("t1tttt_sr7.txt");
|
112 |
|
113 |
//-----------------------------------------------
|
114 |
// Loop over signal regions and load the arrays from the txt files
|
115 |
//-----------------------------------------------
|
116 |
for (int ireg=0; ireg<7; ireg++) {
|
117 |
nentries[ireg] = 0;
|
118 |
ifstream fromfile(filenames.at(ireg));
|
119 |
|
120 |
// if the file does not exist, go to the next one
|
121 |
if ( !fromfile.good() ) {
|
122 |
cout << "file " << filenames.at(ireg) << " does not exist" << endl ;
|
123 |
continue;
|
124 |
}
|
125 |
// the file exist, read from it
|
126 |
cout << "Reading from " << filenames.at(ireg) << endl;
|
127 |
int j=0;
|
128 |
while (!fromfile.eof()) {
|
129 |
fromfile >> glmass[ireg][j]
|
130 |
>> lspmass[ireg][j]
|
131 |
>> eff[ireg][j]
|
132 |
>> err[ireg][j]
|
133 |
>> obslim[ireg][j]
|
134 |
>> explim[ireg][j];
|
135 |
j++;
|
136 |
nentries[ireg]++;
|
137 |
}
|
138 |
fromfile.close();
|
139 |
// the last entry is junk
|
140 |
nentries[ireg] = nentries[ireg]-1;
|
141 |
}
|
142 |
//-------------------------------------------
|
143 |
// Now check that the line order is OK
|
144 |
// And that the number of lines is consistent
|
145 |
//--------------------------------------------
|
146 |
bool bad=false;
|
147 |
int numbLines = 0;
|
148 |
for (int ii=0; ii<7; ii++) {
|
149 |
for (int jj=ii+1; jj<7; jj++) {
|
150 |
if (nentries[ii]*nentries[jj] !=0 ) {
|
151 |
numbLines = nentries[ii];
|
152 |
if (nentries[ii] != nentries[jj]) {
|
153 |
cout << "Mismatched lines: file " << filenames.at(ii);
|
154 |
cout << " has " << nentries[ii] << " lines" << endl;
|
155 |
cout << " file " << filenames.at(jj);
|
156 |
cout << " has " << nentries[jj] << " lines" << endl;
|
157 |
bad=true;
|
158 |
}
|
159 |
for (int line=0; line<nentries[ii]; line++) {
|
160 |
if (glmass[ii][line] != glmass[jj][line]) {
|
161 |
cout << "Gluino mass on line " << line
|
162 |
<< " of files " << filenames.at(ii) << " and "
|
163 |
<< filenames.at(jj) << " do not match" << endl;
|
164 |
bad=true;
|
165 |
}
|
166 |
if (lspmass[ii][line] != lspmass[jj][line]) {
|
167 |
cout << "LSP mass on line " << line
|
168 |
<< " of files " << filenames.at(ii) << " and "
|
169 |
<< filenames.at(jj) << " do not match" << endl;
|
170 |
bad=true;
|
171 |
}
|
172 |
}
|
173 |
}
|
174 |
}
|
175 |
}
|
176 |
if (bad) return;
|
177 |
cout << " The file formats match...Now fill ntuple" <<endl;
|
178 |
|
179 |
//-----------------------------------------------------
|
180 |
// Now the content of the files is loaded in arrays
|
181 |
// We are going to stick everything in an ntuple
|
182 |
//-----------------------------------------------------
|
183 |
TFile* babyFile_ = TFile::Open(Form("%s", ntFile), "RECREATE");
|
184 |
babyFile_->cd();
|
185 |
babyTree_ = new TTree("tree", "A Baby Ntuple");
|
186 |
babyTree_->SetDirectory(0);
|
187 |
|
188 |
// define the variables
|
189 |
float glmass_; // gluino mass
|
190 |
float lspmass_; // lsp mass
|
191 |
float xsec_; // xsec
|
192 |
float xsecup_; // xsec, one sigma high
|
193 |
float xsecdwn_; // xsec, one sigma down
|
194 |
int bestsr_ = 0; // best signal region
|
195 |
|
196 |
// Here come the efficiencies for the various SR.
|
197 |
// The last one is the "best", based on the best expected limit
|
198 |
float effsr1_ = -1.;
|
199 |
float effsr2_ = -1.;
|
200 |
float effsr3_ = -1.;
|
201 |
float effsr4_ = -1.;
|
202 |
float effsr5_ = -1.;
|
203 |
float effsr6_ = -1.;
|
204 |
float effsr7_ = -1.;
|
205 |
float effsrb_ = -1.; // b = best SR
|
206 |
|
207 |
// same as above but for the efficiency error
|
208 |
float errsr1_ = -1.;
|
209 |
float errsr2_ = -1.;
|
210 |
float errsr3_ = -1.;
|
211 |
float errsr4_ = -1.;
|
212 |
float errsr5_ = -1.;
|
213 |
float errsr6_ = -1.;
|
214 |
float errsr7_ = -1.;
|
215 |
float errsrb_ = -1.; // b = best SR
|
216 |
|
217 |
// same as above but for observed limit
|
218 |
float obslimsr1_ = -1.;
|
219 |
float obslimsr2_ = -1.;
|
220 |
float obslimsr3_ = -1.;
|
221 |
float obslimsr4_ = -1.;
|
222 |
float obslimsr5_ = -1.;
|
223 |
float obslimsr6_ = -1.;
|
224 |
float obslimsr7_ = -1.;
|
225 |
float obslimsrb_ = -1.; // b = best SR
|
226 |
|
227 |
// same as above but for expected limit
|
228 |
float explimsr1_ = -1.;
|
229 |
float explimsr2_ = -1.;
|
230 |
float explimsr3_ = -1.;
|
231 |
float explimsr4_ = -1.;
|
232 |
float explimsr5_ = -1.;
|
233 |
float explimsr6_ = -1.;
|
234 |
float explimsr7_ = -1.;
|
235 |
float explimsrb_ = -1.; // b = best SR
|
236 |
|
237 |
|
238 |
// put the branches in the tree
|
239 |
babyTree_->Branch("glmass", &glmass_);
|
240 |
babyTree_->Branch("lspmass", &lspmass_);
|
241 |
babyTree_->Branch("xsec", &xsec_);
|
242 |
babyTree_->Branch("xsecup", &xsecup_);
|
243 |
babyTree_->Branch("xsecdwn", &xsecdwn_);
|
244 |
babyTree_->Branch("bestsr", &bestsr_);
|
245 |
|
246 |
|
247 |
babyTree_->Branch("effsr1", &effsr1_);
|
248 |
babyTree_->Branch("effsr2", &effsr2_);
|
249 |
babyTree_->Branch("effsr3", &effsr3_);
|
250 |
babyTree_->Branch("effsr4", &effsr4_);
|
251 |
babyTree_->Branch("effsr5", &effsr5_);
|
252 |
babyTree_->Branch("effsr6", &effsr6_);
|
253 |
babyTree_->Branch("effsr7", &effsr7_);
|
254 |
babyTree_->Branch("effsrb", &effsrb_); // b = best SR
|
255 |
|
256 |
babyTree_->Branch("errsr1", &errsr1_);
|
257 |
babyTree_->Branch("errsr2", &errsr2_);
|
258 |
babyTree_->Branch("errsr3", &errsr3_);
|
259 |
babyTree_->Branch("errsr4", &errsr4_);
|
260 |
babyTree_->Branch("errsr5", &errsr5_);
|
261 |
babyTree_->Branch("errsr6", &errsr6_);
|
262 |
babyTree_->Branch("errsr7", &errsr7_);
|
263 |
babyTree_->Branch("errsrb", &errsrb_); // b = best SR
|
264 |
|
265 |
babyTree_->Branch("obslimsr1", &obslimsr1_);
|
266 |
babyTree_->Branch("obslimsr2", &obslimsr2_);
|
267 |
babyTree_->Branch("obslimsr3", &obslimsr3_);
|
268 |
babyTree_->Branch("obslimsr4", &obslimsr4_);
|
269 |
babyTree_->Branch("obslimsr5", &obslimsr5_);
|
270 |
babyTree_->Branch("obslimsr6", &obslimsr6_);
|
271 |
babyTree_->Branch("obslimsr7", &obslimsr7_);
|
272 |
babyTree_->Branch("obslimsrb", &obslimsrb_); // b = best SR
|
273 |
|
274 |
babyTree_->Branch("explimsr1", &explimsr1_);
|
275 |
babyTree_->Branch("explimsr2", &explimsr2_);
|
276 |
babyTree_->Branch("explimsr3", &explimsr3_);
|
277 |
babyTree_->Branch("explimsr4", &explimsr4_);
|
278 |
babyTree_->Branch("explimsr5", &explimsr5_);
|
279 |
babyTree_->Branch("explimsr6", &explimsr6_);
|
280 |
babyTree_->Branch("explimsr7", &explimsr7_);
|
281 |
babyTree_->Branch("explimsrb", &explimsrb_); // b = best SR
|
282 |
|
283 |
// fill the variables.
|
284 |
// Note: the best region is the one with the smallest
|
285 |
// ratio of expected limit and efficiency
|
286 |
cout << "Filling an ntuple with " << numbLines << " entries" << endl;
|
287 |
for (int il=0; il<numbLines; il++) {
|
288 |
float best = 999999.;
|
289 |
int ibest = -1;
|
290 |
if (nentries[0] != 0) {
|
291 |
glmass_ = glmass[0][il];
|
292 |
lspmass_ = lspmass[0][il];
|
293 |
effsr1_ = eff[0][il];
|
294 |
errsr1_ = eff[0][il];
|
295 |
obslimsr1_ = obslim[0][il];
|
296 |
explimsr1_ = explim[0][il];
|
297 |
float temp = explimsr1_/effsr1_;
|
298 |
if (temp < best) {
|
299 |
best = temp;
|
300 |
ibest = 1;
|
301 |
}
|
302 |
}
|
303 |
if (nentries[1] != 0) {
|
304 |
glmass_ = glmass[1][il];
|
305 |
lspmass_ = lspmass[1][il];
|
306 |
effsr2_ = eff[1][il];
|
307 |
errsr2_ = eff[1][il];
|
308 |
obslimsr2_ = obslim[1][il];
|
309 |
explimsr2_ = explim[1][il];
|
310 |
float temp = explimsr2_/effsr2_;
|
311 |
if (temp < best) {
|
312 |
best = temp;
|
313 |
ibest = 2;
|
314 |
}
|
315 |
}
|
316 |
if (nentries[2] != 0) {
|
317 |
glmass_ = glmass[2][il];
|
318 |
lspmass_ = lspmass[2][il];
|
319 |
effsr3_ = eff[2][il];
|
320 |
errsr3_ = eff[2][il];
|
321 |
obslimsr3_ = obslim[2][il];
|
322 |
explimsr3_ = explim[2][il];
|
323 |
float temp = explimsr3_/effsr3_;
|
324 |
if (temp < best) {
|
325 |
best = temp;
|
326 |
ibest = 3;
|
327 |
}
|
328 |
}
|
329 |
if (nentries[3] != 0) {
|
330 |
glmass_ = glmass[3][il];
|
331 |
lspmass_ = lspmass[3][il];
|
332 |
effsr4_ = eff[3][il];
|
333 |
errsr4_ = eff[3][il];
|
334 |
obslimsr4_ = obslim[3][il];
|
335 |
explimsr4_ = explim[3][il];
|
336 |
float temp = explimsr4_/effsr4_;
|
337 |
if (temp < best) {
|
338 |
best = temp;
|
339 |
ibest = 4;
|
340 |
}
|
341 |
}
|
342 |
if (nentries[4] != 0) {
|
343 |
glmass_ = glmass[4][il];
|
344 |
lspmass_ = lspmass[4][il];
|
345 |
effsr5_ = eff[4][il];
|
346 |
errsr5_ = eff[4][il];
|
347 |
obslimsr5_ = obslim[4][il];
|
348 |
explimsr5_ = explim[4][il];
|
349 |
float temp = explimsr5_/effsr5_;
|
350 |
if (temp < best) {
|
351 |
best = temp;
|
352 |
ibest = 5;
|
353 |
}
|
354 |
}
|
355 |
if (nentries[5] != 0) {
|
356 |
glmass_ = glmass[5][il];
|
357 |
lspmass_ = lspmass[5][il];
|
358 |
effsr6_ = eff[5][il];
|
359 |
errsr6_ = eff[5][il];
|
360 |
obslimsr6_ = obslim[5][il];
|
361 |
explimsr6_ = explim[5][il];
|
362 |
float temp = explimsr6_/effsr6_;
|
363 |
if (temp < best) {
|
364 |
best = temp;
|
365 |
ibest = 6;
|
366 |
}
|
367 |
}
|
368 |
if (nentries[6] != 0) {
|
369 |
glmass_ = glmass[6][il];
|
370 |
lspmass_ = lspmass[6][il];
|
371 |
effsr7_ = eff[6][il];
|
372 |
errsr7_ = eff[6][il];
|
373 |
obslimsr7_ = obslim[6][il];
|
374 |
explimsr7_ = explim[6][il];
|
375 |
float temp = explimsr7_/effsr7_;
|
376 |
if (temp < best) {
|
377 |
best = temp;
|
378 |
ibest = 7;
|
379 |
}
|
380 |
}
|
381 |
bestsr_ = ibest;
|
382 |
effsrb_ = eff[ibest-1][il];
|
383 |
errsrb_ = eff[ibest-1][il];
|
384 |
obslimsrb_ = obslim[ibest-1][il];
|
385 |
explimsrb_ = explim[ibest-1][il];
|
386 |
|
387 |
// Fill the cross-section
|
388 |
bool done = false;
|
389 |
for (int xs=0; xs<nxsec-1; xs++) {
|
390 |
if ( (glmass_ >= xsec_mass[xs] && glmass_ < xsec_mass[xs+1]) ||
|
391 |
(glmass_ < xsec_mass[xs] && xs==0) ) {
|
392 |
|
393 |
float dd1 = glmass_ - xsec_mass[xs];
|
394 |
float dd2 = xsec_mass[xs+1] - xsec_mass[xs];
|
395 |
xsec_ = xsec[xs] + (xsec[xs+1] - xsec[xs]) *(dd1/dd2);
|
396 |
xsecup_ = xsecup[xs] + (xsecup[xs+1] - xsecup[xs]) *(dd1/dd2);
|
397 |
xsecdwn_ = xsecdwn[xs] + (xsecdwn[xs+1] - xsecdwn[xs])*(dd1/dd2);
|
398 |
done = true;
|
399 |
break;
|
400 |
}
|
401 |
}
|
402 |
if (!done) {
|
403 |
int xs = nxsec - 1;
|
404 |
float dd1 = glmass_ - xsec_mass[xs];
|
405 |
float dd2 = xsec_mass[xs] - xsec_mass[xs-1];
|
406 |
xsec_ = xsec[xs] + (xsec[xs] - xsec[xs-1]) *(dd1/dd2);
|
407 |
xsecup_ = xsecup[xs] + (xsecup[xs] - xsecup[xs-1]) *(dd1/dd2);
|
408 |
xsecdwn_ = xsecdwn[xs] + (xsecdwn[xs] - xsecdwn[xs-1])*(dd1/dd2);
|
409 |
}
|
410 |
|
411 |
|
412 |
// Now fill the ntuple
|
413 |
babyTree_->Fill();
|
414 |
|
415 |
}
|
416 |
// Now close the file
|
417 |
|
418 |
babyFile_->cd();
|
419 |
babyTree_->Write();
|
420 |
babyFile_->Close();
|
421 |
|
422 |
|
423 |
}
|
424 |
|