ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/cbrown/Development/Plotting/Modules/ZbTools.C
Revision: 1.2
Committed: Tue Jul 10 16:00:42 2012 UTC (12 years, 9 months ago) by buchmann
Content type: text/plain
Branch: MAIN
Changes since 1.1: +202 -23 lines
Log Message:
Z+b analysis update: now producing all plots from M.Musich's presentation

File Contents

# User Rev Content
1 buchmann 1.1 #include <iostream>
2     #include <vector>
3     #include <sys/stat.h>
4    
5     #include <TCut.h>
6     #include <TROOT.h>
7     #include <TCanvas.h>
8     #include <TMath.h>
9     #include <TColor.h>
10     #include <TPaveText.h>
11     #include <TRandom.h>
12     #include <TH1.h>
13     #include <TH2.h>
14     #include <TF1.h>
15     #include <TSQLResult.h>
16     #include <TProfile.h>
17     #include <TLegendEntry.h>
18    
19     using namespace std;
20    
21     using namespace PlottingSetup;
22    
23     namespace ZbData {
24     vector<float> data_over_mc;
25    
26     TCut ZplusBsel("pt1>20&&pt2>20&&mll>60&&mll<120");
27     TCut LeadingB("bTagProbCSVBP[0]>0.679");
28     TCut EtaB("abs(pfJetGoodEta[0])<1.3");
29     TCut PhiZcut("pfJetDphiZ[0]>2.7");
30    
31     }
32    
33     using namespace ZbData;
34    
35     void print_yield(TCut cut) {
36     TH1F *data = allsamples.Draw("data", "mll",1,0,10000000, "m_{ll} [GeV]", "events", cut,data,luminosity);
37     cout << " Yields for " << (const char*) cut << endl;
38     cout << " data: " << data->Integral() << endl;
39 buchmann 1.2 THStack mcm = allsamples.DrawStack("mc", "mll",1,0,10000000, "m_{ll} [GeV]", "events", cut,mc,luminosity);
40 buchmann 1.1 int nhists = mcm.GetHists()->GetEntries();
41     cout << "Going to loop over " << nhists << " histograms " << endl;
42     TFile *test = new TFile("test.root","RECREATE");
43     mcm.Write();
44     test->Close();
45     for (int i = 0; i < nhists; i++) {
46     TH1* hist = static_cast<TH1*>(mcm.GetHists()->At(i));
47     cout << " " << hist->GetName() << ": " << hist->Integral() << endl;
48     }
49     cout << endl << endl;
50     delete data;
51     }
52    
53 buchmann 1.2 void draw_kin_variable(string variable, TCut cut, int nbins, float min, float max, string xlabel, string saveas, bool dologscale, string speciallabel="") {
54 buchmann 1.1
55     TCanvas *ckin = new TCanvas("ckin","kin variable canvas");
56     ckin->SetLogy(dologscale);
57     TH1F *datahisto = allsamples.Draw("datahisto",variable,nbins,min,max, xlabel, "events",cut,data,luminosity);
58     datahisto->SetMarkerSize(DataMarkerSize);
59     THStack mcstack = allsamples.DrawStack("mcstack",variable,nbins,min,max, xlabel, "events",cut,mc,luminosity);
60 buchmann 1.2 if(dologscale) datahisto->SetMaximum(3*datahisto->GetMaximum());
61     else datahisto->SetMaximum(1.5*datahisto->GetMaximum());
62     TText *speciall;
63     if(speciallabel!="") {
64     speciall = new TText(0.2,0.8,speciallabel.c_str());
65     speciall->Draw();
66     }
67 buchmann 1.1 datahisto->Draw("e1");
68     ckin->Update();
69     mcstack.Draw("same");
70     datahisto->Draw("same,e1");
71     TLegend *kinleg = allsamples.allbglegend();
72     kinleg->Draw();
73    
74     TPad *kinpad = new TPad("kinpad","kinpad",0,0,1,1);
75     kinpad->cd();
76     kinpad->SetLogy(dologscale);
77     datahisto->Draw("e1");
78     mcstack.Draw("same");
79     datahisto->Draw("same,e1");
80     datahisto->Draw("same,axis");
81     kinleg->Draw();
82     DrawPrelim();
83     saveas="kin/"+saveas;
84     save_with_ratio(datahisto,mcstack,kinpad->cd(),saveas);
85    
86    
87    
88    
89     ZbData::data_over_mc.push_back(datahisto->Integral()/CollapseStack(mcstack)->Integral());
90    
91     datahisto->Delete();
92     delete ckin;
93     }
94    
95     Value getfrom2Dmap(TH2F *map, int ixbin, int iybin) {
96     Value sum(0,0);
97     for(int iy=1;iy<=iybin;iy++) {
98     sum=sum+Value(map->GetBinContent(ixbin,iy),map->GetBinError(ixbin,iy));
99     }
100     return sum;
101     }
102    
103    
104     void data_mc_agreement_2d() {
105    
106     vector<int> emptyvector;
107     float binX[5]={10,50,100,200,1000};
108     float binY[5]={0.0,0.1,0.15,0.2,0.3};
109     vector<float>binningx;
110     vector<float>binningy;
111     binningx.assign(binX,binX+5);
112     binningy.assign(binY,binY+5);
113    
114     TCut cut(ZplusBsel&&LeadingB&&EtaB&&PhiZcut);
115    
116     TCanvas *agrcan = new TCanvas("agrcan","agrcan");
117    
118     cout << "Drawing 2d for data " << endl;
119     TH2F *alphazptplane_data = allsamples.Draw("alphazptplane_data","alpha:pt",binningx,binningy,"Z p_{T}","#alpha",cut,data,luminosity,emptyvector,false);
120     cout << "Drawing 2d for MC " << endl;
121     TH2F *alphazptplane_mc = allsamples.Draw("alphazptplane_mc", "alpha:pt",binningx,binningy,"Z p_{T}","#alpha",cut,mc, luminosity,emptyvector,false);
122    
123    
124     TFile *f = new TFile("test.root","RECREATE");
125     alphazptplane_data->Write();
126     f->Close();
127    
128     for(int ix=1;ix<binningx.size();ix++) {
129     cout << "pt in " << binningx[ix-1] << " , " << binningx[ix] << endl;
130     TGraphErrors gr;
131     for(int iy=1;iy<binningy.size();iy++) {
132     Value mc = Value((getfrom2Dmap(alphazptplane_mc,ix,iy)).getValue(),0.0);// not taking MC statistics into account
133     Value data = getfrom2Dmap(alphazptplane_data,ix,iy);
134     cout << " For an alpha cut at " << binningy[iy] << " we get : mc (" << mc << ") , data (" << data << ") --> data/mc = " << data/mc << endl;
135     gr.SetPoint(iy-1,binningy[iy],(data/mc).getValue());
136     gr.SetPointError(iy-1,0,(data/mc).getError());
137     }
138    
139     gr.Draw("AP*");
140     gr.GetYaxis()->SetTitle("data / mc");
141     gr.GetYaxis()->CenterTitle();
142     gr.GetXaxis()->SetTitle("#alpha");
143     gr.GetXaxis()->CenterTitle();
144     gr.GetXaxis()->SetRangeUser(0,0.4);
145     gr.GetYaxis()->SetRangeUser(0.5,1.5);
146     gr.Draw("AP*");
147     DrawPrelim();
148     stringstream ptbin;
149     ptbin << binningx[ix-1] << " < p_{T}^{Z} < " << binningx[ix] << " GeV";
150     TText *ptt = write_title(ptbin.str());
151     ptt->SetX(0.7);
152     ptt->SetY(0.7);
153     ptt->SetTextSize(0.03);
154     ptt->Draw();
155     TLine *l = new TLine(0.0,1.0,0.4,1.0);
156     l->SetLineStyle(2);
157     l->SetLineColor(kBlue);
158     l->Draw();
159    
160     stringstream saveas;
161     saveas << "Zplusb_data_over_mc___" << binningx[ix-1] << "to" << binningx[ix];
162     CompleteSave(agrcan,saveas.str());
163     // agrcan->SaveAs("tester.png");
164     }
165    
166     delete agrcan;
167    
168     }
169    
170     void print_all_b_yields() {
171     cout << "Basic selection with a b jet" << endl;
172     print_yield(ZplusBsel&&"pfJetGoodNumBtag>0");
173 buchmann 1.2 cout << "Leading jet is a b-jet " << endl;
174 buchmann 1.1 print_yield(ZplusBsel&&LeadingB);
175     cout << "|#eta_{b}|<1.3 " << endl;
176     print_yield(ZplusBsel&&LeadingB&&EtaB);
177     cout << "|#delta#phi(b<Z)|>2.7" << endl;
178     print_yield(ZplusBsel&&LeadingB&&EtaB&&PhiZcut);
179 buchmann 1.2 cout << "10<ptZ<1000 GeV" << endl;
180 buchmann 1.1 print_yield(ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000");
181     cout << "#alpha < 0.3" << endl;
182     print_yield(ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"alpha<0.3");
183     cout << "#alpha < 0.2" << endl;
184     print_yield(ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"alpha<0.2");
185     cout << "#alpha < 0.15" << endl;
186     print_yield(ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"alpha<0.15");
187     cout << "#alpha < 0.1" << endl;
188     print_yield(ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"alpha<0.1");
189     cout << "#alpha < 0.05" << endl;
190     print_yield(ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"alpha<0.05");
191     }
192    
193     void draw_Zb_kin_vars() {
194     draw_kin_variable("pfBJetDphiZ",ZplusBsel&&LeadingB,100,0,3.2,"#delta#phi (Z,b lead)","DeltaPhiZBlead",1);
195     draw_kin_variable("pt",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000",50,0,200,"Z p_{T}","Zpt",1);
196     draw_kin_variable("pfJetGoodPt[1]",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000",50,0,200,"Sub-Leading Jet Pt","Jet2Pt",1);
197     draw_kin_variable("pfJetGoodPtlBtag[0]",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000",50,0,200,"Leading Jet Pt (B)","LeadingJetPt",1);
198     draw_kin_variable("alpha",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000",20,0,2,"#alpha","alpha",1);
199     }
200    
201     void draw_mpf_vars() {
202     ZbData::data_over_mc.clear();
203 buchmann 1.2 draw_kin_variable("mpf",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"alpha<0.3",20,0,2,"Z+b MPF","mpf_alpha_smaller_0p3",0,"#alpha<0.3, 10 GeV < p_{T}^{Z} < 1000 GeV");
204     draw_kin_variable("mpf",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"alpha<0.2",20,0,2,"Z+b MPF","mpf_alpha_smaller_0p2",0,"#alpha<0.2, 10 GeV < p_{T}^{Z} < 1000 GeV");
205     draw_kin_variable("mpf",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"alpha<0.15",20,0,2,"Z+b MPF","mpf_alpha_smaller_0p15",0,"#alpha<0.15, 10 GeV < p_{T}^{Z} < 1000 GeV");
206     draw_kin_variable("mpf",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"alpha<0.1",20,0,2,"Z+b MPF","mpf_alpha_smaller_0p1",0,"#alpha<0.1, 10 GeV < p_{T}^{Z} < 1000 GeV");
207     draw_kin_variable("mpf",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"alpha<0.05",20,0,2,"Z+b MPF","mpf_alpha_smaller_0p05",0,"#alpha<0.05, 10 GeV < p_{T}^{Z} < 1000 GeV");
208    
209     draw_kin_variable("pfJetGoodPt[0]/pt",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"alpha<0.3",20,0,2,"p_{T} b-jet / p_{T} Z","ptb_over_ptz___alpha_smaller_0p3",0,"#alpha<0.3, 10 GeV < p_{T}^{Z} < 1000 GeV");
210     draw_kin_variable("pfJetGoodPt[0]/pt",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"alpha<0.2",20,0,2,"p_{T} b-jet / p_{T} Z","ptb_over_ptz___alpha_smaller_0p2",0,"#alpha<0.2, 10 GeV < p_{T}^{Z} < 1000 GeV");
211     draw_kin_variable("pfJetGoodPt[0]/pt",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"alpha<0.15",20,0,2,"p_{T} b-jet / p_{T} Z","ptb_over_ptz___alpha_smaller_0p15",0,"#alpha<0.15, 10 GeV < p_{T}^{Z} < 1000 GeV");
212     draw_kin_variable("pfJetGoodPt[0]/pt",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"alpha<0.1",20,0,2,"p_{T} b-jet / p_{T} Z","ptb_over_ptz___alpha_smaller_0p1",0,"#alpha<0.1, 10 GeV < p_{T}^{Z} < 1000 GeV");
213     draw_kin_variable("pfJetGoodPt[0]/pt",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"alpha<0.05",20,0,2,"p_{T} b-jet / p_{T} Z","ptb_over_ptz__mpf_alpha_smaller_0p05",0,"#alpha<0.05, 10 GeV < p_{T}^{Z} < 1000 GeV");
214 buchmann 1.1 }
215    
216 buchmann 1.2 Value get_Zb_data_over_mc(string variable, TCut cut) {
217     // write_warning(__FUNCTION__,"Debugging this function, therefore always returning 3.1415+/-1.010101"); return Value(3.1415,1.010101);
218     TH1F *hdata = allsamples.Draw("hdata",variable,1000,0,30, "#alpha [GeV]", "events", cut,data,luminosity);
219     TH1F *hmc = allsamples.Draw("hmc" , variable,1000,0,30, "#alpha [GeV]", "events", cut,mc, luminosity);
220    
221     float a=hdata->GetMean();
222     float b=hmc->GetMean();
223     float da=hdata->GetMeanError();
224     float db=hmc->GetMeanError();
225     float factor = a / b;
226     float error = TMath::Sqrt( (1/(b*b)) * (da*da) + ((a*a)/(b*b))*db*db);
227     delete hdata;
228     delete hmc;
229     return Value(factor,error);
230     }
231    
232    
233    
234     void new_data_mc_agreement_2d() {
235     TCut basecut(ZplusBsel&&LeadingB&&EtaB&&PhiZcut);
236    
237     const int nalphacuts=5;
238     float alphacuts[nalphacuts] = {0.1,0.15,0.2,0.25,0.3};
239     const int nptcuts=5;
240     float ptcuts[nptcuts]={10,50,100,200,1000};
241    
242     float MPF_Results[nalphacuts][nptcuts];
243     float MPF_Errors[nalphacuts][nptcuts];
244     float RABS_Results[nalphacuts][nptcuts];
245     float RABS_Errors[nalphacuts][nptcuts];
246    
247    
248     for(int ia=0;ia<nalphacuts;ia++) {
249     for(int ipt=0;ipt<nptcuts-1;ipt++) {
250     stringstream specialcut;
251     specialcut << "((pt>" << ptcuts[ipt] << " && pt< " << ptcuts[ipt+1] << ") && (alpha<" << alphacuts[ia] << "))";
252     Value MPF_data_over_mc = get_Zb_data_over_mc("mpf",TCut(basecut && specialcut.str().c_str()));
253     Value RABS_data_over_mc = get_Zb_data_over_mc("pfJetGoodPt[0]/pt",TCut(basecut && specialcut.str().c_str()));
254    
255     MPF_Results[ia][ipt]=MPF_data_over_mc.getValue();
256     MPF_Errors[ia][ipt]=MPF_data_over_mc.getError();
257    
258     RABS_Results[ia][ipt]=RABS_data_over_mc.getValue();
259     RABS_Errors[ia][ipt]=RABS_data_over_mc.getError();
260    
261     cout << "alpha cut at " << alphacuts[ia] << " and pt range " << ptcuts[ipt] << " , " << ptcuts[ipt+1] << " \t : \t" << specialcut.str() << " \t \t --> " << MPF_data_over_mc << " (MPF) " << RABS_data_over_mc << " (RABS)" << endl;
262     }
263     }
264    
265     float MPF_ExtraPolatedResults[nptcuts];
266     float RABS_ExtraPolatedResults[nptcuts];
267    
268     TCanvas *can = new TCanvas("can");
269    
270     TGraphErrors *MPF_FinalGraph = new TGraphErrors();
271     TGraphErrors *RABS_FinalGraph = new TGraphErrors();
272    
273     for(int ipt=0;ipt<nptcuts-1;ipt++) {
274    
275     stringstream filename;
276     filename << "Extrapolation/Zplusb_data_over_mc___" << ptcuts[ipt] << "_to_" << ptcuts[ipt+1];
277     cout << "Going to create histo for " << filename.str() << endl;
278     TGraphErrors *mpf_gr =new TGraphErrors();
279     TGraphErrors *rabs_gr =new TGraphErrors();
280     for(int ia=0;ia<nalphacuts;ia++) {
281     mpf_gr->SetPoint(ia,alphacuts[ia],MPF_Results[ia][ipt]);
282     mpf_gr->SetPointError(ia,0,MPF_Errors[ia][ipt]);
283     rabs_gr->SetPoint(ia,alphacuts[ia],RABS_Results[ia][ipt]);
284     rabs_gr->SetPointError(ia,0,RABS_Errors[ia][ipt]);
285     }
286    
287     mpf_gr->SetMarkerStyle(21);
288     mpf_gr->SetMarkerColor(kBlue);
289     mpf_gr->Draw("AP");
290     mpf_gr->Fit("pol1");
291     TF1 *mpf_pol=(TF1*)mpf_gr->GetFunction("pol1");
292     float mpf_a=mpf_pol->GetParameter(0);
293     float mpf_b=mpf_pol->GetParameter(1);
294     MPF_FinalGraph->SetPoint(ipt,0.5*(ptcuts[ipt]+ptcuts[ipt+1]),mpf_a);
295     MPF_FinalGraph->SetPointError(ipt,0,mpf_pol->GetParError(0));
296    
297     MPF_ExtraPolatedResults[ipt]=mpf_a;
298    
299     stringstream mpf_info;
300     mpf_info << "#splitline{#splitline{" << ptcuts[ipt] << " GeV < p_{T}^{Z} < " << ptcuts[ipt+1] << " GeV }{ (MPF) }}{Extrapolated value: " << std::setprecision(3) << mpf_a << "}";
301     TText *mpf_mark = write_title(mpf_info.str());
302     mpf_mark->SetX(0.75);
303     mpf_mark->SetTextSize(0.03);
304     mpf_mark->SetY(0.75);
305     mpf_mark->Draw();
306    
307     string filenamebkp=filename.str();
308     filename << "__MPF";
309     CompleteSave(can,filename.str());
310     cout << "MPF : " << mpf_a << " + " << mpf_b << " * alpha " << endl;
311    
312     filename.str("");
313    
314     rabs_gr->SetMarkerStyle(21);
315     rabs_gr->SetMarkerColor(kBlue);
316     rabs_gr->Draw("AP");
317     rabs_gr->Fit("pol1");
318     TF1 *rabs_pol=(TF1*)rabs_gr->GetFunction("pol1");
319     float rabs_a=rabs_pol->GetParameter(0);
320     float rabs_b=rabs_pol->GetParameter(1);
321     RABS_FinalGraph->SetPoint(ipt,0.5*(ptcuts[ipt]+ptcuts[ipt+1]),rabs_a);
322     RABS_FinalGraph->SetPointError(ipt,0,rabs_pol->GetParError(0));
323     cout << "!+!+!+!+!+!+!+!+!+ Just added a point to the final plot : " << rabs_a << " +/- " << rabs_pol->GetParError(0) << endl;
324    
325    
326     stringstream rabs_info;
327     rabs_info << "#splitline{#splitline{" << ptcuts[ipt] << " GeV < p_{T}^{Z} < " << ptcuts[ipt+1] << " GeV }{ (R_{abs}) }}{Extrapolated value: " << std::setprecision(3) << rabs_a << "}";
328     TText *rabs_mark = write_title(rabs_info.str());
329     rabs_mark->SetX(0.75);
330     rabs_mark->SetTextSize(0.03);
331     rabs_mark->SetY(0.75);
332     rabs_mark->Draw();
333    
334     RABS_ExtraPolatedResults[ipt]=rabs_a;
335     filename << filenamebkp << "__RABS";
336     CompleteSave(can,filename.str());
337     cout << "RABS : " << rabs_a << " + " << rabs_b << " * alpha " << endl;
338     }
339    
340     MPF_FinalGraph->SetMarkerStyle(21);
341     MPF_FinalGraph->SetMarkerColor(kRed);
342     MPF_FinalGraph->Fit("pol0");
343     TF1 *mpf_pol0=(TF1*)MPF_FinalGraph->GetFunction("pol0");
344     MPF_FinalGraph->Draw("AP");
345     stringstream mpf_result;
346     mpf_result << "C_{abs}= " << mpf_pol0->GetParameter(0) << " #pm " << mpf_pol0->GetParError(0);
347     TText *rmpf_final_mark = write_title(mpf_result.str());
348     rmpf_final_mark->SetX(0.75);
349     rmpf_final_mark->SetY(0.75);
350     rmpf_final_mark->SetTextSize(0.03);
351     rmpf_final_mark->Draw();
352    
353     stringstream filename2;
354     filename2 << "Extrapolation/FINAL_RESULT_MPF";
355     CompleteSave(can,filename2.str());
356    
357    
358     RABS_FinalGraph->SetMarkerStyle(21);
359     RABS_FinalGraph->SetMarkerStyle(21);
360     RABS_FinalGraph->SetMarkerColor(kRed);
361     RABS_FinalGraph->Fit("pol0");
362     RABS_FinalGraph->Draw("AP");
363     TF1 *rabs_pol0=(TF1*)RABS_FinalGraph->GetFunction("pol0");
364     stringstream rabs_result;
365     rabs_result << "C_{abs}= " << rabs_pol0->GetParameter(0) << " #pm " << rabs_pol0->GetParError(0);
366     TText *rabs_final_mark = write_title(rabs_result.str());
367     rabs_final_mark->SetX(0.75);
368     rabs_final_mark->SetY(0.75);
369     rabs_final_mark->SetTextSize(0.03);
370     rabs_final_mark->Draw();
371    
372     filename2.str("");
373     filename2 << "Extrapolation/FINAL_RESULT_RABS";
374     CompleteSave(can,filename2.str());
375    
376    
377     delete can;
378    
379    
380     }
381    
382    
383    
384 buchmann 1.1 void do_basic_ZB_analysis() {
385 buchmann 1.2 TCanvas *zbcanvas = new TCanvas("zbcanvas","zbcanvas");
386     write_warning(__FUNCTION__,"TESTING NEW DATA MC AGREEMENT");
387 buchmann 1.1
388 buchmann 1.2 /*
389     write_warning(__FUNCTION__,"Not doing b yields (already implemented)");
390     // print_all_b_yields();
391     write_warning(__FUNCTION__,"Not drawing mpf variables (already implemented)");
392     // draw_mpf_vars();
393     write_warning(__FUNCTION__,"Not drawing mpf variables (already implemented)");
394     draw_Zb_kin_vars();
395 buchmann 1.1
396     data_mc_agreement_2d();
397 buchmann 1.2 */
398     new_data_mc_agreement_2d();
399 buchmann 1.1
400 buchmann 1.2 delete zbcanvas;
401 buchmann 1.1 }