PWG4/AliCaloTrackMCReader.cxx

   1 /**************************************************************************
   2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
   3  *                                                                        *
   4  * Author: The ALICE Off-line Project.                                    *
   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 /* $Id:  $ */
  16
  17 /* History of cvs commits:
  18  *
  19  * $Log$
  20  *
  21  *
  22  */
  23
  24 //_________________________________________________________________________
  25 // Class for reading data (Kinematics) in order to do prompt gamma
  26 // or other particle identification and correlations
  27 //
  28 //*-- Author: Gustavo Conesa (LNF-INFN)
  29 //////////////////////////////////////////////////////////////////////////////
  30
  31
  32 // --- ROOT system ---
  33
  34 #include <TParticle.h>
  35 #include <TH2.h>
  36 #include <TChain.h>
  37 #include <TRandom.h>
  38 #include <TClonesArray.h>
  39
  40 //---- ANALYSIS system ----
  41 #include "AliCaloTrackMCReader.h"
  42 #include "Riostream.h"
  43 #include "AliLog.h"
  44 #include "AliGenEventHeader.h"
  45 #include "AliStack.h"
  46 #include "AliAODCaloCluster.h"
  47 #include "AliAODTrack.h"
  48
  49 ClassImp(AliCaloTrackMCReader)
  50
  51 //____________________________________________________________________________
  52 AliCaloTrackMCReader::AliCaloTrackMCReader() :
  53   AliCaloTrackReader(), fDecayPi0(0),
  54   fNeutralParticlesArray(0x0),    fChargedParticlesArray(0x0),
  55   fStatusArray(0x0), fKeepAllStatus(0), fClonesArrayType(0)
  56 {
  57   //Ctor
  58
  59   //Initialize parameters
  60   InitParameters();
  61   fDataType = kMC;
  62
  63 }
  64
  65 //____________________________________________________________________________
  66 AliCaloTrackMCReader::AliCaloTrackMCReader(const AliCaloTrackMCReader & g) :
  67   AliCaloTrackReader(g), fDecayPi0(g.fDecayPi0),
  68   fNeutralParticlesArray(g.fNeutralParticlesArray?new TArrayI(*g.fNeutralParticlesArray):0x0),
  69   fChargedParticlesArray(g.fChargedParticlesArray?new TArrayI(*g.fChargedParticlesArray):0x0),
  70   fStatusArray(g.fStatusArray?new TArrayI(*g.fStatusArray):0x0),
  71   fKeepAllStatus(g.fKeepAllStatus), fClonesArrayType(g.fClonesArrayType)
  72 {
  73   // cpy ctor
  74 }
  75
  76 //_________________________________________________________________________
  77 AliCaloTrackMCReader & AliCaloTrackMCReader::operator = (const AliCaloTrackMCReader & source)
  78 {
  79   // assignment operator
  80
  81   if(&source == this) return *this;
  82
  83   fDecayPi0 = source.fDecayPi0;
  84
  85   delete fChargedParticlesArray;
  86   fChargedParticlesArray = source.fChargedParticlesArray?new TArrayI(*source.fChargedParticlesArray):0x0;
  87
  88   delete fNeutralParticlesArray;
  89   fNeutralParticlesArray = source.fNeutralParticlesArray?new TArrayI(*source.fNeutralParticlesArray):0x0;
  90
  91   delete fStatusArray;
  92   fStatusArray = source.fStatusArray?new TArrayI(*source.fStatusArray):0x0;
  93
  94   fKeepAllStatus = source.fKeepAllStatus ;
  95   fClonesArrayType = source.fClonesArrayType ;
  96
  97   return *this;
  98
  99 }
 100
 101 //_________________________________
 102 AliCaloTrackMCReader::~AliCaloTrackMCReader() {
 103   //Dtor
 104
 105   if(fChargedParticlesArray) delete fChargedParticlesArray ;
 106   if(fNeutralParticlesArray) delete fNeutralParticlesArray ;
 107   if(fStatusArray) delete fStatusArray ;
 108
 109 }
 110
 111 //____________________________________________________________________________
 112 void AliCaloTrackMCReader::GetVertex(Double_t  v[3]) {
 113   //Return vertex position
 114
 115   TArrayF pv;
 116   GetGenEventHeader()->PrimaryVertex(pv);
 117   v[0]=pv.At(0);
 118   v[1]=pv.At(1);
 119   v[2]=pv.At(2);
 120
 121 }
 122
 123
 124 //_______________________________________________________________
 125 void AliCaloTrackMCReader::InitParameters()
 126 {
 127
 128   //Initialize the parameters of the analysis.
 129
 130   fDecayPi0 = kFALSE;
 131
 132   fChargedParticlesArray = new TArrayI(1);
 133   fChargedParticlesArray->SetAt(11,0);
 134   //Int_t pdgarray[]={12,14,16};// skip neutrinos
 135   //fNeutralParticlesArray = new TArrayI(3, pdgarray);
 136   fNeutralParticlesArray = new TArrayI(3);
 137   fNeutralParticlesArray->SetAt(12,0); fNeutralParticlesArray->SetAt(14,1); fNeutralParticlesArray->SetAt(16,2);
 138   fStatusArray = new TArrayI(1);
 139   fStatusArray->SetAt(1,0);
 140
 141   fKeepAllStatus = kFALSE;
 142   fClonesArrayType = kAliAOD ;
 143
 144 }
 145
 146 //____________________________________________________________________________
 147 void  AliCaloTrackMCReader::FillCalorimeters(const Int_t iParticle, TParticle* particle, TLorentzVector momentum,
 148                                              Int_t &indexPHOS, Int_t &indexEMCAL){
 149   //Fill AODCaloClusters or TParticles lists of PHOS or EMCAL
 150
 151   //In PHOS
 152   if(fFillPHOS && fFidutialCut->IsInFidutialCut(momentum,"PHOS") && momentum.Pt() > fPHOSPtMin){
 153
 154     if(fClonesArrayType == kTParticle) new((*fAODPHOS)[indexPHOS++])       TParticle(*particle) ;
 155     else{
 156
 157       Char_t ttype= AliAODCluster::kPHOSNeutral;
 158       Int_t labels[] = {iParticle};
 159       Float_t x[] = {momentum.X(), momentum.Y(), momentum.Z()};
 160       AliAODCaloCluster *calo = new((*fAODPHOS)[indexPHOS++])
 161         AliAODCaloCluster(iParticle,1,labels,momentum.E(), x, NULL, ttype, 0);
 162       SetCaloClusterPID(particle->GetPdgCode(),calo) ;
 163     }
 164   }
 165   //In EMCAL
 166   else if(fFillEMCAL && fFidutialCut->IsInFidutialCut(momentum,"EMCAL") && momentum.Pt() > fEMCALPtMin){
 167     //cout<<"Reader : E "<<momentum.E()<<" ; Phi"<<momentum.Phi()<< " ; Eta "<<momentum.Eta()<<endl;
 168     if(fClonesArrayType == kTParticle) new((*fAODEMCAL)[indexEMCAL++])       TParticle(*particle) ;
 169     else{
 170       Char_t ttype= AliAODCluster::kEMCALClusterv1;
 171       Int_t labels[] = {iParticle};
 172       Float_t x[] = {momentum.X(), momentum.Y(), momentum.Z()};
 173       AliAODCaloCluster *calo = new((*fAODEMCAL)[indexEMCAL++])
 174         AliAODCaloCluster(iParticle,1,labels,momentum.E(), x, NULL, ttype, 0);
 175       SetCaloClusterPID(particle->GetPdgCode(),calo) ;
 176     }
 177   }
 178 }
 179
 180 //____________________________________________________________________________
 181 void AliCaloTrackMCReader::FillInputEvent()
 182 {
 183   //Create list of particles from EMCAL, PHOS and CTS.
 184
 185   if(fClonesArrayType == kTParticle){
 186     fAODCTS = new TClonesArray("TParticle",0);
 187     fAODEMCAL = new TClonesArray("TParticle",0);
 188     fAODPHOS = new TClonesArray("TParticle",0);
 189   }
 190   else if(fClonesArrayType == kAliAOD){
 191    fAODCTS = new TClonesArray("AliAODTrack",0);
 192    fAODPHOS = new TClonesArray("AliAODCaloCluster",0);
 193    fAODEMCAL = new TClonesArray("AliAODCaloCluster",0);
 194   }
 195   else {AliFatal("Wrong clones type");}
 196
 197
 198   Int_t indexCh      = 0 ;
 199   Int_t indexEMCAL = 0 ;
 200   Int_t indexPHOS = 0 ;
 201
 202   Int_t iParticle = 0 ;
 203   Double_t charge = 0.;
 204
 205   for (iParticle=0 ; iParticle <  GetStack()->GetNprimary() ; iParticle++) {
 206     TParticle * particle = GetStack()->Particle(iParticle);
 207     TLorentzVector momentum;
 208     Float_t p[3];
 209     Float_t x[3];
 210     Int_t pdg = particle->GetPdgCode();
 211
 212     //Keep particles with a given status
 213     if(KeepParticleWithStatus(particle->GetStatusCode()) && (particle->Pt() > 0) ){
 214
 215       charge = TDatabasePDG::Instance()->GetParticle(pdg)->Charge();
 216       particle->Momentum(momentum);
 217       //---------- Charged particles ----------------------
 218       if((charge != 0) && (momentum.Pt() > fCTSPtMin) && (fFidutialCut->IsInFidutialCut(momentum,"CTS"))){
 219         if(fFillCTS){
 220           //Particles in CTS acceptance
 221           if(fDebug > 3 && momentum.Pt() > 0.2)printf("Fill MC CTS :: Selected tracks E %3.2f, pt %3.2f, phi %3.2f, eta %3.2f\n",
 222                                                       momentum.E(),momentum.Pt(),momentum.Phi()*TMath::RadToDeg(),momentum.Eta());
 223
 224           if(fClonesArrayType == kTParticle) new((*fAODCTS)[indexCh++])       TParticle(*particle) ;
 225           else{
 226             x[0] = particle->Vx(); x[1] = particle->Vy(); x[2] = particle->Vz();
 227             p[0] = particle->Px(); p[1] = particle->Py(); p[2] = particle->Pz();
 228             AliAODTrack *aodTrack = new((*fAODCTS)[indexCh++])
 229               AliAODTrack(0, iParticle, p, kTRUE, x, kFALSE,NULL, 0, 0,
 230                           NULL,
 231                           0x0,//primary,
 232                         kFALSE, // No fit performed
 233                           kFALSE, // No fit performed
 234                           AliAODTrack::kPrimary,
 235                           0);
 236             SetTrackChargeAndPID(pdg, aodTrack);
 237           }
 238         }
 239         //Keep some charged particles in calorimeters lists
 240         if((fFillPHOS || fFillEMCAL) && KeepChargedParticles(pdg)) FillCalorimeters(iParticle, particle, momentum, indexPHOS, indexEMCAL);
 241
 242       }//Charged
 243
 244        //-------------Neutral particles ----------------------
 245       else if(charge == 0 && (fFillPHOS || fFillEMCAL)){
 246         //Skip neutrinos or other neutral particles
 247         if(SkipNeutralParticles(pdg) || particle->IsPrimary()) continue ; // protection added (MG)
 248
 249         //Fill particle/calocluster arrays
 250         if(!fDecayPi0) FillCalorimeters(iParticle, particle, momentum, indexPHOS, indexEMCAL);
 251         else {
 252           //Sometimes pi0 are stable for the generator, if needed decay it by hand
 253           if(pdg == 111 ){
 254             if(momentum.Pt() >  fPHOSPtMin || momentum.Pt() >  fEMCALPtMin){
 255               TLorentzVector lvGamma1, lvGamma2 ;
 256               //Double_t angle = 0;
 257
 258               //Decay
 259               MakePi0Decay(momentum,lvGamma1,lvGamma2);//,angle);
 260
 261               //Check if Pi0 is in the acceptance of the calorimeters, if aperture angle is small, keep it
 262               TParticle * pPhoton1 = new TParticle(22,1,iParticle,0,0,0,lvGamma1.Px(),lvGamma1.Py(),
 263                                                    lvGamma1.Pz(),lvGamma1.E(),0,0,0,0);
 264               TParticle * pPhoton2 = new TParticle(22,1,iParticle,0,0,0,lvGamma2.Px(),lvGamma2.Py(),
 265                                                    lvGamma2.Pz(),lvGamma2.E(),0,0,0,0);
 266               //Fill particle/calocluster arrays
 267               FillCalorimeters(iParticle,pPhoton1,lvGamma1, indexPHOS, indexEMCAL);
 268               FillCalorimeters(iParticle,pPhoton2,lvGamma2, indexPHOS, indexEMCAL);
 269             }//pt cut
 270           }//pi0
 271           else FillCalorimeters(iParticle,particle, momentum, indexPHOS, indexEMCAL); //Add the rest
 272         }
 273       }//neutral particles
 274      } //particle with correct status
 275    }//particle loop
 276
 277 }
 278
 279 //________________________________________________________________
 280 Bool_t AliCaloTrackMCReader::KeepParticleWithStatus(Int_t status) const {
 281   //Check if status is equal to one of the  list
 282   //These particles will be used in analysis.
 283   if(!fKeepAllStatus){
 284     for(Int_t i= 0; i < fStatusArray->GetSize(); i++)
 285       if(status ==  fStatusArray->At(i)) return kTRUE ;
 286
 287     return kFALSE;
 288
 289   }
 290   else
 291     return kTRUE ;
 292 }
 293
 294 //________________________________________________________________
 295 Bool_t AliCaloTrackMCReader::KeepChargedParticles(Int_t pdg) const {
 296   //Check if pdg is equal to one of the charged particles list
 297   //These particles will be added to the calorimeters lists.
 298
 299   for(Int_t i= 0; i < fChargedParticlesArray->GetSize(); i++)
 300     if(TMath::Abs(pdg) ==  fChargedParticlesArray->At(i)) return kTRUE ;
 301
 302   return kFALSE;
 303
 304 }
 305
 306 //________________________________________________________________
 307 void AliCaloTrackMCReader::Print(const Option_t * opt) const
 308 {
 309   //Print some relevant parameters set for the analysis
 310   if(! opt)
 311     return;
 312
 313   Info("**** Print **** ", "%s %s", GetName(), GetTitle() ) ;
 314
 315   printf("Decay Pi0?          : %d\n", fDecayPi0) ;
 316   printf("TClonesArray type   : %d\n", fClonesArrayType) ;
 317   printf("Keep all status?    : %d\n", fKeepAllStatus) ;
 318
 319   if(!fKeepAllStatus) printf("Keep particles with status : ");
 320   for(Int_t i= 0; i < fStatusArray->GetSize(); i++)
 321     printf(" %d ; ", fStatusArray->At(i));
 322   printf("\n");
 323
 324   printf("Skip neutral particles in calo : ");
 325   for(Int_t i= 0; i < fNeutralParticlesArray->GetSize(); i++)
 326     printf(" %d ; ", fNeutralParticlesArray->At(i));
 327   printf("\n");
 328
 329   printf("Keep charged particles in calo : ");
 330   for(Int_t i= 0; i < fChargedParticlesArray->GetSize(); i++)
 331     printf(" %d ; ", fChargedParticlesArray->At(i));
 332   printf("\n");
 333
 334 }
 335
 336 //____________________________________________________________________________
 337 void AliCaloTrackMCReader::MakePi0Decay(TLorentzVector &p0, TLorentzVector &p1,
 338                                 TLorentzVector &p2){//, Double_t &angle) {
 339   // Perform isotropic decay pi0 -> 2 photons
 340   // p0 is pi0 4-momentum (inut)
 341   // p1 and p2 are photon 4-momenta (output)
 342   //  cout<<"Boost vector"<<endl;
 343   Double_t mPi0 = TDatabasePDG::Instance()->GetParticle(111)->Mass();
 344   TVector3 b = p0.BoostVector();
 345   //cout<<"Parameters"<<endl;
 346   //Double_t mPi0   = p0.M();
 347   Double_t phi    = TMath::TwoPi() * gRandom->Rndm();
 348   Double_t cosThe = 2 * gRandom->Rndm() - 1;
 349   Double_t cosPhi = TMath::Cos(phi);
 350   Double_t sinPhi = TMath::Sin(phi);
 351   Double_t sinThe = TMath::Sqrt(1-cosThe*cosThe);
 352   Double_t ePi0   = mPi0/2.;
 353   //cout<<"ePi0 "<<ePi0<<endl;
 354   //cout<<"Components"<<endl;
 355   p1.SetPx(+ePi0*cosPhi*sinThe);
 356   p1.SetPy(+ePi0*sinPhi*sinThe);
 357   p1.SetPz(+ePi0*cosThe);
 358   p1.SetE(ePi0);
 359   //cout<<"p1: "<<p1.Px()<<" "<<p1.Py()<<" "<<p1.Pz()<<" "<<p1.E()<<endl;
 360   //cout<<"p1 Mass: "<<p1.Px()*p1.Px()+p1.Py()*p1.Py()+p1.Pz()*p1.Pz()-p1.E()*p1.E()<<endl;
 361   p2.SetPx(-ePi0*cosPhi*sinThe);
 362   p2.SetPy(-ePi0*sinPhi*sinThe);
 363   p2.SetPz(-ePi0*cosThe);
 364   p2.SetE(ePi0);
 365   //cout<<"p2: "<<p2.Px()<<" "<<p2.Py()<<" "<<p2.Pz()<<" "<<p2.E()<<endl;
 366   //cout<<"p2 Mass: "<<p2.Px()*p2.Px()+p2.Py()*p2.Py()+p2.Pz()*p2.Pz()-p2.E()*p2.E()<<endl;
 367   //cout<<"Boost "<<b.X()<<" "<<b.Y()<<" "<<b.Z()<<endl;
 368   p1.Boost(b);
 369   //cout<<"p1: "<<p1.Px()<<" "<<p1.Py()<<" "<<p1.Pz()<<" "<<p1.E()<<endl;
 370   p2.Boost(b);
 371   //cout<<"p2: "<<p2.Px()<<" "<<p2.Py()<<" "<<p2.Pz()<<" "<<p2.E()<<endl;
 372   //cout<<"angle"<<endl;
 373   //angle = p1.Angle(p2.Vect());
 374   //cout<<angle<<endl;
 375 }
 376
 377 //____________________________________________________________________________
 378 void AliCaloTrackMCReader::SetInputEvent(TObject* /*esd*/, TObject* /*aod*/, TObject* mc) {
 379   // Connect the data pointer
 380   SetMC((AliMCEvent*) mc);
 381 }
 382
 383
 384 //________________________________________________________________
 385 Bool_t AliCaloTrackMCReader::SkipNeutralParticles(Int_t pdg) const {
 386   //Check if pdg is equal to one of the neutral particles list
 387   //These particles will be skipped from analysis.
 388
 389   for(Int_t i= 0; i < fNeutralParticlesArray->GetSize(); i++)
 390     if(TMath::Abs(pdg) ==  fNeutralParticlesArray->At(i)) return kTRUE ;
 391
 392   return kFALSE;
 393
 394 }
 395
 396
 397 //____________________________________________________________________
 398 void AliCaloTrackMCReader::SetTrackChargeAndPID(const Int_t pdgCode, AliAODTrack *track) {
 399
 400   Float_t PID[10] = {0., 0., 0., 0., 0., 0., 0., 0., 0., 0.};
 401
 402   switch (pdgCode) {
 403
 404   case 22: // gamma
 405     track->SetCharge(0);
 406     PID[AliAODTrack::kUnknown] = 1.;
 407     track->SetPID(PID);
 408     break;
 409
 410   case 11: // e-
 411     track->SetCharge(-1);
 412     PID[AliAODTrack::kElectron] = 1.;
 413     track->SetPID(PID);
 414     break;
 415
 416   case -11: // e+
 417     track->SetCharge(+1);
 418     PID[AliAODTrack::kElectron] = 1.;
 419     track->SetPID(PID);
 420     break;
 421
 422   case 13: // mu-
 423     track->SetCharge(-1);
 424     PID[AliAODTrack::kMuon] = 1.;
 425     track->SetPID(PID);
 426     break;
 427
 428   case -13: // mu+
 429     track->SetCharge(+1);
 430     PID[AliAODTrack::kMuon] = 1.;
 431     track->SetPID(PID);
 432     break;
 433
 434   case 111: // pi0
 435     track->SetCharge(0);
 436     PID[AliAODTrack::kUnknown] = 1.;
 437     track->SetPID(PID);
 438     break;
 439
 440   case 211: // pi+
 441     track->SetCharge(+1);
 442     PID[AliAODTrack::kPion] = 1.;
 443     track->SetPID(PID);
 444     break;
 445
 446   case -211: // pi-
 447     track->SetCharge(-1);
 448     PID[AliAODTrack::kPion] = 1.;
 449     track->SetPID(PID);
 450     break;
 451
 452   case 130: // K0L
 453     track->SetCharge(0);
 454     PID[AliAODTrack::kUnknown] = 1.;
 455     track->SetPID(PID);
 456     break;
 457
 458   case 321: // K+
 459     track->SetCharge(+1);
 460     PID[AliAODTrack::kKaon] = 1.;
 461     track->SetPID(PID);
 462     break;
 463
 464   case -321: // K-
 465     track->SetCharge(-1);
 466     PID[AliAODTrack::kKaon] = 1.;
 467     track->SetPID(PID);
 468     break;
 469
 470   case 2112: // n
 471     track->SetCharge(0);
 472     PID[AliAODTrack::kUnknown] = 1.;
 473     track->SetPID(PID);
 474     break;
 475
 476   case 2212: // p
 477     track->SetCharge(+1);
 478     PID[AliAODTrack::kProton] = 1.;
 479     track->SetPID(PID);
 480     break;
 481
 482   case -2212: // anti-p
 483     track->SetCharge(-1);
 484     PID[AliAODTrack::kProton] = 1.;
 485     track->SetPID(PID);
 486     break;
 487
 488   case 310: // K0S
 489     track->SetCharge(0);
 490     PID[AliAODTrack::kUnknown] = 1.;
 491     track->SetPID(PID);
 492     break;
 493
 494   case 311: // K0
 495     track->SetCharge(0);
 496     PID[AliAODTrack::kUnknown] = 1.;
 497     track->SetPID(PID);
 498     break;
 499
 500   case -311: // anti-K0
 501     track->SetCharge(0);
 502     PID[AliAODTrack::kUnknown] = 1.;
 503     track->SetPID(PID);
 504     break;
 505
 506   case 221: // eta
 507     track->SetCharge(0);
 508     PID[AliAODTrack::kUnknown] = 1.;
 509     track->SetPID(PID);
 510     break;
 511
 512   case 3122: // lambda
 513     track->SetCharge(0);
 514     PID[AliAODTrack::kUnknown] = 1.;
 515     track->SetPID(PID);
 516     break;
 517
 518   case 3222: // Sigma+
 519     track->SetCharge(+1);
 520     PID[AliAODTrack::kUnknown] = 1.;
 521     track->SetPID(PID);
 522     break;
 523
 524   case 3212: // Sigma0
 525     track->SetCharge(-1);
 526     PID[AliAODTrack::kUnknown] = 1.;
 527     track->SetPID(PID);
 528     break;
 529
 530   case 3112: // Sigma-
 531     track->SetCharge(-1);
 532     PID[AliAODTrack::kUnknown] = 1.;
 533     track->SetPID(PID);
 534     break;
 535
 536   case 3322: // Xi0
 537     track->SetCharge(0);
 538     PID[AliAODTrack::kUnknown] = 1.;
 539     track->SetPID(PID);
 540     break;
 541
 542   case 3312: // Xi-
 543     track->SetCharge(-1);
 544     PID[AliAODTrack::kUnknown] = 1.;
 545     track->SetPID(PID);
 546     break;
 547
 548   case 3334: // Omega-
 549     track->SetCharge(-1);
 550     PID[AliAODTrack::kUnknown] = 1.;
 551     track->SetPID(PID);
 552     break;
 553
 554   case -2112: // n-bar
 555     track->SetCharge(0);
 556     PID[AliAODTrack::kUnknown] = 1.;
 557     track->SetPID(PID);
 558     break;
 559
 560   case -3122: // anti-Lambda
 561     track->SetCharge(0);
 562     PID[AliAODTrack::kUnknown] = 1.;
 563     track->SetPID(PID);
 564     break;
 565
 566   case -3222: // anti-Sigma-
 567     track->SetCharge(-1);
 568     PID[AliAODTrack::kUnknown] = 1.;
 569     track->SetPID(PID);
 570     break;
 571
 572   case -3212: // anti-Sigma0
 573     track->SetCharge(0);
 574     PID[AliAODTrack::kUnknown] = 1.;
 575     track->SetPID(PID);
 576     break;
 577
 578   case -3112: // anti-Sigma+
 579     track->SetCharge(+1);
 580     PID[AliAODTrack::kUnknown] = 1.;
 581     track->SetPID(PID);
 582     break;
 583
 584   case -3322: // anti-Xi0
 585     track->SetCharge(0);
 586     PID[AliAODTrack::kUnknown] = 1.;
 587     track->SetPID(PID);
 588     break;
 589
 590   case -3312: // anti-Xi+
 591     track->SetCharge(+1);
 592     break;
 593
 594   case -3334: // anti-Omega+
 595     track->SetCharge(+1);
 596     PID[AliAODTrack::kUnknown] = 1.;
 597     track->SetPID(PID);
 598     break;
 599
 600   case 411: // D+
 601     track->SetCharge(+1);
 602     PID[AliAODTrack::kUnknown] = 1.;
 603     track->SetPID(PID);
 604     break;
 605
 606   case -411: // D-
 607     track->SetCharge(-1);
 608     PID[AliAODTrack::kUnknown] = 1.;
 609     track->SetPID(PID);
 610     break;
 611
 612   case 421: // D0
 613     track->SetCharge(0);
 614     PID[AliAODTrack::kUnknown] = 1.;
 615     track->SetPID(PID);
 616     break;
 617
 618   case -421: // anti-D0
 619     track->SetCharge(0);
 620     PID[AliAODTrack::kUnknown] = 1.;
 621     track->SetPID(PID);
 622     break;
 623
 624   default : // unknown
 625     track->SetCharge(-99);
 626     PID[AliAODTrack::kUnknown] = 1.;
 627     track->SetPID(PID);
 628  }
 629
 630   return;
 631 }
 632
 633 //____________________________________________________________________
 634 void AliCaloTrackMCReader::SetCaloClusterPID(const Int_t pdgCode, AliAODCaloCluster *calo) {
 635
 636   Float_t PID[9] = {0., 0., 0., 0., 0., 0., 0., 0., 0.};
 637
 638   switch (pdgCode) {
 639
 640   case 22: // gamma
 641     PID[AliAODCaloCluster::kPhoton] = 1.;
 642     calo->SetPID(PID);
 643     break;
 644
 645   case 11: // e-
 646     PID[AliAODCaloCluster::kElectron] = 1.;
 647     calo->SetPID(PID);
 648     break;
 649
 650   case -11: // e+
 651     PID[AliAODCaloCluster::kElectron] = 1.;
 652     calo->SetPID(PID);
 653     break;
 654
 655   case 13: // mu-
 656     PID[AliAODCaloCluster::kCharged] = 1.;
 657     calo->SetPID(PID);
 658     break;
 659
 660   case -13: // mu+
 661     PID[AliAODCaloCluster::kCharged] = 1.;
 662     calo->SetPID(PID);
 663     break;
 664
 665   case 111: // pi0
 666     PID[AliAODCaloCluster::kPi0] = 1.;
 667     calo->SetPID(PID);
 668     break;
 669
 670   case 211: // pi+
 671     PID[AliAODCaloCluster::kCharged] = 1.;
 672     calo->SetPID(PID);
 673     break;
 674
 675   case -211: // pi-
 676     PID[AliAODCaloCluster::kCharged] = 1.;
 677     calo->SetPID(PID);
 678     break;
 679
 680   case 130: // K0L
 681     PID[AliAODCaloCluster::kKaon0] = 1.;
 682     PID[AliAODCaloCluster::kNeutral] = 1;
 683     calo->SetPID(PID);
 684     break;
 685
 686   case 321: // K+
 687     PID[AliAODCaloCluster::kCharged] = 1.;
 688     calo->SetPID(PID);
 689     break;
 690
 691   case -321: // K-
 692     PID[AliAODCaloCluster::kCharged] = 1.;
 693     calo->SetPID(PID);
 694     break;
 695
 696   case 2112: // n
 697     PID[AliAODCaloCluster::kNeutron] = 1.;
 698     PID[AliAODCaloCluster::kNeutral] = 1.;
 699     calo->SetPID(PID);
 700     break;
 701
 702   case 2212: // p
 703     PID[AliAODCaloCluster::kCharged] = 1.;
 704     calo->SetPID(PID);
 705     break;
 706
 707   case -2212: // anti-p
 708     PID[AliAODCaloCluster::kCharged] = 1.;
 709     calo->SetPID(PID);
 710     break;
 711
 712   case 310: // K0S
 713     PID[AliAODCaloCluster::kKaon0] = 1.;
 714     PID[AliAODCaloCluster::kNeutral] = 1.;
 715     calo->SetPID(PID);
 716     break;
 717
 718   case 311: // K0
 719     PID[AliAODCaloCluster::kKaon0] = 1.;
 720     PID[AliAODCaloCluster::kNeutral] = 1.;
 721     calo->SetPID(PID);
 722     break;
 723
 724   case -311: // anti-K0
 725     PID[AliAODCaloCluster::kKaon0] = 1.;
 726     PID[AliAODCaloCluster::kNeutral] = 1.;
 727     calo->SetPID(PID);
 728     break;
 729
 730   case 221: // eta
 731     PID[AliAODCaloCluster::kNeutral] = 1.;
 732     calo->SetPID(PID);
 733     break;
 734
 735   case 3122: // lambda
 736     PID[AliAODCaloCluster::kUnknown] = 1.;
 737     calo->SetPID(PID);
 738     break;
 739
 740   case 3222: // Sigma+
 741     PID[AliAODCaloCluster::kUnknown] = 1.;
 742     calo->SetPID(PID);
 743     break;
 744
 745   case 3212: // Sigma0
 746     PID[AliAODCaloCluster::kUnknown] = 1.;
 747     calo->SetPID(PID);
 748     break;
 749
 750   case 3112: // Sigma-
 751     PID[AliAODCaloCluster::kUnknown] = 1.;
 752     calo->SetPID(PID);
 753     break;
 754
 755   case 3322: // Xi0
 756     PID[AliAODCaloCluster::kUnknown] = 1.;
 757     calo->SetPID(PID);
 758     break;
 759
 760   case 3312: // Xi-
 761     PID[AliAODCaloCluster::kUnknown] = 1.;
 762     calo->SetPID(PID);
 763     break;
 764
 765   case 3334: // Omega-
 766     PID[AliAODCaloCluster::kUnknown] = 1.;
 767     calo->SetPID(PID);
 768     break;
 769
 770   case -2112: // n-bar
 771     PID[AliAODCaloCluster::kNeutron] = 1.;
 772     PID[AliAODCaloCluster::kNeutral] = 1.;
 773     calo->SetPID(PID);
 774     break;
 775
 776   case -3122: // anti-Lambda
 777     PID[AliAODCaloCluster::kUnknown] = 1.;
 778     calo->SetPID(PID);
 779     break;
 780
 781   case -3222: // anti-Sigma-
 782     PID[AliAODCaloCluster::kUnknown] = 1.;
 783     calo->SetPID(PID);
 784     break;
 785
 786   case -3212: // anti-Sigma0
 787     PID[AliAODCaloCluster::kUnknown] = 1.;
 788     calo->SetPID(PID);
 789     break;
 790
 791   case -3112: // anti-Sigma+
 792     PID[AliAODCaloCluster::kUnknown] = 1.;
 793     calo->SetPID(PID);
 794     break;
 795
 796   case -3322: // anti-Xi0
 797     PID[AliAODCaloCluster::kUnknown] = 1.;
 798     calo->SetPID(PID);
 799     break;
 800
 801   case -3312: // anti-Xi+
 802     PID[AliAODCaloCluster::kUnknown] = 1.;
 803     calo->SetPID(PID);
 804     break;
 805
 806   case -3334: // anti-Omega+
 807     PID[AliAODCaloCluster::kUnknown] = 1.;
 808     calo->SetPID(PID);
 809     break;
 810
 811   case 411: // D+
 812     PID[AliAODCaloCluster::kUnknown] = 1.;
 813     calo->SetPID(PID);
 814     break;
 815
 816   case -411: // D-
 817     PID[AliAODCaloCluster::kUnknown] = 1.;
 818     calo->SetPID(PID);
 819     break;
 820
 821   case 421: // D0
 822     PID[AliAODCaloCluster::kUnknown] = 1.;
 823     calo->SetPID(PID);
 824     break;
 825
 826   case -421: // anti-D0
 827     PID[AliAODCaloCluster::kUnknown] = 1.;
 828     calo->SetPID(PID);
 829     break;
 830
 831   default : // unknown
 832     PID[AliAODCaloCluster::kUnknown] = 1.;
 833     calo->SetPID(PID);
 834  }
 835
 836   return;
 837 }