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 SetRanges(TAxis* axis)
  24 {
  25   if (strcmp(axis->GetTitle(), "#eta") == 0)
  26     axis->SetRangeUser(-1.7999, 1.7999);
  27   if (strcmp(axis->GetTitle(), "p_{T} [GeV/c]") == 0)
  28     axis->SetRangeUser(0, 9.9999);
  29   if (strcmp(axis->GetTitle(), "vtx z [cm]") == 0)
  30     axis->SetRangeUser(-15, 14.9999);
  31   if (strcmp(axis->GetTitle(), "Ntracks") == 0)
  32     axis->SetRangeUser(0, 99.9999);
  33 }
  34
  35 void SetRanges(TH1* hist)
  36 {
  37   SetRanges(hist->GetXaxis());
  38   SetRanges(hist->GetYaxis());
  39   SetRanges(hist->GetZaxis());
  40 }
  41
  42
  43 void Prepare3DPlot(TH3* hist)
  44 {
  45   hist->GetXaxis()->SetTitleOffset(1.5);
  46   hist->GetYaxis()->SetTitleOffset(1.5);
  47   hist->GetZaxis()->SetTitleOffset(1.5);
  48
  49   hist->SetStats(kFALSE);
  50 }
  51
  52 void Prepare2DPlot(TH2* hist)
  53 {
  54   hist->SetStats(kFALSE);
  55   hist->GetYaxis()->SetTitleOffset(1.4);
  56
  57   hist->SetMinimum(0);
  58   hist->SetMaximum(gMax);
  59
  60   SetRanges(hist);
  61 }
  62
  63 void Prepare1DPlot(TH1* hist)
  64 {
  65   hist->SetLineWidth(2);
  66   hist->SetStats(kFALSE);
  67
  68   hist->GetXaxis()->SetLabelOffset(0.02);
  69   hist->GetXaxis()->SetTitleOffset(1.3);
  70   hist->GetYaxis()->SetTitleOffset(1.3);
  71
  72   SetRanges(hist);
  73 }
  74
  75 void InitPad()
  76 {
  77   gPad->Range(0, 0, 1, 1);
  78   gPad->SetLeftMargin(0.15);
  79   //gPad->SetRightMargin(0.05);
  80   //gPad->SetTopMargin(0.13);
  81   gPad->SetBottomMargin(0.12);
  82
  83   gPad->SetGridx();
  84   gPad->SetGridy();
  85 }
  86
  87 void InitPadCOLZ()
  88 {
  89   gPad->Range(0, 0, 1, 1);
  90   gPad->SetRightMargin(0.15);
  91   gPad->SetLeftMargin(0.12);
  92
  93   gPad->SetGridx();
  94   gPad->SetGridy();
  95 }
  96
  97 void dNdEta(Bool_t onlyESD = kFALSE)
  98 {
  99   TFile* file = TFile::Open("analysis_esd.root");
 100   TH1* histESD = (TH1*) file->Get("dndeta/dndeta_dNdEta_corrected_2");
 101   TH1* histESDNoPt = (TH1*) file->Get("dndeta/dndeta_dNdEta_2");
 102   TH1* histESDMB = (TH1*) file->Get("dndeta_mb/dndeta_mb_dNdEta_corrected_2");
 103   TH1* histESDMBVtx = (TH1*) file->Get("dndeta_mbvtx/dndeta_mbvtx_dNdEta_corrected_2");
 104
 105   TCanvas* canvas = new TCanvas("dNdEta1", "dNdEta1", 500, 500);
 106
 107   Prepare1DPlot(histESD);
 108   Prepare1DPlot(histESDNoPt);
 109   Prepare1DPlot(histESDMB);
 110   Prepare1DPlot(histESDMBVtx);
 111
 112   histESD->SetLineColor(0);
 113   histESDMB->SetLineColor(0);
 114   histESDMBVtx->SetLineColor(0);
 115
 116   histESD->SetMarkerColor(kRed);
 117   histESDMB->SetMarkerColor(kBlue);
 118   histESDMBVtx->SetMarkerColor(103);
 119
 120   histESD->SetMarkerStyle(20);
 121   histESDMB->SetMarkerStyle(21);
 122   histESDMBVtx->SetMarkerStyle(22);
 123
 124   TH2F* dummy = new TH2F("dummy", "", 100, -1.5, 1.5, 100, 0, histESDMBVtx->GetMaximum() * 1.1);
 125   Prepare1DPlot(dummy);
 126   dummy->SetStats(kFALSE);
 127   dummy->SetXTitle("#eta");
 128   dummy->SetYTitle("dN_{ch}/d#eta");
 129   dummy->GetYaxis()->SetTitleOffset(1);
 130
 131   Float_t etaLimit = 1.1999;
 132
 133   histESDMBVtx->GetXaxis()->SetRangeUser(-etaLimit, etaLimit);
 134   histESDMB->GetXaxis()->SetRangeUser(-etaLimit, etaLimit);
 135   histESD->GetXaxis()->SetRangeUser(-etaLimit, etaLimit);
 136   histESDNoPt->GetXaxis()->SetRangeUser(-etaLimit, etaLimit);
 137
 138   dummy->DrawCopy();
 139   histESDMBVtx->Draw("SAME");
 140   histESDMB->Draw("SAME");
 141   histESD->Draw("SAME");
 142
 143   canvas->SaveAs("dNdEta1.gif");
 144   canvas->SaveAs("dNdEta1.eps");
 145
 146   if (onlyESD)
 147     return;
 148
 149   TFile* file2 = TFile::Open("analysis_mc.root");
 150   TH1* histMC = (TH1*) file2->Get("dndeta/dndeta_dNdEta_corrected_2")->Clone("cloned");
 151
 152   gSystem->Load("libPWG0base");
 153   dNdEtaAnalysis* fdNdEtaAnalysis = new dNdEtaAnalysis("dndeta", "dndeta");
 154   fdNdEtaAnalysis->LoadHistograms();
 155   fdNdEtaAnalysis->Finish(0, 0.3);
 156   TH1* histMCPtCut = fdNdEtaAnalysis->GetdNdEtaHistogram(2);
 157
 158   TCanvas* canvas2 = new TCanvas("dNdEta2", "dNdEta2", 500, 500);
 159
 160   Prepare1DPlot(histMC);
 161   Prepare1DPlot(histMCPtCut);
 162
 163   histMC->SetLineColor(kBlue);
 164   histMCPtCut->SetLineColor(104);
 165   histESDNoPt->SetLineColor(102);
 166
 167   TH2* dummy2 = (TH2F*) dummy->Clone("dummy2");
 168   Prepare1DPlot(dummy2);
 169   dummy2->GetYaxis()->SetRangeUser(0, histESD->GetMaximum() * 1.1);
 170
 171   dummy2->DrawCopy();
 172   histMC->Draw("SAME");
 173 //  histMC->Draw();
 174   histESD->Draw("SAME");
 175   histESDNoPt->Draw("SAME");
 176   histMCPtCut->Draw("SAME");
 177
 178   canvas2->SaveAs("dNdEta2.gif");
 179   canvas2->SaveAs("dNdEta2.eps");
 180
 181   TCanvas* canvas3 = new TCanvas("dNdEta", "dNdEta", 700, 500);
 182   //InitPad();
 183   gPad->SetBottomMargin(0.12);
 184
 185   TLegend* legend = new TLegend(0.35, 0.2, 0.6, 0.4);
 186   legend->SetFillColor(0);
 187   legend->AddEntry(histESDMBVtx, "triggered, vertex");
 188   legend->AddEntry(histESDMB, "triggered");
 189   legend->AddEntry(histESD, "all events");
 190   legend->AddEntry(histMC, "MC prediction");
 191
 192   dummy->DrawCopy();
 193   histESDMBVtx->Draw("SAME");
 194   histESDMB->Draw("SAME");
 195   histESD->Draw("SAME");
 196   histMC->Draw("SAME");
 197
 198   legend->Draw();
 199
 200   canvas3->SaveAs("dNdEta.gif");
 201   canvas3->SaveAs("dNdEta.eps");
 202
 203   TH1* ratio = (TH1*) histMC->Clone("ratio");
 204   TH1* ratioNoPt = (TH1*) histMCPtCut->Clone("ratioNoPt");
 205
 206   ratio->Divide(histESD);
 207   ratioNoPt->Divide(histESDNoPt);
 208
 209   TCanvas* canvas4 = new TCanvas("ratio", "ratio", 700, 500);
 210
 211   ratio->GetXaxis()->SetRangeUser(-0.7999, 0.7999);
 212
 213   ratio->SetLineColor(1);
 214   ratioNoPt->SetLineColor(2);
 215
 216   ratio->Draw();
 217   ratioNoPt->Draw("SAME");
 218
 219   TLegend* legend = new TLegend(0.6, 0.7, 0.95, 0.9);
 220   legend->SetFillColor(0);
 221   legend->AddEntry(ratio, "mc/esd");
 222   legend->AddEntry(ratioNoPt, "mc/esd, not pt cut off corrected");
 223   legend->Draw();
 224 }
 225
 226 void ptSpectrum()
 227 {
 228   TFile* file = TFile::Open("analysis_esd.root");
 229   TH1* histESD = (TH1*) file->Get("dndeta/dndeta_pt");
 230
 231   TFile* file2 = TFile::Open("analysis_mc.root");
 232   TH1* histMC = (TH1*) file2->Get("dndeta/dndeta_pt");
 233
 234   TCanvas* canvas = new TCanvas("ptSpectrum", "ptSpectrum", 500, 500);
 235   InitPad();
 236   gPad->SetLogy();
 237
 238   Prepare1DPlot(histMC);
 239   Prepare1DPlot(histESD);
 240
 241   histESD->SetTitle("");
 242   histESD->GetXaxis()->SetTitle("p_{T} [GeV/c]");
 243   histESD->GetYaxis()->SetTitle("#frac{dN}{d#eta dp_{T}} [c/GeV]");
 244
 245   histMC->SetLineColor(kBlue);
 246   histESD->SetLineColor(kRed);
 247
 248   histESD->GetYaxis()->SetTitleOffset(1.5);
 249   histESD->GetXaxis()->SetRangeUser(0, 4.9999);
 250
 251   histESD->SetMaximum(TMath::Max(histESD->GetMaximum(), histMC->GetMaximum()) * 2);
 252
 253   histESD->Draw();
 254   histMC->Draw("SAME");
 255
 256   canvas->SaveAs("ptSpectrum.gif");
 257   canvas->SaveAs("ptSpectrum.eps");
 258 }
 259
 260 void ptCutoff()
 261 {
 262   gSystem->Load("libPWG0base");
 263
 264   AlidNdEtaCorrection* dNdEtaCorrection = new AlidNdEtaCorrection("dndeta_correction", "dndeta_correction");
 265   dNdEtaCorrection->LoadHistograms("correction_map.root","dndeta_correction");
 266
 267   dNdEtaCorrection->GetMeasuredFraction(0.3, -100, kTRUE);
 268
 269   TH1* hist = dynamic_cast<TH1*> (gROOT->FindObject("gene_dndeta_correction_nTrackToNPart_pt")->Clone("ptcutoff"));
 270
 271   hist->GetXaxis()->SetRangeUser(0, 0.9999);
 272   hist->SetMinimum(0);
 273
 274   hist->SetTitle("Generated Particles");
 275   Prepare1DPlot(hist);
 276
 277   TCanvas* canvas = new TCanvas("ptCutoff", "ptCutoff", 700, 500);
 278   hist->DrawCopy();
 279
 280   TLine* line = new TLine(0.3, 0 - hist->GetMaximum() * 0, 0.3, hist->GetMaximum() * 1.1);
 281   line->SetLineWidth(3);
 282   line->SetLineColor(kRed);
 283   line->Draw();
 284
 285   canvas->SaveAs("ptCutoff.gif");
 286   canvas->SaveAs("ptCutoff.eps");
 287
 288   TH1F* factor = new TH1F("factor", ";#eta;correction factor", 10, -1, 1.000001);
 289   for (Float_t eta = -0.9; eta<1; eta += 0.2)
 290     factor->Fill(eta, 1.0 / dNdEtaCorrection->GetMeasuredFraction(0.3, eta, kFALSE));
 291
 292   TCanvas* canvas = new TCanvas("ptCutoff_factor", "ptCutoff_factor", 700, 500);
 293   InitPad();
 294
 295   Prepare1DPlot(factor);
 296   factor->GetYaxis()->SetRangeUser(1, 2);
 297   factor->GetYaxis()->SetTitleOffset(1);
 298   factor->Draw();
 299
 300   canvas->SaveAs("ptCutoff_factor.eps");
 301 }
 302
 303 void TriggerBiasVtxRecon(const char* fileName = "correction_map.root", const char* folder = "dndeta_correction")
 304 {
 305   TFile* file = TFile::Open(fileName);
 306
 307   TH2* corrTrigger = dynamic_cast<TH2*> (file->Get(Form("%s/corr_%s_trigger", folder, folder)));
 308   TH2* corrVtx = dynamic_cast<TH2*> (file->Get(Form("%s/corr_%s_vtxReco", folder, folder)));
 309
 310   Prepare2DPlot(corrTrigger);
 311   corrTrigger->SetTitle("a) Trigger bias correction");
 312
 313   Prepare2DPlot(corrVtx);
 314   corrVtx->SetTitle("b) Vertex reconstruction correction");
 315
 316   corrTrigger->GetYaxis()->SetTitle("Multiplicity");
 317   corrVtx->GetYaxis()->SetTitle("Multiplicity");
 318
 319   TCanvas* canvas = new TCanvas("TriggerBiasVtxRecon", "TriggerBiasVtxRecon", 1000, 500);
 320   canvas->Divide(2, 1);
 321
 322   canvas->cd(1);
 323   InitPadCOLZ();
 324   corrTrigger->DrawCopy("COLZ");
 325
 326   canvas->cd(2);
 327   InitPadCOLZ();
 328   corrVtx->DrawCopy("COLZ");
 329
 330   canvas->SaveAs(Form("TriggerBiasVtxRecon_%d.gif", gMax));
 331   canvas->SaveAs(Form("TriggerBiasVtxRecon_%d.eps", gMax));
 332
 333   canvas = new TCanvas("TriggerBiasVtxReconZoom", "TriggerBiasVtxReconZoom", 1000, 500);
 334   canvas->Divide(2, 1);
 335
 336   corrTrigger->GetYaxis()->SetRangeUser(0, 5);
 337   corrVtx->GetYaxis()->SetRangeUser(0, 5);
 338
 339   canvas->cd(1);
 340   InitPadCOLZ();
 341   corrTrigger->DrawCopy("COLZ");
 342
 343   canvas->cd(2);
 344   InitPadCOLZ();
 345   corrVtx->DrawCopy("COLZ");
 346
 347   canvas->SaveAs(Form("TriggerBiasVtxReconZoom_%d.gif", gMax));
 348   canvas->SaveAs(Form("TriggerBiasVtxReconZoom_%d.eps", gMax));
 349 }
 350
 351 void TriggerBias(const char* fileName = "correction_map.root")
 352 {
 353   TFile* file = TFile::Open(fileName);
 354
 355   TH2* corr = dynamic_cast<TH2*> (file->Get("dndeta_correction/corr_dndeta_correction_trigger"));
 356
 357   Prepare2DPlot(corr);
 358   corr->SetTitle("Trigger bias correction");
 359
 360   TCanvas* canvas = new TCanvas("TriggerBias", "TriggerBias", 500, 500);
 361   InitPadCOLZ();
 362   corr->DrawCopy("COLZ");
 363
 364   canvas->SaveAs(Form("TriggerBias_%d.gif", gMax));
 365   canvas->SaveAs(Form("TriggerBias_%d.eps", gMax));
 366
 367   corr->GetYaxis()->SetRangeUser(0, 5);
 368
 369   canvas = new TCanvas("TriggerBiasZoom", "TriggerBiasZoom", 500, 500);
 370   InitPadCOLZ();
 371   corr->DrawCopy("COLZ");
 372
 373   canvas->SaveAs(Form("TriggerBiasZoom_%d.gif", gMax));
 374   canvas->SaveAs(Form("TriggerBiasZoom_%d.eps", gMax));
 375 }
 376
 377 void TriggerBias1D(const char* fileName = "correction_map.root", const char* folderName = "dndeta_correction")
 378 {
 379   gSystem->Load("libPWG0base");
 380
 381   TFile* file = TFile::Open(fileName);
 382   AlidNdEtaCorrection* dNdEtaCorrection = new AlidNdEtaCorrection(folderName, folderName);
 383   dNdEtaCorrection->LoadHistograms(fileName, folderName);
 384
 385   TH1* hist = dNdEtaCorrection->GetTriggerCorrection()->Get1DCorrection("x");
 386
 387   TCanvas* canvas = new TCanvas("TriggerBias1D", "TriggerBias1D", 500, 500);
 388   InitPad();
 389
 390   Prepare1DPlot(hist);
 391   hist->SetTitle("");
 392   hist->GetYaxis()->SetTitle("correction factor");
 393   hist->GetYaxis()->SetRangeUser(1, 1.5);
 394   hist->GetYaxis()->SetTitleOffset(1.6);
 395   hist->Draw();
 396
 397   canvas->SaveAs("TriggerBias1D.eps");
 398 }
 399
 400 void VtxRecon()
 401 {
 402   TFile* file = TFile::Open("correction_map.root");
 403
 404   TH2* corr = dynamic_cast<TH2*> (file->Get("dndeta_correction/corr_dndeta_correction_vtxReco"));
 405
 406   Prepare2DPlot(corr);
 407   corr->SetTitle("Vertex reconstruction correction");
 408
 409   TCanvas* canvas = new TCanvas("VtxRecon", "VtxRecon", 500, 500);
 410   InitPadCOLZ();
 411   corr->DrawCopy("COLZ");
 412
 413   canvas->SaveAs(Form("VtxRecon_%d.eps", gMax));
 414   canvas->SaveAs(Form("VtxRecon_%d.eps", gMax));
 415
 416   corr->GetYaxis()->SetRangeUser(0, 5);
 417
 418   canvas = new TCanvas("VtxReconZoom", "VtxReconZoom", 500, 500);
 419   InitPadCOLZ();
 420   corr->DrawCopy("COLZ");
 421
 422   canvas->SaveAs(Form("VtxReconZoom_%d.gif", gMax));
 423   canvas->SaveAs(Form("VtxReconZoom_%d.eps", gMax));
 424 }
 425
 426 void VtxRecon1D(const char* fileName = "correction_map.root", const char* folderName = "dndeta_correction")
 427 {
 428   gSystem->Load("libPWG0base");
 429
 430   TFile* file = TFile::Open(fileName);
 431   AlidNdEtaCorrection* dNdEtaCorrection = new AlidNdEtaCorrection(folderName, folderName);
 432   dNdEtaCorrection->LoadHistograms(fileName, folderName);
 433
 434   TH1* hist = dNdEtaCorrection->GetVertexRecoCorrection()->Get1DCorrection("x");
 435
 436   TCanvas* canvas = new TCanvas("VtxRecon1D", "VtxRecon1D", 500, 500);
 437   InitPad();
 438
 439   Prepare1DPlot(hist);
 440   hist->SetTitle("");
 441   hist->GetYaxis()->SetTitle("correction factor");
 442   hist->GetYaxis()->SetRangeUser(1, 1.8);
 443   hist->GetYaxis()->SetTitleOffset(1.6);
 444   hist->Draw();
 445
 446   canvas->SaveAs("VtxRecon1D.eps");
 447 }
 448
 449 void Track2ParticleAsNumber(const char* fileName = "correction_map.root")
 450 {
 451   gSystem->Load("libPWG0base");
 452
 453   TFile::Open(fileName);
 454   AlidNdEtaCorrection* dNdEtaCorrection = new AlidNdEtaCorrection("dndeta_correction", "dndeta_correction");
 455   dNdEtaCorrection->LoadHistograms(fileName, "dndeta_correction");
 456
 457   TH3F* gene = dNdEtaCorrection->GetTrack2ParticleCorrection()->GetGeneratedHistogram();
 458   TH3F* meas = dNdEtaCorrection->GetTrack2ParticleCorrection()->GetMeasuredHistogram();
 459
 460   TH3F* gene = dNdEtaCorrection->GetTrack2ParticleCorrection()->GetGeneratedHistogram();
 461   TH3F* meas = dNdEtaCorrection->GetTrack2ParticleCorrection()->GetMeasuredHistogram();
 462
 463   gene->GetYaxis()->SetRangeUser(-0.8, 0.8);
 464   meas->GetYaxis()->SetRangeUser(-0.8, 0.8);
 465   gene->GetXaxis()->SetRangeUser(-10, 10);
 466   meas->GetXaxis()->SetRangeUser(-10, 10);
 467
 468   Float_t eff1 = gene->Integral() / meas->Integral();
 469   Float_t error1 = TMath::Sqrt(gene->Integral()) / meas->Integral();
 470
 471   printf("Correction without pT cut: %f +- %f\n", eff1, error1);
 472
 473   gene->GetZaxis()->SetRangeUser(0.3, 10);
 474   meas->GetZaxis()->SetRangeUser(0.3, 10);
 475
 476   Float_t eff2 = gene->Integral() / meas->Integral();
 477   Float_t error2 = TMath::Sqrt(gene->Integral()) / meas->Integral();
 478
 479   printf("Correction with pT cut: %f +- %f\n", eff2, error2);
 480
 481   gene->GetZaxis()->SetRangeUser(0.3, 1);
 482   meas->GetZaxis()->SetRangeUser(0.3, 1);
 483
 484   Float_t eff3 = gene->Integral() / meas->Integral();
 485   Float_t error3 = TMath::Sqrt(gene->Integral()) / meas->Integral();
 486
 487   printf("Correction with 0.3 < pT < 0.5: %f +- %f\n", eff3, error3);
 488 }
 489
 490 void Track2Particle1DCreatePlots(const char* fileName = "correction_map.root", const char* folderName = "dndeta_correction", Float_t upperPtLimit = 10)
 491 {
 492   TFile::Open(fileName);
 493   AlidNdEtaCorrection* dNdEtaCorrection = new AlidNdEtaCorrection(folderName, folderName);
 494   dNdEtaCorrection->LoadHistograms(fileName, folderName);
 495
 496   TH3F* gene = dNdEtaCorrection->GetTrack2ParticleCorrection()->GetGeneratedHistogram();
 497   TH3F* meas = dNdEtaCorrection->GetTrack2ParticleCorrection()->GetMeasuredHistogram();
 498
 499   gene->GetZaxis()->SetRangeUser(0.3, upperPtLimit);
 500   meas->GetZaxis()->SetRangeUser(0.3, upperPtLimit);
 501   gene->GetYaxis()->SetRangeUser(-0.8, 0.8);
 502   meas->GetYaxis()->SetRangeUser(-0.8, 0.8);
 503   AliPWG0Helper::CreateDividedProjections(gene, meas, "x", kTRUE);
 504   gene->GetYaxis()->SetRange(0, 0);
 505   meas->GetYaxis()->SetRange(0, 0);
 506
 507   gene->GetXaxis()->SetRangeUser(-10, 10);
 508   meas->GetXaxis()->SetRangeUser(-10, 10);
 509   AliPWG0Helper::CreateDividedProjections(gene, meas, "y", kTRUE);
 510   gene->GetZaxis()->SetRange(0, 0);
 511   meas->GetZaxis()->SetRange(0, 0);
 512
 513   gene->GetYaxis()->SetRangeUser(-0.8, 0.8);
 514   meas->GetYaxis()->SetRangeUser(-0.8, 0.8);
 515   AliPWG0Helper::CreateDividedProjections(gene, meas, "z", kTRUE);
 516 }
 517
 518 void Track2Particle1D(const char* fileName = "correction_map.root", const char* folder = "dndeta_correction", Float_t upperPtLimit = 9.9)
 519 {
 520   gSystem->Load("libPWG0base");
 521
 522   Track2Particle1DCreatePlots(fileName, folder, upperPtLimit);
 523
 524   TH1* corrX = dynamic_cast<TH1*> (gROOT->FindObject(Form("gene_%s_nTrackToNPart_x_div_meas_%s_nTrackToNPart_x", folder, folder)));
 525   TH1* corrY = dynamic_cast<TH1*> (gROOT->FindObject(Form("gene_%s_nTrackToNPart_y_div_meas_%s_nTrackToNPart_y", folder, folder)));
 526   TH1* corrZ = dynamic_cast<TH1*> (gROOT->FindObject(Form("gene_%s_nTrackToNPart_z_div_meas_%s_nTrackToNPart_z", folder, folder)));
 527
 528   Prepare1DPlot(corrX);
 529   Prepare1DPlot(corrY);
 530   Prepare1DPlot(corrZ);
 531
 532   //corrX->SetTitle("a) z projection");
 533   corrY->SetTitle("a) #eta projection");
 534   corrZ->SetTitle("b) p_{T} projection");
 535
 536   corrY->GetYaxis()->SetTitle("correction factor");
 537   corrZ->GetYaxis()->SetTitle("correction factor");
 538
 539   corrZ->GetXaxis()->SetRangeUser(0, upperPtLimit);
 540
 541   TString canvasName;
 542   canvasName.Form("Track2Particle1D_%s", folder);
 543   TCanvas* canvas = new TCanvas(canvasName, canvasName, 1200, 400);
 544   canvas->Divide(3, 1);
 545
 546   canvas->cd(1);
 547   InitPad();
 548   corrX->DrawCopy();
 549
 550   canvas->cd(2);
 551   InitPad();
 552   corrY->DrawCopy();
 553
 554   canvas->cd(3);
 555   InitPad();
 556   corrZ->DrawCopy();
 557
 558   canvas->SaveAs(Form("Track2Particle1D_%s_%f.gif", fileName, upperPtLimit));
 559   canvas->SaveAs(Form("Track2Particle1D_%s_%f.eps", fileName, upperPtLimit));
 560
 561   //TPaveText* pave = new TPaveText(-0.4, 1.35, 0.4, 1.45);
 562
 563   canvasName.Form("Track2Particle1D_%s_etapt", folder);
 564   TCanvas* canvas = new TCanvas(canvasName, canvasName, 1000, 500);
 565   canvas->Divide(2, 1);
 566
 567   canvas->cd(1);
 568   InitPad();
 569   corrY->GetXaxis()->SetRangeUser(-0.99, 0.99);
 570   corrY->GetYaxis()->SetRangeUser(1, 1.5);
 571   corrY->GetYaxis()->SetTitleOffset(1.5);
 572   corrY->DrawCopy();
 573   TPaveText* pave = new TPaveText(0.3, 0.7, 0.7, 0.8, "NDC");
 574   pave->AddText("|z| < 10 cm");
 575   pave->AddText("0.3 GeV/c < p_{T} < 10 GeV/c");
 576   pave->Draw();
 577
 578   canvas->cd(2);
 579   InitPad();
 580   gPad->SetLogx();
 581   corrZ->GetYaxis()->SetRangeUser(1, 2.5);
 582   corrZ->GetXaxis()->SetRangeUser(0.101, upperPtLimit);
 583   corrZ->GetYaxis()->SetTitleOffset(1.5);
 584   corrZ->DrawCopy();
 585   pave = new TPaveText(0.5, 0.7, 0.8, 0.8, "NDC");
 586   pave->AddText("|z| < 10 cm");
 587   pave->AddText("|#eta| < 0.8");
 588   pave->Draw();
 589
 590   canvas->SaveAs(Form("Track2Particle1D_etapt_%s_%f.eps", fileName, upperPtLimit));
 591   canvas->SaveAs(Form("Track2Particle1D_etapt_%s_%f.gif", fileName, upperPtLimit));
 592 }
 593
 594 void CompareTrack2Particle1D(Float_t upperPtLimit = 9.9)
 595 {
 596   gSystem->Load("libPWG0base");
 597
 598   Track2Particle1DCreatePlots("correction_maponly-positive.root", "dndeta_correction", upperPtLimit);
 599
 600   TH1* posX = dynamic_cast<TH1*> (gROOT->FindObject("gene_nTrackToNPart_x_div_meas_nTrackToNPart_x")->Clone("pos_x"));
 601   TH1* posY = dynamic_cast<TH1*> (gROOT->FindObject("gene_nTrackToNPart_y_div_meas_nTrackToNPart_y")->Clone("pos_y"));
 602   TH1* posZ = dynamic_cast<TH1*> (gROOT->FindObject("gene_nTrackToNPart_z_div_meas_nTrackToNPart_z")->Clone("pos_z"));
 603
 604   Track2Particle1DCreatePlots("correction_maponly-negative.root", "dndeta_correction", upperPtLimit);
 605
 606   TH1* negX = dynamic_cast<TH1*> (gROOT->FindObject("gene_nTrackToNPart_x_div_meas_nTrackToNPart_x")->Clone("neg_x"));
 607   TH1* negY = dynamic_cast<TH1*> (gROOT->FindObject("gene_nTrackToNPart_y_div_meas_nTrackToNPart_y")->Clone("neg_y"));
 608   TH1* negZ = dynamic_cast<TH1*> (gROOT->FindObject("gene_nTrackToNPart_z_div_meas_nTrackToNPart_z")->Clone("neg_z"));
 609
 610   //printf("%f %f %f %f\n", posX->GetBinContent(20), posX->GetBinError(20), negX->GetBinContent(20), negX->GetBinError(20));
 611
 612   posX->Divide(negX);
 613   posY->Divide(negY);
 614   posZ->Divide(negZ);
 615
 616   //printf("%f %f\n", posX->GetBinContent(20), posX->GetBinError(20));
 617
 618   Prepare1DPlot(posX);
 619   Prepare1DPlot(posY);
 620   Prepare1DPlot(posZ);
 621
 622   Float_t min = 0.8;
 623   Float_t max = 1.2;
 624
 625   posX->SetMinimum(min);
 626   posX->SetMaximum(max);
 627   posY->SetMinimum(min);
 628   posY->SetMaximum(max);
 629   posZ->SetMinimum(min);
 630   posZ->SetMaximum(max);
 631
 632   posZ->GetXaxis()->SetRangeUser(0, upperPtLimit);
 633
 634   posX->GetYaxis()->SetTitleOffset(1.7);
 635   posX->GetYaxis()->SetTitle("C_{+} / C_{-}");
 636   posY->GetYaxis()->SetTitleOffset(1.7);
 637   posY->GetYaxis()->SetTitle("C_{+} / C_{-}");
 638   posZ->GetYaxis()->SetTitleOffset(1.7);
 639   posZ->GetYaxis()->SetTitle("C_{+} / C_{-}");
 640
 641   TCanvas* canvas = new TCanvas("CompareTrack2Particle1D", "CompareTrack2Particle1D", 1200, 400);
 642   canvas->Divide(3, 1);
 643
 644   canvas->cd(1);
 645   InitPad();
 646   posX->Draw();
 647
 648   canvas->cd(2);
 649   InitPad();
 650   posY->Draw();
 651
 652   canvas->cd(3);
 653   InitPad();
 654   posZ->Draw();
 655
 656   canvas->SaveAs(Form("CompareTrack2Particle1D_%f.gif", upperPtLimit));
 657   canvas->SaveAs(Form("CompareTrack2Particle1D_%f.eps", upperPtLimit));
 658 }
 659
 660 void Track2Particle2DCreatePlots(const char* fileName = "correction_map.root")
 661 {
 662   TFile::Open(fileName);
 663   AlidNdEtaCorrection* dNdEtaCorrection = new AlidNdEtaCorrection("dndeta_correction", "dndeta_correction");
 664   dNdEtaCorrection->LoadHistograms(fileName, "dndeta_correction");
 665
 666   TH3F* gene = dNdEtaCorrection->GetTrack2ParticleCorrection()->GetGeneratedHistogram();
 667   TH3F* meas = dNdEtaCorrection->GetTrack2ParticleCorrection()->GetMeasuredHistogram();
 668
 669   gene->GetZaxis()->SetRangeUser(0.3, 10);
 670   meas->GetZaxis()->SetRangeUser(0.3, 10);
 671   AliPWG0Helper::CreateDividedProjections(gene, meas, "yx");
 672   gene->GetZaxis()->SetRange(0, 0);
 673   meas->GetZaxis()->SetRange(0, 0);
 674
 675   gene->GetYaxis()->SetRangeUser(-0.8, 0.8);
 676   meas->GetYaxis()->SetRangeUser(-0.8, 0.8);
 677   AliPWG0Helper::CreateDividedProjections(gene, meas, "zx");
 678   gene->GetYaxis()->SetRange(0, 0);
 679   meas->GetYaxis()->SetRange(0, 0);
 680
 681   gene->GetXaxis()->SetRangeUser(-10, 10);
 682   meas->GetXaxis()->SetRangeUser(-10, 10);
 683   AliPWG0Helper::CreateDividedProjections(gene, meas, "zy");
 684   gene->GetXaxis()->SetRange(0, 0);
 685   meas->GetXaxis()->SetRange(0, 0);
 686 }
 687
 688 void Track2Particle2D(const char* fileName = "correction_map.root", const char* folder = "dndeta_correction")
 689 {
 690   gSystem->Load("libPWG0base");
 691
 692   Track2Particle2DCreatePlots(fileName);
 693
 694   TH2* corrYX = dynamic_cast<TH2*> (gROOT->FindObject(Form("gene_%s_nTrackToNPart_yx_div_meas_%s_nTrackToNPart_yx", folder, folder)));
 695   TH2* corrZX = dynamic_cast<TH2*> (gROOT->FindObject(Form("gene_%s_nTrackToNPart_zx_div_meas_%s_nTrackToNPart_zx", folder, folder)));
 696   TH2* corrZY = dynamic_cast<TH2*> (gROOT->FindObject(Form("gene_%s_nTrackToNPart_zy_div_meas_%s_nTrackToNPart_zy", folder, folder)));
 697
 698   /* this reads them from the file
 699   TH2* corrYX = dynamic_cast<TH2*> (file->Get("dndeta_correction/gene_nTrackToNPart_yx_div_meas_nTrackToNPart_yx"));
 700   TH2* corrZX = dynamic_cast<TH2*> (file->Get("dndeta_correction/gene_nTrackToNPart_zx_div_meas_nTrackToNPart_zx"));
 701   TH2* corrZY = dynamic_cast<TH2*> (file->Get("dndeta_correction/gene_nTrackToNPart_zy_div_meas_nTrackToNPart_zy"));*/
 702
 703   Prepare2DPlot(corrYX);
 704   Prepare2DPlot(corrZX);
 705   Prepare2DPlot(corrZY);
 706
 707   const char* title = "";
 708   corrYX->SetTitle(title);
 709   corrZX->SetTitle(title);
 710   corrZY->SetTitle(title);
 711
 712   TCanvas* canvas = new TCanvas("Track2Particle2D", "Track2Particle2D", 1200, 400);
 713   canvas->Divide(3, 1);
 714
 715   canvas->cd(1);
 716   InitPadCOLZ();
 717   corrYX->Draw("COLZ");
 718
 719   canvas->cd(2);
 720   InitPadCOLZ();
 721   corrZX->Draw("COLZ");
 722
 723   canvas->cd(3);
 724   InitPadCOLZ();
 725   corrZY->Draw("COLZ");
 726
 727   canvas->SaveAs(Form("Track2Particle2D_%d.gif", gMax));
 728   canvas->SaveAs(Form("Track2Particle2D_%d.eps", gMax));
 729 }
 730
 731 void CompareTrack2Particle2D()
 732 {
 733   gSystem->Load("libPWG0base");
 734
 735   Track2Particle2DCreatePlots("correction_maponly-positive.root");
 736
 737   TH2* posYX = dynamic_cast<TH2*> (gROOT->FindObject("gene_nTrackToNPart_yx_div_meas_nTrackToNPart_yx")->Clone("pos_yx"));
 738   TH2* posZX = dynamic_cast<TH2*> (gROOT->FindObject("gene_nTrackToNPart_zx_div_meas_nTrackToNPart_zx")->Clone("pos_zx"));
 739   TH2* posZY = dynamic_cast<TH2*> (gROOT->FindObject("gene_nTrackToNPart_zy_div_meas_nTrackToNPart_zy")->Clone("pos_zy"));
 740
 741   Track2Particle2DCreatePlots("correction_maponly-negative.root");
 742
 743   TH2* negYX = dynamic_cast<TH2*> (gROOT->FindObject("gene_nTrackToNPart_yx_div_meas_nTrackToNPart_yx")->Clone("neg_yx"));
 744   TH2* negZX = dynamic_cast<TH2*> (gROOT->FindObject("gene_nTrackToNPart_zx_div_meas_nTrackToNPart_zx")->Clone("neg_zx"));
 745   TH2* negZY = dynamic_cast<TH2*> (gROOT->FindObject("gene_nTrackToNPart_zy_div_meas_nTrackToNPart_zy")->Clone("neg_zy"));
 746
 747   posYX->Divide(negYX);
 748   posZX->Divide(negZX);
 749   posZY->Divide(negZY);
 750
 751   Prepare2DPlot(posYX);
 752   Prepare2DPlot(posZX);
 753   Prepare2DPlot(posZY);
 754
 755   Float_t min = 0.8;
 756   Float_t max = 1.2;
 757
 758   posYX->SetMinimum(min);
 759   posYX->SetMaximum(max);
 760   posZX->SetMinimum(min);
 761   posZX->SetMaximum(max);
 762   posZY->SetMinimum(min);
 763   posZY->SetMaximum(max);
 764
 765   TCanvas* canvas = new TCanvas("CompareTrack2Particle2D", "CompareTrack2Particle2D", 1200, 400);
 766   canvas->Divide(3, 1);
 767
 768   canvas->cd(1);
 769   InitPadCOLZ();
 770   posYX->Draw("COLZ");
 771
 772   canvas->cd(2);
 773   InitPadCOLZ();
 774   posZX->Draw("COLZ");
 775
 776   canvas->cd(3);
 777   InitPadCOLZ();
 778   posZY->Draw("COLZ");
 779
 780   canvas->SaveAs("CompareTrack2Particle2D.gif");
 781   canvas->SaveAs("CompareTrack2Particle2D.eps");
 782 }
 783
 784 void Track2Particle3D()
 785 {
 786   // get left margin proper
 787
 788   TFile* file = TFile::Open("correction_map.root");
 789
 790   TH3* corr = dynamic_cast<TH3*> (file->Get("dndeta_correction/corr_nTrackToNPart"));
 791
 792   corr->SetTitle("Correction Factor");
 793   SetRanges(corr->GetZaxis());
 794
 795   Prepare3DPlot(corr);
 796
 797   TCanvas* canvas = new TCanvas("Track2Particle3D", "Track2Particle3D", 500, 500);
 798   canvas->SetTheta(29.428);
 799   canvas->SetPhi(16.5726);
 800
 801   corr->Draw();
 802
 803   canvas->SaveAs("Track2Particle3D.gif");
 804   canvas->SaveAs("Track2Particle3D.eps");
 805 }
 806
 807 void Track2Particle3DAll()
 808 {
 809   TFile* file = TFile::Open("correction_map.root");
 810
 811   TH3* gene = dynamic_cast<TH3*> (file->Get("dndeta_correction/gene_nTrackToNPart"));
 812   TH3* meas = dynamic_cast<TH3*> (file->Get("dndeta_correction/meas_nTrackToNPart"));
 813   TH3* corr = dynamic_cast<TH3*> (file->Get("dndeta_correction/corr_nTrackToNPart"));
 814
 815   gene->SetTitle("Generated Particles");
 816   meas->SetTitle("Measured Tracks");
 817   corr->SetTitle("Correction Factor");
 818
 819   Prepare3DPlot(gene);
 820   Prepare3DPlot(meas);
 821   Prepare3DPlot(corr);
 822
 823   TCanvas* canvas = new TCanvas("Track2Particle3DAll", "Track2Particle3DAll", 1200, 400);
 824   canvas->Divide(3, 1);
 825
 826   canvas->cd(1);
 827   InitPad();
 828   gene->Draw();
 829
 830   canvas->cd(2);
 831   meas->Draw();
 832
 833   canvas->cd(3);
 834   corr->Draw();
 835
 836   canvas->SaveAs("Track2Particle3DAll.gif");
 837   canvas->SaveAs("Track2Particle3DAll.eps");
 838 }
 839
 840 void MultiplicityMC(Int_t xRangeMax = 50)
 841 {
 842   TFile* file = TFile::Open("multiplicityMC.root");
 843
 844   if (!file)
 845   {
 846     printf("multiplicityMC.root could not be opened.\n");
 847     return;
 848   }
 849
 850   TH1F* fMultiplicityESD = dynamic_cast<TH1F*> (file->Get("fMultiplicityESD"));
 851   TH1F* fMultiplicityMC = dynamic_cast<TH1F*> (file->Get("fMultiplicityMC"));
 852   TH2F* fCorrelation = dynamic_cast<TH2F*> (file->Get("fCorrelation"));
 853
 854   TH1F* correction = new TH1F("MultiplicityMC_correction", "MultiplicityMC_correction;Ntracks;Npart", 76, -0.5, 75.5);
 855   TH1F* correctionWidth = new TH1F("MultiplicityMC_correctionwidth", "MultiplicityMC_correctionwidth;Ntracks;Npart", 76, -0.5, 75.5);
 856   //fMultiplicityMC->GetNbinsX(), fMultiplicityMC->GetXaxis()->GetXmin(), fMultiplicityMC->GetXaxis()->GetXmax());
 857   for (Int_t i=1; i<=correction->GetNbinsX(); ++i)
 858   {
 859     TH1D* proj = fCorrelation->ProjectionX("_px", i, i+1);
 860     proj->Fit("gaus", "0");
 861     correction->SetBinContent(i, proj->GetFunction("gaus")->GetParameter(1));
 862     correctionWidth->SetBinContent(i, proj->GetFunction("gaus")->GetParameter(2));
 863
 864     continue;
 865
 866     // draw for debugging
 867     new TCanvas;
 868     proj->DrawCopy();
 869     proj->GetFunction("gaus")->DrawCopy("SAME");
 870   }
 871
 872   TH1F* fMultiplicityESDCorrected = new TH1F("fMultiplicityESDCorrected", "fMultiplicityESDCorrected", 2010, -0.5, 200.5);
 873
 874   for (Int_t i=1; i<=correction->GetNbinsX(); ++i)
 875   {
 876     Float_t mean = correction->GetBinContent(i);
 877     Float_t width = correctionWidth->GetBinContent(i);
 878
 879     Int_t fillBegin = fMultiplicityESDCorrected->FindBin(mean - width * 3);
 880     Int_t fillEnd   = fMultiplicityESDCorrected->FindBin(mean + width * 3);
 881     printf("bin %d mean %f width %f, filling from %d to %d\n", i, mean, width, fillBegin, fillEnd);
 882
 883     for (Int_t j=fillBegin; j <= fillEnd; ++j)
 884     {
 885       fMultiplicityESDCorrected->AddBinContent(j, TMath::Gaus(fMultiplicityESDCorrected->GetXaxis()->GetBinCenter(j), mean, width, kTRUE) * fMultiplicityESD->GetBinContent(i));
 886     }
 887   }
 888
 889   TH1F* fMultiplicityESDCorrectedRebinned = dynamic_cast<TH1F*> (fMultiplicityESDCorrected->Clone("fMultiplicityESDCorrectedRebinned"));
 890   fMultiplicityESDCorrectedRebinned->Rebin(10);
 891   fMultiplicityESDCorrectedRebinned->Scale(0.1);
 892
 893   TH1F* ratio = dynamic_cast<TH1F*> (fMultiplicityESD->Clone("multiplicity_ratio"));
 894   ratio->SetTitle("ratio;Ntracks;Nreco/Ngene");
 895   ratio->Divide(fMultiplicityMC);
 896
 897   TH1F* ratio2 = dynamic_cast<TH1F*> (fMultiplicityESDCorrectedRebinned->Clone("multiplicity_ratio_corrected"));
 898   ratio2->Divide(fMultiplicityMC);
 899
 900   TCanvas* canvas = new TCanvas("MultiplicityMC", "MultiplicityMC", 1500, 1000);
 901   canvas->Divide(3, 2);
 902
 903   fMultiplicityESD->GetXaxis()->SetRangeUser(0, xRangeMax);
 904   ratio->GetXaxis()->SetRangeUser(0, xRangeMax);
 905   fCorrelation->GetXaxis()->SetRangeUser(0, xRangeMax);
 906   fCorrelation->GetYaxis()->SetRangeUser(0, xRangeMax);
 907   correction->GetXaxis()->SetRangeUser(0, xRangeMax);
 908   fMultiplicityESDCorrected->GetXaxis()->SetRangeUser(0, xRangeMax);
 909   fMultiplicityESDCorrectedRebinned->GetXaxis()->SetRangeUser(0, xRangeMax);
 910
 911   canvas->cd(1); //InitPad();
 912   fMultiplicityESD->Draw();
 913   fMultiplicityMC->SetLineColor(2);
 914   fMultiplicityMC->Draw("SAME");
 915
 916   TLegend* legend = new TLegend(0.6, 0.7, 0.85, 0.85);
 917   legend->AddEntry(fMultiplicityESD, "ESD");
 918   legend->AddEntry(fMultiplicityMC, "MC");
 919   legend->Draw();
 920
 921   canvas->cd(2);
 922   fCorrelation->Draw("COLZ");
 923
 924   canvas->cd(3);
 925   correction->Draw();
 926   //correction->Fit("pol1");
 927   correctionWidth->SetLineColor(2);
 928   correctionWidth->Draw("SAME");
 929
 930   legend = new TLegend(0.2, 0.7, 0.45, 0.85);
 931   legend->AddEntry(correction, "#bar{x}");
 932   legend->AddEntry(correctionWidth, "#sigma");
 933   legend->Draw();
 934
 935   canvas->cd(4);
 936   ratio->Draw();
 937
 938   ratio2->SetLineColor(2);
 939   ratio2->Draw("SAME");
 940
 941   legend = new TLegend(0.6, 0.7, 0.85, 0.85);
 942   legend->AddEntry(ratio, "uncorrected");
 943   legend->AddEntry(ratio2, "corrected");
 944   legend->Draw();
 945
 946   canvas->cd(5);
 947   fMultiplicityESDCorrected->SetLineColor(kBlue);
 948   fMultiplicityESDCorrected->Draw();
 949   fMultiplicityMC->Draw("SAME");
 950   fMultiplicityESD->Draw("SAME");
 951
 952   legend = new TLegend(0.6, 0.7, 0.85, 0.85);
 953   legend->AddEntry(fMultiplicityESDCorrected, "ESD corrected");
 954   legend->AddEntry(fMultiplicityMC, "MC");
 955   legend->AddEntry(fMultiplicityESD, "ESD");
 956   legend->Draw();
 957
 958   canvas->cd(6);
 959   fMultiplicityESDCorrectedRebinned->SetLineColor(kBlue);
 960   fMultiplicityESDCorrectedRebinned->Draw();
 961   fMultiplicityMC->Draw("SAME");
 962
 963   legend = new TLegend(0.6, 0.7, 0.85, 0.85);
 964   legend->AddEntry(fMultiplicityESDCorrectedRebinned, "ESD corrected");
 965   legend->AddEntry(fMultiplicityMC, "MC");
 966   legend->Draw();
 967
 968   canvas->SaveAs("MultiplicityMC.gif");
 969 }
 970
 971 void MultiplicityESD()
 972 {
 973   TFile* file = TFile::Open("multiplicityESD.root");
 974
 975   if (!file)
 976   {
 977     printf("multiplicityESD.root could not be opened.\n");
 978     return;
 979   }
 980
 981   TH1F* fMultiplicityESD = dynamic_cast<TH1F*> (file->Get("fMultiplicity"));
 982
 983   TCanvas* canvas = new TCanvas("MultiplicityESD", "MultiplicityESD", 500, 500);
 984
 985   fMultiplicityESD->Draw();
 986 }
 987
 988 void drawPlots(Int_t max)
 989 {
 990   gMax = max;
 991
 992   ptCutoff();
 993   TriggerBias();
 994   VtxRecon();
 995   Track2Particle2D();
 996   Track2Particle3D();
 997   ptSpectrum();
 998   dNdEta();
 999 }
1000
1001 void drawPlots()
1002 {
1003   drawPlots(5);
1004   drawPlots(2);
1005 }