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