1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
16 //_________________________________________________________________________
17 // Class for analysis utils for MC data
18 // stored in stack or event header.
20 // - method to check the origin of a given track/cluster
21 // - method to obtain the generated jets
23 //*-- Author: Gustavo Conesa (LNF-INFN)
24 //////////////////////////////////////////////////////////////////////////////
27 // --- ROOT system ---
30 #include "TParticle.h"
31 #include "TDatabasePDG.h"
34 //---- ANALYSIS system ----
35 #include "AliMCAnalysisUtils.h"
36 #include "AliCaloTrackReader.h"
38 #include "AliGenPythiaEventHeader.h"
39 #include "AliAODMCParticle.h"
41 ClassImp(AliMCAnalysisUtils)
43 //________________________________________
44 AliMCAnalysisUtils::AliMCAnalysisUtils() :
49 fMCGenerator("PYTHIA")
54 //_______________________________________
55 AliMCAnalysisUtils::~AliMCAnalysisUtils()
57 // Remove all pointers.
65 //_____________________________________________________________________________________________
66 Int_t AliMCAnalysisUtils::CheckCommonAncestor(Int_t index1, Int_t index2,
67 const AliCaloTrackReader* reader,
68 Int_t & ancPDG, Int_t & ancStatus,
69 TLorentzVector & momentum, TVector3 & prodVertex)
71 //Check the first common ancestor of 2 clusters, given the most likely labels of the primaries generating such clusters.
79 if(label1[0]==label2[0])
81 //printf("AliMCAnalysisUtils::CheckCommonAncestor() - Already the same label: %d\n",label1[0]);
87 if(reader->ReadAODMCParticles())
89 TClonesArray * mcparticles = reader->GetAODMCParticles();
91 Int_t label=label1[0];
92 if(label < 0) return -1;
94 while(label > -1 && counter1 < 99)
97 AliAODMCParticle * mom = (AliAODMCParticle *) mcparticles->At(label);
100 label = mom->GetMother() ;
101 label1[counter1]=label;
103 //printf("\t counter %d, label %d\n", counter1,label);
106 //printf("Org label2=%d,\n",label2[0]);
108 if(label < 0) return -1;
110 while(label > -1 && counter2 < 99)
113 AliAODMCParticle * mom = (AliAODMCParticle *) mcparticles->At(label);
116 label = mom->GetMother() ;
117 label2[counter2]=label;
119 //printf("\t counter %d, label %d\n", counter2,label);
123 { //Kine stack from ESDs
124 AliStack * stack = reader->GetStack();
126 Int_t label=label1[0];
127 while(label > -1 && counter1 < 99)
130 TParticle * mom = stack->Particle(label);
133 label = mom->GetFirstMother() ;
134 label1[counter1]=label;
136 //printf("\t counter %d, label %d\n", counter1,label);
139 //printf("Org label2=%d,\n",label2[0]);
142 while(label > -1 && counter2 < 99)
145 TParticle * mom = stack->Particle(label);
148 label = mom->GetFirstMother() ;
149 label2[counter2]=label;
151 //printf("\t counter %d, label %d\n", counter2,label);
153 }// Kine stack from ESDs
154 }//First labels not the same
156 if((counter1==99 || counter2==99) && fDebug >=0)
157 printf("AliMCAnalysisUtils::CheckCommonAncestor() - Genealogy too large c1: %d, c2= %d\n", counter1, counter2);
158 //printf("CheckAncestor:\n");
160 Int_t commonparents = 0;
162 //printf("counters %d %d \n",counter1, counter2);
163 for (Int_t c1 = 0; c1 < counter1; c1++)
165 for (Int_t c2 = 0; c2 < counter2; c2++)
167 if(label1[c1]==label2[c2] && label1[c1]>-1)
169 ancLabel = label1[c1];
172 if(reader->ReadAODMCParticles())
174 AliAODMCParticle * mom = (AliAODMCParticle *) reader->GetAODMCParticles()->At(label1[c1]);
178 ancPDG = mom->GetPdgCode();
179 ancStatus = mom->GetStatus();
180 momentum.SetPxPyPzE(mom->Px(),mom->Py(),mom->Pz(),mom->E());
181 prodVertex.SetXYZ(mom->Xv(),mom->Yv(),mom->Zv());
186 TParticle * mom = (reader->GetStack())->Particle(label1[c1]);
190 ancPDG = mom->GetPdgCode();
191 ancStatus = mom->GetStatusCode();
192 mom->Momentum(momentum);
193 prodVertex.SetXYZ(mom->Vx(),mom->Vy(),mom->Vz());
197 //First ancestor found, end the loops
201 }//second cluster loop
202 }//first cluster loop
208 momentum.SetXYZT(0,0,0,0);
209 prodVertex.SetXYZ(-10,-10,-10);
215 //________________________________________________________________________________________
216 Int_t AliMCAnalysisUtils::CheckOrigin(const Int_t * label, Int_t nlabels,
217 const AliCaloTrackReader* reader, TString calorimeter)
219 //Play with the montecarlo particles if available
223 printf("AliMCAnalysisUtils::CheckOrigin(nlabel<=0) - No MC labels available, please check!!!\n");
227 TObjArray* arrayCluster = 0;
228 if ( calorimeter == "EMCAL" ) arrayCluster = reader->GetEMCALClusters();
229 else if ( calorimeter == "PHOS" ) arrayCluster= reader->GetPHOSClusters();
231 //Select where the information is, ESD-galice stack or AOD mcparticles branch
232 if(reader->ReadStack()){
233 tag = CheckOriginInStack(label, nlabels, reader->GetStack(), arrayCluster);
235 else if(reader->ReadAODMCParticles()){
236 tag = CheckOriginInAOD(label, nlabels, reader->GetAODMCParticles(),arrayCluster);
242 //____________________________________________________________________________________________________
243 Int_t AliMCAnalysisUtils::CheckOrigin(Int_t label, const AliCaloTrackReader* reader, TString calorimeter)
245 //Play with the montecarlo particles if available
249 printf("AliMCAnalysisUtils::CheckOrigin(label<0) - No MC labels available, please check!!!\n");
253 TObjArray* arrayCluster = 0;
254 if ( calorimeter == "EMCAL" ) arrayCluster = reader->GetEMCALClusters();
255 else if ( calorimeter == "PHOS" ) arrayCluster= reader->GetPHOSClusters();
257 Int_t labels[]={label};
259 //Select where the information is, ESD-galice stack or AOD mcparticles branch
260 if(reader->ReadStack()){
261 tag = CheckOriginInStack(labels, 1,reader->GetStack(),arrayCluster);
263 else if(reader->ReadAODMCParticles()){
264 tag = CheckOriginInAOD(labels, 1,reader->GetAODMCParticles(),arrayCluster);
270 //__________________________________________________________________________________________
271 Int_t AliMCAnalysisUtils::CheckOriginInStack(const Int_t *labels, Int_t nlabels,
272 AliStack* stack, const TObjArray* arrayCluster)
274 // Play with the MC stack if available. Tag particles depending on their origin.
275 // Do same things as in CheckOriginInAOD but different input.
277 //generally speaking, label is the MC label of a reconstructed
278 //entity (track, cluster, etc) for which we want to know something
279 //about its heritage, but one can also use it directly with stack
280 //particles not connected to reconstructed entities
285 printf("AliMCAnalysisUtils::CheckOriginInStack() - Stack is not available, check analysis settings in configuration file, STOP!!\n");
290 Int_t label=labels[0];//Most significant particle contributing to the cluster
292 if(label >= 0 && label < stack->GetNtrack())
294 //MC particle of interest is the "mom" of the entity
295 TParticle * mom = stack->Particle(label);
297 Int_t mPdgSign = mom->GetPdgCode();
298 Int_t mPdg = TMath::Abs(mPdgSign);
299 Int_t mStatus = mom->GetStatusCode() ;
300 Int_t iParent = mom->GetFirstMother() ;
302 if(fDebug > 0 && label < 8 && fMCGenerator!="") printf("AliMCAnalysisUtils::CheckOriginInStack() - Mother is parton %d\n",iParent);
304 //GrandParent of the entity
305 TParticle * parent = NULL;
310 parent = stack->Particle(iParent);
314 pPdg = TMath::Abs(parent->GetPdgCode());
315 pStatus = parent->GetStatusCode();
318 else if(fDebug > 0 ) printf("AliMCAnalysisUtils::CheckOriginInStack() - Parent with label %d\n",iParent);
322 printf("AliMCAnalysisUtils::CheckOriginInStack() - Cluster most contributing mother and its parent: \n");
323 printf("\t Mother label %d, pdg %d, status %d\n",iMom, mPdg, mStatus);
324 printf("\t Parent label %d, pdg %d, status %d\n",iParent, pPdg, pStatus);
327 //Check if "mother" of entity is converted, if not, get the first non converted mother
328 if((mPdg == 22 || mPdg == 11) && (pPdg == 22 || pPdg == 11) && mStatus == 0)
330 SetTagBit(tag,kMCConversion);
332 //Check if the mother is photon or electron with status not stable
333 while ((pPdg == 22 || pPdg == 11) && mStatus != 1)
336 iMom = mom->GetFirstMother();
337 mom = stack->Particle(iMom);
338 mPdgSign = mom->GetPdgCode();
339 mPdg = TMath::Abs(mPdgSign);
340 mStatus = mom->GetStatusCode() ;
341 iParent = mom->GetFirstMother() ;
343 if(fDebug > 0 && label < 8 ) printf("AliMCAnalysisUtils::CheckOriginInStack() - Mother is parton %d\n",iParent);
348 parent = stack->Particle(iParent);
351 pPdg = TMath::Abs(parent->GetPdgCode());
352 pStatus = parent->GetStatusCode();
355 else {// in case of gun/box simulations
364 printf("AliMCAnalysisUtils::CheckOriginInStack() - Converted photon/electron: \n");
365 printf("\t Mother label %d, pdg %d, status %d\n",iMom, mPdg, mStatus);
366 printf("\t Parent label %d, pdg %d, status %d\n",iParent, pPdg, pStatus);
369 }//mother and parent are electron or photon and have status 0
370 else if((mPdg == 22 || mPdg == 11) && mStatus == 0)
372 //Still a conversion but only one electron/photon generated. Just from hadrons but not decays.
373 if(pPdg == 2112 || pPdg == 211 || pPdg == 321 ||
374 pPdg == 2212 || pPdg == 130 || pPdg == 13 )
376 SetTagBit(tag,kMCConversion);
377 iMom = mom->GetFirstMother();
378 mom = stack->Particle(iMom);
379 mPdgSign = mom->GetPdgCode();
380 mPdg = TMath::Abs(mPdgSign);
384 printf("AliMCAnalysisUtils::CheckOriginInStack() - Converted hadron: \n");
385 printf("\t Mother label %d, pdg %d, status %d\n",iMom, mPdg, mStatus);
389 //Comment for the next lines, we do not check the parent of the hadron for the moment.
390 //iParent = mom->GetFirstMother() ;
391 //if(fDebug > 0 && label < 8 ) printf("AliMCAnalysisUtils::CheckOriginInStack() - Mother is parton %d\n",iParent);
395 // parent = stack->Particle(iParent);
396 // pPdg = TMath::Abs(parent->GetPdgCode());
399 // conversion into electrons/photons checked
401 //first check for typical charged particles
402 if (mPdg == 13) SetTagBit(tag,kMCMuon);
403 else if(mPdg == 211) SetTagBit(tag,kMCPion);
404 else if(mPdg == 321) SetTagBit(tag,kMCKaon);
405 else if(mPdgSign == 2212) SetTagBit(tag,kMCProton);
406 else if(mPdgSign == -2212) SetTagBit(tag,kMCAntiProton);
407 else if(mPdgSign == 2112) SetTagBit(tag,kMCNeutron);
408 else if(mPdgSign == -2112) SetTagBit(tag,kMCAntiNeutron);
410 //check for pi0 and eta (shouldn't happen unless their decays were turned off)
414 SetTagBit(tag,kMCPi0Decay);
417 printf("AliMCAnalysisUtils::CheckOriginInStack() - First mother is directly pi0, not decayed by generator \n");
419 CheckOverlapped2GammaDecay(labels,nlabels, iMom, stack, tag); //set to kMCPi0 if 2 gammas in same cluster
424 SetTagBit(tag,kMCEtaDecay);
427 printf("AliMCAnalysisUtils::CheckOriginInStack() - First mother is directly eta, not decayed by generator \n");
429 CheckOverlapped2GammaDecay(labels,nlabels, iMom, stack, tag); //set to kMCEta if 2 gammas in same cluster
434 SetTagBit(tag,kMCPhoton);
438 SetTagBit(tag,kMCPi0Decay);
440 if(fDebug > 2 ) printf("AliMCAnalysisUtils::CheckOriginInStack() - PYTHIA pi0 decay photon, parent pi0 with status 11 \n");
442 CheckOverlapped2GammaDecay(labels,nlabels, iParent, stack, tag); //set to kMCPi0 if 2 gammas in same cluster
443 // In case it did not merge, check if the decay companion is lost
444 if(!CheckTagBit(tag, kMCPi0) && !CheckTagBit(tag,kMCDecayPairInCalo))
445 CheckLostDecayPair(arrayCluster,iMom, iParent, stack, tag);
447 //printf("Bit set is Merged %d, Pair in calo %d, Lost %d\n",CheckTagBit(tag, kMCPi0),CheckTagBit(tag,kMCDecayPairInCalo),CheckTagBit(tag,kMCDecayPairLost));
449 else if (pPdg == 221)
451 SetTagBit(tag, kMCEtaDecay);
453 if(fDebug > 2 ) printf("AliMCAnalysisUtils::CheckOriginInStack() - PYTHIA eta decay photon, parent pi0 with status 11 \n");
455 CheckOverlapped2GammaDecay(labels,nlabels, iParent, stack, tag);//set to kMCEta if 2 gammas in same cluster
456 // In case it did not merge, check if the decay companion is lost
457 if(!CheckTagBit(tag, kMCEta) && !CheckTagBit(tag,kMCDecayPairInCalo))
458 CheckLostDecayPair(arrayCluster,iMom, iParent, stack, tag);
460 else if(mStatus == 1)
461 { //undecayed particle
462 if(fMCGenerator == "PYTHIA")
464 if(iParent < 8 && iParent > 5)
466 if(pPdg == 22) SetTagBit(tag,kMCPrompt);
467 else SetTagBit(tag,kMCFragmentation);
469 else if(iParent <= 5)
471 SetTagBit(tag, kMCISR); //Initial state radiation
473 else SetTagBit(tag,kMCUnknown);
476 else if(fMCGenerator == "HERWIG")
484 if(parent->GetFirstMother()<=5) break;
485 iParent = parent->GetFirstMother();
486 parent=stack->Particle(iParent);
487 pStatus= parent->GetStatusCode();
488 pPdg = TMath::Abs(parent->GetPdgCode());
490 }//Look for the parton
492 if(iParent < 8 && iParent > 5)
494 if(pPdg == 22) SetTagBit(tag,kMCPrompt);
495 else SetTagBit(tag,kMCFragmentation);
497 else SetTagBit(tag,kMCISR);//Initial state radiation
499 else SetTagBit(tag,kMCUnknown);
502 else SetTagBit(tag,kMCOtherDecay);
506 //Electron check. Where did that electron come from?
507 else if(mPdg == 11){ //electron
508 if(pPdg == 11 && parent)
510 Int_t iGrandma = parent->GetFirstMother();
513 TParticle* gma = (TParticle*)stack->Particle(iGrandma); //get mother
514 Int_t gPdg = TMath::Abs(gma->GetPdgCode());
516 if (gPdg == 23) { SetTagBit(tag,kMCZDecay); } //parent is Z-boson
517 else if (gPdg == 24) { SetTagBit(tag,kMCWDecay); } //parent is W-boson
521 SetTagBit(tag,kMCElectron);
523 if(fDebug > 0) printf("AliMCAnalysisUtils::CheckOriginInStack() - Checking ancestors of electrons\n");
525 if (pPdg == 111) { SetTagBit(tag,kMCPi0Decay); SetTagBit(tag, kMCDecayDalitz) ; } //Pi0 Dalitz decay
526 else if (pPdg == 221) { SetTagBit(tag,kMCEtaDecay); SetTagBit(tag, kMCDecayDalitz) ; } //Eta Dalitz decay
527 else if((499 < pPdg && pPdg < 600)||(4999 < pPdg && pPdg < 6000)) { SetTagBit(tag,kMCEFromB); } //b-->e decay
528 else if((399 < pPdg && pPdg < 500)||(3999 < pPdg && pPdg < 5000)) { //check charm decay
531 Int_t iGrandma = parent->GetFirstMother();
534 TParticle* gma = (TParticle*)stack->Particle(iGrandma); //get mother of charm
535 Int_t gPdg = TMath::Abs(gma->GetPdgCode());
536 if((499 < gPdg && gPdg < 600)||(4999 < gPdg && gPdg < 6000)) SetTagBit(tag,kMCEFromCFromB); //b-->c-->e
537 else SetTagBit(tag,kMCEFromC); //c-->e
539 else SetTagBit(tag,kMCEFromC); //c-->e
544 //if it is not from any of the above, where is it from?
545 if(pPdg > 10000) SetTagBit(tag,kMCUnknown);
547 else SetTagBit(tag,kMCOtherDecay);
549 if(fDebug > 0 && parent) printf("AliMCAnalysisUtils::CheckOriginInStack() - Status %d Electron from other origin: %s (pPdg = %d) %s (mpdg = %d)\n",mStatus,parent->GetName(),pPdg,mom->GetName(),mPdg);
552 //Cluster was made by something else
556 printf("AliMCAnalysisUtils::CheckOriginInStack() - \tSetting kMCUnknown for cluster from %s (pdg = %d, Parent pdg = %d)\n",
557 mom->GetName(),mPdg,pPdg);
559 SetTagBit(tag,kMCUnknown);
564 if(label < 0 && (fDebug >= 0))
565 printf("AliMCAnalysisUtils::CheckOriginInStack() *** bad label or no stack ***: label %d \n", label);
567 if(label >= stack->GetNtrack() && (fDebug >= 0))
568 printf("AliMCAnalysisUtils::CheckOriginInStack() *** large label ***: label %d, n tracks %d \n", label, stack->GetNtrack());
570 SetTagBit(tag,kMCUnknown);
578 //________________________________________________________________________________________________________
579 Int_t AliMCAnalysisUtils::CheckOriginInAOD(const Int_t *labels, Int_t nlabels,
580 const TClonesArray *mcparticles, const TObjArray* arrayCluster)
582 // Play with the MCParticles in AOD if available. Tag particles depending on their origin.
583 // Do same things as in CheckOriginInStack but different input.
587 printf("AliMCAnalysisUtils::CheckOriginInAOD() - AODMCParticles is not available, check analysis settings in configuration file!!\n");
592 Int_t label=labels[0];//Most significant particle contributing to the cluster
594 Int_t nprimaries = mcparticles->GetEntriesFast();
595 if(label >= 0 && label < nprimaries)
598 AliAODMCParticle * mom = (AliAODMCParticle *) mcparticles->At(label);
600 Int_t mPdgSign = mom->GetPdgCode();
601 Int_t mPdg = TMath::Abs(mPdgSign);
602 Int_t iParent = mom->GetMother() ;
604 if(fDebug > 0 && label < 8 && fMCGenerator!="") printf("AliMCAnalysisUtils::CheckOriginInAOD() - Mother is parton %d\n",iParent);
607 AliAODMCParticle * parent = NULL ;
611 parent = (AliAODMCParticle *) mcparticles->At(iParent);
612 pPdg = TMath::Abs(parent->GetPdgCode());
614 else if(fDebug > 0 ) printf("AliMCAnalysisUtils::CheckOriginInAOD() - Parent with label %d\n",iParent);
618 printf("AliMCAnalysisUtils::CheckOriginInAOD() - Cluster most contributing mother and its parent: \n");
620 printf("\t Mother label %d, pdg %d, Primary? %d, Physical Primary? %d\n",
621 iMom, mPdg, mom->IsPrimary(), mom->IsPhysicalPrimary());
624 printf("\t Parent label %d, pdg %d, Primary? %d, Physical Primary? %d\n",
625 iParent, pPdg, parent->IsPrimary(), parent->IsPhysicalPrimary());
628 //Check if mother is converted, if not, get the first non converted mother
629 if((mPdg == 22 || mPdg == 11) && (pPdg == 22 || pPdg == 11) && !mom->IsPrimary())
631 SetTagBit(tag,kMCConversion);
633 //Check if the mother is photon or electron with status not stable
634 while ((pPdg == 22 || pPdg == 11) && !mom->IsPhysicalPrimary())
637 iMom = mom->GetMother();
638 mom = (AliAODMCParticle *) mcparticles->At(iMom);
639 mPdgSign = mom->GetPdgCode();
640 mPdg = TMath::Abs(mPdgSign);
641 iParent = mom->GetMother() ;
642 if(fDebug > 0 && label < 8 ) printf("AliMCAnalysisUtils::CheckOriginInAOD() - Mother is parton %d\n",iParent);
645 if(iParent >= 0 && parent)
647 parent = (AliAODMCParticle *) mcparticles->At(iParent);
648 pPdg = TMath::Abs(parent->GetPdgCode());
650 // printf("\t While Mother label %d, pdg %d, Primary? %d, Physical Primary? %d\n",iMom, mPdg, mom->IsPrimary(), mom->IsPhysicalPrimary());
651 // printf("\t While Parent label %d, pdg %d, Primary? %d, Physical Primary? %d\n",iParent, pPdg, parent->IsPrimary(), parent->IsPhysicalPrimary());
657 printf("AliMCAnalysisUtils::CheckOriginInAOD() - Converted photon/electron : \n");
659 printf("\t Mother label %d, pdg %d, Primary? %d, Physical Primary? %d\n",
660 iMom, mPdg, mom->IsPrimary(), mom->IsPhysicalPrimary());
663 printf("\t Parent label %d, pdg %d, Primary? %d, Physical Primary? %d\n",
664 iParent, pPdg, parent->IsPrimary(), parent->IsPhysicalPrimary());
667 }//mother and parent are electron or photon and have status 0 and parent is photon or electron
668 else if((mPdg == 22 || mPdg == 11) && !mom->IsPrimary())
670 //Still a conversion but only one electron/photon generated. Just from hadrons
671 if(pPdg == 2112 || pPdg == 211 || pPdg == 321 ||
672 pPdg == 2212 || pPdg == 130 || pPdg == 13 )
674 SetTagBit(tag,kMCConversion);
675 iMom = mom->GetMother();
676 mom = (AliAODMCParticle *) mcparticles->At(iMom);
677 mPdgSign = mom->GetPdgCode();
678 mPdg = TMath::Abs(mPdgSign);
682 printf("AliMCAnalysisUtils::CheckOriginInAOD() - Converted hadron : \n");
683 printf("\t Mother label %d, pdg %d, Primary? %d, Physical Primary? %d\n",iMom, mPdg, mom->IsPrimary(), mom->IsPhysicalPrimary());
687 //Comment for next lines, we do not check the parent of the hadron for the moment.
688 //iParent = mom->GetMother() ;
689 //if(fDebug > 0 && label < 8 ) printf("AliMCAnalysisUtils::CheckOriginInAOD() - Mother is parton %d\n",iParent);
693 // parent = (AliAODMCParticle *) mcparticles->At(iParent);
694 // pPdg = TMath::Abs(parent->GetPdgCode());
698 //printf("Final mother mPDG %d\n",mPdg);
700 // conversion into electrons/photons checked
702 //first check for typical charged particles
703 if (mPdg == 13) SetTagBit(tag,kMCMuon);
704 else if(mPdg == 211) SetTagBit(tag,kMCPion);
705 else if(mPdg == 321) SetTagBit(tag,kMCKaon);
706 else if(mPdgSign == 2212) SetTagBit(tag,kMCProton);
707 else if(mPdgSign == 2112) SetTagBit(tag,kMCNeutron);
708 else if(mPdgSign == -2212) SetTagBit(tag,kMCAntiProton);
709 else if(mPdgSign == -2112) SetTagBit(tag,kMCAntiNeutron);
711 //check for pi0 and eta (shouldn't happen unless their decays were turned off)
714 SetTagBit(tag,kMCPi0Decay);
716 if(fDebug > 2 ) printf("AliMCAnalysisUtils::CheckOriginInAOD() - First mother is directly pi0, not decayed by generator \n");
718 CheckOverlapped2GammaDecay(labels,nlabels, iMom, mcparticles, tag); //set to kMCPi0 if 2 gammas in same cluster
722 SetTagBit(tag,kMCEtaDecay);
724 if(fDebug > 2 ) printf("AliMCAnalysisUtils::CheckOriginInAOD() - First mother is directly eta, not decayed by generator \n");
726 CheckOverlapped2GammaDecay(labels,nlabels, iMom, mcparticles, tag); //set to kMCEta if 2 gammas in same cluster
731 SetTagBit(tag,kMCPhoton);
735 SetTagBit(tag,kMCPi0Decay);
737 if(fDebug > 2 ) printf("AliMCAnalysisUtils::CheckOriginInAOD() - Generator pi0 decay photon \n");
739 CheckOverlapped2GammaDecay(labels,nlabels, iParent, mcparticles, tag); //set to kMCPi0 if 2 gammas in same cluster
740 // In case it did not merge, check if the decay companion is lost
741 if(!CheckTagBit(tag, kMCPi0) && !CheckTagBit(tag,kMCDecayPairInCalo) && !CheckTagBit(tag,kMCDecayPairLost))
742 CheckLostDecayPair(arrayCluster,iMom, iParent, mcparticles, tag);
744 //printf("Bit set is Merged %d, Pair in calo %d, Lost %d\n",CheckTagBit(tag, kMCPi0),CheckTagBit(tag,kMCDecayPairInCalo),CheckTagBit(tag,kMCDecayPairLost));
746 else if (pPdg == 221)
748 SetTagBit(tag, kMCEtaDecay);
750 if(fDebug > 2 ) printf("AliMCAnalysisUtils::CheckOriginInAOD() - Generator eta decay photon \n");
752 CheckOverlapped2GammaDecay(labels,nlabels, iParent, mcparticles, tag); //set to kMCEta if 2 gammas in same cluster
753 // In case it did not merge, check if the decay companion is lost
754 if(!CheckTagBit(tag, kMCEta) && !CheckTagBit(tag,kMCDecayPairInCalo))
755 CheckLostDecayPair(arrayCluster,iMom, iParent, mcparticles, tag);
757 else if(mom->IsPhysicalPrimary() && (fMCGenerator=="PYTHIA" || fMCGenerator=="HERWIG")) //undecayed particle
759 if(iParent < 8 && iParent > 5 )
761 if(pPdg == 22) SetTagBit(tag,kMCPrompt);
762 else SetTagBit(tag,kMCFragmentation);
764 else if(iParent <= 5 && (fMCGenerator=="PYTHIA" || fMCGenerator=="HERWIG"))
766 SetTagBit(tag, kMCISR); //Initial state radiation
768 else SetTagBit(tag,kMCUnknown);
770 else SetTagBit(tag,kMCOtherDecay);
774 //Electron check. Where did that electron come from?
777 if(pPdg == 11 && parent)
779 Int_t iGrandma = parent->GetMother();
782 AliAODMCParticle* gma = (AliAODMCParticle*)mcparticles->At(iGrandma);
783 Int_t gPdg = TMath::Abs(gma->GetPdgCode());
785 if (gPdg == 23) { SetTagBit(tag,kMCZDecay); } //parent is Z-boson
786 else if (gPdg == 24) { SetTagBit(tag,kMCWDecay); } //parent is W-boson
790 SetTagBit(tag,kMCElectron);
792 if(fDebug > 0) printf("AliMCAnalysisUtils::CheckOriginInAOD() - Checking ancestors of electrons");
794 if (pPdg == 111) { SetTagBit(tag,kMCPi0Decay); SetTagBit(tag,kMCDecayDalitz);} //Pi0 Dalitz decay
795 else if (pPdg == 221) { SetTagBit(tag,kMCEtaDecay); SetTagBit(tag,kMCDecayDalitz);} //Eta Dalitz decay
796 else if((499 < pPdg && pPdg < 600)||(4999 < pPdg && pPdg < 6000)) { SetTagBit(tag,kMCEFromB);} //b-hadron decay
797 else if((399 < pPdg && pPdg < 500)||(3999 < pPdg && pPdg < 5000))
798 { //c-hadron decay check
801 Int_t iGrandma = parent->GetMother();
804 AliAODMCParticle* gma = (AliAODMCParticle*)mcparticles->At(iGrandma); //charm's mother
805 Int_t gPdg = TMath::Abs(gma->GetPdgCode());
806 if((499 < gPdg && gPdg < 600)||(4999 < gPdg && gPdg < 6000)) SetTagBit(tag,kMCEFromCFromB); //b-->c-->e decay
807 else SetTagBit(tag,kMCEFromC); //c-hadron decay
809 else SetTagBit(tag,kMCEFromC); //c-hadron decay
812 { //prompt or other decay
813 TParticlePDG* foo = TDatabasePDG::Instance()->GetParticle(pPdg);
814 TParticlePDG* foo1 = TDatabasePDG::Instance()->GetParticle(mPdg);
816 if(fDebug > 0) printf("AliMCAnalysisUtils::CheckOriginInAOD() - Electron from other origin: %s (pPdg = %d) %s (mPdg = %d)\n",
817 foo->GetName(), pPdg,foo1->GetName(),mPdg);
819 if(pPdg > 10000) SetTagBit(tag,kMCUnknown);
820 else SetTagBit(tag,kMCOtherDecay);
823 //cluster was made by something else
826 if(fDebug > 0) printf("AliMCAnalysisUtils::CheckOriginInAOD() - \tSetting kMCUnknown for cluster with pdg = %d, Parent pdg = %d\n",
828 SetTagBit(tag,kMCUnknown);
833 if(label < 0 && (fDebug >= 0) )
834 printf("AliMCAnalysisUtils::CheckOriginInAOD() *** bad label or no mcparticles ***: label %d \n", label);
836 if(label >= mcparticles->GetEntriesFast() && (fDebug >= 0) )
837 printf("AliMCAnalysisUtils::CheckOriginInAOD() *** large label ***: label %d, n tracks %d \n", label, mcparticles->GetEntriesFast());
839 SetTagBit(tag,kMCUnknown);
847 //_________________________________________________________________________________________
848 void AliMCAnalysisUtils::CheckOverlapped2GammaDecay(const Int_t *labels, Int_t nlabels,
849 Int_t mesonIndex, AliStack *stack,
852 //Check if cluster is formed from the contribution of 2 decay photons from pi0 or eta. Input in stack
854 if(labels[0] < 0 || labels[0] > stack->GetNtrack() || nlabels <= 1) {
855 if(fDebug > 2) printf("AliMCAnalysisUtils::CheckOverlapped2GammaDecay(ESD) - Exit : label[0] %d, n primaries %d, nlabels %d \n",
856 labels[0],stack->GetNtrack(), nlabels);
860 TParticle * meson = stack->Particle(mesonIndex);
861 Int_t mesonPdg = meson->GetPdgCode();
862 if(mesonPdg!=111 && mesonPdg!=221)
864 printf("AliMCAnalysisUtils::CheckOverlapped2GammaDecay(ESD) - Wrong pi0/eta PDG : %d \n",mesonPdg);
868 if(fDebug > 2) printf("AliMCAnalysisUtils::CheckOverlapped2GammaDecay(ESD) - %s, label %d\n",meson->GetName(), mesonIndex);
870 //Check if meson decayed into 2 daughters or if both were kept.
871 if(meson->GetNDaughters() != 2)
874 printf("AliMCAnalysisUtils::CheckOverlapped2GammaDecay(ESD) - Not overalapped. Number of daughters is %d, not 2 \n",meson->GetNDaughters());
879 Int_t iPhoton0 = meson->GetDaughter(0);
880 Int_t iPhoton1 = meson->GetDaughter(1);
881 TParticle *photon0 = stack->Particle(iPhoton0);
882 TParticle *photon1 = stack->Particle(iPhoton1);
884 //Check if both daughters are photons
885 if(photon0->GetPdgCode() != 22 || photon1->GetPdgCode()!=22)
888 printf("AliMCAnalysisUtils::CheckOverlapped2GammaDecay(ESD) - Not overalapped. PDG: daughter 1 = %d, of daughter 2 = %d \n",photon0->GetPdgCode(),photon1->GetPdgCode());
893 printf("AliMCAnalysisUtils::CheckOverlapped2GammaDecay(ESD) - Daughter labels : photon0 = %d, photon1 = %d \n",iPhoton0,iPhoton1);
895 //Check if both photons contribute to the cluster
896 Bool_t okPhoton0 = kFALSE;
897 Bool_t okPhoton1 = kFALSE;
899 if(fDebug > 3) printf("AliMCAnalysisUtils::CheckOverlapped2GammaDecay(ESD) - Labels loop:\n");
901 Bool_t conversion = kFALSE;
903 for(Int_t i = 0; i < nlabels; i++)
905 if(fDebug > 3) printf("\t at begin:label %d/%d: %d, ok? photon1 %d, photon2 %d\n", i+1, nlabels, labels[i], okPhoton0, okPhoton1);
907 //If we already found both, break the loop
908 if(okPhoton0 && okPhoton1) break;
910 Int_t index = labels[i];
911 if (iPhoton0 == index)
916 else if (iPhoton1 == index)
922 //Trace back the mother in case it was a conversion
924 if(index >= stack->GetNtrack())
926 printf("AliMCAnalysisUtils::CheckOverlapped2GammaDecay(ESD) Particle index %d larger than size of list %d !!\n",
927 index,stack->GetNtrack());
931 TParticle * daught = stack->Particle(index);
932 Int_t tmpindex = daught->GetFirstMother();
933 if(fDebug > 3) printf("\t Conversion? : mother %d\n",tmpindex);
936 //MC particle of interest is the mother
937 if(fDebug > 3) printf("\t \t parent index %d\n",tmpindex);
938 daught = stack->Particle(tmpindex);
939 if (iPhoton0 == tmpindex)
945 else if (iPhoton1 == tmpindex)
952 tmpindex = daught->GetFirstMother();
954 }//While to check if pi0/eta daughter was one of these contributors to the cluster
956 if(i == 0 && (!okPhoton0 && !okPhoton1) && fDebug>=0)
957 printf("AliMCAnalysisUtils::CheckOverlapped2GammaDecay(ESD) - Something happens, first label should be from a photon decay!\n");
959 }//loop on list of labels
961 //If both photons contribute tag as the corresponding meson.
962 if(okPhoton0 && okPhoton1)
965 printf("AliMCAnalysisUtils::CheckOverlapped2GammaDecay(ESD) - %s OVERLAPPED DECAY \n", meson->GetName());
967 if(!CheckTagBit(tag,kMCConversion) && conversion) SetTagBit(tag,kMCConversion) ;
969 if(mesonPdg == 111) SetTagBit(tag,kMCPi0);
970 else SetTagBit(tag,kMCEta);
975 //________________________________________________________________________________________________________
976 void AliMCAnalysisUtils::CheckOverlapped2GammaDecay(const Int_t *labels, Int_t nlabels, Int_t mesonIndex,
977 const TClonesArray *mcparticles, Int_t & tag )
979 //Check if cluster is formed from the contribution of 2 decay photons from pi0 or eta. Input in AODMCParticles
981 if(labels[0] < 0 || labels[0] > mcparticles->GetEntriesFast() || nlabels <= 1)
984 printf("AliMCAnalysisUtils::CheckOverlapped2GammaDecay(AOD) - Exit : label[0] %d, n primaries %d, nlabels %d \n",
985 labels[0],mcparticles->GetEntriesFast(), nlabels);
989 AliAODMCParticle * meson = (AliAODMCParticle *) mcparticles->At(mesonIndex);
990 Int_t mesonPdg = meson->GetPdgCode();
991 if(mesonPdg != 111 && mesonPdg != 221)
993 printf("AliMCAnalysisUtils::CheckOverlapped2GammaDecay(AOD) - Wrong pi0/eta PDG : %d \n",mesonPdg);
997 if(fDebug > 2) printf("AliMCAnalysisUtils::CheckOverlapped2GammaDecay(AOD) - pdg %d, label %d, ndaughters %d\n", mesonPdg, mesonIndex, meson->GetNDaughters());
1001 if(meson->GetNDaughters() != 2)
1004 printf("AliMCAnalysisUtils::CheckOverlapped2GammaDecay(ESD) - Not overalapped. Number of daughters is %d, not 2 \n",meson->GetNDaughters());
1008 Int_t iPhoton0 = meson->GetDaughter(0);
1009 Int_t iPhoton1 = meson->GetDaughter(1);
1010 //if((iPhoton0 == -1) || (iPhoton1 == -1)){
1012 // printf("AliMCAnalysisUtils::CheckOverlapped2GammaDecay(AOD) - Exit : Not overlapped. At least a daughter do not exists : d1 %d, d2 %d \n", iPhoton0, iPhoton1);
1015 AliAODMCParticle *photon0 = (AliAODMCParticle *) mcparticles->At(iPhoton0);
1016 AliAODMCParticle *photon1 = (AliAODMCParticle *) mcparticles->At(iPhoton1);
1018 //Check if both daughters are photons
1019 if(photon0->GetPdgCode() != 22 && photon1->GetPdgCode()!=22)
1021 printf("AliMCAnalysisUtils::CheckOverlapped2GammaDecay(AOD) - Not overlapped. PDG: daughter 1 = %d, of daughter 2 = %d \n",photon0->GetPdgCode(),photon1->GetPdgCode());
1026 printf("AliMCAnalysisUtils::CheckOverlapped2GammaDecay(AOD) - Daughter labels : photon0 = %d, photon1 = %d \n",iPhoton0,iPhoton1);
1028 //Check if both photons contribute to the cluster
1029 Bool_t okPhoton0 = kFALSE;
1030 Bool_t okPhoton1 = kFALSE;
1033 printf("AliMCAnalysisUtils::CheckOverlapped2GammaDecay(AOD) - Labels loop:\n");
1035 Bool_t conversion = kFALSE;
1037 for(Int_t i = 0; i < nlabels; i++)
1040 printf("\t label %d/%d: %d, ok? %d, %d\n", i, nlabels, labels[i], okPhoton0, okPhoton1);
1042 if(labels[i]<0) continue;
1044 //If we already found both, break the loop
1045 if(okPhoton0 && okPhoton1) break;
1047 Int_t index = labels[i];
1048 if (iPhoton0 == index)
1053 else if (iPhoton1 == index)
1059 //Trace back the mother in case it was a conversion
1061 if(index >= mcparticles->GetEntriesFast())
1063 printf("AliMCAnalysisUtils::CheckOverlapped2GammaDecay(AOD) Particle index %d larger than size of list %d !!\n",index,mcparticles->GetEntriesFast());
1067 AliAODMCParticle * daught = (AliAODMCParticle*) mcparticles->At(index);
1068 Int_t tmpindex = daught->GetMother();
1070 printf("AliMCAnalysisUtils::CheckOverlapped2GammaDecay(AOD) - Conversion? : mother %d\n",tmpindex);
1074 //MC particle of interest is the mother
1076 printf("\t parent index %d\n",tmpindex);
1077 daught = (AliAODMCParticle*) mcparticles->At(tmpindex);
1078 //printf("tmpindex %d\n",tmpindex);
1079 if (iPhoton0 == tmpindex)
1085 else if (iPhoton1 == tmpindex)
1092 tmpindex = daught->GetMother();
1094 }//While to check if pi0/eta daughter was one of these contributors to the cluster
1096 if(i == 0 && (!okPhoton0 && !okPhoton1) && fDebug>=0 )
1097 printf("AliMCAnalysisUtils::CheckOverlapped2GammaDecay(AOD) - Something happens, first label should be from a photon decay!\n");
1099 }//loop on list of labels
1101 //If both photons contribute tag as the corresponding meson.
1102 if(okPhoton0 && okPhoton1)
1105 printf("AliMCAnalysisUtils::CheckOverlapped2GammaDecay(AOD) - %s OVERLAPPED DECAY \n",
1106 (TDatabasePDG::Instance()->GetParticle(mesonPdg))->GetName());
1108 if(!CheckTagBit(tag,kMCConversion) && conversion)
1111 printf("AliMCAnalysisUtils::CheckOverlapped2GammaDecay(AOD) Second decay photon produced a conversion \n");
1112 SetTagBit(tag,kMCConversion) ;
1115 if(mesonPdg == 111) SetTagBit(tag,kMCPi0);
1116 else SetTagBit(tag,kMCEta);
1121 //______________________________________________________________________________________________________
1122 void AliMCAnalysisUtils::CheckLostDecayPair(const TObjArray* arrayCluster, Int_t iMom, Int_t iParent,
1123 AliStack * stack, Int_t & tag)
1125 // Check on ESDs if the current decay photon has the second photon companion lost
1127 if(!arrayCluster || iMom < 0 || iParent < 0|| !stack) return;
1129 TParticle * parent= stack->Particle(iParent);
1131 if(parent->GetNDaughters()!=2)
1133 SetTagBit(tag, kMCDecayPairLost);
1137 Int_t pairLabel = -1;
1138 if ( iMom != parent->GetDaughter(0) ) pairLabel = parent->GetDaughter(0);
1139 else if( iMom != parent->GetDaughter(1) ) pairLabel = parent->GetDaughter(1);
1143 SetTagBit(tag, kMCDecayPairLost);
1147 for(Int_t iclus = 0; iclus < arrayCluster->GetEntriesFast(); iclus++)
1149 AliVCluster * cluster = (AliVCluster*) arrayCluster->At(iclus);
1150 for(UInt_t ilab = 0; ilab< cluster->GetNLabels(); ilab++)
1152 Int_t label = cluster->GetLabels()[ilab];
1153 if(label==pairLabel)
1155 SetTagBit(tag, kMCDecayPairInCalo);
1158 else if(label== iParent || label== iMom)
1162 else // check the ancestry
1164 TParticle * mother = stack->Particle(label);
1165 Int_t momPDG = TMath::Abs(mother->GetPdgCode());
1166 if(momPDG!=11 && momPDG!=22) continue;
1168 //Check if "mother" of entity is converted, if not, get the first non converted mother
1169 Int_t iParentClus = mother->GetFirstMother();
1170 if(iParentClus < 0) continue;
1172 TParticle * parentClus = stack->Particle(iParentClus);
1173 if(!parentClus) continue;
1175 Int_t parentClusPDG = TMath::Abs(parentClus->GetPdgCode());
1176 Int_t parentClusStatus = parentClus->GetStatusCode();
1178 if( parentClusPDG != 22 && parentClusPDG != 11 && parentClusStatus != 0) continue;
1180 //printf("Conversion\n");
1182 //Check if the mother is photon or electron with status not stable
1183 while ((parentClusPDG == 22 || parentClusPDG == 11) && parentClusStatus != 1)
1186 label = iParentClus;
1187 momPDG = parentClusPDG;
1189 iParentClus = parentClus->GetFirstMother();
1190 if(iParentClus < 0) break;
1192 parentClus = stack->Particle(iParentClus);
1193 if(!parentClus) break;
1195 parentClusPDG = TMath::Abs(parentClus->GetPdgCode());
1196 parentClusStatus = parentClus->GetStatusCode() ;
1199 if((momPDG == 22 || parentClusPDG ==22) && (label==pairLabel || iParentClus == pairLabel))
1201 SetTagBit(tag, kMCDecayPairInCalo);
1202 //printf("Conversion is paired\n");
1211 SetTagBit(tag, kMCDecayPairLost);
1215 //______________________________________________________________________________________________________
1216 void AliMCAnalysisUtils::CheckLostDecayPair(const TObjArray * arrayCluster,Int_t iMom, Int_t iParent,
1217 const TClonesArray* mcparticles, Int_t & tag)
1219 // Check on AODs if the current decay photon has the second photon companion lost
1221 if(!arrayCluster || iMom < 0 || iParent < 0|| !mcparticles) return;
1223 AliAODMCParticle * parent = (AliAODMCParticle*) mcparticles->At(iParent);
1225 //printf("*** Check label %d with parent %d\n",iMom, iParent);
1227 if(parent->GetNDaughters()!=2)
1229 SetTagBit(tag, kMCDecayPairLost);
1230 //printf("\t ndaugh = %d\n",parent->GetNDaughters());
1234 Int_t pairLabel = -1;
1235 if ( iMom != parent->GetDaughter(0) ) pairLabel = parent->GetDaughter(0);
1236 else if( iMom != parent->GetDaughter(1) ) pairLabel = parent->GetDaughter(1);
1240 //printf("\t pair Label not found = %d\n",pairLabel);
1241 SetTagBit(tag, kMCDecayPairLost);
1245 //printf("\t *** find pair %d\n",pairLabel);
1247 for(Int_t iclus = 0; iclus < arrayCluster->GetEntriesFast(); iclus++)
1249 AliVCluster * cluster = (AliVCluster*) arrayCluster->At(iclus);
1250 //printf("\t \t ** Cluster %d, nlabels %d\n",iclus,cluster->GetNLabels());
1251 for(UInt_t ilab = 0; ilab< cluster->GetNLabels(); ilab++)
1253 Int_t label = cluster->GetLabels()[ilab];
1255 //printf("\t \t label %d\n",label);
1257 if(label==pairLabel)
1259 //printf("\t \t Pair found\n");
1260 SetTagBit(tag, kMCDecayPairInCalo);
1263 else if(label== iParent || label== iMom)
1265 //printf("\t \t skip\n");
1268 else // check the ancestry
1270 AliAODMCParticle * mother = (AliAODMCParticle*) mcparticles->At(label);
1271 Int_t momPDG = TMath::Abs(mother->GetPdgCode());
1272 if(momPDG!=11 && momPDG!=22)
1274 //printf("\t \t skip, pdg %d\n",momPDG);
1278 //Check if "mother" of entity is converted, if not, get the first non converted mother
1279 Int_t iParentClus = mother->GetMother();
1280 if(iParentClus < 0) continue;
1282 AliAODMCParticle * parentClus = (AliAODMCParticle*) mcparticles->At(iParentClus);
1283 if(!parentClus) continue;
1285 Int_t parentClusPDG = TMath::Abs(parentClus->GetPdgCode());
1286 Int_t parentClusStatus = parentClus->GetStatus();
1288 if( parentClusPDG != 22 && parentClusPDG != 11 && parentClusStatus != 0)
1290 //printf("\t \t skip, not a conversion, parent: pdg %d, status %d\n",parentClusPDG,parentClusStatus);
1294 //printf("\t \t Conversion\n");
1296 //Check if the mother is photon or electron with status not stable
1297 while ((parentClusPDG == 22 || parentClusPDG == 11) && parentClusStatus != 1)
1300 label = iParentClus;
1301 momPDG = parentClusPDG;
1303 iParentClus = parentClus->GetMother();
1304 if(iParentClus < 0) break;
1306 parentClus = (AliAODMCParticle*) mcparticles->At(iParentClus);
1307 if(!parentClus) break;
1309 parentClusPDG = TMath::Abs(parentClus->GetPdgCode());
1310 parentClusStatus = parentClus->GetStatus() ;
1313 if((momPDG == 22 || parentClusPDG ==22) && (label==pairLabel || iParentClus == pairLabel))
1315 SetTagBit(tag, kMCDecayPairInCalo);
1316 //printf("\t \t Conversion is paired: mom %d, parent %d\n",label,iParentClus);
1321 //printf("\t \t Skip, finally label %d, pdg %d, parent label %d, pdg %d, status %d\n",label,momPDG,iParentClus,parentClusPDG,parentClusStatus);
1330 SetTagBit(tag, kMCDecayPairLost);
1335 //_____________________________________________________________________
1336 TList * AliMCAnalysisUtils::GetJets(const AliCaloTrackReader * reader)
1338 //Return list of jets (TParticles) and index of most likely parton that originated it.
1339 AliStack * stack = reader->GetStack();
1340 Int_t iEvent = reader->GetEventNumber();
1341 AliGenEventHeader * geh = reader->GetGenEventHeader();
1342 if(fCurrentEvent!=iEvent){
1343 fCurrentEvent = iEvent;
1344 fJetsList = new TList;
1345 Int_t nTriggerJets = 0;
1346 Float_t tmpjet[]={0,0,0,0};
1348 //printf("Event %d %d\n",fCurrentEvent,iEvent);
1349 //Get outgoing partons
1350 if(stack->GetNtrack() < 8) return fJetsList;
1351 TParticle * parton1 = stack->Particle(6);
1352 TParticle * parton2 = stack->Particle(7);
1354 printf("AliMCAnalysisUtils::GetJets() - parton 6 : %s, pt %2.2f,E %2.2f, phi %2.2f, eta %2.2f \n",
1355 parton1->GetName(),parton1->Pt(),parton1->Energy(),parton1->Phi()*TMath::RadToDeg(),parton1->Eta());
1356 printf("AliMCAnalysisUtils::GetJets() - parton 7 : %s, pt %2.2f,E %2.2f, phi %2.2f, eta %2.2f \n",
1357 parton2->GetName(),parton2->Pt(),parton2->Energy(),parton2->Phi()*TMath::RadToDeg(),parton2->Eta());
1359 // //Trace the jet from the mother parton
1366 // TParticle * tmptmp = new TParticle;
1367 // for(Int_t i = 0; i< stack->GetNprimary(); i++){
1368 // tmptmp = stack->Particle(i);
1370 // if(tmptmp->GetStatusCode() == 1){
1371 // pt = tmptmp->Pt();
1372 // e = tmptmp->Energy();
1373 // Int_t imom = tmptmp->GetFirstMother();
1375 // //printf("1st imom %d\n",imom);
1378 // tmptmp = stack->Particle(imom);
1379 // imom = tmptmp->GetFirstMother();
1380 // //printf("imom %d \n",imom);
1382 // //printf("Last imom %d %d\n",imom1, imom);
1387 // else if (imom1 == 7){
1394 // printf("JET 1, pt %2.2f, e %2.2f; JET 2, pt %2.2f, e %2.2f \n",pt1,e1,pt2,e2);
1396 //Get the jet, different way for different generator
1398 if(fMCGenerator == "PYTHIA")
1400 TParticle * jet = 0x0;
1401 AliGenPythiaEventHeader* pygeh= (AliGenPythiaEventHeader*) geh;
1402 nTriggerJets = pygeh->NTriggerJets();
1404 printf("AliMCAnalysisUtils::GetJets() - PythiaEventHeader: Njets: %d\n",nTriggerJets);
1406 for(Int_t i = 0; i< nTriggerJets; i++)
1408 pygeh->TriggerJet(i, tmpjet);
1409 jet = new TParticle(94, 21, -1, -1, -1, -1, tmpjet[0],tmpjet[1],tmpjet[2],tmpjet[3], 0,0,0,0);
1410 //Assign an outgoing parton as mother
1411 Float_t phidiff1 = TMath::Abs(jet->Phi()-parton1->Phi());
1412 Float_t phidiff2 = TMath::Abs(jet->Phi()-parton2->Phi());
1413 if(phidiff1 > phidiff2) jet->SetFirstMother(7);
1414 else jet->SetFirstMother(6);
1417 printf("AliMCAnalysisUtils::GetJets() - PYTHIA Jet %d: mother %d, pt %2.2f,E %2.2f, phi %2.2f, eta %2.2f \n",
1418 i, jet->GetFirstMother(),jet->Pt(),jet->Energy(),jet->Phi()*TMath::RadToDeg(),jet->Eta());
1419 fJetsList->Add(jet);
1421 }//Pythia triggered jets
1423 else if (fMCGenerator=="HERWIG")
1427 TParticle * tmp = parton1;
1428 if(parton1->GetPdgCode()!=22)
1431 if(tmp->GetFirstDaughter()==-1) return fJetsList;
1432 tmp = stack->Particle(tmp->GetFirstDaughter());
1433 pdg = tmp->GetPdgCode();
1436 //Add found jet to list
1437 TParticle *jet1 = new TParticle(*tmp);
1438 jet1->SetFirstMother(6);
1439 fJetsList->Add(jet1);
1440 //printf("jet 1: first daughter %d, last daughter %d\n", tmp->GetFirstDaughter(), tmp->GetLastDaughter());
1441 //tmp = stack->Particle(tmp->GetFirstDaughter());
1445 printf("AliMCAnalysisUtils::GetJets() - HERWIG Jet 1: mother %d, status %d, pt %2.2f,E %2.2f, phi %2.2f, eta %2.2f \n",
1446 jet1->GetFirstMother(),jet1->GetStatusCode(),jet1->Pt(),jet1->Energy(),jet1->Phi()*TMath::RadToDeg(),jet1->Eta());
1452 if(parton2->GetPdgCode()!=22)
1456 if(tmp->GetFirstDaughter()==-1) return fJetsList;
1457 tmp = stack->Particle(tmp->GetFirstDaughter());
1458 pdg = tmp->GetPdgCode();
1461 //Add found jet to list
1462 TParticle *jet2 = new TParticle(*tmp);
1463 jet2->SetFirstMother(7);
1464 fJetsList->Add(jet2);
1467 printf("AliMCAnalysisUtils::GetJets() - HERWIG Jet 2: mother %d, status %d, pt %2.2f,E %2.2f, phi %2.2f, eta %2.2f \n",
1468 jet2->GetFirstMother(),jet2->GetStatusCode(),jet2->Pt(),jet2->Energy(),jet2->Phi()*TMath::RadToDeg(),jet2->Eta());
1469 //Int_t first = tmp->GetFirstDaughter();
1470 //Int_t last = tmp->GetLastDaughter();
1471 //printf("jet 2: first daughter %d, last daughter %d, pdg %d\n",first, last, tmp->GetPdgCode());
1472 // for(Int_t d = first ; d < last+1; d++){
1473 // tmp = stack->Particle(d);
1474 // if(i == tmp->GetFirstMother())
1475 // printf("Daughter n %d, Mother %d, name %s, status %d, pT %2.2f,E %2.2f, phi %2.2f, eta %2.2f \n",
1476 // d,tmp->GetFirstMother(), tmp->GetName(), tmp->GetStatusCode(),tmp->Pt(),tmp->Energy(),tmp->Phi()*TMath::RadToDeg(),tmp->Eta());
1480 }//Herwig generated jets
1487 //__________________________________________________________________________________________________________
1488 TLorentzVector AliMCAnalysisUtils::GetDaughter(Int_t idaugh, Int_t label,
1489 const AliCaloTrackReader* reader,
1490 Int_t & pdg, Int_t & status, Bool_t & ok, Int_t & daughlabel)
1492 //Return the kinematics of the particle that generated the signal, its pdg and its status and its label mother
1494 TLorentzVector daughter(0,0,0,0);
1496 if(reader->ReadStack())
1498 if(!reader->GetStack())
1501 printf("AliMCAnalysisUtils::GetDaughter() - Stack is not available, check analysis settings in configuration file, STOP!!\n");
1507 Int_t nprimaries = reader->GetStack()->GetNtrack();
1508 if(label >= 0 && label < nprimaries)
1510 TParticle * momP = reader->GetStack()->Particle(label);
1511 daughlabel = momP->GetDaughter(idaugh);
1513 if(daughlabel < 0 || daughlabel >= nprimaries)
1519 TParticle * daughP = reader->GetStack()->Particle(daughlabel);
1520 daughP->Momentum(daughter);
1521 pdg = daughP->GetPdgCode();
1522 status = daughP->GetStatusCode();
1530 else if(reader->ReadAODMCParticles())
1532 TClonesArray* mcparticles = reader->GetAODMCParticles();
1536 printf("AliMCAnalysisUtils::GetDaughter() - AODMCParticles is not available, check analysis settings in configuration file!!\n");
1542 Int_t nprimaries = mcparticles->GetEntriesFast();
1543 if(label >= 0 && label < nprimaries)
1545 AliAODMCParticle * momP = (AliAODMCParticle *) mcparticles->At(label);
1546 daughlabel = momP->GetDaughter(idaugh);
1548 if(daughlabel < 0 || daughlabel >= nprimaries)
1554 AliAODMCParticle * daughP = (AliAODMCParticle *) mcparticles->At(daughlabel);
1555 daughter.SetPxPyPzE(daughP->Px(),daughP->Py(),daughP->Pz(),daughP->E());
1556 pdg = daughP->GetPdgCode();
1557 status = daughP->GetStatus();
1571 //______________________________________________________________________________________________________
1572 TLorentzVector AliMCAnalysisUtils::GetMother(Int_t label, const AliCaloTrackReader* reader, Bool_t & ok)
1574 //Return the kinematics of the particle that generated the signal
1576 Int_t pdg = -1; Int_t status = -1; Int_t momlabel = -1;
1577 return GetMother(label,reader,pdg,status, ok,momlabel);
1580 //_________________________________________________________________________________________
1581 TLorentzVector AliMCAnalysisUtils::GetMother(Int_t label, const AliCaloTrackReader* reader,
1582 Int_t & pdg, Int_t & status, Bool_t & ok)
1584 //Return the kinematics of the particle that generated the signal
1586 Int_t momlabel = -1;
1587 return GetMother(label,reader,pdg,status, ok,momlabel);
1590 //______________________________________________________________________________________________________
1591 TLorentzVector AliMCAnalysisUtils::GetMother(Int_t label, const AliCaloTrackReader* reader,
1592 Int_t & pdg, Int_t & status, Bool_t & ok, Int_t & momlabel)
1594 //Return the kinematics of the particle that generated the signal, its pdg and its status and its label mother
1596 TLorentzVector mom(0,0,0,0);
1598 if(reader->ReadStack())
1600 if(!reader->GetStack())
1603 printf("AliMCAnalysisUtils::GetMother() - Stack is not available, check analysis settings in configuration file, STOP!!\n");
1608 if(label >= 0 && label < reader->GetStack()->GetNtrack())
1610 TParticle * momP = reader->GetStack()->Particle(label);
1611 momP->Momentum(mom);
1612 pdg = momP->GetPdgCode();
1613 status = momP->GetStatusCode();
1614 momlabel = momP->GetFirstMother();
1622 else if(reader->ReadAODMCParticles())
1624 TClonesArray* mcparticles = reader->GetAODMCParticles();
1628 printf("AliMCAnalysisUtils::GetMother() - AODMCParticles is not available, check analysis settings in configuration file!!\n");
1634 Int_t nprimaries = mcparticles->GetEntriesFast();
1635 if(label >= 0 && label < nprimaries)
1637 AliAODMCParticle * momP = (AliAODMCParticle *) mcparticles->At(label);
1638 mom.SetPxPyPzE(momP->Px(),momP->Py(),momP->Pz(),momP->E());
1639 pdg = momP->GetPdgCode();
1640 status = momP->GetStatus();
1641 momlabel = momP->GetMother();
1655 //___________________________________________________________________________________
1656 TLorentzVector AliMCAnalysisUtils::GetMotherWithPDG(Int_t label, Int_t pdg,
1657 const AliCaloTrackReader* reader,
1658 Bool_t & ok, Int_t & momlabel)
1660 //Return the kinematics of the particle that generated the signal
1662 TLorentzVector grandmom(0,0,0,0);
1665 if(reader->ReadStack())
1667 if(!reader->GetStack())
1670 printf("AliMCAnalysisUtils::GetMotherWithPDG() - Stack is not available, check analysis settings in configuration file, STOP!!\n");
1676 if(label >= 0 && label < reader->GetStack()->GetNtrack())
1678 TParticle * momP = reader->GetStack()->Particle(label);
1680 Int_t grandmomLabel = momP->GetFirstMother();
1681 Int_t grandmomPDG = -1;
1682 TParticle * grandmomP = 0x0;
1683 while (grandmomLabel >=0 )
1685 grandmomP = reader->GetStack()->Particle(grandmomLabel);
1686 grandmomPDG = grandmomP->GetPdgCode();
1687 if(grandmomPDG==pdg)
1689 momlabel = grandmomLabel;
1690 grandmom.SetPxPyPzE(grandmomP->Px(),grandmomP->Py(),grandmomP->Pz(),grandmomP->Energy());
1694 grandmomLabel = grandmomP->GetFirstMother();
1698 if(grandmomPDG!=pdg) printf("AliMCAnalysisUtils::GetMotherWithPDG(ESD) - mother with PDG %d, not found! \n",pdg);
1701 else if(reader->ReadAODMCParticles())
1703 TClonesArray* mcparticles = reader->GetAODMCParticles();
1707 printf("AliMCAnalysisUtils::GetMotherWithPDG() - AODMCParticles is not available, check analysis settings in configuration file!!\n");
1713 Int_t nprimaries = mcparticles->GetEntriesFast();
1714 if(label >= 0 && label < nprimaries)
1716 AliAODMCParticle * momP = (AliAODMCParticle *) mcparticles->At(label);
1718 Int_t grandmomLabel = momP->GetMother();
1719 Int_t grandmomPDG = -1;
1720 AliAODMCParticle * grandmomP = 0x0;
1721 while (grandmomLabel >=0 )
1723 grandmomP = (AliAODMCParticle *) mcparticles->At(grandmomLabel);
1724 grandmomPDG = grandmomP->GetPdgCode();
1725 if(grandmomPDG==pdg)
1727 //printf("AliMCAnalysisUtils::GetMotherWithPDG(AOD) - mother with PDG %d FOUND! \n",pdg);
1728 momlabel = grandmomLabel;
1729 grandmom.SetPxPyPzE(grandmomP->Px(),grandmomP->Py(),grandmomP->Pz(),grandmomP->E());
1733 grandmomLabel = grandmomP->GetMother();
1737 if(grandmomPDG!=pdg) printf("AliMCAnalysisUtils::GetMotherWithPDG(AOD) - mother with PDG %d, NOT found! \n",pdg);
1747 //______________________________________________________________________________________________
1748 TLorentzVector AliMCAnalysisUtils::GetGrandMother(Int_t label, const AliCaloTrackReader* reader,
1749 Int_t & pdg, Int_t & status, Bool_t & ok,
1750 Int_t & grandMomLabel, Int_t & greatMomLabel)
1752 //Return the kinematics of the particle that generated the signal
1754 TLorentzVector grandmom(0,0,0,0);
1756 if(reader->ReadStack())
1758 if(!reader->GetStack())
1761 printf("AliMCAnalysisUtils::GetGrandMother() - Stack is not available, check analysis settings in configuration file, STOP!!\n");
1767 if(label >= 0 && label < reader->GetStack()->GetNtrack())
1769 TParticle * momP = reader->GetStack()->Particle(label);
1771 grandMomLabel = momP->GetFirstMother();
1773 TParticle * grandmomP = 0x0;
1775 if (grandMomLabel >=0 )
1777 grandmomP = reader->GetStack()->Particle(grandMomLabel);
1778 pdg = grandmomP->GetPdgCode();
1779 status = grandmomP->GetStatusCode();
1781 grandmom.SetPxPyPzE(grandmomP->Px(),grandmomP->Py(),grandmomP->Pz(),grandmomP->Energy());
1782 greatMomLabel = grandmomP->GetFirstMother();
1787 else if(reader->ReadAODMCParticles())
1789 TClonesArray* mcparticles = reader->GetAODMCParticles();
1793 printf("AliMCAnalysisUtils::GetGrandMother() - AODMCParticles is not available, check analysis settings in configuration file!!\n");
1799 Int_t nprimaries = mcparticles->GetEntriesFast();
1800 if(label >= 0 && label < nprimaries)
1802 AliAODMCParticle * momP = (AliAODMCParticle *) mcparticles->At(label);
1804 grandMomLabel = momP->GetMother();
1806 AliAODMCParticle * grandmomP = 0x0;
1808 if(grandMomLabel >=0 )
1810 grandmomP = (AliAODMCParticle *) mcparticles->At(grandMomLabel);
1811 pdg = grandmomP->GetPdgCode();
1812 status = grandmomP->GetStatus();
1814 grandmom.SetPxPyPzE(grandmomP->Px(),grandmomP->Py(),grandmomP->Pz(),grandmomP->E());
1815 greatMomLabel = grandmomP->GetMother();
1826 //_______________________________________________________________________________________________________________
1827 void AliMCAnalysisUtils::GetMCDecayAsymmetryAngleForPDG(Int_t label, Int_t pdg, const AliCaloTrackReader* reader,
1828 Float_t & asym, Float_t & angle, Bool_t & ok)
1830 //In case of an eta or pi0 decay into 2 photons, get the asymmetry in the energy of the photons
1833 if(reader->ReadStack())
1835 if(!reader->GetStack())
1838 printf("AliMCAnalysisUtils::GetMCDecayAsymmetryForPDG() - Stack is not available, check analysis settings in configuration file, STOP!!\n");
1842 if(label >= 0 && label < reader->GetStack()->GetNtrack())
1844 TParticle * momP = reader->GetStack()->Particle(label);
1846 Int_t grandmomLabel = momP->GetFirstMother();
1847 Int_t grandmomPDG = -1;
1848 TParticle * grandmomP = 0x0;
1849 while (grandmomLabel >=0 )
1851 grandmomP = reader->GetStack()->Particle(grandmomLabel);
1852 grandmomPDG = grandmomP->GetPdgCode();
1854 if(grandmomPDG==pdg) break;
1856 grandmomLabel = grandmomP->GetFirstMother();
1860 if(grandmomPDG==pdg && grandmomP->GetNDaughters()==2)
1862 TParticle * d1 = reader->GetStack()->Particle(grandmomP->GetDaughter(0));
1863 TParticle * d2 = reader->GetStack()->Particle(grandmomP->GetDaughter(1));
1864 if(d1->GetPdgCode() == 22 && d1->GetPdgCode() == 22)
1866 asym = (d1->Energy()-d2->Energy())/grandmomP->Energy();
1867 TLorentzVector lvd1, lvd2;
1870 angle = lvd1.Angle(lvd2.Vect());
1876 printf("AliMCAnalysisUtils::GetMCDecayAsymmetryForPDG(ESD) - mother with PDG %d, not found! \n",pdg);
1881 else if(reader->ReadAODMCParticles())
1883 TClonesArray* mcparticles = reader->GetAODMCParticles();
1887 printf("AliMCAnalysisUtils::GetMCDecayAsymmetryForPDG() - AODMCParticles is not available, check analysis settings in configuration file!!\n");
1893 Int_t nprimaries = mcparticles->GetEntriesFast();
1894 if(label >= 0 && label < nprimaries)
1896 AliAODMCParticle * momP = (AliAODMCParticle *) mcparticles->At(label);
1898 Int_t grandmomLabel = momP->GetMother();
1899 Int_t grandmomPDG = -1;
1900 AliAODMCParticle * grandmomP = 0x0;
1901 while (grandmomLabel >=0 )
1903 grandmomP = (AliAODMCParticle *) mcparticles->At(grandmomLabel);
1904 grandmomPDG = grandmomP->GetPdgCode();
1906 if(grandmomPDG==pdg) break;
1908 grandmomLabel = grandmomP->GetMother();
1912 if(grandmomPDG==pdg && grandmomP->GetNDaughters()==2)
1914 AliAODMCParticle * d1 = (AliAODMCParticle *) mcparticles->At(grandmomP->GetDaughter(0));
1915 AliAODMCParticle * d2 = (AliAODMCParticle *) mcparticles->At(grandmomP->GetDaughter(1));
1916 if(d1->GetPdgCode() == 22 && d1->GetPdgCode() == 22)
1918 asym = (d1->E()-d2->E())/grandmomP->E();
1919 TLorentzVector lvd1, lvd2;
1920 lvd1.SetPxPyPzE(d1->Px(),d1->Py(),d1->Pz(),d1->E());
1921 lvd2.SetPxPyPzE(d2->Px(),d2->Py(),d2->Pz(),d2->E());
1922 angle = lvd1.Angle(lvd2.Vect());
1928 printf("AliMCAnalysisUtils::GetMCDecayAsymmetryForPDG(AOD) - mother with PDG %d, not found! \n",pdg);
1939 //_________________________________________________________________________________________________
1940 Int_t AliMCAnalysisUtils::GetNDaughters(Int_t label, const AliCaloTrackReader* reader, Bool_t & ok)
1942 // Return the the number of daughters of a given MC particle
1945 if(reader->ReadStack())
1947 if(!reader->GetStack())
1950 printf("AliMCAnalysisUtils::GetNDaughters() - Stack is not available, check analysis settings in configuration file, STOP!!\n");
1955 if(label >= 0 && label < reader->GetStack()->GetNtrack())
1957 TParticle * momP = reader->GetStack()->Particle(label);
1959 return momP->GetNDaughters();
1967 else if(reader->ReadAODMCParticles())
1969 TClonesArray* mcparticles = reader->GetAODMCParticles();
1973 printf("AliMCAnalysisUtils::GetNDaughters() - AODMCParticles is not available, check analysis settings in configuration file!!\n");
1979 Int_t nprimaries = mcparticles->GetEntriesFast();
1980 if(label >= 0 && label < nprimaries)
1982 AliAODMCParticle * momP = (AliAODMCParticle *) mcparticles->At(label);
1984 return momP->GetNDaughters();
1998 //_________________________________________________________________________________
1999 Int_t AliMCAnalysisUtils::GetNOverlaps(const Int_t * label, UInt_t nlabels,
2000 Int_t mctag, Int_t mesonLabel,
2001 AliCaloTrackReader * reader, Int_t *overpdg)
2003 // Compare the primary depositing more energy with the rest,
2004 // if no photon/electron (conversion) or neutral meson as comon ancestor, consider it as other particle contributing
2005 // Give as input the meson label in case it was a pi0 or eta merged cluster
2006 // Init overpdg with nlabels
2008 Int_t ancPDG = 0, ancStatus = -1;
2009 TLorentzVector momentum; TVector3 prodVertex;
2011 Int_t noverlaps = 0;
2014 for (UInt_t ilab = 1; ilab < nlabels; ilab++ )
2016 ancLabel = CheckCommonAncestor(label[0],label[ilab],reader,ancPDG,ancStatus,momentum,prodVertex);
2018 //printf("Overlaps, i %d: Main Label %d, second label %d, ancestor: Label %d, pdg %d - tag %d \n",
2019 // ilab,label[0],label[ilab],ancLabel,ancPDG, mctag);
2021 Bool_t overlap = kFALSE;
2026 //printf("\t \t \t No Label = %d\n",ancLabel);
2028 else if( ( ancPDG==111 || ancPDG==221 ) && ( CheckTagBit(mctag,kMCPi0) || CheckTagBit(mctag,kMCEta)) && mesonLabel != ancLabel)
2030 //printf("\t \t meson Label %d, ancestor Label %d\n",mesonLabel,ancLabel);
2033 else if( ancPDG!=22 && TMath::Abs(ancPDG)!=11 && ancPDG != 111 && ancPDG != 221 )
2035 //printf("\t \t \t Non EM PDG = %d\n",ancPDG);
2039 if( !overlap ) continue ;
2041 // We have at least one overlap
2043 //printf("Overlap!!!!!!!!!!!!!!\n");
2047 // What is the origin of the overlap?
2048 Bool_t mOK = 0, gOK = 0;
2049 Int_t mpdg = -999999, gpdg = -1;
2050 Int_t mstatus = -1, gstatus = -1;
2051 Int_t gLabel = -1, ggLabel = -1;
2052 TLorentzVector mother = GetMother (label[ilab],reader,mpdg,mstatus,mOK);
2053 TLorentzVector grandmother = GetGrandMother(label[ilab],reader,gpdg,gstatus,gOK, gLabel,ggLabel);
2055 //printf("\t Overlap!, mother pdg %d; grand mother pdg %d",mpdg,gpdg);
2057 if( ( mpdg == 22 || TMath::Abs(mpdg==11) ) &&
2058 ( gpdg == 22 || TMath::Abs(gpdg==11) ) &&
2061 Int_t labeltmp = gLabel;
2062 while( ( gpdg == 22 || TMath::Abs(gpdg==11) ) && gLabel >=0 )
2065 grandmother = GetGrandMother(labeltmp,reader,gpdg,gstatus,ok, gLabel,ggLabel);
2069 overpdg[noverlaps-1] = mpdg;
2076 //________________________________________________________
2077 void AliMCAnalysisUtils::Print(const Option_t * opt) const
2079 //Print some relevant parameters set for the analysis
2084 printf("***** Print: %s %s ******\n", GetName(), GetTitle() ) ;
2086 printf("Debug level = %d\n",fDebug);
2087 printf("MC Generator = %s\n",fMCGenerator.Data());
2092 //__________________________________________________
2093 void AliMCAnalysisUtils::PrintMCTag(Int_t tag) const
2095 // print the assigned origins to this particle
2097 printf("AliMCAnalysisUtils::PrintMCTag() - tag %d \n photon %d, conv %d, prompt %d, frag %d, isr %d, \n pi0 decay %d, eta decay %d, other decay %d pi0 %d, eta %d \n electron %d, muon %d,pion %d, proton %d, neutron %d, \n kaon %d, a-proton %d, a-neutron %d, unk %d, bad %d\n",
2099 CheckTagBit(tag,kMCPhoton),
2100 CheckTagBit(tag,kMCConversion),
2101 CheckTagBit(tag,kMCPrompt),
2102 CheckTagBit(tag,kMCFragmentation),
2103 CheckTagBit(tag,kMCISR),
2104 CheckTagBit(tag,kMCPi0Decay),
2105 CheckTagBit(tag,kMCEtaDecay),
2106 CheckTagBit(tag,kMCOtherDecay),
2107 CheckTagBit(tag,kMCPi0),
2108 CheckTagBit(tag,kMCEta),
2109 CheckTagBit(tag,kMCElectron),
2110 CheckTagBit(tag,kMCMuon),
2111 CheckTagBit(tag,kMCPion),
2112 CheckTagBit(tag,kMCProton),
2113 CheckTagBit(tag,kMCAntiNeutron),
2114 CheckTagBit(tag,kMCKaon),
2115 CheckTagBit(tag,kMCAntiProton),
2116 CheckTagBit(tag,kMCAntiNeutron),
2117 CheckTagBit(tag,kMCUnknown),
2118 CheckTagBit(tag,kMCBadLabel)