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

Comparing UserCode/cbrown/Development/Plotting/Modules/WZStudy.C (file contents):
Revision 1.2 by buchmann, Thu Oct 11 09:57:19 2012 UTC vs.
Revision 1.4 by buchmann, Wed Nov 7 10:52:31 2012 UTC

# Line 18 | Line 18
18  
19   using namespace std;
20  
21 < void CarryOutWZStudy() {
22 <    //we do not want the
21 > namespace TriLeptons {
22 >    float Rsfof_CorrectionFactor=-1;
23 >    float Rsfof_CorrectionFactor_Error=-1;
24 > }
25 >
26 > string IdentifierForUnderlinedEntry(string Entry) {
27 >  if(Contains(Entry,"t_bar_t__")) return "ttbar";
28 >  if(Contains(Entry,"Single_top__")) return "SingleTop";
29 >  if(Contains(Entry,"W_Jets__")) return "W+Jets";
30 >  if(Contains(Entry,"WZ__")) return "WZ";
31 >  if(Contains(Entry,"Z_Jets")) return "Z+Jets";
32 >  if(Contains(Entry,"DATA")) return "Data";
33 >  if(Contains(Entry,"2l2q_")) return "ZZ->2l2q";
34 >  if(Contains(Entry,"4l_")) return "ZZ->4l";
35 >  return Entry;
36 > }
37 >
38 > void WZ_kin_plot(string variable, int nbins, float min, float high, string xlabel, TCut cut, string saveas) {
39 >  TCanvas *can = new TCanvas("can","can");
40 >  TPad *pad = new TPad("pad","pad",0,0,1,1);
41 >  pad->cd();
42 >
43 >  THStack stack = allsamples.DrawStack("stack",variable,nbins,min,high,xlabel, "events", cut,mc, luminosity);
44 >  TH1F *hdata = allsamples.Draw("hdata",variable,nbins,min,high,xlabel, "events", cut,data, luminosity);
45 >
46 >  float StackMax=CollapseStack(stack)->GetMaximum();
47 >  if(StackMax>hdata->GetMaximum()) hdata->SetMaximum(1.2*StackMax);
48 >  hdata->Draw("e1");
49 >  stack.Draw("same");
50 >  hdata->Draw("e1,same");
51 >
52 >  TLegend *kinleg = allsamples.allbglegend();
53 >  kinleg->Draw();
54 >  DrawPrelim();
55    
56 +  float total=0;
57 +  float wz=0;
58 +  cout << "Numbers with variable " << variable << " and cut " << cut << ": " << endl;
59 +  TIter nextSF = TIter(stack.GetHists());
60 +  TH1F *h;
61 +  stringstream spurity;
62 +  int bracecounter=0;
63 +
64 +
65 +  while ( h = (TH1F*)nextSF() ) {
66 +    cout << IdentifierForUnderlinedEntry(h->GetName()) << " : " << h->Integral() << endl;
67 +    if(IdentifierForUnderlinedEntry(h->GetName())=="WZ") {
68 +      wz+=h->Integral();
69 +      total+=h->Integral();
70 +    } else {
71 +      if(IdentifierForUnderlinedEntry(h->GetName())!="Data") total+=h->Integral();
72 +    }
73 +    bracecounter++;
74 +    spurity << "#splitline{" << IdentifierForUnderlinedEntry(h->GetName()) << " : " << h->Integral() << "}{";
75 +  }
76 +  
77 +  spurity << "#splitline{Purity: " << wz/total << "}{N(events): " << total << "}";
78 +  for(int i=0;i<bracecounter;i++) spurity << "}";
79    
80 < //   string Zmass="((leptonId[0]==leptonId[1])*(leptonCharge[0]*leptonCharge[1]<0)*leptonPairMass[0]   +   ((leptonId[0]!=leptonId[1])||(leptonCharge[0]*leptonCharge[1]>0))*(leptonId[0]==leptonId[2])*(leptonCharge[0]*leptonCharge[2]<0)*leptonPairMass[1]   +    ((((leptonId[0]!=leptonId[1])||(leptonCharge[0]*leptonCharge[1]>0)))*(((leptonId[0]!=leptonId[2])||(leptonCharge[0]*leptonCharge[2]>0))))*(leptonId[1]==leptonId[2])*(leptonCharge[1]*leptonCharge[2]<0)*leptonPairMass[2])";
80 >  TPaveText *PurityPaveText = new TPaveText(0.60, 0.40,0.89, 0.6,"blNDC");
81 >  PurityPaveText->SetFillStyle(4000);
82 >  PurityPaveText->SetBorderSize(0);
83 >  PurityPaveText->SetFillColor(kWhite);
84 >  PurityPaveText->SetTextFont(42);
85 >  PurityPaveText->SetTextSize(0.021);
86 >  PurityPaveText->AddText(spurity.str().c_str());
87 >  PurityPaveText->Draw();
88  
27  string Zmass="((leptonId[0]==leptonId[1])*(leptonCharge[0]*leptonCharge[1]<0)*leptonPairMass[0]   +   ((leptonId[0]!=leptonId[1])||(leptonCharge[0]*leptonCharge[1]>0))*(leptonId[0]==leptonId[2])*(leptonCharge[0]*leptonCharge[2]<0)*leptonPairMass[1]";
89    
90 + //  TText *purity = new TText(0.8,0.6,spurity.str().c_str());
91 + //  purity->Draw();
92 +  cout << "DATA : " << hdata->Integral() << endl;
93    
94 <  string WZPt1="((leptonId[0]==leptonId[1])*(leptonCharge[0]*leptonCharge[1]<0)*leptonPt[0]   +   ((leptonId[0]!=leptonId[1])||(leptonCharge[0]*leptonCharge[1]>0))*(leptonId[0]==leptonId[2])*(leptonCharge[0]*leptonCharge[2]<0)*leptonPt[0]   +    ((((leptonId[0]!=leptonId[1])||(leptonCharge[0]*leptonCharge[1]>0)))*(((leptonId[0]!=leptonId[2])||(leptonCharge[0]*leptonCharge[2]>0))))*(leptonId[1]==leptonId[2])*(leptonCharge[1]*leptonCharge[2]<0)*leptonPt[1])";
95 <  //this is the pt of the first lepton of the Z
94 >  save_with_ratio(hdata,stack,pad->cd(),"WZ_study/"+saveas);
95 >  delete hdata;
96 >  delete can;
97 > }
98 >
99 > void SF_vs_OF_plots(string variable, int nbins, float min, float high, string xlabel, TCut cut, TCut sf, TCut of, int mcordata, string saveas) {
100 >    TCanvas *can = new TCanvas("can","can");
101 >    TPad *pad = new TPad("pad","pad",0,0,1,1);
102 >    pad->cd();
103 >    
104 >    TLegend *kinleg = make_legend();
105 >    if(mcordata==data) kinleg->SetHeader("Data");
106 >    else kinleg->SetHeader("MC");
107 >    
108 >    TH1F *sfdata = allsamples.Draw("sfdata",variable,nbins,min,high,xlabel, "events", cut&&sf,mcordata, luminosity);
109 >    TH1F *ofdata = allsamples.Draw("ofdata",variable,nbins,min,high,xlabel, "events", cut&&of,mcordata, luminosity);
110 >    
111 >    float OFmax=ofdata->GetMaximum();
112 >    if(OFmax>sfdata->GetMaximum()) sfdata->SetMaximum(1.2*OFmax);
113 >    ofdata->SetLineColor(kRed);
114 >    sfdata->Draw("e1");
115 >    ofdata->Draw("histo,same");
116 >    sfdata->Draw("e1,same");
117 >    
118 >    kinleg->AddEntry(sfdata,"Same flavor","p");
119 >    kinleg->AddEntry(ofdata,"Opposite flavor","l");
120 >    kinleg->Draw();
121 >    if(mcordata==data) DrawPrelim();
122 >    else DrawMCPrelim();
123 >    
124 >    save_with_ratio(sfdata,ofdata,pad->cd(),"WZ_study/SF_vs_OF_plots/"+saveas);
125 >    delete sfdata;
126 >    delete ofdata;
127 >    delete can;
128 > }
129 >
130 > void write_report(TCut BasicCut, TCut cut, ofstream &report) {
131 >  THStack stack = allsamples.DrawStack("stack","tri_mll",10,0,1000,"m_{ll}", "events", BasicCut&&cut,mc, luminosity);
132 >  float total=0;
133 >  float wz=0;
134 >  TIter nextSF = TIter(stack.GetHists());
135 >  TH1F *h;
136 >  while ( h = (TH1F*)nextSF() ) {
137 >    cout << IdentifierForUnderlinedEntry(h->GetName()) << " : " << h->Integral() << endl;
138 >    if(IdentifierForUnderlinedEntry(h->GetName())=="WZ") {
139 >      wz+=h->Integral();
140 >      total+=h->Integral();
141 >    } else {
142 >      if(IdentifierForUnderlinedEntry(h->GetName())!="Data") total+=h->Integral();
143 >    }
144 >    report << h->Integral() << ";";
145 >  }
146 >  report << wz/total << ";\n";
147 >  report << flush;
148 >
149 > }
150 >
151 > vector<string> producecuts(string variablename, float cut1, float cut2, float cut3=-1, float cut4=-1, float cut5=-1, float cut6=-1, float cut7=-1, float cut8=-1, float cut9=-1, float cut10=-1, float cut11=-1 ) {
152 >  float allcuts[12] = {cut1,cut2,cut3,cut4,cut5,cut6,cut7,cut8,cut9,cut10,cut11};
153 >  vector<string> cuts;
154 >  for(int i = 0; i < 11; ++i ) {
155 >    if(allcuts[i]<0) continue;
156 >    stringstream nowcut;
157 >    nowcut << variablename << allcuts[i];
158 >    cuts.push_back(nowcut.str());
159 >  }
160 >  return cuts;
161 > }
162 >    
163 >
164 > void OptimizeSelection(TCut BasicCut) {
165 >  ofstream report;
166 >  report.open ("report3.txt");
167 >  
168 >  vector<string> mt   = producecuts("tri_mT>",0,10,15,20,25,30,35,40,45,50);
169 >  vector<string> et   = producecuts("met[4]>",0,10,20,30,40,50);
170 >  vector<string> mll  = producecuts("abs(tri_mll-91.2)<",10000,5,10,15,20);
171 >  vector<string> mlll = producecuts("abs(tri_mlll-91)>",0,5,10,15,20,25,30,40,50);
172 >  
173 >  for(int imt=0;imt<mt.size();imt++) {
174 >    for(int iet=0;iet<et.size();iet++) {
175 >      for(int imll=0;imll<mll.size();imll++) {
176 >        for(int imlll=0;imlll<mlll.size();imlll++) {
177 >          stringstream cut;
178 >          cut << mt[imt] << " && " << et[iet] << " && " << mll[imll] << " && " << mlll[imlll];
179 >          report << mt[imt] << ";" << et[iet] << ";" << mll[imll] << ";" << mlll[imlll] <<";";
180 > //        cout << "Cut is : " << cut.str() << endl;
181 >          write_report(BasicCut, TCut(cut.str().c_str()),report);
182 >        }
183 >      }
184 >    }
185 >  }
186 >  report.close();
187 > }
188 >
189 > void SF_OF_comparison(TCut BasicCut, TCut Selection) {
190 >  int nbins=38;
191 >  float mll_low=10;
192 >  float mll_hi=200;
193 >  
194 >  // all MC samples
195 >  TH1F *sfmc   = allsamples.Draw("sfmc"  ,"tri_submll",nbins,mll_low,mll_hi,"m_{l_{2}l_{3}} [GeV]", "events", BasicCut&&Selection&&TriLepSF,mc  , luminosity);
196 >  TH1F *ofmc   = allsamples.Draw("ofmc"  ,"tri_submll",nbins,mll_low,mll_hi,"m_{l_{2}l_{3}} [GeV]", "events", BasicCut&&Selection&&TriLepOF,mc  , luminosity);
197    
198 +  TCanvas *can = new TCanvas("can","can");
199 +  TPad *pad = new TPad("pad","pad",0,0,1,1);
200 +  pad->cd();
201 +  if(ofmc->GetMaximum()>sfmc->GetMaximum()) sfmc->SetMaximum(1.3*ofmc->GetMaximum());
202 +  ofmc->SetLineColor(kRed);
203 +  TLegend *leg = make_legend();
204 +  leg->AddEntry(sfmc,"same flavor","P");
205 +  leg->AddEntry(ofmc,"opposite flavor","L");
206 +  leg->SetHeader("MC:");
207 +  sfmc->Draw("e1");
208 +  ofmc->Draw("histo,same");
209 +  sfmc->Draw("e1,same");
210 +  leg->Draw();
211 +  DrawMCPrelim();
212 +
213 +  save_with_ratio(sfmc,ofmc,pad->cd(),"WZ_study/OFSFcomparison_MC");
214    
215 <  string WZPt2="((leptonId[0]==leptonId[1])*(leptonCharge[0]*leptonCharge[1]<0)*leptonPt[1]   +   ((leptonId[0]!=leptonId[1])||(leptonCharge[0]*leptonCharge[1]>0))*(leptonId[0]==leptonId[2])*(leptonCharge[0]*leptonCharge[2]<0)*leptonPt[2]   +    ((((leptonId[0]!=leptonId[1])||(leptonCharge[0]*leptonCharge[1]>0)))*(((leptonId[0]!=leptonId[2])||(leptonCharge[0]*leptonCharge[2]>0))))*(leptonId[1]==leptonId[2])*(leptonCharge[1]*leptonCharge[2]<0)*leptonPt[2])";
216 <  //this is the pt of the second lepton of the Z
215 >  // WZ only
216 >  TH1F *sfwz   = allsamples.Draw("sfwz"  ,"tri_submll",nbins,mll_low,mll_hi,"m_{l_{2}l_{3}} [GeV]", "events", BasicCut&&Selection&&TriLepSF,mc  , luminosity,allsamples.FindSample("WZJetsTo3"));
217 >  TH1F *ofwz   = allsamples.Draw("ofwz"  ,"tri_submll",nbins,mll_low,mll_hi,"m_{l_{2}l_{3}} [GeV]", "events", BasicCut&&Selection&&TriLepOF,mc  , luminosity,allsamples.FindSample("WZJetsTo3"));
218    
219 <  string WZPt3="((leptonId[0]==leptonId[1])*(leptonCharge[0]*leptonCharge[1]<0)*leptonPt[2]   +   ((leptonId[0]!=leptonId[1])||(leptonCharge[0]*leptonCharge[1]>0))*(leptonId[0]==leptonId[2])*(leptonCharge[0]*leptonCharge[2]<0)*leptonPt[1]   +    ((((leptonId[0]!=leptonId[1])||(leptonCharge[0]*leptonCharge[1]>0)))*(((leptonId[0]!=leptonId[2])||(leptonCharge[0]*leptonCharge[2]>0))))*(leptonId[1]==leptonId[2])*(leptonCharge[1]*leptonCharge[2]<0)*leptonPt[0])";
220 <  //this is the pt of the lepton pertaining to the W
219 >  TPad *pad2 = new TPad("pad","pad",0,0,1,1);
220 >  pad2->cd();
221 >  if(ofwz->GetMaximum()>sfwz->GetMaximum()) sfwz->SetMaximum(1.3*ofwz->GetMaximum());
222 >  ofwz->SetLineColor(kRed);
223 >  sfwz->Draw("e1");
224 >  ofwz->Draw("histo,same");
225 >  sfwz->Draw("e1,same");
226 >  leg->SetHeader("WZ MC:");
227 >  leg->Draw();
228 >  DrawMCPrelim();
229 >  
230 >  save_with_ratio(sfwz,ofwz,pad2->cd(),"WZ_study/OFSFcomparison_WZ");
231 >  
232 >  // data
233 >  TH1F *sfdata = allsamples.Draw("sfdata","tri_submll",nbins,mll_low,mll_hi,"m_{l_{2}l_{3}} [GeV]", "events", BasicCut&&Selection&&TriLepSF,data, luminosity);
234 >  TH1F *ofdata = allsamples.Draw("ofdata","tri_submll",nbins,mll_low,mll_hi,"m_{l_{2}l_{3}} [GeV]", "events", BasicCut&&Selection&&TriLepOF,data, luminosity);
235 >  
236 >  TPad *pad3 = new TPad("pad","pad",0,0,1,1);
237 >  pad3->cd();
238 >  if(ofdata->GetMaximum()>sfdata->GetMaximum()) sfdata->SetMaximum(1.3*ofdata->GetMaximum());
239 >  ofdata->SetLineColor(kRed);
240 >  sfdata->Draw("e1");
241 >  ofdata->Draw("histo,same");
242 >  sfdata->Draw("e1,same");
243 >  leg->SetHeader("Data:");
244 >  leg->Draw();
245 >  DrawPrelim();
246    
247 <  TCut BasicCut("leptonPt[1]>20 && leptonPt[2]>20 && leptonPt[0]>30 && leptonNum>=3 &&   (((leptonId[0]==leptonId[1])||(leptonCharge[0]*leptonCharge[1]<0))||((leptonId[0]==leptonId[2])*(leptonCharge[0]*leptonCharge[2]<0)))   ");
247 >  save_with_ratio(sfdata,ofdata,pad3->cd(),"WZ_study/OFSFcomparison_data");
248 >
249 >  delete can;
250    
251 <  TCut ZmassWindow("abs((((leptonId[0]==leptonId[1])*(leptonCharge[0]*leptonCharge[1]<0)*leptonPairMass[0]   +   ((leptonId[0]!=leptonId[1])||(leptonCharge[0]*leptonCharge[1]>0))*(leptonId[0]==leptonId[2])*(leptonCharge[0]*leptonCharge[2]<0)*leptonPairMass[1]))-91)<10");
251 >  dout << "SF events: " << sfdata->Integral() << "   (mc : " << sfmc->Integral() << " , wz : " << sfwz->Integral() << " )   WZ purity: " << 100*sfwz->Integral()/sfmc->Integral() << " %" << endl;
252 >  dout << "OF events: " << ofdata->Integral() << "   (mc : " << ofmc->Integral() << " , wz : " << ofwz->Integral() << " )   WZ purity: " << 100*ofwz->Integral()/ofmc->Integral() << " %" << endl;
253    
254 +  dout << "RATIO : SF / OF = " << sfdata->Integral() / ofdata->Integral() << " +/- " << (sfdata->Integral()/ofdata->Integral()) * TMath::Sqrt(1/sfdata->Integral()+1/ofdata->Integral()) << " (stat) " << endl;
255 +    
256 +  TriLeptons::Rsfof_CorrectionFactor=ofdata->Integral() / sfdata->Integral();
257 +  TriLeptons::Rsfof_CorrectionFactor_Error=(ofdata->Integral()/sfdata->Integral()) * TMath::Sqrt(1/sfdata->Integral()+1/ofdata->Integral());
258    
259 + }
260 +
261 + void DrawMatchGenDistribution(string variable, int nbins, float low, float high, bool uselog, string xlabel, string SaveAs, TCut basecut) {
262 +  TLegend *leg = make_legend();
263 +  leg->SetY1(0.7);
264    TCanvas *can = new TCanvas("can","can");
265    TPad *pad = new TPad("pad","pad",0,0,1,1);
266    pad->cd();
267 +  pad->SetLogy(uselog);
268 +  
269 +  TCut IsMatch("(gentri_GoodWMatch&&gentri_GoodZMatch)");
270 +  TCut NoMatch("!(gentri_GoodWMatch&&gentri_GoodZMatch)");
271 +  
272 +  TH1F *goodMatch = allsamples.Draw("goodMatch",variable,nbins,low,high,xlabel,"events",basecut&&IsMatch,mc,luminosity,allsamples.FindSample("WZJetsTo3"));
273 +  TH1F *noMatch =   allsamples.Draw("noMatch",  variable,nbins,low,high,xlabel,"events",basecut&&NoMatch,mc,luminosity,allsamples.FindSample("WZJetsTo3"));
274 +  
275 +  noMatch->SetLineColor(TColor::GetColor("#FF0000"));//an evil red
276 +  goodMatch->SetLineColor(TColor::GetColor("#04B404")); // a nice green
277 +  
278 + //   noMatch->Scale(1.0/noMatch->Integral());
279 + //   goodMatch->Scale(1.0/goodMatch->Integral());
280 +  
281 +  if(noMatch->GetMaximum()>goodMatch->GetMaximum()) goodMatch->SetMaximum(1.3*noMatch->GetMaximum());
282 +  
283 +  goodMatch->Draw("histo");
284 +  noMatch->Draw("histo,same");
285 +  
286 +  leg->AddEntry(goodMatch,"correct match","l");
287 +  leg->AddEntry(noMatch,"incorrect match","l");
288 +  leg->Draw();
289 +  
290 +  save_with_ratio(goodMatch,noMatch,pad->cd(),"WZ_study/GenMatchComparison/"+SaveAs);
291 +  
292 +  delete goodMatch;
293 +  delete noMatch;
294 +  delete can;
295 + }
296  
297 <  THStack m_zmass = allsamples.DrawStack("m_zmass",Zmass,40,0,200, "m_{Z} [GeV]", "events", BasicCut,mc, luminosity);
298 <  TH1F *d_zmass = allsamples.Draw("d_zmass",Zmass,40,0,200, "m_{Z} [GeV]", "events", BasicCut,data, luminosity);
297 > void DrawMatchRecoDistribution(string variable, int nbins, float low, float high, bool uselog, string xlabel, string SaveAs, TCut basecut) {
298 >  TLegend *leg = make_legend();
299 >  leg->SetY1(0.7);
300 >  TCanvas *can = new TCanvas("can","can");
301 >  TPad *pad = new TPad("pad","pad",0,0,1,1);
302 >  pad->cd();
303 >  pad->SetLogy(uselog);
304    
305 <  d_zmass->Draw("e1");
306 <  m_zmass.Draw("same");
307 <  d_zmass->Draw("e1,same");
308 <  TLegend *kinleg = allsamples.allbglegend();
309 <  kinleg->Draw();
305 >  TCut IsMatch("(tri_GoodWMatch&&tri_GoodZMatch)");
306 >  TCut NoMatch("!(tri_GoodWMatch&&tri_GoodZMatch)");
307 >  
308 >  TH1F *goodMatch = allsamples.Draw("goodMatch",variable,nbins,low,high,xlabel,"events",basecut&&IsMatch,mc,luminosity,allsamples.FindSample("WZJetsTo3"));
309 >  TH1F *noMatch =   allsamples.Draw("noMatch",  variable,nbins,low,high,xlabel,"events",basecut&&NoMatch,mc,luminosity,allsamples.FindSample("WZJetsTo3"));
310 >  
311 >  noMatch->SetLineColor(TColor::GetColor("#FF0000"));//an evil red
312 >  goodMatch->SetLineColor(TColor::GetColor("#04B404")); // a nice green
313 >  
314 > //   noMatch->Scale(1.0/noMatch->Integral());
315 > //   goodMatch->Scale(1.0/goodMatch->Integral());
316 >  
317 >  if(noMatch->GetMaximum()>goodMatch->GetMaximum()) goodMatch->SetMaximum(1.3*noMatch->GetMaximum());
318 >  
319 >  goodMatch->Draw("histo");
320 >  noMatch->Draw("histo,same");
321 >  
322 >  leg->AddEntry(goodMatch,"correct match","l");
323 >  leg->AddEntry(noMatch,"incorrect match","l");
324 >  leg->Draw();
325    DrawPrelim();
58  save_with_ratio(d_zmass,m_zmass,pad->cd(),"WZ_study/Zmass");
326    
327 <  cout << "Numbers with Z mass:" << endl;
328 <  TIter nextSF = TIter(m_zmass.GetHists());
327 >  save_with_ratio(goodMatch,noMatch,pad->cd(),"WZ_study/RecoMatchComparison/"+SaveAs);
328 >  
329 >  delete goodMatch;
330 >  delete noMatch;
331 >  delete can;
332 > }
333 >
334 > void GenSFOFCalculator(string title, TCut Cut) {
335 >  TH1F *histoSF = allsamples.Draw("histoSF","gentri_mll",1,0,1000,"processing","events",Cut&&TCut("abs(gentri_id2)==abs(gentri_id3)"),mc,luminosity,allsamples.FindSample("WZJetsTo3"));
336 >  TH1F *histoOF = allsamples.Draw("histoOF","gentri_mll",1,0,1000,"processing","events",Cut&&TCut("abs(gentri_id2)!=abs(gentri_id3)"),mc,luminosity,allsamples.FindSample("WZJetsTo3"));
337 >  cout << title << endl;
338 >  cout << "     " << (const char*) Cut << endl;
339 >  cout << "     SF:    \t" << histoSF->Integral() << " +/- " << sqrt(histoSF->GetEntries())*(histoSF->Integral()/histoSF->GetEntries()) << endl;
340 >  cout << "     OF:    \t" << histoOF->Integral() << " +/- " << sqrt(histoOF->Integral())*(histoOF->Integral()/histoOF->GetEntries()) << endl;
341 >  cout << "     SF/OF: \t" << histoSF->Integral() / histoOF->Integral() << " +/- " << histoSF->GetEntries() / histoOF->GetEntries() * sqrt(1.0/histoSF->GetEntries() + 1.0 / histoOF->GetEntries()) << endl;
342 >  
343 >  delete histoSF;
344 >  delete histoOF;
345 > }
346 >  
347 > void reco_SF_OF_comparison(TCut Base, TCut Selection) {
348 >  
349 >  bool dolog=true;
350 >  bool nolog=false;
351 >  
352 >  DrawMatchRecoDistribution("tri_mlll",60,0,300,dolog,"m_{lll}","TrileptonMass",Base);
353 >  DrawMatchRecoDistribution("tri_mll",60,0,300,dolog,"m_{l_{1}l_{2}}","LeadingDileptonMass",Base);
354 >  DrawMatchRecoDistribution("tri_mT",50,0,100,nolog,"m_{T}","MT",Base);
355 >  DrawMatchRecoDistribution("tri_submll",60,0,300,nolog,"m_{l_{2}l_{3}}","SubleadingDileptonMass",Base);
356 >  DrawMatchRecoDistribution("tri_submll",60,0,300,nolog,"m_{l_{2}l_{3}}","SubleadingDileptonMass_SF",Base&&TCut("abs(tri_id2)==abs(tri_id3)"));
357 >  DrawMatchRecoDistribution("tri_submll",60,0,300,nolog,"m_{l_{2}l_{3}}","SubleadingDileptonMass_OF",Base&&TCut("abs(tri_id2)!=abs(tri_id3)"));
358 >  DrawMatchRecoDistribution("genMET",60,0,300,dolog,"gen MET","MET",Base);
359 >  DrawMatchRecoDistribution("abs(tri_id1)",50,-0.5,49.5,nolog,"id(l_{1})","genID1",Base);
360 >  DrawMatchRecoDistribution("abs(tri_id2)",50,-0.5,49.5,nolog,"id(l_{2})","genID2",Base);
361 >  DrawMatchRecoDistribution("abs(tri_id3)",50,-0.5,49.5,nolog,"id(l_{3})","genID3",Base);
362 >  DrawMatchRecoDistribution("abs(tri_genMID1)",50,-0.5,49.5,nolog,"mother id(l_{1})","genMID1",Base);
363 >  DrawMatchRecoDistribution("abs(tri_genMID2)",50,-0.5,49.5,nolog,"mother id(l_{2})","genMID2",Base);
364 >  DrawMatchRecoDistribution("abs(tri_genMID3)",50,-0.5,49.5,nolog,"mother id(l_{3})","genMID3",Base);
365 >  
366 >  DrawMatchRecoDistribution("tri_mlll",60,0,300,dolog,"m_{lll}","FullSelection/TrileptonMass",Base&&Selection);
367 >  DrawMatchRecoDistribution("tri_mll",60,0,300,dolog,"m_{l_{1}l_{2}}","FullSelection/LeadingDileptonMass",Base&&Selection);
368 >  DrawMatchRecoDistribution("tri_mT",50,0,100,nolog,"m_{T}","FullSelection/MT",Base&&Selection);
369 >  DrawMatchRecoDistribution("tri_submll",60,0,300,nolog,"m_{l_{2}l_{3}}","FullSelection/SubleadingDileptonMass",Base&&Selection);
370 >  DrawMatchRecoDistribution("tri_submll",60,0,300,nolog,"m_{l_{2}l_{3}}","FullSelection/SubleadingDileptonMass_SF",Base&&TCut("abs(tri_id2)==abs(tri_id3)")&&Selection);
371 >  DrawMatchRecoDistribution("tri_submll",60,0,300,nolog,"m_{l_{2}l_{3}}","FullSelection/SubleadingDileptonMass_OF",Base&&TCut("abs(tri_id2)!=abs(tri_id3)")&&Selection);
372 >  DrawMatchRecoDistribution("genMET",60,0,300,dolog,"gen MET","FullSelection/MET",Base&&Selection);
373 >  DrawMatchRecoDistribution("abs(tri_id1)",50,-0.5,49.5,nolog,"id(l_{1})","FullSelection/genID1",Base&&Selection);
374 >  DrawMatchRecoDistribution("abs(tri_id2)",50,-0.5,49.5,nolog,"id(l_{2})","FullSelection/genID2",Base&&Selection);
375 >  DrawMatchRecoDistribution("abs(tri_id3)",50,-0.5,49.5,nolog,"id(l_{3})","FullSelection/genID3",Base&&Selection);
376 >  DrawMatchRecoDistribution("abs(tri_genMID1)",50,-0.5,49.5,nolog,"mother id(l_{1})","FullSelection/genMID1",Base&&Selection);
377 >  DrawMatchRecoDistribution("abs(tri_genMID2)",50,-0.5,49.5,nolog,"mother id(l_{2})","FullSelection/genMID2",Base&&Selection);
378 >  DrawMatchRecoDistribution("abs(tri_genMID3)",50,-0.5,49.5,nolog,"mother id(l_{3})","FullSelection/genMID3",Base&&Selection);
379 > }
380 >  
381 > void generator_SF_OF_comparison() {
382 >  /*we want to analyze the cuts at generator level (and their impact), plus understand how bad matching can affect the selection.
383 >    we want to know the distributions for correctly and incorrectly matched W/Z for the following variables:
384 >      - M_{lll}
385 >      - M_{ll}
386 >      - M_{T}
387 >      - M_{l_{2}l_{3}}
388 >      - MET
389 >  */
390 >  TCut genBase("genNleptons>=3 && gentri_pt1>30 && gentri_pt2>20 && gentri_pt3>20");
391 >  TCut genParticleBase("genNleptons>=3 && pgentri_pt1>30 && pgentri_pt2>20 && pgentri_pt3>20");
392 >  
393 >  TCut Selection("abs(gentri_mlll-91)>30 && genMET>30 && abs(gentri_mll-91)<10 && gentri_mT > 40");
394 >  TCut genParticleSelection("abs(pgentri_mlll-91)>30 && genMET>30 && abs(pgentri_mll-91)<10 && abs(pgentri_submll-91)>20 && pgentri_mT > 40");
395 >  
396 >  bool dolog=true;
397 >  bool nolog=false;
398 >  
399 >  DrawMatchGenDistribution("gentri_mlll",60,0,300,dolog,"m_{lll}","TrileptonMass",genBase);
400 >  DrawMatchGenDistribution("gentri_mll",60,0,300,dolog,"m_{l_{1}l_{2}}","LeadingDileptonMass",genBase);
401 >  DrawMatchGenDistribution("gentri_mT",50,0,100,nolog,"m_{T}","MT",genBase);
402 >  DrawMatchGenDistribution("gentri_submll",60,0,300,nolog,"m_{l_{2}l_{3}}","SubleadingDileptonMass",genBase);
403 >  DrawMatchGenDistribution("gentri_submll",60,0,300,nolog,"m_{l_{2}l_{3}}","SubleadingDileptonMass_SF",genBase&&TCut("abs(gentri_id2)==abs(gentri_id3)"));
404 >  DrawMatchGenDistribution("gentri_submll",60,0,300,nolog,"m_{l_{2}l_{3}}","SubleadingDileptonMass_OF",genBase&&TCut("abs(gentri_id2)!=abs(gentri_id3)"));
405 >  DrawMatchGenDistribution("genMET",60,0,300,dolog,"gen MET","MET",genBase);
406 >  DrawMatchGenDistribution("abs(gentri_id1)",50,-0.5,49.5,nolog,"id(l_{1})","genID1",genBase);
407 >  DrawMatchGenDistribution("abs(gentri_id2)",50,-0.5,49.5,nolog,"id(l_{2})","genID2",genBase);
408 >  DrawMatchGenDistribution("abs(gentri_id3)",50,-0.5,49.5,nolog,"id(l_{3})","genID3",genBase);
409 >  DrawMatchGenDistribution("abs(gentri_genMID1)",50,-0.5,49.5,nolog,"mother id(l_{1})","genMID1",genBase);
410 >  DrawMatchGenDistribution("abs(gentri_genMID2)",50,-0.5,49.5,nolog,"mother id(l_{2})","genMID2",genBase);
411 >  DrawMatchGenDistribution("abs(gentri_genMID3)",50,-0.5,49.5,nolog,"mother id(l_{3})","genMID3",genBase);
412 >  
413 >  DrawMatchGenDistribution("gentri_mlll",60,0,300,dolog,"m_{lll}","FullSelection/TrileptonMass",genBase&&Selection);
414 >  DrawMatchGenDistribution("gentri_mll",60,0,300,dolog,"m_{l_{1}l_{2}}","FullSelection/LeadingDileptonMass",genBase&&Selection);
415 >  DrawMatchGenDistribution("gentri_mT",50,0,100,nolog,"m_{T}","FullSelection/MT",genBase&&Selection);
416 >  DrawMatchGenDistribution("gentri_submll",60,0,300,nolog,"m_{l_{2}l_{3}}","FullSelection/SubleadingDileptonMass",genBase&&Selection);
417 >  DrawMatchGenDistribution("gentri_submll",60,0,300,nolog,"m_{l_{2}l_{3}}","FullSelection/SubleadingDileptonMass_SF",genBase&&TCut("abs(gentri_id2)==abs(gentri_id3)")&&Selection);
418 >  DrawMatchGenDistribution("gentri_submll",60,0,300,nolog,"m_{l_{2}l_{3}}","FullSelection/SubleadingDileptonMass_OF",genBase&&TCut("abs(gentri_id2)!=abs(gentri_id3)")&&Selection);
419 >  DrawMatchGenDistribution("genMET",60,0,300,dolog,"gen MET","FullSelection/MET",genBase&&Selection);
420 >  DrawMatchGenDistribution("abs(gentri_id1)",50,-0.5,49.5,nolog,"id(l_{1})","FullSelection/genID1",genBase&&Selection);
421 >  DrawMatchGenDistribution("abs(gentri_id2)",50,-0.5,49.5,nolog,"id(l_{2})","FullSelection/genID2",genBase&&Selection);
422 >  DrawMatchGenDistribution("abs(gentri_id3)",50,-0.5,49.5,nolog,"id(l_{3})","FullSelection/genID3",genBase&&Selection);
423 >  DrawMatchGenDistribution("abs(gentri_genMID1)",50,-0.5,49.5,nolog,"mother id(l_{1})","FullSelection/genMID1",genBase&&Selection);
424 >  DrawMatchGenDistribution("abs(gentri_genMID2)",50,-0.5,49.5,nolog,"mother id(l_{2})","FullSelection/genMID2",genBase&&Selection);
425 >  DrawMatchGenDistribution("abs(gentri_genMID3)",50,-0.5,49.5,nolog,"mother id(l_{3})","FullSelection/genMID3",genBase&&Selection);
426 >  
427 > return;
428 >
429 >
430 >  GenSFOFCalculator("Basic selection",genBase);
431 >  GenSFOFCalculator("Basic selection, correct Z",genBase&&TCut("gentri_GoodZMatch"));
432 >  GenSFOFCalculator("Basic selection, wrong Z",genBase&&TCut("!gentri_GoodZMatch"));
433 >  GenSFOFCalculator("Outside |m_{23}-Z|<20",genBase&&TCut("abs(gentri_submll-91)>20"));
434 >  GenSFOFCalculator("Outside |m_{23}-Z|<20, Correct Z", genBase&&TCut("abs(gentri_submll-91)>20&&gentri_GoodZMatch"));
435 >  GenSFOFCalculator("Outside |m_{23}-Z|<20, Wrong Z", genBase&&TCut("abs(gentri_submll-91)>20&&!gentri_GoodZMatch"));
436 >  GenSFOFCalculator("Outside |m_{23}-Z|<20, Correct W", genBase&&TCut("abs(gentri_submll-91)>20&&gentri_GoodWMatch"));
437 >  GenSFOFCalculator("Outside |m_{23}-Z|<20, Wrong W", genBase&&TCut("abs(gentri_submll-91)>20&&!gentri_GoodWMatch"));
438 >  GenSFOFCalculator("Outside |m_{23}-Z|<20, Correct Z, correct W", genBase&&TCut("abs(gentri_submll-91)>20&&gentri_GoodZMatch&&gentri_GoodWMatch"));
439 >  GenSFOFCalculator("Outside |m_{23}-Z|<20, Wrong Z or W", genBase&&TCut("abs(gentri_submll-91)>20&&!(gentri_GoodZMatch&&gentri_GoodWMatch)"));
440 >  GenSFOFCalculator("Full selection, Correct Z and W", genBase&&Selection&&TCut("abs(gentri_submll-91)>20&&(gentri_GoodZMatch&&gentri_GoodWMatch)"));
441 >  GenSFOFCalculator("Full selection, Wrong Z or W", genBase&&Selection&&TCut("abs(gentri_submll-91)>20&&!(gentri_GoodZMatch&&gentri_GoodWMatch)"));
442 > //  GenSFOFCalculator("Full selection", genBase&&Selection);
443 > //  GenSFOFCalculator("Full selection with genParticles", genBase&&genParticleSelection&&TCut("abs(gentri_submll-91)>20"));
444 >  
445 >  
446 > }
447 >
448 >
449 > void AlgorithmTest(TCut BasicCut, TCut Selection) {
450 >  TCut Pure("tri_GoodWMatch && tri_GoodZMatch");
451 >  
452 >  int nbins=38;
453 >  float mll_low=10;
454 >  float mll_hi=200;
455 >  
456 >  TCanvas *can = new TCanvas("can","can");
457 >  TH1F *sfwz   = allsamples.Draw("sfwz"  ,"tri_submll",nbins,mll_low,mll_hi,"m_{l_{2}l_{3}} [GeV]", "events", BasicCut&&Selection&&TriLepSF,mc,luminosity,allsamples.FindSample("WZJetsTo3"));
458 >  TH1F *ofwz   = allsamples.Draw("ofwz"  ,"tri_submll",nbins,mll_low,mll_hi,"m_{l_{2}l_{3}} [GeV]", "events", BasicCut&&Selection&&TriLepOF,mc,luminosity,allsamples.FindSample("WZJetsTo3"));
459 >  
460 >  TH1F *sfwzpure = allsamples.Draw("sfwzpure","tri_submll",nbins,mll_low,mll_hi,"m_{l_{2}l_{3}} [GeV]", "events", BasicCut&&Selection&&TriLepSF&&Pure,mc,luminosity,allsamples.FindSample("WZJetsTo3"));
461 >  TH1F *ofwzpure = allsamples.Draw("ofwzpure","tri_submll",nbins,mll_low,mll_hi,"m_{l_{2}l_{3}} [GeV]", "events", BasicCut&&Selection&&TriLepOF&&Pure,mc,luminosity,allsamples.FindSample("WZJetsTo3"));
462 >  
463 >  delete can;
464 >  
465 >  dout << "SF events: " << sfwz->Integral() << "   (pure: " << sfwzpure->Integral() << ")" << endl;
466 >  dout << "OF events: " << ofwz->Integral() << "   (pure: " << ofwzpure->Integral() << ")" << endl;
467 >  
468 >  dout << "RATIO : SF / OF = " << sfwz->Integral() / ofwz->Integral() << " +/- " << sfwz->Integral()/ofwz->Integral() * sqrt(1/sfwz->GetEntries() + 1/ofwz->GetEntries()) << endl;
469 >  dout << "  PURE: SF / OF = " << sfwzpure->Integral() / ofwzpure->Integral() << " +/- " << sfwzpure->Integral()/ofwzpure->Integral() * sqrt(1/sfwzpure->GetEntries() + 1/ofwzpure->GetEntries()) << endl;
470 >  dout << "PURITY: " << endl;
471 >  dout << "   SF : " << 100.0*sfwzpure->Integral()/sfwz->Integral() << " +/- " << 100.0*sfwzpure->Integral()/sfwz->Integral() * sqrt(1/sfwzpure->GetEntries() + 1/sfwz->GetEntries()) <<endl;
472 >  dout << "   OF : " << 100.0*ofwzpure->Integral()/ofwz->Integral() << " +/- " << 100.0*ofwzpure->Integral()/ofwz->Integral() * sqrt(1/ofwzpure->GetEntries() + 1/ofwz->GetEntries()) <<endl;
473 >  
474 >
475 >  delete sfwz;
476 >  delete ofwz;
477 >  delete sfwzpure;
478 >  delete ofwzpure;
479 > }
480 >
481 > void GetResult(TCut BasicCut, TCut Selection) {
482 >    if(TriLeptons::Rsfof_CorrectionFactor<0 || TriLeptons::Rsfof_CorrectionFactor_Error<0) SF_OF_comparison(BasicCut,Selection);
483 >    int nbins=38;
484 >    float mll_low=10;
485 >    float mll_hi=200;
486 >    
487 >    TCanvas *can = new TCanvas("can","can");
488 >    TH1F *sfwz   = allsamples.Draw("sfwz"  ,"tri_submll",nbins,mll_low,mll_hi,"m_{l_{2}l_{3}} [GeV]", "events", BasicCut&&Selection&&TriLepSF,data,luminosity);
489 >    TH1F *ofwz   = allsamples.Draw("ofwz"  ,"tri_submll",nbins,mll_low,mll_hi,"m_{l_{2}l_{3}} [GeV]", "events", BasicCut&&Selection&&TriLepOF,data,luminosity);
490 >    
491 >    delete can;
492 >    
493 >    dout << " ****** FINAL RESULT ****** " << endl;
494 >    dout << "SF events: " << sfwz->Integral() << " +/- " << sqrt(sfwz->Integral()) << endl;
495 >    dout << "OF events: " << ofwz->Integral() << " +/- " << sqrt(ofwz->Integral()) << endl;
496 >    
497 >    float rawstatuncert=sfwz->Integral()/ofwz->Integral() * sqrt(1/sfwz->GetEntries() + 1/ofwz->GetEntries());
498 >    float rawresult=sfwz->Integral() / ofwz->Integral();
499 >    dout << "UNCORRECTED RATIO : SF / OF = " << rawresult << " +/- " << rawstatuncert << endl;
500 >    dout << "CORRECTED RATIO : SF / OF = " << TriLeptons::Rsfof_CorrectionFactor * rawresult << " +/- " << (TriLeptons::Rsfof_CorrectionFactor * rawresult) * sqrt(rawstatuncert*rawstatuncert/(rawresult*rawresult) + ( TriLeptons::Rsfof_CorrectionFactor_Error * TriLeptons::Rsfof_CorrectionFactor_Error )/(TriLeptons::Rsfof_CorrectionFactor * TriLeptons::Rsfof_CorrectionFactor));
501 >    
502 >    
503 >    delete sfwz;
504 >    delete ofwz;
505 > //    delete can;
506 > }
507 >
508 >
509 > void PredictionFailurePlot(string truthvariable, string predictionvariable, int nbins, float min, float max, TCut PredictionCut, TCut TruthCut, string saveas, string xlabel) {
510 >  TLegend *leg = make_legend();
511 >  leg->SetY1(0.7);
512 >  TCanvas *can = new TCanvas("can","can");
513 >  TPad *pad = new TPad("pad","pad",0,0,1,1);
514 >  pad->cd();
515 > //    pad->SetLogy(1);
516 >  
517 >  TH1F *Prediction = allsamples.Draw("Prediction",predictionvariable,nbins,min,max,xlabel,"events",PredictionCut,mc,luminosity,allsamples.FindSample("WZJetsTo3"));
518 >  TH1F *MCTruth    = allsamples.Draw("MCTruth",   truthvariable,     nbins,min,max,xlabel,"events",TruthCut     ,mc,luminosity,allsamples.FindSample("WZJetsTo3"));
519 >  
520 >  Prediction->SetLineColor(TColor::GetColor("#FF0000"));//an evil red
521 >  MCTruth->SetLineColor(TColor::GetColor("#04B404")); // a nice green
522 >  
523 > //  if(Prediction->GetMaximum()>MCTruth->GetMaximum()) MCTruth->SetMaximum(1.3*Prediction->GetMaximum());
524 >  
525 >  if(MCTruth->GetMaximum()>Prediction->GetMaximum()) {
526 >    MCTruth->SetMaximum(MCTruth->GetMaximum()*1.2);
527 >    MCTruth->DrawNormalized("histo");
528 >    Prediction->DrawNormalized("histo,same");
529 >  } else {
530 >    Prediction->SetMaximum(Prediction->GetMaximum()*1.2);
531 >    Prediction->DrawNormalized("histo");
532 >    MCTruth->DrawNormalized("histo,same");
533 >  }
534 >    
535 >  
536 >  leg->AddEntry(Prediction,"prediction (norm.)","l");
537 >  leg->AddEntry(MCTruth,"MC truth (norm.)","l");
538 >  leg->Draw();
539 >  DrawMCPrelim();
540 >  save_with_ratio(MCTruth,Prediction,pad->cd(),"WZ_study/PredictionMethodFailure/"+saveas);
541 >  
542 >  delete Prediction;
543 >  delete MCTruth;
544 >  delete can;
545 >  
546 >  
547 >  
548 > }
549 >
550 > void WhyDoesThePredictionMethodFail(TCut PredictionCut, TCut TruthCut) {
551 >  PredictionFailurePlot("tri_submll","tri_badsubmll",50,0,200,PredictionCut,TruthCut,"SubMll","m_{23}");
552 >  PredictionFailurePlot("tri_mT","tri_badmT",50,0,200,PredictionCut,TruthCut,"MT","m_{T}");
553 >  PredictionFailurePlot("tri_mll","tri_badmll",25,0,200,PredictionCut,TruthCut,"Mll","m_{12}");
554 > }
555 >  
556 >  
557 > void ComputeNewMiscombinationRate(TCut origSelection, TCut sBasicCut) {
558 >  TCanvas *can = new TCanvas("can","can");
559 >  
560 >  string sSelection= (const char*)origSelection;
561 >  sSelection= ReplaceAll(sSelection,"tri_","tri_bad");
562 >  sSelection= ReplaceAll(sSelection,"tri_badmlll","tri_mlll");
563 >  cout << "modified " << (const char*)origSelection << " to " << sSelection << endl;
564 >  
565 >  TCut FlavorBlindSelection((sSelection+(string)"&&(abs(tri_id2)!=abs(tri_id3))").c_str()); // making sure we select OF events
566 >  TCut CorrectSelection    (origSelection && TCut("(abs(tri_id2)!=abs(tri_id3))"));        // making sure we select OF events
567 >  TCut OFZ                 ("tri_mll!=tri_badmll");
568 >  
569 >  TH1F *WZ_OF_OFZ = allsamples.Draw("WZ_OFZ"    ,"tri_badsubmll",1,0.0,1000.0,"m_{l_{2}l_{3}} [GeV]", "events", FlavorBlindSelection&&OFZ ,mc,luminosity,allsamples.FindSample("WZJetsTo3"));
570 >  TH1F *WZ_OF_any = allsamples.Draw("WZ_OF_any" ,"tri_badsubmll",1,0.0,1000.0,"m_{l_{2}l_{3}} [GeV]", "events", FlavorBlindSelection      ,mc,luminosity,allsamples.FindSample("WZJetsTo3"));
571 >  
572 >  dout << "Miscombination rate in WZ sample with Selection in place : " << 100 * WZ_OF_OFZ->Integral() / WZ_OF_any->Integral() << " %    (" << WZ_OF_OFZ->Integral() << " / " << WZ_OF_any->Integral() << ")" << endl;
573 > //   dout << "This corresponds to an estimate in same flavor of        : " <<
574 >  
575 >  delete WZ_OF_any;
576 >  delete WZ_OF_OFZ;
577 >  delete can;
578 > }
579 >
580 > void ComputeMiscombinationRate(TCut origSelection, TCut sBasicCut) {
581 > /*  TCanvas *can = new TCanvas("can","can");
582 >  //idea: use the flavor-blind version of the algorithm to obtain the yields in our signal region.
583 >  string sSelection= (const char*)origSelection;
584 >  sSelection= ReplaceAll(sSelection,"tri_","tri_bad");
585 >  sSelection= ReplaceAll(sSelection,"tri_badmlll","tri_mlll");
586 >  cout << "modified " << (const char*)origSelection << " to " << sSelection << endl;
587 >  
588 >  TCut BadSelection((sSelection+(string)"&&abs(tri_id2)!=abs(tri_id3)").c_str()); // making sure we select OF events
589 >  TCut GoodSelection(origSelection && TCut("abs(tri_id2)!=abs(tri_id3)")); // making sure we select OF events
590 >  
591 >  TCut BasicCut(sBasicCut&&TCut("tri_id2!=tri_id3"));
592 >  
593 >  TCut TruthBad("!(tri_GoodZMatch&&tri_GoodZMatch)");
594 >  TCut MethodBad("!(tri_id2==tri_badid2&&tri_id3==tri_badid3)");
595 >  
596 >  for(int icut=0;icut<2;icut++) {
597 >    TCut Goodused;
598 >    TCut Badused;
599 >    if(icut==0) {
600 >      Goodused=BasicCut;
601 >      Badused =BasicCut;
602 >      cout << "Preselection only: " << endl;
603 >      cout << "      used : " << (const char*) Goodused << endl;
604 >    } else {
605 >      Goodused=BasicCut&&GoodSelection;
606 >      Badused =BasicCut&&BadSelection;
607 >      cout << "Full selection: " << endl;
608 >      cout << "      Good used : " << (const char*) Goodused << endl;
609 >      cout << "      Bad used : " << (const char*) Goodused << endl;
610 >    }
611 >    
612 >    TH1F *Truth_badcomb = allsamples.Draw("Truth_badcomb"    ,"tri_submll",1,0.0,1000.0,"m_{l_{2}l_{3}} [GeV]", "events", Goodused           ,mc,luminosity,allsamples.FindSample("WZJetsTo3"));
613 >    TH1F *Truth_all     = allsamples.Draw("Truth_all"        ,"tri_submll",1,0.0,1000.0,"m_{l_{2}l_{3}} [GeV]", "events", Goodused           ,mc,luminosity,allsamples.FindSample("WZJetsTo3"));
614 >    
615 >    TH1F *Events_WZ    = allsamples.Draw("Events_WZ"    ,"tri_submll",1,0.0,1000.0,"m_{l_{2}l_{3}} [GeV]", "events", Goodused           ,mc,luminosity,allsamples.FindSample("WZJetsTo3"));
616 >    TH1F *bad_WZ       = allsamples.Draw("bad_WZ"       ,"tri_submll",1,0.0,1000.0,"m_{l_{2}l_{3}} [GeV]", "events", Goodused&&TruthBad ,mc,luminosity,allsamples.FindSample("WZJetsTo3"));
617 >    TH1F *MethodBad_WZ = allsamples.Draw("MethodBad_WZ" ,"tri_submll",1,0.0,1000.0,"m_{l_{2}l_{3}} [GeV]", "events", Badused&&MethodBad,mc,luminosity,allsamples.FindSample("WZJetsTo3"));
618 >    TH1F *MethodBad_WZg= allsamples.Draw("MethodBad_WZg","tri_submll",1,0.0,1000.0,"m_{l_{2}l_{3}} [GeV]", "events", Badused,mc,luminosity,allsamples.FindSample("WZJetsTo3"));
619 >
620 >    cout << "    ********** WZ **********" << endl;
621 >    cout << "    N events :          " << Events_WZ->Integral() << "\t (mc events: " << Events_WZ->GetEntries() << ")" << endl;
622 >    cout << "    Bad combinations:   " << bad_WZ->Integral() << "\t (mc events: " << bad_WZ->GetEntries() << ")" << endl;
623 >    cout << "    Predicted bad com.: " << MethodBad_WZ->Integral() << " / " << MethodBad_WZg->Integral() << "\t (mc events: " << MethodBad_WZ->GetEntries() << ")" << endl;
624 >
625 >
626 >    delete Events_WZ;
627 >    delete bad_WZ;
628 >    delete MethodBad_WZ;
629 >    delete MethodBad_WZg;
630 >  }
631 >  delete can;
632 >  
633 >  WhyDoesThePredictionMethodFail(BasicCut&&MethodBad,BasicCut&&TruthBad);
634 > */  
635 > }
636 >  
637 >  
638 >
639 >
640 > void CarryOutWZStudy() {
641 >  TCut BasicCut("leptonNum>2&&tri_pt1>30&&tri_pt2>20&&tri_pt3>20");
642 >  TCut Selection("abs(tri_mlll-91)>30 && tri_mT>40 && met[4] > 30 && abs(tri_mll-91)<10 && tri_submll > 10 && abs(tri_submll-91)>5");
643 >
644 >  //OptimizeSelection(BasicCut);
645 >  
646 >  //ComputeMiscombinationRate(Selection, BasicCut);
647 >  
648 >  ComputeNewMiscombinationRate(Selection, BasicCut);
649 >  return;
650 >  
651 >  WZ_kin_plot("tri_mll",40,0,200,"m_{Z} [GeV]",BasicCut,"Inclusive/Zmass");
652 >  WZ_kin_plot("tri_mlll",40,0,400,"m_{ll} [GeV]",BasicCut,"Inclusive/TrileptonMass");
653 >  WZ_kin_plot("tri_mT",40,0,200,"m_{T} [GeV]",BasicCut,"Inclusive/mT");
654 >  WZ_kin_plot("met[4]", 40,0,400,"raw PFMET [GeV]",BasicCut,"Inclusive/MET");
655 >  WZ_kin_plot("tri_pt1",40,0,400,"p_{T}^{1} [GeV]",BasicCut,"Inclusive/Pt1");
656 >  WZ_kin_plot("tri_pt2",40,0,400,"p_{T}^{2} [GeV]",BasicCut,"Inclusive/Pt2");
657 >  WZ_kin_plot("tri_pt3",40,0,400,"p_{T}^{3} [GeV]",BasicCut,"Inclusive/Pt3");
658 >
659 >  WZ_kin_plot("tri_mll",40,0,200,"m_{Z} [GeV]",BasicCut&&TCut("abs(tri_mlll-91)>30"),"Mlll_cut/Zmass");
660 >  WZ_kin_plot("tri_mlll",40,0,400,"m_{ll} [GeV]",BasicCut&&TCut("abs(tri_mlll-91)>30"),"Mlll_cut/TrileptonMass");
661 >  WZ_kin_plot("tri_mT",40,0,200,"m_{T} [GeV]",BasicCut&&TCut("abs(tri_mlll-91)>30"),"Mlll_cut/mT");
662 >  WZ_kin_plot("met[4]", 40,0,400,"raw PFMET [GeV]",BasicCut&&TCut("abs(tri_mlll-91)>30"),"Mlll_cut/MET");
663 >  
664 >  WZ_kin_plot("tri_mll",40,0,200,"m_{Z} [GeV]",BasicCut&&TCut("abs(tri_mll-91)<10"),"Zwindow/Zmass");
665 >  WZ_kin_plot("tri_mlll",40,0,400,"m_{ll} [GeV]",BasicCut&&TCut("abs(tri_mll-91)<10"),"Zwindow/TrileptonMass");
666 >  WZ_kin_plot("tri_mT",40,0,200,"m_{T} [GeV]",BasicCut&&TCut("abs(tri_mll-91)<10"),"Zwindow/mT");
667 >  WZ_kin_plot("met[4]", 40,0,400,"raw PFMET [GeV]",BasicCut&&TCut("abs(tri_mll-91)<10"),"Zwindow/MET");
668 >
669 >  WZ_kin_plot("tri_mll",40,0,200,"m_{Z} [GeV]",BasicCut,"CutFlow/Zmass__0_CutOnNothing");
670 >  WZ_kin_plot("tri_mll",40,0,200,"m_{Z} [GeV]",BasicCut&&TCut("tri_mT>40"),"CutFlow/Zmass__1_CutOn_MT");
671 >  WZ_kin_plot("tri_mll",40,0,200,"m_{Z} [GeV]",BasicCut&&TCut("tri_mT>40&&abs(tri_mlll-91)>30"),"CutFlow/Zmass__2_CutOn_MT_n_Mll");
672 >  WZ_kin_plot("tri_mll",40,0,200,"m_{Z} [GeV]",BasicCut&&TCut("tri_mT>40&&abs(tri_mlll-91)>30&&met[4]>30"),"CutFlow/Zmass__3_CutOn_MT_n_Mll_n_Met");
673 >  WZ_kin_plot("tri_mll",40,0,200,"m_{Z} [GeV]",BasicCut&&TCut("tri_mT>40&&abs(tri_mlll-91)>30&&met[4]>30&&abs(tri_mll-91)<10"),"CutFlow/Zmass__4_CutOn_MT_n_Mll_n_Met_MllWindow");
674 >  
675 >  WZ_kin_plot("tri_mT",40,0,200,"m_{T} [GeV]",BasicCut,"CutFlow/MT__0_CutOnNothing");
676 >  WZ_kin_plot("tri_mT",40,0,200,"m_{T} [GeV]",BasicCut&&TCut("tri_mT>40"),"CutFlow/MT__1_CutOn_MT");
677 >  WZ_kin_plot("tri_mT",40,0,200,"m_{T} [GeV]",BasicCut&&TCut("tri_mT>40&&abs(tri_mlll-91)>30"),"CutFlow/MT__2_CutOn_MT_n_Mll");
678 >  WZ_kin_plot("tri_mT",40,0,200,"m_{T} [GeV]",BasicCut&&TCut("tri_mT>40&&abs(tri_mlll-91)>30&&met[4]>30"),"CutFlow/MT__3_CutOn_MT_n_Mll_n_Met");
679 >  WZ_kin_plot("tri_mT",40,0,200,"m_{T} [GeV]",BasicCut&&TCut("tri_mT>40&&abs(tri_mlll-91)>30&&met[4]>30&&abs(tri_mll-91)<10"),"CutFlow/MT__4_CutOn_MT_n_Mll_n_Met_MllWindow");
680 >  
681 >  WZ_kin_plot("met[4]",40,0,200,"MET [GeV]",BasicCut,"CutFlow/MET__0_CutOnNothing");
682 >  WZ_kin_plot("met[4]",40,0,200,"MET [GeV]",BasicCut&&TCut("tri_mT>40"),"CutFlow/MET__1_CutOn_MT");
683 >  WZ_kin_plot("met[4]",40,0,200,"MET [GeV]",BasicCut&&TCut("tri_mT>40&&abs(tri_mlll-91)>30"),"CutFlow/MET__2_CutOn_MT_n_Mll");
684 >  WZ_kin_plot("met[4]",40,0,200,"MET [GeV]",BasicCut&&TCut("tri_mT>40&&abs(tri_mlll-91)>30&&met[4]>30"),"CutFlow/MET__3_CutOn_MT_n_Mll_n_Met");
685 >  WZ_kin_plot("met[4]",40,0,200,"MET [GeV]",BasicCut&&TCut("tri_mT>40&&abs(tri_mlll-91)>30&&met[4]>30&&abs(tri_mll-91)<10"),"CutFlow/MET__4_CutOn_MT_n_Mll_n_Met_MllWindow");
686 >
687 >  
688 >  WZ_kin_plot("tri_mll",40,0,200,"m_{Z} [GeV]",BasicCut&&Selection,"Zmass");
689 >  WZ_kin_plot("tri_mlll",40,0,400,"m_{lll} [GeV]",BasicCut&&Selection,"TrileptonMass");
690 >  WZ_kin_plot("tri_mT",40,0,200,"m_{T} [GeV]",BasicCut&&Selection,"mT");
691 >  WZ_kin_plot("met[4]", 40,0,400,"raw PFMET [GeV]",BasicCut&&Selection,"MET");
692 >  WZ_kin_plot("tri_pt1",40,0,400,"p_{T}^{1} [GeV]",BasicCut&&Selection,"Pt1");
693 >  WZ_kin_plot("tri_pt2",40,0,400,"p_{T}^{2} [GeV]",BasicCut&&Selection,"Pt2");
694 >  WZ_kin_plot("tri_pt3",40,0,400,"p_{T}^{3} [GeV]",BasicCut&&Selection,"Pt3");
695 >  
696 >
697 >  WZ_kin_plot("tri_submll",40,0,200,"m_{l_{2}l_{3}} [GeV]",BasicCut&&Selection&&TriLepSF,"mllsub_SF");
698 >  WZ_kin_plot("tri_submll",40,0,200,"m_{l_{2}l_{3}} [GeV]",BasicCut&&Selection&&TriLepOF,"mllsub_OF");
699 >  WZ_kin_plot("tri_submll",40,0,200,"m_{l_{2}l_{3}} [GeV]",BasicCut&&Selection&&TriLepSF&&TCut("tri_id2==0"),"mllsub_SF_ee");
700 >  WZ_kin_plot("tri_submll",40,0,200,"m_{l_{2}l_{3}} [GeV]",BasicCut&&Selection&&TriLepSF&&TCut("tri_id2==1"),"mllsub_SF_mm");
701 >  WZ_kin_plot("tri_submll",40,0,200,"m_{l_{2}l_{3}} [GeV]",BasicCut&&Selection&&TriLepOF&&TCut("tri_id2==0"),"mllsub_OF_em");
702 >  WZ_kin_plot("tri_submll",40,0,200,"m_{l_{2}l_{3}} [GeV]",BasicCut&&Selection&&TriLepOF&&TCut("tri_id2==1"),"mllsub_OF_me");
703 >
704 >  
705 >  SF_OF_comparison(BasicCut,Selection);  
706 > //  generator_SF_OF_comparison();
707 > //  reco_SF_OF_comparison(BasicCut,Selection);
708 >
709 >  SF_vs_OF_plots("tri_submll",40,0,200,"m_{l_{2}l_{3}} [GeV]",BasicCut&&Selection,TriLepSF,TriLepOF,mc,"mllsub_MC");
710 >  SF_vs_OF_plots("tri_submll",40,0,200,"m_{l_{2}l_{3}} [GeV]",BasicCut&&Selection,TriLepSF,TriLepOF,data,"mllsub_data");
711 >  AlgorithmTest(BasicCut,Selection);
712 >  GetResult(BasicCut,Selection);
713 >  
714 >
715 > }
716 >
717 > void Make_Chi2_Resolution_Minimization(TGraph *gr) {
718 >  TCanvas *can2 = new TCanvas("can2","can2");
719 >  can2->cd();
720 >  gStyle->SetOptFit(0);
721 >  gr->Draw("AP*");
722 >  gr->GetXaxis()->SetTitle("MET scale factor");
723 >  gr->GetYaxis()->SetTitle("#chi^{2}");
724 >  gr->GetXaxis()->CenterTitle();
725 >  gr->GetYaxis()->CenterTitle();
726 >  gr->Draw("AP*");
727 >  gr->Fit("pol3");
728 >  TF1 *fit = (TF1*)gr->GetFunction("pol3");
729 >  fit->SetLineWidth(2);
730 >  fit->SetLineColor(kBlue);
731 >  gr->Draw("AP*");
732 >  DrawPrelim();
733 >  CompleteSave(can2,"WZ_study/DYResolution/Result_Chi2");      
734 >  delete can2;
735 > }
736 >
737 > void ExtractDYMissingResolution() {
738 >  
739 >  TCanvas *can = new TCanvas("can","can");
740 >  
741 >  TH1F *histoSF = allsamples.Draw("histoSF","met[4]",50,0,100,"MET [GeV]","events",TCut("id1==id2&&mll>0&&pfJetGoodNum40>=2"),data,luminosity);
742 >  TH1F *histoOF = allsamples.Draw("histoOF","met[4]",50,0,100,"MET [GeV]","events",TCut("id1!=id2&&mll>0&&pfJetGoodNum40>=2"),data,luminosity);
743 >  THStack histoMC = allsamples.DrawStack("histoMC","met[4]",50,0,100,"MET [GeV]","events",TCut("id1==id2&&mll>0&&pfJetGoodNum40>=2"),mc,luminosity);
744 >  TFile *f = new TFile("file.root","RECREATE");
745 >  
746 >  TH1F *MC_wo_DY = CollapseStack(histoMC);
747 >  TIter nextSF = TIter(histoMC.GetHists());
748    TH1F *h;
749    while ( h = (TH1F*)nextSF() ) {
750 <    cout << h->GetName() << " : " << h->Integral() << endl;
750 >    if(Contains(h->GetName(),"Z_Jets")) {
751 >      MC_wo_DY->Add(h,-1);
752 >    }
753    }
754 <  cout << "DATA : " << d_zmass->Integral() << endl;
755 <
756 <  TPad *pad2 = new TPad("pad","pad",0,0,1,1);
757 <  pad2->cd();
754 >  MC_wo_DY->Write();
755 >    
756 >  
757 >  int ntests=150;
758 >  TH1F *DYres[ntests];
759 >  
760 >  float min_chi2=-1;
761 >  int best_index=-1;
762 >  stringstream allresults;
763 >  TGraph *gr = new TGraph();
764 >  for(int i=0;i<ntests;i++) {
765 >    
766 >    float smearby=0.15*i/((float)ntests);
767 >    
768 >    stringstream drawvariable;
769 >    drawvariable << "((" << 1+smearby << ")*met[4])";
770 >    cout << "Going to use the variable " << drawvariable.str() << endl;
771 >    stringstream name;
772 >    name << "histo_DY_" << i;
773 >    DYres[i] = allsamples.Draw(name.str().c_str(),drawvariable.str(),50,0,100,"MET [GeV]","events",TCut("id1==id2&&mll>0&&pfJetGoodNum40>=2"),mc,luminosity, allsamples.FindSample("DYJets"));
774 >    DYres[i]->Write();
775 >    
776 >    TPad *pad = new TPad("pad","pad",0,0,1,1);
777 >    pad->cd();
778 >    
779 >    DrawPrelim();
780 >    MC_wo_DY->SetFillColor(kWhite);
781 >    MC_wo_DY->SetLineColor(kBlue);
782 >    TLegend *leg = make_legend();
783 >    leg->SetY1(0.8);
784 >    leg->AddEntry(histoSF,"data","pl");
785 >    leg->AddEntry(MC_wo_DY,"MC","l");
786 >    
787 >    TH1F *sum = (TH1F*)MC_wo_DY->Clone("sum");
788 >    sum->Add(DYres[i]);
789 >    sum->Draw("histo");
790 >    histoSF->Draw("e1,same");
791 >    stringstream Result;
792 >    
793 >    Double_t chi2;
794 >    Int_t ndf,igood;
795 >    Double_t res=0.0;
796 >    Double_t chi2prob = MarcosChi2TestX(histoSF,sum,chi2,ndf,igood,"UW");
797 >    Result << "#splitline{Smeared by : " << std::setprecision(3) << 100*smearby << " %}{#splitline{#chi^{2} / NDF : " << std::setprecision(3) << chi2 << " / " << ndf << "}{ #chi^{2} prob : " << std::setprecision(3)  << chi2prob << "}}";
798 >    
799 >    allresults << "\n" << smearby << ";" << chi2;
800 >    gr->SetPoint(i,smearby,chi2);
801 >    
802 >    if(min_chi2<0||chi2<min_chi2) {
803 >      min_chi2=chi2;
804 >      best_index=i;
805 >    }
806 >    
807 >    TPaveText *PurityPaveText = new TPaveText(0.60, 0.40,0.89, 0.6,"blNDC");
808 >    PurityPaveText->SetFillStyle(4000);
809 >    PurityPaveText->SetBorderSize(0);
810 >    PurityPaveText->SetFillColor(kWhite);
811 >    PurityPaveText->SetTextFont(42);
812 >    PurityPaveText->SetTextSize(0.042);
813 >    PurityPaveText->AddText(Result.str().c_str());
814 >    PurityPaveText->Draw();
815 >    DrawPrelim();
816 >    leg->Draw();
817 >    
818 >    
819 >    stringstream savethis;
820 >    savethis << "WZ_study/DYResolution/Response" << 100*smearby << "Percent";
821 >    
822 >    save_with_ratio(histoSF,sum,pad->cd(),savethis.str());
823 >  }
824 >  
825 >  gr->SetTitle("Chi2Diagram");
826 >  gr->SetName("Chi2Diagram");
827    
71  THStack m_met = allsamples.DrawStack("m_met","met[4]",40,0,400, "raw PFMET [GeV]", "events", ZmassWindow&&BasicCut,mc, luminosity);
72  TH1F *d_met = allsamples.Draw("d_met","met[4]",40,0,400, "raw PFMET [GeV]", "events", ZmassWindow&&BasicCut,data, luminosity);
828    
829 <  cout << "And with MET , using a 10 GeV window around the Z mass : " << endl;
830 <  nextSF = TIter(m_met.GetHists());
829 >  dout << "Result summary: " << endl;
830 >  dout << allresults.str() << endl;
831 >  
832 >  nextSF = TIter(histoMC.GetHists());
833 >  THStack stack;
834    while ( h = (TH1F*)nextSF() ) {
835 <    cout << h->GetName() << " : " << h->Integral() << endl;
835 >    if(!Contains(h->GetName(),"Z_Jets")) {
836 >      stack.Add(h);
837 >    }
838    }
79  cout << "DATA : " << d_met->Integral() << endl;
839    
840    
841 +  TPad *pad = new TPad("pad","pad",0,0,1,1);
842 +  pad->cd();
843 +
844 +  DYres[best_index]->SetFillColor(kYellow);
845 +  DYres[best_index]->SetLineColor(kYellow);
846    
847 <  d_met->Draw("e1");
848 <  m_met.Draw("same");
849 <  d_met->Draw("e1,same");
86 <  kinleg->Draw();
847 >  stack.Add(DYres[best_index]);
848 >  stack.Draw();
849 >  histoSF->Draw("e1,same");
850    DrawPrelim();
851 <  save_with_ratio(d_met,m_met,pad2->cd(),"WZ_study/MET");
851 >  save_with_ratio(histoSF,CollapseStack(stack),pad->cd(),"WZ_study/DYResolution/Result");
852 >  
853 >  TFile *fr = new TFile("rumba.root","RECREATE");
854 >  gr->Write();
855 >  fr->Close();
856    
90
857    delete can;
858 +  
859 +  Make_Chi2_Resolution_Minimization(gr);
860 +
861   }
862  
863   void WZstudy() {
864 + //  ExtractDYMissingResolution();
865    
866    switch_overunderflow(true);
867    cout << "Essential cut is " << (const char*) essentialcut << endl;
868    cout << "Going to set essential cut to new triggers" << endl;
869    TCut essential_bkp = essentialcut;
870    write_warning(__FUNCTION__,"Need to define trigger requirement for WZ!");
871 <  essentialcut=TCut("");
871 >  essentialcut=TCut("mll>5||mll<6");
872    cout << "Essential cut is now " << (const char*) essentialcut << endl;
873    
874    CarryOutWZStudy();
# Line 108 | Line 878 | void WZstudy() {
878    essentialcut = essential_bkp;
879    cout << "Essential cut is now " << (const char*) essentialcut << "   (restored)" << endl;
880    switch_overunderflow(false);
881 +  
882   }
883 +

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines