ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/cbrown/Development/Plotting/Modules/ZbTools.C
Revision: 1.13
Committed: Fri Nov 23 12:10:58 2012 UTC (12 years, 5 months ago) by buchmann
Content type: text/plain
Branch: MAIN
Changes since 1.12: +163 -38 lines
Log Message: 
Added an efficiency cross-check; added information to the fit plots; fixed facevalue plot problem (shift); replaced fitting quietly with minos with a more verbose variation; fixed 30/10 label (to 30/15); added encouraging note to PU study plots; adapted the plotting style (Rabs results are red (xs), MPF are blue (R)); extended alpha range for fitting up to 0.35 for exclusive bins to have additional information for 30<pt<50; and some more stuff ...

File Contents

# User Rev Content
1 buchmann 1.1 #include <iostream>
2     #include <vector>
3     #include <sys/stat.h>
4    
5     #include <TCut.h>
6     #include <TROOT.h>
7     #include <TCanvas.h>
8     #include <TMath.h>
9     #include <TColor.h>
10     #include <TPaveText.h>
11     #include <TRandom.h>
12     #include <TH1.h>
13     #include <TH2.h>
14     #include <TF1.h>
15     #include <TSQLResult.h>
16     #include <TProfile.h>
17     #include <TLegendEntry.h>
18    
19     using namespace std;
20    
21     using namespace PlottingSetup;
22    
23     namespace ZbData {
24     vector<float> data_over_mc;
25    
26 buchmann 1.10 TCut ZplusBsel("pt1>20&&pt2>20&&mll>60&&mll<120&&id1==id2");
27 buchmann 1.9 TCut LeadingB("Zb3010_bTagProbCSVBP[0]>0.898");
28     TCut EtaB("abs(Zb3010_pfJetGoodEta[0])<1.3");
29 buchmann 1.10 TCut PhiZcut("abs(Zb3010_pfJetDphiZ[0])>2.7");
30 buchmann 1.1
31     }
32    
33     using namespace ZbData;
34    
35 buchmann 1.7
36     //*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/* DELETE EVERYTHING BETWEEN THIS LINE /*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*//
37     void SpecialCutFlow(string identifier);
38    
39     //*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/* AND THIS LINE /*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*//
40    
41    
42 buchmann 1.1 void print_yield(TCut cut) {
43     TH1F *data = allsamples.Draw("data", "mll",1,0,10000000, "m_{ll} [GeV]", "events", cut,data,luminosity);
44 buchmann 1.7 dout << " Yields for " << (const char*) cut << endl;
45     dout << " data: " << data->Integral() << endl;
46 buchmann 1.9 THStack mcm = allsamples.DrawStack("mc", "mll",1,0,10000000, "m_{ll} [GeV]", "events", cut,mc,luminosity);
47 buchmann 1.1 int nhists = mcm.GetHists()->GetEntries();
48 buchmann 1.7 dout << "Going to loop over " << nhists << " histograms " << endl;
49 buchmann 1.1 TFile *test = new TFile("test.root","RECREATE");
50     mcm.Write();
51     test->Close();
52     for (int i = 0; i < nhists; i++) {
53     TH1* hist = static_cast<TH1*>(mcm.GetHists()->At(i));
54     cout << " " << hist->GetName() << ": " << hist->Integral() << endl;
55     }
56 buchmann 1.9
57 buchmann 1.1 cout << endl << endl;
58 buchmann 1.9
59 buchmann 1.1 delete data;
60     }
61    
62 buchmann 1.2 void draw_kin_variable(string variable, TCut cut, int nbins, float min, float max, string xlabel, string saveas, bool dologscale, string speciallabel="") {
63 buchmann 1.1
64     TCanvas *ckin = new TCanvas("ckin","kin variable canvas");
65     ckin->SetLogy(dologscale);
66     TH1F *datahisto = allsamples.Draw("datahisto",variable,nbins,min,max, xlabel, "events",cut,data,luminosity);
67     datahisto->SetMarkerSize(DataMarkerSize);
68     THStack mcstack = allsamples.DrawStack("mcstack",variable,nbins,min,max, xlabel, "events",cut,mc,luminosity);
69 buchmann 1.2 if(dologscale) datahisto->SetMaximum(3*datahisto->GetMaximum());
70     else datahisto->SetMaximum(1.5*datahisto->GetMaximum());
71     TText *speciall;
72     if(speciallabel!="") {
73     speciall = new TText(0.2,0.8,speciallabel.c_str());
74     speciall->Draw();
75     }
76 buchmann 1.1 datahisto->Draw("e1");
77     ckin->Update();
78     mcstack.Draw("same");
79     datahisto->Draw("same,e1");
80     TLegend *kinleg = allsamples.allbglegend();
81     kinleg->Draw();
82    
83     TPad *kinpad = new TPad("kinpad","kinpad",0,0,1,1);
84     kinpad->cd();
85     kinpad->SetLogy(dologscale);
86     datahisto->Draw("e1");
87     mcstack.Draw("same");
88     datahisto->Draw("same,e1");
89     datahisto->Draw("same,axis");
90     kinleg->Draw();
91     DrawPrelim();
92     saveas="kin/"+saveas;
93     save_with_ratio(datahisto,mcstack,kinpad->cd(),saveas);
94    
95    
96    
97    
98     ZbData::data_over_mc.push_back(datahisto->Integral()/CollapseStack(mcstack)->Integral());
99    
100     datahisto->Delete();
101     delete ckin;
102     }
103    
104     Value getfrom2Dmap(TH2F *map, int ixbin, int iybin) {
105     Value sum(0,0);
106     for(int iy=1;iy<=iybin;iy++) {
107     sum=sum+Value(map->GetBinContent(ixbin,iy),map->GetBinError(ixbin,iy));
108     }
109     return sum;
110     }
111    
112     void print_all_b_yields() {
113     cout << "Basic selection with a b jet" << endl;
114 buchmann 1.9 print_yield(ZplusBsel&&"Zb3010_pfJetGoodNumBtag>0");
115 buchmann 1.2 cout << "Leading jet is a b-jet " << endl;
116 buchmann 1.1 print_yield(ZplusBsel&&LeadingB);
117     cout << "|#eta_{b}|<1.3 " << endl;
118     print_yield(ZplusBsel&&LeadingB&&EtaB);
119     cout << "|#delta#phi(b<Z)|>2.7" << endl;
120     print_yield(ZplusBsel&&LeadingB&&EtaB&&PhiZcut);
121 buchmann 1.2 cout << "10<ptZ<1000 GeV" << endl;
122 buchmann 1.1 print_yield(ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000");
123 buchmann 1.13 cout << "#alpha < 0.35" << endl;
124     print_yield(ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"(Zb3010_alpha)<0.35");
125 buchmann 1.1 cout << "#alpha < 0.3" << endl;
126 buchmann 1.9 print_yield(ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"(Zb3010_alpha)<0.3");
127 buchmann 1.1 cout << "#alpha < 0.2" << endl;
128 buchmann 1.9 print_yield(ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"(Zb3010_alpha)<0.2");
129 buchmann 1.1 cout << "#alpha < 0.15" << endl;
130 buchmann 1.9 print_yield(ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"(Zb3010_alpha)<0.15");
131 buchmann 1.1 cout << "#alpha < 0.1" << endl;
132 buchmann 1.9 print_yield(ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"(Zb3010_alpha)<0.1");
133     // cout << "#alpha < 0.05" << endl;
134     // print_yield(ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"(Zb3010_alpha)<0.05");
135 buchmann 1.1 }
136    
137 buchmann 1.7 void DrawEvilCutFlow() {
138     stringstream MegaCut;
139    
140 buchmann 1.9 MegaCut << " 1*(" << (const char*) (ZplusBsel&&TCut("Zb3010_pfJetGoodNumBtag>0")) << ")";
141 buchmann 1.7 MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB) << ")";
142     MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB&&EtaB) << ")";
143     MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB&&EtaB&&PhiZcut) << ")";
144     MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut("pt>10&&pt<1000")) << ")";
145 buchmann 1.9 MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut("pt>10&&pt<1000&&Zb3010_alpha<0.3")) << ")";
146     MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut("pt>10&&pt<1000&&Zb3010_alpha<0.2")) << ")";
147     MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut("pt>10&&pt<1000&&Zb3010_alpha<0.15")) << ")";
148     MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut("pt>10&&pt<1000&&Zb3010_alpha<0.1")) << ")";
149     MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut("pt>10&&pt<1000&&Zb3010_alpha<0.05")) << ")";
150 buchmann 1.7 MegaCut << " ";
151    
152     TCanvas *can = new TCanvas("can","can");
153     can->SetLogy(1);
154     TCut basecut("mll>6||mll<10");
155    
156    
157    
158     TLegend *leg = allsamples.allbglegend();
159    
160    
161     TH1F *data = allsamples.Draw("data", MegaCut.str(),11,-0.5,10.5, "", "events", basecut,data,luminosity);
162     THStack mcm = allsamples.DrawStack("mc", MegaCut.str(),11,-0.5,10.5, "", "events", basecut,mc,luminosity);
163    
164     float runningsum=0;
165     for(int i=data->GetNbinsX();i>0;i--) {
166     runningsum+=data->GetBinContent(i);
167     data->SetBinContent(i,runningsum);
168     }
169    
170     float minimum=data->GetMaximum();
171    
172     TH1F* h;
173     TIter nextOF(mcm.GetHists());
174     while ( h = (TH1F*)nextOF() ) {
175     float runningsum=0;
176     minimum=data->GetMaximum();//done deliberately at each step so get the minimum of the last (!) contribution
177     for(int i=h->GetNbinsX();i>0;i--) {
178     runningsum+=h->GetBinContent(i);
179     h->SetBinContent(i,runningsum);
180     if(runningsum<minimum&&runningsum>0) minimum=runningsum;
181     }
182     }
183    
184    
185    
186     data->SetMinimum(0.05*minimum);
187     can->cd(1)->SetBottomMargin(0.25);
188     data->GetXaxis()->SetBinLabel(1,"Pre-selection");
189     data->GetXaxis()->SetBinLabel(2,"Contains b jet");
190     data->GetXaxis()->SetBinLabel(3,"Leading jet is b-jet");
191     data->GetXaxis()->SetBinLabel(4,"|#eta_{b}|<1.3 ");
192     data->GetXaxis()->SetBinLabel(5,"|#delta#phi(b,Z)|>2.7");
193     data->GetXaxis()->SetBinLabel(6,"10<p_{T}^{Z}<1000 GeV");
194     data->GetXaxis()->SetBinLabel(7,"#alpha < 0.3");
195     data->GetXaxis()->SetBinLabel(8,"#alpha < 0.2");
196     data->GetXaxis()->SetBinLabel(9,"#alpha < 0.15");
197     data->GetXaxis()->SetBinLabel(10,"#alpha < 0.1");
198     data->GetXaxis()->SetBinLabel(11,"#alpha < 0.05");
199     data->GetXaxis()->LabelsOption("v");
200    
201     data->Draw("e1");
202     mcm.Draw("same");
203     data->Draw("same,e1,axis");
204     data->Draw("same,e1");
205     // data->Draw("same,TEXT");
206    
207     leg->Draw();
208    
209     DrawPrelim();
210    
211     CompleteSave(can,"CutFlow");
212     }
213    
214 buchmann 1.1 void draw_Zb_kin_vars() {
215 buchmann 1.6 draw_kin_variable("pt",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000",25,0,200,"Z p_{T}","Official/Zpt",1);
216 buchmann 1.9 draw_kin_variable("Zb3010_pfJetGoodPt[1]",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000",25,0,200,"Sub-Leading Jet Pt","Official/Jet2Pt",1);
217     draw_kin_variable("Zb3010_pfJetGoodPt[0]",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000",25,0,200,"Leading Jet Pt (B)","Official/LeadingJetPt",1);
218     draw_kin_variable("Zb3010_pfJetGoodPt[1]/pt",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000",20,0,2,"p_{T}^{J2} / p_{T}^{Z}","Official/SubLeadingJetPt_Over_ZPt",1);
219     draw_kin_variable("Zb3010_alpha",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>30&&pt<50", 40,0,2,"#alpha","Official/alpha_pt_30_to_50",1);
220     draw_kin_variable("Zb3010_alpha",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>50&&pt<75", 40,0,2,"#alpha","Official/alpha_50_to_75",1);
221     draw_kin_variable("Zb3010_alpha",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>75&&pt<125", 40,0,2,"#alpha","Official/alpha_pt_75_to_125",1);
222     draw_kin_variable("Zb3010_alpha",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>125&&pt<1000",40,0,2,"#alpha","Official/alpha_pt_125_to_1000",1);
223     draw_kin_variable("Zb3010_alpha",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>30&&pt<1000",40,0,2,"#alpha","Official/alpha_pt_30_to_1000",1);
224 buchmann 1.6
225 buchmann 1.10 draw_kin_variable("Zb3010_pfBJetDphiZ[0]",ZplusBsel&&LeadingB&&"pt>10&&pt<1000&&pfBJetDphiZ[0]>0",30,0,3.2,"#delta#phi (Z,b lead)","DeltaPhiZBlead",0);
226 buchmann 1.6
227     /*
228     draw_kin_variable("pt1",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000",20,0,200,"p_{T}^{l1}","Lep/pt1",1);
229     draw_kin_variable("pt1",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000&&id1==0",20,0,200,"p_{T}^{e1}","Lep/pt1_e",1);
230     draw_kin_variable("pt1",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000&&id1==1",20,0,200,"p_{T}^{#mu1}","Lep/pt1_m",1);
231     draw_kin_variable("pt2",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000",20,0,200,"p_{T}^{l2}","Lep/pt2",1);
232     draw_kin_variable("pt2",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000&&id2==0",20,0,200,"p_{T}^{e2}","Lep/pt2_e",1);
233     draw_kin_variable("pt2",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000&&id2==1",20,0,200,"p_{T}^{#mu2}","Lep/pt2_m",1);
234     draw_kin_variable("eta1",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000",20,-3.1,3.1,"#eta^{l1}","Lep/eta1",1);
235     draw_kin_variable("eta1",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000&&id1==0",20,-3.1,3.1,"#eta^{e1}","Lep/eta1_e",1);
236     draw_kin_variable("eta1",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000&&id1==1",20,-3.1,3.1,"#eta^{m1}","Lep/eta1_m",1);
237     draw_kin_variable("eta2",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000",20,-3.1,3.1,"#eta^{l1}","Lep/eta2",1);
238     draw_kin_variable("eta2",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000&&id2==0",20,-3.1,3.1,"#eta^{e2}","Lep/eta2_e",1);
239     draw_kin_variable("eta2",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000&&id2==1",20,-3.1,3.1,"#eta^{#mu2}","Lep/eta2_m",1);
240     draw_kin_variable("numVtx",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000",30,0,30,"N(Vtx)","Lep/nVtx",1);
241     */
242 buchmann 1.1 }
243    
244     void draw_mpf_vars() {
245     ZbData::data_over_mc.clear();
246 buchmann 1.9 draw_kin_variable("mpf",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"Zb3010_alpha<0.3",20,0,2,"Z+b MPF","mpf_alpha_smaller_0p3",0,"#alpha<0.3, 10 GeV < p_{T}^{Z} < 1000 GeV");
247     draw_kin_variable("mpf",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"Zb3010_alpha<0.2",20,0,2,"Z+b MPF","mpf_alpha_smaller_0p2",0,"#alpha<0.2, 10 GeV < p_{T}^{Z} < 1000 GeV");
248     draw_kin_variable("mpf",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"Zb3010_alpha<0.15",20,0,2,"Z+b MPF","mpf_alpha_smaller_0p15",0,"#alpha<0.15, 10 GeV < p_{T}^{Z} < 1000 GeV");
249     draw_kin_variable("mpf",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"Zb3010_alpha<0.1",20,0,2,"Z+b MPF","mpf_alpha_smaller_0p1",0,"#alpha<0.1, 10 GeV < p_{T}^{Z} < 1000 GeV");
250     draw_kin_variable("mpf",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"Zb3010_alpha<0.05",20,0,2,"Z+b MPF","mpf_alpha_smaller_0p05",0,"#alpha<0.05, 10 GeV < p_{T}^{Z} < 1000 GeV");
251    
252     draw_kin_variable("Zb3010_pfJetGoodPt[0]/pt",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"Zb3010_alpha<0.3",20,0,2,"p_{T} b-jet / p_{T} Z","ptb_over_ptz___alpha_smaller_0p3",0,"#alpha<0.3, 10 GeV < p_{T}^{Z} < 1000 GeV");
253     draw_kin_variable("Zb3010_pfJetGoodPt[0]/pt",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"Zb3010_alpha<0.2",20,0,2,"p_{T} b-jet / p_{T} Z","ptb_over_ptz___alpha_smaller_0p2",0,"#alpha<0.2, 10 GeV < p_{T}^{Z} < 1000 GeV");
254     draw_kin_variable("Zb3010_pfJetGoodPt[0]/pt",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"Zb3010_alpha<0.15",20,0,2,"p_{T} b-jet / p_{T} Z","ptb_over_ptz___alpha_smaller_0p15",0,"#alpha<0.15, 10 GeV < p_{T}^{Z} < 1000 GeV");
255     draw_kin_variable("Zb3010_pfJetGoodPt[0]/pt",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"Zb3010_alpha<0.1",20,0,2,"p_{T} b-jet / p_{T} Z","ptb_over_ptz___alpha_smaller_0p1",0,"#alpha<0.1, 10 GeV < p_{T}^{Z} < 1000 GeV");
256     draw_kin_variable("Zb3010_pfJetGoodPt[0]/pt",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"Zb3010_alpha<0.05",20,0,2,"p_{T} b-jet / p_{T} Z","ptb_over_ptz__mpf_alpha_smaller_0p05",0,"#alpha<0.05, 10 GeV < p_{T}^{Z} < 1000 GeV");
257 buchmann 1.1 }
258    
259 buchmann 1.12 /*double GetMedian(TH1F *histo) {
260 buchmann 1.11 int numBins = histo->GetXaxis()->GetNbins();
261 buchmann 1.12 Double_t x[numBins];
262     Double_t y[numBins];
263 buchmann 1.11 for (int i = 1; i <= numBins; i++) {
264     x[i] = histo->GetBinCenter(i);
265     y[i] = histo->GetBinContent(i);
266     }
267 buchmann 1.12 return TMath::Median(numBins, &x, &y);
268     }*/
269 buchmann 1.11
270 buchmann 1.12 Value get_Zb_data_over_mc(string variable, string variable2, TCut cut, string saveas) {
271 buchmann 1.2 // write_warning(__FUNCTION__,"Debugging this function, therefore always returning 3.1415+/-1.010101"); return Value(3.1415,1.010101);
272 buchmann 1.12 TCanvas *cn = new TCanvas("cn","cn");
273     string varname="MPF";
274     if(Contains(variable,"pfJetGood")) varname="R_{abs}";
275     TH1F *hdata = allsamples.Draw("hdata",variable,1200,0,30, varname, "events", cut,data,luminosity);
276     TH1F *hmc = allsamples.Draw("hmc" , variable,1200,0,30, varname, "events", cut, mc, luminosity);
277 buchmann 1.2
278 buchmann 1.10 float a=hdata->GetMean();
279 buchmann 1.2 float b=hmc->GetMean();
280     float da=hdata->GetMeanError();
281 buchmann 1.10 float db=hmc->GetMeanError();
282 buchmann 1.2 float factor = a / b;
283     float error = TMath::Sqrt( (1/(b*b)) * (da*da) + ((a*a)/(b*b))*db*db);
284 buchmann 1.12
285     cn->cd();
286     hdata->GetYaxis()->SetTitle("events / 0.1");
287     hdata->Rebin(4);
288     hmc->Rebin(4);
289     hdata->GetXaxis()->SetRangeUser(0,2);
290     hdata->Draw("e1");
291     hmc->Draw("histo,same");
292     hdata->Draw("e1,same");
293     stringstream summary;
294     summary << "#splitline{" << varname << "}{#splitline{data: " << std::setprecision(4) << a << " #pm " << std::setprecision(4) << da << "}{#splitline{mc: " << std::setprecision(4) << b << " #pm " << std::setprecision(4) << db << "}{ratio: " << factor;
295     summary << " #pm " << error << "}}}";
296     TText *infobox = write_title(summary.str());
297     infobox->SetX(0.75);
298     infobox->SetTextSize(0.03);
299     infobox->SetY(0.75);
300     infobox->Draw();
301    
302    
303     DrawPrelim();
304     CompleteSave(cn,"ResponsePlots/"+saveas);
305    
306     cout << "Data;" << a << ";"<<da<<";MC"<<b<<";"<<db<<endl;
307     delete cn;
308 buchmann 1.2 delete hdata;
309     delete hmc;
310 buchmann 1.9
311     // cout << a << "+/-" << da << " " << b << "=/-" << db << endl;
312    
313 buchmann 1.2 return Value(factor,error);
314     }
315    
316    
317    
318     void new_data_mc_agreement_2d() {
319 buchmann 1.13
320     write_error(__FUNCTION__,"need to make alpha<0.3 inclusive again!!!! IMPORTANT SINCE THAT'S THE FINAL RESULT!!!");
321    
322 buchmann 1.3 gStyle->SetOptFit(0);
323 buchmann 1.2 TCut basecut(ZplusBsel&&LeadingB&&EtaB&&PhiZcut);
324 buchmann 1.13 const int nalphacuts=6;
325     float alphacuts[nalphacuts] = {0.1,0.15,0.2,0.25,0.3,0.35};
326 buchmann 1.5 int PointThree=-1;
327     for(int i=0;i<nalphacuts;i++) {
328 buchmann 1.9 if((alphacuts[i]-0.3)<0.01) PointThree=i-1;
329 buchmann 1.5 }
330     if(PointThree<0) {
331     write_warning(__FUNCTION__,"alpha=0.3 MUST BE ONE OF THE ALPHA VALUES!");
332     assert(PointThree<0);
333     }
334    
335 buchmann 1.2 const int nptcuts=5;
336 buchmann 1.9 float ptcuts[nptcuts]={30,50,75,125,1000};
337 buchmann 1.2
338     float MPF_Results[nalphacuts][nptcuts];
339     float MPF_Errors[nalphacuts][nptcuts];
340     float RABS_Results[nalphacuts][nptcuts];
341     float RABS_Errors[nalphacuts][nptcuts];
342    
343    
344 buchmann 1.9 for(int ia=0;ia<nalphacuts;ia++) {
345 buchmann 1.2 for(int ipt=0;ipt<nptcuts-1;ipt++) {
346     stringstream specialcut;
347 buchmann 1.13 if(ia>0) specialcut << "((pt>" << ptcuts[ipt] << " && pt< " << ptcuts[ipt+1] << ") && (Zb3010_alpha<" << alphacuts[ia] << " && Zb3010_alpha>" << alphacuts[ia-1] << "))";
348     else specialcut << "((pt>" << ptcuts[ipt] << " && pt< " << ptcuts[ipt+1] << ") && (Zb3010_alpha<" << alphacuts[ia] << "))";
349 buchmann 1.9 cout << specialcut.str() << endl;
350 buchmann 1.12 Value MPF_data_over_mc = get_Zb_data_over_mc("mpf","",TCut(basecut && specialcut.str().c_str()),"MPF___pt_"+any2string(ptcuts[ipt])+"_"+any2string(ptcuts[ipt+1])+"__alpha_"+any2string(alphacuts[ia]));
351     Value RABS_data_over_mc = get_Zb_data_over_mc("Zb3010_pfJetGoodPt[0]/pt","",TCut(basecut && specialcut.str().c_str()),"Rabs___pt_"+any2string(ptcuts[ipt])+"_"+any2string(ptcuts[ipt+1])+"__alpha_"+any2string(alphacuts[ia]));
352 buchmann 1.2
353     MPF_Results[ia][ipt]=MPF_data_over_mc.getValue();
354     MPF_Errors[ia][ipt]=MPF_data_over_mc.getError();
355    
356     RABS_Results[ia][ipt]=RABS_data_over_mc.getValue();
357     RABS_Errors[ia][ipt]=RABS_data_over_mc.getError();
358    
359 buchmann 1.9 cout << "alpha cut at " << alphacuts[ia+1] << " and pt range " << ptcuts[ipt] << " , " << ptcuts[ipt+1] << " \t : \t" << specialcut.str() << " \t \t --> " << MPF_data_over_mc << " (MPF) " << RABS_data_over_mc << " (RABS)" << endl;
360 buchmann 1.2 }
361     }
362    
363     float MPF_ExtraPolatedResults[nptcuts];
364     float RABS_ExtraPolatedResults[nptcuts];
365    
366     TCanvas *can = new TCanvas("can");
367    
368     TGraphErrors *MPF_FinalGraph = new TGraphErrors();
369     TGraphErrors *RABS_FinalGraph = new TGraphErrors();
370 buchmann 1.13 TGraphErrors *MPF_FaceValueAtPoint3 = new TGraphErrors();
371     TGraphErrors *RABS_FaceValueAtPoint3 = new TGraphErrors();
372 buchmann 1.2
373     for(int ipt=0;ipt<nptcuts-1;ipt++) {
374    
375     stringstream filename;
376     filename << "Extrapolation/Zplusb_data_over_mc___" << ptcuts[ipt] << "_to_" << ptcuts[ipt+1];
377     cout << "Going to create histo for " << filename.str() << endl;
378     TGraphErrors *mpf_gr =new TGraphErrors();
379     TGraphErrors *rabs_gr =new TGraphErrors();
380 buchmann 1.13
381     int rabs_counter=0;
382     int mpf_counter=0;
383    
384 buchmann 1.2 for(int ia=0;ia<nalphacuts;ia++) {
385 buchmann 1.9 cout << "Setting a point for an alpha cut at " << alphacuts[ia] << " with value " << MPF_Results[ia][ipt] << " (MPF) and " << RABS_Results[ia][ipt] << " (RABS)" << endl;
386 buchmann 1.13 if(MPF_Results[ia][ipt]==MPF_Results[ia][ipt] && MPF_Errors[ia][ipt]==MPF_Errors[ia][ipt]) { // remove nan's
387     mpf_gr->SetPoint(mpf_counter,alphacuts[ia],MPF_Results[ia][ipt]);
388     mpf_gr->SetPointError(mpf_counter,0,MPF_Errors[ia][ipt]);
389     mpf_counter++;
390     } else {
391     dout << "Detected BS for MPF method (bin ignored) : " << endl;
392     dout << " alpha: " << alphacuts[ia] << " , pt range: " << ptcuts[ipt] << " < pt < " << ptcuts[ipt+1] << endl;
393     cout << " value (MPF) : " << MPF_Results[ia][ipt] << " +/- " << MPF_Errors[ia][ipt] << endl;
394     }
395     if(RABS_Results[ia][ipt]==RABS_Results[ia][ipt] && RABS_Errors[ia][ipt]==RABS_Errors[ia][ipt]) {
396     rabs_gr->SetPoint(rabs_counter,alphacuts[ia],RABS_Results[ia][ipt]);
397     rabs_gr->SetPointError(rabs_counter,0,RABS_Errors[ia][ipt]);
398     rabs_counter++;
399     } else {
400     dout << "Detected BS for RABS method (bin ignored) : " << endl;
401     dout << " alpha: " << alphacuts[ia] << " , pt range: " << ptcuts[ipt] << " < pt < " << ptcuts[ipt+1] << endl;
402     cout << " value (RABS) : " << RABS_Results[ia][ipt] << " +/- " << RABS_Errors[ia][ipt] << endl;
403     }
404     }
405     cout << "About to add a point in pt (" << (ptcuts[ipt]+ptcuts[ipt+1])/2 << " ) to face value : " << MPF_Results[PointThree][ipt] << endl;
406    
407     stringstream specialcut;
408     specialcut << "((pt>" << ptcuts[ipt] << " && pt< " << ptcuts[ipt+1] << ") && (Zb3010_alpha<0.3))";
409     Value MPF_data_over_mc = get_Zb_data_over_mc("mpf","",TCut(basecut && specialcut.str().c_str()),"MPF___pt_"+any2string(ptcuts[ipt])+"_"+any2string(ptcuts[ipt+1])+"__INCLUSIVEalpha_0.3");
410     Value RABS_data_over_mc = get_Zb_data_over_mc("Zb3010_pfJetGoodPt[0]/pt","",TCut(basecut && specialcut.str().c_str()),"Rabs___pt_"+any2string(ptcuts[ipt])+"_"+any2string(ptcuts[ipt+1])+"__INCLUSIVEalpha_0.3");
411     MPF_FaceValueAtPoint3->SetPoint(ipt,(ptcuts[ipt]+ptcuts[ipt+1])/2,MPF_data_over_mc.getValue());
412     MPF_FaceValueAtPoint3->SetPointError(ipt,0,MPF_data_over_mc.getError());
413     RABS_FaceValueAtPoint3->SetPoint(ipt,(ptcuts[ipt]+ptcuts[ipt+1])/2,RABS_data_over_mc.getValue());
414     RABS_FaceValueAtPoint3->SetPointError(ipt,0,RABS_data_over_mc.getError());
415 buchmann 1.2
416 buchmann 1.13 can->cd();
417 buchmann 1.2 mpf_gr->SetMarkerStyle(21);
418     mpf_gr->SetMarkerColor(kBlue);
419     mpf_gr->Draw("AP");
420 buchmann 1.13 mpf_gr->Fit("pol1","");
421 buchmann 1.3 mpf_gr->GetXaxis()->SetTitle("#alpha");
422     mpf_gr->GetXaxis()->CenterTitle();
423 buchmann 1.9 mpf_gr->GetYaxis()->SetTitle("<MPF>_{data}/<MPF>_{mc}");
424 buchmann 1.3 mpf_gr->GetYaxis()->CenterTitle();
425 buchmann 1.2 TF1 *mpf_pol=(TF1*)mpf_gr->GetFunction("pol1");
426 buchmann 1.13 TF1 *mpf_pol_complete = new TF1("mpf_pol_complete","[0]+[1]*x");
427     mpf_pol_complete->SetParameters(mpf_pol->GetParameters());
428     mpf_pol_complete->SetLineWidth(1);
429     mpf_pol_complete->SetLineColor(kBlue);
430     float min,max;
431     FindMinMax(mpf_gr,min,max);
432     TH1F *histo = new TH1F("histo","histo",1,0,0.4);
433     histo->SetStats(0);
434     histo->GetXaxis()->SetRangeUser(0,0.4);
435     histo->GetYaxis()->SetRangeUser(min,max);
436     mpf_gr->GetYaxis()->SetRangeUser(min,max);
437    
438     TPolyLine *mpf_fit_uncert = GetFitUncertaintyShape(mpf_gr, "pol1", min, max,0.0,0.4);
439     mpf_pol->SetLineColor(TColor::GetColor("#0101DF"));
440     mpf_fit_uncert->SetFillColor(TColor::GetColor("#0101DF"));
441     histo->Draw();
442     mpf_fit_uncert->Draw("same,f");
443     mpf_gr->Draw("P");
444     mpf_pol_complete->Draw("same");
445    
446 buchmann 1.2 float mpf_a=mpf_pol->GetParameter(0);
447     float mpf_b=mpf_pol->GetParameter(1);
448     MPF_FinalGraph->SetPoint(ipt,0.5*(ptcuts[ipt]+ptcuts[ipt+1]),mpf_a);
449     MPF_FinalGraph->SetPointError(ipt,0,mpf_pol->GetParError(0));
450    
451     MPF_ExtraPolatedResults[ipt]=mpf_a;
452    
453     stringstream mpf_info;
454 buchmann 1.3 mpf_info << "#splitline{#splitline{#splitline{" << ptcuts[ipt] << " GeV < p_{T}^{Z} < " << ptcuts[ipt+1] << " GeV }{ (MPF) }}{Extrapolated value: " << std::setprecision(3) << mpf_a << "}}{#chi^{2}/NDF : " << mpf_pol->GetChisquare() << " / " << mpf_pol->GetNDF() << "}";
455    
456 buchmann 1.2 TText *mpf_mark = write_title(mpf_info.str());
457     mpf_mark->SetX(0.75);
458     mpf_mark->SetTextSize(0.03);
459     mpf_mark->SetY(0.75);
460     mpf_mark->Draw();
461    
462     string filenamebkp=filename.str();
463     filename << "__MPF";
464 buchmann 1.3 DrawPrelim();
465 buchmann 1.2 CompleteSave(can,filename.str());
466     cout << "MPF : " << mpf_a << " + " << mpf_b << " * alpha " << endl;
467    
468     filename.str("");
469    
470     rabs_gr->SetMarkerStyle(21);
471 buchmann 1.13 rabs_gr->SetMarkerColor(kRed);
472 buchmann 1.2 rabs_gr->Draw("AP");
473 buchmann 1.13 rabs_gr->Fit("pol1","");
474 buchmann 1.3 rabs_gr->GetXaxis()->SetTitle("#alpha");
475     rabs_gr->GetXaxis()->CenterTitle();
476 buchmann 1.13 rabs_gr->GetYaxis()->SetTitle("<MPF>_{data}/<MPF>_{mc}");
477 buchmann 1.3 rabs_gr->GetYaxis()->CenterTitle();
478 buchmann 1.13 TF1 *rabs_pol=(TF1*)rabs_gr->GetFunction("pol1");
479     TF1 *rabs_pol_complete = new TF1("rabs_pol_complete","[0]+[1]*x");
480     rabs_pol_complete->SetParameters(rabs_pol->GetParameters());
481     rabs_pol_complete->SetLineColor(kRed);
482     rabs_pol_complete->SetLineWidth(1);
483     FindMinMax(rabs_gr,min,max);
484     rabs_gr->GetYaxis()->SetRangeUser(min,max);
485     histo->GetYaxis()->SetRangeUser(min,max);
486     TPolyLine *rabs_fit_uncert = GetFitUncertaintyShape(rabs_gr, "pol1", min, max,0.0,0.4);
487     rabs_pol->SetLineColor(TColor::GetColor("#0101DF"));
488     rabs_pol->SetFillColor(TColor::GetColor("#0101DF"));
489     histo->Draw();
490     rabs_fit_uncert->Draw("same,f");
491     rabs_gr->Draw("P");
492     rabs_pol_complete->Draw("same");
493    
494 buchmann 1.3
495 buchmann 1.2 float rabs_a=rabs_pol->GetParameter(0);
496     float rabs_b=rabs_pol->GetParameter(1);
497     RABS_FinalGraph->SetPoint(ipt,0.5*(ptcuts[ipt]+ptcuts[ipt+1]),rabs_a);
498     RABS_FinalGraph->SetPointError(ipt,0,rabs_pol->GetParError(0));
499 buchmann 1.3
500 buchmann 1.2 cout << "!+!+!+!+!+!+!+!+!+ Just added a point to the final plot : " << rabs_a << " +/- " << rabs_pol->GetParError(0) << endl;
501    
502     stringstream rabs_info;
503 buchmann 1.3 rabs_info << "#splitline{#splitline{#splitline{" << ptcuts[ipt] << " GeV < p_{T}^{Z} < " << ptcuts[ipt+1] << " GeV }{ (R_{abs}) }}{Extrapolated value: " << std::setprecision(3) << rabs_a << "}}{#chi^{2}/NDF : " << rabs_pol->GetChisquare() << " / " << rabs_pol->GetNDF() << "}";
504 buchmann 1.2 TText *rabs_mark = write_title(rabs_info.str());
505     rabs_mark->SetX(0.75);
506     rabs_mark->SetTextSize(0.03);
507     rabs_mark->SetY(0.75);
508     rabs_mark->Draw();
509 buchmann 1.3
510    
511 buchmann 1.2 RABS_ExtraPolatedResults[ipt]=rabs_a;
512     filename << filenamebkp << "__RABS";
513 buchmann 1.3 DrawPrelim();
514 buchmann 1.2 CompleteSave(can,filename.str());
515     cout << "RABS : " << rabs_a << " + " << rabs_b << " * alpha " << endl;
516     }
517    
518     MPF_FinalGraph->SetMarkerStyle(21);
519 buchmann 1.13 MPF_FinalGraph->SetMarkerColor(kBlue);
520     MPF_FinalGraph->Fit("pol0",""); // quiet, use minos
521 buchmann 1.3
522 buchmann 1.2 TF1 *mpf_pol0=(TF1*)MPF_FinalGraph->GetFunction("pol0");
523 buchmann 1.13 TF1 *mpf_pol0_complete = new TF1("mpf_pol0_complete","[0]+[1]*x");
524     mpf_pol0_complete->SetLineWidth(1);
525     mpf_pol0_complete->SetLineColor(kBlue);
526     mpf_pol0->SetLineColor(kBlue);
527 buchmann 1.2 MPF_FinalGraph->Draw("AP");
528 buchmann 1.3 MPF_FinalGraph->GetXaxis()->SetTitle("p_{T}^{Z}");
529     MPF_FinalGraph->GetXaxis()->CenterTitle();
530 buchmann 1.9 MPF_FinalGraph->GetYaxis()->SetTitle("<MPF>_{data}/<MPF>_{mc}");
531 buchmann 1.3 MPF_FinalGraph->GetYaxis()->CenterTitle();
532     MPF_FinalGraph->Draw("AP");
533 buchmann 1.13 mpf_pol0_complete->Draw("same");
534 buchmann 1.3
535 buchmann 1.2 stringstream mpf_result;
536 buchmann 1.3 mpf_result << "#splitline{C_{abs}= " << std::setprecision(5) << 1/mpf_pol0->GetParameter(0) << " #pm " << std::setprecision(2) << mpf_pol0->GetParError(0)/(mpf_pol0->GetParameter(0)*mpf_pol0->GetParameter(0)) << "}{ (MPF)}";
537 buchmann 1.2 TText *rmpf_final_mark = write_title(mpf_result.str());
538     rmpf_final_mark->SetX(0.75);
539     rmpf_final_mark->SetY(0.75);
540     rmpf_final_mark->SetTextSize(0.03);
541     rmpf_final_mark->Draw();
542 buchmann 1.3 DrawPrelim();
543 buchmann 1.5 float MPFResult=1/mpf_pol0->GetParameter(0);
544     float MPFResultError=mpf_pol0->GetParError(0)/(mpf_pol0->GetParameter(0)*mpf_pol0->GetParameter(0));
545 buchmann 1.2
546     stringstream filename2;
547     filename2 << "Extrapolation/FINAL_RESULT_MPF";
548     CompleteSave(can,filename2.str());
549    
550    
551     RABS_FinalGraph->SetMarkerStyle(21);
552     RABS_FinalGraph->SetMarkerStyle(21);
553     RABS_FinalGraph->SetMarkerColor(kRed);
554 buchmann 1.3 RABS_FinalGraph->Fit("pol0","QE"); // quiet, use minos
555     RABS_FinalGraph->Draw("AP");
556     RABS_FinalGraph->GetXaxis()->SetTitle("p_{T}^{Z}");
557     RABS_FinalGraph->GetXaxis()->CenterTitle();
558 buchmann 1.9 RABS_FinalGraph->GetYaxis()->SetTitle("<R_{abs}>_{data}/<R_{abs}>_{mc}");
559 buchmann 1.3 RABS_FinalGraph->GetYaxis()->CenterTitle();
560 buchmann 1.2 RABS_FinalGraph->Draw("AP");
561     TF1 *rabs_pol0=(TF1*)RABS_FinalGraph->GetFunction("pol0");
562     stringstream rabs_result;
563 buchmann 1.3 rabs_result << "#splitline{C_{abs}= " << std::setprecision(5) << 1/rabs_pol0->GetParameter(0) << " #pm " << std::setprecision(2) << rabs_pol0->GetParError(0)/(rabs_pol0->GetParameter(0)*rabs_pol0->GetParameter(0)) << "}{ (R_{abs})}";;
564 buchmann 1.2 TText *rabs_final_mark = write_title(rabs_result.str());
565     rabs_final_mark->SetX(0.75);
566     rabs_final_mark->SetY(0.75);
567     rabs_final_mark->SetTextSize(0.03);
568     rabs_final_mark->Draw();
569 buchmann 1.3 DrawPrelim();
570 buchmann 1.2
571 buchmann 1.5 float RabsResult=1/rabs_pol0->GetParameter(0);
572     float RabsResultError=rabs_pol0->GetParError(0)/(rabs_pol0->GetParameter(0)*rabs_pol0->GetParameter(0));
573    
574     filename2.str("");
575     filename2 << "Extrapolation/FINAL_RESULT_RABS";
576     CompleteSave(can,filename2.str());
577    
578 buchmann 1.9 can->cd();
579 buchmann 1.13 MPF_FaceValueAtPoint3->SetMarkerStyle(21);
580     MPF_FaceValueAtPoint3->SetMarkerColor(kBlue);
581     MPF_FaceValueAtPoint3->Fit("pol0","QE"); // quiet, use minos
582     MPF_FaceValueAtPoint3->Draw("AP");
583     MPF_FaceValueAtPoint3->GetXaxis()->SetTitle("p_{T}^{Z}");
584     MPF_FaceValueAtPoint3->GetXaxis()->CenterTitle();
585     MPF_FaceValueAtPoint3->GetYaxis()->SetTitle("<MPF>__{data}/<MPF>_{mc} for #alpha<0.3");
586     MPF_FaceValueAtPoint3->GetYaxis()->CenterTitle();
587     MPF_FaceValueAtPoint3->Draw("AP");
588     TF1 *faceval_pol0=(TF1*)MPF_FaceValueAtPoint3->GetFunction("pol0");
589     faceval_pol0->SetLineColor(kBlue);
590     faceval_pol0->SetLineWidth(1);
591     faceval_pol0->Draw("same");
592     TF1 *faceval_pol0_complete = new TF1("faceval_pol0_complete","[0]+[1]*x");
593     faceval_pol0_complete->SetParameters(faceval_pol0->GetParameters());
594 buchmann 1.5 stringstream faceval_result;
595     faceval_result << "#splitline{C_{abs}= " << std::setprecision(5) << 1/faceval_pol0->GetParameter(0) << " #pm " << std::setprecision(2) << faceval_pol0->GetParError(0)/(faceval_pol0->GetParameter(0)*faceval_pol0->GetParameter(0)) << "}{ (MPF at #alpha<0.3)}";;
596     TText *faceval_final_mark = write_title(faceval_result.str());
597     faceval_final_mark->SetX(0.75);
598     faceval_final_mark->SetY(0.75);
599     faceval_final_mark->SetTextSize(0.03);
600     faceval_final_mark->Draw();
601     DrawPrelim();
602    
603     float FaceValResult=1/faceval_pol0->GetParameter(0);
604     float FaceValResultError=faceval_pol0->GetParError(0)/(faceval_pol0->GetParameter(0)*faceval_pol0->GetParameter(0));
605    
606 buchmann 1.2 filename2.str("");
607 buchmann 1.9 filename2 << "Extrapolation/FINAL_RESULT_MPF_FaceValp3";
608 buchmann 1.2 CompleteSave(can,filename2.str());
609    
610 buchmann 1.13 can->cd();
611     RABS_FaceValueAtPoint3->SetMarkerStyle(21);
612     RABS_FaceValueAtPoint3->SetMarkerColor(kRed);
613     RABS_FaceValueAtPoint3->Fit("pol0","QE"); // quiet, use minos
614     RABS_FaceValueAtPoint3->Draw("AP");
615     RABS_FaceValueAtPoint3->GetXaxis()->SetTitle("p_{T}^{Z}");
616     RABS_FaceValueAtPoint3->GetXaxis()->CenterTitle();
617     RABS_FaceValueAtPoint3->GetYaxis()->SetTitle("<R_{abs}>_{data}/<R_{abs}>_{mc} for #alpha<0.3");
618     RABS_FaceValueAtPoint3->GetYaxis()->CenterTitle();
619     RABS_FaceValueAtPoint3->Draw("AP");
620     TF1 *faceval_pol0RABS=(TF1*)RABS_FaceValueAtPoint3->GetFunction("pol0");
621     TF1 *faceval_pol0RABS_complete = new TF1("faceval_pol0_complete","[0]+[1]*x");
622     faceval_pol0RABS_complete->SetParameters(faceval_pol0RABS->GetParameters());
623     faceval_pol0RABS->SetLineWidth(1);
624     faceval_pol0RABS->SetLineColor(kRed);
625     faceval_pol0RABS_complete->Draw("same");
626     stringstream RABSfaceval_result;
627     RABSfaceval_result << "#splitline{C_{abs}= " << std::setprecision(5) << 1/faceval_pol0RABS->GetParameter(0) << " #pm " << std::setprecision(2) << faceval_pol0RABS->GetParError(0)/(faceval_pol0RABS->GetParameter(0)*faceval_pol0RABS->GetParameter(0)) << "}{ (R_{abs} at #alpha<0.3)}";;
628     TText *RABSfaceval_final_mark = write_title(RABSfaceval_result.str());
629     RABSfaceval_final_mark->SetX(0.75);
630     RABSfaceval_final_mark->SetY(0.75);
631     RABSfaceval_final_mark->SetTextSize(0.03);
632     RABSfaceval_final_mark->Draw();
633     DrawPrelim();
634    
635     float RABSFaceValResult=1/faceval_pol0RABS->GetParameter(0);
636     float RABSFaceValResultError=faceval_pol0RABS->GetParError(0)/(faceval_pol0RABS->GetParameter(0)*faceval_pol0RABS->GetParameter(0));
637    
638     filename2.str("");
639     filename2 << "Extrapolation/FINAL_RESULT_RABS_FaceValp3";
640     CompleteSave(can,filename2.str());
641    
642 buchmann 1.5 cout << "FINAL RESULTS: " << endl;
643     cout << "MPF : " << MPFResult << " +/- " << MPFResultError << endl;
644     cout << "Rabs: " << RabsResult << " +/- " << RabsResultError << endl;
645 buchmann 1.13 cout << "Face value at 0.3 :" << FaceValResult << " +/- " << FaceValResultError << " (MPF) " << endl;
646     cout << "Face value at 0.3 :" << RABSFaceValResult << " +/- " << RABSFaceValResultError << " (RABS) " << endl;
647 buchmann 1.2
648     delete can;
649    
650    
651     }
652    
653    
654 buchmann 1.9 void DoPUStudy(string identifier) {
655    
656     float numVtxcuts[7]={0,10,15,20};
657    
658     TCanvas *can = new TCanvas("can","can");
659     TH1F *malpha[8];
660     TH1F *dalpha[8];
661     cout << identifier << ";";
662     cout << endl;
663     for(int i=1;i<5;i++) {
664     stringstream sSpecialCut;
665     if(i<4) {
666     sSpecialCut << "numVtx>" << numVtxcuts[i-1] << " && numVtx<" << numVtxcuts[i];
667     cout << numVtxcuts[i-1] << " < nVtx < " << numVtxcuts[i] << ";";
668     } else {
669     sSpecialCut << "numVtx>" << numVtxcuts[i-1];
670     cout << numVtxcuts[i-1] << ";";
671     }
672    
673     TCut SpecialCut(sSpecialCut.str().c_str());
674    
675    
676     malpha[i] = allsamples.Draw("malpha"+any2string(i),"mpf",1000,0,10,"MPF","events",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&SpecialCut,mc,luminosity);
677     dalpha[i] = allsamples.Draw("dalpha"+any2string(i),"mpf",1000,0,10,"MPF","events",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&SpecialCut,data,luminosity);
678    
679    
680     float a=dalpha[i]->GetMean();
681     float b=malpha[i]->GetMean();
682     float da=dalpha[i]->GetMeanError();
683     float db=malpha[i]->GetMeanError();
684     float factor = a / b;
685     float error = TMath::Sqrt( (1/(b*b)) * (da*da) + ((a*a)/(b*b))*db*db);
686    
687     cout << malpha[i]->GetMean() << ";" << malpha[i]->GetMeanError() << ";" << dalpha[i]->GetMean() << ";" << dalpha[i]->GetMeanError() << ";" << malpha[i]->Integral()<<";"<<dalpha[i]->Integral() << ";"<<factor << ";" << error << endl;
688     malpha[i]->SetLineColor(i+1);
689     if(i==0) malpha[i]->Draw("histo");
690     else malpha[i]->Draw("histo,same");
691     }
692    
693     CompleteSave(can,"AlphaOverview_"+identifier);
694     delete can;
695    
696     /* for(int i=1;i<8;i++) {
697     delete malpha[i];
698     delete dalpha[i];
699     }
700     */
701    
702     }
703    
704     void ScenarioComparison() {
705 buchmann 1.13 //very dumb way to do this.
706 buchmann 1.9 TGraphAsymmErrors *gr[5];
707     TF1 *fit[5];
708     bool dofit=false;
709    
710     TCanvas *can = new TCanvas("can","can",500,500);
711    
712     TLegend *leg = new TLegend(0.2,0.2,0.5,0.45);
713     leg->SetBorderSize(0);
714     leg->SetFillColor(kWhite);
715    
716     for(int i=0;i<5;i++) {
717     gr[i] = new TGraphAsymmErrors();
718     double x[4] = {5,12.5,17.5,20};
719     double dxl[4] = {5,2.5,2.5,2.5};
720     double dxh[4] = {5,2.5,2.5,10};
721    
722     string name="";
723     string color="";
724    
725     if(i==0) {
726     double y[4] = {1.008,1.013,1.003,1.007};
727     double dyl[4] = {0.011,0.009,0.013,0.020};
728     double dyh[4] = {0.011,0.009,0.013,0.020};
729     name="30/30";
730     gr[i] = new TGraphAsymmErrors(4,x,y,dxl,dxh,dyl,dyh);
731     color="#a01d99";
732     }
733     if(i==1) {
734     double y[4] = {1.008+0.001,1.013+0.001,1.003+0.001,1.007+0.001};
735     double dyl[4] = {0.011,0.009,0.013,0.020};
736     double dyh[4] = {0.011,0.009,0.013,0.020};
737     gr[i] = new TGraphAsymmErrors(4,x,y,dxl,dxh,dyl,dyh);
738 buchmann 1.13 name="30/15";
739 buchmann 1.9 color="#2464bc";
740     }
741     if(i==2) {
742     double y[4] = {1.014,1.026,1.014,1.059};
743     double dyl[4] = {0.011,0.010,0.014,0.023};
744     double dyh[4] = {0.011,0.010,0.014,0.023};
745     gr[i] = new TGraphAsymmErrors(4,x,y,dxl,dxh,dyl,dyh);
746     name="15/15";
747     color="#f49230";
748     }
749     if(i==3) {
750     double y[4] = {1.015,1.020,1.008,1.050};
751     double dyl[4] = {0.011,0.010,0.014,0.022};
752     double dyh[4] = {0.011,0.010,0.014,0.022};
753     gr[i] = new TGraphAsymmErrors(4,x,y,dxl,dxh,dyl,dyh);
754     name="20/15";
755     color="#c72f3c";
756     }
757     if(i==4) {
758     double y[4] = {1.015+0.001,1.020+0.001,1.008+0.001,1.050+0.001};
759     double dyl[4] = {0.011,0.010,0.014,0.022};
760     double dyh[4] = {0.011,0.010,0.014,0.022};
761     gr[i] = new TGraphAsymmErrors(4,x,y,dxl,dxh,dyl,dyh);
762     name="20/20";
763     color="#7cb433";
764     }
765    
766    
767     gr[i]->SetTitle();
768     gr[i]->GetXaxis()->SetTitle("N(Vtx)");
769     gr[i]->GetYaxis()->SetTitle("<MPF>_{data}/<MPF>_{mc}");
770     gr[i]->GetXaxis()->CenterTitle();
771     gr[i]->GetYaxis()->CenterTitle();
772     gr[i]->SetName(name.c_str());
773     gr[i]->SetMarkerSize(0.00001);
774     gr[i]->GetYaxis()->SetRangeUser(0.90,1.10);
775     gr[i]->SetLineColor(TColor::GetColor(color.c_str()));
776     gr[i]->SetMarkerColor(TColor::GetColor(color.c_str()));
777     if(dofit) {
778     gr[i]->Fit("pol1");
779     fit[i] = (TF1*)gr[i]->GetFunction("pol1");
780     fit[i]->SetLineColor(TColor::GetColor(color.c_str()));
781     }
782     leg->AddEntry(gr[i],name.c_str(),"l");
783     if(i==0) gr[i]->Draw("AP*");
784     else gr[i]->Draw("P");
785     }
786     leg->Draw();
787     DrawPrelim();
788    
789     CompleteSave(can,"ScenarioComparison");
790     delete can;
791    
792     }
793    
794     void ScenarioComparisonInclusive() {
795 buchmann 1.13 //dumbest way ever to do this. but ok we only need to do it once.
796 buchmann 1.9 TGraphAsymmErrors *gr[5];
797     TF1 *fit[5];
798     bool dofit=true;
799    
800     TCanvas *can = new TCanvas("can","can",500,500);
801    
802     TLegend *leg = new TLegend(0.2,0.2,0.5,0.45);
803     leg->SetBorderSize(0);
804     leg->SetFillColor(kWhite);
805    
806     for(int i=0;i<5;i++) {
807     gr[i] = new TGraphAsymmErrors();
808     double x[4] = {5,12.5,17.5,20};
809     double dxl[4] = {5,2.5,2.5,2.5};
810     double dxh[4] = {5,2.5,2.5,10};
811    
812     string name="";
813     string color="";
814    
815     if(i==0) {
816     double y[4] = {1.01103,1.01496,1.01486,1.02785};
817     double dyl[4] = {0.00235044,0.00225264,0.00313574,0.00541047};
818     double dyh[4] = {0.00235044,0.00225264,0.00313574,0.00541047};
819     name="30/30";
820     gr[i] = new TGraphAsymmErrors(4,x,y,dxl,dxh,dyl,dyh);
821     color="#a01d99";
822     }
823     if(i==1) {
824     double y[4] = {1.01103+0.001,1.01496+0.001,1.01486+0.001,1.02785+0.001};
825     double dyl[4] = {0.00235044,0.00225264,0.00313574,0.00541047};
826     double dyh[4] = {0.00235044,0.00225264,0.00313574,0.00541047};
827     gr[i] = new TGraphAsymmErrors(4,x,y,dxl,dxh,dyl,dyh);
828 buchmann 1.13 name="30/15";
829 buchmann 1.9 color="#2464bc";
830     }
831     if(i==2) {
832     double y[4] = {1.01575,1.02879,1.0325,1.05456};
833     double dyl[4] = {0.00288281,0.00278536,0.00367451,0.00590421};
834     double dyh[4] = {0.00288281,0.00278536,0.00367451,0.00590421};
835     gr[i] = new TGraphAsymmErrors(4,x,y,dxl,dxh,dyl,dyh);
836     name="15/15";
837     color="#f49230";
838     }
839     if(i==3) {
840     double y[4] = {1.01221,1.01922,1.02396,1.04677};
841     double dyl[4] = {0.00263821,0.00260853,0.00353673,0.00584728};
842     double dyh[4] = {0.00263821,0.00260853,0.00353673,0.00584728};
843     gr[i] = new TGraphAsymmErrors(4,x,y,dxl,dxh,dyl,dyh);
844     name="20/15";
845     color="#c72f3c";
846     }
847     if(i==4) {
848     double y[4] = {1.01221+0.001,1.01922+0.001,1.02396+0.001,1.04677+0.001};
849     double dyl[4] = {0.00263821,0.00260853,0.00353673,0.00584728};
850     double dyh[4] = {0.00263821,0.00260853,0.00353673,0.00584728};
851     gr[i] = new TGraphAsymmErrors(4,x,y,dxl,dxh,dyl,dyh);
852     name="20/20";
853     color="#7cb433";
854     }
855    
856    
857     gr[i]->SetTitle();
858     gr[i]->GetXaxis()->SetTitle("N(Vtx)");
859     gr[i]->GetYaxis()->SetTitle("<MPF>_{data}/<MPF>_{mc}");
860     gr[i]->GetXaxis()->CenterTitle();
861     gr[i]->GetYaxis()->CenterTitle();
862     gr[i]->SetName(name.c_str());
863     gr[i]->SetMarkerSize(0.00001);
864     gr[i]->GetYaxis()->SetRangeUser(0.90,1.10);
865     gr[i]->SetLineColor(TColor::GetColor(color.c_str()));
866     gr[i]->SetMarkerColor(TColor::GetColor(color.c_str()));
867     if(dofit) {
868     gr[i]->Fit("pol1");
869     fit[i] = (TF1*)gr[i]->GetFunction("pol1");
870     fit[i]->SetLineColor(TColor::GetColor(color.c_str()));
871     }
872     leg->AddEntry(gr[i],name.c_str(),"l");
873     if(i==0) gr[i]->Draw("AP*");
874     else gr[i]->Draw("P");
875     }
876     leg->Draw();
877     DrawPrelim();
878    
879     CompleteSave(can,"ScenarioComparisonInclusive");
880     delete can;
881    
882     }
883 buchmann 1.2
884 buchmann 1.7 void compare_selection(string identifier) {
885     bool recognized_scenario=false;
886    
887     cout << "Running with identifier " << identifier << endl;
888 buchmann 1.9 if(identifier=="Zb3010") {
889     LeadingB=TCut ("Zb3010_bTagProbCSVBP[0]>0.898");
890     EtaB=TCut("abs(Zb3010_pfJetGoodEta[0])<1.3");
891     PhiZcut=TCut("abs(Zb3010_pfJetDphiZ[0])>2.7");
892     recognized_scenario=true;
893     }
894    
895     if(identifier=="Zb3010") {
896     LeadingB=TCut ("Zb3010_bTagProbCSVBP[0]>0.898");
897     EtaB=TCut("abs(Zb3010_pfJetGoodEta[0])<1.3");
898     PhiZcut=TCut("abs(Zb3010_pfJetDphiZ[0])>2.7");
899 buchmann 1.7 recognized_scenario=true;
900     }
901    
902 buchmann 1.9 if(identifier=="Zb40") {
903     LeadingB=TCut ("Zb40_bTagProbCSVBP[0]>0.898");
904     EtaB=TCut("abs(Zb40_pfJetGoodEta[0])<1.3");
905     PhiZcut=TCut("abs(Zb40_pfJetDphiZ[0])>2.7");
906     recognized_scenario=true;
907     }
908    
909     if(identifier=="Zb301030") {
910     LeadingB=TCut ("Zb301030_bTagProbCSVBP[0]>0.898");
911     EtaB=TCut("abs(Zb301030_pfJetGoodEta[0])<1.3");
912     PhiZcut=TCut("abs(Zb301030_pfJetDphiZ[0])>2.7");
913 buchmann 1.7 recognized_scenario=true;
914     }
915    
916     if(identifier=="Zb1530") {
917 buchmann 1.9 LeadingB=TCut ("Zb1530_bTagProbCSVBP[0]>0.898");
918 buchmann 1.7 EtaB=TCut("abs(Zb1530_pfJetGoodEta[0])<1.3");
919     PhiZcut=TCut("abs(Zb1530_pfJetDphiZ[0])>2.7");
920     recognized_scenario=true;
921     }
922    
923 buchmann 1.9 if(identifier=="Zb301010") {
924     LeadingB=TCut ("Zb301010_bTagProbCSVBP[0]>0.898");
925     EtaB=TCut("abs(Zb301010_pfJetGoodEta[0])<1.3");
926     PhiZcut=TCut("abs(Zb301010_pfJetDphiZ[0])>2.7");
927     recognized_scenario=true;
928     }
929    
930     if(identifier=="Zb1510") {
931     LeadingB=TCut ("Zb1510_bTagProbCSVBP[0]>0.898");
932     EtaB=TCut("abs(Zb1510_pfJetGoodEta[0])<1.3");
933     PhiZcut=TCut("abs(Zb1510_pfJetDphiZ[0])>2.7");
934     recognized_scenario=true;
935     }
936    
937     if(identifier=="Zb301010") {
938     LeadingB=TCut ("Zb301010_bTagProbCSVBP[0]>0.898");
939     EtaB=TCut("abs(Zb301010_pfJetGoodEta[0])<1.3");
940     PhiZcut=TCut("abs(Zb301010_pfJetDphiZ[0])>2.7");
941     recognized_scenario=true;
942     }
943    
944     if(identifier=="ZbMikko") {
945     LeadingB=TCut ("ZbMikko_bTagProbCSVBP[0]>0.898 && ZbMikko__IsClean");
946     EtaB=TCut("abs(ZbMikko_pfJetGoodEta[0])<1.3 && ZbMikko__IsClean");
947     PhiZcut=TCut("abs(ZbMikko_pfJetDphiZ[0])>2.7 && ZbMikko__IsClean");
948 buchmann 1.7 recognized_scenario=true;
949     }
950    
951 buchmann 1.9 if(identifier=="Zb3010_SecEta3") {
952     LeadingB=TCut ("Zb3010_SecEta3_bTagProbCSVBP[0]>0.898");
953     EtaB=TCut("abs(Zb3010_SecEta3_pfJetGoodEta[0])<1.3");
954     PhiZcut=TCut("abs(Zb3010_SecEta3_pfJetDphiZ[0])>2.7");
955 buchmann 1.7 recognized_scenario=true;
956     }
957    
958 buchmann 1.9 if(identifier=="Zb3010_SecEta5") {
959     LeadingB=TCut ("Zb3010_SecEta5_bTagProbCSVBP[0]>0.898");
960     EtaB=TCut("abs(Zb3010_SecEta5_pfJetGoodEta[0])<1.3");
961     PhiZcut=TCut("abs(Zb3010_SecEta5_pfJetDphiZ[0])>2.7");
962 buchmann 1.7 recognized_scenario=true;
963     }
964    
965 buchmann 1.9 if(identifier=="Zb3010_p5Clean") {
966     LeadingB=TCut ("Zb3010_p5Clean_bTagProbCSVBP[0]>0.898 && Zb3010_p5Clean_IsClean");
967     EtaB=TCut("abs(Zb3010_p5Clean_pfJetGoodEta[0])<1.3 && Zb3010_p5Clean_IsClean");
968     PhiZcut=TCut("abs(Zb3010_p5Clean_pfJetDphiZ[0])>2.7 && Zb3010_p5Clean_IsClean");
969 buchmann 1.7 recognized_scenario=true;
970     }
971    
972 buchmann 1.9 if(identifier=="Zb3010CHS") {
973     LeadingB=TCut ("Zb3010CHS_bTagProbCSVBP[0]>0.898");
974     EtaB=TCut("abs(Zb3010CHS_pfJetGoodEta[0])<1.3");
975     PhiZcut=TCut("abs(Zb3010CHS_pfJetDphiZ[0])>2.7");
976 buchmann 1.7 recognized_scenario=true;
977     }
978    
979     assert(recognized_scenario);
980    
981     cout << "Cuts: " << endl;
982     cout << " " << (const char*) LeadingB << endl;
983     cout << " " << (const char*) EtaB << endl;
984     cout << " " << (const char*) PhiZcut << endl;
985    
986     cout << endl << endl;
987    
988 buchmann 1.9 // ScenarioComparison();
989     // ScenarioComparisonInclusive();
990    
991    
992     // DoPUStudy(identifier);
993    
994     // SpecialCutFlow(identifier);
995     // draw_kin_variable("pt",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000",17000,30,200,"Z p_{T}","Official/"+identifier+"/Zpt__"+identifier,1);
996     // draw_kin_variable(identifier+"_pfJetGoodPt[0]",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000",40,0,200,"Leading Jet Pt (B)","Official/"+identifier+"/LeadingJetPt__"+identifier,1);
997     // draw_kin_variable(identifier+"_pfJetGoodPt[1]",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000",40,0,200,"Sub-Leading Jet Pt [GeV]","Official/"+identifier+"/Jet2Pt__"+identifier,1);
998     // draw_kin_variable(identifier+"_pfJetGoodPt[1]/pt",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000",20,0,2,"p_{T}^{J2} / p_{T}^{Z}","Official/"+identifier+"/SubLeadingJetPt_Over_ZPt__"+identifier,1);
999     draw_kin_variable(identifier+"_alpha",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut("pt>30&&pt<50"),20,0,2,"#alpha","Official/"+identifier+"/alpha_pt_30_to_50__"+identifier,1);
1000     draw_kin_variable(identifier+"_alpha",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut("pt>50&&pt<75"),20,0,2,"#alpha","Official/"+identifier+"/alpha_50_to_75__"+identifier,1);
1001     draw_kin_variable(identifier+"_alpha",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut("pt>75&&pt<125"),20,0,2,"#alpha","Official/"+identifier+"/alpha_pt_75_to_125__"+identifier,1);
1002     draw_kin_variable(identifier+"_alpha",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut("pt>125&&pt<1000"),20,0,2,"#alpha","Official/"+identifier+"/alpha_pt_125_to_1000__"+identifier,1);
1003     draw_kin_variable(identifier+"_alpha",ZplusBsel&&LeadingB&&EtaB&&PhiZcut,20,0,2,"#alpha","Official/"+identifier+"/alpha_full__"+identifier,1);
1004     // draw_kin_variable(identifier+"_pfBJetDphiZ[0]",ZplusBsel&&LeadingB&&"pt>10&&pt<1000&&pfBJetDphiZ[0]>0",30,0,3.2,"#delta#phi (Z,b lead)","Official/"+identifier+"/DeltaPhiZBlead__"+identifier,0);
1005 buchmann 1.7
1006 buchmann 1.9 // draw_kin_variable("mpf",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut(((string)"pt>10&&pt<1000&&"+identifier+"_alpha<0.3").c_str()),20,0,2,"Z+b MPF","Official/"+identifier+"/mpf_alpha_smaller_0p3___"+identifier,0,"#alpha<0.3, 10 GeV < p_{T}^{Z} < 1000 GeV");
1007     // draw_kin_variable(identifier+"_pfJetGoodPt[0]/pt",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut(((string)"pt>10&&pt<1000&&"+identifier+"_alpha<0.3").c_str()),20,0,2,"p_{T} b-jet / p_{T} Z","Official/"+identifier+"/ptb_over_ptz___alpha_smaller_0p3______"+identifier,0,"#alpha<0.3, 10 GeV < p_{T}^{Z} < 1000 GeV");
1008 buchmann 1.7
1009    
1010 buchmann 1.9 // */
1011 buchmann 1.7
1012     }
1013    
1014 buchmann 1.13 void GetNumberEventsInsideOutsideAlphaWindow() {
1015    
1016     TCut cut = ZplusBsel&&LeadingB&&EtaB&&PhiZcut;
1017     TCut in = TCut("Zb3010_alpha<=0.3");
1018     TCut out = TCut("Zb3010_alpha>0.3");
1019    
1020     TH1F *data_IN = allsamples.Draw("data_IN", "mll",1,0,150, "nothing [GeV]", "events", cut&&in, data,luminosity);
1021     TH1F *data_OUT = allsamples.Draw("data_OUT", "mll",1,0,150, "nothing [GeV]", "events", cut&&out,data,luminosity);
1022     TH1F *mc_IN = allsamples.Draw("mc_IN", "mll",1,0,150, "nothing [GeV]", "events", cut&&in, mc ,luminosity);
1023     TH1F *mc_OUT = allsamples.Draw("mc_OUT", "mll",1,0,150, "nothing [GeV]", "events", cut&&out,mc ,luminosity);
1024    
1025     dout << " Data: In " << data_IN->Integral() << " vs out " << data_OUT->Integral() << endl;
1026     dout << " MC : In " << mc_IN->Integral() << " vs out " << mc_OUT->Integral() << endl;
1027     }
1028    
1029    
1030 buchmann 1.12 void do_basic_ZB_analysis(int do_inclusive) {
1031    
1032     //https://twiki.cern.ch/twiki/bin/viewauth/CMS/BtagPOG#2011_Data_and_MC
1033     //https://twiki.cern.ch/twiki/pub/CMS/BtagPOG/SFb-ttbar_payload.txt
1034 buchmann 1.9
1035 buchmann 1.12 if(do_inclusive) {
1036     LeadingB=TCut("Zb3010_bTagProbCSVBP[0]>=-10000");
1037     cutWeight=TCut("(weight*(weight<1000)*(is_data+(!is_data)*((id1==id2&&id1==0)*0.95+(id1==id2&&id1==1)*0.88+(id1!=id2)*0.92)))");
1038     write_info(__FUNCTION__,"Inclusive selection - no btag requirement, no btag efficiency/mistag correction (weight)");
1039     } else {
1040     LeadingB=TCut("Zb3010_bTagProbCSVBP[0]>0.898");
1041     cutWeight=TCut("(weight*(weight<1000)*(is_data+(!is_data)*((id1==id2&&id1==0)*0.95+(id1==id2&&id1==1)*0.88+(id1!=id2)*0.92)*Zb3010_BTagWgt))");
1042     write_info(__FUNCTION__,"Exclusive selection - switched on btag requirement and btag efficiency/mistag correction (weight)");
1043     }
1044    
1045     bool doquick=true;
1046 buchmann 1.9
1047 buchmann 1.2 TCanvas *zbcanvas = new TCanvas("zbcanvas","zbcanvas");
1048 buchmann 1.1
1049 buchmann 1.9 // compare_selection("Zb3010_p5Clean");
1050     // compare_selection("Zb3010CHS");
1051     // compare_selection("Zb3010");
1052     // compare_selection("Zb301030");
1053     // compare_selection("Zb1530");
1054     // compare_selection("Zb301010");
1055     // compare_selection("Zb3010_SecEta3");
1056     // compare_selection("Zb3010_SecEta5");
1057     // compare_selection("Zb1510");
1058     // compare_selection("Zb301010");
1059     // compare_selection("ZbMikko");
1060     // compare_selection("Zb3010");
1061     // compare_selection("Zb40");
1062 buchmann 1.8
1063 buchmann 1.11 if(!doquick) {
1064     print_all_b_yields();
1065     draw_mpf_vars();
1066     // DrawEvilCutFlow();
1067    
1068     draw_Zb_kin_vars();
1069    
1070     }
1071 buchmann 1.7
1072 buchmann 1.13 GetNumberEventsInsideOutsideAlphaWindow();
1073     // new_data_mc_agreement_2d();
1074 buchmann 1.1
1075 buchmann 1.12 write_error(__FUNCTION__,"Need to implement switch for inclusive result; inclusive result needs to have btag efficiency reweighting deactivated!");
1076 buchmann 1.2 delete zbcanvas;
1077 buchmann 1.1 }
1078 buchmann 1.7
1079    
1080    
1081    
1082     //*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/* DELETE EVERYTHING BELOW THIS LINE /*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*//
1083    
1084    
1085     const char* concatenate(string bla1, string bla2) {
1086     stringstream bla;
1087     bla << bla1 << bla2;
1088     return bla.str().c_str();
1089     }
1090    
1091     void SpecialCutFlow(string identifier) {
1092     stringstream MegaCut;
1093    
1094    
1095    
1096     MegaCut << " 1*(" << (const char*) (ZplusBsel&&TCut(concatenate(identifier,"_pfJetGoodNumBtag>0"))) << ")";
1097     MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB) << ")";
1098     MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB&&EtaB) << ")";
1099     MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB&&EtaB&&PhiZcut) << ")";
1100 buchmann 1.9 MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut(((string)"pt>10&&pt<1000").c_str())) << ")";
1101     MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut(((string)"pt>10&&pt<1000&&"+identifier+"_alpha<0.3").c_str())) << ")";
1102     MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut(((string)"pt>10&&pt<1000&&"+identifier+"_alpha<0.2").c_str())) << ")";
1103     MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut(((string)"pt>10&&pt<1000&&"+identifier+"_alpha<0.15").c_str())) << ")";
1104     MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut(((string)"pt>10&&pt<1000&&"+identifier+"_alpha<0.1").c_str())) << ")";
1105     MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut(((string)"pt>10&&pt<1000&&"+identifier+"_alpha<0.05").c_str())) << ")";
1106 buchmann 1.7 MegaCut << " ";
1107    
1108     TCanvas *can = new TCanvas("can","can");
1109     can->SetLogy(1);
1110 buchmann 1.9 TCut basecut(ZplusBsel);
1111 buchmann 1.7
1112    
1113    
1114     TLegend *leg = allsamples.allbglegend();
1115    
1116    
1117     TH1F *data = allsamples.Draw("data", MegaCut.str(),11,-0.5,10.5, "", "events", basecut,data,luminosity);
1118     THStack mcm = allsamples.DrawStack("mc", MegaCut.str(),11,-0.5,10.5, "", "events", basecut,mc,luminosity);
1119    
1120     float runningsum=0;
1121     for(int i=data->GetNbinsX();i>0;i--) {
1122     runningsum+=data->GetBinContent(i);
1123     data->SetBinContent(i,runningsum);
1124     }
1125    
1126     float minimum=data->GetMaximum();
1127    
1128     TH1F* h;
1129     TIter nextOF(mcm.GetHists());
1130    
1131     float mcsum=0;
1132     float zb=0;
1133     while ( h = (TH1F*)nextOF() ) {
1134     float runningsum=0;
1135     minimum=data->GetMaximum();//done deliberately at each step so get the minimum of the last (!) contribution
1136     for(int i=h->GetNbinsX();i>0;i--) {
1137     runningsum+=h->GetBinContent(i);
1138     h->SetBinContent(i,runningsum);
1139     if(runningsum<minimum&&runningsum>0) minimum=runningsum;
1140     }
1141     if(Contains(h->GetName(),"Z_b_")) {
1142     zb=h->GetBinContent(h->GetNbinsX());
1143     }
1144     mcsum+=h->GetBinContent(h->GetNbinsX());
1145     }
1146    
1147    
1148    
1149     data->SetMinimum(0.05*minimum);
1150     can->cd(1)->SetBottomMargin(0.25);
1151     data->GetXaxis()->SetBinLabel(1,"Pre-selection");
1152     data->GetXaxis()->SetBinLabel(2,"Contains b jet");
1153     data->GetXaxis()->SetBinLabel(3,"Leading jet is b-jet");
1154     data->GetXaxis()->SetBinLabel(4,"|#eta_{b}|<1.3 ");
1155     data->GetXaxis()->SetBinLabel(5,"|#delta#phi(b,Z)|>2.7");
1156     data->GetXaxis()->SetBinLabel(6,"10<p_{T}^{Z}<1000 GeV");
1157     data->GetXaxis()->SetBinLabel(7,"#alpha < 0.3");
1158     data->GetXaxis()->SetBinLabel(8,"#alpha < 0.2");
1159     data->GetXaxis()->SetBinLabel(9,"#alpha < 0.15");
1160     data->GetXaxis()->SetBinLabel(10,"#alpha < 0.1");
1161     data->GetXaxis()->SetBinLabel(11,"#alpha < 0.05");
1162     data->GetXaxis()->LabelsOption("v");
1163    
1164     stringstream purityinfo;
1165     purityinfo << "Purity: " << std::setprecision(3) << 100*zb/mcsum << " %";
1166     stringstream neventsinfo;
1167     neventsinfo << "Nevents: " << data->GetBinContent(data->GetNbinsX());
1168     TH1F *crap = new TH1F("crap","",1,0,1);
1169     crap->SetLineColor(kWhite);
1170     leg->AddEntry(crap,purityinfo.str().c_str(),"l");
1171     leg->AddEntry(crap,neventsinfo.str().c_str(),"l");
1172    
1173     data->Draw("e1");
1174     mcm.Draw("same");
1175     data->Draw("same,e1,axis");
1176     data->Draw("same,e1");
1177     // data->Draw("same,TEXT");
1178    
1179     leg->Draw();
1180    
1181     DrawPrelim();
1182    
1183     CompleteSave(can,"CutFlow___"+identifier);
1184    
1185     delete crap;
1186     delete data;
1187    
1188    
1189     }