]>
Commit | Line | Data |
---|---|---|
bdcfac30 | 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 * | |
4b707925 | 12 | * about the suitability of this software for any purpose. It is * |
bdcfac30 | 13 | * provided "as is" without express or implied warranty. * |
14 | **************************************************************************/ | |
15 | /* $Id$ */ | |
16 | ||
17 | /* History of cvs commits: | |
18 | * | |
19 | * $Log$ | |
2021085e | 20 | * Revision 1.4 2007/10/29 13:48:42 gustavo |
21 | * Corrected coding violations | |
22 | * | |
bdcfac30 | 23 | * |
24 | */ | |
25 | ||
26 | //_________________________________________________________________________ | |
27 | // Class for reading data (Kinematics) in order to do prompt gamma correlations | |
28 | // Class created from old AliPHOSGammaJet | |
29 | // (see AliRoot versions previous Release 4-09) | |
30 | // | |
31 | //*-- Author: Gustavo Conesa (LNF-INFN) | |
32 | ////////////////////////////////////////////////////////////////////////////// | |
33 | ||
34 | ||
35 | // --- ROOT system --- | |
36 | ||
37 | #include <TFile.h> | |
38 | #include <TParticle.h> | |
39 | #include <TH2.h> | |
40 | #include <TChain.h> | |
41 | #include <TRandom.h> | |
3bb2c538 | 42 | #include <TArrayI.h> |
43 | #include <TClonesArray.h> | |
bdcfac30 | 44 | |
45 | //---- ANALYSIS system ---- | |
46 | #include "AliGammaMCReader.h" | |
47 | #include "Riostream.h" | |
48 | #include "AliLog.h" | |
3bb2c538 | 49 | #include "AliStack.h" |
bdcfac30 | 50 | |
51 | ClassImp(AliGammaMCReader) | |
52 | ||
53 | //____________________________________________________________________________ | |
54 | AliGammaMCReader::AliGammaMCReader() : | |
3bb2c538 | 55 | AliGammaReader(), fDecayPi0(0), fCheckOverlapping(kFALSE), fNeutralParticlesArray(0x0) |
bdcfac30 | 56 | { |
57 | //Ctor | |
58 | ||
59 | //Initialize parameters | |
bdcfac30 | 60 | InitParameters(); |
4b707925 | 61 | fDataType = kMC; |
bdcfac30 | 62 | |
63 | } | |
64 | ||
65 | //____________________________________________________________________________ | |
66 | AliGammaMCReader::AliGammaMCReader(const AliGammaMCReader & g) : | |
3bb2c538 | 67 | AliGammaReader(g), fDecayPi0(g.fDecayPi0), fCheckOverlapping(g.fCheckOverlapping), |
68 | fNeutralParticlesArray(g.fNeutralParticlesArray?new TArrayI(*g.fNeutralParticlesArray):0x0) | |
bdcfac30 | 69 | { |
70 | // cpy ctor | |
71 | } | |
72 | ||
73 | //_________________________________________________________________________ | |
74 | AliGammaMCReader & AliGammaMCReader::operator = (const AliGammaMCReader & source) | |
75 | { | |
76 | // assignment operator | |
77 | ||
78 | if(&source == this) return *this; | |
79 | ||
4b707925 | 80 | fDecayPi0 = source.fDecayPi0; |
81 | fCheckOverlapping = source.fCheckOverlapping; | |
3bb2c538 | 82 | delete fNeutralParticlesArray; |
83 | fNeutralParticlesArray = source.fNeutralParticlesArray?new TArrayI(*source.fNeutralParticlesArray):0x0; | |
bdcfac30 | 84 | return *this; |
85 | ||
86 | } | |
87 | ||
3bb2c538 | 88 | //_________________________________ |
89 | AliGammaMCReader::~AliGammaMCReader() { | |
90 | //Dtor | |
91 | ||
92 | delete fNeutralParticlesArray ; | |
93 | ||
94 | } | |
95 | ||
4b707925 | 96 | //___________________________________________________________________________ |
97 | void AliGammaMCReader::CaseDecayGamma(Int_t index, TParticle * particle, AliStack * sta, | |
98 | TClonesArray * plEMCAL, Int_t &indexEMCAL, | |
99 | TClonesArray * plPHOS, Int_t &indexPHOS){ | |
100 | //In case pi0 are decayed by pythia, check if mother is pi0 and in such case look if | |
101 | //there is overlapping. Send particle=pi0 if there is overlapping. | |
3bb2c538 | 102 | |
4b707925 | 103 | TParticle * pmother =sta->Particle(particle->GetFirstMother()); |
104 | if(pmother->GetPdgCode() == 111 && pmother->GetNDaughters() == 2) {//Do not consider 3 particle decay case | |
105 | Int_t idaug0 = pmother->GetDaughter(0); | |
106 | Int_t idaug1 = pmother->GetDaughter(1); | |
107 | TParticle * pdaug0 = sta -> Particle(idaug0); | |
108 | TParticle * pdaug1 = sta -> Particle(idaug1); | |
3bb2c538 | 109 | |
4b707925 | 110 | if((index == idaug0 && pdaug0->Pt() > fNeutralPtCut) || |
3bb2c538 | 111 | (index == idaug1 && pdaug0->Pt() <= fNeutralPtCut)){//Check decay when first daughter arrives, do nothing with second. |
112 | ||
4b707925 | 113 | FillListWithDecayGammaOrPi0(pmother, pdaug0, pdaug1, plEMCAL, indexEMCAL, plPHOS, indexPHOS); |
114 | ||
3bb2c538 | 115 | } |
4b707925 | 116 | }//mother is a pi0 with 2 daughters |
117 | else{ | |
3bb2c538 | 118 | |
119 | if(IsInPHOS(particle->Phi(),particle->Eta())) { | |
120 | TParticle * pmother =sta->Particle(particle->GetFirstMother()); | |
121 | SetPhotonStatus(particle,pmother); | |
4b707925 | 122 | new((*plPHOS)[indexPHOS++]) TParticle(*particle) ; |
3bb2c538 | 123 | } |
124 | else if(IsInEMCAL(particle->Phi(),particle->Eta())){ | |
125 | TParticle * pmother =sta->Particle(particle->GetFirstMother()); | |
126 | SetPhotonStatus(particle,pmother); | |
4b707925 | 127 | new((*plEMCAL)[indexEMCAL++]) TParticle(*particle) ; |
3bb2c538 | 128 | } |
129 | } | |
4b707925 | 130 | } |
131 | ||
132 | //___________________________________________________________________________ | |
133 | void AliGammaMCReader::CaseGeantDecay(TParticle * particle, AliStack * stack, | |
134 | TClonesArray * plEMCAL, Int_t &indexEMCAL, | |
135 | TClonesArray * plPHOS, Int_t &indexPHOS){ | |
136 | //Find decay gamma from pi0, decayed by GEANT and put them in the list. | |
137 | ||
138 | Int_t ndaug = particle->GetNDaughters() ; | |
139 | TParticle * d1 = new TParticle(); | |
140 | TParticle * d2 = new TParticle(); | |
141 | ||
142 | if(ndaug > 0){//At least 1 daugther | |
143 | d1 = stack->Particle(particle->GetDaughter(0)); | |
144 | if (ndaug > 1 ) | |
145 | d2 = stack->Particle(particle->GetDaughter(1)); | |
146 | ||
147 | FillListWithDecayGammaOrPi0(particle, d1, d2, plEMCAL, indexEMCAL, plPHOS, indexPHOS); | |
148 | ||
149 | }// ndaugh > 0 | |
150 | } | |
151 | ||
152 | //___________________________________________________________________________ | |
153 | void AliGammaMCReader::CasePi0Decay(TParticle * pPi0, | |
154 | TClonesArray * plEMCAL, Int_t &indexEMCAL, | |
155 | TClonesArray * plPHOS, Int_t &indexPHOS){ | |
156 | ||
157 | //Decays pi0, see if aperture angle is small and then add the pi0 or the 2 gamma | |
158 | ||
159 | TLorentzVector lvPi0, lvGamma1, lvGamma2 ; | |
160 | //Double_t angle = 0; | |
161 | ||
162 | lvPi0.SetPxPyPzE(pPi0->Px(),pPi0->Py(),pPi0->Pz(),pPi0->Energy()); | |
163 | ||
164 | //Decay | |
165 | MakePi0Decay(lvPi0,lvGamma1,lvGamma2);//,angle); | |
166 | ||
167 | //Check if Pi0 is in the acceptance of the calorimeters, if aperture angle is small, keep it | |
168 | TParticle * pPhoton1 = new TParticle(22,1,0,0,0,0,lvGamma1.Px(),lvGamma1.Py(), | |
169 | lvGamma1.Pz(),lvGamma1.E(),0,0,0,0); | |
170 | TParticle * pPhoton2 = new TParticle(22,1,0,0,0,0,lvGamma2.Px(),lvGamma2.Py(), | |
171 | lvGamma2.Pz(),lvGamma2.E(),0,0,0,0); | |
172 | ||
173 | FillListWithDecayGammaOrPi0(pPi0, pPhoton1, pPhoton2, plEMCAL, indexEMCAL, plPHOS, indexPHOS); | |
174 | ||
175 | } | |
176 | ||
177 | //_______________________________________________________________ | |
178 | void AliGammaMCReader::InitParameters() | |
179 | { | |
180 | ||
181 | //Initialize the parameters of the analysis. | |
182 | ||
183 | fDecayPi0 = kGeantDecay; | |
184 | fCheckOverlapping = kTRUE ; | |
3bb2c538 | 185 | |
186 | Int_t pdgarray[]={12,14,16};// skip neutrinos | |
187 | fNeutralParticlesArray = new TArrayI(3, pdgarray); | |
188 | ||
4b707925 | 189 | } |
bdcfac30 | 190 | |
191 | //____________________________________________________________________________ | |
192 | void AliGammaMCReader::CreateParticleList(TObject * data, TObject *, | |
4b707925 | 193 | TClonesArray * plCh, |
194 | TClonesArray * plEMCAL, | |
bdcfac30 | 195 | TClonesArray * plPHOS, |
4b707925 | 196 | TClonesArray *plParton,TClonesArray *,TClonesArray *) |
bdcfac30 | 197 | { |
198 | //Create list of particles from EMCAL, PHOS and CTS. | |
199 | AliStack * stack = (AliStack *) data ; | |
200 | Int_t indexCh = plCh->GetEntries() ; | |
201 | Int_t indexEMCAL = plEMCAL->GetEntries() ; | |
202 | Int_t indexPHOS = plPHOS->GetEntries() ; | |
203 | Int_t indexParton = plParton->GetEntries() ; | |
204 | Int_t iParticle = 0 ; | |
205 | Double_t charge = 0.; | |
206 | ||
207 | for (iParticle=0 ; iParticle < stack->GetNprimary() ; iParticle++) { | |
208 | TParticle * particle = stack->Particle(iParticle); | |
4b707925 | 209 | |
bdcfac30 | 210 | //Keep partons |
211 | if(particle->GetStatusCode() == 21 && iParticle>=2){//All partons, not nucleus | |
212 | new((*plParton)[indexParton++]) TParticle(*particle) ; | |
213 | } | |
214 | ||
215 | //Keep Stable particles | |
4b707925 | 216 | if((particle->GetStatusCode() == 1) && (particle->Pt() > 0)){ |
bdcfac30 | 217 | |
218 | charge = TDatabasePDG::Instance()->GetParticle(particle->GetPdgCode())->Charge(); | |
219 | ||
220 | //---------- Charged particles ---------------------- | |
3bb2c538 | 221 | if( fSwitchOnCTS && (charge != 0) && (particle->Pt() > fChargedPtCut)){ |
bdcfac30 | 222 | //Particles in CTS acceptance |
223 | if(TMath::Abs(particle->Eta())<fCTSEtaCut){ | |
224 | //Fill lists | |
225 | new((*plCh)[indexCh++]) TParticle(*particle) ; | |
226 | } | |
227 | } | |
3bb2c538 | 228 | |
bdcfac30 | 229 | //-------------Neutral particles ---------------------- |
3bb2c538 | 230 | else if((charge == 0) && particle->Pt() > fNeutralPtCut ){ //&& |
231 | //TMath::Abs(particle->GetPdgCode())>16){//Avoid neutrinos | |
bdcfac30 | 232 | |
233 | if(particle->GetPdgCode()!=111){ | |
4b707925 | 234 | //In case that we force PYTHIA to decay pi0, and we want to check the overlapping of |
235 | // the decay gamma. | |
236 | if(particle->GetPdgCode() == 22 && fDecayPi0==kDecayGamma){ | |
237 | CaseDecayGamma(iParticle,particle, stack,plEMCAL, indexEMCAL, plPHOS, indexPHOS); //If pythia decays pi0 | |
238 | } | |
239 | else{ | |
3bb2c538 | 240 | //Take out some particles like neutrinos |
2021085e | 241 | if(!SkipNeutralParticles(particle->GetPdgCode()) && !particle->IsPrimary()){ // protection added (MG) |
3bb2c538 | 242 | TParticle * pmother =stack->Particle(particle->GetFirstMother()); |
243 | if(IsInPHOS(particle->Phi(),particle->Eta())){ | |
244 | if(particle->GetPdgCode()==22)SetPhotonStatus(particle,pmother); | |
245 | new((*plPHOS)[indexPHOS++]) TParticle(*particle) ; | |
246 | } | |
247 | else if(IsInEMCAL(particle->Phi(),particle->Eta())){ | |
248 | if(particle->GetPdgCode()==22) SetPhotonStatus(particle,pmother); | |
249 | new((*plEMCAL)[indexEMCAL++]) TParticle(*particle) ; | |
250 | } | |
251 | }//skip neutral particles | |
252 | }//Case kDecayGamma | |
bdcfac30 | 253 | }//no pi0 |
254 | else{ | |
255 | if(fDecayPi0 == kNoDecay){//keep the pi0 do not decay | |
256 | if(IsInPHOS(particle->Phi(),particle->Eta())) | |
257 | new((*plPHOS)[indexPHOS++]) TParticle(*particle) ; | |
258 | else if(IsInEMCAL(particle->Phi(),particle->Eta())) | |
259 | new((*plEMCAL)[indexEMCAL++]) TParticle(*particle) ; | |
3bb2c538 | 260 | }//nodecay |
bdcfac30 | 261 | else if(fDecayPi0 == kDecay) |
4b707925 | 262 | CasePi0Decay(particle,plEMCAL,indexEMCAL,plPHOS, indexPHOS); |
bdcfac30 | 263 | else if(fDecayPi0 == kGeantDecay) |
4b707925 | 264 | CaseGeantDecay(particle, stack,plEMCAL, indexEMCAL, plPHOS, indexPHOS); |
bdcfac30 | 265 | }//pi0 |
266 | }//neutral particle | |
267 | }//stable particle | |
268 | }//particle loop | |
269 | } | |
270 | ||
4b707925 | 271 | |
272 | //___________________________________________________________________________ | |
273 | void AliGammaMCReader::FillListWithDecayGammaOrPi0(TParticle * pPi0, TParticle * pdaug0, TParticle * pdaug1, | |
274 | TClonesArray * plEMCAL, Int_t &indexEMCAL, | |
275 | TClonesArray * plPHOS, Int_t &indexPHOS){ | |
276 | ||
277 | //Check if decay gamma overlapp in calorimeter, in such case keep the pi0, if not keep both photons. | |
4b707925 | 278 | Bool_t overlap = kFALSE ; |
279 | TLorentzVector lv1 , lv2 ; | |
280 | pdaug0->Momentum(lv1); | |
281 | pdaug1->Momentum(lv2); | |
282 | Double_t angle = lv1.Angle(lv2.Vect()); | |
3bb2c538 | 283 | |
284 | if(IsInEMCAL(pPi0->Phi(), pPi0->Eta())) | |
285 | { | |
4b707925 | 286 | if (angle < fEMCALMinAngle){ |
287 | new((*plEMCAL)[indexEMCAL++]) TParticle(*pPi0) ; | |
288 | overlap = kTRUE; | |
289 | } | |
290 | } | |
3bb2c538 | 291 | |
292 | else if(IsInPHOS(pPi0->Phi(), pPi0->Eta())){ | |
293 | if (angle < fPHOSMinAngle){ | |
294 | new((*plPHOS)[indexPHOS++]) TParticle(*pPi0) ; | |
295 | overlap = kTRUE; | |
4b707925 | 296 | } |
297 | }//fCheckOverlapping | |
3bb2c538 | 298 | |
4b707925 | 299 | //Fill with gammas if not overlapp |
300 | if(!overlap){ | |
301 | if(pdaug0->GetPdgCode() == 22 || TMath::Abs(pdaug0->GetPdgCode() ) == 11 ){ | |
3bb2c538 | 302 | pdaug0->SetStatusCode(kPi0DecayPhoton); |
4b707925 | 303 | if(IsInEMCAL(pdaug0->Phi(),pdaug0->Eta()) && pdaug0->Pt() > fNeutralPtCut) |
304 | new((*plEMCAL)[indexEMCAL++]) TParticle(*pdaug0) ; | |
305 | else if(IsInPHOS(pdaug0->Phi(),pdaug0->Eta()) && pdaug0->Pt() > fNeutralPtCut) | |
306 | new((*plPHOS)[indexPHOS++]) TParticle(*pdaug0) ; | |
307 | } | |
308 | ||
309 | if(pdaug1->GetPdgCode() == 22 || TMath::Abs(pdaug1->GetPdgCode() ) == 11 ){ | |
3bb2c538 | 310 | pdaug1->SetStatusCode(kPi0DecayPhoton); |
4b707925 | 311 | if(IsInEMCAL(pdaug1->Phi(),pdaug1->Eta()) && pdaug1->Pt() > fNeutralPtCut) |
312 | new((*plEMCAL)[indexEMCAL++]) TParticle(*pdaug1) ; | |
313 | else if(IsInPHOS(pdaug1->Phi(),pdaug1->Eta()) && pdaug1->Pt() > fNeutralPtCut) | |
314 | new((*plPHOS)[indexPHOS++]) TParticle(*pdaug1) ; | |
315 | } | |
316 | }// overlap? | |
3bb2c538 | 317 | |
4b707925 | 318 | } |
319 | ||
320 | ||
321 | ||
bdcfac30 | 322 | //___________________________________________________________________________ |
3bb2c538 | 323 | Bool_t AliGammaMCReader::IsInEMCAL(Double_t phi, Double_t eta) const { |
bdcfac30 | 324 | //Check if particle is in EMCAL acceptance |
3bb2c538 | 325 | |
326 | if(fSwitchOnEMCAL){ | |
327 | if(phi<0) | |
328 | phi+=TMath::TwoPi(); | |
329 | if( phi > fPhiEMCALCut[0] && phi < fPhiEMCALCut[1] && | |
330 | TMath::Abs(eta)<fEMCALEtaCut) return kTRUE ; | |
331 | else return kFALSE; | |
332 | } | |
bdcfac30 | 333 | return kFALSE ; |
334 | } | |
335 | ||
336 | //___________________________________________________________________________ | |
3bb2c538 | 337 | Bool_t AliGammaMCReader::IsInPHOS(Double_t phi, Double_t eta) const { |
338 | //Check if particle is in PHOS acceptance | |
339 | ||
340 | if(fSwitchOnPHOS){ | |
341 | if(phi<0) | |
342 | phi+=TMath::TwoPi(); | |
343 | if( phi > fPhiPHOSCut[0] && phi < fPhiPHOSCut[1] && | |
344 | TMath::Abs(eta)<fPHOSEtaCut) return kTRUE ; | |
345 | else return kFALSE; | |
346 | } | |
bdcfac30 | 347 | |
348 | return kFALSE ; | |
349 | } | |
350 | ||
bdcfac30 | 351 | //____________________________________________________________________________ |
4b707925 | 352 | void AliGammaMCReader::MakePi0Decay(TLorentzVector &p0, TLorentzVector &p1, |
353 | TLorentzVector &p2){//, Double_t &angle) { | |
bdcfac30 | 354 | // Perform isotropic decay pi0 -> 2 photons |
355 | // p0 is pi0 4-momentum (inut) | |
356 | // p1 and p2 are photon 4-momenta (output) | |
357 | // cout<<"Boost vector"<<endl; | |
358 | Double_t mPi0 = TDatabasePDG::Instance()->GetParticle(111)->Mass(); | |
359 | TVector3 b = p0.BoostVector(); | |
360 | //cout<<"Parameters"<<endl; | |
361 | //Double_t mPi0 = p0.M(); | |
362 | Double_t phi = TMath::TwoPi() * gRandom->Rndm(); | |
363 | Double_t cosThe = 2 * gRandom->Rndm() - 1; | |
364 | Double_t cosPhi = TMath::Cos(phi); | |
365 | Double_t sinPhi = TMath::Sin(phi); | |
366 | Double_t sinThe = TMath::Sqrt(1-cosThe*cosThe); | |
367 | Double_t ePi0 = mPi0/2.; | |
368 | //cout<<"ePi0 "<<ePi0<<endl; | |
369 | //cout<<"Components"<<endl; | |
370 | p1.SetPx(+ePi0*cosPhi*sinThe); | |
371 | p1.SetPy(+ePi0*sinPhi*sinThe); | |
372 | p1.SetPz(+ePi0*cosThe); | |
373 | p1.SetE(ePi0); | |
374 | //cout<<"p1: "<<p1.Px()<<" "<<p1.Py()<<" "<<p1.Pz()<<" "<<p1.E()<<endl; | |
375 | //cout<<"p1 Mass: "<<p1.Px()*p1.Px()+p1.Py()*p1.Py()+p1.Pz()*p1.Pz()-p1.E()*p1.E()<<endl; | |
376 | p2.SetPx(-ePi0*cosPhi*sinThe); | |
377 | p2.SetPy(-ePi0*sinPhi*sinThe); | |
378 | p2.SetPz(-ePi0*cosThe); | |
379 | p2.SetE(ePi0); | |
380 | //cout<<"p2: "<<p2.Px()<<" "<<p2.Py()<<" "<<p2.Pz()<<" "<<p2.E()<<endl; | |
381 | //cout<<"p2 Mass: "<<p2.Px()*p2.Px()+p2.Py()*p2.Py()+p2.Pz()*p2.Pz()-p2.E()*p2.E()<<endl; | |
382 | //cout<<"Boost "<<b.X()<<" "<<b.Y()<<" "<<b.Z()<<endl; | |
383 | p1.Boost(b); | |
384 | //cout<<"p1: "<<p1.Px()<<" "<<p1.Py()<<" "<<p1.Pz()<<" "<<p1.E()<<endl; | |
385 | p2.Boost(b); | |
386 | //cout<<"p2: "<<p2.Px()<<" "<<p2.Py()<<" "<<p2.Pz()<<" "<<p2.E()<<endl; | |
387 | //cout<<"angle"<<endl; | |
4b707925 | 388 | //angle = p1.Angle(p2.Vect()); |
bdcfac30 | 389 | //cout<<angle<<endl; |
390 | } | |
391 | ||
392 | //________________________________________________________________ | |
393 | void AliGammaMCReader::Print(const Option_t * opt) const | |
394 | { | |
395 | ||
396 | //Print some relevant parameters set for the analysis | |
397 | if(! opt) | |
398 | return; | |
399 | ||
400 | Info("Print", "%s %s", GetName(), GetTitle() ) ; | |
401 | ||
bdcfac30 | 402 | printf("Decay Pi0? : %d\n", fDecayPi0) ; |
4b707925 | 403 | printf("Check Overlapping? : %d\n", fCheckOverlapping) ; |
404 | ||
bdcfac30 | 405 | } |
4b707925 | 406 | |
407 | ||
408 | ||
3bb2c538 | 409 | //________________________________________________________________ |
410 | void AliGammaMCReader::SetPhotonStatus(TParticle* pphoton, TParticle* pmother) | |
411 | { | |
4b707925 | 412 | |
3bb2c538 | 413 | //Check the origin of the photon and tag it as decay from pi0, from eta, from other mesons, or prompt or fragmentation. |
414 | ||
415 | if(pmother->GetStatusCode() != 21 ){ | |
416 | if(pmother->GetStatusCode() ==11){ | |
417 | if(pmother->GetPdgCode()==111) pphoton->SetStatusCode(kPi0DecayPhoton);//decay gamma from pi0 | |
418 | if(pmother->GetPdgCode()==221) pphoton->SetStatusCode(kEtaDecayPhoton);//decay gamma from eta | |
419 | else pphoton->SetStatusCode(kOtherDecayPhoton);// decay gamma from other pphotons | |
420 | } | |
421 | else pphoton->SetStatusCode(kUnknown); | |
422 | } | |
423 | else if(pmother->GetPdgCode()==22) pphoton->SetStatusCode(kPromptPhoton);//Prompt photon | |
424 | else pphoton->SetStatusCode(kFragmentPhoton); //Fragmentation photon | |
425 | ||
426 | } | |
427 | ||
428 | //________________________________________________________________ | |
429 | Bool_t AliGammaMCReader::SkipNeutralParticles(Int_t pdg) const { | |
430 | //Check if pdg is equal to one of the neutral particles list | |
431 | //These particles will be skipped from analysis. | |
432 | ||
433 | for(Int_t i= 0; i < fNeutralParticlesArray->GetSize(); i++) | |
434 | if(TMath::Abs(pdg) == fNeutralParticlesArray->At(i)) return kTRUE ; | |
435 | ||
436 | return kFALSE; | |
437 | ||
438 | } |