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