ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/cbrown/Development/Plotting/Modules/MetPlotting.C
Revision: 1.35
Committed: Wed Sep 19 17:14:23 2012 UTC (12 years, 7 months ago) by buchmann
Content type: text/plain
Branch: MAIN
Changes since 1.34: +9 -9 lines
Log Message:
Changed MET prediction mode to high purity (i.e. 10 GeV window instead of 20 GeV window); added HT requirement for Aachen selection again for JZB-based MET Z estimate (what a word...)

File Contents

# User Rev Content
1 buchmann 1.1 #include <iostream>
2    
3     using namespace std;
4    
5 buchmann 1.17 float Get_Met_Z_Prediction(TCut JetCut, float MetCut, int isdata, bool isDYonly);
6    
7 buchmann 1.24 namespace MetPlotsSpace {
8 buchmann 1.28 float Zprediction_Uncertainty=0.2;
9 buchmann 1.30 float Zestimate__data=-1;
10     float Zestimate__data_sys=-1;
11     float Zestimate__data_stat=-1;
12     float Zestimate__mc=-1;
13     float Zestimate__mc_sys=-1;
14     float Zestimate__mc_stat=-1;
15     float Zestimate__dy=-1;
16     float Zestimate__dy_sys=-1;
17     float Zestimate__dy_stat=-1;
18 buchmann 1.24 }
19 buchmann 1.17
20 buchmann 1.30 void ExperimentalMetPrediction(bool QuickRun);
21    
22     void ExtractScaleFactor(TH1F *mllSF,TH1F *mllOF, THStack* mcMllSF, THStack* mcMllOF, TH1F *prediction, TLegend *leg, string saveasSig, TBox *srbox) {
23     Int_t minbin = mllSF->FindBin(20.);
24     Int_t maxbin = mllSF->FindBin(70.-1);
25    
26     // Get yields in OF region
27     Float_t iDataOF = mllOF->Integral();
28     Float_t iDataOFSR = mllOF->Integral(minbin,maxbin);
29     Float_t iMCOF = 0.0;
30     Float_t iMCOFSR = 0.0;
31     TIter nextOF(mcMllOF->GetHists());
32     TH1F* h;
33     while ( h = (TH1F*)nextOF() ) {
34     iMCOF += h->Integral();
35     iMCOFSR += h->Integral(minbin,maxbin);
36     }
37     Float_t scale = iDataOF/iMCOF;
38    
39     // Re-scale OF
40     nextOF = TIter(mcMllOF->GetHists());
41    
42     while ( h = (TH1F*)nextOF() ) {
43     h->Scale(scale);
44     }
45    
46     nextOF = TIter(mcMllOF->GetHists());
47    
48     // Rescale SF and count in signal region
49     std::cout << "Integrating from " << mllSF->GetBinLowEdge(minbin) << " to " << mllSF->GetBinLowEdge(maxbin)+mllSF->GetBinWidth(maxbin) << std::endl;
50    
51     Float_t iDataSFSR = mllSF->Integral(minbin,maxbin);
52     Float_t iMCSFSR = 0.0;
53     TIter nextSF = TIter(mcMllSF->GetHists());
54     while ( h = (TH1F*)nextSF() ) {
55     h->Scale(scale);
56     iMCSFSR += h->Integral(minbin,maxbin);
57     }
58    
59     nextSF = TIter(mcMllSF->GetHists());
60     while ( h = (TH1F*)nextSF() ) {
61     iMCSFSR += h->Integral(minbin,maxbin);
62     }
63     mcMllSF->Modified();
64    
65     TPad* rcan2 = new TPad("rcan2","rcan2",0,0,1,1);
66     rcan2->cd();
67     mllSF->Draw();
68     mcMllSF->Draw("same");
69     prediction->Draw("histo,same");
70     mllSF->Draw("same");
71     DrawPrelim();
72     stringstream leghead;
73     leghead << "MC scaled by " << std::setprecision(2) << scale << "";
74     cout << "SCALE: " << scale << endl;
75     TH1F *histo = new TH1F("histo","histo",1,0,1);histo->SetLineColor(kWhite);
76     leg->AddEntry(histo,leghead.str().c_str(),"l");
77     leg->Draw();
78     srbox->Draw();
79     stringstream saveasSig2;
80     saveasSig2 << saveasSig << "__mcScaled";
81     rcan2->Update();
82     save_with_ratio( mllSF, *mcMllSF, rcan2, saveasSig2.str() );
83 buchmann 1.32
84     // restore original stacks
85     nextOF = TIter(mcMllOF->GetHists());
86    
87     while ( h = (TH1F*)nextOF() ) {
88     h->Scale(1/scale);
89     }
90    
91     nextSF = TIter(mcMllSF->GetHists());
92     while ( h = (TH1F*)nextSF() ) {
93     h->Scale(1/scale);
94     }
95     mcMllSF->Modified();
96     mcMllOF->Modified();
97    
98 buchmann 1.30 }
99    
100    
101    
102    
103    
104    
105    
106    
107    
108 buchmann 1.1 TGraphErrors* MakeErrorGraph(TH1F *histo) {
109    
110     float dx[histo->GetNbinsX()];
111     float dy[histo->GetNbinsX()];
112     float x[histo->GetNbinsX()];
113     float y[histo->GetNbinsX()];
114     for(int i=1;i<=histo->GetNbinsX();i++) {
115     x[i-1]=histo->GetBinCenter(i);
116     y[i-1]=histo->GetBinContent(i);
117     if(i>1) dx[i-1]=(histo->GetBinCenter(i)-histo->GetBinCenter(i-1))/2.0;
118     else dx[i-1]=(histo->GetBinCenter(i+1)-histo->GetBinCenter(i))/2.0;
119     dy[i-1]=histo->GetBinError(i);
120     }
121    
122     TGraphErrors *gr = new TGraphErrors(histo->GetNbinsX(),x,y,dx,dy);
123     gr->SetFillColor(TColor::GetColor("#2E9AFE"));
124     return gr;
125     }
126    
127 fronga 1.22 void ProduceMetPlotsWithCut(TCut cut, string name, float cutat, int njets, bool doMC = false, float ymax = 80 ) {
128 buchmann 1.30
129     bool UseSpecialZprediction=false;
130    
131     if(cutat==100 && name=="") {
132     UseSpecialZprediction=true;
133     bool ReRunEstimate=false;
134     //need to check if the results have already been stored; if not, need to get the estimate!
135     if(MetPlotsSpace::Zestimate__data<0) ReRunEstimate=true;
136     if(MetPlotsSpace::Zestimate__data_stat<0) ReRunEstimate=true;
137     if(MetPlotsSpace::Zestimate__data_sys<0) ReRunEstimate=true;
138     if(MetPlotsSpace::Zestimate__mc<0) ReRunEstimate=true;
139     if(MetPlotsSpace::Zestimate__mc_stat<0) ReRunEstimate=true;
140     if(MetPlotsSpace::Zestimate__mc_sys<0) ReRunEstimate=true;
141     if(MetPlotsSpace::Zestimate__dy<0) ReRunEstimate=true;
142     if(MetPlotsSpace::Zestimate__dy_stat<0) ReRunEstimate=true;
143     if(MetPlotsSpace::Zestimate__dy_sys<0) ReRunEstimate=true;
144     cout << "****************** About to do Z prediction " << endl;
145     if(ReRunEstimate) ExperimentalMetPrediction(true);//doing quick run (i.e. only data)
146     cout << "****************** Done predicting the Z " << endl;
147     }
148 fronga 1.16
149 buchmann 1.1 TCanvas *tcan = new TCanvas("tcan","tcan");
150     cout << "Doing met plots" << endl;
151 buchmann 1.2 stringstream MetBaseCuts;
152 fronga 1.7 MetBaseCuts << "met[4]>" << cutat << "&&" << cut.GetTitle();
153     stringstream snjets;
154     snjets << njets;
155 buchmann 1.2 TCut MetBaseCut(MetBaseCuts.str().c_str());
156 fronga 1.10 TCut nJetsSignal(PlottingSetup::basicqualitycut&&("pfJetGoodNum40>="+snjets.str()).c_str());
157 fronga 1.22 TCut nJetsControl(PlottingSetup::basiccut&&("met[4]>100&&met[4]<150&&pfJetGoodID[0]!=0&&pfJetGoodNum40=="+snjets.str()+"-1").c_str()); // Only njets vary
158 fronga 1.16
159     if ( !PlottingSetup::openBox ) {
160 buchmann 1.34 nJetsSignal += TCut("mll>70");
161 fronga 1.16 }
162 buchmann 1.1
163     //compute SF / OF rate in (CR1+CR2), should give 0.941 +/- 0.05
164 fronga 1.7
165     // Create histograms
166 fronga 1.9 //int nbins = 30;
167     int nbins = 60;
168 buchmann 1.17 float xmin=15., xmax = 315.;
169     TH1F *mllsigEE = allsamples.Draw("mllsigEE","mll",nbins,xmin,xmax,"m_{ee} [GeV]", "events",TCut(cutOSSF&&MetBaseCut&&nJetsSignal&&"id1==0"),data,PlottingSetup::luminosity);
170     TH1F *mllsigMM = allsamples.Draw("mllsigMM","mll",nbins,xmin,xmax,"m_{#mu#mu} [GeV]","events",TCut(cutOSSF&&MetBaseCut&&nJetsSignal&&"id1==1"),data,PlottingSetup::luminosity);
171 fronga 1.22 TH1F *mllscon = allsamples.Draw("mllscon","mll",nbins,xmin,xmax,"m_{ll} [GeV]", "events",TCut(cutOSSF&&cut&&nJetsControl),data,PlottingSetup::luminosity);
172 fronga 1.7 TH1F *mllOsig = allsamples.Draw("mllOsig", "mll",nbins,xmin,xmax,"m_{ll} [GeV]","events",TCut(cutOSOF&&MetBaseCut&&nJetsSignal),data,PlottingSetup::luminosity);
173 fronga 1.22 TH1F *mllOscon = allsamples.Draw("mllOscon","mll",nbins,xmin,xmax,"m_{ll} [GeV]","events",TCut(cutOSOF&&cut&&nJetsControl),data,PlottingSetup::luminosity);
174 fronga 1.7
175 fronga 1.8 TH1F* mllsig = (TH1F*)mllsigEE->Clone("mllsig");
176     mllsig->Add(mllsigMM);
177     mllsig->GetXaxis()->SetTitle("m_{ll} [GeV]");
178    
179     THStack *mcMllsig, *mcMllsigEE,*mcMllsigMM,*mcMllscon,*mcMllsconEE,*mcMllsconMM, *mcMllOsig, *mcMllOscon;
180 fronga 1.7 if ( doMC ) {
181     name += "_mc";
182 fronga 1.8 mcMllsig = new THStack(allsamples.DrawStack("mcMllsig","mll",nbins,xmin,xmax,"m_{ll} [GeV]","events",TCut(cutOSSF&&MetBaseCut&&nJetsSignal),mc,PlottingSetup::luminosity));
183     mcMllsigEE = new THStack(allsamples.DrawStack("mcMllsigEE","mll",nbins,xmin,xmax,"m_{ee} [GeV]","events",TCut(cutOSSF&&MetBaseCut&&nJetsSignal&&"id1==0"),mc,PlottingSetup::luminosity));
184     mcMllsigMM = new THStack(allsamples.DrawStack("mcMllsigMM","mll",nbins,xmin,xmax,"m_{#mu#mu} [GeV]","events",TCut(cutOSSF&&MetBaseCut&&nJetsSignal&&"id1==1"),mc,PlottingSetup::luminosity));
185 fronga 1.22 mcMllscon = new THStack(allsamples.DrawStack("mcMllscon","mll",nbins,xmin,xmax,"m_{ll} [GeV]","events",TCut(cutOSSF&&cut&&nJetsControl),mc,PlottingSetup::luminosity));
186 fronga 1.8 mcMllOsig = new THStack(allsamples.DrawStack("mcMllOsig","mll",nbins,xmin,xmax,"m_{ll} [GeV]","events",TCut(cutOSOF&&MetBaseCut&&nJetsSignal),mc,PlottingSetup::luminosity));
187 fronga 1.22 mcMllOscon= new THStack(allsamples.DrawStack("mcMllOscon","mll",nbins,xmin,xmax,"m_{ll} [GeV]","events",TCut(cutOSOF&&cut&&nJetsControl),mc,PlottingSetup::luminosity));
188 fronga 1.8
189 fronga 1.7 }
190 fronga 1.8
191 buchmann 1.1 mllOsig->SetLineColor(kRed);
192     mllOscon->SetLineColor(kRed);
193    
194 buchmann 1.28 TH1F *zlineshape = allsamples.Draw("zlineshape","mll",nbins,xmin,xmax,"m_{ll} (GeV)","events",cutOSSF&&TCut("pfJetGoodNum40==2")&&cut,data,PlottingSetup::luminosity);
195     TH1F *Ozlineshape = allsamples.Draw("Ozlineshape","mll",nbins,xmin,xmax,"m_{ll} (GeV)","events",cutOSOF&&TCut("pfJetGoodNum40==2")&&cut,data,PlottingSetup::luminosity);
196     zlineshape->Add(Ozlineshape,-1);
197 buchmann 1.13 TH1F *zlineshapeControl = (TH1F*)zlineshape->Clone("zlineshapeControl");
198 buchmann 1.21 // TH1F *zlineshapeFINE = allsamples.Draw("zlineshapeFINE","mll",50*nbins,xmin,xmax,"m_{ll} (GeV)","events",cutOSSF&&TCut("pfJetGoodNum40==1")&&cut,data,PlottingSetup::luminosity);
199     //
200     // float scalefactor = Get_Met_Z_Prediction(nJetsSignal,cutat, data, false) / (zlineshapeFINE->Integral(zlineshapeFINE->FindBin(91.1-20),zlineshapeFINE->FindBin(91.1+20)));
201     // float scalefactor_Control = Get_Met_Z_Prediction(nJetsControl,cutat, data, false) / (zlineshapeFINE->Integral(zlineshapeFINE->FindBin(91.1-20),zlineshapeFINE->FindBin(91.1+20)));
202     // delete zlineshapeFINE;
203    
204    
205 buchmann 1.30
206    
207 buchmann 1.21 Int_t scaleBinLow = mllsig->FindBin(86);
208     Int_t scaleBinHigh = mllsig->FindBin(94);
209     float scalefactor = (mllsig->Integral(scaleBinLow,scaleBinHigh)-mllOsig->Integral(scaleBinLow,scaleBinHigh));
210     scalefactor /= zlineshape->Integral(scaleBinLow,scaleBinHigh);
211 buchmann 1.30
212 buchmann 1.21 float scalefactor_Control = (mllscon->Integral(scaleBinLow,scaleBinHigh)-mllOscon->Integral(scaleBinLow,scaleBinHigh));
213     scalefactor_Control /= zlineshapeControl->Integral(scaleBinLow,scaleBinHigh);
214    
215     cout << "Bins for scaling : " << scaleBinLow << " : " << scaleBinHigh << endl;
216 buchmann 1.30
217     if(UseSpecialZprediction) {
218 buchmann 1.33 scaleBinLow = mllsig->FindBin(81);
219     scaleBinHigh = mllsig->FindBin(101);
220 buchmann 1.30 scalefactor = MetPlotsSpace::Zestimate__data/ (zlineshape->Integral(scaleBinLow,scaleBinHigh));
221     cout << "Dividing: " << MetPlotsSpace::Zestimate__data << " by " << (zlineshape->Integral(scaleBinLow,scaleBinHigh)) << endl;
222     write_warning(__FUNCTION__,"Not using JZB prediction for control region!");
223     }
224    
225 buchmann 1.21 cout << "Scale factors : " << scalefactor << " : " << scalefactor_Control << endl;
226 buchmann 1.30 if(UseSpecialZprediction) cout << " NOTE: Used JZB prediction for scaling! (Bins )" << scaleBinLow << " to " << scaleBinHigh << endl;
227 buchmann 1.17
228 buchmann 1.3 zlineshape->Scale(scalefactor);
229 buchmann 1.30
230 buchmann 1.3 zlineshape->SetLineColor(kBlue);
231     zlineshape->SetLineStyle(2);
232    
233 buchmann 1.30 if(UseSpecialZprediction) {
234     //need to update each bin with correct stat uncert
235     float relDYerr = MetPlotsSpace::Zestimate__data_stat/MetPlotsSpace::Zestimate__data;
236     for(int iz=1;iz<=zlineshape->GetNbinsX();iz++) {
237     float bincontent=zlineshape->GetBinContent(iz);
238     float binerror=zlineshape->GetBinError(iz);
239     float finalerr=0;
240     if(bincontent>0) finalerr+= (binerror/bincontent) * (binerror/bincontent);
241     if(MetPlotsSpace::Zestimate__data>0) finalerr+= relDYerr*relDYerr;
242     finalerr=bincontent * TMath::Sqrt(finalerr);
243     zlineshape->SetBinError(iz,finalerr);
244     }
245     }
246    
247 buchmann 1.13 zlineshapeControl->Scale(scalefactor_Control);
248     zlineshapeControl->SetLineColor(kBlue);
249     zlineshapeControl->SetLineStyle(2);
250    
251 buchmann 1.3 TH1F *subtracted = (TH1F*)mllsig->Clone("subtracted");
252 buchmann 1.13 TH1F *subtractedControl = (TH1F*)mllscon->Clone("subtractedControl");
253 buchmann 1.3 TH1F *baseline = (TH1F*)mllOsig->Clone("baseline");
254 buchmann 1.13 TH1F *baselineControl = (TH1F*)mllOscon->Clone("baselineControl");
255 buchmann 1.3 for(int i=1;i<=subtracted->GetNbinsX();i++) {
256     subtracted->SetBinContent(i,mllsig->GetBinContent(i)-mllOsig->GetBinContent(i));
257 buchmann 1.13 subtractedControl->SetBinContent(i,mllscon->GetBinContent(i)-mllOscon->GetBinContent(i));
258 buchmann 1.3 baseline->SetBinContent(i,0);
259 buchmann 1.13 baselineControl->SetBinContent(i,0);
260 buchmann 1.3 }
261    
262     TH1F *prediction = (TH1F*)mllOsig->Clone("prediction");
263     prediction->Add(zlineshape);
264     prediction->SetLineColor(TColor::GetColor("#CF35CA"));
265 fronga 1.7
266     TH1F *control_prediction = (TH1F*)mllOscon->Clone("control_prediction");
267     control_prediction->SetLineColor(TColor::GetColor("#CF35CA"));
268    
269     // FIX Y RANGE TO EASE COMPARISON
270 fronga 1.22 mllsig->SetMaximum(ymax);
271 buchmann 1.30 float PreviousMinimum=mllsig->GetMinimum();
272     mllsig->SetMinimum(0);
273 fronga 1.22 mllsigEE->SetMaximum(ymax);
274     mllsigMM->SetMaximum(ymax);
275     mllOsig->SetMaximum(ymax);
276     mllOscon->SetMaximum(ymax);
277     subtracted->SetMaximum(60);
278     subtracted->SetMinimum(-30);
279     subtractedControl->SetMaximum(65);
280     subtractedControl->SetMinimum(-30);
281 fronga 1.7
282 buchmann 1.3
283 fronga 1.7 // 1.- Signal region comparison
284 buchmann 1.30 TBox *srbox = new TBox(20,0,70,mllsig->GetMaximum());
285 buchmann 1.1 srbox->SetFillStyle(0);
286     srbox->SetLineColor(TColor::GetColor("#298A08"));
287     srbox->SetLineWidth(3);
288 buchmann 1.23
289 fronga 1.10
290 buchmann 1.2 stringstream MetHeader;
291 fronga 1.19 MetHeader << "N_{j}#geq" << snjets.str() << ", MET>" << cutat << " GeV";
292     stringstream MetHeaderCon;
293 fronga 1.22 MetHeaderCon << "N_{j}=" << njets-1 << ", 100<MET<150 GeV";
294 fronga 1.10 stringstream saveasSig;
295     saveasSig << "MetPlots/mll_sig" << cutat << "__" << name;
296    
297     TLegend* leg;
298     if ( !doMC ) {
299     tcan->cd();
300     //mllsig->GetYaxis()->SetRangeUser(0,mllsig->GetMaximum()*1.2);
301     mllsig->Draw();
302     TGraphErrors *stat3j = MakeErrorGraph(prediction);
303     stat3j->Draw("2,same");
304     mllOsig->Draw("histo,same");
305     zlineshape->Draw("histo,same");
306     prediction->Draw("histo,same");
307     mllsig->Draw("same");
308     DrawPrelim();
309     leg = make_legend();
310 fronga 1.19 leg->SetX1(0.52);
311     leg->SetHeader(MetHeader.str().c_str());
312 fronga 1.10 leg->AddEntry(mllsig,"Data","PL");
313     leg->AddEntry(prediction,"All bg prediction","L");
314 fronga 1.7 leg->AddEntry(mllOsig,"bg without Z","L");
315 buchmann 1.30 if(!UseSpecialZprediction) leg->AddEntry(zlineshape,"Z lineshape","L");
316     else leg->AddEntry(zlineshape,"bg with Z (JZB)","L");
317 fronga 1.7 leg->AddEntry(stat3j,"stat. uncert.","F");
318 fronga 1.10 leg->AddEntry(srbox,"SR","F");
319     leg->Draw();
320     srbox->Draw();
321     CompleteSave(tcan,saveasSig.str());
322     } else {
323 fronga 1.19
324 fronga 1.10 TPad* rcan = new TPad("rcan","rcan",0,0,1,1);
325     rcan->cd();
326     mllsig->Draw();
327     mcMllsig->Draw("same");
328     prediction->Draw("histo,same");
329     mllsig->Draw("same");
330     DrawPrelim();
331     leg = allsamples.allbglegend();
332 fronga 1.19 leg->SetHeader(MetHeader.str().c_str());
333     leg->SetX1(0.52);
334 fronga 1.10 leg->AddEntry(prediction,"All bg prediction","L");
335     leg->AddEntry(srbox,"SR","F");
336     leg->Draw();
337     srbox->Draw();
338     save_with_ratio( mllsig, *mcMllsig, rcan, saveasSig.str() );
339 buchmann 1.30
340     ExtractScaleFactor(mllsig,mllOsig,mcMllsig,mcMllOsig,prediction,leg,saveasSig.str(),srbox);
341 fronga 1.7 }
342 buchmann 1.1
343 fronga 1.8 // 1b. MC: split ee and mumu
344     if ( doMC ) {
345 fronga 1.10 TPad* rcan = new TPad("rcan","rcan",0,0,1,1);
346     rcan->cd();
347 fronga 1.8 mllsigEE->Draw();
348     mcMllsigEE->Draw("same");
349     mllsigEE->Draw("same");
350     DrawPrelim();
351     leg->Draw();
352     srbox->Draw();
353 fronga 1.10 save_with_ratio( mllsigEE, *mcMllsigEE,rcan->cd(),saveasSig.str()+"_ee" );
354 fronga 1.8
355 fronga 1.10 rcan = new TPad("rcan","rcan",0,0,1,1);
356     rcan->cd();
357 fronga 1.8 mllsigMM->Draw();
358     mcMllsigMM->Draw("same");
359     mllsigMM->Draw("same");
360     DrawPrelim();
361     leg->Draw();
362     srbox->Draw();
363 fronga 1.10 save_with_ratio( mllsigMM, *mcMllsigMM,rcan,saveasSig.str()+"_mm" );
364 fronga 1.8 }
365 fronga 1.14
366     // 1c. MC: compare of and sf
367     if ( doMC ) {
368     TH1F* hMcMllsig = CollapseStack( *mcMllsig);
369 fronga 1.15 leg = allsamples.allbglegend("");
370 fronga 1.19 leg->SetHeader(MetHeader.str().c_str());
371 fronga 1.15 // Change "Data" label by hand
372     ((TLegendEntry*)leg->GetListOfPrimitives()->At(0))->SetLabel("Same-flavor (MC)");
373 fronga 1.14 TPad* rcan = new TPad("rcan","rcan",0,0,1,1);
374     rcan->cd();
375 fronga 1.22 hMcMllsig->SetMaximum(ymax);
376 fronga 1.14 hMcMllsig->Draw("E");
377     mcMllOsig->Draw("same,hist");
378     hMcMllsig->Draw("same,E");
379     DrawMCPrelim();
380 fronga 1.19 leg->SetX1(0.52);
381 fronga 1.14 leg->AddEntry(srbox,"SR","F");
382     leg->Draw();
383     srbox->Draw();
384     save_with_ratio( hMcMllsig, *mcMllOsig, rcan, saveasSig.str()+"_mconly");
385    
386     }
387 buchmann 1.1
388 fronga 1.7 // 2.- Signal region comparison - LOG scale
389 fronga 1.10 if ( !doMC ) {
390     tcan->cd();
391     mllsig->SetMinimum(0.2); // FIX Y RANGE TO EASE COMPARISON
392     //mllsig->SetMaximum(mllsig->GetMaximum()*4.0);
393     srbox->SetY2(mllsig->GetMaximum());
394     tcan->SetLogy(1);
395     stringstream saveasSig2;
396     saveasSig2 << "MetPlots/mll_sig_ZLINESHAPE_" << cutat << "__" << name;
397    
398     CompleteSave(tcan,saveasSig2.str());
399     tcan->SetLogy(0);
400     }
401 buchmann 1.1
402 fronga 1.7
403     // 3.- Signal region, background subtracted
404 fronga 1.8 if ( !doMC ) {
405 fronga 1.10 tcan->cd();
406 buchmann 1.13 for(int i=1;i<=subtracted->GetNbinsX();i++) {
407 fronga 1.8 subtracted->SetBinContent(i,subtracted->GetBinContent(i)-zlineshape->GetBinContent(i));
408 buchmann 1.13 subtractedControl->SetBinContent(i,subtractedControl->GetBinContent(i)-zlineshapeControl->GetBinContent(i));
409 fronga 1.8 }
410 buchmann 1.3
411 fronga 1.8 TGraphErrors *subtrerr = MakeErrorGraph(baseline);
412     subtracted->Draw();
413     subtrerr->Draw("2,same");
414     subtracted->Draw("same");
415     DrawPrelim();
416     TLegend *DiffLeg = make_legend();
417 fronga 1.22 DiffLeg->SetX1(0.4);
418 fronga 1.8 DiffLeg->SetFillStyle(0);
419 fronga 1.19 DiffLeg->SetHeader(MetHeader.str().c_str());
420 fronga 1.8 DiffLeg->AddEntry(subtracted,"observed - predicted","PL");
421 fronga 1.19 DiffLeg->AddEntry(subtrerr,"stat. uncert","F");
422 fronga 1.8 DiffLeg->AddEntry((TObject*)0,"","");
423     DiffLeg->AddEntry((TObject*)0,"","");
424     DiffLeg->Draw();
425    
426     stringstream saveasSigSub;
427     saveasSigSub << "MetPlots/mll_sig_SUBTRACTED_" << cutat << "__" << name;
428    
429 fronga 1.22 //CompleteSave(tcan,saveasSigSub.str());
430 buchmann 1.13
431     // 3a.- Control region, background subtracted
432     TGraphErrors *subtrerrControl = MakeErrorGraph(baselineControl);
433     subtractedControl->Draw();
434     subtrerrControl->Draw("2,same");
435     subtractedControl->Draw("same");
436     DrawPrelim();
437 fronga 1.19 DiffLeg->SetHeader(MetHeaderCon.str().c_str());
438 buchmann 1.13 DiffLeg->Draw();
439     saveasSigSub.str("");
440 fronga 1.16 saveasSigSub << "MetPlots/mll_con_SUBTRACTED_" << cutat << "__" << name;
441 fronga 1.22 //CompleteSave(tcan,saveasSigSub.str());
442 buchmann 1.13
443    
444    
445 fronga 1.8 // 4.- Signal region, background subtracted, errors added in quadrature
446     TGraphErrors *subtrerr2 = (TGraphErrors*)subtrerr->Clone("subtrerr2");
447 buchmann 1.13 for(int i=1;i<=subtrerr2->GetN();i++) {
448     subtrerr2->SetPoint(i-1,subtracted->GetBinCenter(i),subtracted->GetBinContent(i));
449     subtrerr2->SetPointError(i-1,subtrerr2->GetErrorX(i),TMath::Sqrt(subtrerr2->GetErrorY(i)*subtrerr2->GetErrorY(i)+subtracted->GetBinError(i)*subtracted->GetBinError(i)));
450 fronga 1.8 }
451     TLine* l = new TLine(subtracted->GetBinLowEdge(1),0.,subtracted->GetBinLowEdge(subtracted->GetNbinsX()-1)+subtracted->GetBinWidth(1),0.);
452     l->SetLineWidth(subtracted->GetLineWidth());
453     subtracted->Draw();
454     subtrerr2->Draw("2,same");
455     l->Draw("same");
456     subtracted->Draw("same");
457     DrawPrelim();
458     TLegend *DiffLeg2 = make_legend();
459 fronga 1.22 DiffLeg2->SetX1(0.4);
460 fronga 1.19 DiffLeg2->SetHeader(MetHeader.str().c_str());
461 fronga 1.8 DiffLeg2->SetFillStyle(0);
462     DiffLeg2->AddEntry(subtracted,"observed - predicted","PL");
463 fronga 1.19 DiffLeg2->AddEntry(subtrerr2,"stat. uncert","F");
464 fronga 1.8 DiffLeg2->AddEntry((TObject*)0,"","");
465     DiffLeg2->AddEntry((TObject*)0,"","");
466     DiffLeg2->Draw();
467    
468     stringstream saveasSigSub2;
469     saveasSigSub2 << "MetPlots/mll_sig_SUBTRACTED_quadr_" << cutat << "__" << name;
470 fronga 1.5
471 fronga 1.8 CompleteSave(tcan,saveasSigSub2.str());
472 fronga 1.5
473 buchmann 1.13
474    
475     //4a.- Control region, background subtracted, errors added in quadrature
476     TGraphErrors *subtrerr2Control = (TGraphErrors*)subtrerrControl->Clone("subtrerr2Control");
477     for(int i=1;i<=subtrerr2Control->GetN();i++) {
478     subtrerr2Control->SetPoint(i-1,subtractedControl->GetBinCenter(i),subtractedControl->GetBinContent(i));
479     float width=subtrerr2Control->GetErrorX(i);
480     if(i==subtrerr2Control->GetN()) width=subtrerr2Control->GetErrorX(i-1);
481     subtrerr2Control->SetPointError(i-1,width,TMath::Sqrt(subtrerr2Control->GetErrorY(i)*subtrerr2Control->GetErrorY(i)+subtractedControl->GetBinError(i)*subtractedControl->GetBinError(i)));
482     }
483     TLine* lControl = new TLine(subtractedControl->GetBinLowEdge(1),0.,subtractedControl->GetBinLowEdge(subtractedControl->GetNbinsX()-1)+subtractedControl->GetBinWidth(1),0.);
484     lControl->SetLineWidth(subtractedControl->GetLineWidth());
485     subtractedControl->Draw();
486     subtrerr2Control->Draw("2,same");
487     lControl->Draw("same");
488     subtractedControl->Draw("same");
489     DrawPrelim();
490 fronga 1.22 DiffLeg2->SetHeader(MetHeaderCon.str().c_str());
491 buchmann 1.13 DiffLeg2->Draw();
492    
493     saveasSigSub2.str("");
494 fronga 1.16 saveasSigSub2 << "MetPlots/mll_con_SUBTRACTED_quadr_" << cutat << "__" << name;
495 buchmann 1.13
496     CompleteSave(tcan,saveasSigSub2.str());
497    
498 fronga 1.8 delete DiffLeg;
499     delete DiffLeg2;
500 buchmann 1.13
501 fronga 1.8 } // !doMC
502 buchmann 1.3
503    
504 fronga 1.7 // 5.- Control region comparison
505 fronga 1.16 // scalefactor = (mllscon->Integral(scaleBinLow,scaleBinHigh)-mllOscon->Integral(scaleBinLow,scaleBinHigh));
506     // scalefactor /= zlineshape->Integral(scaleBinLow,scaleBinHigh);
507     // zlineshape->Scale(scalefactor);
508     control_prediction->Add(zlineshapeControl);
509    
510 fronga 1.22 control_prediction->SetMaximum(ymax); // FIX MAXIMUM TO EASE COMPARISON
511 fronga 1.7
512 fronga 1.14 TBox *cr1box = new TBox(xmin,0,70,control_prediction->GetMaximum());
513 buchmann 1.1 cr1box->SetFillStyle(0);
514     cr1box->SetLineColor(TColor::GetColor("#0404B4"));
515     cr1box->SetLineWidth(3);
516    
517 fronga 1.14 TBox *cr2box = new TBox(120,0,xmax,control_prediction->GetMaximum());
518 buchmann 1.1 cr2box->SetFillStyle(0);
519     cr2box->SetLineColor(TColor::GetColor("#0404B4"));
520     cr2box->SetLineWidth(3);
521     cr2box->SetLineStyle(2);
522    
523 fronga 1.10 stringstream saveasCon;
524     saveasCon << "MetPlots/mll_con" << cutat << "__" << name;
525    
526 fronga 1.7 TLegend *legc;
527 fronga 1.10 //control_prediction->GetYaxis()->SetRangeUser(0,control_prediction->GetMaximum()*1.3);
528     if ( !doMC ) {
529     tcan->cd();
530     control_prediction->Draw("hist");
531     TGraphErrors *stat2j = MakeErrorGraph(control_prediction);
532     stat2j->Draw("2,same");
533     mllOscon->Draw("same,hist");
534 fronga 1.16 zlineshapeControl->Draw("histo,same");
535 fronga 1.10 control_prediction->Draw("histo,same");
536     mllscon->Draw("same");
537     DrawPrelim();
538 fronga 1.7 legc = make_legend();
539 fronga 1.19 legc->SetX1(0.52);
540     legc->SetHeader(MetHeaderCon.str().c_str());
541 fronga 1.10 legc->AddEntry(mllscon,"Data","PL");
542     legc->AddEntry(control_prediction,"All bg","L");
543 fronga 1.16 legc->AddEntry(zlineshapeControl,"Z lineshape","L");
544 fronga 1.7 legc->AddEntry(mllOscon,"bg without Z","L");
545     legc->AddEntry(stat2j,"stat. uncert.","F");
546 fronga 1.10 legc->AddEntry(cr1box,"CR1","F");
547     legc->AddEntry(cr2box,"CR2","F");
548     legc->Draw();
549     cr1box->Draw();
550     cr2box->Draw();
551     CompleteSave(tcan,saveasCon.str());
552     } else {
553     TPad* rcan = new TPad("rcan","rcan",0,0,1,1);
554     rcan->cd();
555     control_prediction->Draw("hist");
556     mcMllscon->Draw("same");
557     control_prediction->Draw("histo,same");
558     mllscon->Draw("same");
559     DrawPrelim();
560     legc = allsamples.allbglegend();
561 fronga 1.19 legc->SetX1(0.52);
562     legc->SetHeader(MetHeaderCon.str().c_str());
563 fronga 1.10 legc->AddEntry(control_prediction,"All bg","L");
564     legc->AddEntry(cr1box,"CR1","F");
565     legc->AddEntry(cr2box,"CR2","F");
566     legc->Draw();
567     cr1box->Draw();
568     cr2box->Draw();
569     save_with_ratio( mllscon, *mcMllscon, rcan, saveasCon.str());
570     }
571 buchmann 1.1
572 fronga 1.7 // 6. - Opposite-flavour data/MC comparison
573     if ( doMC ) {
574 fronga 1.10 TPad* rcan = new TPad("rcan","rcan",0,0,1,1);
575     rcan->cd();
576 fronga 1.7 mllOsig->SetLineColor(kBlack);
577     mllOsig->Draw();
578     mcMllOsig->Draw("same");
579     mllOsig->Draw("same");
580     TLegend *legsdm = allsamples.allbglegend();
581 fronga 1.19 legsdm->SetHeader((MetHeader.str()+", OF").c_str());
582     legsdm->SetX1(0.52);
583 fronga 1.7 legsdm->Draw();
584     stringstream saveasSigOF;
585     saveasSigOF << "MetPlots/mll_sig_of_" << cutat << "__" << name;
586 fronga 1.10 save_with_ratio( mllOsig, *mcMllOsig, rcan, saveasSigOF.str());
587 fronga 1.7
588 fronga 1.10 rcan = new TPad("rcan","rcan",0,0,1,1);
589     rcan->cd();
590 fronga 1.7 mllOscon->SetLineColor(kBlack);
591     mllOscon->Draw();
592     mcMllOscon->Draw("same");
593     mllOscon->Draw("same");
594     TLegend *legcdm = allsamples.allbglegend();
595 fronga 1.19 legcdm->SetHeader((MetHeaderCon.str()+", OF").c_str());
596     legcdm->SetX1(0.52);
597 fronga 1.7 legcdm->Draw();
598     stringstream saveasConOF;
599     saveasConOF << "MetPlots/mll_con_of_" << cutat << "__" << name;
600 fronga 1.10 save_with_ratio( mllOscon, *mcMllOscon, rcan, saveasConOF.str());
601    
602 fronga 1.7 delete legsdm;
603     delete legcdm;
604 fronga 1.10 }
605 fronga 1.7
606 fronga 1.10 // Memory clean-up
607     if (doMC) {
608 fronga 1.7 delete mcMllscon;
609     delete mcMllOscon;
610     delete mcMllsig;
611 fronga 1.8 delete mcMllsigEE;
612     delete mcMllsigMM;
613 fronga 1.7 delete mcMllOsig;
614     }
615 buchmann 1.1
616     delete cr1box;
617     delete cr2box;
618     delete srbox;
619     delete legc;
620     delete leg;
621 fronga 1.7
622 buchmann 1.1 delete mllscon;
623     delete mllOscon;
624     delete mllsig;
625 fronga 1.8 delete mllsigEE;
626     delete mllsigMM;
627 buchmann 1.1 delete mllOsig;
628 fronga 1.7 delete zlineshape;
629 fronga 1.16 delete zlineshapeControl;
630 buchmann 1.6 delete tcan;
631 buchmann 1.1 }
632    
633 fronga 1.14 void DoMetPlots(string datajzb, string mcjzb) {
634 buchmann 1.27 switch_overunderflow(true);
635 fronga 1.7 float metCuts[] = { 100., 150. };
636 fronga 1.22 float ymax[] = { 90., 80. };
637 fronga 1.7 int jetCuts[] = { 3, 2 };
638 fronga 1.11 //int jetCuts[] = { 3, 3 };
639 buchmann 1.35 string leptCuts[] = { "pt1>20&&pt2>20", "pt1>20&&pt2>10&&pfTightHT>100" };
640 fronga 1.7 bool nomc(0),domc(1);
641     for ( int i=0; i<2; ++i ) {
642 buchmann 1.33 // ProduceMetPlotsWithCut(TCut(("mll>15&&"+leptCuts[i]).c_str()),"",metCuts[i],jetCuts[i],nomc,ymax[i]);
643 fronga 1.22 ProduceMetPlotsWithCut(TCut(("mll>15&&"+leptCuts[i]).c_str()),"",metCuts[i],jetCuts[i],domc,ymax[i]);
644     continue;
645 fronga 1.19 // stringstream jzbcut;
646     // jzbcut << "((is_data&&("<<datajzb<<")>0)||(!is_data&&("<<mcjzb<<")>0))";
647     // ProduceMetPlotsWithCut(TCut((jzbcut.str()+"&&mll>20&&"+leptCuts[i]).c_str()),"JZBpos",metCuts[i], jetCuts[i],domc);
648 fronga 1.11 //continue;
649 fronga 1.7 if (!PlottingSetup::is53reco) { // Old 5_2
650 fronga 1.22 ProduceMetPlotsWithCut(TCut(("mll>15&&pfJetGoodNumBtag==0&&"+leptCuts[i]).c_str()),"bTagVeto30",metCuts[i], jetCuts[i], nomc,ymax[i]);
651     ProduceMetPlotsWithCut(TCut(("mll>15&&pfJetGoodNumBtag>0&&"+leptCuts[i]).c_str()),"AtLeastOneBJet30",metCuts[i],jetCuts[i], nomc,ymax[i]);
652     ProduceMetPlotsWithCut(TCut(("mll>15&&pfJetGoodNumBtag==0&&"+leptCuts[i]).c_str()),"bTagVeto30",metCuts[i], jetCuts[i], domc,ymax[i]);
653     ProduceMetPlotsWithCut(TCut(("mll>15&&pfJetGoodNumBtag>0&&"+leptCuts[i]).c_str()),"AtLeastOneBJet30",metCuts[i],jetCuts[i], domc,ymax[i]);
654 fronga 1.7 } else {
655 fronga 1.22 ProduceMetPlotsWithCut(TCut(("mll>15&&pfJetGoodNumBtag30==0&&"+leptCuts[i]).c_str()),"bTagVeto30",metCuts[i], jetCuts[i],nomc,ymax[i]);
656     ProduceMetPlotsWithCut(TCut(("mll>15&&pfJetGoodNumBtag30>0&&"+leptCuts[i]).c_str()),"AtLeastOneBJet30",metCuts[i],jetCuts[i],nomc,ymax[i]);
657     ProduceMetPlotsWithCut(TCut(("mll>15&&pfJetGoodNumBtag30==0&&"+leptCuts[i]).c_str()),"bTagVeto30",metCuts[i], jetCuts[i],domc,ymax[i]);
658     ProduceMetPlotsWithCut(TCut(("mll>15&&pfJetGoodNumBtag30>0&&"+leptCuts[i]).c_str()),"AtLeastOneBJet30",metCuts[i], jetCuts[i],domc,ymax[i]);
659 fronga 1.7 }
660     }
661 buchmann 1.27 switch_overunderflow(false);
662 buchmann 1.1 }
663 buchmann 1.12
664 buchmann 1.13 void compute_r_in_out() {
665     int nbins=10;
666     float xmin=10;
667     float xmax=10;
668 buchmann 1.33 TCut InCut("mll>81&&mll<101&&pfJetGoodNum40>0");
669 buchmann 1.13 TCut OutCut("mll>20&&mll<70&&pfJetGoodNum40>0");
670    
671     TH1F *mllIN =allsamples.Draw("mllIN", "mll",nbins,xmin,xmax,"m_{ee} [GeV]", "events",TCut(cutOSSF&&InCut),mc,PlottingSetup::luminosity,allsamples.FindSample("DYJetsToLLNoTau"));
672     TH1F *mllOUT =allsamples.Draw("mllOUT","mll",nbins,xmin,xmax,"m_{ee} [GeV]", "events",TCut(cutOSSF&&OutCut),mc,PlottingSetup::luminosity,allsamples.FindSample("DYJetsToLLNoTau"));
673    
674     cout << "IN: " << mllIN->Integral() << endl;
675     cout << "OUT: " << mllOUT->Integral() << endl;
676     cout << "r_oi: " <<((float) mllOUT->Integral()) / mllIN->Integral()<< endl;
677    
678     delete mllIN;
679     delete mllOUT;
680    
681     }
682    
683 buchmann 1.17
684    
685     void LabelHisto(TH1 *MET_ratio,string titlex, string titley) {
686     MET_ratio->GetXaxis()->SetTitle(titlex.c_str());
687     MET_ratio->GetXaxis()->CenterTitle();
688     MET_ratio->GetYaxis()->SetTitle(titley.c_str());
689     MET_ratio->GetYaxis()->CenterTitle();
690     }
691    
692 buchmann 1.23 TH1F* GetPredictedAndObservedMetShapes(TCut JetCut, string sPositiveCut,string sNegativeCut,string CorrectedMet,string ObservedMet, string JZBPosvar, string JZBNegvar, float MetCut, int is_data, bool isDYonly, bool isAachen) {
693 buchmann 1.17
694     //Steps:
695     // 1) Prepare samples and histo definition (with "optimal" binning for MET cut)
696     // 2) Fill MET histograms
697     // 3) Fill JZB histograms
698     // 4) Draw them and store them
699     // 5) return predicted MET distribution as is (i.e. not scaled by factor of 2!)
700    
701 buchmann 1.24 cout << "*************************************" << endl;
702 buchmann 1.26 // cout << "** SUMMARY BEFORE STARTING DRAWING **" << endl;
703     // cout << "MET variable: " << ObservedMet << endl;
704     // cout << "Corr. MET var:" << CorrectedMet << endl;
705     // cout << "JZB pos. var: " << JZBPosvar << endl;
706     // cout << "JZB neg. var: " << JZBNegvar << endl;
707     // cout << "JZB pos cut : " << sPositiveCut << endl;
708     // cout << "JZB neg cut : " << sNegativeCut << endl;
709 buchmann 1.30
710    
711     if(isAachen) MetPlotsSpace::Zprediction_Uncertainty=0.3;
712 buchmann 1.17 //Step 1: Prepare samples and histo definition
713     vector<int> SelectedSamples;
714     if(is_data==mc&&isDYonly) {
715     SelectedSamples=allsamples.FindSample("Z_em_DYJetsToL");
716     if(SelectedSamples.size()==0) {
717     write_error(__FUNCTION__,"Cannot continue, there seems to be no DY sample without Taus - goodbye!");
718     assert(SelectedSamples.size()>0);
719     }
720     }
721    
722     float DisplayedBinSize=10.0; // this is the bin size that we use for plotting
723    
724 buchmann 1.21 float BinWidth=1.0;
725 buchmann 1.17 float xmin=0;
726 buchmann 1.23 float xmax=110;
727     if(isAachen) xmax=160;
728 buchmann 1.21 if(MetCut>=xmax) xmax=MetCut+10;
729 buchmann 1.17 int nbins=int((xmax-xmin)/BinWidth);
730 buchmann 1.23
731     float pt2cut=20;
732     if(isAachen)pt2cut=10;
733    
734 buchmann 1.17 stringstream basiccut;
735 buchmann 1.23 basiccut << (const char*) JetCut << "&&" << (const char*) Restrmasscut << "&&" << (const char*) leptoncut << "&&pt1>20&&pt2>" << pt2cut;
736 buchmann 1.17
737    
738     stringstream cMET_observed;
739     cMET_observed << "(" << basiccut.str() << "&&(" << sPositiveCut << ")&&" << (const char*) cutOSSF << ")";
740     stringstream cMET_ttbar_pred;
741     cMET_ttbar_pred << "(" << basiccut.str() << "&&(" << sPositiveCut << ")&&" << (const char*) cutOSOF << ")";
742     stringstream cMET_osof_pred;
743     cMET_osof_pred << "(" << basiccut.str() << "&&(" << sNegativeCut << ")&&" << (const char*) cutOSOF << ")";
744     stringstream cMET_ossf_pred;
745     cMET_ossf_pred << "(" << basiccut.str() << "&&(" << sNegativeCut << ")&&" << (const char*) cutOSSF << ")";
746    
747     write_warning(__FUNCTION__,"Once the rush is over you might want to define the potential sidebands ... ");
748 buchmann 1.28
749 buchmann 1.17 //Step 2: Fill Met histograms
750 buchmann 1.28 float bottommargin=gStyle->GetPadBottomMargin();
751     float canvas_height=gStyle->GetCanvasDefH();
752     float canvas_width=gStyle->GetCanvasDefW();
753     float ratiospace=0.25;// space the ratio should take up (relative to original pad)
754    
755     float ratiobottommargin=0.3;
756     float ratiotopmargin=0.1;
757    
758     float xstretchfactor=((1-ratiospace)*(1-gStyle->GetPadTopMargin()))/((1)*ratiospace);
759    
760     TCanvas *main_canvas = new TCanvas("main_canvas","main_canvas",(Int_t)canvas_width,(Int_t)(canvas_height*(1+ratiospace)));
761     TPad *mainpad = new TPad("mainpad","mainpad",0,1-(1.0/(1+ratiospace)),1,1);//top (main) pad
762     TPad *coverpad = new TPad("coverpad","coverpad",gStyle->GetPadLeftMargin()-0.008,1-(1.0/(1+ratiospace))-0.04,1,1-(1.0/(1+ratiospace))+0.103);//pad covering up the x scale
763     TPad *bottompad = new TPad("bottompad", "Ratio Pad",0,0,1,(1-(1-bottommargin)/(1+ratiospace))-0.015); //bottom pad
764    
765     main_canvas->Range(0,0,1,1);
766     main_canvas->SetBorderSize(0);
767     main_canvas->SetFrameFillColor(0);
768    
769     mainpad->Draw();
770     mainpad->cd();
771     mainpad->SetLogy(1);
772     mainpad->Range(0,0,1,1);
773     mainpad->SetFillColor(kWhite);
774     mainpad->SetBorderSize(0);
775     mainpad->SetFrameFillColor(0);
776    
777    
778    
779    
780 buchmann 1.17 TH1F *MET_observed = allsamples.Draw("MET_observed",ObservedMet,nbins,xmin,xmax,"MET [GeV]","events",
781 buchmann 1.26 TCut(cMET_observed.str().c_str()),is_data,PlottingSetup::luminosity,SelectedSamples);
782 buchmann 1.17 TH1F *MET_ossf_pred = allsamples.Draw("MET_ossf_pred",CorrectedMet,nbins,xmin,xmax,"MET [GeV]","events",
783 buchmann 1.26 TCut(cMET_ossf_pred.str().c_str()),is_data,PlottingSetup::luminosity,SelectedSamples);
784 buchmann 1.17 TH1F *MET_osof_pred = allsamples.Draw("MET_osof_pred",CorrectedMet,nbins,xmin,xmax,"MET [GeV]","events",
785 buchmann 1.26 TCut(cMET_osof_pred.str().c_str()),is_data,PlottingSetup::luminosity,SelectedSamples);
786 buchmann 1.17 TH1F *MET_ttbar_pred= allsamples.Draw("MET_ttbar_pred",ObservedMet,nbins,xmin,xmax,"MET [GeV]","events",
787 buchmann 1.26 TCut(cMET_ttbar_pred.str().c_str()),is_data,PlottingSetup::luminosity,SelectedSamples);
788 buchmann 1.17
789 buchmann 1.25
790     if(isDYonly && is_data==mc) {
791     TH1F *MET_truth = allsamples.Draw("MET_truth",ObservedMet,1,MetCut,10000,"MET [GeV]","events",TCut(((string)"met[4]>"+any2string(MetCut)).c_str())&&cutOSSF&&TCut(basiccut.str().c_str()),is_data,PlottingSetup::luminosity,SelectedSamples);
792     write_info(__FUNCTION__,"DY Truth is : "+any2string(MET_truth->Integral()));
793     delete MET_truth;
794     }
795    
796    
797 buchmann 1.17 TH1F *MET_predicted=(TH1F*)MET_ossf_pred->Clone("MET_predicted");
798     MET_predicted->Add(MET_osof_pred,-1);
799     MET_predicted->Add(MET_ttbar_pred);
800     MET_predicted->SetLineColor(kRed);
801     MET_observed->SetLineColor(kBlack);
802    
803     TH1F *MET_Z_prediction=(TH1F*)MET_ossf_pred->Clone("MET_Z_prediction");
804     MET_Z_prediction->Add(MET_osof_pred,-1);
805     MET_Z_prediction->SetLineColor(kBlue);
806    
807     LabelHisto(MET_observed,"MET (GeV)","events");
808    
809     //Step 3: Fill JZB histograms
810 buchmann 1.25
811 buchmann 1.17 TH1F *JZB_observed = allsamples.Draw("JZB_observed",JZBPosvar,nbins,xmin,xmax,"JZB [GeV]","events",
812     TCut(cMET_observed.str().c_str()),is_data,PlottingSetup::luminosity,SelectedSamples);
813     TH1F *JZB_ossf_pred = allsamples.Draw("JZB_ossf_pred",JZBNegvar,nbins,xmin,xmax,"JZB [GeV]","events",
814     TCut(cMET_ossf_pred.str().c_str()),is_data,PlottingSetup::luminosity,SelectedSamples);
815     TH1F *JZB_osof_pred = allsamples.Draw("JZB_osof_pred",JZBNegvar,nbins,xmin,xmax,"JZB [GeV]","events",
816     TCut(cMET_osof_pred.str().c_str()),is_data,PlottingSetup::luminosity,SelectedSamples);
817     TH1F *JZB_ttbar_pred= allsamples.Draw("JZB_ttbar_pred",JZBPosvar,nbins,xmin,xmax,"JZB [GeV]","events",
818     TCut(cMET_ttbar_pred.str().c_str()),is_data,PlottingSetup::luminosity,SelectedSamples);
819    
820     TH1F *JZB_predicted=(TH1F*)JZB_ossf_pred->Clone("JZB_predicted");
821     JZB_predicted->Add(JZB_osof_pred,-1);
822     JZB_predicted->Add(JZB_ttbar_pred);
823     JZB_predicted->SetLineColor(kRed);
824     JZB_observed->SetLineColor(kBlack);
825    
826     TH1F *JZB_Z_prediction=(TH1F*)JZB_ossf_pred->Clone("JZB_Z_prediction");
827     JZB_Z_prediction->Add(JZB_osof_pred,-1);
828     MET_Z_prediction->SetLineColor(kBlue);
829    
830     LabelHisto(JZB_observed,"JZB (GeV)","events");
831    
832     // Step 4: Draw them and store them
833    
834     TLegend *legend = new TLegend(0.6,0.6,0.89,0.89);
835    
836     MET_ttbar_pred->SetLineColor(TColor::GetColor("#005C00"));
837     JZB_ttbar_pred->SetLineColor(TColor::GetColor("#005C00"));
838    
839     legend->SetFillColor(kWhite);
840     legend->SetBorderSize(0);
841     legend->AddEntry(MET_predicted,"prediction","l");
842     legend->AddEntry(MET_observed,"observed","p");
843     legend->AddEntry(MET_Z_prediction,"predicted Z","l");
844     legend->AddEntry(MET_ttbar_pred,"OF-based prediction","l");
845    
846     if(is_data==mc) legend->SetHeader("Simulation:");
847     if(is_data==mc&&isDYonly) legend->SetHeader("DY #rightarrow ee,#mu#mu only:");
848     if(is_data==data) legend->SetHeader("Data:");
849    
850     stringstream SaveJZBname;
851     stringstream SaveMETname;
852     if(is_data==data) {
853     SaveJZBname << "MetPrediction/JZBdistribution_Data_METCutAt" << MetCut;
854     SaveMETname << "MetPrediction/METdistribution_Data_METCutAt" << MetCut;
855     }
856     if(is_data==mc&&!isDYonly) {
857     SaveJZBname << "MetPrediction/JZBdistribution_FullMC_METCutAt" << MetCut;
858     SaveMETname << "MetPrediction/METdistribution_FullMC_METCutAt" << MetCut;
859     }
860     if(is_data==mc&&isDYonly) {
861     SaveJZBname << "MetPrediction/JZBdistribution_DYMC_METCutAt" << MetCut;
862     SaveMETname << "MetPrediction/METdistribution_DYMC_METCutAt" << MetCut;
863     }
864    
865 buchmann 1.26 dout << "Shape summary (MET>50) for ";
866     if(is_data==data) cout << "data";
867     if(is_data==mc&&isDYonly) cout<< "DY ";
868     if(is_data==mc&&!isDYonly) cout << " Full MC";
869     cout << " : " << endl;
870    
871 buchmann 1.24 dout << " observed : " << MET_observed->Integral(MET_observed->FindBin(50),MET_observed->FindBin(xmax)) << endl;
872     dout << " predicted : " << MET_predicted->Integral(MET_predicted->FindBin(50),MET_predicted->FindBin(xmax)) << endl;
873     dout << " Z pred : " << MET_Z_prediction->Integral(MET_Z_prediction->FindBin(50),MET_Z_prediction->FindBin(xmax)) << endl;
874     dout << " ttbar : " << MET_ttbar_pred->Integral(MET_ttbar_pred->FindBin(50),MET_ttbar_pred->FindBin(xmax)) << endl;
875    
876    
877 buchmann 1.17 TH1F *ZpredClone = (TH1F*)MET_Z_prediction->Clone("ZpredClone");
878     ZpredClone->SetLineColor(kBlue);
879     MET_observed->Rebin(int(DisplayedBinSize/BinWidth));
880     ZpredClone->Rebin(int(DisplayedBinSize/BinWidth));
881     MET_predicted->Rebin(int(DisplayedBinSize/BinWidth));
882     MET_ttbar_pred->Rebin(int(DisplayedBinSize/BinWidth));
883    
884     TH1F *JZBZpredClone = (TH1F*)JZB_Z_prediction->Clone("ZpredClone");
885     JZBZpredClone->SetLineColor(kBlue);
886     JZB_observed->Rebin(int(DisplayedBinSize/BinWidth));
887     JZBZpredClone->Rebin(int(DisplayedBinSize/BinWidth));
888     JZB_predicted->Rebin(int(DisplayedBinSize/BinWidth));
889     JZB_ttbar_pred->Rebin(int(DisplayedBinSize/BinWidth));
890    
891     TH1F *JZB_ratio = (TH1F*)JZB_observed->Clone("JZB_ratio");
892     JZB_ratio->Divide(JZB_predicted);
893     LabelHisto(JZB_ratio,"JZB (GeV)","obs/pred");
894     TH1F *MET_ratio = (TH1F*)MET_observed->Clone("MET_ratio");
895     MET_ratio->Divide(MET_predicted);
896 buchmann 1.28 MET_observed->SetMaximum(5*MET_observed->GetMaximum());
897     JZB_observed->SetMaximum(5*JZB_observed->GetMaximum());
898     MET_observed->SetMinimum(0.5);
899     JZB_observed->SetMinimum(0.5);
900 buchmann 1.17 LabelHisto(MET_ratio,"MET (GeV)","obs/pred");
901 buchmann 1.24 TBox *sysenvelope = new TBox(xmin,1.0-MetPlotsSpace::Zprediction_Uncertainty,xmax,1.0+MetPlotsSpace::Zprediction_Uncertainty);
902 buchmann 1.17 sysenvelope->SetFillColor(TColor::GetColor("#82FA58")); // light green
903     sysenvelope->SetLineWidth(0);
904 buchmann 1.24 TBox *dsysenvelope = new TBox(xmin,1.0-2*MetPlotsSpace::Zprediction_Uncertainty,xmax,1.0+2*MetPlotsSpace::Zprediction_Uncertainty);
905     dsysenvelope->SetFillColor(TColor::GetColor("#F3F781")); // light yellow
906     dsysenvelope->SetLineWidth(0);
907 buchmann 1.28
908     MET_ratio->GetYaxis()->SetNdivisions(502,false);
909     JZB_ratio->GetYaxis()->SetNdivisions(502,false);
910    
911 buchmann 1.17
912     MET_observed->Draw("e1");
913     MET_ttbar_pred->Draw("histo,same");
914     ZpredClone->Draw("histo,same");
915     MET_predicted->Draw("histo,same");
916     MET_observed->Draw("e1,same");
917     legend->Draw();
918     if(is_data==data) DrawPrelim();
919     else DrawMCPrelim();
920    
921 buchmann 1.28 mainpad->Modified();
922     main_canvas->cd();
923     coverpad->Draw();
924     coverpad->cd();
925     coverpad->Range(0,0,1,1);
926     coverpad->SetFillColor(kWhite);
927     coverpad->SetBorderSize(0);
928     coverpad->SetFrameFillColor(0);
929     coverpad->Modified();
930     main_canvas->cd();
931     bottompad->SetTopMargin(ratiotopmargin);
932     bottompad->SetBottomMargin(ratiobottommargin);
933     bottompad->Draw();
934 buchmann 1.17 bottompad->cd();
935 buchmann 1.28 bottompad->Range(0,0,1,1);
936     bottompad->SetFillColor(kWhite);
937    
938 buchmann 1.17 MET_ratio->GetYaxis()->SetRangeUser(0,2);
939 buchmann 1.28 MET_ratio->GetXaxis()->SetLabelSize(xstretchfactor*MET_ratio->GetXaxis()->GetLabelSize());
940     MET_ratio->GetYaxis()->SetLabelSize(xstretchfactor*MET_ratio->GetYaxis()->GetLabelSize());
941     MET_ratio->GetXaxis()->SetTitleSize(xstretchfactor*gStyle->GetTitleSize());
942     MET_ratio->GetYaxis()->SetTitleSize(xstretchfactor*gStyle->GetTitleSize());
943    
944 buchmann 1.17 MET_ratio->Draw("e1");
945 buchmann 1.28 // dsysenvelope->Draw();
946 buchmann 1.17 sysenvelope->Draw();
947     MET_ratio->Draw("AXIS,same");
948     MET_ratio->Draw("e1,same");
949     TLine *metl = new TLine(xmin,1.0,xmax,1.0);
950     metl->SetLineColor(kBlue);
951     metl->Draw();
952 buchmann 1.28 CompleteSave(main_canvas,SaveMETname.str());
953    
954     mainpad->cd();
955 buchmann 1.17
956     JZB_observed->Draw("e1");
957     JZB_ttbar_pred->Draw("histo,same");
958     JZBZpredClone->Draw("histo,same");
959     JZB_predicted->Draw("histo,same");
960     JZB_observed->Draw("e1,same");
961     legend->Draw();
962     if(is_data==data) DrawPrelim();
963     else DrawMCPrelim();
964    
965 buchmann 1.28 main_canvas->cd();
966     coverpad->Draw();
967     main_canvas->cd();
968     bottompad->Draw();
969 buchmann 1.17 bottompad->cd();
970     JZB_ratio->GetYaxis()->SetRangeUser(0,2);
971 buchmann 1.28
972     JZB_ratio->GetXaxis()->SetLabelSize(xstretchfactor*JZB_ratio->GetXaxis()->GetLabelSize());
973     JZB_ratio->GetYaxis()->SetLabelSize(xstretchfactor*JZB_ratio->GetYaxis()->GetLabelSize());
974     JZB_ratio->GetXaxis()->SetTitleSize(xstretchfactor*gStyle->GetTitleSize());
975     JZB_ratio->GetYaxis()->SetTitleSize(xstretchfactor*gStyle->GetTitleSize());
976    
977 buchmann 1.17 JZB_ratio->Draw("e1");
978 buchmann 1.28 // dsysenvelope->Draw();
979 buchmann 1.17 sysenvelope->Draw();
980     JZB_ratio->Draw("AXIS,same");
981     JZB_ratio->Draw("e1,same");
982     metl->Draw();
983    
984 buchmann 1.28 CompleteSave(main_canvas,SaveJZBname.str());
985 buchmann 1.17
986 buchmann 1.28 delete main_canvas;
987 buchmann 1.17 delete MET_observed;
988     delete MET_predicted;
989     //do NOT delete MET_Z_prediction (it's the return value)
990     delete MET_osof_pred;
991     delete MET_ossf_pred;
992     delete MET_ttbar_pred;
993    
994     delete JZB_observed;
995     delete JZB_predicted;
996     delete JZB_osof_pred;
997     delete JZB_ossf_pred;
998     delete JZB_Z_prediction;
999     delete JZB_ttbar_pred;
1000    
1001     return MET_Z_prediction;
1002     }
1003    
1004 buchmann 1.25 float extract_correction(string jzbvariable) {
1005     int position = (int)jzbvariable.find("[1]");
1006     if(position==-1) return 0.0;
1007     string correction=jzbvariable.substr(position+3,jzbvariable.length()-position-3);
1008     position = (int)correction.find("*");
1009     if(position==-1) return 0.0;
1010     correction=correction.substr(0,position);
1011     float correctionvalue=atof(correction.c_str());
1012     assert(correctionvalue<1&&correctionvalue>0);
1013     return correctionvalue;
1014     }
1015    
1016 buchmann 1.23 float Get_Met_Z_Prediction(TCut JetCut, float MetCut, int isdata, bool isDYonly, bool isAachen=false) {
1017 buchmann 1.17 dout << "Going to compute Z region prediction for a MET cut at " << MetCut << " GeV" << endl;
1018     // Steps:
1019 buchmann 1.25 // 1) Get peak
1020     // 2) use the peak and pt correction for sample splitting
1021     // and for MET distribution shifting
1022 buchmann 1.17 // 3) compute the estimate for MET>MetCut
1023    
1024     // do this for data if isdata==data, otherwise for MC (full closure if isDYonly==false, otherwise use only DY sample)
1025    
1026 buchmann 1.28 // Step 0 : If we're dealing with DY, we need to make sure PURW is off!
1027     // string bkpcutweight = (const char*) cutWeight;
1028     // if(isdata==mc && isDYonly) cutWeight=TCut("1.0");
1029    
1030    
1031 buchmann 1.25 // Step 1) Get peak
1032 buchmann 1.17 float MCPeakNoPtCorr=0,MCPeakErrorNoPtCorr=0,DataPeakNoPtCorr=0,DataPeakErrorNoPtCorr=0,MCSigma=0,DataSigma=0;
1033     stringstream resultsNoPtCorr;
1034     stringstream NoPtCorrdatajzb;
1035     stringstream NoPtCorrmcjzb;
1036    
1037 buchmann 1.24 if(isAachen) {
1038     //need to make sure that none of the typical basic cuts contain problematic selections!
1039     string Sleptoncut = (const char*) leptoncut;
1040     if((int)Sleptoncut.find("pt2>20")>-1) {
1041     write_error(__FUNCTION__,"You're trying to compute the Aachen estimate but are requiring pt2>20 ... please check your config.");
1042     assert((int)Sleptoncut.find("pt2>20")==-1);
1043     }
1044     } else {
1045     string Sleptoncut = (const char*) leptoncut;
1046     if((int)Sleptoncut.find("pt2>10")>-1) {
1047     write_error(__FUNCTION__,"You're trying to compute the ETH estimate but are requiring pt2>10 ... please check your config.");
1048     assert((int)Sleptoncut.find("pt2>10")==-1);
1049     }
1050     }
1051    
1052    
1053 buchmann 1.25 float Ptcorrection=0.0;
1054 buchmann 1.24
1055 buchmann 1.25 if(isdata==data) Ptcorrection=extract_correction(PlottingSetup::jzbvariabledata);
1056     else Ptcorrection=extract_correction(PlottingSetup::jzbvariablemc);
1057 buchmann 1.24
1058 buchmann 1.29 bool OverFlowStatus=addoverunderflowbins;
1059    
1060 buchmann 1.25 find_peaks(MCPeakNoPtCorr,MCPeakErrorNoPtCorr, DataPeakNoPtCorr,DataPeakErrorNoPtCorr,resultsNoPtCorr,true,NoPtCorrdatajzb,NoPtCorrmcjzb,(const char*) JetCut, true);
1061 buchmann 1.17
1062 buchmann 1.29 switch_overunderflow(OverFlowStatus);
1063 buchmann 1.28
1064 buchmann 1.25 float PeakPosition=0.0;
1065     string jzbvariable;
1066 buchmann 1.17 if(isdata==data) {
1067     PeakPosition=DataPeakNoPtCorr;
1068 buchmann 1.25 jzbvariable=jzbvariabledata;
1069 buchmann 1.17 dout << "Found peak in data at " << DataPeakNoPtCorr << " +/- " << DataPeakErrorNoPtCorr << " ; will use this result (" << PeakPosition << ")" << endl;
1070     } else {
1071     PeakPosition=MCPeakNoPtCorr;
1072 buchmann 1.25 jzbvariable=jzbvariablemc;
1073 buchmann 1.17 dout << "Found peak in mc at " << MCPeakNoPtCorr << " +/- " << MCPeakErrorNoPtCorr << " ; will use this result (" << PeakPosition << ")" << endl;
1074     }
1075    
1076     // Step 2: Use peak for sample splitting and MET shifting
1077 buchmann 1.25 string CorrectedMet="met[4]-"+any2string(Ptcorrection)+"*pt +"+any2string(abs(1.0*(PeakPosition)));
1078     if(2*(PeakPosition)<0) CorrectedMet="met[4]-"+any2string(Ptcorrection)+"*pt -"+any2string(abs(1.0*(PeakPosition)));
1079 buchmann 1.17
1080     stringstream sPositiveCut;
1081 buchmann 1.25 if(PeakPosition>0) sPositiveCut << "((" << jzbvariable << "-" << PeakPosition << ")>0)";
1082     else sPositiveCut << "( " << jzbvariable << "+" << abs(PeakPosition) << ")>0)";
1083 buchmann 1.17
1084     stringstream sNegativeCut;
1085 buchmann 1.25 if(PeakPosition<0) sNegativeCut << "((" << jzbvariable << "+" << abs(PeakPosition) << ")<0)";
1086     else sNegativeCut << "(( " << jzbvariable << "-" << abs(PeakPosition) << ")<0)";
1087 buchmann 1.17
1088     string ObservedMet="met[4]";
1089    
1090     stringstream JZBPosvar;
1091 buchmann 1.25 JZBPosvar<<jzbvariable;
1092     if(PeakPosition>0) JZBPosvar << "-" << PeakPosition;
1093     else JZBPosvar << "+" << abs(PeakPosition);
1094    
1095 buchmann 1.17 stringstream JZBNegvar;
1096 buchmann 1.25 JZBNegvar<<"-(" << jzbvariable;
1097     if(PeakPosition>0) JZBNegvar << "-" << PeakPosition << ")";
1098     else JZBNegvar << "+" << abs(PeakPosition) << ")";
1099    
1100 buchmann 1.17
1101     // Step 3: Compute estimate
1102 buchmann 1.23 TH1F *predicted = GetPredictedAndObservedMetShapes(JetCut, sPositiveCut.str(),sNegativeCut.str(),CorrectedMet,ObservedMet,JZBPosvar.str(),JZBNegvar.str(), MetCut, isdata, isDYonly, isAachen);
1103 buchmann 1.17 float ZregionZestimate=0;
1104     for(int ibin=1;ibin<=(int)predicted->GetNbinsX();ibin++) {
1105     if(predicted->GetBinLowEdge(ibin)+predicted->GetBinWidth(ibin)>MetCut) {
1106     ZregionZestimate+=2*(predicted->GetBinContent(ibin));
1107     }
1108     }
1109    
1110 buchmann 1.25 cout << " Z region estimate in MET>" << MetCut << " for this sample: " << ZregionZestimate << endl;
1111 buchmann 1.30 if(isdata==data) {
1112     MetPlotsSpace::Zestimate__data=ZregionZestimate;
1113     MetPlotsSpace::Zestimate__data_stat=2*TMath::Sqrt(ZregionZestimate/2);
1114     MetPlotsSpace::Zestimate__data_sys=ZregionZestimate*MetPlotsSpace::Zprediction_Uncertainty;
1115     }
1116     if(isdata==mc && isDYonly) {
1117     MetPlotsSpace::Zestimate__dy=ZregionZestimate;
1118     MetPlotsSpace::Zestimate__dy_stat=2*TMath::Sqrt(ZregionZestimate/2);
1119     MetPlotsSpace::Zestimate__dy_sys=ZregionZestimate*MetPlotsSpace::Zprediction_Uncertainty;
1120     }
1121     if(isdata==mc && !isDYonly) {
1122     MetPlotsSpace::Zestimate__mc=ZregionZestimate;
1123     MetPlotsSpace::Zestimate__mc_stat=2*TMath::Sqrt(ZregionZestimate/2);
1124     MetPlotsSpace::Zestimate__mc_sys=ZregionZestimate*MetPlotsSpace::Zprediction_Uncertainty;
1125     }
1126    
1127    
1128 buchmann 1.25
1129 buchmann 1.28 // if(isdata==mc && isDYonly) cutWeight=TCut(bkpcutweight.c_str());
1130    
1131 buchmann 1.17 return ZregionZestimate;
1132     }
1133    
1134 buchmann 1.30 void ExperimentalMetPrediction(bool QuickRun=false) {
1135 buchmann 1.17
1136 buchmann 1.28 switch_overunderflow(true);
1137 buchmann 1.23 bool isAachen=false;
1138 buchmann 1.24
1139 buchmann 1.35 bool HighPurityMode=false; // High Purity = |mll-91|<10 GeV , else <20
1140 buchmann 1.24
1141 buchmann 1.30 if(QuickRun) {
1142     isAachen=false;
1143 buchmann 1.32 HighPurityMode=true;
1144 buchmann 1.30 }
1145    
1146 buchmann 1.24 string restrmasscutbkp=(const char*) PlottingSetup::Restrmasscut;
1147    
1148 buchmann 1.32 if(HighPurityMode) PlottingSetup::Restrmasscut=TCut("abs(mll-91)<10");
1149     else PlottingSetup::Restrmasscut= TCut("abs(mll-91)<20");
1150 buchmann 1.24
1151    
1152    
1153 buchmann 1.23 cout << "Aachen mode (20/10, 2 jets) ? " << isAachen << endl;
1154 buchmann 1.24 cout << "High Purity mode? " << HighPurityMode << endl;
1155    
1156 buchmann 1.23
1157     if(isAachen) write_warning(__FUNCTION__,"Please don't forget to adapt the global lepton cut (to 20/10) for Aachen!");
1158 buchmann 1.24 stringstream snjets;
1159     if(isAachen) snjets << 2;
1160     else snjets << 3;
1161     float maxMET=100;
1162     if(isAachen) maxMET=150;
1163 buchmann 1.12
1164 buchmann 1.17 TCut nJetsSignal(PlottingSetup::basicqualitycut&&("pfJetGoodNum40>="+snjets.str()).c_str());
1165 buchmann 1.35 /*
1166 buchmann 1.17 cout << " ***** TESTING Z PREDICTION ***** " << endl;
1167     cout << "Notation (you can copy & paste this to evaluate it further)" << endl;
1168     cout << "Cut;Data;MC;DY;" << endl;
1169 buchmann 1.30 float DataEstimate = Get_Met_Z_Prediction(Restrmasscut&&nJetsSignal,maxMET, data, false, isAachen);
1170     float DYEstimate;
1171     if(!QuickRun) DYEstimate = Get_Met_Z_Prediction(Restrmasscut&&nJetsSignal,maxMET, mc, true, isAachen);
1172     float MCEstimate;
1173     if(!QuickRun) MCEstimate = Get_Met_Z_Prediction(Restrmasscut&&nJetsSignal,maxMET, mc, false, isAachen);
1174    
1175     if(QuickRun) return;
1176     cout << maxMET << ";" << DataEstimate << ";" << MCEstimate << ";" << DYEstimate << endl;
1177     dout << "Found estimate in data of " << DataEstimate << endl;
1178 buchmann 1.35 */
1179 buchmann 1.24 float Diff=20.0;
1180     if(HighPurityMode) Diff=10;
1181     TCut cut("mll>20&&pt1>20&&pt2>20");
1182 buchmann 1.35 if (isAachen) cut = TCut("mll>20&&pt1>20&&pt2>10&&pfTightHT>100");
1183 buchmann 1.28
1184     TCanvas *qcan = new TCanvas("qcan","qcan");
1185 buchmann 1.35 TH1F *zlineshape = allsamples.Draw("zlineshape","mll",int((91+25-18)*5),18,91+25,"m_{ll} (GeV)","events",cutOSSF&&TCut("pfJetGoodNum40==2")&&cut,data,PlottingSetup::luminosity); // bins of 0.2 GeV
1186     TH1F *Ozlineshape = allsamples.Draw("Ozlineshape","mll",int((91+25-18)*5),18,91+25,"m_{ll} (GeV)","events",cutOSOF&&TCut("pfJetGoodNum40==2")&&cut,data,PlottingSetup::luminosity); // bins of 0.2 GeV
1187 buchmann 1.28 zlineshape->Add(Ozlineshape,-1);
1188     delete qcan;
1189 buchmann 1.24 float a = (zlineshape->Integral(zlineshape->FindBin(20),zlineshape->FindBin(70)));
1190 buchmann 1.32 float b = (zlineshape->Integral(zlineshape->FindBin(91-Diff),zlineshape->FindBin(91+Diff)));
1191 buchmann 1.24 float r = a/b;
1192     float dr= (a/b)*TMath::Sqrt(1/a+1/b);
1193 buchmann 1.35 /*
1194 buchmann 1.24 float SysUncertainty = TMath::Sqrt(DataEstimate*DataEstimate*dr*dr + r*r*(DataEstimate*MetPlotsSpace::Zprediction_Uncertainty*DataEstimate*MetPlotsSpace::Zprediction_Uncertainty));
1195     float StatUncertainty = TMath::Sqrt(DataEstimate);
1196    
1197     cout << "Z estimate in peak : " << DataEstimate << " +/- " << DataEstimate*MetPlotsSpace::Zprediction_Uncertainty << " (sys) +/- " << TMath::Sqrt(2*DataEstimate) << " (stat) " << endl;
1198     cout << "Z ESTIMATE IN SR : " << DataEstimate*r << " +/- " << SysUncertainty << " (sys) +/- " << StatUncertainty << " (stat) " << endl;
1199 buchmann 1.35 // cout << endl;*/
1200 buchmann 1.24 cout << "r = " << r << " +/- " << dr << endl;
1201    
1202    
1203    
1204     delete zlineshape;
1205    
1206     PlottingSetup::Restrmasscut=TCut(restrmasscutbkp.c_str());
1207 buchmann 1.28 switch_overunderflow(false);
1208 buchmann 1.24
1209 buchmann 1.12 }
1210 buchmann 1.17