PWGLF/SPECTRA/PiKaPr/TOF/pPb502/macros/TOFmatchMC.C

   1 #include "CommonDefs.C"
   2
   3 const Char_t *destdir = "plots";
   4
   5 enum EHisto_t {
   6   kAcceptedTracks,
   7   kMatchedTracks,
   8   kMismatchedTracks,
   9   kUncorrelatedTracks,
  10   kMatchedCorrelatedTracks,
  11   kMatchedGoodTracks,
  12   kNHistos
  13 };
  14 const Char_t *histoName[kNHistos] = {
  15   "hAcceptedTracks",
  16   "hMatchedTracks",
  17   "hMismatchedTracks",
  18   "hUncorrelatedTracks",
  19   "hMatchedCorrelatedTracks",
  20   "hMatchedGoodTracks"
  21 };
  22
  23 enum EParam_t {
  24   kCentrality,
  25   kPt,
  26   kEta,
  27   kPhi,
  28   kNParams
  29 };
  30
  31 TOFmatchMC(const Char_t *filename, Int_t rapidityPart = -1, Int_t evMax = kMaxInt)
  32 {
  33
  34   /* include path for ACLic */
  35   gSystem->AddIncludePath("-I$ALICE_ROOT/include");
  36   gSystem->AddIncludePath("-I$ALICE_ROOT/TOF");
  37   /* load libraries */
  38   gSystem->Load("libANALYSIS");
  39   gSystem->Load("libANALYSISalice");
  40   /* build analysis task class */
  41   gROOT->LoadMacro("AliAnalysisParticle.cxx+g");
  42   gROOT->LoadMacro("AliAnalysisEvent.cxx+g");
  43   gROOT->LoadMacro("AliAnalysisTrack.cxx+g");
  44
  45   /* open file, get tree and connect */
  46   TFile *filein = TFile::Open(filename);
  47   TTree *treein = (TTree *)filein->Get("aodTree");
  48   printf("got \"aodTree\": %d entries\n", treein->GetEntries());
  49   AliAnalysisEvent *analysisEvent = new AliAnalysisEvent();
  50   TClonesArray *analysisTrackArray = new TClonesArray("AliAnalysisTrack");
  51   AliAnalysisTrack *analysisTrack = NULL;
  52   treein->SetBranchAddress("AnalysisEvent", &analysisEvent);
  53   treein->SetBranchAddress("AnalysisTrack", &analysisTrackArray);
  54
  55   /* open enabled flag map */
  56   TH1F *hEnabledFlag = NULL;
  57   if (enabledChannelsFileName) {
  58     TFile *enabledfile = TFile::Open(enabledChannelsFileName);
  59     hEnabledFlag = (TH1F *)enabledfile->Get("hEnabledFlag");
  60   }
  61
  62   /* binning */
  63   for (Int_t ieta = 0; ieta < NetaBins + 1; ieta++)
  64     etaBin[ieta] = etaMin + ieta * etaStep;
  65   for (Int_t iphi = 0; iphi < NphiBins + 1; iphi++)
  66     phiBin[iphi] = phiMin + iphi * phiStep;
  67   /* THnSparse */
  68   Int_t NparamBins[kNParams] = {NcentralityBins, NptBins, NetaBins, NphiBins};
  69   Double_t *paramBin[kNParams] = {centralityBin, ptBin, etaBin, phiBin};
  70   THnSparseF *hHisto[kNHistos][AliPID::kSPECIES][kNCharges];
  71   for (Int_t ihisto = 0; ihisto < kNHistos; ihisto++)
  72     for (Int_t ipart = 0; ipart < AliPID::kSPECIES; ipart++)
  73       for (Int_t icharge = 0; icharge < kNCharges; icharge++) {
  74         hHisto[ihisto][ipart][icharge] = new THnSparseF(Form("hHisto_%s_%s_%s", histoName[ihisto], AliPID::ParticleName(ipart), chargeName[icharge]), "", kNParams, NparamBins);
  75         for (Int_t iparam = 0; iparam < kNParams; iparam++)
  76           hHisto[ihisto][ipart][icharge]->SetBinEdges(iparam, paramBin[iparam]);
  77       }
  78   /* histos */
  79   TH2F *hFEAMap = new TH2F("hFEAMap", "", 72, 0., 18., 91, 0., 91.);
  80
  81   /* start stopwatch */
  82   TStopwatch timer;
  83   timer.Start();
  84
  85   /* loop over events */
  86   Double_t param[kNParams];
  87   Int_t part, charge;
  88   Int_t index, sector, sectorStrip, padx, fea;
  89   Float_t hitmapx, hitmapy;
  90   AliTOFcalibHisto calib;
  91   calib.LoadCalibHisto();
  92   for (Int_t iev = 0; iev < treein->GetEntries() && iev < evMax; iev++) {
  93     /* get event */
  94     treein->GetEvent(iev);
  95     if (iev % 1000 == 0) printf("iev = %d\n", iev);
  96     /* check event */
  97     if (!analysisEvent->AcceptEvent(acceptEventType)) continue;
  98
  99     /*** ACCEPTED EVENT ***/
 100
 101     /* get centrality */
 102     param[kCentrality] = analysisEvent->GetCentralityPercentile(centralityEstimator);
 103
 104     /* loop over tracks */
 105     for (Int_t itrk = 0; itrk < analysisTrackArray->GetEntries(); itrk++) {
 106       /* get track */
 107       analysisTrack = (AliAnalysisTrack *)analysisTrackArray->At(itrk);
 108       if (!analysisTrack) continue;
 109       /* check charged primary with defined PID */
 110       part = analysisTrack->GetMCPID();
 111       if (!analysisTrack->IsMCPrimary() || part < 0 || analysisTrack->GetSign() == 0.) continue;
 112       /* check accepted track */
 113       if (!analysisTrack->AcceptTrack()) continue;
 114       /* check rapidity */
 115       if (rapidityPart < 0) {
 116         if ((analysisTrack->GetY(AliPID::ParticleMass(part)) - rapidityShift) > rapidityMaxCut ||
 117             (analysisTrack->GetY(AliPID::ParticleMass(part)) - rapidityShift) < rapidityMinCut) continue;
 118       }
 119       else {
 120         if ((analysisTrack->GetY(AliPID::ParticleMass(rapidityPart)) - rapidityShift) > rapidityMaxCut ||
 121             (analysisTrack->GetY(AliPID::ParticleMass(rapidityPart)) - rapidityShift) < rapidityMinCut) continue;
 122       }
 123
 124       /*** ACCEPTED TRACK ***/
 125
 126       /* get track info */
 127       param[kPt] = analysisTrack->GetPt();
 128       param[kEta] = analysisTrack->GetEta();
 129       param[kPhi] = analysisTrack->GetPhi();
 130       charge = analysisTrack->GetSign() > 0. ? kPositive : kNegative;
 131
 132       /* fill accepted tracks histos */
 133       hHisto[kAcceptedTracks][part][charge]->Fill(param);
 134
 135       /* check TOF PID */
 136       if (!analysisTrack->HasTOFPID(hEnabledFlag)) continue;
 137
 138       /*** ACCEPTED TRACK WITH TOF SIGNAL ***/
 139
 140       /* fill FEA map */
 141       sector = calib.GetCalibMap(AliTOFcalibHisto::kSector, index);
 142       sectorStrip = calib.GetCalibMap(AliTOFcalibHisto::kSectorStrip, index);
 143       padx = calib.GetCalibMap(AliTOFcalibHisto::kPadX, index);
 144       fea = padx / 12;
 145       hitmapx = sector + ((Double_t)(3 - fea) + 0.5) / 4.;
 146       hitmapy = sectorStrip;
 147       hFEAMap->Fill(hitmapx, hitmapy);
 148
 149       /* fill matched tracks histos */
 150       hHisto[kMatchedTracks][part][charge]->Fill(param);
 151
 152       /* check mismatch */
 153       if (!analysisTrack->IsMismatchMC()) {
 154         hHisto[kMatchedGoodTracks][part][charge]->Fill(param);
 155         hHisto[kMatchedCorrelatedTracks][part][charge]->Fill(param);
 156         continue;
 157       }
 158
 159       /*** MIS-MATCHED TRACK ***/
 160
 161       /* fill mis-matched tracks histos */
 162       hHisto[kMismatchedTracks][part][charge]->Fill(param);
 163
 164       /* check uncorrelated mismatch */
 165       if (!analysisTrack->IsUncorrelatedMismatchMC()) {
 166         hHisto[kMatchedCorrelatedTracks][part][charge]->Fill(param);
 167         continue;
 168       }
 169
 170       /*** UNCORRELATED MIS-MATCHED TRACK ***/
 171
 172       /* fill uncorrelated mis-matched tracks histos */
 173       hHisto[kUncorrelatedTracks][part][charge]->Fill(param);
 174
 175     }
 176   }
 177
 178   /* stop stopwatch */
 179   timer.Stop();
 180   timer.Print();
 181
 182   /* output */
 183   if (rapidityPart < 0)
 184     TFile *fileout = TFile::Open(Form("TOFmatchMC.%s", filename), "RECREATE");
 185   else
 186     TFile *fileout = TFile::Open(Form("TOFmatchMC.%s.%s", AliPID::ParticleName(rapidityPart), filename), "RECREATE");
 187   hFEAMap->Write();
 188   for (Int_t ihisto = 0; ihisto < kNHistos; ihisto++)
 189     for (Int_t ipart = 0; ipart < AliPID::kSPECIES; ipart++)
 190       for (Int_t icharge = 0; icharge < kNCharges; icharge++)
 191         hHisto[ihisto][ipart][icharge]->Write();
 192   fileout->Close();
 193
 194   if (rapidityPart < 0)
 195     TString str = Form("TOFmatchMC.%s", filename);
 196   else
 197     TString str = Form("TOFmatchMC.%s.%s", AliPID::ParticleName(rapidityPart), filename);
 198   TOFmatchMC_efficiencyPt(str.Data());
 199 }
 200
 201 //_____________________________________________________________________________-
 202
 203 TOFmatchMC_efficiencyPt(const Char_t *filename)
 204 {
 205
 206   /* get data */
 207   TFile *filein = TFile::Open(filename);
 208   THnSparseF *hHisto[kNHistos][AliPID::kSPECIES][kNCharges];
 209   TH1D *hHistoPt_MB[kNHistos][AliPID::kSPECIES][kNCharges], *hHistoPt_centrality[NcentralityBins][kNHistos][AliPID::kSPECIES][kNCharges];
 210   TH1D *hHistoAllPt_MB[kNHistos], *hHistoAllPt_centrality[NcentralityBins][kNHistos];
 211   /* loop over histos */
 212   for (Int_t ihisto = 0; ihisto < kNHistos; ihisto++) {
 213
 214     /* INCLUSIVE */
 215
 216     hHistoAllPt_MB[ihisto] = new TH1D(Form("hHistoAllPt_MB_%s", histoName[ihisto]), "", NptBins, ptBin);
 217     for (Int_t icent = 0; icent < NcentralityBins; icent++)
 218       hHistoAllPt_centrality[icent][ihisto] = new TH1D(Form("hHistoAllPt_centrality%d_%s", icent, histoName[ihisto]), "", NptBins, ptBin);
 219
 220     /* SINGLE PARTICLE */
 221
 222     for (Int_t ipart = 0; ipart < AliPID::kSPECIES; ipart++) {
 223       for (Int_t icharge = 0; icharge < kNCharges; icharge++) {
 224
 225         /* get histo */
 226         hHisto[ihisto][ipart][icharge] = (THnSparseF *)filein->Get(Form("hHisto_%s_%s_%s", histoName[ihisto], AliPID::ParticleName(ipart), chargeName[icharge]));
 227
 228         /* MB projection */
 229         hHistoPt_MB[ihisto][ipart][icharge] = hHisto[ihisto][ipart][icharge]->Projection(kPt);
 230         hHistoPt_MB[ihisto][ipart][icharge]->SetName(Form("hHistoPt_MB_%s_%s_%s", histoName[ihisto], AliPID::ParticleName(ipart), chargeName[icharge]));
 231         hHistoPt_MB[ihisto][ipart][icharge]->Sumw2();
 232         hHistoAllPt_MB[ihisto]->Add(hHistoPt_MB[ihisto][ipart][icharge]);
 233
 234         /* centrality projection */
 235         for (Int_t icent = 0; icent < NcentralityBins; icent++) {
 236           hHisto[ihisto][ipart][icharge]->GetAxis(kCentrality)->SetRange(icent + 1, icent + 1);
 237           hHistoPt_centrality[icent][ihisto][ipart][icharge] = hHisto[ihisto][ipart][icharge]->Projection(kPt);
 238           hHistoPt_centrality[icent][ihisto][ipart][icharge]->SetName(Form("hHistoPt_centrality%d_%s_%s_%s", icent, histoName[ihisto], AliPID::ParticleName(ipart), chargeName[icharge]));
 239           hHistoPt_centrality[icent][ihisto][ipart][icharge]->Sumw2();
 240           hHistoAllPt_centrality[icent][ihisto]->Add(hHistoPt_centrality[icent][ihisto][ipart][icharge]);
 241         }
 242       }
 243     }
 244   }
 245
 246   /* output */
 247   TString str = filename;
 248   str.Insert(str.Length() - TString(".root").Length(), ".efficiencyPt");
 249   TFile *fileout = TFile::Open(str.Data(), "RECREATE");
 250
 251   /* efficiencies/fractions and write */
 252   TH1D *hEfficiencyPt_MB[kNHistos][AliPID::kSPECIES][kNCharges], *hEfficiencyPt_centrality[NcentralityBins][kNHistos][AliPID::kSPECIES][kNCharges], *hEfficiencyPt_ratioMB_centrality[NcentralityBins][kNHistos][AliPID::kSPECIES][kNCharges];
 253   TH1D *hEfficiencyAllPt_MB[kNHistos], *hEfficiencyAllPt_centrality[NcentralityBins][kNHistos], *hEfficiencyAllPt_ratioMB_centrality[NcentralityBins][kNHistos];
 254
 255
 256   TH1D *hFractionPt_MB[kNHistos][AliPID::kSPECIES][kNCharges], *hFractionPt_centrality[NcentralityBins][kNHistos][AliPID::kSPECIES][kNCharges], *hFractionPt_ratioMB_centrality[NcentralityBins][kNHistos][AliPID::kSPECIES][kNCharges];
 257   for (Int_t ihisto = 0; ihisto < kNHistos; ihisto++) {
 258
 259     if (ihisto == kAcceptedTracks) continue;
 260
 261     /* INCLUSIVE */
 262
 263     /* MB efficiency */
 264     hEfficiencyAllPt_MB[ihisto] = new TH1D(*hHistoAllPt_MB[ihisto]);
 265     hEfficiencyAllPt_MB[ihisto]->SetName(Form("hEfficiencyAllPt_MB_%s", histoName[ihisto]));
 266     hEfficiencyAllPt_MB[ihisto]->SetLineWidth(2);
 267     hEfficiencyAllPt_MB[ihisto]->SetLineColor(1);
 268     hEfficiencyAllPt_MB[ihisto]->SetMarkerStyle(20);
 269     hEfficiencyAllPt_MB[ihisto]->SetMarkerColor(1);
 270     hEfficiencyAllPt_MB[ihisto]->Divide(hEfficiencyAllPt_MB[ihisto], hHistoAllPt_MB[kAcceptedTracks], 1, 1, "B");
 271     hEfficiencyAllPt_MB[ihisto]->Write();
 272
 273     /* multiplicity/centrality efficiency */
 274     for (Int_t icent = 0; icent < NcentralityBins; icent++) {
 275       hEfficiencyAllPt_centrality[icent][ihisto] = new TH1D(*hHistoAllPt_centrality[icent][ihisto]);
 276       hEfficiencyAllPt_centrality[icent][ihisto]->SetName(Form("hEfficiencyAllPt_centrality%d_%s", icent, histoName[ihisto]));
 277       hEfficiencyAllPt_centrality[icent][ihisto]->SetLineWidth(2);
 278       hEfficiencyAllPt_centrality[icent][ihisto]->SetLineColor(multcentColor[icent]);
 279       hEfficiencyAllPt_centrality[icent][ihisto]->SetMarkerStyle(20);
 280       hEfficiencyAllPt_centrality[icent][ihisto]->SetMarkerColor(multcentColor[icent]);
 281       hEfficiencyAllPt_centrality[icent][ihisto]->Divide(hEfficiencyAllPt_centrality[icent][ihisto], hHistoAllPt_centrality[icent][kAcceptedTracks], 1, 1, "B");
 282       hEfficiencyAllPt_centrality[icent][ihisto]->Write();
 283
 284       /* ratio wrt. MB */
 285       hEfficiencyAllPt_ratioMB_centrality[icent][ihisto] = new TH1D(*hEfficiencyAllPt_centrality[icent][ihisto]);
 286       hEfficiencyAllPt_ratioMB_centrality[icent][ihisto]->SetName(Form("hEfficiencyAllPt_ratioMB_centrality%d_%s", icent, histoName[ihisto]));
 287       hEfficiencyAllPt_ratioMB_centrality[icent][ihisto]->SetLineWidth(2);
 288       hEfficiencyAllPt_ratioMB_centrality[icent][ihisto]->SetLineColor(multcentColor[icent]);
 289       hEfficiencyAllPt_ratioMB_centrality[icent][ihisto]->SetMarkerStyle(20);
 290       hEfficiencyAllPt_ratioMB_centrality[icent][ihisto]->SetMarkerColor(multcentColor[icent]);
 291       hEfficiencyAllPt_ratioMB_centrality[icent][ihisto]->Divide(hEfficiencyAllPt_MB[ihisto]);
 292       hEfficiencyAllPt_ratioMB_centrality[icent][ihisto]->Write();
 293     }
 294
 295     /* SINGLE PARTICLE */
 296
 297     for (Int_t ipart = 0; ipart < AliPID::kSPECIES; ipart++) {
 298       for (Int_t icharge = 0; icharge < kNCharges; icharge++) {
 299
 300         /* MB efficiency */
 301         hEfficiencyPt_MB[ihisto][ipart][icharge] = new TH1D(*hHistoPt_MB[ihisto][ipart][icharge]);
 302         hEfficiencyPt_MB[ihisto][ipart][icharge]->SetName(Form("hEfficiencyPt_MB_%s_%s_%s", histoName[ihisto], AliPID::ParticleName(ipart), chargeName[icharge]));
 303         hEfficiencyPt_MB[ihisto][ipart][icharge]->SetLineWidth(2);
 304         hEfficiencyPt_MB[ihisto][ipart][icharge]->SetLineColor(1);
 305         hEfficiencyPt_MB[ihisto][ipart][icharge]->SetMarkerStyle(20);
 306         hEfficiencyPt_MB[ihisto][ipart][icharge]->SetMarkerColor(1);
 307         hEfficiencyPt_MB[ihisto][ipart][icharge]->Divide(hEfficiencyPt_MB[ihisto][ipart][icharge], hHistoPt_MB[kAcceptedTracks][ipart][icharge], 1, 1, "B");
 308         hEfficiencyPt_MB[ihisto][ipart][icharge]->Write();
 309
 310         /* multiplicity/centrality efficiency */
 311         for (Int_t icent = 0; icent < NcentralityBins; icent++) {
 312           hEfficiencyPt_centrality[icent][ihisto][ipart][icharge] = new TH1D(*hHistoPt_centrality[icent][ihisto][ipart][icharge]);
 313           hEfficiencyPt_centrality[icent][ihisto][ipart][icharge]->SetName(Form("hEfficiencyPt_centrality%d_%s_%s_%s", icent, histoName[ihisto], AliPID::ParticleName(ipart), chargeName[icharge]));
 314           hEfficiencyPt_centrality[icent][ihisto][ipart][icharge]->SetLineWidth(2);
 315           hEfficiencyPt_centrality[icent][ihisto][ipart][icharge]->SetLineColor(multcentColor[icent]);
 316           hEfficiencyPt_centrality[icent][ihisto][ipart][icharge]->SetMarkerStyle(20);
 317           hEfficiencyPt_centrality[icent][ihisto][ipart][icharge]->SetMarkerColor(multcentColor[icent]);
 318           hEfficiencyPt_centrality[icent][ihisto][ipart][icharge]->Divide(hEfficiencyPt_centrality[icent][ihisto][ipart][icharge], hHistoPt_centrality[icent][kAcceptedTracks][ipart][icharge], 1, 1, "B");
 319           hEfficiencyPt_centrality[icent][ihisto][ipart][icharge]->Write();
 320
 321           /* ratio wrt. MB */
 322           hEfficiencyPt_ratioMB_centrality[icent][ihisto][ipart][icharge] = new TH1D(*hEfficiencyPt_centrality[icent][ihisto][ipart][icharge]);
 323           hEfficiencyPt_ratioMB_centrality[icent][ihisto][ipart][icharge]->SetName(Form("hEfficiencyPt_ratioMB_centrality%d_%s_%s_%s", icent, histoName[ihisto], AliPID::ParticleName(ipart), chargeName[icharge]));
 324           hEfficiencyPt_ratioMB_centrality[icent][ihisto][ipart][icharge]->SetLineWidth(2);
 325           hEfficiencyPt_ratioMB_centrality[icent][ihisto][ipart][icharge]->SetLineColor(multcentColor[icent]);
 326           hEfficiencyPt_ratioMB_centrality[icent][ihisto][ipart][icharge]->SetMarkerStyle(20);
 327           hEfficiencyPt_ratioMB_centrality[icent][ihisto][ipart][icharge]->SetMarkerColor(multcentColor[icent]);
 328           hEfficiencyPt_ratioMB_centrality[icent][ihisto][ipart][icharge]->Divide(hEfficiencyPt_MB[ihisto][ipart][icharge]);
 329           hEfficiencyPt_ratioMB_centrality[icent][ihisto][ipart][icharge]->Write();
 330         }
 331
 332         /* fractions */
 333
 334         if (ihisto == kAcceptedTracks || ihisto == kMatchedTracks) continue;
 335
 336         /* MB fraction */
 337         hFractionPt_MB[ihisto][ipart][icharge] = new TH1D(*hHistoPt_MB[ihisto][ipart][icharge]);
 338         hFractionPt_MB[ihisto][ipart][icharge]->SetName(Form("hFractionPt_MB_%s_%s_%s", histoName[ihisto], AliPID::ParticleName(ipart), chargeName[icharge]));
 339         hFractionPt_MB[ihisto][ipart][icharge]->SetLineWidth(2);
 340         hFractionPt_MB[ihisto][ipart][icharge]->SetLineColor(1);
 341         hFractionPt_MB[ihisto][ipart][icharge]->SetMarkerStyle(20);
 342         hFractionPt_MB[ihisto][ipart][icharge]->SetMarkerColor(1);
 343         hFractionPt_MB[ihisto][ipart][icharge]->Divide(hFractionPt_MB[ihisto][ipart][icharge], hHistoPt_MB[kMatchedTracks][ipart][icharge], 1, 1, "B");
 344         hFractionPt_MB[ihisto][ipart][icharge]->Write();
 345
 346         /* multiplicity/centrality fraction */
 347         for (Int_t icent = 0; icent < NcentralityBins; icent++) {
 348           hFractionPt_centrality[icent][ihisto][ipart][icharge] = new TH1D(*hHistoPt_centrality[icent][ihisto][ipart][icharge]);
 349           hFractionPt_centrality[icent][ihisto][ipart][icharge]->SetName(Form("hFractionPt_centrality%d_%s_%s_%s", icent, histoName[ihisto], AliPID::ParticleName(ipart), chargeName[icharge]));
 350           hFractionPt_centrality[icent][ihisto][ipart][icharge]->SetLineWidth(2);
 351           hFractionPt_centrality[icent][ihisto][ipart][icharge]->SetLineColor(multcentColor[icent]);
 352           hFractionPt_centrality[icent][ihisto][ipart][icharge]->SetMarkerStyle(20);
 353           hFractionPt_centrality[icent][ihisto][ipart][icharge]->SetMarkerColor(multcentColor[icent]);
 354           hFractionPt_centrality[icent][ihisto][ipart][icharge]->Divide(hFractionPt_centrality[icent][ihisto][ipart][icharge], hHistoPt_centrality[icent][kMatchedTracks][ipart][icharge], 1, 1, "B");
 355           hFractionPt_centrality[icent][ihisto][ipart][icharge]->Write();
 356
 357           /* ratio wrt. MB */
 358           hFractionPt_ratioMB_centrality[icent][ihisto][ipart][icharge] = new TH1D(*hFractionPt_centrality[icent][ihisto][ipart][icharge]);
 359           hFractionPt_ratioMB_centrality[icent][ihisto][ipart][icharge]->SetName(Form("hFractionPt_ratioMB_centrality%d_%s_%s_%s", icent, histoName[ihisto], AliPID::ParticleName(ipart), chargeName[icharge]));
 360           hFractionPt_ratioMB_centrality[icent][ihisto][ipart][icharge]->SetLineWidth(2);
 361           hFractionPt_ratioMB_centrality[icent][ihisto][ipart][icharge]->SetLineColor(multcentColor[icent]);
 362           hFractionPt_ratioMB_centrality[icent][ihisto][ipart][icharge]->SetMarkerStyle(20);
 363           hFractionPt_ratioMB_centrality[icent][ihisto][ipart][icharge]->SetMarkerColor(multcentColor[icent]);
 364           hFractionPt_ratioMB_centrality[icent][ihisto][ipart][icharge]->Divide(hFractionPt_MB[ihisto][ipart][icharge]);
 365           hFractionPt_ratioMB_centrality[icent][ihisto][ipart][icharge]->Write();
 366         }
 367       }
 368     }
 369   }
 370
 371   fileout->Close();
 372
 373   TOFmatchMC_matchingEfficiency(str.Data());
 374   TOFmatchMC_matchingEfficiency_mismatchCorrected(str.Data());
 375
 376 }
 377
 378 //_____________________________________________________________________________-
 379
 380 TOFmatchMC_efficiencyEta(const Char_t *filename)
 381 {
 382
 383   /* get data */
 384   TFile *filein = TFile::Open(filename);
 385   THnSparseF *hHisto[kNHistos][AliPID::kSPECIES][kNCharges];
 386   TH1D *hHistoEta_MB[kNHistos][AliPID::kSPECIES][kNCharges], *hHistoEta_centrality[NcentralityBins][kNHistos][AliPID::kSPECIES][kNCharges];
 387   TH1D *hHistoEta_MB_pt[NptsubBins][kNHistos][AliPID::kSPECIES][kNCharges], *hHistoEta_centrality_pt[NcentralityBins][NptsubBins][kNHistos][AliPID::kSPECIES][kNCharges];
 388   /* loop over histos */
 389   for (Int_t ihisto = 0; ihisto < kNHistos; ihisto++)
 390     for (Int_t ipart = 0; ipart < AliPID::kSPECIES; ipart++)
 391       for (Int_t icharge = 0; icharge < kNCharges; icharge++) {
 392
 393         /* get histo */
 394         hHisto[ihisto][ipart][icharge] = (THnSparseF *)filein->Get(Form("hHisto_%s_%s_%s", histoName[ihisto], AliPID::ParticleName(ipart), chargeName[icharge]));
 395
 396         /* MB projection */
 397         hHisto[ihisto][ipart][icharge]->GetAxis(kPt)->SetRange(0, 0);
 398         hHisto[ihisto][ipart][icharge]->GetAxis(kCentrality)->SetRange(0, 0);
 399         hHistoEta_MB[ihisto][ipart][icharge] = hHisto[ihisto][ipart][icharge]->Projection(kEta);
 400         hHistoEta_MB[ihisto][ipart][icharge]->SetName(Form("hHistoEta_MB_%s_%s_%s", histoName[ihisto], AliPID::ParticleName(ipart), chargeName[icharge]));
 401         hHistoEta_MB[ihisto][ipart][icharge]->Sumw2();
 402         /* pt bins */
 403         for (Int_t ipt = 0; ipt < NptsubBins; ipt++) {
 404           hHisto[ihisto][ipart][icharge]->GetAxis(kPt)->SetRange(ptsubBinMin[ipt] + 1, ptsubBinMax[ipt] + 1);
 405           hHisto[ihisto][ipart][icharge]->GetAxis(kCentrality)->SetRange(0, 0);
 406           hHistoEta_MB_pt[ipt][ihisto][ipart][icharge] = hHisto[ihisto][ipart][icharge]->Projection(kEta);
 407           hHistoEta_MB_pt[ipt][ihisto][ipart][icharge]->SetName(Form("hHistoEta_MB_pt%d_%s_%s_%s", ipt, histoName[ihisto], AliPID::ParticleName(ipart), chargeName[icharge]));
 408           hHistoEta_MB_pt[ipt][ihisto][ipart][icharge]->Sumw2();
 409         }
 410
 411         /* centrality projection */
 412         for (Int_t icent = 0; icent < NcentralityBins; icent++) {
 413           hHisto[ihisto][ipart][icharge]->GetAxis(kPt)->SetRange(0, 0);
 414           hHisto[ihisto][ipart][icharge]->GetAxis(kCentrality)->SetRange(icent + 1, icent + 1);
 415           hHistoEta_centrality[icent][ihisto][ipart][icharge] = hHisto[ihisto][ipart][icharge]->Projection(kEta);
 416           hHistoEta_centrality[icent][ihisto][ipart][icharge]->SetName(Form("hHistoEta_centrality%d_%s_%s_%s", icent, histoName[ihisto], AliPID::ParticleName(ipart), chargeName[icharge]));
 417           hHistoEta_centrality[icent][ihisto][ipart][icharge]->Sumw2();
 418           /* pt bins */
 419           for (Int_t ipt = 0; ipt < NptsubBins; ipt++) {
 420             hHisto[ihisto][ipart][icharge]->GetAxis(kPt)->SetRange(ptsubBinMin[ipt] + 1, ptsubBinMax[ipt] + 1);
 421             hHisto[ihisto][ipart][icharge]->GetAxis(kCentrality)->SetRange(icent + 1, icent + 1);
 422             hHistoEta_centrality_pt[icent][ipt][ihisto][ipart][icharge] = hHisto[ihisto][ipart][icharge]->Projection(kEta);
 423             hHistoEta_centrality_pt[icent][ipt][ihisto][ipart][icharge]->SetName(Form("hHistoEta_centrality%d_pt%d_%s_%s_%s", icent, ipt, histoName[ihisto], AliPID::ParticleName(ipart), chargeName[icharge]));
 424             hHistoEta_centrality_pt[icent][ipt][ihisto][ipart][icharge]->Sumw2();
 425           }
 426         }
 427       }
 428
 429   /* output */
 430   TString str = filename;
 431   str.Insert(str.Length() - TString(".root").Length(), ".efficiencyEta");
 432   TFile *fileout = TFile::Open(str.Data(), "RECREATE");
 433
 434   /* efficiencies/fractions and write */
 435   TH1D *hEfficiencyEta_MB[kNHistos][AliPID::kSPECIES][kNCharges], *hEfficiencyEta_centrality[NcentralityBins][kNHistos][AliPID::kSPECIES][kNCharges], *hEfficiencyEta_ratioMB_centrality[NcentralityBins][kNHistos][AliPID::kSPECIES][kNCharges];
 436   TH1D *hEfficiencyEta_MB_pt[NptsubBins][kNHistos][AliPID::kSPECIES][kNCharges], *hEfficiencyEta_centrality_pt[NcentralityBins][NptsubBins][kNHistos][AliPID::kSPECIES][kNCharges], *hEfficiencyEta_ratioMB_centrality_pt[NcentralityBins][NptsubBins][kNHistos][AliPID::kSPECIES][kNCharges];
 437
 438
 439   TH1D *hFractionEta_MB[kNHistos][AliPID::kSPECIES][kNCharges], *hFractionEta_centrality[NcentralityBins][kNHistos][AliPID::kSPECIES][kNCharges], *hFractionEta_ratioMB_centrality[NcentralityBins][kNHistos][AliPID::kSPECIES][kNCharges];
 440   for (Int_t ihisto = 0; ihisto < kNHistos; ihisto++) {
 441     for (Int_t ipart = 0; ipart < AliPID::kSPECIES; ipart++)
 442       for (Int_t icharge = 0; icharge < kNCharges; icharge++) {
 443
 444         if (ihisto == kAcceptedTracks) continue;
 445
 446         /* MB efficiency */
 447         hEfficiencyEta_MB[ihisto][ipart][icharge] = new TH1D(*hHistoEta_MB[ihisto][ipart][icharge]);
 448         hEfficiencyEta_MB[ihisto][ipart][icharge]->SetName(Form("hEfficiencyEta_MB_%s_%s_%s", histoName[ihisto], AliPID::ParticleName(ipart), chargeName[icharge]));
 449         hEfficiencyEta_MB[ihisto][ipart][icharge]->SetLineWidth(2);
 450         hEfficiencyEta_MB[ihisto][ipart][icharge]->SetLineColor(1);
 451         hEfficiencyEta_MB[ihisto][ipart][icharge]->SetMarkerStyle(20);
 452         hEfficiencyEta_MB[ihisto][ipart][icharge]->SetMarkerColor(1);
 453         hEfficiencyEta_MB[ihisto][ipart][icharge]->Divide(hHistoEta_MB[kAcceptedTracks][ipart][icharge]);
 454         hEfficiencyEta_MB[ihisto][ipart][icharge]->Write();
 455         /* pt bins */
 456         for (Int_t ipt = 0; ipt < NptsubBins; ipt++) {
 457           hEfficiencyEta_MB_pt[ipt][ihisto][ipart][icharge] = new TH1D(*hHistoEta_MB_pt[ipt][ihisto][ipart][icharge]);
 458           hEfficiencyEta_MB_pt[ipt][ihisto][ipart][icharge]->SetName(Form("hEfficiencyEta_MB_pt%d_%s_%s_%s", ipt, histoName[ihisto], AliPID::ParticleName(ipart), chargeName[icharge]));
 459           hEfficiencyEta_MB_pt[ipt][ihisto][ipart][icharge]->SetLineWidth(2);
 460           hEfficiencyEta_MB_pt[ipt][ihisto][ipart][icharge]->SetLineColor(1);
 461           hEfficiencyEta_MB_pt[ipt][ihisto][ipart][icharge]->SetMarkerStyle(20);
 462           hEfficiencyEta_MB_pt[ipt][ihisto][ipart][icharge]->SetMarkerColor(1);
 463           hEfficiencyEta_MB_pt[ipt][ihisto][ipart][icharge]->Divide(hHistoEta_MB_pt[ipt][kAcceptedTracks][ipart][icharge]);
 464           hEfficiencyEta_MB_pt[ipt][ihisto][ipart][icharge]->Write();
 465         }
 466
 467         /* multiplicity/centrality efficiency */
 468         for (Int_t icent = 0; icent < NcentralityBins; icent++) {
 469           hEfficiencyEta_centrality[icent][ihisto][ipart][icharge] = new TH1D(*hHistoEta_centrality[icent][ihisto][ipart][icharge]);
 470           hEfficiencyEta_centrality[icent][ihisto][ipart][icharge]->SetName(Form("hEfficiencyEta_centrality%d_%s_%s_%s", icent, histoName[ihisto], AliPID::ParticleName(ipart), chargeName[icharge]));
 471           hEfficiencyEta_centrality[icent][ihisto][ipart][icharge]->SetLineWidth(2);
 472           hEfficiencyEta_centrality[icent][ihisto][ipart][icharge]->SetLineColor(multcentColor[icent]);
 473           hEfficiencyEta_centrality[icent][ihisto][ipart][icharge]->SetMarkerStyle(20);
 474           hEfficiencyEta_centrality[icent][ihisto][ipart][icharge]->SetMarkerColor(multcentColor[icent]);
 475           hEfficiencyEta_centrality[icent][ihisto][ipart][icharge]->Divide(hHistoEta_centrality[icent][kAcceptedTracks][ipart][icharge]);
 476           hEfficiencyEta_centrality[icent][ihisto][ipart][icharge]->Write();
 477
 478           /* ratio wrt. MB */
 479           hEfficiencyEta_ratioMB_centrality[icent][ihisto][ipart][icharge] = new TH1D(*hEfficiencyEta_centrality[icent][ihisto][ipart][icharge]);
 480           hEfficiencyEta_ratioMB_centrality[icent][ihisto][ipart][icharge]->SetName(Form("hEfficiencyEta_ratioMB_centrality%d_%s_%s_%s", icent, histoName[ihisto], AliPID::ParticleName(ipart), chargeName[icharge]));
 481           hEfficiencyEta_ratioMB_centrality[icent][ihisto][ipart][icharge]->SetLineWidth(2);
 482           hEfficiencyEta_ratioMB_centrality[icent][ihisto][ipart][icharge]->SetLineColor(multcentColor[icent]);
 483           hEfficiencyEta_ratioMB_centrality[icent][ihisto][ipart][icharge]->SetMarkerStyle(20);
 484           hEfficiencyEta_ratioMB_centrality[icent][ihisto][ipart][icharge]->SetMarkerColor(multcentColor[icent]);
 485           hEfficiencyEta_ratioMB_centrality[icent][ihisto][ipart][icharge]->Divide(hEfficiencyEta_MB[ihisto][ipart][icharge]);
 486           hEfficiencyEta_ratioMB_centrality[icent][ihisto][ipart][icharge]->Write();
 487         }
 488
 489       }
 490   }
 491
 492   fileout->Close();
 493
 494 }
 495
 496 //_____________________________________________________________________________-
 497
 498 TOFmatchMC_centralityDependence(const Char_t *filename)
 499 {
 500
 501   Double_t fitMin[AliPID::kSPECIES] = {0.5, 0.5, 0.5, 0.5, 0.5};
 502   Double_t fitMax[AliPID::kSPECIES] = {3.0, 3.0, 3.0, 3.0, 5.0};
 503
 504   TF1 *pol0 = (TF1 *)gROOT->GetFunction("pol0");
 505   TF1 *pol1 = (TF1 *)gROOT->GetFunction("pol1");
 506   pol0->SetRange(0., 5.0);
 507   TFile *filein = TFile::Open(filename);
 508
 509   /* output */
 510   TString str = filename;
 511   str.Insert(str.Length() - TString(".root").Length(), ".centralityDependence");
 512   TFile *fileout = TFile::Open(str.Data(), "RECREATE");
 513
 514   TH1D *hEfficiencyPt_ratioMB;
 515   TH1D *hEfficiencyCentrality[kNHistos][AliPID::kSPECIES][kNCharges];
 516   TH1D *hChi2Centrality[kNHistos][AliPID::kSPECIES][kNCharges];
 517   TH1D *hEfficiencyAllCentrality[kNHistos];
 518   TH1F *hHistoEffRatio = new TH1F("hHistoEffRatio", "", 100, 0.5, 1.5);
 519   for (Int_t ihisto = 0; ihisto < kNHistos; ihisto++) {
 520     if (ihisto == kAcceptedTracks) continue;
 521
 522     /* INCLUSIVE */
 523
 524     hEfficiencyAllCentrality[ihisto] = new TH1D(Form("hEfficiencyAllCentrality_ratioMB_%s", histoName[ihisto]), "", NcentralityBins, centralityBin);
 525     hEfficiencyAllCentrality[ihisto]->SetLineWidth(2);
 526     hEfficiencyAllCentrality[ihisto]->SetLineColor(1);
 527     hEfficiencyAllCentrality[ihisto]->SetMarkerStyle(20);
 528     hEfficiencyAllCentrality[ihisto]->SetMarkerColor(1);
 529
 530     for (Int_t icent = 0; icent < NcentralityBins; icent++) {
 531
 532       hEfficiencyPt_ratioMB = (TH1D *)filein->Get(Form("hEfficiencyAllPt_ratioMB_centrality%d_%s", icent, histoName[ihisto]));
 533       hEfficiencyPt_ratioMB->Fit(pol0, "q0", "", 0., 5.0);
 534       hEfficiencyAllCentrality[ihisto]->SetBinContent(icent + 1, pol0->GetParameter(0));
 535       hEfficiencyAllCentrality[ihisto]->SetBinError(icent + 1, pol0->GetParError(0));
 536       hEfficiencyPt_ratioMB->Add(pol0, -1.);
 537       hEfficiencyPt_ratioMB->Write();
 538
 539     }
 540     hEfficiencyAllCentrality[ihisto]->Write();
 541
 542
 543     /* SINGLE PARTICLE */
 544
 545     for (Int_t ipart = 2; ipart < AliPID::kSPECIES; ipart++)
 546       for (Int_t icharge = 0; icharge < kNCharges; icharge++) {
 547
 548         hEfficiencyCentrality[ihisto][ipart][icharge] = new TH1D(Form("hEfficiencyCentrality_ratioMB_%s_%s_%s", histoName[ihisto], AliPID::ParticleName(ipart), chargeName[icharge]), "", NcentralityBins, centralityBin);
 549         hEfficiencyCentrality[ihisto][ipart][icharge]->SetLineWidth(2);
 550         hEfficiencyCentrality[ihisto][ipart][icharge]->SetLineColor(particleColor[ipart]);
 551         hEfficiencyCentrality[ihisto][ipart][icharge]->SetMarkerStyle(chargeMarker[icharge]);
 552         hEfficiencyCentrality[ihisto][ipart][icharge]->SetMarkerColor(particleColor[ipart]);
 553
 554         for (Int_t icent = 0; icent < NcentralityBins; icent++) {
 555
 556           hEfficiencyPt_ratioMB = (TH1D *)filein->Get(Form("hEfficiencyPt_ratioMB_centrality%d_%s_%s_%s", icent, histoName[ihisto], AliPID::ParticleName(ipart), chargeName[icharge]));
 557
 558           hEfficiencyPt_ratioMB->Fit(pol0, "q0", "", fitMin[ipart], fitMax[ipart]);
 559           hEfficiencyCentrality[ihisto][ipart][icharge]->SetBinContent(icent + 1, pol0->GetParameter(0));
 560           hEfficiencyCentrality[ihisto][ipart][icharge]->SetBinError(icent + 1, pol0->GetParError(0));
 561           pol0->SetRange(fitMin[ipart], fitMax[ipart]);
 562           hEfficiencyPt_ratioMB->Add(pol0, -1.);
 563           hEfficiencyPt_ratioMB->Write();
 564
 565         }
 566
 567         hEfficiencyCentrality[ihisto][ipart][icharge]->Write();
 568       }
 569   }
 570
 571   fileout->Close();
 572 }
 573
 574 //_____________________________________________________________________________-
 575
 576 TOFmatchMC_centralityDependenceFit(const Char_t *filename, Int_t ihisto = kMatchedGoodTracks)
 577 {
 578
 579   const Char_t *particleLabel[5][2] = {"", "", "", "", "#pi^{+}", "#pi^{-}", "K^{+}", "K^{-}", "p", "#bar{p}"};
 580
 581   TF1 *fCentralityDependence = new TF1("fCentralityDependence", "[0] + [1] * (1. - TMath::Exp(-x / [2]))", 0., 90.);
 582   fCentralityDependence->SetParameter(0, 0.98);
 583   fCentralityDependence->SetParameter(1, 0.05);
 584   fCentralityDependence->SetParameter(2, 50.);
 585
 586   TFile *filein = TFile::Open(filename);
 587   TH1D *hEfficiencyAllCentrality;
 588   TH1D *hEfficiencyCentrality[AliPID::kSPECIES][kNCharges];
 589   hEfficiencyAllCentrality = (TH1D *)filein->Get(Form("hEfficiencyAllCentrality_ratioMB_%s", histoName[ihisto]));
 590   TCanvas *cFit = new TCanvas("cFit");
 591   hEfficiencyAllCentrality->Fit(fCentralityDependence);
 592   hEfficiencyAllCentrality->SetMaximum(1.05);
 593   hEfficiencyAllCentrality->SetMinimum(0.95);
 594   TCanvas *cRatios = new TCanvas("cRatios");
 595   cRatios->Divide(2, 3);
 596   for (Int_t ipart = 2; ipart < AliPID::kSPECIES; ipart++) {
 597     for (Int_t icharge = 0; icharge < 2; icharge++) {
 598       cRatios->cd(icharge + 1 + 2 * (ipart - 2));
 599       cRatios->cd(icharge + 1 + 2 * (ipart - 2))->SetGridx();
 600       cRatios->cd(icharge + 1 + 2 * (ipart - 2))->SetGridy();
 601       hEfficiencyCentrality[ipart][icharge] = (TH1D *)filein->Get(Form("hEfficiencyCentrality_ratioMB_%s_%s_%s", histoName[ihisto], AliPID::ParticleName(ipart), chargeName[icharge]));
 602       hEfficiencyCentrality[ipart][icharge]->Divide(fCentralityDependence);
 603       hEfficiencyCentrality[ipart][icharge]->SetMaximum(1.05);
 604       hEfficiencyCentrality[ipart][icharge]->SetMinimum(0.95);
 605       hEfficiencyCentrality[ipart][icharge]->SetTitle(Form("%s;centrality percentile;efficiency ratio wrt. fitted centrality dependence", particleLabel[ipart][icharge]));
 606       hEfficiencyCentrality[ipart][icharge]->SetStats(kFALSE);
 607       hEfficiencyCentrality[ipart][icharge]->Draw();
 608     }
 609   }
 610
 611
 612 }
 613
 614 //_____________________________________________________________________________-
 615
 616 TOFmatchMC_efficiencyPt_MB_plot(const Char_t *filename)
 617 {
 618   TOFmatchMC_efficiencyPt_MB_plot(filename, kMatchedGoodTracks, 2, 0, 20, 4);
 619   TOFmatchMC_efficiencyPt_MB_plot(filename, kMatchedGoodTracks, 2, 1, 25, 4);
 620   TOFmatchMC_efficiencyPt_MB_plot(filename, kMatchedGoodTracks, 3, 0, 20, 8);
 621   TOFmatchMC_efficiencyPt_MB_plot(filename, kMatchedGoodTracks, 3, 1, 25, 8);
 622   TOFmatchMC_efficiencyPt_MB_plot(filename, kMatchedGoodTracks, 4, 0, 20, 2);
 623   TOFmatchMC_efficiencyPt_MB_plot(filename, kMatchedGoodTracks, 4, 1, 25, 2);
 624 }
 625
 626 TH1D *
 627 TOFmatchMC_efficiencyPt_MB_plot(const Char_t *filename, Int_t ihisto = kMatchedGoodTracks, Int_t ipart, Int_t icharge, Int_t marker = 20, Int_t color = 2, Option_t *opt = "")
 628 {
 629
 630   TCanvas *cCanvas1 = new TCanvas("cCanvas1");
 631   TCanvas *cCanvas2 = new TCanvas("cCanvas2");
 632
 633   Double_t ptMin[AliPID::kSPECIES] = {0.5, 0.5, 0.5, 0.5, 0.5};
 634   Double_t ptMax[AliPID::kSPECIES] = {3.0, 3.0, 3.0, 3.0, 5.0};
 635
 636   TF1 *fEff = new TF1("fEff", "[0] * TMath::Exp(-TMath::Power(x / [1], [2])) / (1. + TMath::Exp(([3] - x) / [4]))", 0., 5.0);
 637   fEff->SetParameter(0, 0.5);
 638   fEff->SetParameter(1, 0.5);
 639   fEff->SetParameter(2, 1.);
 640   fEff->SetParameter(3, 0.5);
 641   fEff->SetParameter(4, 0.5);
 642   fEff->SetParLimits(4, 0.3, 2.0);
 643
 644   TFile *fileout = TFile::Open(Form("%s/efficiencyPt_MB_%s_%s.root", destdir, AliPID::ParticleName(ipart), chargeName[icharge]), "RECREATE");
 645
 646   TFile *filein = TFile::Open(filename);
 647   TH1D *hEfficiencyPt = (TH1D *)filein->Get(Form("hEfficiencyPt_MB_%s_%s_%s", histoName[ihisto], AliPID::ParticleName(ipart), chargeName[icharge]));
 648
 649   hEfficiencyPt->Fit(fEff, "0q", "IME", 0.5, 5.0);
 650   hEfficiencyPt->SetTitle(Form("%s;p_{T} (GeV/c);acceptance #times efficiency;", partChargeName[ipart][icharge]));
 651   hEfficiencyPt->SetMinimum(0.2);
 652   hEfficiencyPt->SetMaximum(0.7);
 653   hEfficiencyPt->SetMarkerStyle(marker);
 654   hEfficiencyPt->SetMarkerColor(color);
 655   hEfficiencyPt->SetMarkerSize(1.5);
 656   hEfficiencyPt->SetLineWidth(2);
 657   hEfficiencyPt->SetLineColor(color);
 658   hEfficiencyPt->GetXaxis()->SetRangeUser(ptMin[ipart] + 0.001, ptMax[ipart] - 0.001);
 659   hEfficiencyPt->SetStats(kFALSE);
 660   cCanvas1->cd();
 661   hEfficiencyPt->DrawCopy();
 662   fEff->DrawCopy("same");
 663   fileout->cd();
 664   hEfficiencyPt->Write("hEfficiencyPt");
 665   fEff->Write("fEff");
 666
 667   cCanvas1->SetGridx();
 668   cCanvas1->SetGridy();
 669   cCanvas1->SaveAs(Form("%s/efficiencyPt_MB_%s_%s.C", destdir, AliPID::ParticleName(ipart), chargeName[icharge]));
 670   cCanvas1->SaveAs(Form("%s/efficiencyPt_MB_%s_%s.png", destdir, AliPID::ParticleName(ipart), chargeName[icharge]));
 671   cCanvas1->SaveAs(Form("%s/efficiencyPt_MB_%s_%s.eps", destdir, AliPID::ParticleName(ipart), chargeName[icharge]));
 672
 673   TH1D *hRatioPt = new TH1D(*hEfficiencyPt);
 674   hRatioPt->Divide(fEff);
 675   hRatioPt->SetTitle(Form("%s;p_{T} (GeV/c);ratio wrt. fitted dependence;", partChargeName[ipart][icharge]));
 676   hRatioPt->SetMinimum(0.9);
 677   hRatioPt->SetMaximum(1.1);
 678   cCanvas2->cd();
 679   hRatioPt->DrawCopy();
 680   fileout->cd();
 681   hRatioPt->Write("hRatioPt");
 682
 683   cCanvas2->SetGridx();
 684   cCanvas2->SetGridy();
 685   cCanvas2->SaveAs(Form("%s/efficiencyPt_ratioFit_MB_%s_%s.C", destdir, AliPID::ParticleName(ipart), chargeName[icharge]));
 686   cCanvas2->SaveAs(Form("%s/efficiencyPt_ratioFit_MB_%s_%s.png", destdir, AliPID::ParticleName(ipart), chargeName[icharge]));
 687   cCanvas2->SaveAs(Form("%s/efficiencyPt_ratioFit_MB_%s_%s.eps", destdir, AliPID::ParticleName(ipart), chargeName[icharge]));
 688
 689
 690   fileout->Close();
 691
 692   //  hEfficiencyPt->Add(fEff, -1.);
 693   return hEfficiencyPt;
 694 }
 695
 696
 697 //_____________________________________________________________________________-
 698
 699 TOFmatchMC_efficiencyPt_centrality_plot(const Char_t *filename)
 700 {
 701   TOFmatchMC_efficiencyPt_centrality_plot(filename, kMatchedGoodTracks, 2, 0, 20, 4);
 702   TOFmatchMC_efficiencyPt_centrality_plot(filename, kMatchedGoodTracks, 2, 1, 25, 4);
 703   TOFmatchMC_efficiencyPt_centrality_plot(filename, kMatchedGoodTracks, 3, 0, 20, 8);
 704   TOFmatchMC_efficiencyPt_centrality_plot(filename, kMatchedGoodTracks, 3, 1, 25, 8);
 705   TOFmatchMC_efficiencyPt_centrality_plot(filename, kMatchedGoodTracks, 4, 0, 20, 2);
 706   TOFmatchMC_efficiencyPt_centrality_plot(filename, kMatchedGoodTracks, 4, 1, 25, 2);
 707 }
 708
 709 TOFmatchMC_efficiencyPt_fit(TH1 *h, TF1 *f, Int_t param)
 710 {
 711
 712   if (param >= 0) {
 713     for (Int_t ipar = 0; ipar < f->GetNpar(); ipar++)
 714       if (ipar != param)
 715         f->FixParameter(ipar, f->GetParameter(ipar));
 716       else
 717         f->ReleaseParameter(ipar);
 718   }
 719   else {
 720     for (Int_t ipar = 0; ipar < f->GetNpar(); ipar++)
 721       f->ReleaseParameter(ipar);
 722   }
 723   h->Fit(f);
 724
 725 }
 726
 727 TH1D *
 728 TOFmatchMC_efficiencyPt_centrality_plot(const Char_t *filename, Int_t ihisto = kMatchedGoodTracks, Int_t ipart, Int_t icharge, Int_t marker = 20, Int_t color = 2, Option_t *opt = "")
 729 {
 730
 731   TVirtualFitter::SetMaxIterations(1000000);
 732
 733   /* load HistoUtils */
 734   gROOT->LoadMacro("HistoUtils.C");
 735
 736   TCanvas *cCanvas1 = new TCanvas("cCanvas1");
 737   TCanvas *cCanvas2 = new TCanvas("cCanvas2");
 738   TCanvas *cCanvas3 = new TCanvas("cCanvas3");
 739
 740   Double_t ptMin[AliPID::kSPECIES] = {0.5, 0.5, 0.5, 0.5, 0.5};
 741   Double_t ptMax[AliPID::kSPECIES] = {3.0, 3.0, 3.0, 3.0, 5.0};
 742
 743   /* fit minimum-bias efficiency pt */
 744
 745   TF1 *fEff = new TF1("fEff", "[0] * TMath::Exp(-TMath::Power([1] / x, [2])) + [3] * x", ptMin[ipart], ptMax[ipart]);
 746   fEff->SetParameter(0, 0.5);
 747   fEff->SetParameter(1, 0.1);
 748   fEff->SetParameter(2, 2.);
 749   fEff->SetParameter(3, 0.);
 750   fEff->SetParameter(4, 0.);
 751   Int_t nPars = fEff->GetNpar();
 752   Int_t scalePar = 0.;
 753   TF1 *fEffCopy = new TF1(*fEff);
 754
 755   TFile *fileout = TFile::Open(Form("%s/efficiencyPt_centrality_%s_%s.root", destdir, AliPID::ParticleName(ipart), chargeName[icharge]), "RECREATE");
 756
 757   TFile *filein = TFile::Open(filename);
 758   TH1D *hEfficiencyPt = (TH1D *)filein->Get(Form("hEfficiencyPt_MB_%s_%s_%s", histoName[ihisto], AliPID::ParticleName(ipart), chargeName[icharge]));
 759   hEfficiencyPt->Fit(fEff, "0", "", 0.5, 1.0);
 760   hEfficiencyPt->Fit(fEff, "0", "", 0.5, 1.5);
 761   hEfficiencyPt->Fit(fEff, "0", "", 0.5, 2.0);
 762   hEfficiencyPt->Fit(fEff, "0", "IMRE", ptMin[ipart], ptMax[ipart]);
 763   hEfficiencyPt->Fit(fEff, "0", "IMRE", ptMin[ipart], ptMax[ipart]);
 764   hEfficiencyPt->Fit(fEff, "0", "IMRE", ptMin[ipart], ptMax[ipart]);
 765
 766   /* build efficiency profile */
 767   for (Int_t ipar = 0; ipar < nPars; ipar++) {
 768     fEffCopy->SetParameter(ipar, fEff->GetParameter(ipar));
 769     fEffCopy->SetParError(ipar, fEff->GetParError(ipar));
 770   }
 771   TProfile *pEfficiencyPt = new TProfile(Form("pEfficiencyPt_MB_%s_%s_%s", histoName[ihisto], AliPID::ParticleName(ipart), chargeName[icharge]), Form("%s;p_{T} (GeV/c);acceptance #times efficiency;", partChargeName[ipart][icharge]), NptBins, ptBin, "s");
 772   HistoUtils_Function2Profile(fEffCopy, pEfficiencyPt);
 773
 774   /* draw */
 775   cCanvas1->cd();
 776   pEfficiencyPt->SetMarkerStyle(marker);
 777   pEfficiencyPt->SetMarkerColor(color);
 778   pEfficiencyPt->SetMarkerSize(0);
 779   pEfficiencyPt->SetLineWidth(2);
 780   pEfficiencyPt->SetLineColor(color);
 781   pEfficiencyPt->SetFillStyle(3001);
 782   pEfficiencyPt->SetFillColor(color);
 783   pEfficiencyPt->GetXaxis()->SetRangeUser(ptMin[ipart] + 0.001, ptMax[ipart] - 0.001);
 784   pEfficiencyPt->SetStats(kFALSE);
 785   pEfficiencyPt->DrawCopy("E2");
 786
 787   hEfficiencyPt->SetTitle(Form("%s;p_{T} (GeV/c);acceptance #times efficiency;", partChargeName[ipart][icharge]));
 788   hEfficiencyPt->SetMarkerStyle(marker);
 789   hEfficiencyPt->SetMarkerColor(color);
 790   hEfficiencyPt->SetMarkerSize(1.5);
 791   hEfficiencyPt->SetLineWidth(2);
 792   hEfficiencyPt->SetLineColor(color);
 793   hEfficiencyPt->GetXaxis()->SetRangeUser(ptMin[ipart] + 0.001, ptMax[ipart] - 0.001);
 794   hEfficiencyPt->SetStats(kFALSE);
 795   hEfficiencyPt->DrawCopy("same");
 796   fEff->DrawCopy("same");
 797   cCanvas1->Update();
 798
 799   /* write */
 800   fileout->cd();
 801   pEfficiencyPt->Write("pEfficiencyPt");
 802   hEfficiencyPt->Write("hEfficiencyPt");
 803   fEff->Write("fEff");
 804
 805   TH1D *hEfficiencyCent = new TH1D("hEfficiencyCent", Form("%s;centrality percentile;acceptance x efficiency scale factor;", partChargeName[ipart][icharge]), NcentralityBins, centralityBin);
 806   hEfficiencyCent->SetMinimum(0.9);
 807   hEfficiencyCent->SetMaximum(1.1);
 808   hEfficiencyCent->SetMarkerStyle(marker);
 809   hEfficiencyCent->SetMarkerColor(color);
 810   hEfficiencyCent->SetMarkerSize(1.5);
 811   hEfficiencyCent->SetLineWidth(2);
 812   hEfficiencyCent->SetLineColor(color);
 813   hEfficiencyCent->SetStats(kFALSE);
 814
 815   TProfile *pEfficiencyPt_cent[NcentralityBins];
 816   TH1D *hEfficiencyPt_cent[NcentralityBins];
 817   TH1D *hRatioPt_cent[NcentralityBins];
 818
 819   /* fix efficiency shape and release scale factor */
 820   for (Int_t ipar = 0; ipar < nPars; ipar++)
 821     fEff->FixParameter(ipar, fEff->GetParameter(ipar));
 822   fEff->ReleaseParameter(scalePar);
 823
 824   gStyle->SetOptStat(1100);
 825   for (Int_t icent = 0; icent < NcentralityBins; icent++) {
 826
 827     hEfficiencyPt_cent[icent] = (TH1D *)filein->Get(Form("hEfficiencyPt_centrality%d_%s_%s_%s", icent, histoName[ihisto], AliPID::ParticleName(ipart), chargeName[icharge]));
 828     hEfficiencyPt_cent[icent]->Fit(fEff, "", "IME", ptMin[ipart], ptMax[ipart]);
 829
 830     /* build efficiency profile */
 831     fEffCopy->SetParameter(scalePar, fEff->GetParameter(scalePar));
 832     fEffCopy->SetParError(scalePar, fEff->GetParError(scalePar));
 833     pEfficiencyPt_cent[icent] = new TProfile(Form("pEfficiencyPt_centrality%d_%s_%s_%s", icent, histoName[ihisto], AliPID::ParticleName(ipart), chargeName[icharge]), Form("%s;p_{T} (GeV/c);acceptance #times efficiency;", partChargeName[ipart][icharge]), NptBins, ptBin, "s");
 834     HistoUtils_Function2Profile(fEffCopy, pEfficiencyPt_cent[icent]);
 835
 836     /* draw */
 837     cCanvas1->cd();
 838     pEfficiencyPt_cent[icent]->SetMarkerStyle(marker);
 839     pEfficiencyPt_cent[icent]->SetMarkerColor(color);
 840     pEfficiencyPt_cent[icent]->SetMarkerSize(0);
 841     pEfficiencyPt_cent[icent]->SetLineWidth(2);
 842     pEfficiencyPt_cent[icent]->SetLineColor(color);
 843     pEfficiencyPt_cent[icent]->SetFillStyle(3001);
 844     pEfficiencyPt_cent[icent]->SetFillColor(color);
 845     pEfficiencyPt_cent[icent]->GetXaxis()->SetRangeUser(ptMin[ipart] + 0.001, ptMax[ipart] - 0.001);
 846     pEfficiencyPt_cent[icent]->SetStats(kFALSE);
 847     pEfficiencyPt_cent[icent]->DrawCopy("E2");
 848
 849     hEfficiencyPt_cent[icent]->SetTitle(Form("%s (%d-%d\%);p_{T} (GeV/c);acceptance #times efficiency;", partChargeName[ipart][icharge], (Int_t)centralityBin[icent], (Int_t)centralityBin[icent + 1]));
 850     hEfficiencyPt_cent[icent]->SetMarkerStyle(marker);
 851     hEfficiencyPt_cent[icent]->SetMarkerColor(color);
 852     hEfficiencyPt_cent[icent]->SetMarkerSize(1.5);
 853     hEfficiencyPt_cent[icent]->SetLineWidth(2);
 854     hEfficiencyPt_cent[icent]->SetLineColor(color);
 855     hEfficiencyPt_cent[icent]->GetXaxis()->SetRangeUser(ptMin[ipart] + 0.001, ptMax[ipart] - 0.001);
 856     hEfficiencyPt_cent[icent]->SetStats(kFALSE);
 857     hEfficiencyPt_cent[icent]->DrawCopy("same");
 858     fEff->DrawCopy("same");
 859     cCanvas1->Update();
 860
 861     fileout->cd();
 862     pEfficiencyPt_cent[icent]->Write(Form("pEfficiencyPt_cent%d", icent));
 863     hEfficiencyPt_cent[icent]->Write(Form("hEfficiencyPt_cent%d", icent));
 864     fEff->Write(Form("fEff_cent%d", icent));
 865
 866     hEfficiencyCent->SetBinContent(icent + 1, fEff->GetParameter(scalePar));
 867     hEfficiencyCent->SetBinError(icent + 1, fEff->GetParError(scalePar));
 868
 869     cCanvas1->SetGridx();
 870     cCanvas1->SetGridy();
 871     cCanvas1->SaveAs(Form("%s/efficiencyPt_centrality%d_%s_%s.C", destdir, icent, AliPID::ParticleName(ipart), chargeName[icharge]));
 872     cCanvas1->SaveAs(Form("%s/efficiencyPt_centrality%d_%s_%s.png", destdir, icent, AliPID::ParticleName(ipart), chargeName[icharge]));
 873     cCanvas1->SaveAs(Form("%s/efficiencyPt_centrality%d_%s_%s.eps", destdir, icent, AliPID::ParticleName(ipart), chargeName[icharge]));
 874
 875     hRatioPt_cent[icent] = new TH1D(*hEfficiencyPt_cent[icent]);
 876     hRatioPt_cent[icent]->Divide(fEff);
 877     hRatioPt_cent[icent]->SetTitle(Form("%s (%d-%d\%);p_{T} (GeV/c);ratio wrt. fitted dependence;", partChargeName[ipart][icharge], (Int_t)centralityBin[icent], (Int_t)centralityBin[icent + 1]));
 878     hRatioPt_cent[icent]->SetMinimum(0.75);
 879     hRatioPt_cent[icent]->SetMaximum(1.25);
 880     cCanvas2->cd();
 881     hRatioPt_cent[icent]->DrawCopy();
 882     fileout->cd();
 883     hRatioPt_cent[icent]->Write(Form("hRatioPt_cent%d", icent));
 884
 885
 886     cCanvas2->SetGridx();
 887     cCanvas2->SetGridy();
 888     cCanvas2->SaveAs(Form("%s/efficiencyPt_ratioFit_centrality%d_%s_%s.C", destdir, icent, AliPID::ParticleName(ipart), chargeName[icharge]));
 889     cCanvas2->SaveAs(Form("%s/efficiencyPt_ratioFit_centrality%d_%s_%s.png", destdir, icent, AliPID::ParticleName(ipart), chargeName[icharge]));
 890     cCanvas2->SaveAs(Form("%s/efficiencyPt_ratioFit_centrality%d_%s_%s.eps", destdir, icent, AliPID::ParticleName(ipart), chargeName[icharge]));
 891
 892   }
 893
 894   cCanvas3->cd();
 895   hEfficiencyCent->DrawCopy();
 896   fileout->cd();
 897   hEfficiencyCent->Write("hEfficiencyCent");
 898
 899   cCanvas3->SetGridx();
 900   cCanvas3->SetGridy();
 901   cCanvas3->SaveAs(Form("%s/efficiencyCentrality_scaleFactor_%s_%s.C", destdir, AliPID::ParticleName(ipart), chargeName[icharge]));
 902   cCanvas3->SaveAs(Form("%s/efficiencyCentrality_scaleFactor_%s_%s.png", destdir, AliPID::ParticleName(ipart), chargeName[icharge]));
 903   cCanvas3->SaveAs(Form("%s/efficiencyCentrality_scaleFactor_%s_%s.eps", destdir, AliPID::ParticleName(ipart), chargeName[icharge]));
 904
 905   //  hEfficiencyPt->Add(fEff, -1.);
 906   return hEfficiencyCent;
 907 }
 908
 909
 910 //_____________________________________________________________________________-
 911
 912 TH1D *
 913 TOFmatchMC_efficiencyPt_centrality_all_plot(const Char_t *filename, Int_t ihisto = kMatchedGoodTracks, Int_t marker = 20, Int_t color = 1, Option_t *opt = "")
 914 {
 915
 916   /*
 917     this function measures the centrality-dependent scale factor
 918     of the matching efficiency, assuming that the scale factor
 919     is common to all particles species
 920
 921     1.) first the inclusive pt-dependence is fitted for MB events with scale factor = 1.
 922     2.) the shape of the pt-dependence is the fixed and the scale factor released
 923     3.) using a fixed-shape the inclusive pt-dependence is fitted in centrality bins with only the scale factor as a free parameter
 924
 925   */
 926
 927   TCanvas *cCanvas1 = new TCanvas("cCanvas1");
 928   TCanvas *cCanvas2 = new TCanvas("cCanvas2");
 929
 930   /* pt-dependent matching efficiency function with scale factor */
 931   TF1 *fEff = new TF1("fEff", "[5] * [0] * TMath::Exp(-TMath::Power(x / [1], [2])) / (1. + TMath::Exp(([3] - x) / [4]))", 0., 5.0);
 932   fEff->SetParameter(0, 0.5);
 933   fEff->SetParameter(1, 0.5);
 934   fEff->SetParameter(2, 1.);
 935   fEff->SetParameter(3, 0.5);
 936   fEff->SetParameter(4, 0.5);
 937   fEff->SetParLimits(4, 0.3, 2.0);
 938   fEff->FixParameter(5, 1.);
 939
 940   TFile *fileout = TFile::Open(Form("%s/efficiencyPt_MB_all.root", destdir), "RECREATE");
 941
 942   /* fit inclusive pt-dependence in MB events */
 943   TFile *filein = TFile::Open(filename);
 944   TH1D *hEfficiencyPt = (TH1D *)filein->Get(Form("hEfficiencyAllPt_MB_%s", histoName[ihisto]));
 945   hEfficiencyPt->Fit(fEff, "0", "IME", 0.5, 5.0);
 946   hEfficiencyPt->SetTitle("all particles;p_{T} (GeV/c);acceptance #times efficiency;");
 947   hEfficiencyPt->SetMinimum(0.2);
 948   hEfficiencyPt->SetMaximum(0.7);
 949   hEfficiencyPt->SetMarkerStyle(marker);
 950   hEfficiencyPt->SetMarkerColor(color);
 951   hEfficiencyPt->SetMarkerSize(1.5);
 952   hEfficiencyPt->SetLineWidth(2);
 953   hEfficiencyPt->SetLineColor(color);
 954   hEfficiencyPt->GetXaxis()->SetRangeUser(0.5 + 0.001, 5.0 - 0.001);
 955   hEfficiencyPt->SetStats(kFALSE);
 956   cCanvas1->cd();
 957   hEfficiencyPt->DrawCopy();
 958   fEff->DrawCopy("same");
 959   fileout->cd();
 960   hEfficiencyPt->Write("hEfficiencyPt");
 961   fEff->Write("fEff");
 962
 963   cCanvas1->SetGridx();
 964   cCanvas1->SetGridy();
 965   cCanvas1->SaveAs(Form("%s/efficiencyPt_MB_all.C", destdir));
 966   cCanvas1->SaveAs(Form("%s/efficiencyPt_MB_all.png", destdir));
 967   cCanvas1->SaveAs(Form("%s/efficiencyPt_MB_all.eps", destdir));
 968
 969   TH1D *hRatioPt = new TH1D(*hEfficiencyPt);
 970   hRatioPt->Divide(fEff);
 971   hRatioPt->SetTitle("all particles;p_{T} (GeV/c);ratio wrt. fitted dependence;");
 972   hRatioPt->SetMinimum(0.9);
 973   hRatioPt->SetMaximum(1.1);
 974   cCanvas2->cd();
 975   hRatioPt->DrawCopy();
 976   fileout->cd();
 977   hRatioPt->Write("hRatioPt");
 978
 979   cCanvas2->SetGridx();
 980   cCanvas2->SetGridy();
 981   cCanvas2->SaveAs(Form("%s/efficiencyPt_ratioFit_MB_all.C", destdir));
 982   cCanvas2->SaveAs(Form("%s/efficiencyPt_ratioFit_MB_all.png", destdir));
 983   cCanvas2->SaveAs(Form("%s/efficiencyPt_ratioFit_MB_all.eps", destdir));
 984
 985   /* fix efficiency shape and release scale factor */
 986   fEff->FixParameter(0, fEff->GetParameter(0));
 987   fEff->FixParameter(1, fEff->GetParameter(1));
 988   fEff->FixParameter(2, fEff->GetParameter(2));
 989   fEff->FixParameter(3, fEff->GetParameter(3));
 990   fEff->FixParameter(4, fEff->GetParameter(4));
 991   fEff->ReleaseParameter(5);
 992
 993   TH1D *hEfficiencyCent = new TH1D("hEfficiencyCent", "all particles;centrality percentile;acceptance x efficiency scale factor;", NcentralityBins, centralityBin);
 994   hEfficiencyCent->SetMinimum(0.95);
 995   hEfficiencyCent->SetMaximum(1.05);
 996   hEfficiencyCent->SetMarkerStyle(marker);
 997   hEfficiencyCent->SetMarkerColor(color);
 998   hEfficiencyCent->SetMarkerSize(1.5);
 999   hEfficiencyCent->SetLineWidth(2);
1000   hEfficiencyCent->SetLineColor(color);
1001   hEfficiencyCent->SetStats(kFALSE);
1002
1003
1004   /* fit inclusive pt-dependence in centrality bins
1005      with fixed pt-shape and free scale factor */
1006   TH1D *hEfficiencyPt_cent[NcentralityBins];
1007   TH1D *hRatioPt_cent[NcentralityBins];
1008   for (Int_t icent = 0; icent < NcentralityBins; icent++) {
1009
1010     hEfficiencyPt_cent[icent] = (TH1D *)filein->Get(Form("hEfficiencyAllPt_centrality%d_%s", icent, histoName[ihisto]));
1011     hEfficiencyPt_cent[icent]->Fit(fEff, "", "IME", 0.5, 5.0);
1012
1013     hEfficiencyPt_cent[icent]->SetTitle(Form("all particles (%d-%d\%);p_{T} (GeV/c);acceptance #times efficiency;", (Int_t)centralityBin[icent], (Int_t)centralityBin[icent + 1]));
1014     hEfficiencyPt_cent[icent]->SetMinimum(0.2);
1015     hEfficiencyPt_cent[icent]->SetMaximum(0.8);
1016     hEfficiencyPt_cent[icent]->SetMarkerStyle(marker);
1017     hEfficiencyPt_cent[icent]->SetMarkerColor(color);
1018     hEfficiencyPt_cent[icent]->SetMarkerSize(1.5);
1019     hEfficiencyPt_cent[icent]->SetLineWidth(2);
1020     hEfficiencyPt_cent[icent]->SetLineColor(color);
1021     hEfficiencyPt_cent[icent]->GetXaxis()->SetRangeUser(0.5 + 0.001, 5.0 - 0.001);
1022     hEfficiencyPt_cent[icent]->SetStats(kFALSE);
1023     cCanvas1->cd();
1024     hEfficiencyPt_cent[icent]->DrawCopy();
1025     fEff->DrawCopy("same");
1026     hEfficiencyPt_cent[icent]->Write(Form("hEfficiencyPt_cent%d", icent));
1027     fEff->Write(Form("fEff_cent%d", icent));
1028
1029     hEfficiencyCent->SetBinContent(icent + 1, fEff->GetParameter(5));
1030     hEfficiencyCent->SetBinError(icent + 1, fEff->GetParError(5));
1031
1032     cCanvas1->SetGridx();
1033     cCanvas1->SetGridy();
1034     cCanvas1->SaveAs(Form("%s/efficiencyPt_centrality%d_all.C", destdir, icent));
1035     cCanvas1->SaveAs(Form("%s/efficiencyPt_centrality%d_all.png", destdir, icent));
1036     cCanvas1->SaveAs(Form("%s/efficiencyPt_centrality%d_all.eps", destdir, icent));
1037
1038     hRatioPt_cent[icent] = new TH1D(*hEfficiencyPt_cent[icent]);
1039     hRatioPt_cent[icent]->Divide(fEff);
1040     hRatioPt_cent[icent]->SetTitle(Form("all particles (%d-%d\%);p_{T} (GeV/c);ratio wrt. fitted dependence;", (Int_t)centralityBin[icent], (Int_t)centralityBin[icent + 1]));
1041     hRatioPt_cent[icent]->SetMinimum(0.9);
1042     hRatioPt_cent[icent]->SetMaximum(1.1);
1043     cCanvas2->cd();
1044     hRatioPt_cent[icent]->DrawCopy();
1045     fileout->cd();
1046     hRatioPt_cent[icent]->Write(Form("hRatioPt_cent%d", icent));
1047
1048     cCanvas2->SetGridx();
1049     cCanvas2->SetGridy();
1050     cCanvas2->SaveAs(Form("%s/efficiencyPt_ratioFit_centrality%d_all.C", destdir, icent));
1051     cCanvas2->SaveAs(Form("%s/efficiencyPt_ratioFit_centrality%d_all.png", destdir, icent));
1052     cCanvas2->SaveAs(Form("%s/efficiencyPt_ratioFit_centrality%d_all.eps", destdir, icent));
1053
1054
1055   }
1056
1057   TF1 *fEffCent = new TF1("fEffCent", "[0] - [1] * TMath::Exp(-[2] * TMath::Power(x, [3]))", 0., 90.);
1058   fEffCent->SetParameter(0, 1.02);
1059   fEffCent->SetParameter(1, 0.04);
1060   fEffCent->SetParameter(2, 0.001);
1061   fEffCent->SetParameter(3, 2.);
1062   hEfficiencyCent->Fit(fEffCent, "q0", "IME", 0., 90.);
1063
1064   TCanvas *cCanvas3 = new TCanvas("cCanvas3");
1065   hEfficiencyCent->DrawCopy();
1066   fEffCent->DrawCopy("same");
1067   fileout->cd();
1068   hEfficiencyCent->Write("hEfficiencyCent");
1069   fEffCent->Write("fEffCent");
1070
1071   cCanvas3->SetGridx();
1072   cCanvas3->SetGridy();
1073   cCanvas3->SaveAs(Form("%s/efficiencyCentrality_scaleFactor_all.C", destdir));
1074   cCanvas3->SaveAs(Form("%s/efficiencyCentrality_scaleFactor_all.png", destdir));
1075   cCanvas3->SaveAs(Form("%s/efficiencyCentrality_scaleFactor_all.eps", destdir));
1076
1077   TCanvas *cCanvas4 = new TCanvas("cCanvas4");
1078
1079   TH1D *hRatioCent = new TH1D(*hEfficiencyCent);
1080   hRatioCent->Divide(fEffCent);
1081   hRatioCent->SetTitle(Form("all particles;centrality percentile;ratio wrt. fitted dependence;"));
1082   hRatioCent->SetMinimum(0.95);
1083   hRatioCent->SetMaximum(1.05);
1084   cCanvas4->cd();
1085   hRatioCent->DrawCopy();
1086   fileout->cd();
1087   hRatioCent->Write("hRatioCent");
1088
1089   cCanvas4->SetGridx();
1090   cCanvas4->SetGridy();
1091   cCanvas4->SaveAs(Form("%s/efficiencyCentrality_scaleFactor_ratioFit_all.C", destdir));
1092   cCanvas4->SaveAs(Form("%s/efficiencyCentrality_scaleFactor_ratioFit_all.png", destdir));
1093   cCanvas4->SaveAs(Form("%s/efficiencyCentrality_scaleFactor_ratioFit_all.eps", destdir));
1094
1095
1096
1097   //  hEfficiencyPt->Add(fEff, -1.);
1098   return hEfficiencyCent;
1099 }
1100
1101
1102 //_____________________________________________________________________________-
1103
1104 TOFmatchMC_efficiencyPt_MB_smooth_plot(const Char_t *filename)
1105 {
1106   TOFmatchMC_efficiencyPt_MB_smooth_plot(filename, kMatchedGoodTracks, 2, 0, 20, 4);
1107   TOFmatchMC_efficiencyPt_MB_smooth_plot(filename, kMatchedGoodTracks, 2, 1, 25, 4);
1108   TOFmatchMC_efficiencyPt_MB_smooth_plot(filename, kMatchedGoodTracks, 3, 0, 20, 8);
1109   TOFmatchMC_efficiencyPt_MB_smooth_plot(filename, kMatchedGoodTracks, 3, 1, 25, 8);
1110   TOFmatchMC_efficiencyPt_MB_smooth_plot(filename, kMatchedGoodTracks, 4, 0, 20, 2);
1111   TOFmatchMC_efficiencyPt_MB_smooth_plot(filename, kMatchedGoodTracks, 4, 1, 25, 2);
1112 }
1113
1114 TH1D *
1115 TOFmatchMC_efficiencyPt_MB_smooth_plot(const Char_t *filename, Int_t ihisto = kMatchedGoodTracks, Int_t ipart, Int_t icharge, Int_t marker = 20, Int_t color = 2, Option_t *opt = "")
1116 {
1117
1118   TCanvas *cCanvas1 = new TCanvas("cCanvas1");
1119   TCanvas *cCanvas2 = new TCanvas("cCanvas2");
1120
1121   Double_t ptMin[AliPID::kSPECIES] = {0.5, 0.5, 0.5, 0.5, 0.5};
1122   Double_t ptMax[AliPID::kSPECIES] = {3.0, 3.0, 3.0, 3.0, 5.0};
1123
1124   TF1 *fEff = new TF1("fEff", "[0] + [1] * x - [2] * TMath::Exp(-[3] * TMath::Power(x, [4]))", 1.0, 5.0);
1125   fEff->SetParameter(0, 0.5);
1126   fEff->SetParameter(1, 0.);
1127   fEff->SetParameter(2, 0.5);
1128   fEff->SetParameter(3, 1.);
1129   fEff->SetParameter(4, 2.);
1130
1131
1132   TFile *filein = TFile::Open(filename);
1133   TH1D *hEfficiencyPt = (TH1D *)filein->Get(Form("hEfficiencyPt_MB_%s_%s_%s", histoName[ihisto], AliPID::ParticleName(ipart), chargeName[icharge]));
1134   hEfficiencyPt->Fit(fEff, "0", "IME", 1.0, 5.0);
1135   hEfficiencyPt->SetTitle(Form("%s;p_{T} (GeV/c);acceptance #times efficiency;", partChargeName[ipart][icharge]));
1136   hEfficiencyPt->SetMinimum(0.2);
1137   hEfficiencyPt->SetMaximum(0.7);
1138   hEfficiencyPt->SetMarkerStyle(marker);
1139   hEfficiencyPt->SetMarkerColor(color);
1140   hEfficiencyPt->SetMarkerSize(1.5);
1141   hEfficiencyPt->SetLineWidth(2);
1142   hEfficiencyPt->SetLineColor(color);
1143   hEfficiencyPt->GetXaxis()->SetRangeUser(ptMin[ipart] + 0.001, ptMax[ipart] - 0.001);
1144   hEfficiencyPt->SetStats(kFALSE);
1145   cCanvas1->cd();
1146   hEfficiencyPt->Draw(opt);
1147   fEff->Draw("same");
1148
1149   cCanvas1->SetGridx();
1150   cCanvas1->SetGridy();
1151   cCanvas1->SaveAs(Form("%s/efficiencyPt_MB_smooth_%s_%s.C", destdir, AliPID::ParticleName(ipart), chargeName[icharge]));
1152   cCanvas1->SaveAs(Form("%s/efficiencyPt_MB_smooth_%s_%s.root", destdir, AliPID::ParticleName(ipart), chargeName[icharge]));
1153   cCanvas1->SaveAs(Form("%s/efficiencyPt_MB_smooth_%s_%s.png", destdir, AliPID::ParticleName(ipart), chargeName[icharge]));
1154   cCanvas1->SaveAs(Form("%s/efficiencyPt_MB_smooth_%s_%s.eps", destdir, AliPID::ParticleName(ipart), chargeName[icharge]));
1155
1156   TH1D *hRatioPt = new TH1D(*hEfficiencyPt);
1157   hRatioPt->Divide(fEff);
1158   hRatioPt->SetTitle(Form("%s;p_{T} (GeV/c);ratio wrt. fitted dependence;", partChargeName[ipart][icharge]));
1159   hRatioPt->SetMinimum(0.9);
1160   hRatioPt->SetMaximum(1.1);
1161   cCanvas2->cd();
1162   hRatioPt->Draw();
1163
1164   cCanvas2->SetGridx();
1165   cCanvas2->SetGridy();
1166   cCanvas2->SaveAs(Form("%s/efficiencyPt_ratioFit_MB_smooth_%s_%s.C", destdir, AliPID::ParticleName(ipart), chargeName[icharge]));
1167   cCanvas2->SaveAs(Form("%s/efficiencyPt_ratioFit_MB_smooth_%s_%s.root", destdir, AliPID::ParticleName(ipart), chargeName[icharge]));
1168   cCanvas2->SaveAs(Form("%s/efficiencyPt_ratioFit_MB_smooth_%s_%s.png", destdir, AliPID::ParticleName(ipart), chargeName[icharge]));
1169   cCanvas2->SaveAs(Form("%s/efficiencyPt_ratioFit_MB_smooth_%s_%s.eps", destdir, AliPID::ParticleName(ipart), chargeName[icharge]));
1170
1171
1172   /* smoothing step */
1173   TH1D *hEfficiencyPt_smooth = new TH1D(*hEfficiencyPt);
1174   Double_t ratio, value, error;
1175   for (Int_t ipt = 0; ipt < hRatioPt->GetNbinsX(); ipt++) {
1176     if (hRatioPt->GetBinCenter(ipt + 1) < 1.5) continue;
1177     ratio = hRatioPt->GetBinContent(ipt + 1);
1178     value = hEfficiencyPt->GetBinContent(ipt + 1);
1179     error = hEfficiencyPt->GetBinError(ipt + 1);
1180     value /= ratio;
1181     hEfficiencyPt_smooth->SetBinContent(ipt + 1, value);
1182   }
1183   new TCanvas();
1184   hEfficiencyPt_smooth->Draw();
1185
1186   //  hEfficiencyPt->Add(fEff, -1.);
1187   return hEfficiencyPt;
1188 }
1189
1190 TOFmatchMC_matchingEfficiency(const Char_t *filename, Bool_t useGFcorrection = kTRUE)
1191 {
1192
1193   Int_t marker[2] = {20, 21};
1194   Int_t color[AliPID::kSPECIES] = {1, 1, 4, 8, 2};
1195   Char_t *partLatex[AliPID::kSPECIES][2] = {
1196     "", "", "#pi^{+}", "K^{+}", "p",
1197     "", "", "#pi^{-}", "K^{-}", "#bar{p}"
1198   };
1199
1200   TFile *filein = TFile::Open(filename);
1201   if (useGFcorrection)
1202     TFile *fileout = TFile::Open("TOF_matchingEfficiency.root", "RECREATE");
1203   else
1204     TFile *fileout = TFile::Open("TOF_matchingEfficiency_noGF.root", "RECREATE");
1205   TH1D *hEff;
1206   TF1 *fGF;
1207   Char_t title[1024];
1208   for (Int_t icent = -1; icent < NcentralityBins; icent++)
1209     for (Int_t icharge = 0; icharge < kNCharges; icharge++)
1210       for (Int_t ipart = 2; ipart < AliPID::kSPECIES; ipart++) {
1211         if (icent == -1)
1212           hEff = (TH1D *)filein->Get(Form("hEfficiencyPt_MB_hMatchedCorrelatedTracks_%s_%s", AliPID::ParticleName(ipart), chargeName[icharge]));
1213         else
1214           hEff = (TH1D *)filein->Get(Form("hEfficiencyPt_centrality%d_hMatchedCorrelatedTracks_%s_%s", icent, AliPID::ParticleName(ipart), chargeName[icharge]));
1215         /* geant-fluka correction */
1216         fGF = TOFmatchMC_geantflukaCorrection(ipart, icharge);
1217         if (useGFcorrection)
1218           hEff->Divide(fGF);
1219         if (icent == -1)
1220           sprintf(title, "%s (MB);p_{T} (GeV/c);matching efficiency;", partLatex[ipart][icharge]);
1221         else
1222           sprintf(title, "%s (%d-%d%%);p_{T} (GeV/c);matching efficiency;", partLatex[ipart][icharge], (Int_t)centralityBin[icent], (Int_t)centralityBin[icent + 1]);
1223         hEff->SetMarkerStyle(marker[icharge]);
1224         hEff->SetMarkerColor(color[ipart]);
1225         hEff->SetLineColor(1);
1226         hEff->SetLineWidth(1);
1227         hEff->SetTitle(title);
1228         if (icent == -1)
1229           hEff->SetName(Form("hMatchEff_MB_%s_%s", AliPID::ParticleName(ipart), chargeName[icharge]));
1230         else
1231           hEff->SetName(Form("hMatchEff_cent%d_%s_%s", icent, AliPID::ParticleName(ipart), chargeName[icharge]));
1232         fileout->cd();
1233         hEff->Write();
1234       }
1235
1236   fileout->Close();
1237 }
1238
1239 TOFmatchMC_matchingEfficiency_mismatchCorrected(const Char_t *filename, Bool_t useGFcorrection = kTRUE)
1240 {
1241
1242   Int_t marker[2] = {20, 21};
1243   Int_t color[AliPID::kSPECIES] = {1, 1, 4, 8, 2};
1244   Char_t *partLatex[AliPID::kSPECIES][2] = {
1245     "", "", "#pi^{+}", "K^{+}", "p",
1246     "", "", "#pi^{-}", "K^{-}", "#bar{p}"
1247   };
1248
1249   TFile *filein = TFile::Open(filename);
1250   if (useGFcorrection)
1251     TFile *fileout = TFile::Open("TOF_matchingEfficiency_mismatchCorrected.root", "RECREATE");
1252   else
1253     TFile *fileout = TFile::Open("TOF_matchingEfficiency_mismatchCorrected_noGF.root", "RECREATE");
1254   TH1D *hEff;
1255   TF1 *fGF;
1256   Char_t title[1024];
1257   for (Int_t icent = -1; icent < NcentralityBins; icent++)
1258     for (Int_t icharge = 0; icharge < kNCharges; icharge++)
1259       for (Int_t ipart = 2; ipart < AliPID::kSPECIES; ipart++) {
1260         if (icent == -1)
1261           hEff = (TH1D *)filein->Get(Form("hEfficiencyPt_MB_hMatchedTracks_%s_%s", AliPID::ParticleName(ipart), chargeName[icharge]));
1262         else
1263           hEff = (TH1D *)filein->Get(Form("hEfficiencyPt_centrality%d_hMatchedTracks_%s_%s", icent, AliPID::ParticleName(ipart), chargeName[icharge]));
1264         /* geant-fluka correction */
1265         fGF = TOFmatchMC_geantflukaCorrection(ipart, icharge);
1266         if (useGFcorrection)
1267           hEff->Divide(fGF);
1268         if (icent == -1)
1269 sprintf(title, "%s (MB);p_{T} (GeV/c);matching efficiency;", partLatex[ipart][icharge]);
1270         else
1271           sprintf(title, "%s (%d-%d%%);p_{T} (GeV/c);matching efficiency;", partLatex[ipart][icharge], (Int_t)centralityBin[icent], (Int_t)centralityBin[icent + 1]);
1272         hEff->SetMarkerStyle(marker[icharge]);
1273         hEff->SetMarkerColor(color[ipart]);
1274         hEff->SetLineColor(1);
1275         hEff->SetLineWidth(1);
1276         hEff->SetTitle(title);
1277         if (icent == -1)
1278           hEff->SetName(Form("hMatchEff_MB_%s_%s", AliPID::ParticleName(ipart), chargeName[icharge]));
1279         else
1280           hEff->SetName(Form("hMatchEff_cent%d_%s_%s", icent, AliPID::ParticleName(ipart), chargeName[icharge]));
1281         fileout->cd();
1282         hEff->Write();
1283       }
1284
1285   fileout->Close();
1286 }
1287
1288 TOFmatchMC_matchingEfficiency_acceptTRDin(const Char_t *filename, Bool_t useGFcorrection = kTRUE)
1289 {
1290
1291   Int_t marker[2] = {20, 21};
1292   Int_t color[AliPID::kSPECIES] = {1, 1, 4, 8, 2};
1293   Char_t *partLatex[AliPID::kSPECIES][2] = {
1294     "", "", "#pi^{+}", "K^{+}", "p",
1295     "", "", "#pi^{-}", "K^{-}", "#bar{p}"
1296   };
1297
1298   TFile *filein = TFile::Open(filename);
1299   if (useGFcorrection)
1300     TFile *fileout = TFile::Open("TOF_matchingEfficiency_acceptTRDin.root", "RECREATE");
1301   else
1302     TFile *fileout = TFile::Open("TOF_matchingEfficiency_acceptTRDin_noGF.root", "RECREATE");
1303   TH1D *hEff;
1304   TF1 *fGF;
1305   Char_t title[1024];
1306   for (Int_t icent = 0; icent < NcentralityBins; icent++)
1307     for (Int_t icharge = 0; icharge < kNCharges; icharge++)
1308       for (Int_t ipart = 2; ipart < AliPID::kSPECIES; ipart++) {
1309         hEff = (TH1D *)filein->Get(Form("hEfficiencyPt_centrality%d_hMatchedGoodTracks_%s_%s", icent, AliPID::ParticleName(ipart), chargeName[icharge]));
1310         /* geant-fluka correction */
1311         fGF = TOFmatchMC_geantflukaCorrection_acceptTRDin(ipart, icharge);
1312         if (useGFcorrection)
1313           hEff->Divide(fGF);
1314         sprintf(title, "%s (%d-%d%%);p_{T} (GeV/c);matching efficiency;", partLatex[ipart][icharge], (Int_t)centralityBin[icent], (Int_t)centralityBin[icent + 1]);
1315         hEff->SetMarkerStyle(marker[icharge]);
1316         hEff->SetMarkerColor(color[ipart]);
1317         hEff->SetLineColor(1);
1318         hEff->SetLineWidth(1);
1319         hEff->SetTitle(title);
1320         hEff->SetName(Form("hMatchEff_cent%d_%s_%s", icent, AliPID::ParticleName(ipart), chargeName[icharge]));
1321         fileout->cd();
1322         hEff->Write();
1323       }
1324
1325   fileout->Close();
1326 }
1327
1328 TOFmatchMC_matchingEfficiency_rejectTRDin(const Char_t *filename, Bool_t useGFcorrection = kTRUE)
1329 {
1330
1331   Int_t marker[2] = {20, 21};
1332   Int_t color[AliPID::kSPECIES] = {1, 1, 4, 8, 2};
1333   Char_t *partLatex[AliPID::kSPECIES][2] = {
1334     "", "", "#pi^{+}", "K^{+}", "p",
1335     "", "", "#pi^{-}", "K^{-}", "#bar{p}"
1336   };
1337
1338   TFile *filein = TFile::Open(filename);
1339   if (useGFcorrection)
1340     TFile *fileout = TFile::Open("TOF_matchingEfficiency_rejectTRDin.root", "RECREATE");
1341   else
1342     TFile *fileout = TFile::Open("TOF_matchingEfficiency_rejectTRDin_noGF.root", "RECREATE");
1343   TH1D *hEff;
1344   TF1 *fGF;
1345   Char_t title[1024];
1346   for (Int_t icent = 0; icent < NcentralityBins; icent++)
1347     for (Int_t icharge = 0; icharge < kNCharges; icharge++)
1348       for (Int_t ipart = 2; ipart < AliPID::kSPECIES; ipart++) {
1349         hEff = (TH1D *)filein->Get(Form("hEfficiencyPt_centrality%d_hMatchedGoodTracks_%s_%s", icent, AliPID::ParticleName(ipart), chargeName[icharge]));
1350         /* geant-fluka correction */
1351         fGF = TOFmatchMC_geantflukaCorrection_rejectTRDin(ipart, icharge);
1352         if (useGFcorrection)
1353           hEff->Divide(fGF);
1354         sprintf(title, "%s (%d-%d%%);p_{T} (GeV/c);matching efficiency;", partLatex[ipart][icharge], (Int_t)centralityBin[icent], (Int_t)centralityBin[icent + 1]);
1355         hEff->SetMarkerStyle(marker[icharge]);
1356         hEff->SetMarkerColor(color[ipart]);
1357         hEff->SetLineColor(1);
1358         hEff->SetLineWidth(1);
1359         hEff->SetTitle(title);
1360         hEff->SetName(Form("hMatchEff_cent%d_%s_%s", icent, AliPID::ParticleName(ipart), chargeName[icharge]));
1361         fileout->cd();
1362         hEff->Write();
1363       }
1364
1365   fileout->Close();
1366 }
1367
1368 TF1 *
1369 TOFmatchMC_geantflukaCorrection(Int_t ipart, Int_t icharge)
1370 {
1371
1372   if (ipart == 3 && icharge == kNegative) {
1373     TF1 *f = new TF1(Form("fGeantFluka_%s_%s", AliPID::ParticleName(ipart), chargeName[icharge]), "MatchingPtGeantFlukaCorrectionKaMinus(x)", 0., 5.);
1374     return f;
1375   }
1376   //  else if (ipart == 4 && icharge == kNegative) {
1377   //    TF1 *f = new TF1(Form("fGeantFluka_%s_%s", AliPID::ParticleName(ipart), chargeName[icharge]), "MatchingPtGeantFlukaCorrectionPrMinus(x)", 0., 5.);
1378   //  }
1379   else
1380     TF1 *f = new TF1(Form("fGeantFluka_%s_%s", AliPID::ParticleName(ipart), chargeName[icharge]), "MatchingPtGeantFlukaCorrectionNull(x)", 0., 5.);
1381
1382   return f;
1383 }
1384
1385
1386 TF1 *
1387 TOFmatchMC_geantflukaCorrection_acceptTRDin(Int_t ipart, Int_t icharge)
1388 {
1389
1390   if (ipart == 3 && icharge == kNegative) {
1391     TF1 *f = new TF1(Form("fGeantFluka_%s_%s", AliPID::ParticleName(ipart), chargeName[icharge]), "MatchingPtGeantFlukaCorrectionKaMinus_acceptTRDin(x)", 0., 5.);
1392     return f;
1393   }
1394   else if (ipart == 4 && icharge == kNegative) {
1395     TF1 *f = new TF1(Form("fGeantFluka_%s_%s", AliPID::ParticleName(ipart), chargeName[icharge]), "MatchingPtGeantFlukaCorrectionPrMinus_acceptTRDin(x)", 0., 5.);
1396   }
1397   else
1398     TF1 *f = new TF1(Form("fGeantFluka_%s_%s", AliPID::ParticleName(ipart), chargeName[icharge]), "MatchingPtGeantFlukaCorrectionNull(x)", 0., 5.);
1399
1400   return f;
1401 }
1402
1403
1404 TF1 *
1405 TOFmatchMC_geantflukaCorrection_rejectTRDin(Int_t ipart, Int_t icharge)
1406 {
1407
1408   if (ipart == 3 && icharge == kNegative) {
1409     TF1 *f = new TF1(Form("fGeantFluka_%s_%s", AliPID::ParticleName(ipart), chargeName[icharge]), "MatchingPtGeantFlukaCorrectionKaMinus_rejectTRDin(x)", 0., 5.);
1410     return f;
1411   }
1412   else if (ipart == 4 && icharge == kNegative) {
1413     TF1 *f = new TF1(Form("fGeantFluka_%s_%s", AliPID::ParticleName(ipart), chargeName[icharge]), "MatchingPtGeantFlukaCorrectionPrMinus_rejectTRDin(x)", 0., 5.);
1414   }
1415   else
1416     TF1 *f = new TF1(Form("fGeantFluka_%s_%s", AliPID::ParticleName(ipart), chargeName[icharge]), "MatchingPtGeantFlukaCorrectionNull(x)", 0., 5.);
1417
1418   return f;
1419 }
1420
1421
1422 Double_t
1423 MatchingPtGeantFlukaCorrectionNull(Double_t pTmc)
1424 {
1425   return 1.;
1426 }
1427
1428 Double_t
1429 MatchingPtGeantFlukaCorrectionPrMinus(Double_t pTmc, Int_t ntrd = 12)
1430 {
1431   Float_t ptTPCoutP =pTmc*(1-6.81059e-01*TMath::Exp(-pTmc*4.20094));
1432   Float_t scale = (ntrd * 0.14638485 + (18 - ntrd) * 0.02406) / 18.;
1433   return (TMath::Power(1 - 0.129758*TMath::Exp(-ptTPCoutP*0.679612),scale/0.03471));
1434 }
1435
1436 Double_t
1437 MatchingPtGeantFlukaCorrectionPrMinus_acceptTRDin(Double_t pTmc)
1438 {
1439   Float_t ptTPCoutP =pTmc*(1-6.81059e-01*TMath::Exp(-pTmc*4.20094));
1440   return (TMath::Power(1 - 0.129758*TMath::Exp(-ptTPCoutP*0.679612),0.14638485/0.03471));
1441 }
1442
1443 Double_t
1444 MatchingPtGeantFlukaCorrectionPrMinus_rejectTRDin(Double_t pTmc)
1445 {
1446   Float_t ptTPCoutP =pTmc*(1-6.81059e-01*TMath::Exp(-pTmc*4.20094));
1447   return (TMath::Power(1 - 0.129758*TMath::Exp(-ptTPCoutP*0.679612),0.02406/0.03471));
1448 }
1449
1450 Double_t
1451 MatchingPtGeantFlukaCorrectionKaMinus(Double_t pTmc, Int_t ntrd = 12)
1452 {
1453   Float_t ptTPCoutK=pTmc*(1- 3.37297e-03/pTmc/pTmc - 3.26544e-03/pTmc);
1454   Float_t scale = (ntrd * 0.14638485 + (18 - ntrd) * 0.02406) / 18.;
1455   return TMath::Min((TMath::Power(0.972865 + 0.0117093*ptTPCoutK,scale/0.03471)), 1.);
1456 }
1457
1458 Double_t
1459 MatchingPtGeantFlukaCorrectionKaMinus_acceptTRDin(Double_t pTmc)
1460 {
1461   Float_t ptTPCoutK=pTmc*(1- 3.37297e-03/pTmc/pTmc - 3.26544e-03/pTmc);
1462   return TMath::Min((TMath::Power(0.972865 + 0.0117093*ptTPCoutK,0.14638485/0.03471)), 1.);
1463 }
1464
1465 Double_t
1466 MatchingPtGeantFlukaCorrectionKaMinus_rejectTRDin(Double_t pTmc)
1467 {
1468   Float_t ptTPCoutK=pTmc*(1- 3.37297e-03/pTmc/pTmc - 3.26544e-03/pTmc);
1469   return TMath::Min((TMath::Power(0.972865 + 0.0117093*ptTPCoutK,0.02406/0.03471)), 1.);
1470 }