ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/cbrown/Development/Plotting/Modules/ZbTools.C
(Generate patch)

Comparing UserCode/cbrown/Development/Plotting/Modules/ZbTools.C (file contents):
Revision 1.12 by buchmann, Mon Nov 19 16:21:52 2012 UTC vs.
Revision 1.21 by buchmann, Thu Jan 17 21:58:22 2013 UTC

# Line 21 | Line 21 | using namespace std;
21   using namespace PlottingSetup;
22  
23   namespace ZbData {
24 – vector<float> data_over_mc;
25 –
24   TCut ZplusBsel("pt1>20&&pt2>20&&mll>60&&mll<120&&id1==id2");
25 < TCut LeadingB("Zb3010_bTagProbCSVBP[0]>0.898");
26 < TCut EtaB("abs(Zb3010_pfJetGoodEta[0])<1.3");
27 < TCut PhiZcut("abs(Zb3010_pfJetDphiZ[0])>2.7");
25 > TCut LeadingB("ZbCHS3010_bTagProbCSVBP[0]>0.898");
26 > TCut EtaB("abs(ZbCHS3010_pfJetGoodEta[0])<1.3");
27 > TCut PhiZcut("abs(ZbCHS3010_pfJetDphiZ[0])>2.7");
28    
29   }
30  
31 + TCut STDWEIGHT;
32 +
33 + class ZbPointResult {
34 + public:
35 +  ZbPointResult(string var, float ptlow, float pthigh, float alphalow, float alphahigh, Value Data, Value MC);
36 +  float PtLow;
37 +  float PtHigh;
38 +  float AlphaLow;
39 +  float AlphaHigh;
40 +  string variable;
41 +  Value Data;
42 +  Value MC;
43 +  Value Ratio;
44 + };
45 +
46 + ZbPointResult::ZbPointResult(string var, float ptlow, float pthigh, float alphalow, float alphahigh, Value data, Value mc) :
47 +  variable(var), PtLow(ptlow), PtHigh(pthigh), AlphaLow(alphalow), AlphaHigh(alphahigh), Data(data), MC(mc) {
48 +    Ratio = data/mc;
49 + }
50 +  
51 + std::ostream &operator<<(std::ostream &ostr, ZbPointResult &b)
52 + {
53 +  ostr << b.PtLow << ";" << b.PtHigh << ";" << b.AlphaLow << ";" << b.AlphaHigh << ";" << b.Data.getValue() << ";" << b.Data.getError() << ";" << b.MC.getValue() << ";" << b.MC.getError() << ";" << b.Ratio.getValue() << ";" << b.Ratio.getError();
54 +  return ostr;
55 + }
56 +
57 +
58 + class ZbResultContainer {
59 + public:
60 +  string name;  
61 +  vector<ZbPointResult> MPF_Rabs_location;
62 +  Value MPF_Result_FaceValue;
63 +  Value Rabs_Result_FaceValue;
64 +  Value MPF_Result_Extrapolated;
65 +  Value Rabs_Result_Extrapolated;
66 +  
67 +  ZbResultContainer(string name);
68 +  ZbResultContainer(const ZbResultContainer &old);
69 + };
70 +
71 +
72 + ZbResultContainer::ZbResultContainer(string _name) {
73 +  name=_name;
74 +  cout<< "Result Container for \"" << name << "\" has been initialized" << endl;
75 +  
76 + }
77 +
78 + std::ostream &operator<<(std::ostream &ostr, ZbResultContainer &b)
79 + {
80 +  ostr << "Results for " << b.name << endl;
81 +  ostr << " MPF and Rabs locations: " << endl;
82 +  for(int i=0;i<b.MPF_Rabs_location.size();i++) ostr << "   " << b.MPF_Rabs_location[i] << endl;
83 +  ostr << " Results: " << endl;
84 +  ostr << "   Face Value: " << endl;
85 +  ostr << "     MPF " << b.MPF_Result_FaceValue << endl;
86 +  ostr << "     Rabs " << b.Rabs_Result_FaceValue << endl;
87 +  ostr << "   Extrapolated: " << endl;
88 +  ostr << "     MPF " << b.MPF_Result_Extrapolated << endl;
89 +  ostr << "             Rabs " << b.Rabs_Result_Extrapolated << endl;
90 +  
91 +  return ostr;
92 + }
93 +
94 +
95 +
96   using namespace ZbData;
97  
98  
# Line 39 | Line 102 | void SpecialCutFlow(string identifier);
102   //*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*            AND THIS LINE             /*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*//
103  
104  
105 < void print_yield(TCut cut) {
105 > class CutYields{
106 > public:
107 >  float data;
108 >  float Zb;
109 >  float Zc;
110 >  float Zl;
111 >  float tt;
112 >  float SiTo;
113 >  float Dibosons;
114 >  float WJets;
115 >  float Other;
116 >  
117 >  CutYields();
118 > };
119 >
120 > CutYields::CutYields() {
121 >  this->data=0;
122 >  this->Zb=0;
123 >  this->Zc=0;
124 >  this->Zl=0;
125 >  this->tt=0;
126 >  this->SiTo=0;
127 >  this->Dibosons=0;
128 >  this->WJets=0;
129 >  this->Other=0;
130 > }
131 >
132 > std::ostream &operator<<(std::ostream &ostr, CutYields y)
133 > {
134 >  return ostr << y.data << ";" << y.Zb << ";" << y.Zc << ";" << y.Zl << ";" << y.tt << ";" << y.SiTo << ";" << y.Dibosons << ";" << y.WJets << ";" << endl;
135 > }
136 >
137 >
138 > CutYields print_yield(TCut cut) {
139 >  
140 >  CutYields Yield;
141    TH1F *data  = allsamples.Draw("data",    "mll",1,0,10000000, "m_{ll} [GeV]", "events", cut,data,luminosity);
142    dout << "  Yields for " << (const char*) cut << endl;
143 <  dout << "     data: " << data->Integral() << endl;
143 > //  dout << "     data: " << data->Integral() << endl;
144 >  
145 >  Yield.data=data->Integral();
146    THStack mcm = allsamples.DrawStack("mc", "mll",1,0,10000000, "m_{ll} [GeV]", "events", cut,mc,luminosity);
147    int nhists = mcm.GetHists()->GetEntries();
48 –  dout << "Going to loop over " << nhists << " histograms " << endl;
148    TFile *test = new TFile("test.root","RECREATE");
149    mcm.Write();
150    test->Close();
151    for (int i = 0; i < nhists; i++) {
152      TH1* hist = static_cast<TH1*>(mcm.GetHists()->At(i));
153 <    cout << "     " << hist->GetName() << ": " << hist->Integral() << endl;
153 >    string name=hist->GetName();
154 >    if(Contains(name,"t_bar_t")) Yield.tt+=hist->Integral();
155 >    if(Contains(name,"W_Jets")) Yield.WJets+=hist->Integral();
156 >    if(Contains(name,"Single_top")) Yield.SiTo+=hist->Integral();
157 >    if(Contains(name,"Dibosons")) Yield.Dibosons+=hist->Integral();
158 >    if(Contains(name,"Z_l_")) Yield.Zl+=hist->Integral();
159 >    if(Contains(name,"Z_c_")) Yield.Zc+=hist->Integral();
160 >    if(Contains(name,"Z_b_")) Yield.Zb+=hist->Integral();
161    }
162 <    
162 >  
163 >  cout << Yield << endl;
164    cout << endl << endl;
165    
166    delete data;
167 +  CleanLegends();
168 +  DeleteStack(mcm);
169 +
170 +  return Yield;
171   }
172  
173   void draw_kin_variable(string variable, TCut cut, int nbins, float min, float max, string xlabel, string saveas, bool dologscale, string speciallabel="") {
63 –  
174    TCanvas *ckin = new TCanvas("ckin","kin variable canvas");
175 <  ckin->SetLogy(dologscale);
175 >  if(dologscale) ckin->SetLogy(1);
176    TH1F *datahisto = allsamples.Draw("datahisto",variable,nbins,min,max, xlabel, "events",cut,data,luminosity);
177    datahisto->SetMarkerSize(DataMarkerSize);
178    THStack mcstack = allsamples.DrawStack("mcstack",variable,nbins,min,max, xlabel, "events",cut,mc,luminosity);
179    if(dologscale) datahisto->SetMaximum(3*datahisto->GetMaximum());
180    else datahisto->SetMaximum(1.5*datahisto->GetMaximum());
181 <  TText *speciall;
182 <  if(speciallabel!="") {
183 <    speciall = new TText(0.2,0.8,speciallabel.c_str());
184 <    speciall->Draw();
181 >  
182 >  float SumMC=0.0;
183 >  TIter nextH(mcstack.GetHists());
184 >  TH1F* h;
185 >  while ( h = (TH1F*)nextH() ) {
186 >    SumMC += h->Integral();
187    }
188 +  float scale = datahisto->Integral()/SumMC;
189 +
190 +  // Re-scale stack to data integral
191 +  nextH = TIter(mcstack.GetHists());
192 +  
193 +  while ( h = (TH1F*)nextH() ) {
194 +    h->Scale(scale);
195 +  }
196 +  
197 +  TText *speciall;
198    datahisto->Draw("e1");
199    ckin->Update();
200 <  mcstack.Draw("same");
200 >  mcstack.Draw("histo,same");
201    datahisto->Draw("same,e1");
202    TLegend *kinleg = allsamples.allbglegend();
203 +
204    kinleg->Draw();
205  
206    TPad *kinpad = new TPad("kinpad","kinpad",0,0,1,1);
207    kinpad->cd();
208 <  kinpad->SetLogy(dologscale);
208 >  if(dologscale) kinpad->SetLogy(1);
209    datahisto->Draw("e1");
210 <  mcstack.Draw("same");
210 >  mcstack.Draw("histo,same");
211    datahisto->Draw("same,e1");
212    datahisto->Draw("same,axis");
213    kinleg->Draw();
214    DrawPrelim();
215    saveas="kin/"+saveas;
216 <  save_with_ratio(datahisto,mcstack,kinpad->cd(),saveas);
94 <
95 <  
216 >  kinpad->cd();
217    
218 +  if(speciallabel!="") {
219 +    speciall = write_title(speciallabel);
220 +    speciall->SetX(0.18);
221 +    speciall->SetTextAlign(11);
222 +    speciall->SetTextSize(0.03);
223 +    speciall->SetY(0.85);
224 +    speciall->Draw();
225 +  }
226    
227 <  ZbData::data_over_mc.push_back(datahisto->Integral()/CollapseStack(mcstack)->Integral());
227 >  save_with_ratio(datahisto,mcstack,kinpad->cd(),saveas);
228    
229    datahisto->Delete();
230 +  delete kinleg;
231    delete ckin;
232 +  CleanLegends();
233 +  DeleteStack(mcstack);
234   }
235 <  
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 <  
235 >
236   void print_all_b_yields() {
237 +  vector<CutYields> yields;
238    cout << "Basic selection with a b jet" << endl;
239 <  print_yield(ZplusBsel&&"Zb3010_pfJetGoodNumBtag>0");
239 >  yields.push_back(print_yield(ZplusBsel&&"ZbCHS3010_pfJetGoodNumBtag>0"));
240    cout << "Leading jet is a b-jet " << endl;
241 <  print_yield(ZplusBsel&&LeadingB);
241 >  yields.push_back(print_yield(ZplusBsel&&LeadingB));
242    cout << "|#eta_{b}|<1.3 " << endl;
243 <  print_yield(ZplusBsel&&LeadingB&&EtaB);
243 >  yields.push_back(print_yield(ZplusBsel&&LeadingB&&EtaB));
244    cout << "|#delta#phi(b<Z)|>2.7" << endl;
245 <  print_yield(ZplusBsel&&LeadingB&&EtaB&&PhiZcut);
246 <  cout << "10<ptZ<1000 GeV" << endl;
247 <  print_yield(ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000");
245 >  yields.push_back(print_yield(ZplusBsel&&LeadingB&&EtaB&&PhiZcut));
246 >  cout << "30<ptZ<1000 GeV" << endl;
247 >  yields.push_back(print_yield(ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>30&&pt<1000"));
248 >  cout << "#alpha < 0.35" << endl;
249 >  yields.push_back(print_yield(ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>30&&pt<1000"&&"(ZbCHS3010_alpha)<0.35"));
250    cout << "#alpha < 0.3" << endl;
251 <  print_yield(ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"(Zb3010_alpha)<0.3");
252 <  cout << "#alpha < 0.2" << endl;
253 <  print_yield(ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"(Zb3010_alpha)<0.2");
251 >  yields.push_back(print_yield(ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>30&&pt<1000"&&"(ZbCHS3010_alpha)<0.3"));
252 > /*  cout << "#alpha < 0.2" << endl;
253 >  yields.push_back(print_yield(ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>30&&pt<1000"&&"(ZbCHS3010_alpha)<0.2"));
254    cout << "#alpha < 0.15" << endl;
255 <  print_yield(ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"(Zb3010_alpha)<0.15");
255 >  yields.push_back(print_yield(ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>30&&pt<1000"&&"(ZbCHS3010_alpha)<0.15"));
256    cout << "#alpha < 0.1" << endl;
257 <  print_yield(ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"(Zb3010_alpha)<0.1");
258 < //  cout << "#alpha < 0.05" << endl;
259 < //  print_yield(ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"(Zb3010_alpha)<0.05");
257 >  yields.push_back(print_yield(ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>30&&pt<1000"&&"(ZbCHS3010_alpha)<0.1"));*/
258 >
259 >  for(int iy=0;iy<yields.size();iy++) {
260 >    cout << yields[iy] << endl;
261 >  }
262   }  
263  
264   void DrawEvilCutFlow() {
265    stringstream MegaCut;
266    
267 <  MegaCut << "   1*(" << (const char*) (ZplusBsel&&TCut("Zb3010_pfJetGoodNumBtag>0")) << ")";
267 >  MegaCut << "   1*(" << (const char*) (ZplusBsel&&TCut("ZbCHS3010_pfJetGoodNumBtag>0")) << ")";
268    MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB) << ")";
269    MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB&&EtaB) << ")";
270    MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB&&EtaB&&PhiZcut) << ")";
271 <  MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut("pt>10&&pt<1000")) << ")";
272 <  MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut("pt>10&&pt<1000&&Zb3010_alpha<0.3")) << ")";
273 <  MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut("pt>10&&pt<1000&&Zb3010_alpha<0.2")) << ")";
274 <  MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut("pt>10&&pt<1000&&Zb3010_alpha<0.15")) << ")";
275 <  MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut("pt>10&&pt<1000&&Zb3010_alpha<0.1")) << ")";
276 <  MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut("pt>10&&pt<1000&&Zb3010_alpha<0.05")) << ")";
271 >  MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut("pt>30&&pt<1000")) << ")";
272 >  MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut("pt>30&&pt<1000&&ZbCHS3010_alpha<0.3")) << ")";
273 >  MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut("pt>30&&pt<1000&&ZbCHS3010_alpha<0.2")) << ")";
274 >  MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut("pt>30&&pt<1000&&ZbCHS3010_alpha<0.15")) << ")";
275 >  MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut("pt>30&&pt<1000&&ZbCHS3010_alpha<0.1")) << ")";
276 >  MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut("pt>30&&pt<1000&&ZbCHS3010_alpha<0.05")) << ")";
277    MegaCut << " ";
278    
279    TCanvas *can = new TCanvas("can","can");
# Line 157 | Line 286 | void DrawEvilCutFlow() {
286    
287    
288    TH1F *data  = allsamples.Draw("data",    MegaCut.str(),11,-0.5,10.5, "", "events", basecut,data,luminosity);
289 <   THStack mcm = allsamples.DrawStack("mc", MegaCut.str(),11,-0.5,10.5, "", "events", basecut,mc,luminosity);
289 >  THStack mcm = allsamples.DrawStack("mc", MegaCut.str(),11,-0.5,10.5, "", "events", basecut,mc,luminosity);
290  
291    float runningsum=0;
292    for(int i=data->GetNbinsX();i>0;i--) {
# Line 188 | Line 317 | void DrawEvilCutFlow() {
317    data->GetXaxis()->SetBinLabel(3,"Leading jet is b-jet");
318    data->GetXaxis()->SetBinLabel(4,"|#eta_{b}|<1.3 ");
319    data->GetXaxis()->SetBinLabel(5,"|#delta#phi(b,Z)|>2.7");
320 <  data->GetXaxis()->SetBinLabel(6,"10<p_{T}^{Z}<1000 GeV");
320 >  data->GetXaxis()->SetBinLabel(6,"30<p_{T}^{Z}<1000 GeV");
321    data->GetXaxis()->SetBinLabel(7,"#alpha < 0.3");
322    data->GetXaxis()->SetBinLabel(8,"#alpha < 0.2");
323    data->GetXaxis()->SetBinLabel(9,"#alpha < 0.15");
# Line 197 | Line 326 | void DrawEvilCutFlow() {
326    data->GetXaxis()->LabelsOption("v");
327    
328    data->Draw("e1");
329 <  mcm.Draw("same");
329 >  mcm.Draw("histo,same");
330    data->Draw("same,e1,axis");
331    data->Draw("same,e1");
332   //  data->Draw("same,TEXT");
# Line 207 | Line 336 | void DrawEvilCutFlow() {
336    DrawPrelim();
337    
338    CompleteSave(can,"CutFlow");
339 +  CleanLegends();
340 +  DeleteStack(mcm);
341 +
342   }
343  
344   void draw_Zb_kin_vars() {
345 <   draw_kin_variable("pt",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000",25,0,200,"Z p_{T}","Official/Zpt",1);
346 <   draw_kin_variable("Zb3010_pfJetGoodPt[1]",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000",25,0,200,"Sub-Leading Jet Pt","Official/Jet2Pt",1);
347 <   draw_kin_variable("Zb3010_pfJetGoodPt[0]",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000",25,0,200,"Leading Jet Pt (B)","Official/LeadingJetPt",1);
348 <   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);
349 <   draw_kin_variable("Zb3010_alpha",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>30&&pt<50",   40,0,2,"#alpha","Official/alpha_pt_30_to_50",1);
350 <   draw_kin_variable("Zb3010_alpha",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>50&&pt<75",  40,0,2,"#alpha","Official/alpha_50_to_75",1);
351 <   draw_kin_variable("Zb3010_alpha",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>75&&pt<125", 40,0,2,"#alpha","Official/alpha_pt_75_to_125",1);
352 <   draw_kin_variable("Zb3010_alpha",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>125&&pt<1000",40,0,2,"#alpha","Official/alpha_pt_125_to_1000",1);
353 <    draw_kin_variable("Zb3010_alpha",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>30&&pt<1000",40,0,2,"#alpha","Official/alpha_pt_30_to_1000",1);
354 <  
355 <  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);
356 <  
225 < /*
226 <  draw_kin_variable("pt1",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000",20,0,200,"p_{T}^{l1}","Lep/pt1",1);
227 <  draw_kin_variable("pt1",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000&&id1==0",20,0,200,"p_{T}^{e1}","Lep/pt1_e",1);
228 <  draw_kin_variable("pt1",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000&&id1==1",20,0,200,"p_{T}^{#mu1}","Lep/pt1_m",1);
229 <  draw_kin_variable("pt2",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000",20,0,200,"p_{T}^{l2}","Lep/pt2",1);
230 <  draw_kin_variable("pt2",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000&&id2==0",20,0,200,"p_{T}^{e2}","Lep/pt2_e",1);
231 <  draw_kin_variable("pt2",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000&&id2==1",20,0,200,"p_{T}^{#mu2}","Lep/pt2_m",1);
232 <  draw_kin_variable("eta1",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000",20,-3.1,3.1,"#eta^{l1}","Lep/eta1",1);
233 <  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);
234 <  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);
235 <  draw_kin_variable("eta2",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000",20,-3.1,3.1,"#eta^{l1}","Lep/eta2",1);
236 <  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);
237 <  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);
238 <  draw_kin_variable("numVtx",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000",30,0,30,"N(Vtx)","Lep/nVtx",1);
239 <  */
345 >  
346 >   draw_kin_variable("pt",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut("pt>30&&pt<1000"),25,0,200,"Z p_{T}","Official/Zpt",1);
347 >   draw_kin_variable("ZbCHS3010_pfJetGoodPt[1]",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut("pt>30&&pt<1000"),25,0,200,"Sub-Leading Jet Pt","Official/Jet2Pt",1);
348 >   draw_kin_variable("ZbCHS3010_pfJetGoodPt[0]",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut("pt>30&&pt<1000"),25,0,200,"Leading Jet Pt (B)","Official/LeadingJetPt",1);
349 >   draw_kin_variable("ZbCHS3010_pfJetGoodPt[1]/pt",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut("pt>30&&pt<1000"),20,0,2,"p_{T}^{J2} / p_{T}^{Z}","Official/SubLeadingJetPt_Over_ZPt",1);
350 >   draw_kin_variable("ZbCHS3010_alpha",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut("pt>30&&pt<50"),   40,0,2,"#alpha","Official/alpha_pt_30_to_50",1);
351 >   draw_kin_variable("ZbCHS3010_alpha",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut("pt>50&&pt<75"),  40,0,2,"#alpha","Official/alpha_50_to_75",1);
352 >   draw_kin_variable("ZbCHS3010_alpha",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut("pt>75&&pt<125"), 40,0,2,"#alpha","Official/alpha_pt_75_to_125",1);
353 >   draw_kin_variable("ZbCHS3010_alpha",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut("pt>125&&pt<1000"),40,0,2,"#alpha","Official/alpha_pt_125_to_1000",1);
354 >   draw_kin_variable("ZbCHS3010_alpha",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut("pt>30&&pt<1000"),40,0,2,"#alpha","Official/alpha_pt_30_to_1000",1);
355 >  
356 >   draw_kin_variable("ZbCHS3010_pfBJetDphiZ[0]",ZplusBsel&&LeadingB&&TCut("pt>30&&pt<1000&&pfBJetDphiZ[0]>0"),30,0,3.2,"#delta#phi (Z,b lead)","DeltaPhiZBlead",0);
357   }
358  
359   void draw_mpf_vars() {
360 <  ZbData::data_over_mc.clear();
361 <  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");
362 <  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");
363 <  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");
364 <  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");
365 <  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");
366 <  
367 <  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");
368 <  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");
369 <  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");
370 <  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");
254 <  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");
360 >  draw_kin_variable("mpf",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut("pt>30&&pt<1000&&ZbCHS3010_alpha<0.3"), 20,0.0,2.0,"Z+b MPF","mpf_alpha_smaller_0p3",0,"#alpha<0.3, 30 GeV < p_{T}^{Z} < 1000 GeV");
361 >  draw_kin_variable("mpf",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut("pt>30&&pt<1000&&ZbCHS3010_alpha<0.2"), 20,0.0,2.0,"Z+b MPF","mpf_alpha_smaller_0p2",0,"#alpha<0.2, 30 GeV < p_{T}^{Z} < 1000 GeV");
362 >  draw_kin_variable("mpf",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut("pt>30&&pt<1000&&ZbCHS3010_alpha<0.15"),20,0.0,2.0,"Z+b MPF","mpf_alpha_smaller_0p15",0,"#alpha<0.15, 30 GeV < p_{T}^{Z} < 1000 GeV");
363 >  draw_kin_variable("mpf",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut("pt>30&&pt<1000&&ZbCHS3010_alpha<0.1"), 20,0.0,2.0,"Z+b MPF","mpf_alpha_smaller_0p1",0,"#alpha<0.1, 30 GeV < p_{T}^{Z} < 1000 GeV");
364 >  draw_kin_variable("mpf",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut("pt>30&&pt<1000&&ZbCHS3010_alpha<0.05"),20,0.0,2.0,"Z+b MPF","mpf_alpha_smaller_0p05",0,"#alpha<0.05, 30 GeV < p_{T}^{Z} < 1000 GeV");
365 >  
366 >  draw_kin_variable("ZbCHS3010_pfJetGoodPt[0]/pt",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut("pt>30&&pt<1000&&ZbCHS3010_alpha<0.3"), 20,0.0,2.0,"p_{T} b-jet / p_{T} Z","ptb_over_ptz___alpha_smaller_0p3",    0,"#alpha<0.3, 30 GeV < p_{T}^{Z} < 1000 GeV");
367 >  draw_kin_variable("ZbCHS3010_pfJetGoodPt[0]/pt",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut("pt>30&&pt<1000&&ZbCHS3010_alpha<0.2"), 20,0.0,2.0,"p_{T} b-jet / p_{T} Z","ptb_over_ptz___alpha_smaller_0p2",    0,"#alpha<0.2, 30 GeV < p_{T}^{Z} < 1000 GeV");
368 >  draw_kin_variable("ZbCHS3010_pfJetGoodPt[0]/pt",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut("pt>30&&pt<1000&&ZbCHS3010_alpha<0.15"),20,0.0,2.0,"p_{T} b-jet / p_{T} Z","ptb_over_ptz___alpha_smaller_0p15",   0,"#alpha<0.15, 30 GeV < p_{T}^{Z} < 1000 GeV");
369 >  draw_kin_variable("ZbCHS3010_pfJetGoodPt[0]/pt",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut("pt>30&&pt<1000&&ZbCHS3010_alpha<0.1"), 20,0.0,2.0,"p_{T} b-jet / p_{T} Z","ptb_over_ptz___alpha_smaller_0p1",    0,"#alpha<0.1, 30 GeV < p_{T}^{Z} < 1000 GeV");
370 >  draw_kin_variable("ZbCHS3010_pfJetGoodPt[0]/pt",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut("pt>30&&pt<1000&&ZbCHS3010_alpha<0.05"),20,0.0,2.0,"p_{T} b-jet / p_{T} Z","ptb_over_ptz__mpf_alpha_smaller_0p05",0,"#alpha<0.05, 30 GeV < p_{T}^{Z} < 1000 GeV");
371   }  
372  
373 < /*double GetMedian(TH1F *histo) {
258 <    int numBins = histo->GetXaxis()->GetNbins();
259 <    Double_t x[numBins];
260 <    Double_t y[numBins];
261 <    for (int i = 1; i <= numBins; i++) {
262 <        x[i] = histo->GetBinCenter(i);
263 <        y[i] = histo->GetBinContent(i);
264 <    }
265 <    return TMath::Median(numBins, &x, &y);
266 < }*/
267 <
268 < Value get_Zb_data_over_mc(string variable, string variable2, TCut cut, string saveas) {
373 > Value get_Zb_data_over_mc(string ContainerName, ZbResultContainer &res, string variable, string variable2, TCut cut, string saveas, float pt1, float pt2, float a1, float a2) {
374   //  write_warning(__FUNCTION__,"Debugging this function, therefore always returning 3.1415+/-1.010101"); return Value(3.1415,1.010101);
375    TCanvas *cn = new TCanvas("cn","cn");
376    string varname="MPF";
377 +  
378 +  bool DoFit=false;
379 +  bool UseFit=false;
380 +  
381 +  
382    if(Contains(variable,"pfJetGood")) varname="R_{abs}";
383 <  TH1F *hdata  = allsamples.Draw("hdata",variable,1200,0,30, varname, "events",  cut,data,luminosity);
384 <  TH1F *hmc    = allsamples.Draw("hmc" , variable,1200,0,30, varname, "events", cut, mc,  luminosity);
383 > //  TH1F *hdata  = allsamples.Draw("hdata",variable,1200,0,30, varname, "events",  cut,data,luminosity);
384 > //  TH1F *hmc    = allsamples.Draw("hmc" , variable,1200,0,30, varname, "events", cut, mc,  luminosity);
385 >  TH1F *hdata  = allsamples.Draw("hdata",variable,1200,0.5,1.5, varname, "events",  cut,data,luminosity); // making it more precise
386 >  TH1F *hmc    = allsamples.Draw("hmc" , variable,1200,0.5,1.5, varname, "events", cut, mc,  luminosity); // making it more precise
387 >  hmc->SetLineColor(kBlue);
388    
389    float a=hdata->GetMean();
390    float b=hmc->GetMean();
# Line 280 | Line 393 | Value get_Zb_data_over_mc(string variabl
393    float factor = a / b;
394    float error = TMath::Sqrt(  (1/(b*b)) * (da*da) + ((a*a)/(b*b))*db*db);
395    
396 +  TF1 *gaus_data = new TF1("gaus_data","gaus",0.7,1.3);
397 +  TF1 *gaus_mc = new TF1("gaus_mc","gaus",0.7,1.3);
398 +  if(DoFit) {
399 +    gaus_data->SetParameter(0,hdata->GetMaximum());
400 +    gaus_data->SetParameter(2,hdata->GetMean());
401 +    gaus_data->SetParameter(1,hdata->GetRMS());
402 +    hdata->Fit(gaus_data);
403 +    
404 +    gaus_mc->SetParameter(0,hmc->GetMaximum());
405 +    gaus_mc->SetParameter(1,hmc->GetMean());
406 +    gaus_mc->SetParameter(2,hmc->GetRMS());
407 +    hmc->Fit(gaus_mc);
408 +  }
409 +  
410    cn->cd();
411    hdata->GetYaxis()->SetTitle("events / 0.1");
412 <  hdata->Rebin(4);
413 <  hmc->Rebin(4);
412 >  hdata->SetMarkerColor(kRed);
413 >  if(DoFit) {
414 >    gaus_mc->SetLineColor(kBlue);
415 >    gaus_data->SetLineColor(kRed);
416 >  }
417 >  
418    hdata->GetXaxis()->SetRangeUser(0,2);
419 +  hdata->GetYaxis()->SetTitleOffset(1.3);
420    hdata->Draw("e1");
421    hmc->Draw("histo,same");
422    hdata->Draw("e1,same");
423 +  if(DoFit) {
424 +    gaus_data->Draw("same");
425 +    gaus_mc->Draw("same");
426 +  }
427 +  
428    stringstream summary;
429 <  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;
430 <  summary << " #pm " << error << "}}}";
429 >  summary << "#splitline{#bf{" << varname << "}}{#splitline{" << pt1 << " < p_{T}^{Z} < " << pt2 << ", " << a1 << " #leq #alpha < " << a2 << "}{#splitline{data: " << std::setprecision(4) << a << " #pm " << std::setprecision(4) << da << "}{";
430 >  summary << "#splitline{mc:   " << std::setprecision(4) << b << " #pm " << std::setprecision(4) << db << "}{#splitline{ratio: " << factor;
431 >  summary << " #pm " << error;
432 >  if(DoFit) {
433 >    summary << "}{#splitline{}{#splitline{data fit: " << std::setprecision(4) << gaus_data->GetParameter(1) << "+/-";
434 >    summary << std::setprecision(4) << gaus_data->GetParError(1) << "}{#splitline{";
435 >    summary << "mc fit:" << std::setprecision(4) << gaus_mc->GetParameter(1) << "#pm" << std::setprecision(4) << gaus_mc->GetParError(1) << "}{ratio: ";
436 >    summary << gaus_data->GetParameter(1)/gaus_mc->GetParameter(1) << "#pm" << gaus_data->GetParameter(1)/gaus_mc->GetParameter(1) * sqrt((gaus_data->GetParError(1) * gaus_data->GetParError(1))/(gaus_data->GetParameter(1)*gaus_data->GetParameter(1)) + (gaus_mc->GetParError(1) * gaus_mc->GetParError(1))/(gaus_mc->GetParameter(1)*gaus_mc->GetParameter(1)) ) << "}}}}";
437 >    summary << "}}}}";
438 >  } else {
439 >    summary << "}{}}}}}";
440 >  }
441 >  
442 >  if(UseFit) {
443 >    write_error(__FUNCTION__,"Requested using fit result for final result - not implemented yet!");
444 >    assert(0);
445 >  }
446 >    
447 >  
448    TText *infobox = write_title(summary.str());
449    infobox->SetX(0.75);
450    infobox->SetTextSize(0.03);
# Line 299 | Line 453 | Value get_Zb_data_over_mc(string variabl
453      
454      
455    DrawPrelim();
456 <  CompleteSave(cn,"ResponsePlots/"+saveas);
456 >  CompleteSave(cn,ContainerName+"/ResponsePlots/"+saveas);
457    
458 <  cout << "Data;" << a << ";"<<da<<";MC"<<b<<";"<<db<<endl;
458 >  res.MPF_Rabs_location.push_back(ZbPointResult(variable,pt1,pt2,a1,a2,Value(a,da),Value(b,db)));
459 >  cout << "Data;" << a << ";"<<da<<";MC;"<<b<<";"<<db<<endl;
460    delete cn;
461    delete hdata;
462    delete hmc;
463  
309 – //    cout << a << "+/-" << da << "        " << b << "=/-" << db << endl;
310 –
464    return Value(factor,error);
465   }
466    
467 + int nFakes=1;
468 +
469 + float PurityInclusive=0.05;
470 + float PurityLoose=0.20;
471 + float PurityMedium=0.54;
472 + float PurityTight=0.82;
473 +
474 +
475 + void DeterminePurity(TCut basecut, int nptcuts, float ptcuts[], string AlphaVariable, string ContainerName, float FaceValueToBeTakenAt) {
476 +  
477 +  bool ComputePurity=false;
478 +  
479 +  if(ContainerName=="inclusive") ComputePurity=true;  
480 +  if(ContainerName=="loose") ComputePurity=true;
481 +  if(ContainerName=="medium") ComputePurity=true;
482 +  if(ContainerName=="reference") ComputePurity=true;
483 +  
484 +  if(!ComputePurity) return;
485 +  
486 +  stringstream specialcut;
487 +  specialcut << "((pt>30&& pt< 1000) && (" << AlphaVariable << "<" << FaceValueToBeTakenAt << "))";
488 +  
489 +  THStack mcm    = allsamples.DrawStack("mcm" , "mpf",1200,0,30, "MPF", "events", TCut(basecut && specialcut.str().c_str()), mc,  luminosity);
490 +  int nhists = mcm.GetHists()->GetEntries();
491 +  
492 +  float total=0;
493 +  float Zb=0;
494 +  
495 +  for (int i = 0; i < nhists; i++) {
496 +    TH1* hist = static_cast<TH1*>(mcm.GetHists()->At(i));
497 +    string name=hist->GetName();
498 +    if(Contains(name,"t_bar_t")) total+=hist->Integral();
499 +    if(Contains(name,"W_Jets")) total+=hist->Integral();
500 +    if(Contains(name,"Single_top")) total+=hist->Integral();
501 +    if(Contains(name,"Dibosons")) total+=hist->Integral();
502 +    if(Contains(name,"Z_l_")) total+=hist->Integral();
503 +    if(Contains(name,"Z_c_")) total+=hist->Integral();
504 +    if(Contains(name,"Z_b_")) {
505 +      total+=hist->Integral();
506 +      Zb+=hist->Integral();
507 +    }
508 +  }
509    
510 +  float purity=Zb/total;
511 +  dout << "Determined a purity of " << 100*purity << " % for the container " << ContainerName << endl;
512 +  
513 +  if(ContainerName=="inclusive") PurityInclusive=purity;
514 +  if(ContainerName=="loose") PurityLoose=purity;
515 +  if(ContainerName=="medium") PurityMedium=purity;
516 +  if(ContainerName=="reference") PurityTight=purity;
517 +  
518 +  DeleteStack(mcm);
519 + }
520  
521 < void new_data_mc_agreement_2d() {
521 >
522 > ZbResultContainer new_data_mc_agreement_2d(bool gofast, string AlphaVariable, string ContainerName, TCut BTagCut, TCut BTagWgt, float FaceValueToBeTakenAt) {
523 >  
524 >  
525 >  if(gofast) write_info(__FUNCTION__,"Fast mode activated for "+ContainerName);
526 >  ZbResultContainer results(ContainerName);
527 >  
528 >  cutWeight=STDWEIGHT*BTagWgt;
529 >  LeadingB=BTagCut;
530 >  
531 >  
532    gStyle->SetOptFit(0);
533    TCut basecut(ZplusBsel&&LeadingB&&EtaB&&PhiZcut);
534 <  const int nalphacuts=5;
535 <  float alphacuts[nalphacuts] = {0.1,0.15,0.2,0.25,0.3};
536 <  int PointThree=-1;
322 <  for(int i=0;i<nalphacuts;i++) {
323 <    if((alphacuts[i]-0.3)<0.01) PointThree=i-1;
324 <  }
325 <  if(PointThree<0) {
326 <    write_warning(__FUNCTION__,"alpha=0.3 MUST BE ONE OF THE ALPHA VALUES!");
327 <    assert(PointThree<0);
328 <  }
534 >  const int nalphacuts=6;
535 >  float alphacuts[nalphacuts] = {0.1,0.15,0.2,0.25,0.3,0.35};
536 >  
537    
538    const int nptcuts=5;
539    float ptcuts[nptcuts]={30,50,75,125,1000};
540    
541 +  DeterminePurity(basecut, nptcuts,ptcuts, AlphaVariable, ContainerName, FaceValueToBeTakenAt);
542 +  
543    float MPF_Results[nalphacuts][nptcuts];
544    float MPF_Errors[nalphacuts][nptcuts];
545    float RABS_Results[nalphacuts][nptcuts];
546    float RABS_Errors[nalphacuts][nptcuts];
547    
548    
549 <  for(int ia=0;ia<nalphacuts;ia++) {
549 >  for(int ia=0;ia<nalphacuts && !gofast;ia++) {
550      for(int ipt=0;ipt<nptcuts-1;ipt++) {
551        stringstream specialcut;
552 <      specialcut << "((pt>" << ptcuts[ipt] << " && pt< " << ptcuts[ipt+1] << ") && (Zb3010_alpha<" << alphacuts[ia] << "))";
552 >      float alow,ahigh;
553 >      if(ia>0) {
554 >        specialcut << "((pt>" << ptcuts[ipt] << " && pt< " << ptcuts[ipt+1] << ") && (" << AlphaVariable << "<" << alphacuts[ia] << " && " << AlphaVariable << ">" << alphacuts[ia-1] << "))";
555 >        alow=alphacuts[ia-1];
556 >        ahigh=alphacuts[ia];
557 >      } else {
558 >        specialcut << "((pt>" << ptcuts[ipt] << " && pt< " << ptcuts[ipt+1] << ") && (" << AlphaVariable << "<" << alphacuts[ia] << "))";
559 >        alow=0;
560 >        ahigh=alphacuts[ia];
561 >      }
562        cout << specialcut.str() << endl;
563 <      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]));
564 <      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]));
563 >      Value MPF_data_over_mc = get_Zb_data_over_mc(ContainerName,results,"mpf","",TCut(basecut && specialcut.str().c_str()),"MPF___pt_"+any2string(ptcuts[ipt])+"_"+any2string(ptcuts[ipt+1])+"__alpha_"+any2string(alphacuts[ia]), ptcuts[ipt], ptcuts[ipt+1],alow,ahigh);
564 >      Value RABS_data_over_mc = get_Zb_data_over_mc(ContainerName,results,"ZbCHS3010_pfJetGoodPt[0]/pt","",TCut(basecut && specialcut.str().c_str()),"Rabs___pt_"+any2string(ptcuts[ipt])+"_"+any2string(ptcuts[ipt+1])+"__alpha_"+any2string(alphacuts[ia]), ptcuts[ipt], ptcuts[ipt+1],alow,ahigh);
565        
566        MPF_Results[ia][ipt]=MPF_data_over_mc.getValue();
567        MPF_Errors[ia][ipt]=MPF_data_over_mc.getError();
# Line 358 | Line 577 | void new_data_mc_agreement_2d() {
577    float RABS_ExtraPolatedResults[nptcuts];
578    
579    TCanvas *can = new TCanvas("can");
580 +  can->cd(1)->SetLeftMargin(0.15);
581    
582    TGraphErrors *MPF_FinalGraph = new TGraphErrors();
583    TGraphErrors *RABS_FinalGraph = new TGraphErrors();
584 <  TGraphErrors *FaceValueAtPoint3 = new TGraphErrors();
584 >  TGraphErrors *MPF_FaceValueAtPoint3 = new TGraphErrors();
585 >  TGraphErrors *RABS_FaceValueAtPoint3 = new TGraphErrors();
586    
587    for(int ipt=0;ipt<nptcuts-1;ipt++) {
588      
# Line 370 | Line 591 | void new_data_mc_agreement_2d() {
591      cout << "Going to create histo for " << filename.str() << endl;
592      TGraphErrors *mpf_gr =new TGraphErrors();
593      TGraphErrors *rabs_gr =new TGraphErrors();
373 –    for(int ia=0;ia<nalphacuts;ia++) {
374 –      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;
375 –      mpf_gr->SetPoint(ia,alphacuts[ia],MPF_Results[ia][ipt]);     // avoiding empty bins (nan's)
376 –      mpf_gr->SetPointError(ia,0,MPF_Errors[ia][ipt]);             // avoiding empty bins (nan's)
377 –      rabs_gr->SetPoint(ia,alphacuts[ia],RABS_Results[ia][ipt]); // avoiding empty bins (nan's)
378 –      rabs_gr->SetPointError(ia,0,RABS_Errors[ia][ipt]);         // avoiding empty bins (nan's)
379 –    }
380 –    cout << "About to add a point in pt (" << ptcuts[ipt]+ptcuts[ipt+1] << " ) to face value : " << MPF_Results[PointThree][ipt] << endl;
381 –    FaceValueAtPoint3->SetPoint(ipt,(ptcuts[ipt]+ptcuts[ipt+1]),MPF_Results[PointThree][ipt]);
382 –    FaceValueAtPoint3->SetPointError(ipt,0,MPF_Errors[PointThree][ipt]);
383 –    
384 –    mpf_gr->SetMarkerStyle(21);
385 –    mpf_gr->SetMarkerColor(kBlue);
386 –    mpf_gr->Draw("AP");
387 –    mpf_gr->Fit("pol1","QE"); // quiet, use minos
388 –    mpf_gr->GetXaxis()->SetTitle("#alpha");
389 –    mpf_gr->GetXaxis()->CenterTitle();
390 –    mpf_gr->GetYaxis()->SetTitle("<MPF>_{data}/<MPF>_{mc}");
391 –    mpf_gr->GetYaxis()->CenterTitle();
392 –    TF1 *mpf_pol=(TF1*)mpf_gr->GetFunction("pol1");
393 –    float mpf_a=mpf_pol->GetParameter(0);
394 –    float mpf_b=mpf_pol->GetParameter(1);
395 –    MPF_FinalGraph->SetPoint(ipt,0.5*(ptcuts[ipt]+ptcuts[ipt+1]),mpf_a);
396 –    MPF_FinalGraph->SetPointError(ipt,0,mpf_pol->GetParError(0));
397 –    
398 –    MPF_ExtraPolatedResults[ipt]=mpf_a;
399 –    
400 –    stringstream mpf_info;
401 –    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() << "}";
402 –    
403 –    TText *mpf_mark = write_title(mpf_info.str());
404 –    mpf_mark->SetX(0.75);
405 –    mpf_mark->SetTextSize(0.03);
406 –    mpf_mark->SetY(0.75);
407 –    mpf_mark->Draw();
594      
595 <    string filenamebkp=filename.str();
596 <    filename << "__MPF";
411 <    DrawPrelim();
412 <    CompleteSave(can,filename.str());
413 <    cout << "MPF : " << mpf_a << " + " << mpf_b << " * alpha " << endl;
595 >    int rabs_counter=0;
596 >    int mpf_counter=0;
597      
598 <    filename.str("");
598 >    for(int ia=0;ia<nalphacuts && !gofast;ia++) {
599 >      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;
600 >      if(MPF_Results[ia][ipt]==MPF_Results[ia][ipt] && MPF_Errors[ia][ipt]==MPF_Errors[ia][ipt]) { // remove nan's
601 >       mpf_gr->SetPoint(mpf_counter,alphacuts[ia],MPF_Results[ia][ipt]);
602 >       mpf_gr->SetPointError(mpf_counter,0,MPF_Errors[ia][ipt]);
603 >       mpf_counter++;
604 >      } else {
605 >        dout << "Detected BS for MPF method (bin ignored) : " << endl;
606 >        dout << "       alpha: " << alphacuts[ia] << " , pt range: " << ptcuts[ipt] << " < pt < " << ptcuts[ipt+1] << endl;
607 >        cout << "       value (MPF) : " << MPF_Results[ia][ipt] << " +/- " << MPF_Errors[ia][ipt] << endl;
608 >      }
609 >      if(RABS_Results[ia][ipt]==RABS_Results[ia][ipt] && RABS_Errors[ia][ipt]==RABS_Errors[ia][ipt]) {
610 >       rabs_gr->SetPoint(rabs_counter,alphacuts[ia],RABS_Results[ia][ipt]);
611 >       rabs_gr->SetPointError(rabs_counter,0,RABS_Errors[ia][ipt]);
612 >       rabs_counter++;
613 >      } else {
614 >        dout << "Detected BS for RABS method (bin ignored) : " << endl;
615 >        dout << "       alpha: " << alphacuts[ia] << " , pt range: " << ptcuts[ipt] << " < pt < " << ptcuts[ipt+1] << endl;
616 >        cout << "       value (RABS) : " << RABS_Results[ia][ipt] << " +/- " << RABS_Errors[ia][ipt] << endl;
617 >      }
618 >    }
619      
620 <    rabs_gr->SetMarkerStyle(21);
621 <    rabs_gr->SetMarkerColor(kBlue);
622 <    rabs_gr->Fit("pol1","QE"); // quiet, use minos
623 <    rabs_gr->Draw("AP");
624 <    rabs_gr->GetXaxis()->SetTitle("#alpha");
625 <    rabs_gr->GetXaxis()->CenterTitle();
626 <    rabs_gr->GetYaxis()->SetTitle("<R_{abs}>_{data}/<R_{abs}>_{mc}");
627 <    rabs_gr->GetYaxis()->CenterTitle();
620 >    stringstream specialcut;
621 >    specialcut << "((pt>" << ptcuts[ipt] << " && pt< " << ptcuts[ipt+1] << ") && (" << AlphaVariable << "<" << FaceValueToBeTakenAt << "))";
622 >    Value MPF_data_over_mc = get_Zb_data_over_mc(ContainerName,results,"mpf","",TCut(basecut && specialcut.str().c_str()),"MPF___pt_"+any2string(ptcuts[ipt])+"_"+any2string(ptcuts[ipt+1])+"__INCLUSIVEalpha_"+any2string(FaceValueToBeTakenAt),ptcuts[ipt],ptcuts[ipt+1],0,FaceValueToBeTakenAt);
623 >    Value RABS_data_over_mc = get_Zb_data_over_mc(ContainerName,results,"ZbCHS3010_pfJetGoodPt[0]/pt","",TCut(basecut && specialcut.str().c_str()),"Rabs___pt_"+any2string(ptcuts[ipt])+"_"+any2string(ptcuts[ipt+1])+"__INCLUSIVEalpha_"+any2string(FaceValueToBeTakenAt),ptcuts[ipt],ptcuts[ipt+1],0,FaceValueToBeTakenAt);
624 >    cout << "About to add a point in pt (" << (ptcuts[ipt]+ptcuts[ipt+1])/2 << " ) to face value : " << MPF_data_over_mc << " (MPF) , " <<  RABS_data_over_mc << "(RABS)" << endl;
625 >    MPF_FaceValueAtPoint3->SetPoint(ipt,(ptcuts[ipt]+ptcuts[ipt+1])/2,MPF_data_over_mc.getValue());
626 >    MPF_FaceValueAtPoint3->SetPointError(ipt,0,MPF_data_over_mc.getError());
627 >    RABS_FaceValueAtPoint3->SetPoint(ipt,(ptcuts[ipt]+ptcuts[ipt+1])/2,RABS_data_over_mc.getValue());
628 >    RABS_FaceValueAtPoint3->SetPointError(ipt,0,RABS_data_over_mc.getError());
629      
426 –    TF1 *rabs_pol=(TF1*)rabs_gr->GetFunction("pol1");
427 –    float rabs_a=rabs_pol->GetParameter(0);
428 –    float rabs_b=rabs_pol->GetParameter(1);
429 –    RABS_FinalGraph->SetPoint(ipt,0.5*(ptcuts[ipt]+ptcuts[ipt+1]),rabs_a);
430 –    RABS_FinalGraph->SetPointError(ipt,0,rabs_pol->GetParError(0));
630      
432 –    cout << "!+!+!+!+!+!+!+!+!+ Just added a point to the final plot : " << rabs_a << " +/- " << rabs_pol->GetParError(0) << endl;
631      
632 +    if(!gofast) {
633 +      can->cd();
634 +      mpf_gr->SetMarkerStyle(21);
635 +      mpf_gr->SetMarkerColor(kBlue);
636 +      mpf_gr->Draw("AP");
637 +      mpf_gr->Fit("pol1","");
638 +      mpf_gr->GetXaxis()->SetTitle("#alpha");
639 +      mpf_gr->GetXaxis()->CenterTitle();
640 +      mpf_gr->GetYaxis()->SetTitle("<MPF>_{data}/<MPF>_{mc}");
641 +      mpf_gr->GetYaxis()->CenterTitle();
642 +      TF1 *mpf_pol=(TF1*)mpf_gr->GetFunction("pol1");
643 +      mpf_pol->SetLineWidth(0);
644 +      TF1 *mpf_pol_complete = new TF1("mpf_pol_complete","[0]+[1]*x");
645 +      mpf_pol_complete->SetParameters(mpf_pol->GetParameters());
646 +      mpf_pol_complete->SetLineWidth(1);
647 +      mpf_pol_complete->SetLineColor(kBlue);
648 +      
649 +      float min,max;
650 +      FindMinMax(mpf_gr,min,max);
651 +      TH1F *histo = new TH1F("histo","histo",1,0,0.4);
652 +      histo->GetXaxis()->SetTitle("#alpha");
653 +      histo->GetYaxis()->SetTitle("<MPF>_{data}/<MPF>_{mc}");
654 +      histo->GetXaxis()->CenterTitle();
655 +      histo->GetYaxis()->CenterTitle();
656 +      histo->GetYaxis()->SetTitleOffset(1.3);
657 +      histo->SetStats(0);
658 +      histo->GetXaxis()->SetRangeUser(0,0.4);
659 +      histo->GetYaxis()->SetRangeUser(min,max);
660 +      mpf_gr->GetYaxis()->SetRangeUser(min,max);
661 +      
662 +      TPolyLine *mpf_fit_uncert = GetFitUncertaintyShape(mpf_gr, "pol1", min, max,0.0,0.4);
663 +      mpf_pol->SetLineColor(TColor::GetColor("#0101DF"));
664 +      mpf_fit_uncert->SetFillColor(TColor::GetColor("#5353FC"));
665 +      histo->GetYaxis()->SetTitle("<MPF>_{data}/<MPF>_{mc}");
666 +      histo->Draw();
667 +      mpf_fit_uncert->Draw("F");
668 +      mpf_fit_uncert->Draw("");
669 +      mpf_gr->Draw("P");
670 +      mpf_pol_complete->Draw("same");
671 +      
672 +      float mpf_a=mpf_pol->GetParameter(0);
673 +      float mpf_b=mpf_pol->GetParameter(1);
674 +      MPF_FinalGraph->SetPoint(ipt,0.5*(ptcuts[ipt]+ptcuts[ipt+1]),mpf_a);
675 +      MPF_FinalGraph->SetPointError(ipt,0,mpf_pol->GetParError(0));
676 +      
677 +      MPF_ExtraPolatedResults[ipt]=mpf_a;
678 +      
679 +      stringstream mpf_info;
680 +      mpf_info <<  "#splitline{#splitline{#splitline{" << ptcuts[ipt] << " GeV < p_{T}^{Z} < " << ptcuts[ipt+1] << " GeV }{ (MPF) }}{Extrapolated value: " << std::setprecision(4) << mpf_a << "}}{#chi^{2}/NDF : " << mpf_pol->GetChisquare() << " / " << mpf_pol->GetNDF() << "}";
681 +      
682 +      TText *mpf_mark = write_title(mpf_info.str());
683 +      mpf_mark->SetX(0.75);
684 +      mpf_mark->SetTextSize(0.03);
685 +      mpf_mark->SetY(0.75);
686 +      mpf_mark->Draw();
687 +      
688 +      string filenamebkp=filename.str();
689 +      filename << "__MPF";
690 +      DrawPrelim();
691 +      CompleteSave(can,ContainerName+"/"+filename.str());
692 +      cout << "MPF : " << mpf_a << " + " << mpf_b << " * alpha " << endl;
693 +      
694 +      filename.str("");
695 +      
696 +      rabs_gr->SetMarkerStyle(21);
697 +      rabs_gr->SetMarkerColor(kRed);
698 +      rabs_gr->Draw("AP");
699 +      rabs_gr->Fit("pol1","");
700 +      rabs_gr->GetXaxis()->SetTitle("#alpha");
701 +      rabs_gr->GetXaxis()->CenterTitle();
702 +      rabs_gr->GetYaxis()->SetTitle("<R_{abs}>_{data}/<R_{abs}>_{mc}");
703 +      rabs_gr->GetYaxis()->CenterTitle();
704 +      TF1 *rabs_pol=(TF1*)rabs_gr->GetFunction("pol1");
705 +      rabs_pol->SetLineWidth(0);
706 +      TF1 *rabs_pol_complete = new TF1("rabs_pol_complete","[0]+[1]*x");
707 +      rabs_pol_complete->SetParameters(rabs_pol->GetParameters());
708 +      rabs_pol_complete->SetLineColor(kRed);
709 +      rabs_pol_complete->SetLineWidth(1);
710 +      FindMinMax(rabs_gr,min,max);
711 +      rabs_gr->GetYaxis()->SetRangeUser(min,max);
712 +      histo->GetYaxis()->SetRangeUser(min,max);
713 +      TPolyLine *rabs_fit_uncert = GetFitUncertaintyShape(rabs_gr, "pol1", min, max,0.0,0.4);
714 +      rabs_pol->SetLineColor(TColor::GetColor("#FA5858"));
715 +      rabs_pol->SetFillColor(TColor::GetColor("#FA5858"));
716 +      rabs_fit_uncert->SetLineColor(TColor::GetColor("#FA5858"));
717 +      rabs_fit_uncert->SetFillColor(TColor::GetColor("#FA5858"));
718 +      histo->GetYaxis()->SetTitle("<R_{abs}>_{data}/<R_{abs}>_{mc}");
719 +      histo->Draw();
720 +      rabs_fit_uncert->Draw("F");
721 +      rabs_fit_uncert->Draw("");
722 +      rabs_gr->Draw("P");
723 +      rabs_pol_complete->Draw("same");
724 +      
725 +      float rabs_a=rabs_pol->GetParameter(0);
726 +      float rabs_b=rabs_pol->GetParameter(1);
727 +      RABS_FinalGraph->SetPoint(ipt,0.5*(ptcuts[ipt]+ptcuts[ipt+1]),rabs_a);
728 +      RABS_FinalGraph->SetPointError(ipt,0,rabs_pol->GetParError(0));
729 +      
730 +      cout << "!+!+!+!+!+!+!+!+!+ Just added a point to the final plot : " << rabs_a << " +/- " << rabs_pol->GetParError(0) << endl;
731 +      
732 +      stringstream rabs_info;
733 +      rabs_info << "#splitline{#splitline{#splitline{" << ptcuts[ipt] << " GeV < p_{T}^{Z} < " << ptcuts[ipt+1] << " GeV }{ (R_{abs}) }}{Extrapolated value: " << std::setprecision(4) << rabs_a << "}}{#chi^{2}/NDF : " << rabs_pol->GetChisquare() << " / " << rabs_pol->GetNDF() << "}";
734 +      TText *rabs_mark = write_title(rabs_info.str());
735 +      rabs_mark->SetX(0.75);
736 +      rabs_mark->SetTextSize(0.03);
737 +      rabs_mark->SetY(0.75);
738 +      rabs_mark->Draw();
739 +      
740 +      RABS_ExtraPolatedResults[ipt]=rabs_a;
741 +      filename << filenamebkp << "__RABS";
742 +      DrawPrelim();
743 +      CompleteSave(can,ContainerName+"/"+filename.str());
744 +      cout << "RABS : " << rabs_a << " + " << rabs_b << " * alpha " << endl;
745 +    }
746 +  }
747 +  
748 +  float MPFResult;
749 +  float MPFResultError;
750 +  
751 +  float RabsResult;
752 +  float RabsResultError;
753 +
754 +  if(!gofast) {
755 +    MPF_FinalGraph->SetMarkerStyle(21);
756 +    MPF_FinalGraph->SetMarkerColor(kBlue);
757 +    MPF_FinalGraph->Fit("pol0",""); // quiet, use minos
758 +    
759 +    TF1 *mpf_pol0=(TF1*)MPF_FinalGraph->GetFunction("pol0");
760 +    TF1 *mpf_pol0_complete = new TF1("mpf_pol0_complete","[0]+[1]*x");
761 +    mpf_pol0_complete->SetLineWidth(1);
762 +    mpf_pol0_complete->SetLineColor(kBlue);
763 +    mpf_pol0->SetLineColor(kBlue);
764 +    MPF_FinalGraph->Draw("AP");
765 +    MPF_FinalGraph->GetXaxis()->SetTitle("p_{T}^{Z}");
766 +    MPF_FinalGraph->GetXaxis()->CenterTitle();
767 +    MPF_FinalGraph->GetYaxis()->SetTitle("<MPF>_{data}/<MPF>_{mc}");
768 +    MPF_FinalGraph->GetYaxis()->CenterTitle();
769 +    MPF_FinalGraph->Draw("AP");
770 +    mpf_pol0_complete->Draw("same");
771 +    
772 +    stringstream mpf_result;
773 +    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)}";
774 +    TText *rmpf_final_mark = write_title(mpf_result.str());
775 +    rmpf_final_mark->SetX(0.75);
776 +    rmpf_final_mark->SetY(0.75);
777 +    rmpf_final_mark->SetTextSize(0.03);
778 +    rmpf_final_mark->Draw();
779 +    DrawPrelim();
780 +    MPFResult=1/mpf_pol0->GetParameter(0);
781 +    MPFResultError=mpf_pol0->GetParError(0)/(mpf_pol0->GetParameter(0)*mpf_pol0->GetParameter(0));
782      
783 <    stringstream rabs_info;
784 <    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() << "}";
785 <    TText *rabs_mark = write_title(rabs_info.str());
786 <    rabs_mark->SetX(0.75);
787 <    rabs_mark->SetTextSize(0.03);
788 <    rabs_mark->SetY(0.75);
789 <    rabs_mark->Draw();
783 >    stringstream filename2;
784 >    filename2 << "Extrapolation/FINAL_RESULT_MPF";
785 >    CompleteSave(can,ContainerName+"/"+filename2.str());
786 >    
787 >    RABS_FinalGraph->SetMarkerStyle(21);
788 >    RABS_FinalGraph->SetMarkerStyle(21);
789 >    RABS_FinalGraph->SetMarkerColor(kRed);
790 >    RABS_FinalGraph->Fit("pol0","QE"); // quiet, use minos
791 >    RABS_FinalGraph->Draw("AP");
792 >    RABS_FinalGraph->GetXaxis()->SetTitle("p_{T}^{Z}");
793 >    RABS_FinalGraph->GetXaxis()->CenterTitle();
794 >    RABS_FinalGraph->GetYaxis()->SetTitle("<R_{abs}>_{data}/<R_{abs}>_{mc}");
795 >    RABS_FinalGraph->GetYaxis()->CenterTitle();
796 >    RABS_FinalGraph->Draw("AP");
797 >    TF1 *rabs_pol0=(TF1*)RABS_FinalGraph->GetFunction("pol0");
798 >    stringstream rabs_result;
799 >    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})}";;
800 >    TText *rabs_final_mark = write_title(rabs_result.str());
801 >    rabs_final_mark->SetX(0.75);
802 >    rabs_final_mark->SetY(0.75);
803 >    rabs_final_mark->SetTextSize(0.03);
804 >    rabs_final_mark->Draw();
805 >    DrawPrelim();
806      
807 +    RabsResult=1/rabs_pol0->GetParameter(0);
808 +    RabsResultError=rabs_pol0->GetParError(0)/(rabs_pol0->GetParameter(0)*rabs_pol0->GetParameter(0));
809      
810 <
811 <    RABS_ExtraPolatedResults[ipt]=rabs_a;
812 <    filename << filenamebkp << "__RABS";
447 <    DrawPrelim();
448 <    CompleteSave(can,filename.str());
449 <    cout << "RABS : " << rabs_a << " + " << rabs_b << " * alpha " << endl;
810 >    filename2.str("");
811 >    filename2 << "Extrapolation/FINAL_RESULT_RABS";
812 >    CompleteSave(can,ContainerName+"/"+filename2.str());
813    }
814 <  
452 <  MPF_FinalGraph->SetMarkerStyle(21);
453 <  MPF_FinalGraph->SetMarkerColor(kRed);
454 <  MPF_FinalGraph->Fit("pol0","QE"); // quiet, use minos
455 <  
456 <  TF1 *mpf_pol0=(TF1*)MPF_FinalGraph->GetFunction("pol0");
457 <  MPF_FinalGraph->Draw("AP");
458 <  MPF_FinalGraph->GetXaxis()->SetTitle("p_{T}^{Z}");
459 <  MPF_FinalGraph->GetXaxis()->CenterTitle();
460 <  MPF_FinalGraph->GetYaxis()->SetTitle("<MPF>_{data}/<MPF>_{mc}");
461 <  MPF_FinalGraph->GetYaxis()->CenterTitle();
462 <  MPF_FinalGraph->Draw("AP");
463 <  
464 <  stringstream mpf_result;
465 <  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)}";
466 <  TText *rmpf_final_mark = write_title(mpf_result.str());
467 <  rmpf_final_mark->SetX(0.75);
468 <  rmpf_final_mark->SetY(0.75);
469 <  rmpf_final_mark->SetTextSize(0.03);
470 <  rmpf_final_mark->Draw();
471 <  DrawPrelim();
472 <  float MPFResult=1/mpf_pol0->GetParameter(0);
473 <  float MPFResultError=mpf_pol0->GetParError(0)/(mpf_pol0->GetParameter(0)*mpf_pol0->GetParameter(0));
474 <  
475 <  stringstream filename2;
476 <  filename2 << "Extrapolation/FINAL_RESULT_MPF";
477 <  CompleteSave(can,filename2.str());
478 <  
479 <  
480 <  RABS_FinalGraph->SetMarkerStyle(21);
481 <  RABS_FinalGraph->SetMarkerStyle(21);
482 <  RABS_FinalGraph->SetMarkerColor(kRed);
483 <  RABS_FinalGraph->Fit("pol0","QE"); // quiet, use minos
484 <  RABS_FinalGraph->Draw("AP");
485 <  RABS_FinalGraph->GetXaxis()->SetTitle("p_{T}^{Z}");
486 <  RABS_FinalGraph->GetXaxis()->CenterTitle();
487 <  RABS_FinalGraph->GetYaxis()->SetTitle("<R_{abs}>_{data}/<R_{abs}>_{mc}");
488 <  RABS_FinalGraph->GetYaxis()->CenterTitle();
489 <  RABS_FinalGraph->Draw("AP");
490 <  TF1 *rabs_pol0=(TF1*)RABS_FinalGraph->GetFunction("pol0");
491 <  stringstream rabs_result;
492 <  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})}";;
493 <  TText *rabs_final_mark = write_title(rabs_result.str());
494 <  rabs_final_mark->SetX(0.75);
495 <  rabs_final_mark->SetY(0.75);
496 <  rabs_final_mark->SetTextSize(0.03);
497 <  rabs_final_mark->Draw();
498 <  DrawPrelim();
499 <  
500 <  float RabsResult=1/rabs_pol0->GetParameter(0);
501 <  float RabsResultError=rabs_pol0->GetParError(0)/(rabs_pol0->GetParameter(0)*rabs_pol0->GetParameter(0));
502 <  
503 <  filename2.str("");
504 <  filename2 << "Extrapolation/FINAL_RESULT_RABS";
505 <  CompleteSave(can,filename2.str());
506 <  
814 >
815    can->cd();
816 <  TCanvas *FaceValueCanvas = new TCanvas("FaceValueCanvas","FaceValueCanvas");
817 <  FaceValueAtPoint3->SetMarkerStyle(21);
818 <  FaceValueAtPoint3->SetMarkerStyle(21);
819 <  FaceValueAtPoint3->Fit("pol0","QE"); // quiet, use minos
820 <  FaceValueAtPoint3->Draw("AP");
821 <  FaceValueAtPoint3->GetXaxis()->SetTitle("p_{T}^{Z}");
822 <  FaceValueAtPoint3->GetXaxis()->CenterTitle();
823 <  FaceValueAtPoint3->GetYaxis()->SetTitle("data/mc for #alpha<0.3");
824 <  FaceValueAtPoint3->GetYaxis()->CenterTitle();
825 <  FaceValueAtPoint3->Draw("AP");
826 <  TF1 *faceval_pol0=(TF1*)FaceValueAtPoint3->GetFunction("pol0");
816 >  MPF_FaceValueAtPoint3->SetMarkerStyle(21);
817 >  MPF_FaceValueAtPoint3->SetMarkerColor(kBlue);
818 >  MPF_FaceValueAtPoint3->Fit("pol0","QE"); // quiet, use minos
819 >  MPF_FaceValueAtPoint3->Draw("AP");
820 >  MPF_FaceValueAtPoint3->GetXaxis()->SetTitle("p_{T}^{Z}");
821 >  MPF_FaceValueAtPoint3->GetXaxis()->CenterTitle();
822 >  MPF_FaceValueAtPoint3->GetYaxis()->SetTitle(((string)"<MPF>__{data}/<MPF>_{mc} for #alpha<"+any2string(FaceValueToBeTakenAt)).c_str());
823 >  MPF_FaceValueAtPoint3->GetYaxis()->CenterTitle();
824 >  MPF_FaceValueAtPoint3->Draw("AP");
825 >  TF1 *faceval_pol0=(TF1*)MPF_FaceValueAtPoint3->GetFunction("pol0");
826 >  faceval_pol0->SetLineColor(kBlue);
827 >  faceval_pol0->SetLineWidth(1);
828 >  faceval_pol0->Draw("same");
829 >  TF1 *faceval_pol0_complete = new TF1("faceval_pol0_complete","[0]+[1]*x");
830 >  faceval_pol0_complete->SetParameters(faceval_pol0->GetParameters());
831    stringstream faceval_result;
832 <  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)}";;
832 >  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<" << FaceValueToBeTakenAt << ")}";;
833    TText *faceval_final_mark = write_title(faceval_result.str());
834    faceval_final_mark->SetX(0.75);
835    faceval_final_mark->SetY(0.75);
# Line 528 | Line 840 | void new_data_mc_agreement_2d() {
840    float FaceValResult=1/faceval_pol0->GetParameter(0);
841    float FaceValResultError=faceval_pol0->GetParError(0)/(faceval_pol0->GetParameter(0)*faceval_pol0->GetParameter(0));
842    
843 +  stringstream filename2;
844    filename2.str("");
845    filename2 << "Extrapolation/FINAL_RESULT_MPF_FaceValp3";
846 <  CompleteSave(can,filename2.str());
846 >  CompleteSave(can,ContainerName+"/"+filename2.str());
847 >  
848 >  can->cd();
849 >  RABS_FaceValueAtPoint3->SetMarkerStyle(21);
850 >  RABS_FaceValueAtPoint3->SetMarkerColor(kRed);
851 >  RABS_FaceValueAtPoint3->Fit("pol0","QE"); // quiet, use minos
852 >  RABS_FaceValueAtPoint3->Draw("AP");
853 >  RABS_FaceValueAtPoint3->GetXaxis()->SetTitle("p_{T}^{Z}");
854 >  RABS_FaceValueAtPoint3->GetXaxis()->CenterTitle();
855 >  RABS_FaceValueAtPoint3->GetYaxis()->SetTitle(((string)"<R_{abs}>_{data}/<R_{abs}>_{mc} for #alpha<"+any2string(FaceValueToBeTakenAt)).c_str());
856 >  RABS_FaceValueAtPoint3->GetYaxis()->CenterTitle();
857 >  RABS_FaceValueAtPoint3->Draw("AP");
858 >  TF1 *faceval_pol0RABS=(TF1*)RABS_FaceValueAtPoint3->GetFunction("pol0");
859 >  TF1 *faceval_pol0RABS_complete = new TF1("faceval_pol0_complete","[0]+[1]*x");
860 >  faceval_pol0RABS_complete->SetParameters(faceval_pol0RABS->GetParameters());
861 >  faceval_pol0RABS->SetLineWidth(1);
862 >  faceval_pol0RABS->SetLineColor(kRed);
863 >  faceval_pol0RABS_complete->Draw("same");
864 >  stringstream RABSfaceval_result;
865 >  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<" << FaceValueToBeTakenAt << ")}";;
866 >  TText *RABSfaceval_final_mark = write_title(RABSfaceval_result.str());
867 >  RABSfaceval_final_mark->SetX(0.75);
868 >  RABSfaceval_final_mark->SetY(0.75);
869 >  RABSfaceval_final_mark->SetTextSize(0.03);
870 >  RABSfaceval_final_mark->Draw();
871 >  DrawPrelim();
872 >  
873 >  float RABSFaceValResult=1/faceval_pol0RABS->GetParameter(0);
874 >  float RABSFaceValResultError=faceval_pol0RABS->GetParError(0)/(faceval_pol0RABS->GetParameter(0)*faceval_pol0RABS->GetParameter(0));
875 >  
876 >  filename2.str("");
877 >  filename2 << "Extrapolation/FINAL_RESULT_RABS_FaceValp3";
878 >  CompleteSave(can,ContainerName+"/"+filename2.str());
879    
880    cout << "FINAL RESULTS: " << endl;
881 <  cout << "MPF : " << MPFResult  << " +/- " << MPFResultError  << endl;
882 <  cout << "Rabs: " << RabsResult << " +/- " << RabsResultError << endl;
883 <  cout << "Face value at 0.3 :" << FaceValResult << " +/- " << FaceValResultError << endl;
881 >  if(!gofast) cout << "MPF : " << MPFResult  << " +/- " << MPFResultError  << endl;
882 >  if(!gofast) cout << "Rabs: " << RabsResult << " +/- " << RabsResultError << endl;
883 >  cout << "Face value at " << FaceValueToBeTakenAt << " :" << FaceValResult << " +/- " << FaceValResultError << "  (MPF) " << endl;
884 >  cout << "Face value at " << FaceValueToBeTakenAt << " :" << RABSFaceValResult << " +/- " << RABSFaceValResultError << "  (RABS) " << endl;
885    
886    delete can;
887    
888 <  
888 >  results.MPF_Result_FaceValue=Value(FaceValResult,FaceValResultError);
889 >  results.Rabs_Result_FaceValue=Value(RABSFaceValResult,RABSFaceValResultError);
890 >  if(!gofast) results.MPF_Result_Extrapolated=Value(MPFResult,MPFResultError);
891 >  if(!gofast) results.Rabs_Result_Extrapolated=Value(RabsResult,RabsResultError);
892 >  cout << results << endl;
893 >  return results;
894   }
895  
896  
# Line 562 | Line 913 | void DoPUStudy(string identifier) {
913        cout << numVtxcuts[i-1] << ";";
914      }
915      
916 +    sSpecialCut << " && " << identifier << "_alpha<0.3";
917      TCut SpecialCut(sSpecialCut.str().c_str());
918      
919      
920 <    malpha[i] = allsamples.Draw("malpha"+any2string(i),"mpf",1000,0,10,"MPF","events",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&SpecialCut,mc,luminosity);
921 <    dalpha[i] = allsamples.Draw("dalpha"+any2string(i),"mpf",1000,0,10,"MPF","events",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&SpecialCut,data,luminosity);
920 >    malpha[i] = allsamples.Draw("malpha"+any2string(i),"mpf",1000,0,2,"MPF","events",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&SpecialCut,mc,luminosity);
921 >    dalpha[i] = allsamples.Draw("dalpha"+any2string(i),"mpf",1000,0,2,"MPF","events",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&SpecialCut,data,luminosity);
922      
923      
924      float a=dalpha[i]->GetMean();
# Line 577 | Line 929 | void DoPUStudy(string identifier) {
929      float error = TMath::Sqrt(  (1/(b*b)) * (da*da) + ((a*a)/(b*b))*db*db);
930    
931      cout << malpha[i]->GetMean() << ";" << malpha[i]->GetMeanError() << ";" << dalpha[i]->GetMean() << ";" << dalpha[i]->GetMeanError() << ";" << malpha[i]->Integral()<<";"<<dalpha[i]->Integral() << ";"<<factor << ";" << error <<  endl;
580 –    malpha[i]->SetLineColor(i+1);
581 –    if(i==0) malpha[i]->Draw("histo");
582 –    else malpha[i]->Draw("histo,same");
583 –  }
584 –  
585 –  CompleteSave(can,"AlphaOverview_"+identifier);
586 –  delete can;
587 –  
588 – /*  for(int i=1;i<8;i++) {
932      delete malpha[i];
933      delete dalpha[i];
934    }
935 <  */
935 >  
936 >  delete can;
937    
938   }
939  
940   void ScenarioComparison() {
941 <  
941 > //very dumb way to do this.
942 >
943    TGraphAsymmErrors *gr[5];
944    TF1 *fit[5];
945 <  bool dofit=false;
945 >  bool dofit=true;
946    
947    TCanvas *can = new TCanvas("can","can",500,500);
948    
# Line 605 | Line 950 | void ScenarioComparison() {
950    leg->SetBorderSize(0);
951    leg->SetFillColor(kWhite);
952    
953 <  for(int i=0;i<5;i++) {
953 >  for(int i=0;i<6;i++) {
954      gr[i] = new TGraphAsymmErrors();
955 <    double x[4]   = {5,12.5,17.5,20};
956 <    double dxl[4] = {5,2.5,2.5,2.5};
955 >    double x[4]   = {5,12.5,17.5,25};
956 >    double dxl[4] = {5,2.5,2.5,5};
957      double dxh[4] = {5,2.5,2.5,10};
958      
959      string name="";
960      string color="";
961      
962      if(i==0) {
963 <      double y[4]   = {1.008,1.013,1.003,1.007};
964 <      double dyl[4] = {0.011,0.009,0.013,0.020};
965 <      double dyh[4] = {0.011,0.009,0.013,0.020};
963 >      double y[4]   = {1.01935,1.01558,0.996438,1.01901};
964 >      double dyl[4] = {0.0132973,0.0131661,0.017629,0.0319767};
965 >      double dyh[4] = {0.0132973,0.0131661,0.017629,0.0319767};
966        name="30/30";
967        gr[i] = new TGraphAsymmErrors(4,x,y,dxl,dxh,dyl,dyh);
968        color="#a01d99";
969      }
970      if(i==1) {
971 <      double y[4]   = {1.008+0.001,1.013+0.001,1.003+0.001,1.007+0.001};
972 <      double dyl[4] = {0.011,0.009,0.013,0.020};
973 <      double dyh[4] = {0.011,0.009,0.013,0.020};
971 >      double y[4]   = {1.00577,1.00706,0.991029,0.954098};
972 >      double dyl[4] = {0.0135575,0.0136548,0.0170778,0.0263919};
973 >      double dyh[4] = {0.0135575,0.0136548,0.0170778,0.0263919};
974        gr[i] = new TGraphAsymmErrors(4,x,y,dxl,dxh,dyl,dyh);
975 <      name="30/10";
975 >      name="30/15";
976        color="#2464bc";
977      }
978      if(i==2) {
979 <      double y[4]   = {1.014,1.026,1.014,1.059};
980 <      double dyl[4] = {0.011,0.010,0.014,0.023};
981 <      double dyh[4] = {0.011,0.010,0.014,0.023};
979 >      double y[4]   = {1.023,1.00648,1.02217,1.03181};
980 >      double dyl[4] = {0.0157372,0.0187778,0.0261512,0.0422295};
981 >      double dyh[4] = {0.0157372,0.0187778,0.0261512,0.0422295};
982        gr[i] = new TGraphAsymmErrors(4,x,y,dxl,dxh,dyl,dyh);
983        name="15/15";
984        color="#f49230";
985      }
986      if(i==3) {
987 <      double y[4]   = {1.015,1.020,1.008,1.050};
988 <      double dyl[4] = {0.011,0.010,0.014,0.022};
989 <      double dyh[4] = {0.011,0.010,0.014,0.022};
987 >      double y[4]   = {1.02029,1.01338,0.999769,1.00357};
988 >      double dyl[4] = {0.0129509,0.0120127,0.0160537,0.0268972};
989 >      double dyh[4] = {0.0129509,0.0120127,0.0160537,0.0268972};
990        gr[i] = new TGraphAsymmErrors(4,x,y,dxl,dxh,dyl,dyh);
991 <      name="20/15";
991 >      name="CHS 30/30";
992        color="#c72f3c";
993      }
994      if(i==4) {
995 <      double y[4]   = {1.015+0.001,1.020+0.001,1.008+0.001,1.050+0.001};
996 <      double dyl[4] = {0.011,0.010,0.014,0.022};
997 <      double dyh[4] = {0.011,0.010,0.014,0.022};
995 >      double y[4]   = {1.00274,1.01056,0.993214,0.997963};
996 >      double dyl[4] = {0.0127179,0.0140344,0.0162814,0.0278914};
997 >      double dyh[4] = {0.0127179,0.0140344,0.0162814,0.0278914};
998        gr[i] = new TGraphAsymmErrors(4,x,y,dxl,dxh,dyl,dyh);
999 <      name="20/20";
999 >      name="CHS 30/15";
1000        color="#7cb433";
1001      }
1002 <    
1003 <    
1002 >    if(i==5) {
1003 >      double y[4]   = {1.04103,1.03726,1.04649,1.05383};
1004 >      double dyl[4] = {0.0156465,0.0168194,0.0231982,0.0376114};
1005 >      double dyh[4] = {0.0156465,0.0168194,0.0231982,0.0376114};
1006 >      gr[i] = new TGraphAsymmErrors(4,x,y,dxl,dxh,dyl,dyh);
1007 >      name="CHS 15/15";
1008 >      color="#DF013A";
1009 >    }
1010 >    
1011 >    
1012      gr[i]->SetTitle();
1013      gr[i]->GetXaxis()->SetTitle("N(Vtx)");
1014      gr[i]->GetYaxis()->SetTitle("<MPF>_{data}/<MPF>_{mc}");
# Line 684 | Line 1037 | void ScenarioComparison() {
1037   }
1038  
1039   void ScenarioComparisonInclusive() {
1040 <  
1040 >  //dumbest way ever to do this. but ok we only need to do it once.
1041    TGraphAsymmErrors *gr[5];
1042    TF1 *fit[5];
1043    bool dofit=true;
1044    
1045    TCanvas *can = new TCanvas("can","can",500,500);
1046    
1047 <  TLegend *leg = new TLegend(0.2,0.2,0.5,0.45);
1047 >  TLegend *leg = new TLegend(0.2,0.2,0.65,0.45);
1048    leg->SetBorderSize(0);
1049    leg->SetFillColor(kWhite);
1050    
1051 <  for(int i=0;i<5;i++) {
1051 >  for(int i=0;i<6;i++) {
1052      gr[i] = new TGraphAsymmErrors();
1053 <    double x[4]   = {5,12.5,17.5,20};
1054 <    double dxl[4] = {5,2.5,2.5,2.5};
1055 <    double dxh[4] = {5,2.5,2.5,10};
1053 >    double x[4]   = {5,12.5,17.5,25};
1054 >    double dxl[4] = {5,2.5,2.5,5};
1055 >    double dxh[4] = {5,2.5,2.5,5};
1056      
1057      string name="";
1058      string color="";
1059      
1060      if(i==0) {
1061 <      double y[4]   = {1.01103,1.01496,1.01486,1.02785};
1062 <      double dyl[4] = {0.00235044,0.00225264,0.00313574,0.00541047};
1063 <      double dyh[4] = {0.00235044,0.00225264,0.00313574,0.00541047};
1061 >      double y[4]   = {1.01093,1.01379,1.0147,1.02766};
1062 >      double dyl[4] = {0.00243596,0.00242534,0.00341941,0.00600663};
1063 >      double dyh[4] = {0.00243596,0.00242534,0.00341941,0.00600663};
1064        name="30/30";
1065        gr[i] = new TGraphAsymmErrors(4,x,y,dxl,dxh,dyl,dyh);
1066        color="#a01d99";
1067      }
1068      if(i==1) {
1069 <      double y[4]   = {1.01103+0.001,1.01496+0.001,1.01486+0.001,1.02785+0.001};
1070 <      double dyl[4] = {0.00235044,0.00225264,0.00313574,0.00541047};
1071 <      double dyh[4] = {0.00235044,0.00225264,0.00313574,0.00541047};
1069 >      double y[4]   = {1.0068,1.00779,0.999692,0.994484};
1070 >      double dyl[4] = {0.00249233,0.00257169,0.00371399,0.00631138};
1071 >      double dyh[4] = {0.00249233,0.00257169,0.00371399,0.00631138};
1072        gr[i] = new TGraphAsymmErrors(4,x,y,dxl,dxh,dyl,dyh);
1073 <      name="30/10";
1073 >      name="30/15";
1074        color="#2464bc";
1075      }
1076      if(i==2) {
1077 <      double y[4]   = {1.01575,1.02879,1.0325,1.05456};
1078 <      double dyl[4] = {0.00288281,0.00278536,0.00367451,0.00590421};
1079 <      double dyh[4] = {0.00288281,0.00278536,0.00367451,0.00590421};
1077 >      double y[4]   = {1.01208,1.01991,1.02185,1.03479};
1078 >      double dyl[4] = {0.00366745,0.00400001,0.00584054,0.0104661};
1079 >      double dyh[4] = {0.00366745,0.00400001,0.00584054,0.0104661};
1080        gr[i] = new TGraphAsymmErrors(4,x,y,dxl,dxh,dyl,dyh);
1081        name="15/15";
1082        color="#f49230";
1083      }
1084      if(i==3) {
1085 <      double y[4]   = {1.01221,1.01922,1.02396,1.04677};
1086 <      double dyl[4] = {0.00263821,0.00260853,0.00353673,0.00584728};
1087 <      double dyh[4] = {0.00263821,0.00260853,0.00353673,0.00584728};
1085 >      double y[4]   = {1.01177,1.01568,1.01616,1.03211};
1086 >      double dyl[4] = {0.00217317,0.00203014,0.00276838,0.00465643};
1087 >      double dyh[4] = {0.00217317,0.00203014,0.00276838,0.00465643};
1088        gr[i] = new TGraphAsymmErrors(4,x,y,dxl,dxh,dyl,dyh);
1089 <      name="20/15";
1089 >      name="CHS 30/30";
1090        color="#c72f3c";
1091      }
1092      if(i==4) {
1093 <      double y[4]   = {1.01221+0.001,1.01922+0.001,1.02396+0.001,1.04677+0.001};
1094 <      double dyl[4] = {0.00263821,0.00260853,0.00353673,0.00584728};
1095 <      double dyh[4] = {0.00263821,0.00260853,0.00353673,0.00584728};
1093 >      double y[4]   = {1.00751,1.00964,1.00804,1.01073};
1094 >      double dyl[4] = {0.00224851,0.00215563,0.0029592,0.00499805};
1095 >      double dyh[4] = {0.00224851,0.00215563,0.0029592,0.00499805};
1096        gr[i] = new TGraphAsymmErrors(4,x,y,dxl,dxh,dyl,dyh);
1097 <      name="20/20";
1097 >      name="CHS 30/15";
1098        color="#7cb433";
1099      }
1100 +    if(i==5) {
1101 +      double y[4]   = {1.01557,1.0318,1.03535,1.06331};
1102 +      double dyl[4] = {0.00300783,0.00302486,0.00421064,0.00739416};
1103 +      double dyh[4] = {0.00300783,0.00302486,0.00421064,0.00739416};
1104 +      gr[i] = new TGraphAsymmErrors(4,x,y,dxl,dxh,dyl,dyh);
1105 +      name="CHS 15/15";
1106 +      color="#DF013A";
1107 +    }
1108      
1109 <    
1109 >    gStyle->SetOptFit(0);
1110      gr[i]->SetTitle();
1111      gr[i]->GetXaxis()->SetTitle("N(Vtx)");
1112      gr[i]->GetYaxis()->SetTitle("<MPF>_{data}/<MPF>_{mc}");
# Line 761 | Line 1122 | void ScenarioComparisonInclusive() {
1122        fit[i] = (TF1*)gr[i]->GetFunction("pol1");
1123        fit[i]->SetLineColor(TColor::GetColor(color.c_str()));
1124      }
1125 <    leg->AddEntry(gr[i],name.c_str(),"l");
1125 >    
1126 >    stringstream nameplus;
1127 >    nameplus << name << " (slope: " << std::setprecision(3) << fit[i]->GetParameter(1) << ")";
1128 >    
1129 >    leg->AddEntry(gr[i],nameplus.str().c_str(),"l");
1130      if(i==0) gr[i]->Draw("AP*");
1131      else gr[i]->Draw("P");
1132    }
# Line 777 | Line 1142 | void compare_selection(string identifier
1142    bool recognized_scenario=false;
1143    
1144    cout << "Running with identifier " << identifier << endl;
1145 +  if(identifier=="Zb30") {
1146 +      LeadingB=TCut ("Zb30_bTagProbCSVBP[0]>0.898");
1147 +      EtaB=TCut("abs(Zb30_pfJetGoodEta[0])<1.3");
1148 +      PhiZcut=TCut("abs(Zb30_pfJetDphiZ[0])>2.7");
1149 +      recognized_scenario=true;
1150 +  }
1151    if(identifier=="Zb3010") {
1152        LeadingB=TCut ("Zb3010_bTagProbCSVBP[0]>0.898");
1153        EtaB=TCut("abs(Zb3010_pfJetGoodEta[0])<1.3");
# Line 784 | Line 1155 | void compare_selection(string identifier
1155        recognized_scenario=true;
1156    }
1157    
1158 <  if(identifier=="Zb3010") {
1159 <      LeadingB=TCut ("Zb3010_bTagProbCSVBP[0]>0.898");
1160 <      EtaB=TCut("abs(Zb3010_pfJetGoodEta[0])<1.3");
1161 <      PhiZcut=TCut("abs(Zb3010_pfJetDphiZ[0])>2.7");
1158 >  if(identifier=="ZbCHS3010") {
1159 >      LeadingB=TCut ("ZbCHS1510_bTagProbCSVBP[0]>0.898");
1160 >      EtaB=TCut("abs(ZbCHS1510_pfJetGoodEta[0])<1.3");
1161 >      PhiZcut=TCut("abs(ZbCHS1510_pfJetDphiZ[0])>2.7");
1162        recognized_scenario=true;
1163    }
1164    
# Line 798 | Line 1169 | void compare_selection(string identifier
1169        recognized_scenario=true;
1170    }
1171    
1172 <  if(identifier=="Zb301030") {
1173 <      LeadingB=TCut ("Zb301030_bTagProbCSVBP[0]>0.898");
1174 <      EtaB=TCut("abs(Zb301030_pfJetGoodEta[0])<1.3");
1175 <      PhiZcut=TCut("abs(Zb301030_pfJetDphiZ[0])>2.7");
1172 >  if(identifier=="ZbCHS301030") {
1173 >      LeadingB=TCut ("ZbCHS301030_bTagProbCSVBP[0]>0.898");
1174 >      EtaB=TCut("abs(ZbCHS301030_pfJetGoodEta[0])<1.3");
1175 >      PhiZcut=TCut("abs(ZbCHS301030_pfJetDphiZ[0])>2.7");
1176        recognized_scenario=true;
1177    }
1178    
# Line 812 | Line 1183 | void compare_selection(string identifier
1183        recognized_scenario=true;
1184    }
1185    
1186 <  if(identifier=="Zb301010") {
1187 <      LeadingB=TCut ("Zb301010_bTagProbCSVBP[0]>0.898");
1188 <      EtaB=TCut("abs(Zb301010_pfJetGoodEta[0])<1.3");
1189 <      PhiZcut=TCut("abs(Zb301010_pfJetDphiZ[0])>2.7");
1186 >  if(identifier=="ZbCHS301010") {
1187 >      LeadingB=TCut ("ZbCHS301010_bTagProbCSVBP[0]>0.898");
1188 >      EtaB=TCut("abs(ZbCHS301010_pfJetGoodEta[0])<1.3");
1189 >      PhiZcut=TCut("abs(ZbCHS301010_pfJetDphiZ[0])>2.7");
1190        recognized_scenario=true;
1191    }
1192    
# Line 826 | Line 1197 | void compare_selection(string identifier
1197        recognized_scenario=true;
1198    }
1199    
1200 <  if(identifier=="Zb301010") {
1201 <      LeadingB=TCut ("Zb301010_bTagProbCSVBP[0]>0.898");
1202 <      EtaB=TCut("abs(Zb301010_pfJetGoodEta[0])<1.3");
1203 <      PhiZcut=TCut("abs(Zb301010_pfJetDphiZ[0])>2.7");
1200 >  if(identifier=="ZbCHS301010") {
1201 >      LeadingB=TCut ("ZbCHS301010_bTagProbCSVBP[0]>0.898");
1202 >      EtaB=TCut("abs(ZbCHS301010_pfJetGoodEta[0])<1.3");
1203 >      PhiZcut=TCut("abs(ZbCHS301010_pfJetDphiZ[0])>2.7");
1204        recognized_scenario=true;
1205    }
1206    
# Line 840 | Line 1211 | void compare_selection(string identifier
1211        recognized_scenario=true;
1212    }
1213    
1214 <  if(identifier=="Zb3010_SecEta3") {
1215 <      LeadingB=TCut ("Zb3010_SecEta3_bTagProbCSVBP[0]>0.898");
1216 <      EtaB=TCut("abs(Zb3010_SecEta3_pfJetGoodEta[0])<1.3");
1217 <      PhiZcut=TCut("abs(Zb3010_SecEta3_pfJetDphiZ[0])>2.7");
1214 >  if(identifier=="ZbCHS3010_SecEta3") {
1215 >      LeadingB=TCut ("ZbCHS3010_SecEta3_bTagProbCSVBP[0]>0.898");
1216 >      EtaB=TCut("abs(ZbCHS3010_SecEta3_pfJetGoodEta[0])<1.3");
1217 >      PhiZcut=TCut("abs(ZbCHS3010_SecEta3_pfJetDphiZ[0])>2.7");
1218        recognized_scenario=true;
1219    }
1220    
1221 <  if(identifier=="Zb3010_SecEta5") {
1222 <      LeadingB=TCut ("Zb3010_SecEta5_bTagProbCSVBP[0]>0.898");
1223 <      EtaB=TCut("abs(Zb3010_SecEta5_pfJetGoodEta[0])<1.3");
1224 <      PhiZcut=TCut("abs(Zb3010_SecEta5_pfJetDphiZ[0])>2.7");
1221 >  if(identifier=="ZbCHS3010_SecEta5") {
1222 >      LeadingB=TCut ("ZbCHS3010_SecEta5_bTagProbCSVBP[0]>0.898");
1223 >      EtaB=TCut("abs(ZbCHS3010_SecEta5_pfJetGoodEta[0])<1.3");
1224 >      PhiZcut=TCut("abs(ZbCHS3010_SecEta5_pfJetDphiZ[0])>2.7");
1225        recognized_scenario=true;
1226    }
1227    
1228 <  if(identifier=="Zb3010_p5Clean") {
1229 <      LeadingB=TCut ("Zb3010_p5Clean_bTagProbCSVBP[0]>0.898 && Zb3010_p5Clean_IsClean");
1230 <      EtaB=TCut("abs(Zb3010_p5Clean_pfJetGoodEta[0])<1.3 && Zb3010_p5Clean_IsClean");
1231 <      PhiZcut=TCut("abs(Zb3010_p5Clean_pfJetDphiZ[0])>2.7 && Zb3010_p5Clean_IsClean");
1228 >  if(identifier=="ZbCHS3010_p5Clean") {
1229 >      LeadingB=TCut ("ZbCHS3010_p5Clean_bTagProbCSVBP[0]>0.898 && ZbCHS3010_p5Clean_IsClean");
1230 >      EtaB=TCut("abs(ZbCHS3010_p5Clean_pfJetGoodEta[0])<1.3 && ZbCHS3010_p5Clean_IsClean");
1231 >      PhiZcut=TCut("abs(ZbCHS3010_p5Clean_pfJetDphiZ[0])>2.7 && ZbCHS3010_p5Clean_IsClean");
1232        recognized_scenario=true;
1233    }
1234    
1235 <  if(identifier=="Zb3010CHS") {
1236 <      LeadingB=TCut ("Zb3010CHS_bTagProbCSVBP[0]>0.898");
1237 <      EtaB=TCut("abs(Zb3010CHS_pfJetGoodEta[0])<1.3");
1238 <      PhiZcut=TCut("abs(Zb3010CHS_pfJetDphiZ[0])>2.7");
1235 >  if(identifier=="ZbCHS3010") {
1236 >      LeadingB=TCut ("ZbCHS3010_bTagProbCSVBP[0]>0.898");
1237 >      EtaB=TCut("abs(ZbCHS3010_pfJetGoodEta[0])<1.3");
1238 >      PhiZcut=TCut("abs(ZbCHS3010_pfJetDphiZ[0])>2.7");
1239 >      recognized_scenario=true;
1240 >  }
1241 >  if(identifier=="ZbCHS3010") {
1242 >      LeadingB=TCut ("ZbCHS3010_bTagProbCSVBP[0]>0.898");
1243 >      EtaB=TCut("abs(ZbCHS3010_pfJetGoodEta[0])<1.3");
1244 >      PhiZcut=TCut("abs(ZbCHS3010_pfJetDphiZ[0])>2.7");
1245 >      recognized_scenario=true;
1246 >  }
1247 >  if(identifier=="ZbCHS30") {
1248 >      LeadingB=TCut ("ZbCHS30_bTagProbCSVBP[0]>0.898");
1249 >      EtaB=TCut("abs(ZbCHS30_pfJetGoodEta[0])<1.3");
1250 >      PhiZcut=TCut("abs(ZbCHS30_pfJetDphiZ[0])>2.7");
1251        recognized_scenario=true;
1252    }
1253    
# Line 876 | Line 1259 | void compare_selection(string identifier
1259    cout << "   " << (const char*) PhiZcut << endl;
1260    
1261    cout << endl << endl;
1262 +
1263 + //  ScenarioComparison();
1264 + //  ScenarioComparisonInclusive();
1265 +
1266 + //    DoPUStudy(identifier);
1267 + //    SpecialCutFlow(identifier);
1268    
1269 < //    ScenarioComparison();
1270 < //    ScenarioComparisonInclusive();
1271 <  
1272 <  
884 < //     DoPUStudy(identifier);
885 <  
886 < //   SpecialCutFlow(identifier);
887 < //    draw_kin_variable("pt",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000",17000,30,200,"Z p_{T}","Official/"+identifier+"/Zpt__"+identifier,1);
888 < //   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);
889 < //   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);
890 < //   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);
1269 > //    draw_kin_variable("pt",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>30&&pt<1000",17000,30,200,"Z p_{T}","Official/"+identifier+"/Zpt__"+identifier,1);
1270 > //    draw_kin_variable(identifier+"_pfJetGoodPt[0]",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>30&&pt<1000",40,0,200,"Leading Jet Pt (B)","Official/"+identifier+"/LeadingJetPt__"+identifier,1);
1271 > //    draw_kin_variable(identifier+"_pfJetGoodPt[1]",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>30&&pt<1000",40,0,200,"Sub-Leading Jet Pt [GeV]","Official/"+identifier+"/Jet2Pt__"+identifier,1);
1272 > //    draw_kin_variable(identifier+"_pfJetGoodPt[1]/pt",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>30&&pt<1000",20,0,2,"p_{T}^{J2} / p_{T}^{Z}","Official/"+identifier+"/SubLeadingJetPt_Over_ZPt__"+identifier,1);
1273    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);
1274    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);
1275    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);
1276    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);
1277    draw_kin_variable(identifier+"_alpha",ZplusBsel&&LeadingB&&EtaB&&PhiZcut,20,0,2,"#alpha","Official/"+identifier+"/alpha_full__"+identifier,1);
1278 < //   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);
1278 >   draw_kin_variable(identifier+"_pfBJetDphiZ[0]",ZplusBsel&&LeadingB&&"pt>30&&pt<1000&&pfBJetDphiZ[0]>0",30,0,3.2,"#delta#phi (Z,b lead)","Official/"+identifier+"/DeltaPhiZBlead__"+identifier,0);
1279    
1280 < //   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");
1281 < //   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");
1282 <
1280 >   draw_kin_variable("mpf",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut(((string)"pt>30&&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, 30 GeV < p_{T}^{Z} < 1000 GeV");
1281 >   draw_kin_variable(identifier+"_pfJetGoodPt[0]/pt",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut(((string)"pt>30&&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, 30 GeV < p_{T}^{Z} < 1000 GeV");
1282 >   draw_kin_variable(identifier+"_pfJetGoodPt[0]/pt",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut(((string)"pt>30&&pt<50&&"+identifier+"_alpha<0.3").c_str()),20,0,2,"p_{T} b-jet / p_{T} Z","Official/"+identifier+"/ptb_over_ptz___alpha_smaller_0p3___pt3050___"+identifier,0,"#alpha<0.3, 30 GeV < p_{T}^{Z} < 1000 GeV}{30<p_{T}^{Z}<50 GeV}");
1283    
1284 < //   */
1284 > }
1285 >
1286 > void GetNumberEventsInsideOutsideAlphaWindow(TCut specialcut, string title) {
1287    
1288 + //  TCut cut =  ZplusBsel&&LeadingB&&PhiZcut&&specialcut&&TCut("ZbCHS3010_alpha<0.3");
1289 +  TCut cut = specialcut;
1290 +  TCut in  =  EtaB;
1291 +  TCut out =  TCut("abs(ZbCHS3010_pfJetGoodEta[0])>1.3");
1292 +
1293 +  TH1F *data_IN  = allsamples.Draw("data_IN",    "mll",1,0,150, "nothing [GeV]", "events", cut&&in, data,luminosity);
1294 +  TH1F *data_OUT = allsamples.Draw("data_OUT",   "mll",1,0,150, "nothing [GeV]", "events", cut&&out,data,luminosity);
1295 +  TH1F *mc_IN    = allsamples.Draw("mc_IN",      "mll",1,0,150, "nothing [GeV]", "events", cut&&in, mc  ,luminosity);
1296 +  TH1F *mc_OUT   = allsamples.Draw("mc_OUT",     "mll",1,0,150, "nothing [GeV]", "events", cut&&out,mc  ,luminosity);
1297 +  
1298 +  dout << title << " : " << endl;
1299 +  dout << "  Data: In " << data_IN->Integral() << " vs out " << data_OUT->Integral() << endl;
1300 +  dout << "  MC  : In " << mc_IN->Integral()   << " vs out " << mc_OUT->Integral()   << endl;  
1301 +  dout << "       ratio in  : " << data_IN->Integral()/mc_IN->Integral() << " +/- " << (data_IN->Integral()/mc_IN->Integral()) * sqrt(1.0/mc_IN->Integral()+ 1.0/data_IN->Integral()) << endl;
1302 +  dout << "       ratio out : " << data_OUT->Integral()/mc_OUT->Integral() << " +/- " << (data_OUT->Integral()/mc_OUT->Integral()) * sqrt(1.0/mc_OUT->Integral()+ 1.0/data_OUT->Integral()) << endl;
1303 +  
1304 +  delete data_IN;
1305 +  delete data_OUT;
1306 +  delete mc_IN;
1307 +  delete mc_OUT;
1308   }
1309  
1310 < void do_basic_ZB_analysis(int do_inclusive) {
1310 > void TEST_GetNumberEventsInsideOutsideAlphaWindow() {
1311 >  TCanvas *ca = new TCanvas("ca","ca");
1312 > //   GetNumberEventsInsideOutsideAlphaWindow(TCut("id1==id2&&id1==0"), "ee");
1313 > //   GetNumberEventsInsideOutsideAlphaWindow(TCut("id1==id2&&id1==1"), "mm");
1314 > //   GetNumberEventsInsideOutsideAlphaWindow(TCut("id1==id2"), "ee/mm");
1315 >  
1316 >  cout << "BASE SELECTION" << endl;
1317 >  GetNumberEventsInsideOutsideAlphaWindow(TCut("id1==id2&&id1==0"), "ee");
1318 >  GetNumberEventsInsideOutsideAlphaWindow(TCut("id1==id2&&id1==1"), "mm");
1319 >  GetNumberEventsInsideOutsideAlphaWindow(TCut("id1==id2"), "ee/mm");
1320 >
1321 >  cout << "BASE SELECTION: Z+b" << endl;
1322 >  GetNumberEventsInsideOutsideAlphaWindow(ZplusBsel&&TCut("id1==id2&&id1==0"), "ee");
1323 >  GetNumberEventsInsideOutsideAlphaWindow(ZplusBsel&&TCut("id1==id2&&id1==1"), "mm");
1324 >  GetNumberEventsInsideOutsideAlphaWindow(ZplusBsel&&TCut("id1==id2"), "ee/mm");
1325 >  
1326 >  cout << "Leading B" << endl;
1327 >  GetNumberEventsInsideOutsideAlphaWindow(LeadingB&&ZplusBsel&&TCut("id1==id2&&id1==0"), "ee");
1328 >  GetNumberEventsInsideOutsideAlphaWindow(LeadingB&&ZplusBsel&&TCut("id1==id2&&id1==1"), "mm");
1329 >  GetNumberEventsInsideOutsideAlphaWindow(LeadingB&&ZplusBsel&&TCut("id1==id2"), "ee/mm");
1330 >  
1331 >  cout << "PhiZcut" << endl;
1332 >  GetNumberEventsInsideOutsideAlphaWindow(PhiZcut&&LeadingB&&ZplusBsel&&TCut("id1==id2&&id1==0"), "ee");
1333 >  GetNumberEventsInsideOutsideAlphaWindow(PhiZcut&&LeadingB&&ZplusBsel&&TCut("id1==id2&&id1==1"), "mm");
1334 >  GetNumberEventsInsideOutsideAlphaWindow(PhiZcut&&LeadingB&&ZplusBsel&&TCut("id1==id2"), "ee/mm");
1335 >  
1336 >  cout << "alpha cut" << endl;
1337 >  GetNumberEventsInsideOutsideAlphaWindow(PhiZcut&&LeadingB&&ZplusBsel&&TCut("id1==id2&&id1==0")&&TCut("ZbCHS3010_alpha<0.3"), "ee");
1338 >  GetNumberEventsInsideOutsideAlphaWindow(PhiZcut&&LeadingB&&ZplusBsel&&TCut("id1==id2&&id1==1")&&TCut("ZbCHS3010_alpha<0.3"), "mm");
1339 >  GetNumberEventsInsideOutsideAlphaWindow(PhiZcut&&LeadingB&&ZplusBsel&&TCut("id1==id2")&&TCut("ZbCHS3010_alpha<0.3"), "ee/mm");
1340    
1341 <  //https://twiki.cern.ch/twiki/bin/viewauth/CMS/BtagPOG#2011_Data_and_MC
909 <  //https://twiki.cern.ch/twiki/pub/CMS/BtagPOG/SFb-ttbar_payload.txt
1341 >  delete ca;
1342    
1343 <  if(do_inclusive) {
1344 <    LeadingB=TCut("Zb3010_bTagProbCSVBP[0]>=-10000");
1345 <    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)))");
1346 <    write_info(__FUNCTION__,"Inclusive selection - no btag requirement, no btag efficiency/mistag correction (weight)");
1347 <  } else {
1348 <    LeadingB=TCut("Zb3010_bTagProbCSVBP[0]>0.898");
1349 <    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))");
1350 <    write_info(__FUNCTION__,"Exclusive selection - switched on btag requirement and btag efficiency/mistag correction (weight)");
1351 <  }
1352 <  
1353 <  bool doquick=true;
1354 <  
1355 <  TCanvas *zbcanvas = new TCanvas("zbcanvas","zbcanvas");
1356 <  
1357 < //     compare_selection("Zb3010_p5Clean");
1358 < //     compare_selection("Zb3010CHS");
1359 < //     compare_selection("Zb3010");
1360 < //     compare_selection("Zb301030");
1361 < //     compare_selection("Zb1530");
1362 < //     compare_selection("Zb301010");
1363 < //     compare_selection("Zb3010_SecEta3");
1364 < //     compare_selection("Zb3010_SecEta5");
1365 < //     compare_selection("Zb1510");
1366 < //     compare_selection("Zb301010");
1367 < //     compare_selection("ZbMikko");
1368 < //     compare_selection("Zb3010");
1369 < //     compare_selection("Zb40");
1343 >  
1344 > }
1345 >  
1346 >
1347 >
1348 >    
1349 > void do_basic_ZB_analysis(bool DoCompleteAnalysis) {
1350 >  STDWEIGHT=TCut(cutWeight);
1351 >  cutWeight=TCut(STDWEIGHT*TCut("ZbCHS3010_BTagWgtT"));
1352 >  cout << "The lepton requirement is " << (const char*) leptoncut << endl;
1353 >  bool doquick=false;
1354 >  
1355 >   TCanvas *zbcanvas = new TCanvas("zbcanvas","zbcanvas");
1356 >  
1357 > //   ScenarioComparison();
1358 > //   ScenarioComparisonInclusive();
1359 >
1360 > //   compare_selection("ZbCHS3010_p5Clean");
1361 > //   compare_selection("ZbCHS3010");
1362 > //   compare_selection("Zb30");
1363 > //   compare_selection("Zb3010");
1364 > //   compare_selection("Zb1510");
1365 > //   compare_selection("ZbCHS30");
1366 > //   compare_selection("ZbCHS3010");
1367 > //   compare_selection("ZbCHS3010");
1368 > //   compare_selection("ZbCHS301010");
1369 > //   compare_selection("ZbCHS3010_SecEta3");
1370 > //   compare_selection("ZbCHS3010_SecEta5");
1371 > //   compare_selection("Zb1510");
1372 > //   compare_selection("ZbCHS301010");
1373 > //   compare_selection("ZbMikko");
1374 > //   compare_selection("ZbCHS3010");
1375 > //   compare_selection("Zb40");
1376  
1377    if(!doquick) {
1378 <    print_all_b_yields();
1378 >    //draw_kin_variable("mll","",28,60.0,200.0,"m_{ll}","mll",1,"");//nice and quick test to see if the normalization is ok (i.e. all data is processed etc.)
1379 >    //print_all_b_yields();
1380      draw_mpf_vars();
1381 <    //  DrawEvilCutFlow();
1381 >    //DrawEvilCutFlow();
1382          
1383      draw_Zb_kin_vars();
1384          
1385     }
1386    
1387 <  new_data_mc_agreement_2d();
1387 > //  GetNumberEventsInsideOutsideAlphaWindow();
1388 >  
1389 >  if(DoCompleteAnalysis) {
1390 >    /* need to run the following scenarios:
1391 >     * - normal
1392 >     * - inclusive
1393 >     * - medium WP
1394 >     * - alpha up
1395 >     * - alpha down
1396 >     */
1397 >    
1398 >    bool fast=true;bool slow=false;//don't change these
1399 >    
1400 >    ZbResultContainer inclusive = new_data_mc_agreement_2d(fast,"ZbCHS3010_alpha","inclusive",TCut("mll>0"),TCut("1.0"),0.3);
1401 >    ZbResultContainer Reference = new_data_mc_agreement_2d(slow,"ZbCHS3010_alpha","reference",TCut("ZbCHS3010_bTagProbCSVBP[0]>0.898"),TCut("ZbCHS3010_BTagWgtT"),0.3);
1402 >    
1403 >    ZbResultContainer NeutrinoQ = new_data_mc_agreement_2d(fast,"ZbCHS3010_alpha","NeutrinoQ",TCut("ZbCHS3010_bTagProbCSVBP[0]>0.898"),TCut("ZbCHS3010_BTagWgtTUp*(1.0*softMuon)"),0.3);
1404 >    ZbResultContainer ANeutrinoQ = new_data_mc_agreement_2d(fast,"ZbCHS3010_alpha","ANeutrinoQ",TCut("ZbCHS3010_bTagProbCSVBP[0]>0.898"),TCut("ZbCHS3010_BTagWgtTUp*(1.0*(1.0-softMuon))"),0.3);
1405 >    
1406 >    ZbResultContainer effmistagUp = new_data_mc_agreement_2d(fast,"ZbCHS3010_alpha","effmistagUp",TCut("ZbCHS3010_bTagProbCSVBP[0]>0.898"),TCut("ZbCHS3010_BTagWgtTUp"),0.3);
1407 >    ZbResultContainer effmistagDown = new_data_mc_agreement_2d(fast,"ZbCHS3010_alpha","effmistagDown",TCut("ZbCHS3010_bTagProbCSVBP[0]>0.898"),TCut("ZbCHS3010_BTagWgtTDown"),0.3);
1408 >    ZbResultContainer medium = new_data_mc_agreement_2d(fast,"ZbCHS3010_alpha","medium",TCut("ZbCHS3010_bTagProbCSVBP[0]>0.679"),TCut("ZbCHS3010_BTagWgtM"),0.3);
1409 >    ZbResultContainer loose = new_data_mc_agreement_2d(fast,"ZbCHS3010_alpha","loose",TCut("ZbCHS3010_bTagProbCSVBP[0]>0.244"),TCut("ZbCHS3010_BTagWgtL"),0.3);
1410 >    ZbResultContainer AlphaUp = new_data_mc_agreement_2d(fast,"ZbCHS3010_alphaUp","AlphaUp",TCut("ZbCHS3010_bTagProbCSVBP[0]>0.898"),TCut("ZbCHS3010_BTagWgtT"),0.3);
1411 >    ZbResultContainer AlphaDown = new_data_mc_agreement_2d(fast,"ZbCHS3010_alphaDown","AlphaDown",TCut("ZbCHS3010_bTagProbCSVBP[0]>0.898"),TCut("ZbCHS3010_BTagWgtT"),0.3);
1412 >    
1413 >    ZbResultContainer L5inclusive = new_data_mc_agreement_2d(fast,"ZbCHS3010_alphaL5","L5inclusive",TCut("mll>0"),TCut("1.0"),0.3);
1414 >    ZbResultContainer L5Reference = new_data_mc_agreement_2d(fast,"ZbCHS3010_alphaL5","L5reference",TCut("ZbCHS3010_bTagProbCSVBP[0]>0.898"),TCut("ZbCHS3010_BTagWgtT"),0.3);
1415 >    
1416 >    ZbResultContainer AlphaThresholdVariationDown10 =    new_data_mc_agreement_2d(fast,"ZbCHS3010_alpha","AlphaThresholdVariationDown10",   TCut("ZbCHS3010_bTagProbCSVBP[0]>0.898"),   TCut("ZbCHS3010_BTagWgtT"),0.3-0.1*0.3);
1417 >    ZbResultContainer AlphaThresholdVariationUp10 =      new_data_mc_agreement_2d(fast,"ZbCHS3010_alpha","AlphaThresholdVariationUp10",     TCut("ZbCHS3010_bTagProbCSVBP[0]>0.898"),   TCut("ZbCHS3010_BTagWgtT"),0.3+0.1*0.3);
1418 >    ZbResultContainer AlphaThresholdVariationDown30 =    new_data_mc_agreement_2d(fast,"ZbCHS3010_alpha","AlphaThresholdVariationDown30",   TCut("ZbCHS3010_bTagProbCSVBP[0]>0.898"),   TCut("ZbCHS3010_BTagWgtT"),0.3-0.3*0.3);
1419 >    ZbResultContainer AlphaThresholdVariationUp30 =      new_data_mc_agreement_2d(fast,"ZbCHS3010_alpha","AlphaThresholdVariationUp30",     TCut("ZbCHS3010_bTagProbCSVBP[0]>0.898"),   TCut("ZbCHS3010_BTagWgtT"),0.3+0.3*0.3);
1420 >    ZbResultContainer AlphaThresholdVariationDown10inc = new_data_mc_agreement_2d(fast,"ZbCHS3010_alpha","AlphaThresholdVariationDown10inc",TCut("mll>0"),                              TCut("1.0"),0.3-0.1*0.3);
1421 >    ZbResultContainer AlphaThresholdVariationUp10inc =   new_data_mc_agreement_2d(fast,"ZbCHS3010_alpha","AlphaThresholdVariationUp10inc",  TCut("mll>0"),                              TCut("1.0"),0.3+0.1*0.3);
1422 >    ZbResultContainer AlphaThresholdVariationDown30inc = new_data_mc_agreement_2d(fast,"ZbCHS3010_alpha","AlphaThresholdVariationDown30inc",TCut("mll>0"),                              TCut("1.0"),0.3-0.3*0.3);
1423 >    ZbResultContainer AlphaThresholdVariationUp30inc =   new_data_mc_agreement_2d(fast,"ZbCHS3010_alpha","AlphaThresholdVariationUp30inc",  TCut("mll>0"),                              TCut("1.0"),0.3+0.3*0.3);
1424 >
1425 >    //and now let's compute some systematics!
1426 >    //1 alpha variation
1427 >    float SysMPF_Alpha_down = Reference.MPF_Result_FaceValue.getValue()-AlphaDown.MPF_Result_FaceValue.getValue();
1428 >    float SysMPF_Alpha_up = Reference.MPF_Result_FaceValue.getValue()-AlphaUp.MPF_Result_FaceValue.getValue();
1429 >    
1430 >    //2 btagging efficiency/mistag correction
1431 >    float SysMPF_EFF_down = Reference.MPF_Result_FaceValue.getValue()-effmistagDown.MPF_Result_FaceValue.getValue();
1432 >    float SysMPF_EFF_up = Reference.MPF_Result_FaceValue.getValue()-effmistagUp.MPF_Result_FaceValue.getValue();
1433 >    
1434 >    //3 Presence of soft muons (->neutrino question)
1435 >    float SysMPF_Neutrino_up = Reference.MPF_Result_FaceValue.getValue()-NeutrinoQ.MPF_Result_FaceValue.getValue();
1436 >    float SysMPF_Neutrino_down = Reference.MPF_Result_FaceValue.getValue()-ANeutrinoQ.MPF_Result_FaceValue.getValue();
1437 >    
1438 >    //4 purity
1439 >    zbcanvas->cd();
1440 >    gStyle->SetOptFit(0);
1441 >    write_info(__FUNCTION__,"Don't know the exact purity at the moment, still needs to be determined!");
1442 >    TH1F *purityhisto = new TH1F("purityhisto","purityhisto",1,0,1);
1443 >    purityhisto->GetXaxis()->SetTitle("Purity");
1444 >    purityhisto->GetYaxis()->SetTitle("C_{abs}");
1445 >    purityhisto->GetXaxis()->CenterTitle();
1446 >    purityhisto->GetYaxis()->CenterTitle();
1447 >    
1448 >    TGraphErrors *gSysMPF_Purity = new TGraphErrors();
1449 >    
1450 >    gSysMPF_Purity->SetPoint(0,PurityInclusive,inclusive.MPF_Result_FaceValue.getValue());
1451 >    gSysMPF_Purity->SetPointError(0,0.0,inclusive.MPF_Result_FaceValue.getError());
1452 >    
1453 >    gSysMPF_Purity->SetPoint(1,PurityLoose,loose.MPF_Result_FaceValue.getValue());//x value?
1454 >    gSysMPF_Purity->SetPointError(1,0.0,loose.MPF_Result_FaceValue.getError());
1455 >    gSysMPF_Purity->SetPoint(2,PurityMedium,medium.MPF_Result_FaceValue.getValue());//x value?
1456 >    gSysMPF_Purity->SetPointError(2,0.0,medium.MPF_Result_FaceValue.getError());
1457 >    gSysMPF_Purity->SetPoint(3,PurityTight,Reference.MPF_Result_FaceValue.getValue());//x value?
1458 >    gSysMPF_Purity->SetPointError(3,0.0,Reference.MPF_Result_FaceValue.getError());
1459 >    
1460 >    
1461 >    float min,max;
1462 >    FindMinMax(gSysMPF_Purity,min,max);
1463 >    purityhisto->SetStats(0);
1464 >    purityhisto->GetYaxis()->SetRangeUser(min,max);
1465 >    
1466 >    gSysMPF_Purity->Fit("pol1");
1467 >    TPolyLine *purity_fit_uncert = GetFitUncertaintyShape(gSysMPF_Purity, "pol1", min, max,0.0,1.0);
1468 >    TF1* pufit = (TF1*)gSysMPF_Purity->GetFunction("pol1");
1469 >    
1470 >    TF1* pufit_c = new TF1("pufit_c","[0]+[1]*x");
1471 >    pufit_c->SetParameters(pufit->GetParameters());
1472 >    pufit_c->SetLineColor(TColor::GetColor("#0101DF"));
1473 >    purity_fit_uncert->SetFillColor(TColor::GetColor("#5353FC"));
1474 >    pufit->SetLineColor(TColor::GetColor("#5353FC"));
1475 >    
1476 >    purityhisto->Draw();
1477 >    purity_fit_uncert->Draw("f");
1478 >    gSysMPF_Purity->Draw("P*");
1479 >    pufit_c->Draw("same");
1480 >    DrawPrelim();
1481 >    
1482 >    float ExtrapolatedResult=pufit->GetParameter(0)+1.0*pufit->GetParameter(1);
1483 >    float dExtrapolatedResult=sqrt(pufit->GetParError(0)*pufit->GetParError(0)+1.0*1.0*pufit->GetParError(1)*pufit->GetParError(1));
1484 >    
1485 >    stringstream summary;
1486 >    summary << "#splitline{Reference: " << std::setprecision(4) << Reference.MPF_Result_FaceValue.getValue() << "+/-" << std::setprecision(4) << Reference.MPF_Result_FaceValue.getError() << "}{";
1487 >    summary << "Extrapolation: " << std::setprecision(4) << ExtrapolatedResult << " +/- " << std::setprecision(4) << dExtrapolatedResult << "}";
1488 >    
1489 >    dout << "Have found an extrapolated result at 1.0 of " << ExtrapolatedResult << " ;  will use the difference between this result and the face value as ";
1490 >    dout << "purity-related systematic uncertainty" << endl;
1491 >    
1492 >    TText *infobox = write_title(summary.str());
1493 >    infobox->SetX(0.75);
1494 >    infobox->SetTextSize(0.03);
1495 >    infobox->SetY(0.75);
1496 >    infobox->Draw();
1497 >    
1498 >    CompleteSave(zbcanvas,"Systematics/Purity");
1499 >    
1500 >    
1501 >
1502 >    float SysMPF_Purity = abs(ExtrapolatedResult-Reference.MPF_Result_FaceValue.getValue());
1503 >    
1504 >    float sys_alpha    = abs(SysMPF_Alpha_down)>abs(SysMPF_Alpha_up)?abs(SysMPF_Alpha_down):abs(SysMPF_Alpha_up);
1505 >    float sys_eff      = abs(SysMPF_EFF_down)>abs(SysMPF_EFF_up)?abs(SysMPF_EFF_down):abs(SysMPF_EFF_up);
1506 >    float sys_purity   = abs(SysMPF_Purity);
1507 >    float sys_neutrino = abs(SysMPF_Neutrino_down)>abs(SysMPF_Neutrino_up)?abs(SysMPF_Neutrino_down):abs(SysMPF_Neutrino_up);
1508 >    float sys          = sqrt(sys_alpha*sys_alpha+sys_eff*sys_eff+sys_purity*sys_purity+sys_neutrino*sys_neutrino);
1509 >    
1510 >    dout << "MPF METHOD:" << endl;
1511 >    dout << "   Systematics :    down    up     (selected)" << endl;
1512 >    dout << "      Oversmearing         : " << SysMPF_Alpha_down << "   ,  " << SysMPF_Alpha_up << "   ( " << sys_alpha << ") " << endl;
1513 >    dout << "      eff/mistag variation : " << SysMPF_EFF_down << "   ,  " << SysMPF_EFF_up << "   ( " << sys_eff << ") " << endl;
1514 >    dout << "      purity               : " << SysMPF_Purity << "   ( " << sys_purity << " )" << endl;
1515 >    dout << "      neutrino             : " << SysMPF_Neutrino_down << "   ,  " << SysMPF_Neutrino_up << "   ( " << sys_neutrino << ") " << endl;
1516 >    dout << "      total                : " << sys << "     (note that this is the sys on C_{abs}^{b}, not yet on C_{corr} - there will also be an additional source!)" << endl;
1517 >    dout << endl << endl;
1518 >    dout << "   FINAL RESULTS : " << endl;
1519 >    dout << "     C_{abs}^{b}   = " << Reference.MPF_Result_FaceValue << " (stat) +/- " << sys << " (sys) " << endl;
1520 >    dout << "     C_{abs}^{inc} = " << inclusive.MPF_Result_FaceValue << " (stat) " << endl;
1521 >    
1522 >    float sysfinal = Reference.MPF_Result_FaceValue.getValue()/inclusive.MPF_Result_FaceValue.getValue();
1523 >    sysfinal *= sqrt(sys*sys/(Reference.MPF_Result_FaceValue.getValue()*Reference.MPF_Result_FaceValue.getValue()));
1524 >    //yes, this could be simplified but it would look like a bug.
1525 >    
1526 >
1527 >    
1528 >    
1529 >    dout << "Note: Varying the alpha cut (->face value) we get the following results: " << endl;
1530 >    dout << "-30%: " << AlphaThresholdVariationDown30.MPF_Result_FaceValue/AlphaThresholdVariationDown30inc.MPF_Result_FaceValue << " (stat)" << endl;
1531 >    dout << "-10%: " << AlphaThresholdVariationDown10.MPF_Result_FaceValue/AlphaThresholdVariationDown10inc.MPF_Result_FaceValue << " (stat)" << endl;
1532 >    dout << "  0%: " << Reference.MPF_Result_FaceValue/inclusive.MPF_Result_FaceValue << " (stat) +/- " << endl;
1533 >    dout << "+10%: " << AlphaThresholdVariationUp10.MPF_Result_FaceValue/AlphaThresholdVariationUp10inc.MPF_Result_FaceValue << " (stat)" << endl;
1534 >    dout << "+30%: " << AlphaThresholdVariationUp30.MPF_Result_FaceValue/AlphaThresholdVariationUp30inc.MPF_Result_FaceValue << " (stat)" << endl;
1535 >    
1536 >    dout << "          XCheck for Rabs, varying the alpha cut (->face value) we get the following results: " << endl;
1537 >    dout << "             -30%: " << AlphaThresholdVariationDown30.Rabs_Result_FaceValue/AlphaThresholdVariationDown30inc.Rabs_Result_FaceValue << " (stat)" << endl;
1538 >    dout << "             -10%: " << AlphaThresholdVariationDown10.Rabs_Result_FaceValue/AlphaThresholdVariationDown10inc.Rabs_Result_FaceValue << " (stat)" << endl;
1539 >    dout << "               0%: " << Reference.Rabs_Result_FaceValue/inclusive.Rabs_Result_FaceValue << " (stat) +/- " << endl;
1540 >    dout << "             +10%: " << AlphaThresholdVariationUp10.Rabs_Result_FaceValue/AlphaThresholdVariationUp10inc.Rabs_Result_FaceValue << " (stat)" << endl;
1541 >    dout << "             +30%: " << AlphaThresholdVariationUp30.Rabs_Result_FaceValue/AlphaThresholdVariationUp30inc.Rabs_Result_FaceValue << " (stat)" << endl;
1542 >
1543 >    
1544 >    float sys_alphavar_up   = abs(AlphaThresholdVariationUp30.MPF_Result_FaceValue.getValue()/AlphaThresholdVariationUp30inc.MPF_Result_FaceValue.getValue() - Reference.MPF_Result_FaceValue.getValue()/inclusive.MPF_Result_FaceValue.getValue());
1545 >    float sys_alphavar_down = abs(AlphaThresholdVariationDown30.MPF_Result_FaceValue.getValue()/AlphaThresholdVariationDown30inc.MPF_Result_FaceValue.getValue() - Reference.MPF_Result_FaceValue.getValue()/inclusive.MPF_Result_FaceValue.getValue());
1546 >    float sys_alphavar = sys_alphavar_up>sys_alphavar_down?sys_alphavar_up:sys_alphavar_down;
1547 >    
1548 >    dout << " -> Additional systematic (from 30%) : " << sys_alphavar << endl;
1549 >    sysfinal=sqrt(sysfinal*sysfinal+sys_alphavar*sys_alphavar);
1550 >    dout << "Final systematic is : " << sysfinal << endl;
1551 >    
1552 >    dout << "   ******************************************** <  FINAL RESULT > ******************************************** " << endl;
1553 >    dout << "     C_{corr}      = " << Reference.MPF_Result_FaceValue/inclusive.MPF_Result_FaceValue << " (stat) +/- " << sysfinal << " (sys) " << endl;
1554 >    dout << "   ******************************************** < /FINAL RESULT > ******************************************** " << endl;
1555 >    
1556 >    dout << "Note that if L5 corrections are applied, we get " << L5Reference.MPF_Result_FaceValue/L5inclusive.MPF_Result_FaceValue << " (stat) " << endl;
1557 >    
1558 >    LeadingB=TCut("ZbCHS3010_bTagProbCSVBP[0]>0.898");
1559 >    
1560 >  } else {
1561 >    ZbResultContainer Reference = new_data_mc_agreement_2d(false,"ZbCHS3010_alpha","reference",TCut("ZbCHS3010_bTagProbCSVBP[0]>0.898"),TCut("ZbCHS3010_BTagWgtT"),0.3);
1562 >  }
1563 >  
1564    
950 –  write_error(__FUNCTION__,"Need to implement switch for inclusive result; inclusive result needs to have btag efficiency reweighting deactivated!");
1565    delete zbcanvas;
1566   }
1567  
# Line 972 | Line 1586 | void SpecialCutFlow(string identifier) {
1586    MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB) << ")";
1587    MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB&&EtaB) << ")";
1588    MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB&&EtaB&&PhiZcut) << ")";
1589 <  MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut(((string)"pt>10&&pt<1000").c_str())) << ")";
1590 <  MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut(((string)"pt>10&&pt<1000&&"+identifier+"_alpha<0.3").c_str())) << ")";
1591 <  MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut(((string)"pt>10&&pt<1000&&"+identifier+"_alpha<0.2").c_str())) << ")";
1592 <  MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut(((string)"pt>10&&pt<1000&&"+identifier+"_alpha<0.15").c_str())) << ")";
1593 <  MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut(((string)"pt>10&&pt<1000&&"+identifier+"_alpha<0.1").c_str())) << ")";
1594 <  MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut(((string)"pt>10&&pt<1000&&"+identifier+"_alpha<0.05").c_str())) << ")";
1589 >  MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut(((string)"pt>30&&pt<1000").c_str())) << ")";
1590 >  MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut(((string)"pt>30&&pt<1000&&"+identifier+"_alpha<0.3").c_str())) << ")";
1591 >  MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut(((string)"pt>30&&pt<1000&&"+identifier+"_alpha<0.2").c_str())) << ")";
1592 >  MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut(((string)"pt>30&&pt<1000&&"+identifier+"_alpha<0.15").c_str())) << ")";
1593 >  MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut(((string)"pt>30&&pt<1000&&"+identifier+"_alpha<0.1").c_str())) << ")";
1594 >  MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut(((string)"pt>30&&pt<1000&&"+identifier+"_alpha<0.05").c_str())) << ")";
1595    MegaCut << " ";
1596    
1597    TCanvas *can = new TCanvas("can","can");
# Line 990 | Line 1604 | void SpecialCutFlow(string identifier) {
1604    
1605    
1606    TH1F *data  = allsamples.Draw("data",    MegaCut.str(),11,-0.5,10.5, "", "events", basecut,data,luminosity);
1607 <   THStack mcm = allsamples.DrawStack("mc", MegaCut.str(),11,-0.5,10.5, "", "events", basecut,mc,luminosity);
1607 >  THStack mcm = allsamples.DrawStack("mc", MegaCut.str(),11,-0.5,10.5, "", "events", basecut,mc,luminosity);
1608  
1609    float runningsum=0;
1610    for(int i=data->GetNbinsX();i>0;i--) {
# Line 1028 | Line 1642 | void SpecialCutFlow(string identifier) {
1642    data->GetXaxis()->SetBinLabel(3,"Leading jet is b-jet");
1643    data->GetXaxis()->SetBinLabel(4,"|#eta_{b}|<1.3 ");
1644    data->GetXaxis()->SetBinLabel(5,"|#delta#phi(b,Z)|>2.7");
1645 <  data->GetXaxis()->SetBinLabel(6,"10<p_{T}^{Z}<1000 GeV");
1645 >  data->GetXaxis()->SetBinLabel(6,"30<p_{T}^{Z}<1000 GeV");
1646    data->GetXaxis()->SetBinLabel(7,"#alpha < 0.3");
1647    data->GetXaxis()->SetBinLabel(8,"#alpha < 0.2");
1648    data->GetXaxis()->SetBinLabel(9,"#alpha < 0.15");
# Line 1046 | Line 1660 | void SpecialCutFlow(string identifier) {
1660    leg->AddEntry(crap,neventsinfo.str().c_str(),"l");
1661    
1662    data->Draw("e1");
1663 <  mcm.Draw("same");
1663 >  mcm.Draw("histo,same");
1664    data->Draw("same,e1,axis");
1665    data->Draw("same,e1");
1666   //  data->Draw("same,TEXT");
# Line 1059 | Line 1673 | void SpecialCutFlow(string identifier) {
1673    
1674    delete crap;
1675    delete data;
1676 <  
1676 >  CleanLegends();
1677 >  DeleteStack(mcm);
1678 >
1679    
1680   }

Diff Legend

– Removed lines
+ Added lines
< Changed lines
> Changed lines