ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/SchieferD/jetcalib/mctopmass_x.cpp
Revision: 1.6
Committed: Sun Dec 9 18:12:17 2007 UTC (17 years, 4 months ago) by schiefer
Branch: MAIN
CVS Tags: HEAD
Changes since 1.5: +107 -30 lines
Error occurred while calculating annotation data.
Log Message:
*** empty log message ***

File Contents

# Content
1 ////////////////////////////////////////////////////////////////////////////////
2 //
3 // mctopmass_x
4 // -----------
5 //
6 // 11/24/2007 Philipp Schieferdecker <philipp.schieferdecker@cern.ch>
7 ////////////////////////////////////////////////////////////////////////////////
8
9 #include "utils/cmdline.h"
10 #include "jetcalib/jetcorr.h"
11
12 #include <TROOT.h>
13 #include <TStyle.h>
14 #include <TRint.h>
15 #include <TFile.h>
16 #include <TTree.h>
17 #include <TKey.h>
18 #include <TCanvas.h>
19 #include <TLegend.h>
20 #include <TLine.h>
21 #include <TText.h>
22 #include <TPaveText.h>
23 #include <TLatex.h>
24 #include <TH1F.h>
25 #include <THStack.h>
26 #include <TF1.h>
27 #include <TMultiGraph.h>
28 #include <TGraphErrors.h>
29 #include <TLorentzVector.h>
30
31 #include <iostream>
32 #include <iomanip>
33 #include <string>
34 #include <sstream>
35 #include <vector>
36 #include <set>
37 #include <cmath>
38
39
40 //
41 // define collection of histograms recorded per algorithm
42 //
43 class histograms
44 {
45 public:
46 histograms(const std::string& alg_name,float pt_min,float dr_max)
47 : algname(alg_name)
48 , ptmin(pt_min)
49 , drmax(dr_max)
50 , ndilep(0)
51 , nlepjets(0)
52 , nalljets(0)
53 , nlj_njt4(0)
54 , naj_njt6(0)
55 , ntop_lj(0)
56 , ntop_aj(0)
57 , nttbar_aj(0)
58 {
59 TDirectory* currentDir = gDirectory;
60 gROOT->cd();
61
62 hmt[0] =new TH1F(("mt_mc_" +algname).c_str(),"",100,174.9,175.1);
63 hmw[0] =new TH1F(("mw_mc_" +algname).c_str(),"",100, 80.3, 80.5);
64 hmt_lj[0]=new TH1F(("mt_mc_lj_"+algname).c_str(),"",100,174.9,175.1);
65 hmw_lj[0]=new TH1F(("mw_mc_lj_"+algname).c_str(),"",100, 80.3, 80.5);
66 hmt_aj[0]=new TH1F(("mt_mc_aj_"+algname).c_str(),"",100,174.9,175.1);
67 hmw_aj[0]=new TH1F(("mw_mc_aj_"+algname).c_str(),"",100, 80.3, 80.5);
68
69 hmt[1] =new TH1F(("mt_gen_" +algname).c_str(),"",100,0,300);
70 hmw[1] =new TH1F(("mw_gen_" +algname).c_str(),"",100,0,150);
71 hmt_lj[1]=new TH1F(("mt_gen_lj_"+algname).c_str(),"",100,0,300);
72 hmw_lj[1]=new TH1F(("mw_gen_lj_"+algname).c_str(),"",100,0,150);
73 hmt_aj[1]=new TH1F(("mt_gen_aj_"+algname).c_str(),"",100,0,300);
74 hmw_aj[1]=new TH1F(("mw_gen_aj_"+algname).c_str(),"",100,0,150);
75
76 hmt[2] =new TH1F(("mt_rec_" +algname).c_str(),"",100,0,300);
77 hmw[2] =new TH1F(("mw_rec_" +algname).c_str(),"",100,0,150);
78 hmt_lj[2]=new TH1F(("mt_rec_lj_"+algname).c_str(),"",100,0,300);
79 hmw_lj[2]=new TH1F(("mw_rec_lj_"+algname).c_str(),"",100,0,150);
80 hmt_aj[2]=new TH1F(("mt_rec_aj_"+algname).c_str(),"",100,0,300);
81 hmw_aj[2]=new TH1F(("mw_rec_aj_"+algname).c_str(),"",100,0,150);
82
83 hmt[3] =new TH1F(("mt_corr_" +algname).c_str(),"",100,0,300);
84 hmw[3] =new TH1F(("mw_corr_" +algname).c_str(),"",100,0,150);
85 hmt_lj[3]=new TH1F(("mt_corr_lj_"+algname).c_str(),"",100,0,300);
86 hmw_lj[3]=new TH1F(("mw_corr_lj_"+algname).c_str(),"",100,0,150);
87 hmt_aj[3]=new TH1F(("mt_corr_aj_"+algname).c_str(),"",100,0,300);
88 hmw_aj[3]=new TH1F(("mw_corr_aj_"+algname).c_str(),"",100,0,150);
89
90 hmt[4] =new TH1F(("mt_l5_" +algname).c_str(),"",100,0,300);
91 hmw[4] =new TH1F(("mw_l5_" +algname).c_str(),"",100,0,150);
92 hmt_lj[4]=new TH1F(("mt_l5_lj_"+algname).c_str(),"",100,0,300);
93 hmw_lj[4]=new TH1F(("mw_l5_lj_"+algname).c_str(),"",100,0,150);
94 hmt_aj[4]=new TH1F(("mt_l5_aj_"+algname).c_str(),"",100,0,300);
95 hmw_aj[4]=new TH1F(("mw_l5_aj_"+algname).c_str(),"",100,0,150);
96
97 hjetmult = new TH1F(("jetmult_" +algname).c_str(),"",31,-0.5,30.5);
98 hjetmult_lj = new TH1F(("jetmult_lj_"+algname).c_str(),"",31,-0.5,30.5);
99 hjetmult_aj = new TH1F(("jetmult_aj_"+algname).c_str(),"",31,-0.5,30.5);
100
101 hseleff_lj = new TH1F(("seleff_lj_"+algname).c_str(),"",46,5,50);
102 hseleff_aj = new TH1F(("seleff_aj_"+algname).c_str(),"",46,5,50);
103
104 hmw[0] ->SetXTitle("m_{q#bar{q}'} [GeV]");
105 hmw_lj[0]->SetXTitle("m_{q#bar{q}'} [GeV]");
106 hmw_aj[0]->SetXTitle("m_{q#bar{q}'} [GeV]");
107 hmt[0] ->SetXTitle("m_{bq#bar{q}'} [GeV]");
108 hmt_lj[0]->SetXTitle("m_{bq#bar{q}'} [GeV]");
109 hmt_aj[0]->SetXTitle("m_{bq#bar{q}'} [GeV]");
110
111 for (int i=1;i<5;i++) {
112 hmw[i] ->SetXTitle("m_{W} [GeV]");
113 hmw_lj[i]->SetXTitle("m_{W} [GeV]");
114 hmw_aj[i]->SetXTitle("m_{W} [GeV]");
115 hmt[i] ->SetXTitle("m_{top} [GeV]");
116 hmt_lj[i]->SetXTitle("m_{top} [GeV]");
117 hmt_aj[i]->SetXTitle("m_{top} [GeV]");
118 }
119 std::stringstream ssxjetmult;
120 ssxjetmult<<"jet multiplicity (p_{T}>"<<ptmin<<")";
121 hjetmult ->SetXTitle(ssxjetmult.str().c_str());
122 hjetmult_lj->SetXTitle(ssxjetmult.str().c_str());
123 hjetmult_aj->SetXTitle(ssxjetmult.str().c_str());
124 hseleff_lj->SetXTitle("p_{T}^{min} [GeV]");
125 hseleff_aj->SetXTitle("p_{T}^{min} [GeV]");
126 hseleff_lj->SetYTitle("#epsilon(Njet#geq4)");
127 hseleff_aj->SetYTitle("#epsilon(Njet#geq6)");
128
129 currentDir->cd();
130 }
131
132 histograms(TDirectory* d,float pt_min,float dr_max)
133 : algname(d->GetName())
134 , ptmin(pt_min)
135 , drmax(dr_max)
136 {
137 TH1I* hcounts = (TH1I*)d->Get(("counts_"+algname).c_str());
138 ndilep =(unsigned int)hcounts->GetBinContent(1);
139 nlepjets =(unsigned int)hcounts->GetBinContent(2);
140 nalljets =(unsigned int)hcounts->GetBinContent(3);
141 nlj_njt4 =(unsigned int)hcounts->GetBinContent(4);
142 naj_njt6 =(unsigned int)hcounts->GetBinContent(5);
143 ntop_lj =(unsigned int)hcounts->GetBinContent(6);
144 ntop_aj =(unsigned int)hcounts->GetBinContent(7);
145 nttbar_aj=(unsigned int)hcounts->GetBinContent(8);
146
147 hmt[0] =(TH1F*)d->Get(("mt_mc_" +algname).c_str());
148 hmw[0] =(TH1F*)d->Get(("mw_mc_" +algname).c_str());
149 hmt_lj[0]=(TH1F*)d->Get(("mt_mc_lj_"+algname).c_str());
150 hmw_lj[0]=(TH1F*)d->Get(("mw_mc_lj_"+algname).c_str());
151 hmt_aj[0]=(TH1F*)d->Get(("mt_mc_aj_"+algname).c_str());
152 hmw_aj[0]=(TH1F*)d->Get(("mw_mc_aj_"+algname).c_str());
153
154 hmt[1] =(TH1F*)d->Get(("mt_gen_" +algname).c_str());
155 hmw[1] =(TH1F*)d->Get(("mw_gen_" +algname).c_str());
156 hmt_lj[1]=(TH1F*)d->Get(("mt_gen_lj_"+algname).c_str());
157 hmw_lj[1]=(TH1F*)d->Get(("mw_gen_lj_"+algname).c_str());
158 hmt_aj[1]=(TH1F*)d->Get(("mt_gen_aj_"+algname).c_str());
159 hmw_aj[1]=(TH1F*)d->Get(("mw_gen_aj_"+algname).c_str());
160
161 hmt[2] =(TH1F*)d->Get(("mt_rec_" +algname).c_str());
162 hmw[2] =(TH1F*)d->Get(("mw_rec_" +algname).c_str());
163 hmt_lj[2]=(TH1F*)d->Get(("mt_rec_lj_"+algname).c_str());
164 hmw_lj[2]=(TH1F*)d->Get(("mw_rec_lj_"+algname).c_str());
165 hmt_aj[2]=(TH1F*)d->Get(("mt_rec_aj_"+algname).c_str());
166 hmw_aj[2]=(TH1F*)d->Get(("mw_rec_aj_"+algname).c_str());
167
168 hmt[3] =(TH1F*)d->Get(("mt_corr_" +algname).c_str());
169 hmw[3] =(TH1F*)d->Get(("mw_corr_" +algname).c_str());
170 hmt_lj[3]=(TH1F*)d->Get(("mt_corr_lj_"+algname).c_str());
171 hmw_lj[3]=(TH1F*)d->Get(("mw_corr_lj_"+algname).c_str());
172 hmt_aj[3]=(TH1F*)d->Get(("mt_corr_aj_"+algname).c_str());
173 hmw_aj[3]=(TH1F*)d->Get(("mw_corr_aj_"+algname).c_str());
174
175 hmt[4] =(TH1F*)d->Get(("mt_l5_" +algname).c_str());
176 hmw[4] =(TH1F*)d->Get(("mw_l5_" +algname).c_str());
177 hmt_lj[4]=(TH1F*)d->Get(("mt_l5_lj_"+algname).c_str());
178 hmw_lj[4]=(TH1F*)d->Get(("mw_l5_lj_"+algname).c_str());
179 hmt_aj[4]=(TH1F*)d->Get(("mt_l5_aj_"+algname).c_str());
180 hmw_aj[4]=(TH1F*)d->Get(("mw_l5_aj_"+algname).c_str());
181
182 hjetmult =(TH1F*)d->Get(("jetmult_" +algname).c_str());
183 hjetmult_lj=(TH1F*)d->Get(("jetmult_lj_"+algname).c_str());
184 hjetmult_aj=(TH1F*)d->Get(("jetmult_aj_"+algname).c_str());
185
186 hseleff_lj=(TH1F*)d->Get(("seleff_lj_"+algname).c_str());
187 hseleff_aj=(TH1F*)d->Get(("seleff_aj_"+algname).c_str());
188
189 set_directory(gROOT);
190 }
191
192 bool save_to_file()
193 {
194 std::stringstream ssfilename;
195 ssfilename<<"topmass_pt"<<ptmin<<"_dr"<<drmax<<".root";
196 TFile* f = new TFile(ssfilename.str().c_str(),"UPDATE");
197 TDirectory* d = (TDirectory*)f->GetDirectory(algname.c_str());
198 if (d==0) d=f->mkdir(algname.c_str());
199 d->DeleteAll();
200 d->cd();
201 TH1I* hcounts = new TH1I(("counts_"+algname).c_str(),"",8,0,7);
202 hcounts->SetBinContent(1,ndilep);
203 hcounts->SetBinContent(2,nlepjets);
204 hcounts->SetBinContent(3,nalljets);
205 hcounts->SetBinContent(4,nlj_njt4);
206 hcounts->SetBinContent(5,naj_njt6);
207 hcounts->SetBinContent(6,ntop_lj);
208 hcounts->SetBinContent(7,ntop_aj);
209 hcounts->SetBinContent(8,nttbar_aj);
210 set_directory(d);
211 f->Write();
212 set_directory(gROOT);
213 f->Close();
214 delete f;
215 gROOT->cd();
216 return true;
217 }
218
219 void set_directory(TDirectory* d)
220 {
221 for (int i=0;i<5;i++) {
222 hmw[i] ->SetDirectory(d);
223 hmt[i] ->SetDirectory(d);
224 hmw_lj[i]->SetDirectory(d);
225 hmt_lj[i]->SetDirectory(d);
226 hmw_aj[i]->SetDirectory(d);
227 hmt_aj[i]->SetDirectory(d);
228 }
229 hjetmult ->SetDirectory(d);
230 hjetmult_lj->SetDirectory(d);
231 hjetmult_aj->SetDirectory(d);
232 hseleff_lj ->SetDirectory(d);
233 hseleff_aj ->SetDirectory(d);
234 }
235
236 std::string algname;
237 float ptmin;
238 float drmax;
239 unsigned int ndilep;
240 unsigned int nlepjets;
241 unsigned int nalljets;
242 unsigned int nlj_njt4;
243 unsigned int naj_njt6;
244 unsigned int ntop_lj;
245 unsigned int ntop_aj;
246 unsigned int nttbar_aj;
247
248 TH1F* hmt[5];
249 TH1F* hmw[5];
250 TH1F* hmt_lj[5];
251 TH1F* hmw_lj[5];
252 TH1F* hmt_aj[5];
253 TH1F* hmw_aj[5];
254
255 TH1F* hjetmult;
256 TH1F* hjetmult_lj;
257 TH1F* hjetmult_aj;
258 TH1F* hseleff_lj;
259 TH1F* hseleff_aj;
260 };
261
262
263 //
264 // DEFINE GLOBAL FUNCTIONS
265 //
266 histograms* analyze_tree(const std::string& filename,float ptmin,float drmax);
267
268 void plot_mcmass(histograms* histos,const std::string& channel="");
269 void plot_topmass(histograms* histos,const std::string& channel="");
270 void plot_Wmass(histograms* histos,const std::string& channel="");
271 void plot_jetmult(histograms* histos,const std::string& channel="");
272
273 void table_seleff(std::vector<histograms*>vec_histos,const std::string& channel);
274 void plot_seleff(std::vector<histograms*>vec_histos,const std::string& channel);
275 void plot_seleff_vspt(std::vector<histograms*>vec_histos,const std::string& channel);
276
277 void plot_Wmass_mean(std::vector<histograms*> vec_histos,
278 const std::string& channel="");
279 void plot_topmass_mean(std::vector<histograms*> vec_histos,
280 const std::string& channel="");
281 void plot_Wmass_width(std::vector<histograms*> vec_histos,
282 const std::string& channel="");
283 void plot_topmass_width(std::vector<histograms*> vec_histos,
284 const std::string& channel="");
285
286 TGraphErrors* Wmass_mean(std::vector<histograms*> vec_histos,int level,
287 const std::string& channel="");
288 TGraphErrors* topmass_mean(std::vector<histograms*> vec_histos,int level,
289 const std::string& channel="");
290 TGraphErrors* Wmass_width(std::vector<histograms*> vec_histos,int level,
291 const std::string& channel="");
292 TGraphErrors* topmass_width(std::vector<histograms*> vec_histos,int level,
293 const std::string& channel="");
294
295 bool parse_filename(const std::string& filename,float& ptmin,float& drmax);
296
297 std::string alg_label(const std::string& alg_name);
298 Color_t alg_color(const std::string& alg_name);
299 Style_t alg_style(const std::string& alg_name);
300
301
302 using namespace std;
303
304
305 //
306 // MAIN
307 //
308
309 //______________________________________________________________________________
310 int main(int argc,char**argv)
311 {
312 cmdline cl;
313 cl.parse(argc,argv);
314
315 vector<string> input = cl.get_vector<string>("input");
316 float ptmin = cl.get_value<float> ("ptmin", 5.0);
317 float drmax = cl.get_value<float> ("drmax", 0.5);
318 string channel = cl.get_value<string> ("channel","");
319
320 bool show_all = cl.get_value<bool>("show_all", false);
321 bool show_mcmass = cl.get_value<bool>("show_mcmass", false);
322 bool show_Wmass = cl.get_value<bool>("show_Wmass", false);
323 bool show_topmass = cl.get_value<bool>("show_topmass", false);
324 bool show_jetmult = cl.get_value<bool>("show_jetmult", false);
325 bool show_seleff_table = cl.get_value<bool>("show_seleff_table", false);
326 bool show_seleff = cl.get_value<bool>("show_seleff", false);
327 bool show_seleff_vspt = cl.get_value<bool>("show_seleff_vspt", false);
328 bool show_Wmass_mean = cl.get_value<bool>("show_Wmass_mean", false);
329 bool show_topmass_mean = cl.get_value<bool>("show_topmass_mean", false);
330 bool show_Wmass_width = cl.get_value<bool>("show_Wmass_width", false);
331 bool show_topmass_width = cl.get_value<bool>("show_topmass_width",false);
332
333 if (!channel.empty()&&channel!="lj"&&channel!="aj") {
334 cout<<"invalid channel '"<<channel<<"'"<<endl;
335 return 0;
336 }
337
338 if(!cl.check()) return 0;
339 argc=1;
340 TRint* app = new TRint("mctopmass_x",&argc,argv);
341 cl.print();
342
343
344 vector<histograms*> vec_histos;
345
346 if (input[0].find("topmass")!=string::npos) {
347 if (!parse_filename(input[0],ptmin,drmax)) return 0;
348 set<string> jetalgs;
349 for (unsigned int i=1;i<input.size();i++) {
350 unsigned int pos;
351 string algname(input[i]);
352 while ((pos=algname.find('/')!=string::npos)) algname=algname.substr(pos+1);
353 if ((pos=algname.find(".root"))!=string::npos) algname=algname.substr(0,pos);
354 jetalgs.insert(algname);
355 }
356 TFile* f = new TFile(input[0].c_str(),"READ");
357 TIter next(f->GetListOfKeys());
358 TKey* key;
359 while ((key = (TKey*)next())!=0) {
360 string classname = key->GetClassName();
361 if (classname!="TDirectoryFile") continue;
362 TDirectoryFile* d = (TDirectoryFile*)key->ReadObj();
363 string dirname = d->GetName();
364 if (jetalgs.empty()||jetalgs.find(dirname)!=jetalgs.end())
365 vec_histos.push_back(new histograms(d,ptmin,drmax));
366 }
367 }
368 else {
369 for (vector<string>::const_iterator it=input.begin();it!=input.end();++it) {
370 histograms* histos = analyze_tree(*it,ptmin,drmax);
371 if (0!=histos) {
372 vec_histos.push_back(histos);
373 histos->save_to_file();
374 }
375 }
376 }
377
378 if (vec_histos.size()==0) return 0;
379
380 gStyle->SetOptStat(0);
381
382 vector<histograms*>::iterator it;
383 for (it=vec_histos.begin();it!=vec_histos.end();++it) {
384 if (show_mcmass ||show_all) plot_mcmass(*it,channel);
385 if (show_Wmass ||show_all) plot_Wmass(*it,channel);
386 if (show_topmass||show_all) plot_topmass(*it,channel);
387 if (show_jetmult||show_all) plot_jetmult(*it,channel);
388 }
389 if (show_seleff_table ||show_all) table_seleff(vec_histos,channel);
390 if (show_seleff ||show_all) plot_seleff(vec_histos,channel);
391 if (show_seleff_vspt ||show_all) plot_seleff_vspt(vec_histos,channel);
392 if (show_Wmass_mean ||show_all) plot_Wmass_mean(vec_histos,channel);
393 if (show_topmass_mean ||show_all) plot_topmass_mean(vec_histos,channel);
394 if (show_Wmass_width ||show_all) plot_Wmass_width(vec_histos,channel);
395 if (show_topmass_width||show_all) plot_topmass_width(vec_histos,channel);
396
397 app->Run();
398
399 return 0;
400 }
401
402
403 //
404 // IMPLEMENTATION OF GLOBAL FUNCTIONS
405 //
406
407 //______________________________________________________________________________
408 histograms* analyze_tree(const string& filename,float ptmin,float drmax)
409 {
410 TFile* f = new TFile(filename.c_str(),"READ"); if (!f->IsOpen()) return 0;
411 TTree* t = (TTree*)f->Get("t"); if (0==t) return 0;
412
413 char njt; t->SetBranchAddress("njt", &njt);
414
415 float jte[100]; t->SetBranchAddress("jte", jte);
416 float jtpt[100]; t->SetBranchAddress("jtpt", jtpt);
417 float jteta[100]; t->SetBranchAddress("jteta", jteta);
418 float jtphi[100]; t->SetBranchAddress("jtphi", jtphi);
419 float jtjes[100][1]; t->SetBranchAddress("jtjes", jtjes);
420
421 short jtgenflv[100]; t->SetBranchAddress("jtgenflv",jtgenflv);
422 float jtgene[100]; t->SetBranchAddress("jtgene", jtgene);
423 float jtgenpt[100]; t->SetBranchAddress("jtgenpt", jtgenpt);
424 float jtgeneta[100]; t->SetBranchAddress("jtgeneta",jtgeneta);
425 float jtgenphi[100]; t->SetBranchAddress("jtgenphi",jtgenphi);
426 float jtgendr[100]; t->SetBranchAddress("jtgendr", jtgendr);
427
428 char decay; t->SetBranchAddress("decay", &decay);
429 float jtmce[100]; t->SetBranchAddress("jtmce", jtmce);
430 float jtmcpt[100]; t->SetBranchAddress("jtmcpt", jtmcpt);
431 float jtmceta[100]; t->SetBranchAddress("jtmceta", jtmceta);
432 float jtmcphi[100]; t->SetBranchAddress("jtmcphi", jtmcphi);
433 float jtmcdr[100]; t->SetBranchAddress("jtmcdr", jtmcdr);
434
435 jetcorr l5_b;
436 jetcorr l5_c;
437 jetcorr l5_uds;
438
439 string path = getenv("CMSANA_DIR");
440 if (!l5_b .initialize(path+"/data/jetcorr/b.level5")) return 0;
441 if (!l5_c .initialize(path+"/data/jetcorr/c.level5")) return 0;
442 if (!l5_uds.initialize(path+"/data/jetcorr/uds.level5")) return 0;
443
444 string algname = filename;
445 unsigned int pos;
446 while ((pos = algname.find("/"))!=string::npos) algname = algname.substr(pos);
447 algname = algname.substr(0,algname.find("."));
448
449 histograms* histos = new histograms(algname,ptmin,drmax);
450
451 int nevts=t->GetEntries();
452
453 int ndl=0;
454 int nlj=0;
455 int naj=0;
456
457 int nlj_njt4 =0;
458 int naj_njt6 =0;
459
460 int ntop_lj =0;
461 int ntop_aj =0;
462 int nttbar_aj=0;
463
464 for (int ievt=0;ievt<nevts;ievt++) {
465
466 t->GetEntry(ievt);
467
468 if (decay==0) { naj++; }
469 if (decay==1) { nlj++; }
470 if (decay==2) { ndl++; continue; }
471
472 if (decay==1&&njt>=4) {
473 float pt = 5.0;
474 while (pt<=50.0&&pt<=jtpt[3]) { histos->hseleff_lj->Fill(pt); pt+=1.0; }
475 }
476
477 if (decay==0&&njt>=6) {
478 float pt = 5.0;
479 while (pt<=50.0&&pt<=jtpt[5]) { histos->hseleff_aj->Fill(pt); pt+=1.0; }
480 }
481
482 TLorentzVector* b[5] ={0,0,0,0,0};
483 TLorentzVector* bbar[5] ={0,0,0,0,0};
484 TLorentzVector* q1[5] ={0,0,0,0,0};
485 TLorentzVector* qbar1[5]={0,0,0,0,0};
486 TLorentzVector* q2[5] ={0,0,0,0,0};
487 TLorentzVector* qbar2[5]={0,0,0,0,0};
488
489 int jetmult(0);
490
491 for (int ijt=0;ijt<njt;ijt++) {
492
493 if (jtpt[ijt]<ptmin) continue;
494 jetmult++;
495
496 if (jtgenflv[ijt]==0||std::abs(jtgenflv[ijt])>5) continue;
497 if (jtgendr[ijt]>drmax) continue;
498
499 TLorentzVector* lvmc = new TLorentzVector();
500 lvmc->SetPtEtaPhiE(jtmcpt[ijt],jtmceta[ijt],jtmcphi[ijt],jtmce[ijt]);
501 TLorentzVector* lvgen = new TLorentzVector();
502 lvgen->SetPtEtaPhiE(jtgenpt[ijt],jtgeneta[ijt],jtgenphi[ijt],jtgene[ijt]);
503 TLorentzVector* lvrec = new TLorentzVector();
504 lvrec->SetPtEtaPhiE(jtpt[ijt],jteta[ijt],jtphi[ijt],jte[ijt]);
505 TLorentzVector* lvcorr = new TLorentzVector(*lvrec);
506 *lvcorr *= jtjes[ijt][0];
507 TLorentzVector* lvl5 = new TLorentzVector(*lvcorr);
508
509 if (jtgenflv[ijt]==+5) {
510 if (b[0]==0) {
511 *lvl5 *= l5_b.eval(lvl5->Pt(),lvl5->Eta());
512 b[0]=lvmc;
513 b[1]=lvgen;
514 b[2]=lvrec;
515 b[3]=lvcorr;
516 b[4]=lvl5;
517 }
518 }
519 else if (jtgenflv[ijt]==-5) {
520 if (bbar[0]==0) {
521 *lvl5 *= l5_b.eval(lvl5->Pt(),lvl5->Eta());
522 bbar[0]=lvmc;
523 bbar[1]=lvgen;
524 bbar[2]=lvrec;
525 bbar[3]=lvcorr;
526 bbar[4]=lvl5;
527 }
528 }
529 else if (jtgenflv[ijt]==+1) {
530 *lvl5 *= l5_uds.eval(lvl5->Pt(),lvl5->Eta());
531 if (q2[0]==0) {
532 q2[0]=lvmc;
533 q2[1]=lvgen;
534 q2[2]=lvrec;
535 q2[3]=lvcorr;
536 q2[4]=lvl5;
537 }
538 }
539 else if (jtgenflv[ijt]==-1) {
540 if (qbar1[0]==0) {
541 *lvl5 *= l5_uds.eval(lvl5->Pt(),lvl5->Eta());
542 qbar1[0]=lvmc;
543 qbar1[1]=lvgen;
544 qbar1[2]=lvrec;
545 qbar1[3]=lvcorr;
546 qbar1[4]=lvl5;
547 }
548 }
549 else if (jtgenflv[ijt]==+2) {
550 if (q1[0]==0) {
551 *lvl5 *= l5_uds.eval(lvl5->Pt(),lvl5->Eta());
552 q1[0]=lvmc;
553 q1[1]=lvgen;
554 q1[2]=lvrec;
555 q1[3]=lvcorr;
556 q1[4]=lvl5;
557 }
558 }
559 else if (jtgenflv[ijt]==-2) {
560 if (qbar2[0]==0) {
561 *lvl5 *= l5_uds.eval(lvl5->Pt(),lvl5->Eta());
562 qbar2[0]=lvmc;
563 qbar2[1]=lvgen;
564 qbar2[2]=lvrec;
565 qbar2[3]=lvcorr;
566 qbar2[4]=lvl5;
567 }
568 }
569 else if (jtgenflv[ijt]==+3) {
570 if (q2[0]==0) {
571 *lvl5 *= l5_uds.eval(lvl5->Pt(),lvl5->Eta());
572 q2[0]=lvmc;
573 q2[1]=lvgen;
574 q2[2]=lvrec;
575 q2[3]=lvcorr;
576 q2[4]=lvl5;
577 }
578 }
579 else if (jtgenflv[ijt]==-3) {
580 if (qbar1[0]==0) {
581 *lvl5 *= l5_uds.eval(lvl5->Pt(),lvl5->Eta());
582 qbar1[0]=lvmc;
583 qbar1[1]=lvgen;
584 qbar1[2]=lvrec;
585 qbar1[3]=lvcorr;
586 qbar1[4]=lvl5;
587 }
588 }
589 else if (jtgenflv[ijt]==+4) {
590 if (q1[0]==0) {
591 *lvl5 *= l5_c.eval(lvl5->Pt(),lvl5->Eta());
592 q1[0]=lvmc;
593 q1[1]=lvgen;
594 q1[2]=lvrec;
595 q1[3]=lvcorr;
596 q1[4]=lvl5;
597 }
598 }
599 else if (jtgenflv[ijt]==-4) {
600 if (qbar2[0]==0) {
601 *lvl5 *= l5_c.eval(lvl5->Pt(),lvl5->Eta());
602 qbar2[0]=lvmc;
603 qbar2[1]=lvgen;
604 qbar2[2]=lvrec;
605 qbar2[3]=lvcorr;
606 qbar2[4]=lvl5;
607 }
608 }
609 else assert(0);
610 }
611
612 bool has_t = (b[0] !=0&&q1[0]!=0&&qbar1[0]!=0);
613 bool has_tbar = (bbar[0]!=0&&q2[0]!=0&&qbar2[0]!=0);
614
615 if (decay==0&&has_t&&has_tbar) nttbar_aj++;
616
617 assert(decay==0||(!has_t||!has_tbar));
618
619 histos->hjetmult->Fill(jetmult);
620 if (decay==1) histos->hjetmult_lj->Fill(jetmult);
621 if (decay==0) histos->hjetmult_aj->Fill(jetmult);
622
623 if (decay==1&&jetmult>=4) nlj_njt4++;
624 if (decay==0&&jetmult>=6) naj_njt6++;
625
626 if (has_t) {
627 if (decay==1) ntop_lj++;
628 if (decay==0) ntop_aj++;
629 for (int i=0;i<5;i++) {
630 TLorentzVector Wplus = *q1[i] + *qbar1[i];
631 TLorentzVector t = *b[i] + Wplus;
632 histos->hmw[i]->Fill(Wplus.M());
633 histos->hmt[i]->Fill(t.M());
634 if (decay==1) {
635 histos->hmw_lj[i]->Fill(Wplus.M());
636 histos->hmt_lj[i]->Fill(t.M());
637 }
638 else if (decay==0) {
639 histos->hmw_aj[i]->Fill(Wplus.M());
640 histos->hmt_aj[i]->Fill(t.M());
641 }
642 }
643 }
644
645 if (has_tbar) {
646 if (decay==1) ntop_lj++;
647 if (decay==0) ntop_aj++;
648 for (int i=0;i<5;i++) {
649 TLorentzVector Wminus = *q2[i] + *qbar2[i];
650 TLorentzVector tbar = *bbar[i] + Wminus;
651 histos->hmw[i]->Fill(Wminus.M());
652 histos->hmt[i]->Fill(tbar.M());
653 if (decay==1) {
654 histos->hmw_lj[i]->Fill(Wminus.M());
655 histos->hmt_lj[i]->Fill(tbar.M());
656 }
657 else if (decay==0) {
658 histos->hmw_aj[i]->Fill(Wminus.M());
659 histos->hmt_aj[i]->Fill(tbar.M());
660 }
661 }
662 }
663
664 if (b[0]) { for (int i=0;i<5;i++) delete b[i]; }
665 if (bbar[0]) { for (int i=0;i<5;i++) delete bbar[i]; }
666 if (q1[0]) { for (int i=0;i<5;i++) delete q1[i]; }
667 if (qbar1[0]){ for (int i=0;i<5;i++) delete qbar1[i];}
668 if (q2[0]) { for (int i=0;i<5;i++) delete q2[i]; }
669 if (qbar2[0]){ for (int i=0;i<5;i++) delete qbar2[i];}
670
671 }
672
673 for (int i=1;i<=histos->hseleff_lj->GetNbinsX();i++) {
674 float nsel=histos->hseleff_lj->GetBinContent(i);
675 float nall=nlj;
676 float eff = (nall>0) ? nsel/nall : 1.0;
677 float eeff = (nall>0) ? std::sqrt(nsel/nall/nall*(1-nsel/nall)):0.0;
678 histos->hseleff_lj->SetBinContent(i,eff);
679 histos->hseleff_lj->SetBinError(i,eeff);
680 }
681
682 for (int i=1;i<=histos->hseleff_aj->GetNbinsX();i++) {
683 float nsel=histos->hseleff_aj->GetBinContent(i);
684 float nall=naj;
685 float eff = (nall>0) ? nsel/nall : 1.0;
686 float eeff = (nall>0) ? std::sqrt(nsel/nall/nall*(1-nsel/nall)):0.0;
687 histos->hseleff_aj->SetBinContent(i,eff);
688 histos->hseleff_aj->SetBinError(i,eeff);
689 }
690
691 histos->ndilep = ndl;
692 histos->nlepjets = nlj;
693 histos->nalljets = naj;
694
695 histos->nlj_njt4 = nlj_njt4;
696 histos->naj_njt6 = naj_njt6;
697 histos->ntop_lj = ntop_lj;
698 histos->ntop_aj = ntop_aj;
699 histos->nttbar_aj = nttbar_aj;
700
701 return histos;
702 }
703
704 //______________________________________________________________________________
705 void plot_mcmass(histograms* histos,const string& channel)
706 {
707 stringstream sscname;
708 string ch = (channel.empty()) ? channel : "_"+channel;
709 sscname<<"mcmass"<<ch<<"_"<<histos->algname
710 <<"_pt"<<histos->ptmin<<"_dr0"<<histos->drmax*10;
711 TCanvas* cmcmass = new TCanvas(sscname.str().c_str(),
712 sscname.str().c_str(),
713 1000,500);
714 cmcmass->Divide(2,1);
715
716 TH1F* hmw =
717 (channel.empty()) ? histos->hmw[0] :
718 (channel=="lj") ? histos->hmw_lj[0] : histos->hmw_aj[0];
719 TH1F* hmt =
720 (channel.empty()) ? histos->hmt[0] :
721 (channel=="lj") ? histos->hmt_lj[0] : histos->hmt_aj[0];
722
723 cmcmass->cd(1);
724 gPad->SetLogy();
725 hmw->Draw();
726 hmw->SetTitle("W mass (parton-level)");
727 hmw->SetFillStyle(1001);
728 hmw->SetFillColor(kBlue);
729 stringstream meanW; meanW<<"Mean: "<<setprecision(4)<<hmw->GetMean();
730 stringstream rmsW; rmsW <<"RMS : "<<setprecision(1)<<hmw->GetRMS();
731 TText txtW; txtW.SetTextColor(kBlue);
732 txtW.DrawTextNDC(0.23,0.8, meanW.str().c_str());
733 txtW.DrawTextNDC(0.23,0.74,rmsW.str().c_str());
734
735 cmcmass->cd(2);
736 gPad->SetLogy();
737 hmt->Draw();
738 hmt->SetTitle("top mass (parton-level)");
739 hmt->SetFillStyle(1001);
740 hmt->SetFillColor(kRed);
741 stringstream meantop; meantop<<"Mean: "<<setprecision(3)<<hmt->GetMean();
742 stringstream rmstop; rmstop <<"RMS : "<<setprecision(1)<<hmt->GetRMS();
743 TText txttop; txttop.SetTextColor(kRed);
744 txttop.DrawTextNDC(0.23,0.8, meantop.str().c_str());
745 txttop.DrawTextNDC(0.23,0.74,rmstop.str().c_str());
746
747 cmcmass->Modified();
748 }
749
750
751 //______________________________________________________________________________
752 void plot_topmass(histograms* histos,const string& channel)
753 {
754 stringstream sscname;
755 string ch = (channel.empty()) ? channel : "_"+channel;
756 sscname<<"topmass"<<ch<<"_"<<histos->algname
757 <<"_pt"<<histos->ptmin<<"_dr0"<<histos->drmax*10;
758 TCanvas* ctopmass = new TCanvas(sscname.str().c_str(),
759 sscname.str().c_str(),
760 700,700);
761 ctopmass->SetLeftMargin(0.15);
762 ctopmass->SetRightMargin(0.05);
763
764 TH1F* hgen =
765 (channel.empty()) ? histos->hmt[1] :
766 (channel=="lj") ? histos->hmt_lj[1] : histos->hmt_aj[1];
767 TH1F* hrec =
768 (channel.empty()) ? histos->hmt[2] :
769 (channel=="lj") ? histos->hmt_lj[2] : histos->hmt_aj[2];
770 TH1F* hcorr =
771 (channel.empty()) ? histos->hmt[3] :
772 (channel=="lj") ? histos->hmt_lj[3] : histos->hmt_aj[3];
773 TH1F* hl5 =
774 (channel.empty()) ? histos->hmt[4] :
775 (channel=="lj") ? histos->hmt_lj[4] : histos->hmt_aj[4];
776
777 hgen->Draw();
778 hrec->Draw("SAME");
779 hcorr->Draw("SAME");
780 hl5->Draw("SAME");
781
782 hgen->SetTitle("top mass: GEN-REC-CORR-L5");
783 hgen ->SetLineStyle(4);
784 hrec ->SetLineStyle(3);
785 hcorr->SetLineStyle(2);
786 hl5 ->SetLineStyle(1);
787 hgen ->SetFillStyle(1001);
788 hgen ->SetFillColor(kGreen);
789 hcorr ->SetFillStyle(1001);
790 hcorr ->SetFillColor(kOrange+1);
791 hl5 ->SetFillStyle(1001);
792 hl5 ->SetFillColor(kRed);
793 stringstream genmean; genmean <<"Mean: "<<setprecision(4)<<hgen->GetMean();
794 stringstream genrms; genrms <<"RMS : "<<setprecision(3)<<hgen->GetRMS();
795 stringstream recmean; recmean <<"Mean: "<<setprecision(4)<<hrec->GetMean();
796 stringstream recrms; recrms <<"RMS : "<<setprecision(3)<<hrec->GetRMS();
797 stringstream corrmean;corrmean<<"Mean: "<<setprecision(4)<<hcorr->GetMean();
798 stringstream corrrms; corrrms <<"RMS : "<<setprecision(3)<<hcorr->GetRMS();
799 stringstream l5mean; l5mean <<"Mean: "<<setprecision(4)<<hl5->GetMean();
800 stringstream l5rms; l5rms <<"RMS : "<<setprecision(3)<<hl5->GetRMS();
801
802 TText gentxt; gentxt.SetTextColor(kGreen);
803 TText rectxt; rectxt.SetTextColor(kBlack);
804 TText corrtxt;corrtxt.SetTextColor(kOrange+1);
805 TText l5txt; l5txt.SetTextColor(kRed);
806
807 gentxt .DrawTextNDC(0.35,0.80,"GEN");
808 gentxt .DrawTextNDC(0.35,0.75,genmean.str().c_str());
809 gentxt .DrawTextNDC(0.35,0.70,genrms.str().c_str());
810 rectxt .DrawTextNDC(0.17,0.50,"REC");
811 rectxt .DrawTextNDC(0.17,0.45,recmean.str().c_str());
812 rectxt .DrawTextNDC(0.17,0.40,recrms.str().c_str());
813 corrtxt.DrawTextNDC(0.73,0.35,"CORR");
814 corrtxt.DrawTextNDC(0.73,0.30,corrmean.str().c_str());
815 corrtxt.DrawTextNDC(0.73,0.25,corrrms.str().c_str());
816 l5txt .DrawTextNDC(0.66,0.53,"L5(CORR+FLV)");
817 l5txt .DrawTextNDC(0.66,0.48,l5mean.str().c_str());
818 l5txt .DrawTextNDC(0.66,0.43,l5rms.str().c_str());
819
820 TLatex algtxt;
821 algtxt.SetNDC();
822 algtxt.DrawLatex(0.66,0.885,alg_label(histos->algname).c_str());
823
824 hgen->SetMaximum(1.1*hgen->GetMaximum());
825 TLine* mtopline = new TLine(175,hgen->GetMinimum(),175,hgen->GetMaximum());
826 mtopline->SetLineWidth(2);
827 mtopline->Draw("SAME");
828
829 ctopmass->Modified();
830 }
831
832
833 //______________________________________________________________________________
834 void plot_Wmass(histograms* histos,const string& channel)
835 {
836 stringstream sscname;
837 string ch = (channel.empty()) ? channel : "_"+channel;
838 sscname<<"Wmass"<<ch<<"_"<<histos->algname
839 <<"_pt"<<histos->ptmin<<"_dr0"<<histos->drmax*10;
840 TCanvas* cWmass = new TCanvas(sscname.str().c_str(),
841 sscname.str().c_str(),
842 700,700);
843 cWmass->SetLeftMargin(0.15);
844 cWmass->SetRightMargin(0.05);
845
846 TH1F* hgen =
847 (channel.empty()) ? histos->hmw[1] :
848 (channel=="lj") ? histos->hmw_lj[1] : histos->hmw_aj[1];
849 TH1F* hrec =
850 (channel.empty()) ? histos->hmw[2] :
851 (channel=="lj") ? histos->hmw_lj[2] : histos->hmw_aj[2];
852 TH1F* hcorr =
853 (channel.empty()) ? histos->hmw[3] :
854 (channel=="lj") ? histos->hmw_lj[3] : histos->hmw_aj[3];
855 TH1F* hl5 =
856 (channel.empty()) ? histos->hmw[4] :
857 (channel=="lj") ? histos->hmw_lj[4] : histos->hmw_aj[4];
858
859 hgen ->Draw();
860 hrec ->Draw("SAME");
861 hcorr->Draw("SAME");
862 hl5 ->Draw("SAME");
863
864 hgen->SetTitle("W mass: GEN-REC-CORR-L5");
865 hgen ->SetLineStyle(4);
866 hrec ->SetLineStyle(3);
867 hcorr->SetLineStyle(2);
868 hl5 ->SetLineStyle(1);
869 hgen ->SetFillStyle(1001);
870 hgen ->SetFillColor(kGreen);
871 hcorr ->SetFillStyle(1001);
872 hcorr ->SetFillColor(kAzure+1);
873 hl5 ->SetFillStyle(1001);
874 hl5 ->SetFillColor(kBlue);
875 stringstream genmean; genmean <<"Mean: "<<setprecision(3)<<hgen->GetMean();
876 stringstream genrms; genrms <<"RMS : "<<setprecision(3)<<hgen->GetRMS();
877 stringstream recmean; recmean <<"Mean: "<<setprecision(3)<<hrec->GetMean();
878 stringstream recrms; recrms <<"RMS : "<<setprecision(3)<<hrec->GetRMS();
879 stringstream corrmean;corrmean<<"Mean: "<<setprecision(3)<<hcorr->GetMean();
880 stringstream corrrms; corrrms <<"RMS : "<<setprecision(3)<<hcorr->GetRMS();
881 stringstream l5mean; l5mean <<"Mean: "<<setprecision(3)<<hl5->GetMean();
882 stringstream l5rms; l5rms <<"RMS : "<<setprecision(3)<<hl5->GetRMS();
883
884 TText gentxt; gentxt.SetTextColor(kGreen);
885 TText rectxt; rectxt.SetTextColor(kBlack);
886 TText corrtxt; corrtxt.SetTextColor(kAzure+1);
887 TText l5txt; l5txt.SetTextColor(kBlue);
888
889 gentxt .DrawTextNDC(0.30,0.80,"GEN");
890 gentxt .DrawTextNDC(0.30,0.75,genmean.str().c_str());
891 gentxt .DrawTextNDC(0.30,0.70,genrms.str().c_str());
892 rectxt .DrawTextNDC(0.17,0.50,"REC");
893 rectxt .DrawTextNDC(0.17,0.45,recmean.str().c_str());
894 rectxt .DrawTextNDC(0.17,0.40,recrms.str().c_str());
895 corrtxt.DrawTextNDC(0.72,0.35,"CORR");
896 corrtxt.DrawTextNDC(0.72,0.30,corrmean.str().c_str());
897 corrtxt.DrawTextNDC(0.72,0.25,corrrms.str().c_str());
898 l5txt .DrawTextNDC(0.61,0.53,"L5(CORR+FLV)");
899 l5txt .DrawTextNDC(0.61,0.48,l5mean.str().c_str());
900 l5txt .DrawTextNDC(0.61,0.43,l5rms.str().c_str());
901
902 TLatex algtxt;
903 algtxt.SetNDC();
904 algtxt.DrawLatex(0.66,0.885,alg_label(histos->algname).c_str());
905
906 hgen->SetMaximum(1.1*hgen->GetMaximum());
907 TLine* mWline = new TLine(80.42,hgen->GetMinimum(),80.42,hgen->GetMaximum());
908 mWline->SetLineWidth(2);
909 mWline->Draw("SAME");
910
911 cWmass->Modified();
912 }
913
914
915 //______________________________________________________________________________
916 void plot_jetmult(histograms* histos,const string& channel)
917 {
918 stringstream sscname;
919 string ch = (channel.empty()) ? channel : "_"+channel;
920 sscname<<"jetmult"<<ch<<"_"<<histos->algname
921 <<"_pt"<<histos->ptmin<<"_dr0"<<histos->drmax*10;
922 TCanvas* cjetmult = new TCanvas(sscname.str().c_str(),
923 sscname.str().c_str(),
924 700,700);
925
926 TH1F* h(0);
927 if (channel.empty()) h=histos->hjetmult;
928 else if (channel=="lj") h=histos->hjetmult_lj;
929 else if (channel=="aj") h=histos->hjetmult_aj;
930 h->Draw("BAR");
931 h->SetTitle("jet multiplicity");
932
933 stringstream ssmean; ssmean<<"Mean: "<<setprecision(2)<<h->GetMean();
934 TText meantxt; meantxt.DrawTextNDC(0.6,0.6,ssmean.str().c_str());
935
936 TLatex algtxt;
937 algtxt.SetNDC();
938 algtxt.DrawLatex(0.62,0.885,alg_label(histos->algname).c_str());
939 if (!channel.empty()) {
940 string ch = (channel=="lj") ? "lepton+jets" : "alljets";
941 Color_t clr = (channel=="lj") ? kYellow-9 : kGreen-9;
942 TPaveText* ptch = new TPaveText(0.6,0.7,0.85,0.8,"brNDC");
943 ptch->SetFillColor(clr);
944 ptch->AddText(ch.c_str());
945 ptch->Draw();
946 }
947 cjetmult->Modified();
948 }
949
950
951 //______________________________________________________________________________
952 void table_seleff(vector<histograms*>vec_histos,const string& channel)
953 {
954 if (channel.empty()||channel=="lj") {
955 float nlj = vec_histos[0]->nlepjets;
956 vector<float> eff_njt;
957 vector<float> eff_bqq;
958 for (vector<histograms*>::iterator it=vec_histos.begin();
959 it!=vec_histos.end();++it) {
960 eff_njt.push_back(100.*(*it)->nlj_njt4/nlj);
961 eff_bqq.push_back(100.*(*it)->ntop_lj/(float)((*it)->nlj_njt4));
962 }
963 cout<<"lepton+jets, eff(njt>=4):"<<endl;
964 for (vector<float>::iterator it=eff_njt.begin();it!=eff_njt.end();++it)
965 cout<<" & "<<setprecision(3)<<(*it)<<flush;
966 cout<<" \\\\\n"<<endl;
967 cout<<"lepton+jets, eff(t->bqq'):"<<endl;
968 for (vector<float>::iterator it=eff_bqq.begin();it!=eff_bqq.end();++it)
969 cout<<" & "<<setprecision(3)<<(*it)<<flush;
970 cout<<" \\\\\n"<<endl;
971 }
972 if (channel.empty()||channel=="aj") {
973 float naj = vec_histos[0]->nalljets;
974 vector<float> eff_njt;
975 vector<float> eff_bqq;
976 vector<float> eff_bbqqqq;
977 for (vector<histograms*>::iterator it=vec_histos.begin();
978 it!=vec_histos.end();++it) {
979 eff_njt.push_back(100.*(*it)->naj_njt6/naj);
980 eff_bqq.push_back(100.*(*it)->ntop_aj/2./(float)((*it)->naj_njt6));
981 eff_bbqqqq.push_back(100.*(*it)->nttbar_aj/(float)((*it)->naj_njt6));
982 }
983 cout<<"alljets, eff(njt>=6):"<<endl;
984 for (vector<float>::iterator it=eff_njt.begin();it!=eff_njt.end();++it)
985 cout<<" & "<<setprecision(3)<<(*it)<<flush;
986 cout<<" \\\\\n"<<endl;
987 cout<<"alljets, eff(t->bqq'):"<<endl;
988 for (vector<float>::iterator it=eff_bqq.begin();it!=eff_bqq.end();++it)
989 cout<<" & "<<setprecision(3)<<(*it)<<flush;
990 cout<<" \\\\\n"<<endl;
991 cout<<"alljets, eff(ttbar->bbqq'qq'):"<<endl;
992 for (vector<float>::iterator it=eff_bbqqqq.begin();it!=eff_bbqqqq.end();++it)
993 cout<<" & "<<setprecision(3)<<(*it)<<flush;
994 cout<<" \\\\\n"<<endl;
995 }
996
997 }
998
999
1000 //______________________________________________________________________________
1001 void plot_seleff(vector<histograms*> vec_histos,const string& channel)
1002 {
1003 stringstream sscname_lj;
1004 sscname_lj<<"seleff_lj_pt"<<vec_histos[0]->ptmin<<"_dr0"<<vec_histos[0]->drmax*10;
1005 stringstream sscname_aj;
1006 sscname_aj<<"seleff_aj_pt"<<vec_histos[0]->ptmin<<"_dr0"<<vec_histos[0]->drmax*10;
1007
1008 if (channel.empty()||channel=="lj") {
1009 TCanvas* cseleff_lj = new TCanvas(sscname_lj.str().c_str(),
1010 sscname_lj.str().c_str(),
1011 1000,350);
1012 cseleff_lj->SetLogy();
1013 cseleff_lj->SetGridy();
1014 cseleff_lj->SetLeftMargin(0.12);
1015 cseleff_lj->SetRightMargin(0.05);
1016 cseleff_lj->SetTopMargin(0.23);
1017
1018 TLegend* leg_lj = new TLegend(0.50,0.785,0.95,0.99);
1019 leg_lj->SetFillColor(10);
1020
1021 int nbins = vec_histos.size()*2+1;
1022 TH1F* hnsel_lj=new TH1F("nsel_lj","",nbins,0,nbins+1);
1023 TH1F* hntop_lj=new TH1F("ntop_lj","",nbins,0,nbins+1);
1024 THStack* stack_lj=new THStack("seleff_lj","#epsilon(t#rightarrow bq#bar{q}')");
1025 stack_lj->Add(hntop_lj);
1026 stack_lj->Add(hnsel_lj);
1027
1028 hnsel_lj->SetFillStyle(1001);
1029 hntop_lj->SetFillStyle(1001);
1030 hnsel_lj->SetFillColor(kBlue+2);
1031 hntop_lj->SetFillColor(kYellow);
1032
1033 for (unsigned int i=0;i<vec_histos.size();i++) {
1034 hnsel_lj->SetBinContent(i*2+2,vec_histos[i]->nlj_njt4);
1035 hntop_lj->SetBinContent(i*2+2,vec_histos[i]->ntop_lj);
1036 }
1037 stack_lj->Draw();
1038
1039 float ptmin = vec_histos[0]->ptmin;
1040 stringstream ssnsel;ssnsel<<"N(jet)#geq4 (jet pT#geq"<<ptmin<<")";
1041 stringstream ssntop;ssntop<<"N(t#rightarrow bq#bar{q}') "
1042 <<"matched (jet pT#geq"<<ptmin<<")";
1043 leg_lj->AddEntry(hnsel_lj,ssnsel.str().c_str(),"F");
1044 leg_lj->AddEntry(hntop_lj,ssntop.str().c_str(),"F");
1045 leg_lj->Draw();
1046
1047 for (unsigned int i=0;i<vec_histos.size();i++)
1048 stack_lj->GetHistogram()->GetXaxis()
1049 ->SetBinLabel(i*2+2,alg_label(vec_histos[i]->algname).c_str());
1050 stack_lj->GetHistogram()->LabelsOption("u");
1051 stack_lj->GetHistogram()->GetXaxis()->SetLabelOffset(0.01);
1052 stack_lj->GetHistogram()->GetXaxis()->SetLabelSize(0.065);
1053
1054 TPaveText* pt_lj = new TPaveText(0.12,0.82,0.3,0.90,"brNDC");
1055 pt_lj->SetFillColor(kYellow-9);
1056 pt_lj->AddText("lepton+jets");
1057 pt_lj->Draw();
1058
1059 cseleff_lj->Modified();
1060 }
1061 if (channel.empty()||channel=="aj") {
1062 TCanvas* cseleff_aj = new TCanvas(sscname_aj.str().c_str(),
1063 sscname_aj.str().c_str(),
1064 1000,350);
1065 cseleff_aj->SetLogy();
1066 cseleff_aj->SetGridy();
1067 cseleff_aj->SetLeftMargin(0.12);
1068 cseleff_aj->SetRightMargin(0.05);
1069 cseleff_aj->SetTopMargin(0.23);
1070
1071 TLegend* leg_aj = new TLegend(0.50,0.785,0.95,0.99);
1072 leg_aj->SetFillColor(10);
1073
1074 int nbins = vec_histos.size()*2+1;
1075 TH1F* hnsel_aj =new TH1F("nsel_aj", "",nbins,0,nbins+1);
1076 TH1F* hntop_aj =new TH1F("ntop_aj", "",nbins,0,nbins+1);
1077 TH1F* hnttbar_aj=new TH1F("nttbar_aj","",nbins,0,nbins+1);
1078 THStack* stack_aj =new THStack("seleff_aj","#epsilon(t#rightarrow bq#bar{q}')");
1079 stack_aj->Add(hnttbar_aj);
1080 stack_aj->Add(hntop_aj);
1081 stack_aj->Add(hnsel_aj);
1082
1083 hnsel_aj ->SetFillStyle(1001);
1084 hntop_aj ->SetFillStyle(1001);
1085 hnttbar_aj->SetFillStyle(1001);
1086 hnsel_aj ->SetFillColor(kBlue+2);
1087 hntop_aj ->SetFillColor(kYellow);
1088 hnttbar_aj->SetFillColor(kRed);
1089
1090 for (unsigned int i=0;i<vec_histos.size();i++) {
1091 hnsel_aj ->SetBinContent(i*2+2,vec_histos[i]->nlj_njt4);
1092 hntop_aj ->SetBinContent(i*2+2,vec_histos[i]->ntop_aj);
1093 hnttbar_aj->SetBinContent(i*2+2,vec_histos[i]->nttbar_aj);
1094 }
1095 stack_aj->Draw();
1096
1097 float ptmin = vec_histos[0]->ptmin;
1098 stringstream ssnsel; ssnsel<<"N(jet)#geq6 (jet pT#geq"<<ptmin<<")";
1099 stringstream ssntop; ssntop<<"N(t#rightarrow bq#bar{q}') "
1100 <<"matched (jet pT#geq"<<ptmin<<")";
1101 stringstream ssnttbar;ssnttbar<<"N(t#bar{t}#rightarrow "
1102 <<"b#bar{b}q#bar{q}'q#bar{q}') "
1103 <<"matched (jet pT#geq"<<ptmin<<")";
1104 leg_aj->AddEntry(hnsel_aj, ssnsel .str().c_str(),"F");
1105 leg_aj->AddEntry(hntop_aj, ssntop .str().c_str(),"F");
1106 leg_aj->AddEntry(hnttbar_aj,ssnttbar.str().c_str(),"F");
1107 leg_aj->Draw();
1108
1109 for (unsigned int i=0;i<vec_histos.size();i++)
1110 stack_aj->GetHistogram()->GetXaxis()
1111 ->SetBinLabel(i*2+2,alg_label(vec_histos[i]->algname).c_str());
1112 stack_aj->GetHistogram()->LabelsOption("u");
1113 stack_aj->GetHistogram()->GetXaxis()->SetLabelOffset(0.01);
1114 stack_aj->GetHistogram()->GetXaxis()->SetLabelSize(0.065);
1115
1116 TPaveText* pt_aj = new TPaveText(0.12,0.82,0.3,0.90,"brNDC");
1117 pt_aj->SetFillColor(kGreen-9);
1118 pt_aj->AddText("alljets");
1119 pt_aj->Draw();
1120
1121 cseleff_aj->Modified();
1122 }
1123 }
1124
1125
1126 //______________________________________________________________________________
1127 void plot_seleff_vspt(vector<histograms*> vec_histos,const string& channel)
1128 {
1129 stringstream sscname_lj;
1130 sscname_lj<<"seleff_vspt_lj_dr0"<<vec_histos[0]->drmax*10;
1131 stringstream sscname_aj;
1132 sscname_aj<<"seleff_vspt_aj_dr0"<<vec_histos[0]->drmax*10;
1133
1134 std::vector<histograms*>::iterator it;
1135
1136 if (channel.empty()||channel=="lj") {
1137 TCanvas* cseleff_lj = new TCanvas(sscname_lj.str().c_str(),
1138 sscname_lj.str().c_str(),
1139 700,700);
1140 TLegend* leg_lj = new TLegend(0.55,0.85,0.85,0.85-vec_histos.size()*0.045);
1141 leg_lj->SetFillColor(10);
1142 string draw_option="PEX0";
1143 for (it=vec_histos.begin();it!=vec_histos.end();++it) {
1144 TH1F* h = (*it)->hseleff_lj;
1145 h->SetTitle("Jet Selection Effificency");
1146 h->SetLineColor(alg_color((*it)->algname));
1147 h->SetMarkerColor(alg_color((*it)->algname));
1148 h->SetMarkerStyle(alg_style((*it)->algname));
1149 h->Draw(draw_option.c_str());
1150 if (draw_option=="PEX0") draw_option+="SAME";
1151 leg_lj->AddEntry(h,alg_label((*it)->algname).c_str());
1152 }
1153 leg_lj->Draw();
1154
1155 TPaveText* pt_lj = new TPaveText(0.67,0.9,0.9,0.98,"brNDC");
1156 pt_lj->SetFillColor(kYellow-9);
1157 pt_lj->AddText("lepton+jets");
1158 pt_lj->Draw();
1159
1160 cseleff_lj->Modified();
1161 }
1162 if (channel.empty()||channel=="aj") {
1163 TCanvas* cseleff_aj = new TCanvas(sscname_aj.str().c_str(),
1164 sscname_aj.str().c_str(),
1165 700,700);
1166 TLegend* leg_aj = new TLegend(0.55,0.85,0.85,0.85-vec_histos.size()*0.045);
1167 leg_aj->SetFillColor(10);
1168 string draw_option="PEX0";
1169 for (it=vec_histos.begin();it!=vec_histos.end();++it) {
1170 TH1F* h = (*it)->hseleff_aj;
1171 h->SetTitle("Jet Selection Efficiency");
1172 h->SetLineColor(alg_color((*it)->algname));
1173 h->SetMarkerColor(alg_color((*it)->algname));
1174 h->SetMarkerStyle(alg_style((*it)->algname));
1175 h->Draw(draw_option.c_str());
1176 if (draw_option=="PEX0") draw_option+="SAME";
1177 leg_aj->AddEntry(h,alg_label((*it)->algname).c_str());
1178 }
1179 leg_aj->Draw();
1180
1181 TPaveText* pt_aj = new TPaveText(0.7,0.9,0.9,0.98,"brNDC");
1182 pt_aj->SetFillColor(kGreen-9);
1183 pt_aj->AddText("alljets");
1184 pt_aj->Draw();
1185
1186 cseleff_aj->Modified();
1187 }
1188 }
1189
1190
1191 //______________________________________________________________________________
1192 void plot_Wmass_mean(std::vector<histograms*> vec_histos,
1193 const std::string& channel)
1194 {
1195 stringstream sscname;
1196 string ch = (channel.empty()) ? channel : "_"+channel;
1197 sscname<<"WmassMean"<<ch<<"_pt"<<vec_histos[0]->ptmin
1198 <<"_dr0"<<vec_histos[0]->drmax*10;
1199 TCanvas* cWmean = new TCanvas(sscname.str().c_str(),
1200 sscname.str().c_str(),
1201 1000,350);
1202 cWmean->SetTopMargin(0.09);
1203 cWmean->SetLeftMargin(0.1);
1204 cWmean->SetRightMargin(0.15);
1205
1206 TMultiGraph* mg = new TMultiGraph();
1207 TLegend* leg = new TLegend(0.87,0.17,0.98,0.91);
1208 leg->SetFillColor(10);
1209
1210 TGraphErrors* ggen = Wmass_mean(vec_histos,1,channel);
1211 TGraphErrors* grec = Wmass_mean(vec_histos,2,channel);
1212 TGraphErrors* gcorr= Wmass_mean(vec_histos,3,channel);
1213 TGraphErrors* gl5 = Wmass_mean(vec_histos,4,channel);
1214
1215 mg->Add(ggen);
1216 mg->Add(grec);
1217 mg->Add(gcorr);
1218 mg->Add(gl5);
1219
1220 ggen ->SetLineWidth(2);
1221 grec ->SetLineWidth(2);
1222 gcorr->SetLineWidth(2);
1223 gl5 ->SetLineWidth(2);
1224
1225 ggen->SetLineColor(kGreen+1);
1226 ggen->SetMarkerColor(kGreen+1);
1227 grec->SetLineColor(kBlue+2);
1228 grec->SetMarkerColor(kBlue+2);
1229 gcorr->SetLineColor(kMagenta);
1230 gcorr->SetMarkerColor(kMagenta);
1231 gl5->SetLineColor(kRed);
1232 gl5->SetMarkerColor(kRed);
1233
1234 ggen->SetMarkerStyle(kFullStar);
1235 ggen->SetMarkerSize(1.5);
1236 grec->SetMarkerStyle(kFullSquare);
1237 gcorr->SetMarkerStyle(kOpenCircle);
1238 gl5->SetMarkerStyle(kFullCircle);
1239
1240 leg->AddEntry(gcorr,"CORR","lp");
1241 leg->AddEntry(gl5, "L5", "lp");
1242 leg->AddEntry(ggen, "GEN", "lp");
1243 leg->AddEntry(grec, "REC", "lp");
1244
1245 mg->Draw("AP");
1246 mg->GetHistogram()->SetTitle("W mass (mean)");
1247 mg->GetHistogram()->SetYTitle("m_{W} [GeV]");
1248
1249 for (unsigned int i=0;i<vec_histos.size();i++) {
1250 int ibin=0; while ((i+1)>mg->GetHistogram()->GetBinLowEdge(ibin)) { ibin++;}
1251 mg->GetHistogram()->GetXaxis()
1252 ->SetBinLabel(ibin,alg_label(vec_histos[i]->algname).c_str());
1253 }
1254 mg->GetHistogram()->GetXaxis()->LabelsOption("u");
1255 mg->GetHistogram()->GetXaxis()->SetLabelOffset(0.01);
1256 mg->GetHistogram()->GetXaxis()->SetLabelSize(0.065);
1257 mg->GetHistogram()->GetYaxis()->SetTitleOffset(0.75);
1258 leg->Draw();
1259
1260 TLatex mWtxt;
1261 mWtxt.DrawLatex(0.5,81.0,"m_{W}=80.42 GeV");
1262
1263 TLatex pttxt;
1264 stringstream sspt; sspt<<"jet p_{T}#geq "<<vec_histos[0]->ptmin<<" GeV";
1265 pttxt.SetNDC();
1266 pttxt.DrawLatex(0.7,0.92,sspt.str().c_str());
1267
1268 TAxis* xaxis = mg->GetHistogram()->GetXaxis();
1269 TLine* mWline=new TLine(xaxis->GetXmin(),80.42,xaxis->GetXmax(),80.42);
1270 mWline->SetLineWidth(2);
1271 mWline->Draw("SAME");
1272
1273 cWmean->Modified();
1274 }
1275
1276
1277 //______________________________________________________________________________
1278 void plot_topmass_mean(std::vector<histograms*> vec_histos,
1279 const std::string& channel)
1280 {
1281 stringstream sscname;
1282 string ch = (channel.empty()) ? channel : "_"+channel;
1283 sscname<<"topmassMean"<<ch<<"_pt"<<vec_histos[0]->ptmin
1284 <<"_dr0"<<vec_histos[0]->drmax*10;
1285 TCanvas* ctopmean = new TCanvas(sscname.str().c_str(),
1286 sscname.str().c_str(),
1287 1000,350);
1288
1289 ctopmean->SetTopMargin(0.09);
1290 ctopmean->SetLeftMargin(0.1);
1291 ctopmean->SetRightMargin(0.15);
1292
1293 TMultiGraph* mg = new TMultiGraph();
1294 TLegend* leg = new TLegend(0.87,0.17,0.98,0.91);
1295 leg->SetFillColor(10);
1296
1297 TGraphErrors* ggen = topmass_mean(vec_histos,1,channel);
1298 TGraphErrors* grec = topmass_mean(vec_histos,2,channel);
1299 TGraphErrors* gcorr= topmass_mean(vec_histos,3,channel);
1300 TGraphErrors* gl5 = topmass_mean(vec_histos,4,channel);
1301
1302 mg->Add(ggen);
1303 mg->Add(grec);
1304 mg->Add(gcorr);
1305 mg->Add(gl5);
1306
1307 ggen ->SetLineWidth(2);
1308 grec ->SetLineWidth(2);
1309 gcorr->SetLineWidth(2);
1310 gl5 ->SetLineWidth(2);
1311
1312 ggen->SetLineColor(kGreen+1);
1313 ggen->SetMarkerColor(kGreen+1);
1314 grec->SetLineColor(kBlue+2);
1315 grec->SetMarkerColor(kBlue+2);
1316 gcorr->SetLineColor(kMagenta);
1317 gcorr->SetMarkerColor(kMagenta);
1318 gl5->SetLineColor(kRed);
1319 gl5->SetMarkerColor(kRed);
1320
1321 ggen->SetMarkerStyle(kFullStar);
1322 ggen->SetMarkerSize(1.5);
1323 grec->SetMarkerStyle(kFullSquare);
1324 gcorr->SetMarkerStyle(kOpenCircle);
1325 gl5->SetMarkerStyle(kFullCircle);
1326
1327 leg->AddEntry(gcorr,"CORR","lp");
1328 leg->AddEntry(gl5, "L5", "lp");
1329 leg->AddEntry(ggen, "GEN", "lp");
1330 leg->AddEntry(grec, "REC", "lp");
1331
1332 mg->Draw("AP");
1333 mg->GetHistogram()->SetTitle("top mass (mean)");
1334 mg->GetHistogram()->SetYTitle("m_{top} [GeV]");
1335
1336 for (unsigned int i=0;i<vec_histos.size();i++) {
1337 int ibin=0; while ((i+1)>mg->GetHistogram()->GetBinLowEdge(ibin)) { ibin++;}
1338 mg->GetHistogram()->GetXaxis()
1339 ->SetBinLabel(ibin,alg_label(vec_histos[i]->algname).c_str());
1340 }
1341 mg->GetHistogram()->GetXaxis()->LabelsOption("u");
1342 mg->GetHistogram()->GetXaxis()->SetLabelOffset(0.01);
1343 mg->GetHistogram()->GetXaxis()->SetLabelSize(0.065);
1344 mg->GetHistogram()->GetYaxis()->SetTitleOffset(0.75);
1345 leg->Draw();
1346
1347 TLatex mtoptxt;
1348 mtoptxt.DrawLatex(0.5,178.6,"m_{top}=175 GeV");
1349
1350 TLatex pttxt;
1351 stringstream sspt; sspt<<"jet p_{T}#geq "<<vec_histos[0]->ptmin<<" GeV";
1352 pttxt.SetNDC();
1353 pttxt.DrawLatex(0.7,0.92,sspt.str().c_str());
1354
1355 TAxis* xaxis = mg->GetHistogram()->GetXaxis();
1356 TLine* mtopline = new TLine(xaxis->GetXmin(),175.,xaxis->GetXmax(),175.);
1357 mtopline->SetLineWidth(2);
1358 mtopline->Draw("SAME");
1359
1360 ctopmean->Modified();
1361 }
1362
1363
1364 //______________________________________________________________________________
1365 void plot_Wmass_width(std::vector<histograms*> vec_histos,
1366 const std::string& channel)
1367 {
1368 stringstream sscname;
1369 string ch = (channel.empty()) ? channel : "_"+channel;
1370 sscname<<"WmassWidth"<<ch<<"_pt"<<vec_histos[0]->ptmin
1371 <<"_dr0"<<vec_histos[0]->drmax*10;
1372 TCanvas* cWwidth = new TCanvas(sscname.str().c_str(),
1373 sscname.str().c_str(),
1374 1000,350);
1375
1376 cWwidth->SetTopMargin(0.09);
1377 cWwidth->SetLeftMargin(0.1);
1378 cWwidth->SetRightMargin(0.15);
1379
1380 TMultiGraph* mg = new TMultiGraph();
1381 TLegend* leg = new TLegend(0.87,0.17,0.98,0.91);
1382 leg->SetFillColor(10);
1383
1384 TGraphErrors* ggen = Wmass_width(vec_histos,1,channel);
1385 TGraphErrors* grec = Wmass_width(vec_histos,2,channel);
1386 TGraphErrors* gcorr= Wmass_width(vec_histos,3,channel);
1387 TGraphErrors* gl5 = Wmass_width(vec_histos,4,channel);
1388
1389 mg->Add(ggen);
1390 mg->Add(grec);
1391 mg->Add(gcorr);
1392 mg->Add(gl5);
1393
1394 ggen->SetLineColor(kGreen+1);
1395 ggen->SetMarkerColor(kGreen+1);
1396 grec->SetLineColor(kBlue+2);
1397 grec->SetMarkerColor(kBlue+2);
1398 gcorr->SetLineColor(kMagenta);
1399 gcorr->SetMarkerColor(kMagenta);
1400 gl5->SetLineColor(kRed);
1401 gl5->SetMarkerColor(kRed);
1402
1403 ggen->SetMarkerStyle(kFullStar);
1404 ggen->SetMarkerSize(1.5);
1405 grec->SetMarkerStyle(kFullSquare);
1406 gcorr->SetMarkerStyle(kOpenCircle);
1407 gl5->SetMarkerStyle(kFullCircle);
1408
1409 leg->AddEntry(gcorr,"CORR","lp");
1410 leg->AddEntry(gl5, "L5", "lp");
1411 leg->AddEntry(ggen, "GEN", "lp");
1412 leg->AddEntry(grec, "REC", "lp");
1413
1414 mg->Draw("AP");
1415 mg->GetHistogram()->SetTitle("W mass width (RMS)");
1416 mg->GetHistogram()->SetYTitle("RMS(m_{W}) [GeV]");
1417
1418 for (unsigned int i=0;i<vec_histos.size();i++) {
1419 int ibin=0; while ((i+1)>mg->GetHistogram()->GetBinLowEdge(ibin)) { ibin++;}
1420 mg->GetHistogram()->GetXaxis()
1421 ->SetBinLabel(ibin,alg_label(vec_histos[i]->algname).c_str());
1422 }
1423 mg->GetHistogram()->GetXaxis()->LabelsOption("u");
1424 mg->GetHistogram()->GetXaxis()->SetLabelOffset(0.01);
1425 mg->GetHistogram()->GetXaxis()->SetLabelSize(0.065);
1426 mg->GetHistogram()->GetYaxis()->SetTitleOffset(0.75);
1427 leg->Draw();
1428
1429 TLatex pttxt;
1430 stringstream sspt; sspt<<"jet p_{T}#geq "<<vec_histos[0]->ptmin<<" GeV";
1431 pttxt.SetNDC();
1432 pttxt.DrawLatex(0.7,0.92,sspt.str().c_str());
1433
1434 cWwidth->Modified();
1435 }
1436
1437
1438 //______________________________________________________________________________
1439 void plot_topmass_width(std::vector<histograms*> vec_histos,
1440 const std::string& channel)
1441 {
1442 stringstream sscname;
1443 string ch = (channel.empty()) ? channel : "_"+channel;
1444 sscname<<"topmassWidth"<<ch<<"_pt"<<vec_histos[0]->ptmin
1445 <<"_dr0"<<vec_histos[0]->drmax*10;
1446 TCanvas* ctopwidth = new TCanvas(sscname.str().c_str(),
1447 sscname.str().c_str(),
1448 1000,350);
1449
1450 ctopwidth->SetTopMargin(0.09);
1451 ctopwidth->SetLeftMargin(0.1);
1452 ctopwidth->SetRightMargin(0.15);
1453
1454 TMultiGraph* mg = new TMultiGraph();
1455 TLegend* leg = new TLegend(0.87,0.17,0.98,0.91);
1456 leg->SetFillColor(10);
1457
1458 TGraphErrors* ggen = topmass_width(vec_histos,1,channel);
1459 TGraphErrors* grec = topmass_width(vec_histos,2,channel);
1460 TGraphErrors* gcorr= topmass_width(vec_histos,3,channel);
1461 TGraphErrors* gl5 = topmass_width(vec_histos,4,channel);
1462
1463 mg->Add(ggen);
1464 mg->Add(grec);
1465 mg->Add(gcorr);
1466 mg->Add(gl5);
1467
1468 ggen->SetLineColor(kGreen+1);
1469 ggen->SetMarkerColor(kGreen+1);
1470 grec->SetLineColor(kBlue+2);
1471 grec->SetMarkerColor(kBlue+2);
1472 gcorr->SetLineColor(kMagenta);
1473 gcorr->SetMarkerColor(kMagenta);
1474 gl5->SetLineColor(kRed);
1475 gl5->SetMarkerColor(kRed);
1476
1477 ggen->SetMarkerStyle(kFullStar);
1478 ggen->SetMarkerSize(1.5);
1479 grec->SetMarkerStyle(kFullSquare);
1480 gcorr->SetMarkerStyle(kOpenCircle);
1481 gl5->SetMarkerStyle(kFullCircle);
1482
1483 leg->AddEntry(gcorr,"CORR","lp");
1484 leg->AddEntry(gl5, "L5", "lp");
1485 leg->AddEntry(ggen, "GEN", "lp");
1486 leg->AddEntry(grec, "REC", "lp");
1487
1488 mg->Draw("AP");
1489 mg->GetHistogram()->SetTitle("top mass width (RMS)");
1490 mg->GetHistogram()->SetYTitle("RMS(m_{top}) [GeV]");
1491
1492 for (unsigned int i=0;i<vec_histos.size();i++) {
1493 int ibin=0; while ((i+1)>mg->GetHistogram()->GetBinLowEdge(ibin)) { ibin++;}
1494 mg->GetHistogram()->GetXaxis()
1495 ->SetBinLabel(ibin,alg_label(vec_histos[i]->algname).c_str());
1496 }
1497 mg->GetHistogram()->GetXaxis()->LabelsOption("u");
1498 mg->GetHistogram()->GetXaxis()->SetLabelOffset(0.01);
1499 mg->GetHistogram()->GetXaxis()->SetLabelSize(0.065);
1500 mg->GetHistogram()->GetYaxis()->SetTitleOffset(0.75);
1501 leg->Draw();
1502
1503 TLatex pttxt;
1504 stringstream sspt; sspt<<"jet p_{T}#geq "<<vec_histos[0]->ptmin<<" GeV";
1505 pttxt.SetNDC();
1506 pttxt.DrawLatex(0.7,0.92,sspt.str().c_str());
1507
1508 ctopwidth->Modified();
1509 }
1510
1511
1512 //______________________________________________________________________________
1513 TGraphErrors* Wmass_mean(std::vector<histograms*> vec_histos,int level,
1514 const std::string& channel)
1515 {
1516 int n = (int)vec_histos.size();
1517 double *x = new double[n];
1518 double *y = new double[n];
1519 double *ex = new double[n];
1520 double *ey = new double[n];
1521 for (int i=0;i<n;i++) {
1522 x[i] = i+1;
1523 ex[i] = 0.0;
1524 if (channel.empty()) {
1525 y[i] =vec_histos[i]->hmw[level]->GetMean();
1526 ey[i]=vec_histos[i]->hmw[level]->GetMeanError();
1527 }
1528 else if (channel=="lj") {
1529 y[i] =vec_histos[i]->hmw_lj[level]->GetMean();
1530 ey[i]=vec_histos[i]->hmw_lj[level]->GetMeanError();
1531 }
1532 else if (channel=="aj") {
1533 y[i] =vec_histos[i]->hmw_aj[level]->GetMean();
1534 ey[i]=vec_histos[i]->hmw_aj[level]->GetMeanError();
1535 }
1536 else assert(0);
1537 }
1538
1539 return new TGraphErrors(n,x,y,ex,ey);
1540 }
1541
1542
1543 //______________________________________________________________________________
1544 TGraphErrors* topmass_mean(std::vector<histograms*> vec_histos,int level,
1545 const std::string& channel)
1546 {
1547 int n = (int)vec_histos.size();
1548 double *x = new double[n];
1549 double *y = new double[n];
1550 double *ex = new double[n];
1551 double *ey = new double[n];
1552 for (int i=0;i<n;i++) {
1553 x[i] = i+1;
1554 ex[i] = 0.0;
1555 if (channel.empty()) {
1556 y[i] =vec_histos[i]->hmt[level]->GetMean();
1557 ey[i]=vec_histos[i]->hmt[level]->GetMeanError();
1558 }
1559 else if (channel=="lj") {
1560 y[i] =vec_histos[i]->hmt_lj[level]->GetMean();
1561 ey[i]=vec_histos[i]->hmt_lj[level]->GetMeanError();
1562 }
1563 else if (channel=="aj") {
1564 y[i] =vec_histos[i]->hmt_aj[level]->GetMean();
1565 ey[i]=vec_histos[i]->hmt_aj[level]->GetMeanError();
1566 }
1567 else assert(0);
1568 }
1569
1570 return new TGraphErrors(n,x,y,ex,ey);
1571 }
1572
1573
1574 //______________________________________________________________________________
1575 TGraphErrors* Wmass_width(std::vector<histograms*> vec_histos,int level,
1576 const std::string& channel)
1577 {
1578 int n = (int)vec_histos.size();
1579 double *x = new double[n];
1580 double *y = new double[n];
1581 double *ex = new double[n];
1582 double *ey = new double[n];
1583 for (int i=0;i<n;i++) {
1584 x[i] = i+1;
1585 ex[i] = 0.0;
1586 if (channel.empty()) {
1587 y[i] =vec_histos[i]->hmw[level]->GetRMS();
1588 ey[i]=vec_histos[i]->hmw[level]->GetRMSError();
1589 }
1590 else if (channel=="lj") {
1591 y[i] =vec_histos[i]->hmw_lj[level]->GetRMS();
1592 ey[i]=vec_histos[i]->hmw_lj[level]->GetRMSError();
1593 }
1594 else if (channel=="aj") {
1595 y[i] =vec_histos[i]->hmw_aj[level]->GetRMS();
1596 ey[i]=vec_histos[i]->hmw_aj[level]->GetRMSError();
1597 }
1598 else assert(0);
1599 }
1600
1601 return new TGraphErrors(n,x,y,ex,ey);
1602 }
1603
1604
1605 //______________________________________________________________________________
1606 TGraphErrors* topmass_width(std::vector<histograms*> vec_histos,int level,
1607 const std::string& channel)
1608 {
1609 int n = (int)vec_histos.size();
1610 double *x = new double[n];
1611 double *y = new double[n];
1612 double *ex = new double[n];
1613 double *ey = new double[n];
1614 for (int i=0;i<n;i++) {
1615 x[i] = i+1;
1616 ex[i] = 0.0;
1617 if (channel.empty()) {
1618 y[i] =vec_histos[i]->hmt[level]->GetRMS();
1619 ey[i]=vec_histos[i]->hmt[level]->GetRMSError();
1620 }
1621 else if (channel=="lj") {
1622 y[i] =vec_histos[i]->hmt_lj[level]->GetRMS();
1623 ey[i]=vec_histos[i]->hmt_lj[level]->GetRMSError();
1624 }
1625 else if (channel=="aj") {
1626 y[i] =vec_histos[i]->hmt_aj[level]->GetRMS();
1627 ey[i]=vec_histos[i]->hmt_aj[level]->GetRMSError();
1628 }
1629 else assert(0);
1630 }
1631
1632 return new TGraphErrors(n,x,y,ex,ey);
1633 }
1634
1635
1636 //______________________________________________________________________________
1637 bool parse_filename(const std::string& filename,float& ptmin,float& drmax)
1638 {
1639 unsigned int pos;
1640 string tmp=filename;
1641 while ((pos=tmp.find("/"))!=string::npos) tmp=tmp.substr(pos+1);
1642 pos = tmp.find("topmass_pt");
1643 if (pos==string::npos) {
1644 cout<<"ERROR parsing "<<filename<<endl;
1645 return false;
1646 }
1647 tmp = tmp.substr(pos+10);
1648
1649 pos = tmp.find("_dr");
1650 if (pos==string::npos) {
1651 cout<<"ERROR parsing "<<filename<<endl;
1652 return false;
1653 }
1654 stringstream sspt; sspt<<tmp.substr(0,pos);tmp=tmp.substr(pos+3);
1655 sspt>>ptmin;
1656 cout<<"ptmin="<<ptmin<<endl;
1657
1658 pos = tmp.find(".root");
1659 if (pos==string::npos) {
1660 cout<<"ERROR parsing "<<filename<<endl;
1661 return false;
1662 }
1663 stringstream ssdr; ssdr<<tmp.substr(0,pos);
1664 ssdr>>drmax;
1665 cout<<"drmax="<<drmax<<endl;
1666
1667 return true;
1668 }
1669
1670
1671 //______________________________________________________________________________
1672 string alg_label(const string& alg_name)
1673 {
1674 string alg=alg_name.substr(0,2);
1675 string par=alg_name.substr(2,alg_name.find('.'));
1676 if (alg=="kt") return "Fast k_{T}, D=0."+par;
1677 else if (alg=="sc") return "SISCone, R=0."+par;
1678 else if (alg=="mc") return "MidCone, R=0."+par;
1679 else if (alg=="ic") return "ItCone, R=0."+par;
1680 else return "UNKNOWN";
1681
1682 }
1683
1684
1685 //______________________________________________________________________________
1686 Color_t alg_color(const string& alg_name)
1687 {
1688 if (alg_name=="kt3") return kBlue-9;
1689 if (alg_name=="kt4") return kBlue;
1690 if (alg_name=="kt6") return kBlue+3;
1691 if (alg_name=="sc4") return kMagenta-9;
1692 if (alg_name=="sc5") return kMagenta;
1693 if (alg_name=="sc7") return kMagenta+4;
1694 if (alg_name=="mc4") return kCyan-9;
1695 if (alg_name=="mc5") return kCyan;
1696 if (alg_name=="mc7") return kCyan+4;
1697 if (alg_name=="ic4") return kGreen-9;
1698 if (alg_name=="ic5") return kGreen;
1699 if (alg_name=="ic7") return kGreen+4;
1700 return kBlack;
1701 }
1702
1703
1704
1705 //______________________________________________________________________________
1706 Style_t alg_style(const string& alg_name)
1707 {
1708 if (alg_name=="kt3") return kFullTriangleDown;
1709 if (alg_name=="kt4") return kFullCircle;
1710 if (alg_name=="kt6") return kFullTriangleUp;
1711 if (alg_name=="sc4") return kFullTriangleDown;
1712 if (alg_name=="sc5") return kFullCircle;
1713 if (alg_name=="sc7") return kFullTriangleUp;
1714 if (alg_name=="mc4") return kFullTriangleDown;
1715 if (alg_name=="mc5") return kFullCircle;
1716 if (alg_name=="mc7") return kFullTriangleUp;
1717 if (alg_name=="ic4") return kFullTriangleDown;
1718 if (alg_name=="ic5") return kFullCircle;
1719 if (alg_name=="ic7") return kFullTriangleUp;
1720 return kOpenCircle;
1721 }