PWG0/dNdEta/drawPlots.C

   1 /* $Id$ */
   2
   3 #if !defined(__CINT__) || defined(__MAKECINT__)
   4
   5 #include "AliPWG0Helper.h"
   6 #include "dNdEtaAnalysis.h"
   7 #include "AlidNdEtaCorrection.h"
   8
   9 #include <TCanvas.h>
  10 #include <TFile.h>
  11 #include <TH1.h>
  12 #include <TH2F.h>
  13 #include <TH3F.h>
  14 #include <TLine.h>
  15 #include <TSystem.h>
  16
  17 #endif
  18
  19 Int_t gMax = 5;
  20
  21 extern TSystem* gSystem;
  22
  23 void loadlibs()
  24 {
  25   gSystem->Load("libANALYSIS");
  26   gSystem->Load("libPWG0base");
  27 }
  28
  29 void SetRanges(TAxis* axis)
  30 {
  31   if (strcmp(axis->GetTitle(), "#eta") == 0)
  32     axis->SetRangeUser(-1.7999, 1.7999);
  33   if (strcmp(axis->GetTitle(), "p_{T} [GeV/c]") == 0)
  34     axis->SetRangeUser(0, 9.9999);
  35   if (strcmp(axis->GetTitle(), "vtx z [cm]") == 0)
  36     axis->SetRangeUser(-15, 14.9999);
  37   if (strcmp(axis->GetTitle(), "Ntracks") == 0)
  38     axis->SetRangeUser(0, 99.9999);
  39 }
  40
  41 void SetRanges(TH1* hist)
  42 {
  43   SetRanges(hist->GetXaxis());
  44   SetRanges(hist->GetYaxis());
  45   SetRanges(hist->GetZaxis());
  46 }
  47
  48
  49 void Prepare3DPlot(TH3* hist)
  50 {
  51   hist->GetXaxis()->SetTitleOffset(1.5);
  52   hist->GetYaxis()->SetTitleOffset(1.5);
  53   hist->GetZaxis()->SetTitleOffset(1.5);
  54
  55   hist->SetStats(kFALSE);
  56 }
  57
  58 void Prepare2DPlot(TH2* hist)
  59 {
  60   hist->SetStats(kFALSE);
  61   hist->GetYaxis()->SetTitleOffset(1.4);
  62
  63   hist->SetMinimum(0);
  64   hist->SetMaximum(gMax);
  65
  66   SetRanges(hist);
  67 }
  68
  69 void Prepare1DPlot(TH1* hist)
  70 {
  71   hist->SetLineWidth(2);
  72   hist->SetStats(kFALSE);
  73
  74   hist->GetXaxis()->SetLabelOffset(0.02);
  75   hist->GetXaxis()->SetTitleOffset(1.3);
  76   hist->GetYaxis()->SetTitleOffset(1.3);
  77
  78   SetRanges(hist);
  79 }
  80
  81 void InitPad()
  82 {
  83   gPad->Range(0, 0, 1, 1);
  84   gPad->SetLeftMargin(0.15);
  85   //gPad->SetRightMargin(0.05);
  86   //gPad->SetTopMargin(0.13);
  87   gPad->SetBottomMargin(0.12);
  88
  89   gPad->SetGridx();
  90   gPad->SetGridy();
  91 }
  92
  93 void InitPadCOLZ()
  94 {
  95   gPad->Range(0, 0, 1, 1);
  96   gPad->SetRightMargin(0.15);
  97   gPad->SetLeftMargin(0.12);
  98
  99   gPad->SetGridx();
 100   gPad->SetGridy();
 101 }
 102
 103 // --- end of helpers --- begin functions ---
 104
 105 void DrawOverview(const char* fileName = "correction_map.root", const char* dirName = "dndeta_correction")
 106 {
 107   loadlibs();
 108   if (!TFile::Open(fileName))
 109     return;
 110
 111   AlidNdEtaCorrection* dNdEtaCorrection = new AlidNdEtaCorrection(dirName, dirName);
 112   if (!dNdEtaCorrection->LoadHistograms())
 113     return;
 114
 115   dNdEtaCorrection->Finish();
 116
 117   dNdEtaCorrection->DrawOverview();
 118 }
 119
 120 void ComparedNdEta(const char* ESDfolder = "dndeta", const char* MCfolder = "dndeta", const char* esdFile = "analysis_esd.root", const char* mcFile = "analysis_mc.root")
 121 {
 122   gSystem->Load("libPWG0base");
 123
 124   TFile::Open(esdFile);
 125   dNdEtaAnalysis* fdNdEtaAnalysisESD = new dNdEtaAnalysis(ESDfolder, ESDfolder);
 126   fdNdEtaAnalysisESD->LoadHistograms();
 127
 128   TFile::Open(mcFile);
 129   dNdEtaAnalysis* fdNdEtaAnalysisMC = new dNdEtaAnalysis(MCfolder, MCfolder);
 130   fdNdEtaAnalysisMC->LoadHistograms();
 131   //fdNdEtaAnalysisMC->Finish(0, 0.3, AlidNdEtaCorrection::kNone);
 132
 133   for (Int_t i=0; i<dNdEtaAnalysis::kVertexBinning; ++i)
 134     fdNdEtaAnalysisESD->GetdNdEtaPtCutOffCorrectedHistogram(i)->Divide(fdNdEtaAnalysisMC->GetdNdEtaPtCutOffCorrectedHistogram(i));
 135
 136   fdNdEtaAnalysisESD->DrawHistograms();
 137 }
 138
 139 void CompareVertexDist(Int_t plot = 0, const char* correctionMapFile = "correction_map.root", const char* correctionMapFolder = "dndeta_correction")
 140 {
 141   gSystem->Load("libPWG0base");
 142
 143   const char* ESDfolder = 0;
 144
 145   if (plot == 0) // all
 146     ESDfolder = "dndeta";
 147   else if (plot == 1) // mb
 148     ESDfolder = "dndeta_mb";
 149   else if (plot == 2) // mb vtx
 150     ESDfolder = "dndeta_mbvtx";
 151
 152   TFile::Open("analysis_esd.root");
 153   dNdEtaAnalysis* fdNdEtaAnalysisESD = new dNdEtaAnalysis(ESDfolder, ESDfolder);
 154   fdNdEtaAnalysisESD->LoadHistograms();
 155
 156   AlidNdEtaCorrection* dNdEtaCorrection = new AlidNdEtaCorrection(correctionMapFolder, correctionMapFolder);
 157   dNdEtaCorrection->LoadHistograms(correctionMapFile, correctionMapFolder);
 158
 159   TH2F* hist = 0;
 160
 161   if (plot == 0) // all
 162     hist = dNdEtaCorrection->GetTriggerBiasCorrection()->GetGeneratedHistogram();
 163   else if (plot == 1) // mb
 164     hist = dNdEtaCorrection->GetTriggerBiasCorrection()->GetMeasuredHistogram();
 165   else if (plot == 2) // mb vtx
 166     hist = dNdEtaCorrection->GetVertexRecoCorrection()->GetMeasuredHistogram();
 167
 168   TH1* proj = hist->ProjectionX();
 169
 170   TH1* vertex = fdNdEtaAnalysisESD->GetVtxHistogram();
 171   for (Int_t i=1; i<=vertex->GetNbinsX(); ++i)
 172   {
 173     Float_t value = proj->GetBinContent(proj->FindBin(vertex->GetBinCenter(i)));
 174     if (value != 0)
 175     {
 176       printf("vtx = %f, esd = %f, corr = %f, ratio = %f\n", vertex->GetBinCenter(i), vertex->GetBinContent(i), value, vertex->GetBinContent(i) / value);
 177       vertex->SetBinContent(i, vertex->GetBinContent(i) / value);
 178     }
 179   }
 180
 181   new TCanvas;
 182   vertex->DrawCopy();
 183 }
 184
 185 void CompareTrack2ParticleWithAnalysisData(const char* correctionMapFile = "correction_map.root", const char* correctionMapFolder = "dndeta_correction")
 186 {
 187   gSystem->Load("libPWG0base");
 188
 189   TFile::Open("analysis_esd.root");
 190   dNdEtaAnalysis* fdNdEtaAnalysisESD = new dNdEtaAnalysis("dndeta_mbvtx", "dndeta_mbvtx");
 191   fdNdEtaAnalysisESD->LoadHistograms();
 192
 193   AlidNdEtaCorrection* dNdEtaCorrection = new AlidNdEtaCorrection(correctionMapFolder, correctionMapFolder);
 194   dNdEtaCorrection->LoadHistograms(correctionMapFile, correctionMapFolder);
 195
 196   //TH1* histESD = fdNdEtaAnalysisESD->GetUncorrectedHistogram();
 197   //TH1* histCorr = dNdEtaCorrection->GetTrack2ParticleCorrection()->GetMeasuredHistogram();
 198
 199   TH1* histESD = fdNdEtaAnalysisESD->GetHistogram();
 200   TH1* histCorr = dNdEtaCorrection->GetTrack2ParticleCorrection()->GetGeneratedHistogram();
 201
 202   TH1F* diff = new TH1F("diff", "diff", 100, 0.98, 1.02);
 203
 204   new TCanvas;
 205   histESD->Draw();
 206
 207   new TCanvas;
 208   histCorr->Draw();
 209
 210   for (Int_t x=1; x<=histESD->GetNbinsX(); ++x)
 211     for (Int_t y=1; y<=histESD->GetNbinsY(); ++y)
 212       for (Int_t z=1; z<=histESD->GetNbinsZ(); ++z)
 213       {
 214         Float_t value1 = histESD->GetBinContent(x, y, z);
 215         Float_t value2 = histCorr->GetBinContent(histCorr->FindBin(histESD->GetXaxis()->GetBinCenter(x), histESD->GetYaxis()->GetBinCenter(y), histESD->GetZaxis()->GetBinCenter(z)));
 216
 217         if (value2 > 0 && value1 > 0)
 218         {
 219           printf("%f %f %f\n", value1, value2, value1 / value2);
 220           diff->Fill(value1 / value2);
 221         }
 222       }
 223
 224   new TCanvas;
 225   diff->Draw();
 226 }
 227
 228 Double_t PrintIntegratedDeviation(TH1* histMC, TH1* histESD, const char* desc = "")
 229 {
 230   Double_t avgMC = 0;
 231   Double_t avgESD = 0;
 232   for (Int_t bin = histMC->FindBin(-0.7999); bin <= histMC->FindBin(0.7999); bin++)
 233   {
 234     avgMC += histMC->GetBinContent(bin);
 235     avgESD += histESD->GetBinContent(bin);
 236   }
 237   Int_t nBins = histMC->FindBin(0.7999) - histMC->FindBin(-0.7999) + 1;
 238
 239   avgMC /= nBins;
 240   avgESD /= nBins;
 241
 242   // deviation when integrate in |eta| < 0.8 between mc and esd
 243   Double_t diffFullRange = (avgMC - avgESD) / avgMC;
 244
 245   Printf("%s: Integrated deviation in |eta| < 0.8 is %.2f %%", desc, diffFullRange * 100);
 246
 247   return diffFullRange;
 248 }
 249
 250 void dNdEta(Bool_t onlyESD = kFALSE, Bool_t save = kTRUE)
 251 {
 252   TFile* file = TFile::Open("analysis_esd.root");
 253   TH1* histESD = (TH1*) file->Get("dndeta/dNdEta_corrected");
 254   TH1* histESDNoPt = (TH1*) file->Get("dndeta/dNdEta");
 255   TH1* histESDMB = (TH1*) file->Get("dndetaTr/dNdEta_corrected");
 256   TH1* histESDMBNoPt = (TH1*) file->Get("dndetaTr/dNdEta");
 257   TH1* histESDMBVtx = (TH1*) file->Get("dndetaTrVtx/dNdEta_corrected");
 258   TH1* histESDMBVtxNoPt = (TH1*) file->Get("dndetaTrVtx/dNdEta");
 259   TH1* histESDMBTracksNoPt = (TH1*) file->Get("dndetaTracks/dNdEta");
 260
 261   Prepare1DPlot(histESD);
 262   Prepare1DPlot(histESDMB);
 263   Prepare1DPlot(histESDMBVtx);
 264
 265   Prepare1DPlot(histESDNoPt);
 266   Prepare1DPlot(histESDMBNoPt);
 267   Prepare1DPlot(histESDMBVtxNoPt);
 268   Prepare1DPlot(histESDMBTracksNoPt);
 269
 270   histESD->SetLineWidth(0);
 271   histESDMB->SetLineWidth(0);
 272   histESDMBVtx->SetLineWidth(0);
 273
 274   histESDNoPt->SetLineWidth(0);
 275   histESDMBNoPt->SetLineWidth(0);
 276   histESDMBVtxNoPt->SetLineWidth(0);
 277
 278   histESD->SetMarkerColor(1);
 279   histESDMB->SetMarkerColor(2);
 280   histESDMBVtx->SetMarkerColor(3);
 281
 282   histESD->SetLineColor(1);
 283   histESDMB->SetLineColor(2);
 284   histESDMBVtx->SetLineColor(3);
 285
 286   histESDNoPt->SetMarkerColor(1);
 287   histESDMBNoPt->SetMarkerColor(2);
 288   histESDMBVtxNoPt->SetMarkerColor(3);
 289   histESDMBTracksNoPt->SetMarkerColor(4);
 290
 291   histESD->SetMarkerStyle(20);
 292   histESDMB->SetMarkerStyle(21);
 293   histESDMBVtx->SetMarkerStyle(22);
 294
 295   histESDNoPt->SetMarkerStyle(20);
 296   histESDMBNoPt->SetMarkerStyle(21);
 297   histESDMBVtxNoPt->SetMarkerStyle(22);
 298   histESDMBTracksNoPt->SetMarkerStyle(23);
 299
 300   TH2F* dummy = new TH2F("dummy", "", 100, -1.5, 1.5, 1000, 0, histESDMBVtx->GetMaximum() * 1.1);
 301   Prepare1DPlot(dummy);
 302   dummy->SetStats(kFALSE);
 303   dummy->SetXTitle("#eta");
 304   dummy->SetYTitle("dN_{ch}/d#eta");
 305   dummy->GetYaxis()->SetTitleOffset(1);
 306
 307   Float_t etaLimit = 1.2999;
 308
 309   histESDMBVtx->GetXaxis()->SetRangeUser(-etaLimit, etaLimit);
 310   histESDMB->GetXaxis()->SetRangeUser(-etaLimit, etaLimit);
 311   histESD->GetXaxis()->SetRangeUser(-etaLimit, etaLimit);
 312
 313   histESDNoPt->GetXaxis()->SetRangeUser(-etaLimit, etaLimit);
 314   histESDMBNoPt->GetXaxis()->SetRangeUser(-etaLimit, etaLimit);
 315   histESDMBVtxNoPt->GetXaxis()->SetRangeUser(-etaLimit, etaLimit);
 316   histESDMBTracksNoPt->GetXaxis()->SetRangeUser(-etaLimit, etaLimit);
 317
 318   TCanvas* canvas = new TCanvas("dNdEta1", "dNdEta1", 500, 500);
 319
 320   dummy->DrawCopy();
 321   histESDMBVtx->Draw("SAME");
 322   histESDMB->Draw("SAME");
 323   histESD->Draw("SAME");
 324
 325   if (save)
 326   {
 327     canvas->SaveAs("dNdEta1.gif");
 328     canvas->SaveAs("dNdEta1.eps");
 329   }
 330
 331   if (onlyESD)
 332     return;
 333
 334   loadlibs();
 335
 336   TFile* file2 = TFile::Open("analysis_mc.root");
 337   TH1* histMC = (TH1*) file2->Get("dndeta/dNdEta_corrected")->Clone("cloned");
 338   TH1* histMCTr = (TH1*) file2->Get("dndetaTr/dNdEta_corrected")->Clone("cloned2");
 339   TH1* histMCTrVtx = (TH1*) file2->Get("dndetaTrVtx/dNdEta_corrected")->Clone("cloned3");
 340
 341   dNdEtaAnalysis* fdNdEtaAnalysis = new dNdEtaAnalysis("dndeta", "dndeta");
 342   fdNdEtaAnalysis->LoadHistograms();
 343   dNdEtaAnalysis* fdNdEtaAnalysis2 = (dNdEtaAnalysis*) fdNdEtaAnalysis->Clone();
 344
 345   fdNdEtaAnalysis->Finish(0, 0.3, AlidNdEtaCorrection::kNone);
 346   TH1* histMCPtCut = (TH1*) fdNdEtaAnalysis->GetdNdEtaHistogram(0)->Clone("histMCPtCut");
 347
 348   fdNdEtaAnalysis = new dNdEtaAnalysis("dndetaTr", "dndetaTr");
 349   fdNdEtaAnalysis->LoadHistograms();
 350   fdNdEtaAnalysis->Finish(0, 0.3, AlidNdEtaCorrection::kNone);
 351   TH1* histMCTrPtCut = fdNdEtaAnalysis->GetdNdEtaHistogram(0);
 352
 353   fdNdEtaAnalysis = new dNdEtaAnalysis("dndetaTrVtx", "dndetaTrVtx");
 354   fdNdEtaAnalysis->LoadHistograms();
 355   fdNdEtaAnalysis->Finish(0, 0.3, AlidNdEtaCorrection::kNone);
 356   TH1* histMCTrVtxPtCut = fdNdEtaAnalysis->GetdNdEtaHistogram(0);
 357
 358   fdNdEtaAnalysis = new dNdEtaAnalysis("dndetaTracks", "dndetaTracks");
 359   fdNdEtaAnalysis->LoadHistograms();
 360   fdNdEtaAnalysis->Finish(0, 0.3, AlidNdEtaCorrection::kNone);
 361   TH1* histMCTracksPtCut = fdNdEtaAnalysis->GetdNdEtaHistogram(0);
 362
 363   TCanvas* canvas2 = new TCanvas("dNdEta2", "dNdEta2", 500, 500);
 364
 365   Prepare1DPlot(histMC);
 366   Prepare1DPlot(histMCTr);
 367   Prepare1DPlot(histMCTrVtx);
 368
 369   Prepare1DPlot(histMCPtCut);
 370   Prepare1DPlot(histMCTrPtCut);
 371   Prepare1DPlot(histMCTrVtxPtCut);
 372   if (histMCTracksPtCut)
 373     Prepare1DPlot(histMCTracksPtCut);
 374
 375   histMC->SetLineColor(1);
 376   histMCTr->SetLineColor(2);
 377   histMCTrVtx->SetLineColor(3);
 378
 379   histMCPtCut->SetLineColor(1);
 380   histMCTrPtCut->SetLineColor(2);
 381   histMCTrVtxPtCut->SetLineColor(3);
 382   if (histMCTracksPtCut)
 383     histMCTracksPtCut->SetLineColor(4);
 384
 385   TH2* dummy2 = (TH2F*) dummy->Clone("dummy2");
 386   Prepare1DPlot(dummy2);
 387   dummy2->GetYaxis()->SetRangeUser(0, histESDMBVtx->GetMaximum() * 1.1);
 388
 389   dummy2->DrawCopy();
 390   histMC->Draw("SAME");
 391   histMCTr->Draw("SAME");
 392   histMCTrVtx->Draw("SAME");
 393   histESD->Draw("SAME");
 394   histESDMB->Draw("SAME");
 395   histESDMBVtx->Draw("SAME");
 396   histESDNoPt->Draw("SAME");
 397   histESDMBNoPt->Draw("SAME");
 398   histESDMBVtxNoPt->Draw("SAME");
 399   histESDMBTracksNoPt->Draw("SAME");
 400   histMCPtCut->Draw("SAME");
 401   histMCTrPtCut->Draw("SAME");
 402   histMCTrVtxPtCut->Draw("SAME");
 403   if (histMCTracksPtCut)
 404     histMCTracksPtCut->Draw("SAME");
 405
 406   if (save)
 407   {
 408     canvas2->SaveAs("dNdEta2.gif");
 409     canvas2->SaveAs("dNdEta2.eps");
 410   }
 411
 412   TH1* ratio = (TH1*) histMC->Clone("ratio");
 413   TH1* ratioNoPt = (TH1*) histMCPtCut->Clone("ratioNoPt");
 414
 415   ratio->Divide(histESD);
 416   ratioNoPt->Divide(histESDNoPt);
 417
 418   ratio->GetXaxis()->SetRangeUser(-etaLimit, etaLimit);
 419
 420   ratio->SetLineColor(1);
 421   ratioNoPt->SetLineColor(2);
 422
 423   Double_t average = 0;       // average deviation from 1 in ratio (depends on the number of bins if statistical)
 424   for (Int_t bin = ratio->FindBin(-0.7999); bin <= ratio->FindBin(0.7999); bin++)
 425     average += TMath::Abs(ratio->GetBinContent(bin) - 1);
 426   Int_t nBins = ratio->FindBin(0.7999) - ratio->FindBin(-0.7999) + 1;
 427   average /= nBins;
 428   Printf("Average deviation in |eta| < 0.8 is %.2f %%", average * 100);
 429
 430   PrintIntegratedDeviation(histMC, histESD, "all events");
 431   PrintIntegratedDeviation(histMCTr, histESDMB, "triggered");
 432   PrintIntegratedDeviation(histMCTrVtx, histESDMBVtx, "trigger, vertex");
 433   PrintIntegratedDeviation(histMCPtCut, histESDNoPt, "all events (no pt corr)");
 434   PrintIntegratedDeviation(histMCTrPtCut, histESDMBNoPt, "triggered (no pt corr)");
 435   PrintIntegratedDeviation(histMCTrVtxPtCut, histESDMBVtxNoPt, "trigger, vertex (no pt corr)");
 436
 437   TCanvas* canvas3 = new TCanvas("dNdEta", "dNdEta", 700, 600);
 438   canvas3->Range(0, 0, 1, 1);
 439   //canvas3->Divide(1, 2, 0, 0);
 440
 441   //canvas3->cd(1);
 442   TPad* pad1 = new TPad("dNdEta_1", "", 0, 0.5, 0.98, 0.98);
 443   pad1->Draw();
 444
 445   TPad* pad2 = new TPad("dNdEta_2", "", 0, 0.02, 0.98, 0.5);
 446   pad2->Draw();
 447
 448   pad1->SetRightMargin(0.05);
 449   pad2->SetRightMargin(0.05);
 450
 451   // no border between them
 452   pad1->SetBottomMargin(0);
 453   pad2->SetTopMargin(0);
 454
 455   pad1->cd();
 456
 457   TLegend* legend = new TLegend(0.4, 0.05, 0.65, 0.3);
 458   legend->SetFillColor(0);
 459   legend->AddEntry(histESDMBVtx, "triggered, vertex");
 460   legend->AddEntry(histESDMB, "triggered");
 461   legend->AddEntry(histESD, "all events");
 462   legend->AddEntry(histMC, "MC prediction");
 463
 464   dummy->GetXaxis()->SetLabelSize(0.06);
 465   dummy->GetYaxis()->SetLabelSize(0.06);
 466   dummy->GetXaxis()->SetTitleSize(0.06);
 467   dummy->GetYaxis()->SetTitleSize(0.06);
 468   dummy->GetYaxis()->SetTitleOffset(0.7);
 469   dummy->DrawCopy();
 470   histESDMBVtx->Draw("SAME");
 471   histESDMB->Draw("SAME");
 472   histESD->Draw("SAME");
 473   histMC->Draw("SAME");
 474
 475   legend->Draw();
 476
 477   pad2->cd();
 478   pad2->SetBottomMargin(0.15);
 479
 480   Float_t min = TMath::Min(0.961, ratio->GetMinimum() * 0.95);
 481   Float_t max = TMath::Max(1.049, ratio->GetMaximum() * 1.05);
 482
 483   TH1F dummy3("dummy3", ";#eta;Ratio: MC / ESD", 1, -1.5, 1.5);
 484   dummy3.SetStats(kFALSE);
 485   dummy3.SetBinContent(1, 1);
 486   dummy3.GetYaxis()->SetRangeUser(min, max);
 487   dummy3.SetLineWidth(2);
 488   dummy3.GetXaxis()->SetLabelSize(0.06);
 489   dummy3.GetYaxis()->SetLabelSize(0.06);
 490   dummy3.GetXaxis()->SetTitleSize(0.06);
 491   dummy3.GetYaxis()->SetTitleSize(0.06);
 492   dummy3.GetYaxis()->SetTitleOffset(0.7);
 493   dummy3.DrawCopy();
 494
 495   ratio->Draw("SAME");
 496
 497   //pad2->Draw();
 498
 499   canvas3->Modified();
 500
 501   if (save)
 502   {
 503     canvas3->SaveAs("dNdEta.gif");
 504     canvas3->SaveAs("dNdEta.eps");
 505   }
 506
 507   TCanvas* canvas4 = new TCanvas("ratio", "ratio", 700, 500);
 508
 509   ratio->Draw();
 510   ratioNoPt->Draw("SAME");
 511
 512   TLegend* legend = new TLegend(0.6, 0.7, 0.95, 0.9);
 513   legend->SetFillColor(0);
 514   legend->AddEntry(ratio, "mc/esd");
 515   legend->AddEntry(ratioNoPt, "mc/esd, not pt cut off corrected");
 516   legend->Draw();
 517 }
 518
 519 void ptSpectrum()
 520 {
 521   TFile* file = TFile::Open("analysis_esd.root");
 522   TH1* histESD = (TH1*) file->Get("dndeta/dndeta_pt");
 523
 524   TFile* file2 = TFile::Open("analysis_mc.root");
 525   TH1* histMC = (TH1*) file2->Get("dndeta/dndeta_pt");
 526
 527   TCanvas* canvas = new TCanvas("ptSpectrum", "ptSpectrum", 500, 500);
 528   InitPad();
 529   gPad->SetLogy();
 530
 531   Prepare1DPlot(histMC);
 532   Prepare1DPlot(histESD);
 533
 534   histESD->SetTitle("");
 535   histESD->GetXaxis()->SetTitle("p_{T} [GeV/c]");
 536   histESD->GetYaxis()->SetTitle("#frac{dN}{d#eta dp_{T}} [c/GeV]");
 537
 538   histMC->SetLineColor(kBlue);
 539   histESD->SetLineColor(kRed);
 540
 541   histESD->GetYaxis()->SetTitleOffset(1.5);
 542   histESD->GetXaxis()->SetRangeUser(0, 4.9999);
 543
 544   histESD->SetMaximum(TMath::Max(histESD->GetMaximum(), histMC->GetMaximum()) * 2);
 545
 546   histESD->Draw();
 547   histMC->Draw("SAME");
 548
 549   canvas->SaveAs("ptSpectrum.gif");
 550   canvas->SaveAs("ptSpectrum.eps");
 551 }
 552
 553 void ptCutoff()
 554 {
 555   gSystem->Load("libPWG0base");
 556
 557   TFile::Open("correction_map.root");
 558   AlidNdEtaCorrection* dNdEtaCorrection = new AlidNdEtaCorrection("dndeta_correction", "dndeta_correction");
 559   dNdEtaCorrection->LoadHistograms();
 560
 561   dNdEtaCorrection->GetMeasuredFraction(AlidNdEtaCorrection::kINEL, 0.3, -100, kTRUE);
 562
 563   TH1* hist = dynamic_cast<TH1*> (gROOT->FindObject("generated_pt")->Clone("ptcutoff"));
 564
 565   hist->GetXaxis()->SetRangeUser(0, 0.9999);
 566   hist->SetMinimum(0);
 567
 568   hist->SetTitle("Generated Particles");
 569   Prepare1DPlot(hist);
 570
 571   TCanvas* canvas = new TCanvas("ptCutoff", "ptCutoff", 700, 500);
 572   hist->DrawCopy();
 573
 574   TLine* line = new TLine(0.3, 0 - hist->GetMaximum() * 0, 0.3, hist->GetMaximum() * 1.1);
 575   line->SetLineWidth(3);
 576   line->SetLineColor(kRed);
 577   line->Draw();
 578
 579   canvas->SaveAs("ptCutoff.gif");
 580   canvas->SaveAs("ptCutoff.eps");
 581
 582   TH1F* factor = new TH1F("factor", ";#eta;correction factor", 20, -1, 1.000001);
 583   factor->SetLineWidth(2);
 584   for (Float_t eta = -0.95; eta<1; eta += 0.1)
 585     factor->Fill(eta, 1.0 / dNdEtaCorrection->GetMeasuredFraction(AlidNdEtaCorrection::kINEL, 0.3, eta, kFALSE));
 586
 587   TCanvas* canvas = new TCanvas("ptCutoff_factor", "ptCutoff_factor", 700, 500);
 588   InitPad();
 589
 590   Prepare1DPlot(factor);
 591   factor->GetYaxis()->SetRangeUser(1, 2);
 592   factor->GetYaxis()->SetTitleOffset(1);
 593   factor->Draw();
 594
 595   canvas->SaveAs("ptCutoff_factor.eps");
 596 }
 597
 598 void TriggerBiasVtxRecon(const char* fileName = "correction_map.root", const char* folder = "dndeta_correction")
 599 {
 600   gSystem->Load("libPWG0base");
 601
 602   TFile::Open(fileName);
 603   AlidNdEtaCorrection* dNdEtaCorrection = new AlidNdEtaCorrection("dndeta_correction", "dndeta_correction");
 604   dNdEtaCorrection->LoadHistograms();
 605
 606   TH2* corrTrigger = dNdEtaCorrection->GetTriggerBiasCorrectionINEL()->GetEventCorrection()->GetCorrectionHistogram();
 607   TH2* corrVtx = dNdEtaCorrection->GetVertexRecoCorrection()->GetEventCorrection()->GetCorrectionHistogram();
 608
 609   Prepare2DPlot(corrTrigger);
 610   corrTrigger->SetTitle("b) Trigger bias correction");
 611
 612   Prepare2DPlot(corrVtx);
 613   corrVtx->SetTitle("a) Vertex reconstruction correction");
 614
 615   corrTrigger->GetYaxis()->SetTitle("Multiplicity");
 616   corrVtx->GetYaxis()->SetTitle("Multiplicity");
 617
 618   TCanvas* canvas = new TCanvas("TriggerBiasVtxRecon", "TriggerBiasVtxRecon", 1000, 500);
 619   canvas->Divide(2, 1);
 620
 621   canvas->cd(1);
 622   InitPadCOLZ();
 623   corrVtx->DrawCopy("COLZ");
 624
 625   canvas->cd(2);
 626   InitPadCOLZ();
 627   corrTrigger->DrawCopy("COLZ");
 628
 629   canvas->SaveAs(Form("TriggerBiasVtxRecon_%d.gif", gMax));
 630   canvas->SaveAs(Form("TriggerBiasVtxRecon_%d.eps", gMax));
 631
 632   canvas = new TCanvas("TriggerBiasVtxReconZoom", "TriggerBiasVtxReconZoom", 1000, 500);
 633   canvas->Divide(2, 1);
 634
 635   corrTrigger->GetYaxis()->SetRangeUser(0, 5);
 636   corrVtx->GetYaxis()->SetRangeUser(0, 5);
 637
 638   canvas->cd(1);
 639   InitPadCOLZ();
 640   corrVtx->DrawCopy("COLZ");
 641
 642   canvas->cd(2);
 643   InitPadCOLZ();
 644   corrTrigger->DrawCopy("COLZ");
 645
 646   canvas->SaveAs(Form("TriggerBiasVtxReconZoom_%d.gif", gMax));
 647   canvas->SaveAs(Form("TriggerBiasVtxReconZoom_%d.eps", gMax));
 648 }
 649
 650 void TriggerBias(const char* fileName = "correction_map.root")
 651 {
 652   TFile* file = TFile::Open(fileName);
 653
 654   TH2* corr = dynamic_cast<TH2*> (file->Get("dndeta_correction/corr_dndeta_correction_trigger"));
 655
 656   Prepare2DPlot(corr);
 657   corr->SetTitle("Trigger bias correction");
 658
 659   TCanvas* canvas = new TCanvas("TriggerBias", "TriggerBias", 500, 500);
 660   InitPadCOLZ();
 661   corr->DrawCopy("COLZ");
 662
 663   canvas->SaveAs(Form("TriggerBias_%d.gif", gMax));
 664   canvas->SaveAs(Form("TriggerBias_%d.eps", gMax));
 665
 666   corr->GetYaxis()->SetRangeUser(0, 5);
 667
 668   canvas = new TCanvas("TriggerBiasZoom", "TriggerBiasZoom", 500, 500);
 669   InitPadCOLZ();
 670   corr->DrawCopy("COLZ");
 671
 672   canvas->SaveAs(Form("TriggerBiasZoom_%d.gif", gMax));
 673   canvas->SaveAs(Form("TriggerBiasZoom_%d.eps", gMax));
 674 }
 675
 676 void TriggerBias1D(const char* fileName = "correction_map.root", const char* folderName = "dndeta_correction")
 677 {
 678   gSystem->Load("libPWG0base");
 679
 680   TFile* file = TFile::Open(fileName);
 681   AlidNdEtaCorrection* dNdEtaCorrection = new AlidNdEtaCorrection(folderName, folderName);
 682   dNdEtaCorrection->LoadHistograms();
 683
 684   TH1* hist = dNdEtaCorrection->GetTriggerBiasCorrectionINEL()->GetEventCorrection()->Get1DCorrection("x");
 685   TH1* hist2 = dNdEtaCorrection->GetTriggerBiasCorrectionINEL()->GetEventCorrection()->Get1DCorrection("y", -10, 10);
 686
 687   TCanvas* canvas = new TCanvas("TriggerBias1D", "TriggerBias1D", 1000, 500);
 688   canvas->Divide(2, 1);
 689
 690   canvas->cd(1);
 691   InitPad();
 692
 693   Prepare1DPlot(hist);
 694   hist->SetTitle("");
 695   hist->GetYaxis()->SetTitle("correction factor");
 696   hist->GetYaxis()->SetRangeUser(1, 1.5);
 697   hist->GetYaxis()->SetTitleOffset(1.6);
 698   hist->Draw();
 699
 700   canvas->cd(2);
 701   InitPad();
 702
 703   Prepare1DPlot(hist2);
 704   hist2->SetTitle("");
 705   hist2->GetYaxis()->SetTitle("correction factor");
 706   hist2->GetXaxis()->SetRangeUser(0, 5);
 707   hist2->GetYaxis()->SetTitleOffset(1.6);
 708   hist2->GetXaxis()->SetTitle("multiplicity");
 709   hist2->Draw();
 710
 711   TPaveText* pave = new TPaveText(0.6, 0.8, 0.8, 0.85, "NDC");
 712   pave->SetFillColor(0);
 713   pave->AddText("|z| < 10 cm");
 714   pave->Draw();
 715
 716   canvas->SaveAs("TriggerBias1D.eps");
 717 }
 718
 719 void VtxRecon()
 720 {
 721   TFile* file = TFile::Open("correction_map.root");
 722
 723   TH2* corr = dynamic_cast<TH2*> (file->Get("dndeta_correction/corr_dndeta_correction_vtxReco"));
 724
 725   Prepare2DPlot(corr);
 726   corr->SetTitle("Vertex reconstruction correction");
 727
 728   TCanvas* canvas = new TCanvas("VtxRecon", "VtxRecon", 500, 500);
 729   InitPadCOLZ();
 730   corr->DrawCopy("COLZ");
 731
 732   canvas->SaveAs(Form("VtxRecon_%d.eps", gMax));
 733   canvas->SaveAs(Form("VtxRecon_%d.eps", gMax));
 734
 735   corr->GetYaxis()->SetRangeUser(0, 5);
 736
 737   canvas = new TCanvas("VtxReconZoom", "VtxReconZoom", 500, 500);
 738   InitPadCOLZ();
 739   corr->DrawCopy("COLZ");
 740
 741   canvas->SaveAs(Form("VtxReconZoom_%d.gif", gMax));
 742   canvas->SaveAs(Form("VtxReconZoom_%d.eps", gMax));
 743 }
 744
 745 void VtxRecon1D(const char* fileName = "correction_map.root", const char* folderName = "dndeta_correction")
 746 {
 747   gSystem->Load("libPWG0base");
 748
 749   TFile* file = TFile::Open(fileName);
 750   AlidNdEtaCorrection* dNdEtaCorrection = new AlidNdEtaCorrection(folderName, folderName);
 751   dNdEtaCorrection->LoadHistograms();
 752
 753   TH1* hist = dNdEtaCorrection->GetVertexRecoCorrection()->GetEventCorrection()->Get1DCorrection("x");
 754   TH1* hist2 = dNdEtaCorrection->GetVertexRecoCorrection()->GetEventCorrection()->Get1DCorrection("y", -10, 10);
 755
 756   TCanvas* canvas = new TCanvas("VtxRecon1D", "VtxRecon1D", 1000, 500);
 757   canvas->Divide(2, 1);
 758
 759   canvas->cd(1);
 760   InitPad();
 761
 762   Prepare1DPlot(hist);
 763   hist->SetTitle("");
 764   hist->GetYaxis()->SetTitle("correction factor");
 765   hist->GetYaxis()->SetRangeUser(1, 1.8);
 766   hist->GetYaxis()->SetTitleOffset(1.6);
 767   hist->DrawCopy();
 768
 769   canvas->cd(2);
 770   InitPad();
 771
 772   Prepare1DPlot(hist2);
 773   hist2->SetTitle("");
 774   hist2->GetYaxis()->SetTitle("correction factor");
 775   hist2->GetXaxis()->SetRangeUser(0, 20);
 776   hist2->GetYaxis()->SetTitleOffset(1.6);
 777   hist2->GetXaxis()->SetTitle("multiplicity");
 778   hist2->Draw();
 779
 780   TPaveText* pave = new TPaveText(0.6, 0.8, 0.8, 0.85, "NDC");
 781   pave->SetFillColor(0);
 782   pave->AddText("|z| < 10 cm");
 783   pave->Draw();
 784
 785   canvas->SaveAs("VtxRecon1D.eps");
 786
 787   Correction1DCreatePlots(fileName, folderName, 9.9, 2);
 788
 789   TH1* corrX = dynamic_cast<TH1*> (gROOT->FindObject("generated_x_div_measured_x"));
 790   TH1* corrZ = dynamic_cast<TH1*> (gROOT->FindObject("generated_z_div_measured_z"));
 791
 792   Prepare1DPlot(corrX);
 793   Prepare1DPlot(corrZ);
 794
 795   corrX->GetYaxis()->SetTitleOffset(1.5);
 796   corrZ->GetYaxis()->SetTitleOffset(1.5);
 797
 798   corrX->SetTitle("a) z projection");
 799   corrZ->SetTitle("b) p_{T} projection");
 800
 801   corrX->GetYaxis()->SetTitle("correction factor");
 802   corrZ->GetYaxis()->SetTitle("correction factor");
 803
 804   corrZ->GetXaxis()->SetRangeUser(0.11, 9.9);
 805
 806   TString canvasName;
 807   canvasName.Form("VtxRecon1D_Track");
 808   TCanvas* canvas = new TCanvas(canvasName, canvasName, 800, 400);
 809   canvas->Divide(2, 1);
 810
 811   canvas->cd(1);
 812   InitPad();
 813   corrX->DrawCopy();
 814
 815   canvas->cd(2);
 816   InitPad();
 817   gPad->SetLogx();
 818   corrZ->Draw();
 819
 820   canvas->SaveAs("VtxRecon1D_Track.eps");
 821   canvas->SaveAs("VtxRecon1D_Track.gif");
 822 }
 823
 824 void Track2ParticleAsNumber(const char* fileName = "correction_map.root")
 825 {
 826   gSystem->Load("libPWG0base");
 827
 828   TFile::Open(fileName);
 829   AlidNdEtaCorrection* dNdEtaCorrection = new AlidNdEtaCorrection("dndeta_correction", "dndeta_correction");
 830   dNdEtaCorrection->LoadHistograms(fileName, "dndeta_correction");
 831
 832   TH3F* gene = dNdEtaCorrection->GetTrack2ParticleCorrection()->GetGeneratedHistogram();
 833   TH3F* meas = dNdEtaCorrection->GetTrack2ParticleCorrection()->GetMeasuredHistogram();
 834
 835   TH3F* gene = dNdEtaCorrection->GetTrack2ParticleCorrection()->GetGeneratedHistogram();
 836   TH3F* meas = dNdEtaCorrection->GetTrack2ParticleCorrection()->GetMeasuredHistogram();
 837
 838   gene->GetYaxis()->SetRangeUser(-0.8, 0.8);
 839   meas->GetYaxis()->SetRangeUser(-0.8, 0.8);
 840   gene->GetXaxis()->SetRangeUser(-10, 10);
 841   meas->GetXaxis()->SetRangeUser(-10, 10);
 842
 843   Float_t eff1 = gene->Integral() / meas->Integral();
 844   Float_t error1 = TMath::Sqrt(gene->Integral()) / meas->Integral();
 845
 846   printf("Correction without pT cut: %f +- %f\n", eff1, error1);
 847
 848   gene->GetZaxis()->SetRangeUser(0.3, 10);
 849   meas->GetZaxis()->SetRangeUser(0.3, 10);
 850
 851   Float_t eff2 = gene->Integral() / meas->Integral();
 852   Float_t error2 = TMath::Sqrt(gene->Integral()) / meas->Integral();
 853
 854   printf("Correction with pT cut: %f +- %f\n", eff2, error2);
 855
 856   gene->GetZaxis()->SetRangeUser(0.3, 1);
 857   meas->GetZaxis()->SetRangeUser(0.3, 1);
 858
 859   Float_t eff3 = gene->Integral() / meas->Integral();
 860   Float_t error3 = TMath::Sqrt(gene->Integral()) / meas->Integral();
 861
 862   printf("Correction with 0.3 < pT < 0.5: %f +- %f\n", eff3, error3);
 863 }
 864
 865 void Correction1DCreatePlots(const char* fileName = "correction_map.root", const char* folderName = "dndeta_correction", Float_t upperPtLimit = 9.9, Int_t correctionType = 0)
 866 {
 867   TFile::Open(fileName);
 868   AlidNdEtaCorrection* dNdEtaCorrection = new AlidNdEtaCorrection(folderName, folderName);
 869   dNdEtaCorrection->LoadHistograms();
 870
 871   TH3F* gene = dNdEtaCorrection->GetCorrection(correctionType)->GetTrackCorrection()->GetGeneratedHistogram();
 872   TH3F* meas = dNdEtaCorrection->GetCorrection(correctionType)->GetTrackCorrection()->GetMeasuredHistogram();
 873
 874   gene->GetZaxis()->SetRangeUser(0.3, upperPtLimit);
 875   meas->GetZaxis()->SetRangeUser(0.3, upperPtLimit);
 876   gene->GetYaxis()->SetRangeUser(-0.8, 0.8);
 877   meas->GetYaxis()->SetRangeUser(-0.8, 0.8);
 878   AliPWG0Helper::CreateDividedProjections(gene, meas, "x", kTRUE);
 879   gene->GetYaxis()->SetRange(0, 0);
 880   meas->GetYaxis()->SetRange(0, 0);
 881
 882   gene->GetXaxis()->SetRangeUser(-10, 10);
 883   meas->GetXaxis()->SetRangeUser(-10, 10);
 884   AliPWG0Helper::CreateDividedProjections(gene, meas, "y", kTRUE);
 885   gene->GetZaxis()->SetRange(0, 0);
 886   meas->GetZaxis()->SetRange(0, 0);
 887
 888   gene->GetYaxis()->SetRangeUser(-0.8, 0.8);
 889   meas->GetYaxis()->SetRangeUser(-0.8, 0.8);
 890   AliPWG0Helper::CreateDividedProjections(gene, meas, "z", kTRUE);
 891 }
 892
 893 void Correction1D(Int_t correctionType = 0, const char* fileName = "correction_map.root", const char* folder = "dndeta_correction", Float_t upperPtLimit = 9.9)
 894 {
 895   gSystem->Load("libPWG0base");
 896
 897   Correction1DCreatePlots(fileName, folder, upperPtLimit, correctionType);
 898
 899   TH1* corrX = dynamic_cast<TH1*> (gROOT->FindObject(Form("generated_x_div_measured_x", folder, folder)));
 900   TH1* corrY = dynamic_cast<TH1*> (gROOT->FindObject(Form("generated_y_div_measured_y", folder, folder)));
 901   TH1* corrZ = dynamic_cast<TH1*> (gROOT->FindObject(Form("generated_z_div_measured_z", folder, folder)));
 902
 903   Prepare1DPlot(corrX);
 904   Prepare1DPlot(corrY);
 905   Prepare1DPlot(corrZ);
 906
 907   corrX->SetTitle("a) z projection");
 908   corrY->SetTitle("b) #eta projection");
 909   corrZ->SetTitle("c) p_{T} projection");
 910
 911   corrX->GetYaxis()->SetTitle("correction factor");
 912   corrY->GetYaxis()->SetTitle("correction factor");
 913   corrZ->GetYaxis()->SetTitle("correction factor");
 914
 915   corrZ->GetXaxis()->SetRangeUser(0, upperPtLimit);
 916
 917   TString canvasName;
 918   canvasName.Form("Correction1D_%s", folder);
 919   TCanvas* canvas = new TCanvas(canvasName, canvasName, 1200, 400);
 920   canvas->Divide(3, 1);
 921
 922   canvas->cd(1);
 923   InitPad();
 924   corrX->DrawCopy();
 925
 926   canvas->cd(2);
 927   InitPad();
 928   corrY->Draw();
 929
 930   canvas->cd(3);
 931   InitPad();
 932   corrZ->Draw();
 933
 934   canvas->SaveAs(Form("Correction1D_%d_%s_%f.gif", correctionType, fileName, upperPtLimit));
 935   canvas->SaveAs(Form("Correction1D_%d_%s_%f.eps", correctionType, fileName, upperPtLimit));
 936 }
 937
 938 void Track2Particle1D(const char* fileName = "correction_map.root", const char* folder = "dndeta_correction", Float_t upperPtLimit = 9.9)
 939 {
 940   gSystem->Load("libPWG0base");
 941
 942   Correction1DCreatePlots(fileName, folder, upperPtLimit, AlidNdEtaCorrection::kTrack2Particle);
 943
 944   TH1* corrX = dynamic_cast<TH1*> (gROOT->FindObject(Form("generated_x_div_measured_x", folder, folder)));
 945   TH1* corrY = dynamic_cast<TH1*> (gROOT->FindObject(Form("generated_y_div_measured_y", folder, folder)));
 946   TH1* corrZ = dynamic_cast<TH1*> (gROOT->FindObject(Form("generated_z_div_measured_z", folder, folder)));
 947
 948   Prepare1DPlot(corrX);
 949   Prepare1DPlot(corrY);
 950   Prepare1DPlot(corrZ);
 951
 952   corrX->SetTitle("a) z projection");
 953   corrY->SetTitle("a) #eta projection");
 954   corrZ->SetTitle("b) p_{T} projection");
 955
 956   corrY->GetYaxis()->SetTitle("correction factor");
 957   corrZ->GetYaxis()->SetTitle("correction factor");
 958
 959   corrZ->GetXaxis()->SetRangeUser(0, upperPtLimit);
 960
 961   TString canvasName;
 962   canvasName.Form("Track2Particle1D_%s", folder);
 963   TCanvas* canvas = new TCanvas(canvasName, canvasName, 1200, 400);
 964   canvas->Divide(3, 1);
 965
 966   canvas->cd(1);
 967   InitPad();
 968   corrX->DrawCopy();
 969
 970   canvas->cd(2);
 971   InitPad();
 972   corrY->Draw();
 973
 974   canvas->cd(3);
 975   InitPad();
 976   corrZ->Draw();
 977
 978   canvas->SaveAs(Form("Track2Particle1D_%s_%f.gif", fileName, upperPtLimit));
 979   canvas->SaveAs(Form("Track2Particle1D_%s_%f.eps", fileName, upperPtLimit));
 980
 981   //TPaveText* pave = new TPaveText(-0.4, 1.35, 0.4, 1.45);
 982
 983   canvasName.Form("Track2Particle1D_%s_etapt", folder);
 984   TCanvas* canvas = new TCanvas(canvasName, canvasName, 1000, 500);
 985   canvas->Divide(2, 1);
 986
 987   canvas->cd(1);
 988   InitPad();
 989   corrY->GetXaxis()->SetRangeUser(-0.99, 0.99);
 990   corrY->GetYaxis()->SetRangeUser(1, 1.5);
 991   corrY->GetYaxis()->SetTitleOffset(1.5);
 992   corrY->DrawCopy();
 993   TPaveText* pave = new TPaveText(0.3, 0.7, 0.7, 0.8, "NDC");
 994   pave->AddText("|z| < 10 cm");
 995   pave->AddText("0.3 GeV/c < p_{T} < 10 GeV/c");
 996   pave->Draw();
 997
 998   canvas->cd(2);
 999   InitPad();
1000   gPad->SetLogx();
1001   corrZ->GetYaxis()->SetRangeUser(1, 2.5);
1002   corrZ->GetXaxis()->SetRangeUser(0.101, upperPtLimit);
1003   corrZ->GetYaxis()->SetTitleOffset(1.5);
1004   corrZ->DrawCopy();
1005   pave = new TPaveText(0.5, 0.7, 0.8, 0.8, "NDC");
1006   pave->AddText("|z| < 10 cm");
1007   pave->AddText("|#eta| < 0.8");
1008   pave->Draw();
1009
1010   canvas->SaveAs(Form("Track2Particle1D_etapt_%s_%f.eps", fileName, upperPtLimit));
1011   canvas->SaveAs(Form("Track2Particle1D_etapt_%s_%f.gif", fileName, upperPtLimit));
1012 }
1013
1014 void CompareTrack2Particle1D(Float_t upperPtLimit = 9.9)
1015 {
1016   gSystem->Load("libPWG0base");
1017
1018   // particle type
1019   for (Int_t particle=0; particle<4; ++particle)
1020   {
1021     TString dirName;
1022     dirName.Form("correction_%d", particle);
1023     Track2Particle1DCreatePlots("systematics-detail-only-positive.root", dirName, upperPtLimit);
1024
1025     TString tmpx, tmpy, tmpz;
1026     tmpx.Form("gene_%s_nTrackToNPart_x_div_meas_%s_nTrackToNPart_x", dirName.Data(), dirName.Data());
1027     tmpy.Form("gene_%s_nTrackToNPart_y_div_meas_%s_nTrackToNPart_y", dirName.Data(), dirName.Data());
1028     tmpz.Form("gene_%s_nTrackToNPart_z_div_meas_%s_nTrackToNPart_z", dirName.Data(), dirName.Data());
1029
1030     TH1* posX = dynamic_cast<TH1*> (gROOT->FindObject(tmpx)->Clone("pos_x"));
1031     TH1* posY = dynamic_cast<TH1*> (gROOT->FindObject(tmpy)->Clone("pos_y"));
1032     TH1* posZ = dynamic_cast<TH1*> (gROOT->FindObject(tmpz)->Clone("pos_z"));
1033
1034     Track2Particle1DCreatePlots("systematics-detail-only-negative.root", dirName, upperPtLimit);
1035
1036     TH1* negX = dynamic_cast<TH1*> (gROOT->FindObject(tmpx)->Clone("neg_x"));
1037     TH1* negY = dynamic_cast<TH1*> (gROOT->FindObject(tmpy)->Clone("neg_y"));
1038     TH1* negZ = dynamic_cast<TH1*> (gROOT->FindObject(tmpz)->Clone("neg_z"));
1039
1040     posX->Divide(negX);
1041     posY->Divide(negY);
1042     posZ->Divide(negZ);
1043
1044     Prepare1DPlot(posX);
1045     Prepare1DPlot(posY);
1046     Prepare1DPlot(posZ);
1047
1048     Float_t min = 0.8;
1049     Float_t max = 1.2;
1050
1051     posX->SetMinimum(min);
1052     posX->SetMaximum(max);
1053     posY->SetMinimum(min);
1054     posY->SetMaximum(max);
1055     posZ->SetMinimum(min);
1056     posZ->SetMaximum(max);
1057
1058     posZ->GetXaxis()->SetRangeUser(0, upperPtLimit);
1059
1060     posX->GetYaxis()->SetTitleOffset(1.7);
1061     posX->GetYaxis()->SetTitle("C_{+} / C_{-}");
1062     posY->GetYaxis()->SetTitleOffset(1.7);
1063     posY->GetYaxis()->SetTitle("C_{+} / C_{-}");
1064     posZ->GetYaxis()->SetTitleOffset(1.7);
1065     posZ->GetYaxis()->SetTitle("C_{+} / C_{-}");
1066
1067     posZ->GetXaxis()->SetRangeUser(0, 1);
1068
1069     TString canvasName;
1070     canvasName.Form("PosNegRatios_%s_%f", ((particle == 0) ? "Pi" : ((particle == 1) ? "K" : ((particle == 2) ? "p" : "other"))), upperPtLimit);
1071
1072     TCanvas* canvas = new TCanvas(canvasName, canvasName, 1200, 400);
1073     canvas->Divide(3, 1);
1074
1075     canvas->cd(1);
1076     InitPad();
1077     posX->DrawCopy();
1078
1079     canvas->cd(2);
1080     InitPad();
1081     posY->DrawCopy();
1082
1083     canvas->cd(3);
1084     InitPad();
1085     posZ->DrawCopy();
1086
1087     canvas->SaveAs(Form("%s.gif", canvas->GetName()));
1088     canvas->SaveAs(Form("%s.eps", canvas->GetName()));
1089   }
1090 }
1091
1092 void Track2Particle2DCreatePlots(const char* fileName = "correction_map.root")
1093 {
1094   TFile::Open(fileName);
1095   AlidNdEtaCorrection* dNdEtaCorrection = new AlidNdEtaCorrection("dndeta_correction", "dndeta_correction");
1096   dNdEtaCorrection->LoadHistograms();
1097
1098   TH3F* gene = dNdEtaCorrection->GetTrack2ParticleCorrection()->GetTrackCorrection()->GetGeneratedHistogram();
1099   TH3F* meas = dNdEtaCorrection->GetTrack2ParticleCorrection()->GetTrackCorrection()->GetMeasuredHistogram();
1100
1101   gene->GetZaxis()->SetRangeUser(0.3, 10);
1102   meas->GetZaxis()->SetRangeUser(0.3, 10);
1103   AliPWG0Helper::CreateDividedProjections(gene, meas, "yx");
1104   gene->GetZaxis()->SetRange(0, 0);
1105   meas->GetZaxis()->SetRange(0, 0);
1106
1107   gene->GetYaxis()->SetRangeUser(-0.8, 0.8);
1108   meas->GetYaxis()->SetRangeUser(-0.8, 0.8);
1109   AliPWG0Helper::CreateDividedProjections(gene, meas, "zx");
1110   gene->GetYaxis()->SetRange(0, 0);
1111   meas->GetYaxis()->SetRange(0, 0);
1112
1113   gene->GetXaxis()->SetRangeUser(-10, 10);
1114   meas->GetXaxis()->SetRangeUser(-10, 10);
1115   AliPWG0Helper::CreateDividedProjections(gene, meas, "zy");
1116   gene->GetXaxis()->SetRange(0, 0);
1117   meas->GetXaxis()->SetRange(0, 0);
1118 }
1119
1120 void Track2Particle2D(const char* fileName = "correction_map.root", const char* folder = "dndeta_correction")
1121 {
1122   gSystem->Load("libPWG0base");
1123
1124   Track2Particle2DCreatePlots(fileName);
1125
1126   TH2* corrYX = dynamic_cast<TH2*> (gROOT->FindObject("generated_yx_div_measured_yx"));
1127   TH2* corrZX = dynamic_cast<TH2*> (gROOT->FindObject("generated_zx_div_measured_zx"));
1128   TH2* corrZY = dynamic_cast<TH2*> (gROOT->FindObject("generated_zy_div_measured_zy"));
1129
1130   /* this reads them from the file
1131   TH2* corrYX = dynamic_cast<TH2*> (file->Get("dndeta_correction/gene_nTrackToNPart_yx_div_meas_nTrackToNPart_yx"));
1132   TH2* corrZX = dynamic_cast<TH2*> (file->Get("dndeta_correction/gene_nTrackToNPart_zx_div_meas_nTrackToNPart_zx"));
1133   TH2* corrZY = dynamic_cast<TH2*> (file->Get("dndeta_correction/gene_nTrackToNPart_zy_div_meas_nTrackToNPart_zy"));*/
1134
1135   Prepare2DPlot(corrYX);
1136   Prepare2DPlot(corrZX);
1137   Prepare2DPlot(corrZY);
1138
1139   const char* title = "";
1140   corrYX->SetTitle(title);
1141   corrZX->SetTitle(title);
1142   corrZY->SetTitle(title);
1143
1144   TCanvas* canvas = new TCanvas("Track2Particle2D", "Track2Particle2D", 1200, 400);
1145   canvas->Divide(3, 1);
1146
1147   canvas->cd(1);
1148   InitPadCOLZ();
1149   corrYX->Draw("COLZ");
1150
1151   canvas->cd(2);
1152   InitPadCOLZ();
1153   corrZX->Draw("COLZ");
1154
1155   canvas->cd(3);
1156   InitPadCOLZ();
1157   corrZY->Draw("COLZ");
1158
1159   canvas->SaveAs(Form("Track2Particle2D_%d.gif", gMax));
1160   canvas->SaveAs(Form("Track2Particle2D_%d.eps", gMax));
1161 }
1162
1163 void CompareTrack2Particle2D()
1164 {
1165   gSystem->Load("libPWG0base");
1166
1167   Track2Particle2DCreatePlots("correction_maponly-positive.root");
1168
1169   TH2* posYX = dynamic_cast<TH2*> (gROOT->FindObject("gene_nTrackToNPart_yx_div_meas_nTrackToNPart_yx")->Clone("pos_yx"));
1170   TH2* posZX = dynamic_cast<TH2*> (gROOT->FindObject("gene_nTrackToNPart_zx_div_meas_nTrackToNPart_zx")->Clone("pos_zx"));
1171   TH2* posZY = dynamic_cast<TH2*> (gROOT->FindObject("gene_nTrackToNPart_zy_div_meas_nTrackToNPart_zy")->Clone("pos_zy"));
1172
1173   Track2Particle2DCreatePlots("correction_maponly-negative.root");
1174
1175   TH2* negYX = dynamic_cast<TH2*> (gROOT->FindObject("gene_nTrackToNPart_yx_div_meas_nTrackToNPart_yx")->Clone("neg_yx"));
1176   TH2* negZX = dynamic_cast<TH2*> (gROOT->FindObject("gene_nTrackToNPart_zx_div_meas_nTrackToNPart_zx")->Clone("neg_zx"));
1177   TH2* negZY = dynamic_cast<TH2*> (gROOT->FindObject("gene_nTrackToNPart_zy_div_meas_nTrackToNPart_zy")->Clone("neg_zy"));
1178
1179   posYX->Divide(negYX);
1180   posZX->Divide(negZX);
1181   posZY->Divide(negZY);
1182
1183   Prepare2DPlot(posYX);
1184   Prepare2DPlot(posZX);
1185   Prepare2DPlot(posZY);
1186
1187   Float_t min = 0.8;
1188   Float_t max = 1.2;
1189
1190   posYX->SetMinimum(min);
1191   posYX->SetMaximum(max);
1192   posZX->SetMinimum(min);
1193   posZX->SetMaximum(max);
1194   posZY->SetMinimum(min);
1195   posZY->SetMaximum(max);
1196
1197   TCanvas* canvas = new TCanvas("CompareTrack2Particle2D", "CompareTrack2Particle2D", 1200, 400);
1198   canvas->Divide(3, 1);
1199
1200   canvas->cd(1);
1201   InitPadCOLZ();
1202   posYX->Draw("COLZ");
1203
1204   canvas->cd(2);
1205   InitPadCOLZ();
1206   posZX->Draw("COLZ");
1207
1208   canvas->cd(3);
1209   InitPadCOLZ();
1210   posZY->Draw("COLZ");
1211
1212   canvas->SaveAs("CompareTrack2Particle2D.gif");
1213   canvas->SaveAs("CompareTrack2Particle2D.eps");
1214 }
1215
1216 void Track2Particle3D()
1217 {
1218   // get left margin proper
1219
1220   TFile* file = TFile::Open("correction_map.root");
1221
1222   TH3* corr = dynamic_cast<TH3*> (file->Get("dndeta_correction/corr_nTrackToNPart"));
1223
1224   corr->SetTitle("Correction Factor");
1225   SetRanges(corr->GetZaxis());
1226
1227   Prepare3DPlot(corr);
1228
1229   TCanvas* canvas = new TCanvas("Track2Particle3D", "Track2Particle3D", 500, 500);
1230   canvas->SetTheta(29.428);
1231   canvas->SetPhi(16.5726);
1232
1233   corr->Draw();
1234
1235   canvas->SaveAs("Track2Particle3D.gif");
1236   canvas->SaveAs("Track2Particle3D.eps");
1237 }
1238
1239 void Track2Particle3DAll()
1240 {
1241   TFile* file = TFile::Open("correction_map.root");
1242
1243   TH3* gene = dynamic_cast<TH3*> (file->Get("dndeta_correction/gene_nTrackToNPart"));
1244   TH3* meas = dynamic_cast<TH3*> (file->Get("dndeta_correction/meas_nTrackToNPart"));
1245   TH3* corr = dynamic_cast<TH3*> (file->Get("dndeta_correction/corr_nTrackToNPart"));
1246
1247   gene->SetTitle("Generated Particles");
1248   meas->SetTitle("Measured Tracks");
1249   corr->SetTitle("Correction Factor");
1250
1251   Prepare3DPlot(gene);
1252   Prepare3DPlot(meas);
1253   Prepare3DPlot(corr);
1254
1255   TCanvas* canvas = new TCanvas("Track2Particle3DAll", "Track2Particle3DAll", 1200, 400);
1256   canvas->Divide(3, 1);
1257
1258   canvas->cd(1);
1259   InitPad();
1260   gene->Draw();
1261
1262   canvas->cd(2);
1263   meas->Draw();
1264
1265   canvas->cd(3);
1266   corr->Draw();
1267
1268   canvas->SaveAs("Track2Particle3DAll.gif");
1269   canvas->SaveAs("Track2Particle3DAll.eps");
1270 }
1271
1272 void MultiplicityMC(Int_t xRangeMax = 50)
1273 {
1274   TFile* file = TFile::Open("multiplicityMC.root");
1275
1276   if (!file)
1277   {
1278     printf("multiplicityMC.root could not be opened.\n");
1279     return;
1280   }
1281
1282   TH1F* fMultiplicityESD = dynamic_cast<TH1F*> (file->Get("fMultiplicityESD"));
1283   TH1F* fMultiplicityMC = dynamic_cast<TH1F*> (file->Get("fMultiplicityMC"));
1284   TH2F* fCorrelation = dynamic_cast<TH2F*> (file->Get("fCorrelation"));
1285
1286   TH1F* correction = new TH1F("MultiplicityMC_correction", "MultiplicityMC_correction;Ntracks;Npart", 76, -0.5, 75.5);
1287   TH1F* correctionWidth = new TH1F("MultiplicityMC_correctionwidth", "MultiplicityMC_correctionwidth;Ntracks;Npart", 76, -0.5, 75.5);
1288   //fMultiplicityMC->GetNbinsX(), fMultiplicityMC->GetXaxis()->GetXmin(), fMultiplicityMC->GetXaxis()->GetXmax());
1289   for (Int_t i=1; i<=correction->GetNbinsX(); ++i)
1290   {
1291     TH1D* proj = fCorrelation->ProjectionX("_px", i, i+1);
1292     proj->Fit("gaus", "0");
1293     correction->SetBinContent(i, proj->GetFunction("gaus")->GetParameter(1));
1294     correctionWidth->SetBinContent(i, proj->GetFunction("gaus")->GetParameter(2));
1295
1296     continue;
1297
1298     // draw for debugging
1299     new TCanvas;
1300     proj->DrawCopy();
1301     proj->GetFunction("gaus")->DrawCopy("SAME");
1302   }
1303
1304   TH1F* fMultiplicityESDCorrected = new TH1F("fMultiplicityESDCorrected", "fMultiplicityESDCorrected", 2010, -0.5, 200.5);
1305
1306   for (Int_t i=1; i<=correction->GetNbinsX(); ++i)
1307   {
1308     Float_t mean = correction->GetBinContent(i);
1309     Float_t width = correctionWidth->GetBinContent(i);
1310
1311     Int_t fillBegin = fMultiplicityESDCorrected->FindBin(mean - width * 3);
1312     Int_t fillEnd   = fMultiplicityESDCorrected->FindBin(mean + width * 3);
1313     printf("bin %d mean %f width %f, filling from %d to %d\n", i, mean, width, fillBegin, fillEnd);
1314
1315     for (Int_t j=fillBegin; j <= fillEnd; ++j)
1316     {
1317       fMultiplicityESDCorrected->AddBinContent(j, TMath::Gaus(fMultiplicityESDCorrected->GetXaxis()->GetBinCenter(j), mean, width, kTRUE) * fMultiplicityESD->GetBinContent(i));
1318     }
1319   }
1320
1321   TH1F* fMultiplicityESDCorrectedRebinned = dynamic_cast<TH1F*> (fMultiplicityESDCorrected->Clone("fMultiplicityESDCorrectedRebinned"));
1322   fMultiplicityESDCorrectedRebinned->Rebin(10);
1323   fMultiplicityESDCorrectedRebinned->Scale(0.1);
1324
1325   TH1F* ratio = dynamic_cast<TH1F*> (fMultiplicityESD->Clone("multiplicity_ratio"));
1326   ratio->SetTitle("ratio;Ntracks;Nreco/Ngene");
1327   ratio->Divide(fMultiplicityMC);
1328
1329   TH1F* ratio2 = dynamic_cast<TH1F*> (fMultiplicityESDCorrectedRebinned->Clone("multiplicity_ratio_corrected"));
1330   ratio2->Divide(fMultiplicityMC);
1331
1332   TCanvas* canvas = new TCanvas("MultiplicityMC", "MultiplicityMC", 1500, 1000);
1333   canvas->Divide(3, 2);
1334
1335   fMultiplicityESD->GetXaxis()->SetRangeUser(0, xRangeMax);
1336   ratio->GetXaxis()->SetRangeUser(0, xRangeMax);
1337   fCorrelation->GetXaxis()->SetRangeUser(0, xRangeMax);
1338   fCorrelation->GetYaxis()->SetRangeUser(0, xRangeMax);
1339   correction->GetXaxis()->SetRangeUser(0, xRangeMax);
1340   fMultiplicityESDCorrected->GetXaxis()->SetRangeUser(0, xRangeMax);
1341   fMultiplicityESDCorrectedRebinned->GetXaxis()->SetRangeUser(0, xRangeMax);
1342
1343   canvas->cd(1); //InitPad();
1344   fMultiplicityESD->Draw();
1345   fMultiplicityMC->SetLineColor(2);
1346   fMultiplicityMC->Draw("SAME");
1347
1348   TLegend* legend = new TLegend(0.6, 0.7, 0.85, 0.85);
1349   legend->AddEntry(fMultiplicityESD, "ESD");
1350   legend->AddEntry(fMultiplicityMC, "MC");
1351   legend->Draw();
1352
1353   canvas->cd(2);
1354   fCorrelation->Draw("COLZ");
1355
1356   canvas->cd(3);
1357   correction->Draw();
1358   //correction->Fit("pol1");
1359   correctionWidth->SetLineColor(2);
1360   correctionWidth->Draw("SAME");
1361
1362   legend = new TLegend(0.2, 0.7, 0.45, 0.85);
1363   legend->AddEntry(correction, "#bar{x}");
1364   legend->AddEntry(correctionWidth, "#sigma");
1365   legend->Draw();
1366
1367   canvas->cd(4);
1368   ratio->Draw();
1369
1370   ratio2->SetLineColor(2);
1371   ratio2->Draw("SAME");
1372
1373   legend = new TLegend(0.6, 0.7, 0.85, 0.85);
1374   legend->AddEntry(ratio, "uncorrected");
1375   legend->AddEntry(ratio2, "corrected");
1376   legend->Draw();
1377
1378   canvas->cd(5);
1379   fMultiplicityESDCorrected->SetLineColor(kBlue);
1380   fMultiplicityESDCorrected->Draw();
1381   fMultiplicityMC->Draw("SAME");
1382   fMultiplicityESD->Draw("SAME");
1383
1384   legend = new TLegend(0.6, 0.7, 0.85, 0.85);
1385   legend->AddEntry(fMultiplicityESDCorrected, "ESD corrected");
1386   legend->AddEntry(fMultiplicityMC, "MC");
1387   legend->AddEntry(fMultiplicityESD, "ESD");
1388   legend->Draw();
1389
1390   canvas->cd(6);
1391   fMultiplicityESDCorrectedRebinned->SetLineColor(kBlue);
1392   fMultiplicityESDCorrectedRebinned->Draw();
1393   fMultiplicityMC->Draw("SAME");
1394
1395   legend = new TLegend(0.6, 0.7, 0.85, 0.85);
1396   legend->AddEntry(fMultiplicityESDCorrectedRebinned, "ESD corrected");
1397   legend->AddEntry(fMultiplicityMC, "MC");
1398   legend->Draw();
1399
1400   canvas->SaveAs("MultiplicityMC.gif");
1401 }
1402
1403 void MultiplicityESD()
1404 {
1405   TFile* file = TFile::Open("multiplicityESD.root");
1406
1407   if (!file)
1408   {
1409     printf("multiplicityESD.root could not be opened.\n");
1410     return;
1411   }
1412
1413   TH1F* fMultiplicityESD = dynamic_cast<TH1F*> (file->Get("fMultiplicity"));
1414
1415   TCanvas* canvas = new TCanvas("MultiplicityESD", "MultiplicityESD", 500, 500);
1416
1417   fMultiplicityESD->Draw();
1418 }
1419
1420 void drawPlots(Int_t max)
1421 {
1422   gMax = max;
1423
1424   ptCutoff();
1425   TriggerBias();
1426   VtxRecon();
1427   Track2Particle2D();
1428   Track2Particle3D();
1429   ptSpectrum();
1430   dNdEta();
1431 }
1432
1433 void drawPlots()
1434 {
1435   drawPlots(5);
1436   drawPlots(2);
1437 }
1438
1439 void CompareCorrection2Measured(const char* dataInput = "analysis_esd_raw.root", const char* correctionMapFile = "correction_map.root", const char* correctionMapFolder = "dndeta_correction")
1440 {
1441   loadlibs();
1442
1443   AlidNdEtaCorrection* dNdEtaCorrection = new AlidNdEtaCorrection(correctionMapFolder, correctionMapFolder);
1444   TFile::Open(correctionMapFile);
1445   dNdEtaCorrection->LoadHistograms();
1446
1447   TFile* file = TFile::Open(dataInput);
1448
1449   if (!file)
1450   {
1451     cout << "Error. File not found" << endl;
1452     return;
1453   }
1454
1455   dNdEtaAnalysis* fdNdEtaAnalysis = new dNdEtaAnalysis("dndeta", "dndeta");
1456   fdNdEtaAnalysis->LoadHistograms("fdNdEtaAnalysisESD");
1457
1458   gROOT->cd();
1459
1460   TH3* hist1 = (TH3*) dNdEtaCorrection->GetTrack2ParticleCorrection()->GetTrackCorrection()->GetMeasuredHistogram()->Clone("mc");
1461   hist1->SetTitle("mc");
1462   Printf("mc contains %f entries", hist1->Integral());
1463   Printf("mc contains %f entries in |vtx-z| < 10, pt > 0.3", hist1->Integral(hist1->GetXaxis()->FindBin(-9.9), hist1->GetXaxis()->FindBin(9.9), 1, hist1->GetNbinsY(), hist1->GetZaxis()->FindBin(0.301), hist1->GetNbinsZ()));
1464
1465   TH3* hist2 = (TH3*) fdNdEtaAnalysis->GetData()->GetTrackCorrection()->GetMeasuredHistogram()->Clone("esd");
1466   hist2->SetTitle("esd");
1467   Printf("esd contains %f entries", hist2->Integral());
1468   Printf("esd contains %f entries in |vtx-z| < 10, pt > 0.3", hist2->Integral(hist2->GetXaxis()->FindBin(-9.9), hist2->GetXaxis()->FindBin(9.9), 1, hist2->GetNbinsY(), hist2->GetZaxis()->FindBin(0.301), hist2->GetNbinsZ()));
1469
1470   AliPWG0Helper::CreateDividedProjections(hist1, hist2);
1471
1472   new TCanvas; gROOT->FindObject("mc_yx_div_esd_yx")->Draw("COLZ");
1473   new TCanvas; gROOT->FindObject("mc_zx_div_esd_zx")->Draw("COLZ");
1474   new TCanvas; gROOT->FindObject("mc_zy_div_esd_zy")->Draw("COLZ");
1475 }
1476
1477 void CompareMeasured2Measured(const char* dataInput = "analysis_esd_raw.root", const char* dataInput2 = "analysis_esd_raw.root")
1478 {
1479   loadlibs();
1480
1481   TFile* file = TFile::Open(dataInput);
1482
1483   if (!file)
1484   {
1485     cout << "Error. File not found" << endl;
1486     return;
1487   }
1488
1489   dNdEtaAnalysis* fdNdEtaAnalysis = new dNdEtaAnalysis("dndeta", "dndeta");
1490   fdNdEtaAnalysis->LoadHistograms("fdNdEtaAnalysisESD");
1491
1492   TFile* file = TFile::Open(dataInput2);
1493
1494   if (!file)
1495   {
1496     cout << "Error. File not found" << endl;
1497     return;
1498   }
1499
1500   dNdEtaAnalysis* fdNdEtaAnalysis2 = new dNdEtaAnalysis("dndeta2", "dndeta2");
1501   fdNdEtaAnalysis2->LoadHistograms("fdNdEtaAnalysisESD");
1502
1503   gROOT->cd();
1504
1505   TH3* hist1 = (TH3*) fdNdEtaAnalysis->GetData()->GetTrackCorrection()->GetMeasuredHistogram()->Clone("esd1");
1506   hist1->SetTitle("esd1");
1507   Printf("esd1 contains %f entries", hist1->GetEntries());
1508   Printf("esd1 contains %f entries in |vtx-z| < 10, pt > 0.3", hist1->Integral(hist1->GetXaxis()->FindBin(-9.9), hist1->GetXaxis()->FindBin(9.9), 1, hist1->GetNbinsY(), hist1->GetZaxis()->FindBin(0.301), hist1->GetNbinsZ()));
1509
1510   TH3* hist2 = (TH3*) fdNdEtaAnalysis2->GetData()->GetTrackCorrection()->GetMeasuredHistogram()->Clone("esd2");
1511   hist2->SetTitle("esd2");
1512   Printf("esd2 contains %f entries", hist2->GetEntries());
1513   Printf("esd2 contains %f entries in |vtx-z| < 10, pt > 0.3", hist2->Integral(hist2->GetXaxis()->FindBin(-9.9), hist2->GetXaxis()->FindBin(9.9), 1, hist2->GetNbinsY(), hist2->GetZaxis()->FindBin(0.301), hist2->GetNbinsZ()));
1514
1515   AliPWG0Helper::CreateDividedProjections(hist1, hist2);
1516
1517   new TCanvas; gROOT->FindObject("esd1_yx_div_esd2_yx")->Draw("COLZ");
1518   new TCanvas; gROOT->FindObject("esd1_zx_div_esd2_zx")->Draw("COLZ");
1519   new TCanvas; gROOT->FindObject("esd1_zy_div_esd2_zy")->Draw("COLZ");
1520
1521   TH2* event1 = (TH2*) fdNdEtaAnalysis->GetData()->GetEventCorrection()->GetMeasuredHistogram()->Clone("event1");
1522   TH2* event2 = (TH2*) fdNdEtaAnalysis2->GetData()->GetEventCorrection()->GetMeasuredHistogram()->Clone("event2");
1523
1524   Printf("event1 contains %f entries", event1->GetEntries());
1525   Printf("event2 contains %f entries", event2->GetEntries());
1526   Printf("event1 integral is %f", event1->Integral());
1527   Printf("event2 integral is %f", event2->Integral());
1528   Printf("event1 contains %f entries in |vtx-z| < 10", event1->Integral(event1->GetXaxis()->FindBin(-9.9), event1->GetXaxis()->FindBin(9.9), 1, event1->GetNbinsY()));
1529   Printf("event2 contains %f entries in |vtx-z| < 10", event2->Integral(event2->GetXaxis()->FindBin(-9.9), event2->GetXaxis()->FindBin(9.9), 1, event2->GetNbinsY()));
1530
1531   projx1 = event1->ProjectionX();
1532   projx2 = event2->ProjectionX();
1533
1534   new TCanvas; projx1->DrawCopy(); projx2->SetLineColor(2); projx2->DrawCopy("SAME");
1535
1536   projx1->Divide(projx2);
1537   new TCanvas; projx1->Draw();
1538
1539   event1->Divide(event2);
1540   new TCanvas; event1->Draw("COLZ");
1541
1542 }
1543