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