PWGCF/EBYE/PIDFluctuation/task/AliEbyEPidRatioEffCont.cxx

   1 /**************************************************************************
   2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
   3  *                                                                        *
   4  * Author: ALICE Offline.                                                 *
   5  * Contributors are mentioned in the code where appropriate.              *
   6  *                                                                        *
   7  * Permission to use, copy, modify and distribute this software and its   *
   8  * documentation strictly for non-commercial purposes is hereby granted   *
   9  * without fee, provided that the above copyright notice appears in all   *
  10  * copies and that both the copyright notice and this permission notice   *
  11  * appear in the supporting documentation. The authors make no claims     *
  12  * about the suitability of this software for any purpose. It is          *
  13  * provided "as is" without express or implied warranty.                  *
  14  **************************************************************************/
  15
  16
  17 //=========================================================================//
  18 //             AliEbyE Analysis for Particle Ratio Fluctuation             //
  19 //                   Deepika Rathee  | Satyajit Jena                       //
  20 //                   drathee@cern.ch | sjena@cern.ch                       //
  21 //                  Date: Wed Jul  9 18:38:30 CEST 2014                    //
  22 //          New approch to find particle ratio to reduce memory            //
  23 //                             (Test Only)                                 //
  24 //=========================================================================//
  25
  26 #include "TMath.h"
  27 #include "TAxis.h"
  28
  29 #include "AliESDEvent.h"
  30 #include "AliStack.h"
  31 #include "AliMCEvent.h"
  32 #include "AliESDtrackCuts.h"
  33 #include "AliAODEvent.h"
  34 #include "AliAODMCParticle.h"
  35
  36 #include "AliEbyEPidRatioEffCont.h"
  37
  38 using namespace std;
  39
  40
  41 ClassImp(AliEbyEPidRatioEffCont)
  42
  43 //________________________________________________________________________
  44 AliEbyEPidRatioEffCont::AliEbyEPidRatioEffCont() :
  45   AliEbyEPidRatioBase("EffCont", "EffCont"),
  46   fLabelsRec(NULL),
  47   fHnEffMc(NULL),
  48   fHnContMc(NULL),
  49   fHnEffRec(NULL),
  50   fHnContRec(NULL) {
  51   // Constructor
  52
  53   AliLog::SetClassDebugLevel("AliEbyEPidRatioEffCont",10);
  54 }
  55
  56 //________________________________________________________________________
  57 AliEbyEPidRatioEffCont::~AliEbyEPidRatioEffCont() {
  58   // Destructor
  59
  60   for (Int_t ii = 0; ii < 2; ++ii) {
  61     for (Int_t kk = 0; kk < 4; ++kk)
  62       if (fLabelsRec[ii][kk]) delete[] fLabelsRec[ii][kk];
  63     if (fLabelsRec[ii]) delete[] fLabelsRec[ii];
  64   }
  65   if (fLabelsRec) delete[] fLabelsRec;
  66
  67
  68 }
  69
  70 //________________________________________________________________________
  71 void AliEbyEPidRatioEffCont::Process() {
  72   // -- Process event
  73
  74   // -- Setup (clean, create and fill) MC labels
  75   FillMCLabels();
  76
  77   // -- Fill  MC histograms for efficiency studies
  78   FillMCEffHist();
  79
  80   return;
  81 }
  82
  83 //________________________________________________________________________
  84 void AliEbyEPidRatioEffCont::Init() {
  85   // -- Init eventwise
  86
  87   fLabelsRec = new Int_t**[2];
  88   for (Int_t ii = 0 ; ii < 2; ++ii) {
  89     fLabelsRec[ii] = new Int_t*[4];
  90     for (Int_t kk = 0 ; kk < 4; ++kk)
  91       fLabelsRec[ii][kk] = NULL;
  92   }
  93 }
  94
  95 //________________________________________________________________________
  96 void AliEbyEPidRatioEffCont::CreateHistograms() {
  97   // Copied from NetParticle class
  98   Int_t    binHnEff[10] = { AliEbyEPidRatioHelper::fgkfHistNBinsCent, 4,
  99                             AliEbyEPidRatioHelper::fgkfHistNBinsSign, 2,      2  ,      2  ,
 100                             AliEbyEPidRatioHelper::fgkfHistNBinsEta,
 101                             AliEbyEPidRatioHelper::fgkfHistNBinsRap,
 102                             AliEbyEPidRatioHelper::fgkfHistNBinsPhi,
 103                             AliEbyEPidRatioHelper::fgkfHistNBinsPt};
 104
 105   Double_t minHnEff[10] = {AliEbyEPidRatioHelper::fgkfHistRangeCent[0], -0.5,
 106                            AliEbyEPidRatioHelper::fgkfHistRangeSign[0], -0.5, -0.5, -0.5,
 107                            AliEbyEPidRatioHelper::fgkfHistRangeEta[0],
 108                            AliEbyEPidRatioHelper::fgkfHistRangeRap[0],
 109                            AliEbyEPidRatioHelper::fgkfHistRangePhi[0],
 110                            AliEbyEPidRatioHelper::fgkfHistRangePt[0]};
 111
 112
 113   Double_t maxHnEff[10] = {AliEbyEPidRatioHelper::fgkfHistRangeCent[1], 3.5,
 114                            AliEbyEPidRatioHelper::fgkfHistRangeSign[1], 1.5,      1.5,      1.5,
 115                            AliEbyEPidRatioHelper::fgkfHistRangeEta[1],
 116                            AliEbyEPidRatioHelper::fgkfHistRangeRap[1],
 117                            AliEbyEPidRatioHelper::fgkfHistRangePhi[1],
 118                            AliEbyEPidRatioHelper::fgkfHistRangePt[1]};
 119
 120   fHnEffMc    = new THnSparseF("hnEffMc", "cent:pid:SignMC:findable:recStatus:pidStatus:etaMC:yMC:phiMC:ptMC", 10, binHnEff, minHnEff, maxHnEff);
 121   fHnEffRec   = new THnSparseF("hnEffRec", "cent:pid:SignMC:findable:recStatus:pidStatus:etaRec:yRec:phiRec:ptRec", 10, binHnEff, minHnEff, maxHnEff);
 122   //fHnEffRecMc = new THnSparseF("hnEffRecMMc", "cent:pid:SignMC:findable:recStatus:pidStatus:deltaEta:deltaY:deltaPhi:deltaPt" 10, binHnEff, minHnEff, maxHnEff);
 123
 124   fHnEffMc->Sumw2();
 125   fHnEffRec->Sumw2();
 126
 127
 128
 129   fHnEffMc->GetAxis(0)->SetTitle("centrality");                   //  0-5|5-10|10-20|20-30|30-40|40-50|50-60|60-70|70-80|80-90 --> 10 bins
 130   fHnEffMc->GetAxis(1)->SetTitle("N_{ch}|N_{#pi}|N_{K}|N_{p}");                //  0 | 1 | 2 | 3
 131   fHnEffMc->GetAxis(2)->SetTitle("sign");                         //  -1 | 0 | +1
 132   fHnEffMc->GetAxis(3)->SetTitle("findable");                     //  0 not findable      |  1 findable
 133   fHnEffMc->GetAxis(4)->SetTitle("recStatus");                    //  0 not reconstructed |  1 reconstructed
 134   fHnEffMc->GetAxis(5)->SetTitle("recPid");                       //  0 not accepted      |  1 accepted
 135   fHnEffMc->GetAxis(6)->SetTitle("#eta_{MC}");                    //  eta  [-0.9, 0.9]
 136   fHnEffMc->GetAxis(7)->SetTitle("#it{y}_{MC}");                  //  rapidity  [-0.5, 0.5]
 137   fHnEffMc->GetAxis(8)->SetTitle("#varphi_{MC} (rad)");           //  phi  [ 0. , 2Pi]
 138   fHnEffMc->GetAxis(9)->SetTitle("#it{p}_{T,MC} (GeV/#it{c})");   //  pT   [ 0.2, 2.3]
 139
 140
 141   fHnEffRec->GetAxis(0)->SetTitle("centrality");                   //  0-5|5-10|10-20|20-30|30-40|40-50|50-60|60-70|70-80|80-90 --> 10 bins
 142   fHnEffRec->GetAxis(1)->SetTitle("N_{ch}|N_{#pi}|N_{K}|N_{p}");                //  0 | 1 | 2 | 3
 143   fHnEffRec->GetAxis(2)->SetTitle("sign");                         //  -1 | 0 | +1
 144   fHnEffRec->GetAxis(3)->SetTitle("findable");                     //  0 not findable      |  1 findable
 145   fHnEffRec->GetAxis(4)->SetTitle("recStatus");                    //  0 not reconstructed |  1 reconstructed
 146   fHnEffRec->GetAxis(5)->SetTitle("recPid");                       //  0 not accepted      |  1 accepted
 147   fHnEffRec->GetAxis(6)->SetTitle("#eta_{Rec}");                   //  eta  [-0.9, 0.9]
 148   fHnEffRec->GetAxis(7)->SetTitle("#it{y}_{Rec}");                 //  rapidity  [-0.5, 0.5]
 149   fHnEffRec->GetAxis(8)->SetTitle("#varphi_{Rec} (rad)");          //  phi  [ 0. , 2Pi]
 150   fHnEffRec->GetAxis(9)->SetTitle("#it{p}_{T,Rec} (GeV/#it{c})");  //  pt   [ 0.2, 2.3]
 151
 152
 153   // fHnEff->GetAxis(0)->SetTitle("centrality");                   //  0-5|5-10|10-20|20-30|30-40|40-50|50-60|60-70|70-80|80-90 --> 10 bins
 154   // fHnEff->GetAxis(1)->SetTitle("N_{ch}|N_{#pi}|N_{K}|N_{p}");                //  0 | 1 | 2 | 3
 155   // fHnEff->GetAxis(2)->SetTitle("sign");                         //  -1 | 0 | +1
 156   // fHnEff->GetAxis(3)->SetTitle("findable");                     //  0 not findable      |  1 findable
 157   // fHnEff->GetAxis(4)->SetTitle("recStatus");                    //  0 not reconstructed |  1 reconstructed
 158   // fHnEff->GetAxis(5)->SetTitle("recPid");                       //  0 not accepted      |  1 accepted
 159   // fHnEff->GetAxis(6)->SetTitle("#eta_{MC}-#eta_{Rec}");                      //  eta  [-0.9, 0.9]
 160   // fHnEff->GetAxis(7)->SetTitle("#it{y}_{MC}-#it{y}_{Rec}");                  //  rapidity  [-0.5, 0.5]
 161   // fHnEff->GetAxis(8)->SetTitle("#varphi_{MC}-#varphi_{Rec} (rad)");          //  phi  [ -2Pi , 2Pi]
 162   // fHnEff->GetAxis(9)->SetTitle("#it{p}_{T,MC}-#it{p}_{T,Rec} (GeV/#it{c})"); //  pt   [ -2.3, 2.3]
 163   //  fHnEff->GetAxis(10)->SetTitle("sign_{MC}-sign_{Rec}");                      //  -2 | 0 | +2
 164
 165
 166
 167
 168   fHelper->BinLogAxis(fHnEffMc, 9);
 169   fHelper->BinLogAxis(fHnEffRec, 9);
 170
 171   /*
 172      >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> Copied from NetParticle
 173
 174      creation -> findable -> reconstructed -> pid_TPC+TOF
 175         (1)         (2)          (3)            (4)
 176                                                                       ||   findable | recStatus | recPid
 177      1) all primary probeParticles_MC                                 ||      -           -         -
 178      2) all findable primary probeParticles_MC                        ||      x           -         -
 179      3) all reconstructed primary probeParticles_MC                   ||      x           x         -
 180      4) all reconstructed primary probeParticles_MC & recPid_TPC+TOF  ||      x           x         x
 181
 182      <<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
 183   */
 184
 185   // ------------------------------------------------------------------
 186   // -- Create THnSparse - Cont
 187   // ------------------------------------------------------------------
 188
 189   Int_t    binHnCont[8] = {AliEbyEPidRatioHelper::fgkfHistNBinsCent, 4,
 190                            AliEbyEPidRatioHelper::fgkfHistNBinsSign, 8,
 191                            AliEbyEPidRatioHelper::fgkfHistNBinsEta,
 192                            AliEbyEPidRatioHelper::fgkfHistNBinsRap,
 193                            AliEbyEPidRatioHelper::fgkfHistNBinsPhi,
 194                            AliEbyEPidRatioHelper::fgkfHistNBinsPt};
 195
 196   Double_t minHnCont[8] = {AliEbyEPidRatioHelper::fgkfHistRangeCent[0], -0.5,
 197                            AliEbyEPidRatioHelper::fgkfHistRangeSign[0],  0.5,
 198                            AliEbyEPidRatioHelper::fgkfHistRangeEta[0],
 199                            AliEbyEPidRatioHelper::fgkfHistRangeRap[0],
 200                            AliEbyEPidRatioHelper::fgkfHistRangePhi[0],
 201                            AliEbyEPidRatioHelper::fgkfHistRangePt[0]};
 202
 203   Double_t maxHnCont[8] = {AliEbyEPidRatioHelper::fgkfHistRangeCent[1], 3.5,
 204                            AliEbyEPidRatioHelper::fgkfHistRangeSign[1], 8.5,
 205                            AliEbyEPidRatioHelper::fgkfHistRangeEta[1],
 206                            AliEbyEPidRatioHelper::fgkfHistRangeRap[1],
 207                            AliEbyEPidRatioHelper::fgkfHistRangePhi[1],
 208                            AliEbyEPidRatioHelper::fgkfHistRangePt[1]};
 209
 210   fHnContMc  = new THnSparseF("hnContMc", "cent:pid:SignMC:contPart:etaMC:yMC:phiMC:ptMC",8, binHnCont, minHnCont, maxHnCont);
 211   fHnContRec = new THnSparseF("hnContRec", "cent:pid:SignRec:contPart:etaRec:yRec:phiRec:ptRec",8, binHnCont, minHnCont, maxHnCont);
 212
 213   fHnContMc->Sumw2();
 214   fHnContRec->Sumw2();
 215
 216   fHnContMc->GetAxis(0)->SetTitle("centrality");                   //  0-5|5-10|10-20|20-30|30-40|40-50|50-60|60-70|70-80|80-90 --> 10 bins
 217   fHnContMc->GetAxis(1)->SetTitle("N_{ch}|N_{#pi}|N_{K}|N_{p}");                //  0 | 1 | 2 | 3
 218   fHnContMc->GetAxis(2)->SetTitle("sign");                         //  -1 | 0 | +1
 219   fHnContMc->GetAxis(3)->SetTitle("contPart");                     //  1 pi | 2 K | 3 p | 4 e | 5 mu | 6 other | 7 p from WeakDecay | 8 p from Material
 220   fHnContMc->GetAxis(4)->SetTitle("#eta_{MC}");                    //  eta  [-0.9,0.9]
 221   fHnContMc->GetAxis(5)->SetTitle("#it{y}_{MC}");                  //  rapidity  [-0.5, 0.5]
 222   fHnContMc->GetAxis(6)->SetTitle("#varphi_{MC} (rad)");           //  phi  [ 0. ,2Pi]
 223   fHnContMc->GetAxis(7)->SetTitle("#it{p}_{T,MC} (GeV/#it{c})");   //  pT   [ 0.2,2.3]
 224
 225
 226   fHnContRec->GetAxis(0)->SetTitle("centrality");                   //  0-5|5-10|10-20|20-30|30-40|40-50|50-60|60-70|70-80|80-90 --> 10 bins
 227   fHnContRec->GetAxis(1)->SetTitle("N_{ch}|N_{#pi}|N_{K}|N_{p}");                //  0 | 1 | 2 | 3
 228   fHnContRec->GetAxis(2)->SetTitle("sign");                         //  -1 | 0 | +1
 229   fHnContRec->GetAxis(3)->SetTitle("contPart");                     //  1 pi | 2 K | 3 p | 4 e | 5 mu | 6 other | 7 p from WeakDecay | 8 p from Material
 230   fHnContRec->GetAxis(4)->SetTitle("#eta_{Rec}");                   //  eta  [-0.9, 0.9]
 231   fHnContRec->GetAxis(5)->SetTitle("#it{y}_{Rec}");                 //  rapidity  [-0.5, 0.5]
 232   fHnContRec->GetAxis(6)->SetTitle("#varphi_{Rec} (rad)");          //  phi  [ 0. , 2Pi]
 233   fHnContRec->GetAxis(7)->SetTitle("#it{p}_{T,Rec} (GeV/#it{c})"); //  pt   [ 0.2, 2.3]
 234
 235
 236   //  fHnCont->GetAxis(12)->SetTitle("#eta_{MC}-#eta_{Rec}");                      //  eta  [-0.9, 0.9]
 237   // fHnCont->GetAxis(13)->SetTitle("#it{y}_{MC}-#it{y}_{Rec}");                  //  rapidity  [-0.5, 0.5]
 238   // fHnCont->GetAxis(14)->SetTitle("#varphi_{MC}-#varphi_{Rec} (rad)");          //  phi  [ -2Pi , 2Pi]
 239   // fHnCont->GetAxis(15)->SetTitle("#it{p}_{T,MC}-#it{p}_{T,Rec} (GeV/#it{c})"); //  pt   [ -2.3, 2.3]
 240   // fHnCont->GetAxis(16)->SetTitle("sign_{MC}-sign_{Rec}");                      //  -2 | 0 | +2
 241   // fHnCont->GetAxis(17)->SetTitle("N_{ch}|N_{#pi}|N_{K}|N_{p}");                //  0 | 1 | 2 | 3
 242
 243   fHelper->BinLogAxis(fHnContMc,  7);
 244   fHelper->BinLogAxis(fHnContRec, 7);
 245
 246   return;
 247 }
 248
 249 //________________________________________________________________________
 250 Int_t AliEbyEPidRatioEffCont::Setup() {
 251   // -- Setup eventwise
 252
 253   // -- Create label arrays
 254   for(Int_t i = 0; i < 4; i++) {
 255   fLabelsRec[0][i] = new Int_t[fNTracks];
 256   if(!fLabelsRec[0][i]) {
 257     AliError("Cannot create fLabelsRec[0]");
 258     return -1;
 259   }
 260
 261   fLabelsRec[1][i] = new Int_t[fNTracks];
 262   if(!fLabelsRec[1][i]) {
 263     AliError("Cannot create fLabelsRec[1] for PID");
 264     return -1;
 265   }
 266
 267   for(Int_t ii = 0; ii < fNTracks; ++ii) {
 268     fLabelsRec[0][i][ii] = 0;
 269     fLabelsRec[1][i][ii] = 0;
 270   }
 271   }
 272   return 0;
 273 }
 274
 275 //________________________________________________________________________
 276 void AliEbyEPidRatioEffCont::Reset() {
 277   // -- Reset eventwise
 278   for(Int_t i = 0; i < 4; i++) {
 279     for (Int_t ii = 0; ii < 2 ; ++ii) {
 280       if (fLabelsRec[ii][i])
 281         delete[] fLabelsRec[ii][i];
 282       fLabelsRec[ii][i] = NULL;
 283     }
 284   }
 285 }
 286
 287 //________________________________________________________________________
 288 void AliEbyEPidRatioEffCont::FillMCLabels() {
 289   Float_t etaRange[2];
 290   fESDTrackCuts->GetEtaRange(etaRange[0],etaRange[1]);
 291
 292   Float_t ptRange[2];
 293   fESDTrackCuts->GetPtRange(ptRange[0],ptRange[1]);
 294
 295   // -- Track Loop
 296   for (Int_t idxTrack = 0; idxTrack < fNTracks; ++idxTrack) {
 297
 298     AliVTrack *track = (fESD) ? static_cast<AliVTrack*>(fESD->GetTrack(idxTrack)) : static_cast<AliVTrack*>(fAOD->GetTrack(idxTrack));
 299
 300     // -- Check if track is accepted for basic parameters
 301     if (!fHelper->IsTrackAcceptedBasicCharged(track))
 302       continue;
 303
 304     // -- Check if accepted - ESD
 305     if (fESD && !fESDTrackCuts->AcceptTrack(dynamic_cast<AliESDtrack*>(track)))
 306       continue;
 307
 308     // -- Check if accepted - AOD
 309     if (fAOD){
 310       AliAODTrack * trackAOD = dynamic_cast<AliAODTrack*>(track);
 311
 312       if (!trackAOD) {
 313         AliError("Pointer to dynamic_cast<AliAODTrack*>(track) = ZERO");
 314         continue;
 315       }
 316       if (!trackAOD->TestFilterBit(fAODtrackCutBit))
 317         continue;
 318
 319       // -- Check if in pT and eta range (is done in ESDTrackCuts for ESDs)
 320       if(!(track->Pt() > ptRange[0] && track->Pt() <= ptRange[1] && TMath::Abs(track->Eta()) <= etaRange[1]))
 321         continue;
 322     }
 323
 324     Int_t gPdgCode = 0;
 325
 326     Int_t iPid = 0;
 327     Double_t pid[3];
 328     if      (fHelper->IsTrackAcceptedPID(track, pid, (AliPID::kPion)))  {  iPid = 1; gPdgCode = 211;}
 329     else if (fHelper->IsTrackAcceptedPID(track, pid, (AliPID::kKaon)))  {  iPid = 2; gPdgCode = 321;}
 330     else if (fHelper->IsTrackAcceptedPID(track, pid, (AliPID::kProton))){  iPid = 3; gPdgCode = 2212;}
 331     else iPid = 0;
 332
 333     //  cout << " --- EFF ---- " << iPid << "  " << gPdgCode << endl;
 334
 335     Double_t yP;
 336     if (!fHelper->IsTrackAcceptedRapidity(track, yP, iPid))
 337       continue;
 338
 339     if (!fHelper->IsTrackAcceptedDCA(track))
 340       continue;
 341
 342     Int_t label  = TMath::Abs(track->GetLabel());
 343
 344     // -- Fill Label of all reconstructed
 345     if(iPid != 0) fLabelsRec[0][0][idxTrack]    = label;
 346     fLabelsRec[0][iPid][idxTrack] = label;
 347
 348     // -- Fill Label of all reconstructed && recPid_TPC+TOF
 349     if(iPid != 0) fLabelsRec[1][0][idxTrack]    = label;
 350     fLabelsRec[1][iPid][idxTrack] = label;
 351
 352     // -- Check for contamination and fill contamination THnSparse
 353     CheckContTrack(track, iPid, gPdgCode);
 354     if(iPid != 0) CheckContTrack(track, 0, 0);
 355
 356   } // for (Int_t idxTrack = 0; idxTrack < fNTracks; ++idxTrack) {
 357
 358   return;
 359 }
 360
 361 //________________________________________________________________________
 362 void AliEbyEPidRatioEffCont::CheckContTrack(AliVTrack *track, Int_t iPid, Int_t gPdgCode) {
 363   Int_t label     = TMath::Abs(track->GetLabel());
 364   Float_t signRec = track->Charge();
 365
 366
 367   AliVParticle* particle = (fESD) ? fMCEvent->GetTrack(label) : static_cast<AliVParticle*>(fArrayMC->At(label));
 368   if (!particle)
 369     return;
 370
 371   Bool_t isPhysicalPrimary = (fESD) ? fStack->IsPhysicalPrimary(label): (static_cast<AliAODMCParticle*>(particle))->IsPhysicalPrimary();
 372   if (iPid == 0) {
 373     if (particle->PdgCode() == (signRec*gPdgCode))
 374       if (isPhysicalPrimary)
 375         return;
 376   }
 377   else {
 378     if (isPhysicalPrimary)
 379       return;
 380   }
 381
 382   // -- Check if secondaries from material or weak decay
 383   Bool_t isSecondaryFromWeakDecay = (fESD) ? fStack->IsSecondaryFromWeakDecay(label) : (static_cast<AliAODMCParticle*>(particle))->IsSecondaryFromWeakDecay();
 384   Bool_t isSecondaryFromMaterial  = (fESD) ? fStack->IsSecondaryFromMaterial(label)  : (static_cast<AliAODMCParticle*>(particle))->IsSecondaryFromMaterial();
 385
 386   // -- Get PDG Charge of contaminating particle
 387   Float_t signMC = 0.;
 388   if      (particle->Charge() == 0.) signMC =  0.;
 389   else if (particle->Charge() <  0.) signMC = -1.;
 390   else if (particle->Charge() >  0.) signMC =  1.;
 391
 392   // -- Get contaminating particle
 393   Float_t contPart = 0;
 394   if        (isSecondaryFromWeakDecay)                contPart = 7; // probeParticle from WeakDecay
 395   else if   (isSecondaryFromMaterial)                 contPart = 8; // probeParticle from Material
 396   else {
 397     if      (TMath::Abs(particle->PdgCode()) ==  211) contPart = 1; // pion
 398     else if (TMath::Abs(particle->PdgCode()) ==  321) contPart = 2; // kaon
 399     else if (TMath::Abs(particle->PdgCode()) == 2212) contPart = 3; // proton
 400     else if (TMath::Abs(particle->PdgCode()) ==   11) contPart = 4; // electron
 401     else if (TMath::Abs(particle->PdgCode()) ==   13) contPart = 5; // muon
 402     else                                              contPart = 6; // other
 403   }
 404
 405   // cout << " --- CONT ---- " << iPid << "  " << gPdgCode << endl;
 406
 407   // -- Get Reconstructed y
 408   //    yRec = y for identified particles | yRec = eta for charged particles
 409   Double_t yRec  = 0.;
 410   fHelper->IsTrackAcceptedRapidity(track, yRec, iPid);
 411
 412   Double_t deltaPhi = particle->Phi()-track->Phi();
 413   if (TMath::Abs(deltaPhi) > TMath::TwoPi()) {
 414     if (deltaPhi < 0)
 415       deltaPhi += TMath::TwoPi();
 416     else
 417       deltaPhi -= TMath::TwoPi();
 418   }
 419
 420   Double_t hnContMc[8]  = {fCentralityBin,iPid,signMC,contPart,particle->Eta(),particle->Y(),particle->Phi(),particle->Pt()};
 421   Double_t hnContRec[8] = {fCentralityBin,iPid,signRec,contPart, track->Eta(),yRec,track->Phi(),track->Pt()};
 422   fHnContMc->Fill(hnContMc);
 423   fHnContRec->Fill(hnContRec);
 424
 425 }
 426
 427 //________________________________________________________________________
 428 void AliEbyEPidRatioEffCont::FillMCEffHist() {
 429   // Fill efficiency THnSparse for ESDs
 430
 431   Float_t etaRange[2];
 432   fESDTrackCuts->GetEtaRange(etaRange[0],etaRange[1]);
 433
 434   Int_t nPart  = (fESD) ? fStack->GetNprimary() : fArrayMC->GetEntriesFast();
 435
 436   for (Int_t idxMC = 0; idxMC < nPart; ++idxMC) {
 437     AliVParticle* particle = (fESD) ? fMCEvent->GetTrack(idxMC) : static_cast<AliVParticle*>(fArrayMC->At(idxMC));
 438
 439     // -- Check basic MC properties -> charged physical primary
 440     if (!fHelper->IsParticleAcceptedBasicCharged(particle, idxMC))
 441       continue;
 442
 443     Int_t iPid = 0;
 444     Int_t gPdgCode = 0;
 445     if ( particle->PdgCode()      ==  211 ) {  iPid = 1; gPdgCode = 211;}
 446     else if ( particle->PdgCode() ==  321 ) {  iPid = 2; gPdgCode = 321;}
 447     else if ( particle->PdgCode() == 2212 ) {  iPid = 3; gPdgCode = 2212;}
 448     else {iPid = 0; gPdgCode = 0;}
 449
 450     // -- Check if accepted in rapidity window -- for identified particles
 451     Double_t yMC;
 452     if (iPid != 0 && !fHelper->IsParticleAcceptedRapidity(particle, yMC, iPid))
 453       continue;
 454
 455     // -- Check if accepted in eta window -- for charged particles
 456     if (iPid == 0 && TMath::Abs(particle->Eta()) > etaRange[1])
 457       continue;
 458
 459     // -- Check if probeParticle / anti-probeParticle
 460     //    > skip check if PID is not required
 461     if (iPid == 0 && TMath::Abs(particle->PdgCode()) != gPdgCode)
 462       continue;
 463
 464     // -- Get sign of particle
 465     Float_t signMC    = (particle->PdgCode() < 0) ? -1. : 1.;
 466
 467     // -- Get if particle is findable --- not availible for AODs yet
 468     Float_t findable  = (fESD) ? Float_t(fHelper->IsParticleFindable(idxMC)) : 1.;
 469
 470     // -- Get recStatus and pidStatus
 471     Float_t recStatus = 0.;
 472     Float_t recPid    = 0.;
 473
 474     // -- Get Reconstructed values
 475     Float_t etaRec  = 0.;
 476     Float_t phiRec  = 0.;
 477     Float_t ptRec   = 0.;
 478     Double_t yRec   = 0.;
 479     Float_t signRec = 0.;
 480
 481     // -- Loop over all labels
 482     for (Int_t idxRec=0; idxRec < fNTracks; ++idxRec) {
 483       if (idxMC == fLabelsRec[0][iPid][idxRec]) {
 484         recStatus = 1.;
 485
 486         if (idxMC == fLabelsRec[1][iPid][idxRec])
 487           recPid = 1.;
 488
 489         AliVTrack *track = NULL;
 490         if(fESD)
 491           track = fESD->GetTrack(idxRec);
 492         else if(fAOD)
 493           track = fAOD->GetTrack(idxRec);
 494
 495         if (track) {
 496           // if no track present (which should not happen)
 497           // -> pt = 0. , which is not in the looked at range
 498
 499           // -- Get Reconstructed values
 500           etaRec  = track->Eta();
 501           phiRec  = track->Phi();
 502           ptRec   = track->Pt();
 503           signRec = track->Charge();
 504           fHelper->IsTrackAcceptedRapidity(track, yRec, iPid); // yRec = y for identified particles | yRec = eta for charged particles
 505         }
 506         break;
 507       }
 508     } // for (Int_t idxRec=0; idxRec < fNTracks; ++idxRec) {
 509
 510     Double_t deltaPhi = particle->Phi()-phiRec;
 511     if (TMath::Abs(deltaPhi) > TMath::TwoPi()) {
 512       if (deltaPhi < 0)
 513         deltaPhi += TMath::TwoPi();
 514       else
 515         deltaPhi -= TMath::TwoPi();
 516     }
 517
 518     if(iPid != 0) {
 519       Double_t hnEffMc[10]  = {fCentralityBin,0,signMC,findable, recStatus,recPid,particle->Eta(), particle->Y(), particle->Phi(),particle->Pt()};
 520       Double_t hnEffRec[10] = {fCentralityBin,0,signRec,findable, recStatus,recPid,etaRec, yRec, phiRec, ptRec};
 521       fHnEffMc->Fill(hnEffMc);
 522       fHnEffRec->Fill(hnEffRec);
 523     }
 524     Double_t hnEffMc[10]  = {fCentralityBin,iPid,signMC,findable, recStatus,recPid,particle->Eta(), particle->Y(), particle->Phi(),particle->Pt()};
 525     Double_t hnEffRec[10] = {fCentralityBin,iPid,signRec,findable, recStatus,recPid,etaRec, yRec, phiRec, ptRec};
 526     fHnEffMc->Fill(hnEffMc);
 527     fHnEffRec->Fill(hnEffRec);
 528
 529
 530
 531   } // for (Int_t idxMC = 0; idxMC < nPart; ++idxMC) {
 532
 533   return;
 534 }