ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/cbrown/AnalysisFramework/miscFR/JZBeff.C
Revision: 1.1
Committed: Fri Jun 24 07:44:49 2011 UTC (13 years, 10 months ago) by fronga
Content type: text/plain
Branch: MAIN
CVS Tags: cbaf_4_98ifb_paper, cbaf_4p7ifb, Honeypot, cbaf_2p1ifb, HEAD
Log Message:
First commit of various scripts used in the first 2011 analyses.

File Contents

# User Rev Content
1 fronga 1.1 //
2     // Test of the JZB selection efficiency
3     //
4    
5     #include <iostream>
6     #include <iomanip>
7    
8     #include "TTree.h"
9     #include "TFile.h"
10     #include "TH1F.h"
11     #include "TCut.h"
12     #include "TLine.h"
13     #include "TMath.h"
14     #include "TCanvas.h"
15     #include "TString.h"
16     #include "TLegend.h"
17     #include "TProfile.h"
18     #include "TF1.h"
19    
20     Int_t nBins = 100;
21     Float_t jzbMin = -207;
22     Float_t jzbMax = 243;
23     Float_t jzbSel = 100;
24    
25    
26     TFile* getFile(TString name) {
27     // TString dcPath("dcap://t3se01.psi.ch:22125/pnfs/psi.ch/cms/trivcat/store/user/buchmann/");
28     TString dcPath("/scratch/fronga/MC_v1.29/");
29     TFile* file = TFile::Open(dcPath+name+".root");
30     if ( file->IsOpen() ) return file;
31    
32     exit(-1);
33    
34     }
35    
36     //____________________________________________________________________________________
37     // Efficiency plot
38     TH1F* plotEff(TTree* events, TCut kbase, int count) {
39    
40     // Define new histogram
41     char hname[30]; sprintf(hname,"hJzbEff%d",count);
42     TH1F* hJzbEff = new TH1F(hname,"JZB selection efficiency ; JZB (GeV/c); Efficiency",
43     nBins,jzbMin,jzbMax);
44     Float_t step = (jzbMax-jzbMin)/static_cast<Float_t>(nBins);
45    
46     events->Draw("jzb[1]","genJZBSel>-400"&&kbase,"goff");
47     Float_t maxEff = events->GetSelectedRows();
48     std::cout << hname << " " << maxEff << std::endl;
49    
50     std::cout << "JZB max = " << jzbMax << std::endl;
51     // Loop over steps to get efficiency curve
52     char cut[256];
53     for ( Int_t iBin = 0; iBin<nBins; ++iBin ) {
54     sprintf(cut,"genJZBSel>%3f",jzbMin+iBin*step);
55     events->Draw("jzb[1]",TCut(cut)&&kbase,"goff");
56     Float_t eff = static_cast<Float_t>(events->GetSelectedRows())/maxEff;
57     // std::cout << "COUCOU " << __LINE__ << std::endl;
58     hJzbEff->SetBinContent(iBin+1,eff);
59     hJzbEff->SetBinError(iBin+1,TMath::Sqrt(eff*(1-eff)/maxEff));
60     }
61     return hJzbEff;
62    
63    
64     }
65    
66     //____________________________________________________________________________________
67     // MC expectation
68     int expected(Float_t jzbCut = 100.) {
69    
70     TFile* fLM4 = getFile("LM4_SUSY_sftsht_7TeV-pythia6");
71     TFile* fLM8 = getFile("LM8_SUSY_sftsht_7TeV-pythia6");
72     TFile* fTTJ = getFile("TTJets_TuneZ2_7TeV-madgraph-tauola");
73     TFile* fZLL = getFile("DYJetsToLL_TuneZ2_M-50_7TeV-madgraph-tauola");
74    
75     // TCanvas* mycanvas = new TCanvas("mycanvas","canvas",600,500);
76     // mycanvas->Draw();
77    
78     // Acceptance cuts
79     char cut[256];
80     sprintf(cut,"jzb[1]>%f",jzbCut);
81     TCut kbase("pfJetGoodNum>2&&abs(mll-91.2)<20&&id1==id2"&&TCut(cut));
82     std::cout << kbase.GetTitle() << std::endl;
83    
84     Float_t lumi = 191.; // [1/pb]
85     Float_t xs[] = { 2.54, 1.029, 157.5, 3048.0 };
86     // Float_t sys[] = { 5/91.9, 5/88.5 }; // 50 GeV
87     Float_t sys[] = { 6.6/90.4, 6.9/85.3 }; // 100 GeV
88     Float_t nevents[] = { 242322.0, 242322.0, 1161621.0, 2313911.0};
89    
90     Long_t nentries = ((TTree*)fLM4->Get("events"))->Draw("jzb[1]",kbase,"goff");
91     int isample = 0;
92     Float_t rescaling = lumi/nevents[isample]*xs[isample];
93     std::cout << "LM4: " << nentries << " " << nentries*rescaling << "+-" << sqrt(nentries)*rescaling << "+-" << nentries*sys[isample]*rescaling << std::endl;
94     std::cout << "LM4 (sigma x BR x acc): " << nentries/nevents[isample]*xs[isample] << std::endl;
95    
96     ++isample;
97     nentries = ((TTree*)fLM8->Get("events"))->Draw("jzb[1]",kbase,"goff");
98     rescaling = lumi/nevents[isample]*xs[isample];
99     std::cout << "LM8: " << nentries << " " << nentries*rescaling << "+-" << sqrt(nentries)*rescaling << "+-" << nentries*sys[isample]*rescaling << std::endl;
100     std::cout << "LM8 (sigma x BR x acc): " << nentries/nevents[isample]*xs[isample] << std::endl;
101    
102     // ++isample;
103     // nentries = ((TTree*)fTTJ->Get("events"))->Draw("jzb[1]",kbase,"goff");
104     // std::cout << "TTJ: " << nentries << " " << nentries*lumi/nevents[isample]*xs[isample] << std::endl;
105    
106     // ++isample;
107     // nentries = ((TTree*)fZLL->Get("events"))->Draw("jzb[1]",kbase,"goff");
108     // std::cout << "ZLL: " << nentries << " " << nentries*lumi/nevents[isample]*xs[isample] << std::endl;
109    
110     return 0;
111    
112     }
113    
114    
115     //____________________________________________________________________________________
116     // JZB selection efficiency
117     int JZBeff(void) {
118    
119     TFile* fLM4 = getFile("LM4_SUSY_sftsht_7TeV-pythia6");
120     TFile* fLM8 = getFile("LM8_SUSY_sftsht_7TeV-pythia6");
121     TFile* fTTJ = getFile("TTJets_TuneZ2_7TeV-madgraph-tauola");
122     TFile* fZLL = getFile("DYJetsToLL_TuneZ2_M-50_7TeV-madgraph-tauola");
123    
124     // TFile* fLM4 = TFile::Open("dcap://t3se01.psi.ch:22125/pnfs/psi.ch/cms/trivcat/store/user/buchmann/MCSpring2011PU/LM4_SUSY_sftsht_7TeV-pythia6.root");
125     // TFile* fLM8 = TFile::Open("dcap://t3se01.psi.ch:22125/pnfs/psi.ch/cms/trivcat/store/user/buchmann/MCSpring2011PU/LM8_SUSY_sftsht_7TeV-pythia6.root");
126     // TFile* fTTJ = TFile::Open("dcap://t3se01.psi.ch:22125/pnfs/psi.ch/cms/trivcat/store/user/buchmann/MCSpring2011PU/TTJets_TuneZ2_7TeV-madgraph-tauola.root");
127     // TFile* fZLL = TFile::Open("dcap://t3se01.psi.ch:22125/pnfs/psi.ch/cms/trivcat/store/user/buchmann/MCSpring2011PU/DYJetsToLL_TuneZ2_M-50_7TeV-madgraph-tauola.root");
128    
129     TCanvas* mycanvas = new TCanvas("mycanvas","canvas",600,500);
130     mycanvas->Draw();
131    
132     // Acceptance cuts
133     TCut kbase("abs(genMllSel-91.2)<20&&pfJetGoodNum>2&&genZPt>0&&abs(mll-91.2)<20&&((id1+1)*(id2+1)*ch1*ch2)!=-2");
134    
135     int iplot=0;
136     TH1F* hLM4 = plotEff((TTree*)fLM4->Get("events"),kbase,iplot++);
137     hLM4->SetTitle("LM4");
138     hLM4->SetMarkerColor(kRed);
139     hLM4->SetMarkerStyle(24);
140     hLM4->SetMinimum(0.);
141    
142     TH1F* hLM8 = plotEff((TTree*)fLM8->Get("events"),kbase,iplot++);
143     hLM8->SetTitle("LM8");
144     hLM8->SetMarkerColor(kRed);
145    
146     TH1F* hZLL = plotEff((TTree*)fZLL->Get("events"),kbase,iplot++);
147     hZLL->SetTitle("Z+jets");
148     hZLL->SetMarkerColor(kBlack);
149     hZLL->SetMarkerStyle(21);
150    
151     // Remove ZPt cut for TTJets
152     kbase = TCut("abs(genMllSel-91.2)<20&&pfJetGoodNum>2&&abs(mll-91.2)<20&&((id1+1)*(id2+1)*ch1*ch2)!=-2");
153     TH1F* hTTJ = plotEff((TTree*)fTTJ->Get("events"),kbase,iplot++);
154     hTTJ->SetTitle("ttbar");
155     hTTJ->SetMarkerColor(kGreen-3);
156     hTTJ->SetMarkerStyle(30);
157    
158     // Add line where we cut
159     TLine* lSelect = new TLine(jzbSel,0.,jzbSel,hLM4->GetMaximum()*1.05);
160     lSelect->SetLineWidth(2);
161     lSelect->SetLineColor(kBlue);
162    
163     // Dump some information
164     Int_t bin = hLM4->FindBin(jzbSel); // To get the error
165     std::cout << "Efficiency at JZB==" << jzbSel << std::endl;
166     std::cout << " LM4: " << hLM4->Interpolate(jzbSel) << "+-" << hLM4->GetBinError(bin) << std::endl;
167     std::cout << " LM8: " << hLM8->Interpolate(jzbSel) << "+-" << hLM8->GetBinError(bin) << std::endl;
168     std::cout << " TTJ: " << hTTJ->Interpolate(jzbSel) << "+-" << hTTJ->GetBinError(bin) << std::endl;
169     std::cout << " ZLL: " << hZLL->Interpolate(jzbSel) << "+-" << hZLL->GetBinError(bin) << std::endl;
170    
171     // Draw it all
172     hLM4->Draw("E1");
173     hLM8->Draw("sameE1");
174     hTTJ->Draw("sameE1");
175     hZLL->Draw("sameE1");
176     lSelect->Draw();
177    
178     // Add legend
179     TLegend* legend = new TLegend(0.7,0.7,0.9,0.9);
180     legend->AddEntry(hZLL,hZLL->GetTitle(),"p");
181     legend->AddEntry(hTTJ,hTTJ->GetTitle(),"p");
182     legend->AddEntry(hLM4,hLM4->GetTitle(),"p");
183     legend->AddEntry(hLM8,hLM8->GetTitle(),"p");
184     legend->SetFillColor(0);
185     legend->SetBorderSize(1);
186     legend->Draw();
187    
188     mycanvas->SaveAs("plots/JZBefficiency.eps");
189     mycanvas->SaveAs("plots/macros/JZBefficiency.C");
190    
191     return 0;
192    
193     }
194    
195    
196     //____________________________________________________________________________________
197     // Total selection efficiency (MC)
198     int MCeff(void) {
199    
200     char jzbSelStr[256]; sprintf(jzbSelStr,"%f",jzbSel);
201     // All acceptance cuts at gen. level
202     TCut kbase("abs(genMllSel-91.2)<20&&pfJetGoodNum>2&&genZPt>0&&genJZB>"+TString(jzbSelStr)+"&&genId1==-genId2");
203     // Corresponding reco. cuts
204     TCut ksel("abs(mll-91.2)<20&&id1==id2&&jzb[1]+0.9>"+TString(jzbSelStr));
205     TFile* fLM4 = getFile("LM4_SUSY_sftsht_7TeV-pythia6");
206     TFile* fLM8 = getFile("LM8_SUSY_sftsht_7TeV-pythia6");
207     TFile* fTTJ = getFile("TTJets_TuneZ2_7TeV-madgraph-tauola");
208     // TFile* fZLL = getFile("MCSpring2011PU/DYJetsToLL_TuneZ2_M-50_7TeV-madgraph-tauola");
209    
210     TTree* tLM4 = ((TTree*)fLM4->Get("events"));
211     tLM4->Draw("jzb[1]",kbase&&ksel,"goff");
212     Float_t selLM4 = tLM4->GetSelectedRows();
213     tLM4->Draw("jzb[1]",kbase,"goff");
214     Float_t totLM4 = tLM4->GetSelectedRows();
215    
216     TTree* tLM8 = ((TTree*)fLM8->Get("events"));
217     tLM8->Draw("jzb[1]",kbase&&ksel,"goff");
218     Float_t selLM8 = tLM8->GetSelectedRows();
219     tLM8->Draw("jzb[1]",kbase,"goff");
220     Float_t totLM8 = tLM8->GetSelectedRows();
221    
222     // TTree* tZLL = ((TTree*)fZLL->Get("events"));
223     // tZLL->Draw("jzb[1]",kbase&&ksel,"goff");
224     // Float_t selZLL = tZLL->GetSelectedRows();
225     // tZLL->Draw("jzb[1]",kbase,"goff");
226     // Float_t totZLL = tZLL->GetSelectedRows();
227    
228     // // TTBar: remove ZPt cut
229     // kbase = TCut("genZPtSel>0&&genPt1Sel>20&&abs(genEta1Sel)<2.4&&abs(genEta2Sel)<2.4&&genPt2Sel>20&&abs(genMllSel-91.2)<20&&pfJetGoodNum>1&&genId1==-genId2&&genJZBSel>50");
230     // TTree* tTTJ = ((TTree*)fTTJ->Get("events"));
231     // tTTJ->Draw("jzb[1]",kbase&&ksel,"goff");
232     // Float_t selTTJ = tTTJ->GetSelectedRows();
233     // tTTJ->Draw("jzb[1]",kbase,"goff");
234     // Float_t totTTJ = tTTJ->GetSelectedRows();
235    
236     std::cout << "LM4: " << selLM4/totLM4 << "+-" << TMath::Sqrt(selLM4/totLM4*(1-selLM4/totLM4)/totLM4) << std::endl;
237     std::cout << "LM8: " << selLM8/totLM8 << "+-" << TMath::Sqrt(selLM8/totLM8*(1-selLM8/totLM8)/totLM8) << std::endl;
238     // std::cout << "TTJ: " << selTTJ/totTTJ << "+-" << TMath::Sqrt(selTTJ/totTTJ*(1-selTTJ/totTTJ)/totTTJ) << std::endl;
239     // std::cout << "ZLL: " << selZLL/totZLL << "+-" << TMath::Sqrt(selZLL/totZLL*(1-selZLL/totZLL)/totZLL) << std::endl;
240    
241    
242    
243     return 0;
244    
245     }
246    
247     //________________________________________________________________________
248     // Effect of energy scale on efficiency
249     int JZBjetScale(TString fileName = "LM4_SUSY_sftsht_7TeV-pythia6", Float_t jzbSelection=-1, TString plotName = "" ) {
250     TFile* file = getFile(fileName);
251     TCut kbase("abs(genMllSel-91.2)<20&&genZPt>0");
252     TCut ksel("abs(mll-91.2)<20&&((id1+1)*(id2+1)*ch1*ch2)!=-2");
253     TCut nJets("pfJetGoodNum>2");
254     TCut nJetsP("pfJetGoodNum33>2");
255     TCut nJetsM("pfJetGoodNum27>2");
256    
257     if ( jzbSelection>0 ) jzbSel = jzbSelection;
258    
259     if ( !(plotName.Length()>1) ) plotName = fileName;
260    
261     TCanvas* mycanvas = new TCanvas("mycanvas","canvas",600,500);
262     mycanvas->Draw();
263    
264     int iplot=0;
265     nBins = 1; jzbMin = jzbSel*0.95; jzbMax = jzbSel*1.05;
266     TH1F* hist = plotEff((TTree*)file->Get("events"),(kbase&&ksel&&nJets),iplot++);
267     hist->SetTitle(plotName+" p_{T}>30 GeV/c");
268     hist->SetMarkerColor(kRed);
269     hist->SetMarkerStyle(24);
270    
271     TH1F* histp = plotEff((TTree*)file->Get("events"),(kbase&&ksel&&nJetsP),iplot++);
272     histp->SetTitle(plotName+" "+nJetsP.GetTitle());
273     histp->SetMarkerColor(kRed);
274     histp->SetMarkerStyle(21);
275    
276     TH1F* histm = plotEff((TTree*)file->Get("events"),(kbase&&ksel&&nJetsM),iplot++);
277     histm->SetTitle(plotName+" "+nJetsM.GetTitle());
278     histm->SetMarkerColor(kRed);
279     histm->SetMarkerStyle(22);
280    
281     // Add line where we cut
282     TLine* lSelect = new TLine(jzbSel,0.,jzbSel,hist->GetMaximum()*1.05);
283     lSelect->SetLineStyle(2);
284     lSelect->SetLineColor(kBlue);
285    
286     // Dump some information
287     Float_t eff = hist->Interpolate(jzbSel);
288     Float_t effp = histp->Interpolate(jzbSel);
289     Float_t effm = histm->Interpolate(jzbSel);
290     std::cout << "Efficiency at JZB==" << jzbSel << std::endl;
291     std::cout << plotName << ": " << eff << std::endl;
292     std::cout << plotName << "p: " << effp << " (" << (effp-eff)/eff*100. << "%)" << std::endl;
293     std::cout << plotName << "m: " << effm << " (" << (effm-eff)/eff*100. << "%)" << std::endl;
294    
295     // Draw it all
296     hist->Draw("E1");
297     histp->Draw("sameE1");
298     histm->Draw("sameE1");
299     lSelect->Draw();
300    
301     // // Add legend
302     // TLegend* legend = new TLegend(0.7,0.7,0.9,0.9);
303     // legend->AddEntry(hist, hist->GetTitle(),"p");
304     // legend->AddEntry(histp,histp->GetTitle(),"p");
305     // legend->AddEntry(histm,histm->GetTitle(),"p");
306     // legend->SetFillColor(0);
307     // legend->SetBorderSize(1);
308     // legend->Draw();
309    
310     return 0;
311    
312     }
313    
314     //________________________________________________________________________
315     // Effect of energy scale on JZB efficiency
316     int JZBscale(TString fileName = "LM4_SUSY_sftsht_7TeV-pythia6", Float_t jzbSelection=-1, TString plotName = "" ) {
317    
318     Float_t syst = 0.1; //0.05; //0.1;
319    
320     TFile* file = new TFile("/scratch/fronga/MC_v1.29/"+fileName+".root");
321     TCut kbase("abs(genMllSel-91.2)<20&&genZPt>0&&pfJetGoodNum>2");
322     TCut ksel("abs(mll-91.2)<20&&((id1+1)*(id2+1)*ch1*ch2)!=-2");
323    
324     TCanvas* mycanvas = new TCanvas("mycanvas","canvas",600,500);
325     mycanvas->Draw();
326    
327     if ( jzbSelection>0 ) jzbSel = jzbSelection;
328    
329     nBins = 50;
330     jzbMin = 0.5*jzbSel;
331     jzbMax = 2.0*jzbSel;
332    
333     int iplot=0;
334     TTree* events = (TTree*)file->Get("events");
335     TH1F* hist = plotEff(events,kbase&&ksel,iplot++);
336     hist->SetTitle(plotName+" p_{T}>30 GeV/c");
337     hist->SetLineColor(kRed);
338     hist->SetMinimum(0.);
339     hist->SetMaximum(1.);
340    
341     // Add lines where we cut
342     TLine* lSelect = new TLine(jzbSel,0.,jzbSel,hist->GetMaximum());
343     lSelect->SetLineColor(kBlue);
344    
345     // Dump some information
346     Float_t eff = hist->Interpolate(jzbSel);
347     Float_t effp = hist->Interpolate(jzbSel*(1.+syst));
348     Float_t effm = hist->Interpolate(jzbSel*(1.-syst));
349     std::cout << fileName << " eff at JZB==" << jzbSel << ": " << eff << std::endl;
350     std::cout << fileName << " eff at JZB==" << jzbSel*(1.+syst) << ": " << effp << " (" << ((effp-eff)/eff)*100. << "%)" << std::endl;
351     std::cout << fileName << " eff at JZB==" << jzbSel*(1.-syst) << ": " << effm << " (" << ((effm-eff)/eff)*100. << "%)" << std::endl;
352    
353     // Draw it all
354     hist->Draw("E1");
355     lSelect->Draw();
356     lSelect->SetLineStyle(2);
357     lSelect->DrawLine(jzbSel*(1+syst),0.,jzbSel*(1+syst),hist->GetMaximum());
358     lSelect->DrawLine(jzbSel*(1-syst),0.,jzbSel*(1-syst),hist->GetMaximum());
359    
360     mycanvas->SaveAs("plots/jzbScale"+fileName+".eps");
361     mycanvas->SaveAs("plots/macros/jzbScale"+fileName+".C");
362    
363     return 0;
364    
365     }
366    
367     //________________________________________________________________________
368     // JZB response (true/reco. vs. true)
369     int JZBresponse(TString fileName = "LM4_SUSY_sftsht_7TeV-pythia6", bool isMET = kFALSE, Float_t myJzbMax = 200., Int_t nPeriods = 9 ) {
370     TFile* file = getFile(fileName);
371    
372     TCanvas* mycanvas = new TCanvas("mycanvas","canvas",600,500);
373     mycanvas->Draw();
374    
375     jzbMin = 20;
376    
377     TCut kbase("abs(genMllSel-91.2)<20&&genZPtSel>0&&pfJetGoodNum>2");
378     TCut ksel("abs(mll-91.2)<20&&((id1+1)*(id2+1)*ch1*ch2)!=-2");
379    
380     TProfile* hJzbResp = new TProfile("hJzbResp","JZB response ; JZB true (GeV/c); JZB reco. / JZB true",
381     nPeriods, jzbMin, myJzbMax, "" );
382    
383     TTree* events = (TTree*)file->Get("events");
384     if (!isMET) events->Project("hJzbResp","jzb[1]/genJZBSel:genJZBSel",kbase&&ksel);
385     else events->Project("hJzbResp","met[4]/genMET:genMET",kbase&&ksel);
386    
387     hJzbResp->SetMaximum(1.2);
388     hJzbResp->SetMinimum(0.2);
389     hJzbResp->Fit("pol0");
390    
391     // Draw line at 1
392     TLine* lSelect = new TLine(hJzbResp->GetXaxis()->GetXmin(),1.,hJzbResp->GetXaxis()->GetXmax(),1.);
393     lSelect->SetLineColor(kBlue);
394     lSelect->SetLineStyle(2);
395     lSelect->Draw();
396    
397     TString name("JZB");
398     if ( isMET ) name = TString("MET");
399     mycanvas->SaveAs("plots/"+name+"response"+fileName+".eps");
400     mycanvas->SaveAs("plots/macros/"+name+"response"+fileName+".C");
401    
402    
403     return 0;
404    
405     }
406    
407     //________________________________________________________________________________________
408     // Pile-up efficiency
409     int pileup(TString fileName = "LM4_SUSY_sftsht_7TeV-pythia6", Float_t myJzbMax = 140. ) {
410     TFile* file = getFile(fileName);
411    
412     nBins = 16;
413     jzbMax = myJzbMax;
414     TCanvas* mycanvas = new TCanvas("mycanvas","canvas",600,500);
415     mycanvas->Draw();
416    
417     // Acceptance cuts
418     TCut kbase("abs(genMllSel-91.2)<20&&pfJetGoodNum>2&&genZPtSel>0&&abs(mll-91.2)<20&&((id1+1)*(id2+1)*ch1*ch2)!=-2");
419     int iplot=0;
420    
421     TH1F* hLM4 = plotEff((TTree*)file->Get("events"),kbase,iplot++);
422     hLM4->SetTitle(fileName);
423     hLM4->SetMarkerColor(kRed);
424     hLM4->SetMarkerStyle(24);
425     hLM4->SetMinimum(0.);
426    
427     // Nominal function
428     TF1* func = new TF1("func","0.5*TMath::Erfc([0]*x-[1])",jzbMin,jzbMax);
429     func->SetParameter(0,0.03);
430     func->SetParameter(1,0.);
431     hLM4->Fit(func);
432    
433     // Pimped-up function
434     TF1* funcUp = (TF1*)func->Clone();
435     funcUp->SetParameter( 0., func->GetParameter(0)/1.1); // 10% systematic error (up in sigma => 0.1 in erfc)
436     funcUp->SetLineColor(kRed);
437     funcUp->SetLineStyle(2);
438     funcUp->Draw("same");
439    
440     std::cout << funcUp->Eval(jzbSel) << " " << func->Eval(jzbSel)
441     << "(" << (funcUp->Eval(jzbSel)-func->Eval(jzbSel))/func->Eval(jzbSel)*100. << "%)" << std::endl;
442    
443     TLine* lSelect = new TLine(jzbSel,0.,jzbSel,hLM4->GetMaximum()*1.05);
444     lSelect->SetLineWidth(2);
445     lSelect->SetLineColor(kBlue);
446     lSelect->Draw();
447    
448    
449     mycanvas->SaveAs("plots/pileup"+fileName+".eps");
450     mycanvas->SaveAs("plots/macros/pileup"+fileName+".C");
451    
452    
453     return 0;
454    
455     }
456    
457    
458     //________________________________________________________________________________________
459     // Run all efficiency tests
460     int doAll(void) {
461    
462     TString samples[] = { "LM4_SUSY_sftsht_7TeV-pythia6",
463     "LM8_SUSY_sftsht_7TeV-pythia6",
464     "TTJets_TuneZ2_7TeV-madgraph-tauola" };
465     size_t nSamples = 3;
466    
467     std::cout << "MC efficiencies: " << std::endl;
468     MCeff();
469    
470     std::cout << "JZB efficiency: " << std::endl;
471     JZBeff();
472    
473     std::cout << "Jet energy scale: " << std::endl;
474     for ( size_t i=0; i<nSamples; ++i ) JZBjetScale( samples[i] );
475    
476     std::cout << "JZB scale: " << std::endl;
477     for ( size_t i=0; i<nSamples; ++i ) JZBscale( samples[i] );
478    
479     std::cout << "JZB response: " << std::endl;
480     for ( size_t i=0; i<nSamples; ++i ) JZBresponse( samples[i] );
481    
482     std::cout << "Pileup: " << std::endl;
483     for ( size_t i=0; i<nSamples; ++i ) pileup( samples[i] );
484    
485    
486     return 0;
487    
488     }
489    
490     int calcEff(void) {
491    
492     float sysEff[][5] = { { 0.02, 0.02, 0.02 }, // Trigger eff. systematic
493     { 0.02, 0.02, 0.02 }, // Lepton eff. systematic
494     { 0.002, 0.012, 0.012 }, // Jet energy scale
495     { 0.08, 0.011, 0.014 }, // JZB systematic
496     { 0.39, 0.054, 0.069 } }; // pileup systematic
497    
498     float eff[] = { 0.389, 0.933, 0.913 };
499     float effErr[] = { 1.1, 0.6, 0.7 };
500    
501     for ( int i=0; i<3; ++i ) {
502     float totSys = 0;
503     for ( int j=0; j<5; ++j )
504     totSys += sysEff[j][i]*sysEff[j][i]*eff[i]*eff[i];
505     std::cout << "$ " << std::setprecision(3) << eff[i]*100.
506     << " \\pm " << effErr[i] << "(\\text{stat})\\pm "
507     << std::setprecision(2) << sqrt(totSys)*100. << "(\\text{syst}) $" << std::endl;
508     }
509    
510     return 0;
511    
512     }
513