// --- ROOT system ---\r
#include <TMath.h>\r
#include <TList.h>\r
+#include "TParticle.h"\r
\r
//---- ANALYSIS system ----\r
#include "AliMCAnalysisUtils.h"\r
+#include "AliCaloTrackReader.h"\r
#include "AliStack.h"\r
-#include "TParticle.h"\r
#include "AliGenPythiaEventHeader.h"\r
+#include "AliAODMCParticle.h"\r
\r
ClassImp(AliMCAnalysisUtils)\r
\r
//________________________________________________\r
AliMCAnalysisUtils::AliMCAnalysisUtils() : \r
TObject(), fCurrentEvent(-1), fDebug(-1), \r
- fJetsList(new TList), fMCGenerator("PYTHIA"), fpTHardpTJetFactor(7)\r
+ fJetsList(new TList), fMCGenerator("PYTHIA")\r
{\r
//Ctor\r
}\r
//____________________________________________________________________________\r
AliMCAnalysisUtils::AliMCAnalysisUtils(const AliMCAnalysisUtils & mcutils) : \r
TObject(mcutils), fCurrentEvent(mcutils.fCurrentEvent), fDebug(mcutils.fDebug),\r
- fJetsList(mcutils.fJetsList), fMCGenerator(mcutils.fMCGenerator), fpTHardpTJetFactor(mcutils.fpTHardpTJetFactor)\r
+ fJetsList(mcutils.fJetsList), fMCGenerator(mcutils.fMCGenerator)\r
{\r
// cpy ctor\r
\r
fDebug = mcutils.fDebug;\r
fJetsList = mcutils.fJetsList;\r
fMCGenerator = mcutils.fMCGenerator;\r
- fpTHardpTJetFactor = mcutils.fpTHardpTJetFactor;\r
\r
return *this; \r
}\r
} \r
}\r
\r
+\r
//_________________________________________________________________________\r
-Int_t AliMCAnalysisUtils::CheckOrigin(const Int_t label, AliStack * stack) const {\r
- //Play with the MC stack if available\r
- //Check origin of the candidates, good for PYTHIA\r
- \r
+Int_t AliMCAnalysisUtils::CheckOrigin(const Int_t label, AliCaloTrackReader* reader, const Int_t input = 0) {\r
+ //Play with the montecarlo particles if available\r
+ Int_t tag = 0;\r
+ \r
+ //Select where the information is, ESD-galice stack or AOD mcparticles branch\r
+ if(reader->ReadStack()){\r
+ tag = CheckOriginInStack(label, reader->GetStack());\r
+ }\r
+ else if(reader->ReadAODMCParticles()){\r
+ tag = CheckOriginInAOD(label, reader->GetAODMCParticles(input));\r
+ }\r
+ \r
+ return tag ;\r
+} \r
+\r
+//_________________________________________________________________________\r
+Int_t AliMCAnalysisUtils::CheckOriginInStack(const Int_t label, AliStack* stack) {\r
+ // Play with the MC stack if available. Tag particles depending on their origin.\r
+ // Do same things as in CheckOriginInAOD but different input.\r
+ \r
if(!stack) {\r
- printf("AliMCAnalysisUtils::CheckOrigin() - Stack is not available, check analysis settings in configuration file, STOP!!\n");\r
+ printf("AliMCAnalysisUtils::CheckOriginInStack() - Stack is not available, check analysis settings in configuration file, STOP!!\n");\r
abort();\r
}\r
// printf("label %d, ntrack %d, nprim %d\n",label, stack->GetNtrack(), stack->GetNprimary());\r
// TParticle *particle = stack->Particle(i);\r
// //particle->Print();\r
// }\r
+ \r
+\r
+ Int_t tag = 0;\r
if(label >= 0 && label < stack->GetNtrack()){\r
//Mother\r
TParticle * mom = stack->Particle(label);\r
Int_t mPdg = TMath::Abs(mom->GetPdgCode());\r
Int_t mStatus = mom->GetStatusCode() ;\r
Int_t iParent = mom->GetFirstMother() ;\r
- if(fDebug > 0 && label < 8 ) printf("AliMCAnalysisUtils::CheckOrigin: Mother is parton %d\n",iParent);\r
+ if(fDebug > 0 && label < 8 ) printf("AliMCAnalysisUtils::CheckOriginInStack() - Mother is parton %d\n",iParent);\r
\r
//GrandParent\r
TParticle * parent = new TParticle ;\r
pPdg = TMath::Abs(parent->GetPdgCode());\r
pStatus = parent->GetStatusCode(); \r
}\r
- else if(fDebug > 0 ) printf("AliMCAnalysisUtils::CheckOrigin: Parent with label %d\n",iParent);\r
+ else if(fDebug > 0 ) printf("AliMCAnalysisUtils::CheckOriginInStack() - Parent with label %d\n",iParent);\r
\r
- //return tag\r
+ //Check if mother is converted, if not, get the first non converted mother\r
+ if((mPdg == 22 || mPdg == 11) && (pPdg == 22 || pPdg == 11) && mStatus == 0){\r
+ SetTagBit(tag,kMCConversion);\r
+ //Check if the mother is photon or electron with status not stable\r
+ while ((pPdg == 22 || pPdg == 11) && mStatus != 1) {\r
+ //Mother\r
+ mom = stack->Particle(mom->GetFirstMother());\r
+ mPdg = TMath::Abs(mom->GetPdgCode());\r
+ mStatus = mom->GetStatusCode() ;\r
+ iParent = mom->GetFirstMother() ;\r
+ if(fDebug > 0 && label < 8 ) printf("AliMCAnalysisUtils::CheckOriginInStack() - Mother is parton %d\n",iParent);\r
+ \r
+ //GrandParent\r
+ if(iParent > 0){\r
+ parent = stack->Particle(iParent);\r
+ pPdg = TMath::Abs(parent->GetPdgCode());\r
+ pStatus = parent->GetStatusCode(); \r
+ }\r
+ }//while \r
+ }//mother and parent are electron or photon and have status 0\r
+ \r
+ // conversion fo electrons/photons checked \r
+ \r
if(mPdg == 22){ //photon\r
+ SetTagBit(tag,kMCPhoton);\r
if(mStatus == 1){ //undecayed particle\r
if(fMCGenerator == "PYTHIA"){\r
if(iParent < 8 && iParent > 5) {//outgoing partons\r
- if(pPdg == 22) return kMCPrompt;\r
- else return kMCFragmentation;\r
+ if(pPdg == 22) SetTagBit(tag,kMCPrompt);\r
+ else SetTagBit(tag,kMCFragmentation);\r
}//Outgoing partons\r
else if(pStatus == 11){//Decay\r
- if(pPdg == 111) return kMCPi0Decay ;\r
- else if (pPdg == 221) return kMCEtaDecay ;\r
- else return kMCOtherDecay ;\r
+ if(pPdg == 111) SetTagBit(tag,kMCPi0Decay);\r
+ else if (pPdg == 221) SetTagBit(tag, kMCEtaDecay);\r
+ else SetTagBit(tag,kMCOtherDecay);\r
}//Decay\r
- else return kMCISR; //Initial state radiation\r
+ else SetTagBit(tag, kMCISR); //Initial state radiation\r
+ \r
}//PYTHIA\r
\r
else if(fMCGenerator == "HERWIG"){ \r
}//Look for the parton\r
\r
if(iParent < 8 && iParent > 5) {\r
- if(pPdg == 22) return kMCPrompt;\r
- else return kMCFragmentation;\r
+ if(pPdg == 22) SetTagBit(tag,kMCPrompt);\r
+ else SetTagBit(tag,kMCFragmentation);\r
}\r
- return kMCISR;//Initial state radiation\r
+ else SetTagBit(tag,kMCISR);//Initial state radiation\r
}//Not decay\r
else{//Decay\r
- if(pPdg == 111) return kMCPi0Decay ;\r
- else if (pPdg == 221) return kMCEtaDecay ;\r
- else return kMCOtherDecay ;\r
+ if(pPdg == 111) SetTagBit(tag,kMCPi0Decay); \r
+ else if (pPdg == 221) SetTagBit(tag,kMCEtaDecay);\r
+ else SetTagBit(tag,kMCOtherDecay);\r
}//Decay\r
+ \r
}//HERWIG\r
- else return kMCUnknown;\r
- }//Status 1 : Pythia generated\r
- else if(mStatus == 0){\r
- if(pPdg ==22 || pPdg ==11|| pPdg == 2112 || pPdg == 211 || \r
+ \r
+ else SetTagBit(tag,kMCUnknown);\r
+ \r
+ }//Status 1 : Pythia generated\r
+ \r
+ else if(mStatus == 0){\r
+ \r
+ if(pPdg ==22 || pPdg ==11 || pPdg == 2112 || pPdg == 211 || \r
pPdg == 321 || pPdg == 2212 || pPdg == 130 || pPdg == 13 ) \r
- return kMCConversion ;\r
- if(pPdg == 111) return kMCPi0Decay ;\r
- else if (pPdg == 221) return kMCEtaDecay ;\r
- else return kMCOtherDecay ;\r
+ SetTagBit(tag,kMCConversion);\r
+ \r
+ \r
+ if(pPdg == 111) SetTagBit(tag,kMCPi0Decay); \r
+ else if (pPdg == 221) SetTagBit(tag,kMCEtaDecay);\r
+ else SetTagBit(tag,kMCOtherDecay);\r
+ \r
}//status 0 : geant generated\r
+ \r
}//Mother Photon\r
- else if(mPdg == 111) return kMCPi0 ;\r
- else if(mPdg == 221) return kMCEta ;\r
+ else if(mPdg == 111) SetTagBit(tag,kMCPi0);\r
+ else if(mPdg == 221) SetTagBit(tag,kMCEta);\r
\r
//cluster's mother is an electron. Where did that electron come from?\r
else if(mPdg == 11){ //electron\r
-\r
- if(fDebug > 0) printf("AliMCAnalysisUtils::CheckOrigin: Checking ancestors of electrons");\r
+ SetTagBit(tag,kMCElectron); \r
+ if(fDebug > 0) printf("AliMCAnalysisUtils::CheckOriginInStack() - Checking ancestors of electrons");\r
\r
//check first for B and C ancestry, then other possibilities.\r
//An electron from a photon parent could have other particles in\r
//its history and we would want to know that, right?\r
\r
- if(mStatus == 1) { //electron from event generator\r
- if (pPdg == -1) return kMCElectron; //no parent\r
- else if (pPdg == 23) return kMCZDecay; //parent is Z-boson\r
- else if (pPdg == 24) return kMCWDecay; //parent is W-boson\r
+ if(mStatus == 1) { //electron from event generator\r
+ //if (pPdg == -1) return kMCElectron; //no parent\r
+ if (pPdg == 23) SetTagBit(tag,kMCZDecay); //parent is Z-boson\r
+ else if (pPdg == 24) SetTagBit(tag,kMCWDecay); //parent is W-boson\r
else { //check the electron's ancestors for B/C contribution\r
Bool_t bAncestor = kFALSE;\r
Bool_t cAncestor = kFALSE;\r
Int_t aPdg = TMath::Abs(ancestors->GetPdgCode());\r
//Int_t aStatus = ancestors->GetStatusCode();\r
Int_t iAncestors = ancestors->GetFirstMother();\r
- if(fDebug > 0) printf("AliMCAnalysisUtils::CheckOrigin: Scaning the decay chain for bottom/charm generated electron");\r
+ if(fDebug > 0) printf("AliMCAnalysisUtils::CheckOriginInStack() - Scaning the decay chain for bottom/charm generated electron");\r
while(ancestors->IsPrimary()){//searching for ancestors \r
if((499 < aPdg && aPdg < 600)||(4999 < aPdg && aPdg < 6000)) bAncestor = kTRUE;\r
if((399 < aPdg && aPdg < 500)||(3999 < aPdg && aPdg < 5000)) cAncestor = kTRUE;\r
ancestors = stack->Particle(iAncestors);\r
aPdg = ancestors->GetPdgCode();\r
}//searching for ancestors\r
- if(bAncestor && cAncestor) return kMCEFromCFromB;//Decay chain has both B and C\r
- else if(bAncestor && !cAncestor) return kMCEFromB;//Decay chain has only B\r
- else if(!bAncestor && cAncestor) return kMCEFromC;//Decay chain has only C \r
+ if(bAncestor && cAncestor) SetTagBit(tag,kMCEFromCFromB);//Decay chain has both B and C\r
+ else if(bAncestor && !cAncestor) SetTagBit(tag,kMCEFromB);//Decay chain has only B\r
+ else if(!bAncestor && cAncestor) SetTagBit(tag,kMCEFromC);//Decay chain has only C \r
}\r
//if it is not from W,Z or B/C ancestor, where is it from?\r
- if (pPdg == 111) return kMCPi0Decay;//Pi0 Dalitz decay\r
- else if(pPdg == 221) return kMCEtaDecay;//Eta Dalitz decay\r
- else return kMCOtherDecay;\r
+ if (pPdg == 111) SetTagBit(tag,kMCPi0Decay);//Pi0 Dalitz decay\r
+ else if(pPdg == 221) SetTagBit(tag,kMCEtaDecay);//Eta Dalitz decay\r
+ else SetTagBit(tag,kMCOtherDecay);\r
\r
- } else if (mStatus == 0) { //electron from GEANT\r
+ } \r
+ else if (mStatus == 0) { //electron from GEANT\r
\r
//Rewind ancestry and check for electron with status == 1\r
//if we find one, we'll assume that this object is from an\r
Int_t aPdg = TMath::Abs(ancestors->GetPdgCode());\r
Int_t aStatus = ancestors->GetStatusCode();\r
Int_t iAncestors = ancestors->GetFirstMother();\r
- if(fDebug > 0) printf("AliMCAnalysisUtils::CheckOrigin: Scaning the decay chain for bottom/charm electrons");\r
+ if(fDebug > 0) printf("AliMCAnalysisUtils::CheckOriginInStack() - Scaning the decay chain for bottom/charm electrons");\r
while(ancestors->IsPrimary()){//searching for ancestors\r
if(aStatus == 1 && aPdg == 11) eleFromEvGen = kTRUE;\r
- if(eleFromEvGen && aPdg == 23) return kMCZDecay;\r
- if(eleFromEvGen && aPdg == 24) return kMCWDecay;\r
+ if(eleFromEvGen && aPdg == 23) SetTagBit(tag,kMCZDecay);\r
+ if(eleFromEvGen && aPdg == 24) SetTagBit(tag,kMCWDecay);\r
if(eleFromEvGen && ((499 < aPdg && aPdg < 600)||(4999 < aPdg && aPdg < 6000))) bAncestor = kTRUE;\r
if(eleFromEvGen && ((399 < aPdg && aPdg < 500)||(3999 < aPdg && aPdg < 5000))) cAncestor = kTRUE;\r
if(bAncestor && cAncestor) break;\r
ancestors = stack->Particle(iAncestors);\r
aPdg = ancestors->GetPdgCode();\r
}//searching for ancestors\r
- if(bAncestor && cAncestor) return kMCEFromCFromB;//Decay chain has both B and C\r
- else if(bAncestor && !cAncestor) return kMCEFromB;//Decay chain has only B\r
- else if(!bAncestor && cAncestor) return kMCEFromC;//Decay chain has only C\r
+ if(bAncestor && cAncestor) SetTagBit(tag,kMCEFromCFromB);//Decay chain has both B and C\r
+ else if(bAncestor && !cAncestor) SetTagBit(tag,kMCEFromB);//Decay chain has only B\r
+ else if(!bAncestor && cAncestor) SetTagBit(tag,kMCEFromC);//Decay chain has only C\r
if(pPdg ==22 || pPdg ==11|| pPdg == 2112 || pPdg == 211 || \r
pPdg == 321 || pPdg == 2212 || pPdg == 130 || pPdg == 13 ) \r
- return kMCConversion ;\r
- if(pPdg == 111) return kMCPi0Decay ;\r
- else if (pPdg == 221) return kMCEtaDecay ;\r
- else return kMCOtherDecay ;\r
+ SetTagBit(tag,kMCConversion);\r
+ if(pPdg == 111) SetTagBit(tag,kMCPi0Decay);\r
+ else if (pPdg == 221) SetTagBit(tag,kMCEtaDecay);\r
+ else SetTagBit(tag,kMCOtherDecay);\r
} //GEANT check\r
}//electron check\r
- else return kMCUnknown;\r
+ else SetTagBit(tag,kMCUnknown);\r
}//Good label value\r
else{\r
- if(label < 0 ) printf("AliMCAnalysisUtils::CheckOrigin: *** bad label or no stack ***: label %d \n", label);\r
- if(label >= stack->GetNtrack()) printf("AliMCAnalysisUtils::CheckOrigin: *** large label ***: label %d, n tracks %d \n", label, stack->GetNtrack());\r
- return kMCUnknown;\r
+ if(label < 0 ) printf("AliMCAnalysisUtils::CheckOriginInStack() *** bad label or no stack ***: label %d \n", label);\r
+ if(label >= stack->GetNtrack()) printf("AliMCAnalysisUtils::CheckOriginInStack() *** large label ***: label %d, n tracks %d \n", label, stack->GetNtrack());\r
+ SetTagBit(tag,kMCUnknown);\r
}//Bad label\r
\r
- return kMCUnknown;\r
+ return tag;\r
\r
}\r
\r
+\r
//_________________________________________________________________________\r
-TList * AliMCAnalysisUtils::GetJets(const Int_t iEvent, AliStack * stack, const AliGenEventHeader * geh) {\r
- //Return list of jets (TParticles) and index of most likely parton that originated it.\r
+Int_t AliMCAnalysisUtils::CheckOriginInAOD(const Int_t label, TClonesArray *mcparticles) {\r
+ // Play with the MCParticles in AOD if available. Tag particles depending on their origin.\r
+ // Do same things as in CheckOriginInStack but different input.\r
+ if(!mcparticles) {\r
+ printf("AliMCAnalysisUtils::CheckOriginInAOD() - AODMCParticles is not available, check analysis settings in configuration file!!\n");\r
+ \r
+ }\r
+\r
+ // printf(">>>>>>>>> Second Event, AODMCParticles, nMC %d\n", (GetAODMCParticles(1))->GetEntriesFast());\r
+ // for(Int_t i = 0; i< (GetAODMCParticles(1))->GetEntriesFast(); i++){\r
+ // AliAODMCParticle *particle = (AliAODMCParticle*) (GetAODMCParticles(1))->At(i);\r
+ // printf("** index %d, PDG %d, Flag %d, Mother %d, Daughter %d and %d, E %2.3f, pT %2.3f, phi %2.3f, eta %2.3f\n", \r
+ // i, particle->GetPdgCode(), particle->GetFlag(), particle->GetMother(), particle->GetDaughter(0),particle->GetDaughter(1),\r
+ // particle->E(),particle->Pt(), particle->Phi()*TMath::RadToDeg(), particle->Eta());\r
+ // if(particle->IsPhysicalPrimary()) printf("Is Physical Primary !\n");\r
+ // if(!particle->IsPrimary()) printf("Not Primary !\n");\r
+ // }\r
+ \r
\r
+ Int_t tag = 0;\r
+ Int_t nprimaries = mcparticles->GetEntriesFast();\r
+ if(label >= 0 && label < nprimaries ){\r
+ //Mother\r
+ AliAODMCParticle * mom = (AliAODMCParticle *) mcparticles->At(label);\r
+ Int_t mPdg = TMath::Abs(mom->GetPdgCode());\r
+ Int_t iParent = mom->GetMother() ;\r
+ if(fDebug > 0 && label < 8 ) printf("AliMCAnalysisUtils::CheckOriginInAOD() - Mother is parton %d\n",iParent);\r
+ \r
+ //GrandParent\r
+ AliAODMCParticle * parent = new AliAODMCParticle ;\r
+ Int_t pPdg = -1;\r
+ if(iParent >= 0){\r
+ parent = (AliAODMCParticle *) mcparticles->At(iParent);\r
+ pPdg = TMath::Abs(parent->GetPdgCode());\r
+ }\r
+ else if(fDebug > 0 ) printf("AliMCAnalysisUtils::CheckOriginInAOD() - Parent with label %d\n",iParent);\r
+ \r
+ //Check if mother is converted, if not, get the first non converted mother\r
+ if((mPdg == 22 || mPdg == 11) && (pPdg == 22 || pPdg == 11) && !mom->IsPrimary()){\r
+ SetTagBit(tag,kMCConversion);\r
+ //Check if the mother is photon or electron with status not stable\r
+ while ((pPdg == 22 || pPdg == 11) && !mom->IsPhysicalPrimary()) {\r
+ //Mother\r
+ mom = (AliAODMCParticle *) mcparticles->At(mom->GetMother());\r
+ mPdg = TMath::Abs(mom->GetPdgCode());\r
+ iParent = mom->GetMother() ;\r
+ if(fDebug > 0 && label < 8 ) printf("AliMCAnalysisUtils::CheckOriginInAOD() - Mother is parton %d\n",iParent);\r
+ \r
+ //GrandParent\r
+ if(iParent >= 0){\r
+ parent = (AliAODMCParticle *) mcparticles->At(iParent);\r
+ pPdg = TMath::Abs(parent->GetPdgCode());\r
+ }\r
+ }//while \r
+ }//mother and parent are electron or photon and have status 0\r
+ \r
+ // conversion fo electrons/photons checked \r
+ \r
+ if(mPdg == 22){ //photon\r
+ SetTagBit(tag,kMCPhoton);\r
+ if(mom->IsPhysicalPrimary()){ //undecayed particle\r
+ if(iParent < 8 && iParent > 5) {//outgoing partons\r
+ if(pPdg == 22) SetTagBit(tag,kMCPrompt);\r
+ else SetTagBit(tag,kMCFragmentation);\r
+ }//Outgoing partons\r
+ else if(parent->IsPrimary() && !parent->IsPhysicalPrimary()){//Decay\r
+ if(pPdg == 111) SetTagBit(tag,kMCPi0Decay);\r
+ else if (pPdg == 221) SetTagBit(tag, kMCEtaDecay);\r
+ else SetTagBit(tag,kMCOtherDecay);\r
+ }//Decay\r
+ else SetTagBit(tag, kMCISR); //Initial state radiation\r
+ }// Pythia generated\r
+ else if(!mom->IsPrimary()){\r
+ \r
+ if(pPdg ==22 || pPdg ==11 || pPdg == 2112 || pPdg == 211 || \r
+ pPdg == 321 || pPdg == 2212 || pPdg == 130 || pPdg == 13 ) \r
+ SetTagBit(tag,kMCConversion);\r
+ \r
+ if(pPdg == 111) SetTagBit(tag,kMCPi0Decay); \r
+ else if (pPdg == 221) SetTagBit(tag,kMCEtaDecay);\r
+ else SetTagBit(tag,kMCOtherDecay);\r
+ \r
+ }//not primary : geant generated\r
+ else {\r
+// printf("UNKNOWN 1, mom pdg %d, primary %d, physical primary %d; parent %d, pdg %d, primary %d, physical primary %d \n",\r
+// mPdg, mom->IsPrimary(), mom->IsPhysicalPrimary(), iParent, pPdg, parent->IsPrimary(), parent->IsPhysicalPrimary());\r
+ SetTagBit(tag,kMCUnknown);\r
+ }\r
+ }//Mother Photon\r
+ else if(mPdg == 111) SetTagBit(tag,kMCPi0);\r
+ else if(mPdg == 221) SetTagBit(tag,kMCEta);\r
+ \r
+ //cluster's mother is an electron. Where did that electron come from?\r
+ else if(mPdg == 11){ //electron\r
+ SetTagBit(tag,kMCElectron); \r
+ if(fDebug > 0) printf("AliMCAnalysisUtils::CheckOriginInAOD() - Checking ancestors of electrons");\r
+ \r
+ //check first for B and C ancestry, then other possibilities.\r
+ //An electron from a photon parent could have other particles in\r
+ //its history and we would want to know that, right?\r
+ \r
+ if(mom->IsPhysicalPrimary()) { //electron from event generator\r
+ //if (pPdg == -1) return kMCElectron; //no parent\r
+ if (pPdg == 23) SetTagBit(tag,kMCZDecay); //parent is Z-boson\r
+ else if (pPdg == 24) SetTagBit(tag,kMCWDecay); //parent is W-boson\r
+ else { //check the electron's ancestors for B/C contribution\r
+ Bool_t bAncestor = kFALSE;\r
+ Bool_t cAncestor = kFALSE;\r
+ AliAODMCParticle * ancestors = (AliAODMCParticle *) mcparticles->At(label);\r
+ Int_t aPdg = TMath::Abs(ancestors->GetPdgCode());\r
+ //Int_t aStatus = ancestors->GetFlag();\r
+ Int_t iAncestors = ancestors->GetMother();\r
+ if(fDebug > 0) printf("AliMCAnalysisUtils::CheckOriginInAOD() - Scaning the decay chain for bottom/charm generated electron");\r
+ while(ancestors->IsPrimary()){//searching for ancestors \r
+ if((499 < aPdg && aPdg < 600)||(4999 < aPdg && aPdg < 6000)) bAncestor = kTRUE;\r
+ if((399 < aPdg && aPdg < 500)||(3999 < aPdg && aPdg < 5000)) cAncestor = kTRUE;\r
+ if(bAncestor && cAncestor) break;\r
+ iAncestors = ancestors->GetMother();\r
+ if(iAncestors ==-1 ) break; // check what happens here.\r
+ ancestors = (AliAODMCParticle *) mcparticles->At(iAncestors);\r
+ aPdg = ancestors->GetPdgCode();\r
+ }//searching for ancestors\r
+ if(bAncestor && cAncestor) SetTagBit(tag,kMCEFromCFromB);//Decay chain has both B and C\r
+ else if(bAncestor && !cAncestor) SetTagBit(tag,kMCEFromB);//Decay chain has only B\r
+ else if(!bAncestor && cAncestor) SetTagBit(tag,kMCEFromC);//Decay chain has only C \r
+ }\r
+ //if it is not from W,Z or B/C ancestor, where is it from?\r
+ if (pPdg == 111) SetTagBit(tag,kMCPi0Decay);//Pi0 Dalitz decay\r
+ else if(pPdg == 221) SetTagBit(tag,kMCEtaDecay);//Eta Dalitz decay\r
+ else SetTagBit(tag,kMCOtherDecay);\r
+ \r
+ } \r
+ else if (!mom->IsPrimary()) { //electron from GEANT\r
+ \r
+ //Rewind ancestry and check for electron with status == 1\r
+ //if we find one, we'll assume that this object is from an\r
+ //electron but that it may have gone through some showering in\r
+ //material before the detector\r
+ \r
+ //Not a double-counting problem because we are only accessing\r
+ //these histories for MC labels connected to a reco object.\r
+ //If you wanted to use this to sort through the kine stack\r
+ //directly, might it be a problem?\r
+ Bool_t eleFromEvGen = kFALSE;\r
+ Bool_t bAncestor = kFALSE;\r
+ Bool_t cAncestor = kFALSE;\r
+ \r
+ AliAODMCParticle * ancestors = (AliAODMCParticle *) mcparticles->At(label);\r
+ Int_t aPdg = TMath::Abs(ancestors->GetPdgCode());\r
+ //Int_t aStatus = ancestors->GetFlag();\r
+ Int_t iAncestors = ancestors->GetMother();\r
+ if(fDebug > 0) printf("AliMCAnalysisUtils::CheckOriginInAOD() - Scaning the decay chain for bottom/charm electrons");\r
+ while(ancestors->IsPrimary()){//searching for ancestors\r
+ if(ancestors->IsPhysicalPrimary() && aPdg == 11) eleFromEvGen = kTRUE;\r
+ if(eleFromEvGen && aPdg == 23) SetTagBit(tag,kMCZDecay);\r
+ if(eleFromEvGen && aPdg == 24) SetTagBit(tag,kMCWDecay);\r
+ if(eleFromEvGen && ((499 < aPdg && aPdg < 600)||(4999 < aPdg && aPdg < 6000))) bAncestor = kTRUE;\r
+ if(eleFromEvGen && ((399 < aPdg && aPdg < 500)||(3999 < aPdg && aPdg < 5000))) cAncestor = kTRUE;\r
+ if(bAncestor && cAncestor) break;\r
+ iAncestors = ancestors->GetMother();\r
+ if(iAncestors ==-1 ) break; // check what happens here.\r
+ ancestors = (AliAODMCParticle *) mcparticles->At(iAncestors);\r
+ aPdg = ancestors->GetPdgCode();\r
+ }//searching for ancestors\r
+ if(bAncestor && cAncestor) SetTagBit(tag,kMCEFromCFromB);//Decay chain has both B and C\r
+ else if(bAncestor && !cAncestor) SetTagBit(tag,kMCEFromB);//Decay chain has only B\r
+ else if(!bAncestor && cAncestor) SetTagBit(tag,kMCEFromC);//Decay chain has only C\r
+ if(pPdg ==22 || pPdg ==11|| pPdg == 2112 || pPdg == 211 || \r
+ pPdg == 321 || pPdg == 2212 || pPdg == 130 || pPdg == 13 ) \r
+ SetTagBit(tag,kMCConversion);\r
+ if(pPdg == 111) SetTagBit(tag,kMCPi0Decay);\r
+ else if (pPdg == 221) SetTagBit(tag,kMCEtaDecay);\r
+ else SetTagBit(tag,kMCOtherDecay);\r
+ } //GEANT check\r
+ }//electron check\r
+ else {\r
+// printf("UNKNOWN 2, mom pdg %d, primary %d, physical primary %d; parent %d, pdg %d, primary %d, physical primary %d \n",\r
+// mPdg, mom->IsPrimary(), mom->IsPhysicalPrimary(), iParent, pPdg, parent->IsPrimary(), parent->IsPhysicalPrimary());\r
+ SetTagBit(tag,kMCUnknown);\r
+ }\r
+ }//Good label value\r
+ else{\r
+ if(label < 0 ) printf("AliMCAnalysisUtils::CheckOriginInAOD() *** bad label or no stack ***: label %d \n", label);\r
+ if(label >= mcparticles->GetEntriesFast()) printf("AliMCAnalysisUtils::CheckOriginInAOD() *** large label ***: label %d, n tracks %d \n", label, mcparticles->GetEntriesFast());\r
+ SetTagBit(tag,kMCUnknown);\r
+ }//Bad label\r
+ \r
+ return tag;\r
+ \r
+}\r
+\r
+\r
+\r
+//_________________________________________________________________________\r
+TList * AliMCAnalysisUtils::GetJets(AliCaloTrackReader * reader){\r
+ //Return list of jets (TParticles) and index of most likely parton that originated it.\r
+ AliStack * stack = reader->GetStack();\r
+ Int_t iEvent = reader->GetEventNumber(); \r
+ AliGenEventHeader * geh = reader->GetGenEventHeader();\r
if(fCurrentEvent!=iEvent){\r
fCurrentEvent = iEvent;\r
fJetsList = new TList;\r
}\r
\r
\r
-//_________________________________________________________________________\r
-Bool_t AliMCAnalysisUtils::ComparePtHardAndJetPt(const AliGenEventHeader * geh){\r
- // Check the event, if the requested ptHard is much larger than the jet pT, then there is a problem.\r
- // Only for PYTHIA.\r
- \r
- if(fMCGenerator == "PYTHIA"){\r
- TParticle * jet = new TParticle;\r
- AliGenPythiaEventHeader* pygeh= (AliGenPythiaEventHeader*) geh;\r
- Int_t nTriggerJets = pygeh->NTriggerJets();\r
- Float_t ptHard = pygeh->GetPtHard();\r
-\r
- //if(fDebug > 1) printf("AliMCAnalysisUtils::PythiaEventHeader: Njets: %d, pT Hard %f\n",nTriggerJets, ptHard);\r
- Float_t tmpjet[]={0,0,0,0};\r
- for(Int_t ijet = 0; ijet< nTriggerJets; ijet++){\r
- pygeh->TriggerJet(ijet, tmpjet);\r
- jet = new TParticle(94, 21, -1, -1, -1, -1, tmpjet[0],tmpjet[1],tmpjet[2],tmpjet[3], 0,0,0,0);\r
- //Compare jet pT and pt Hard\r
- //if(fDebug > 1) printf("AliMCAnalysisUtils:: %d pycell jet pT %f\n",ijet, jet->Pt());\r
- if(jet->Pt() > fpTHardpTJetFactor * ptHard) {\r
- printf("AliMCAnalysisUtils::PythiaEventHeader: Njets: %d, pT Hard %2.2f, pycell jet pT %2.2f, rejection factor %1.1f\n",\r
- nTriggerJets, ptHard, jet->Pt(), fpTHardpTJetFactor);\r
- return kFALSE;\r
- }\r
- }\r
- } \r
- \r
- return kTRUE ;\r
-\r
-}\r
-\r
//________________________________________________________________\r
void AliMCAnalysisUtils::Print(const Option_t * opt) const\r
{\r
\r
printf("Debug level = %d\n",fDebug);\r
printf("MC Generator = %s\n",fMCGenerator.Data());\r
- printf("fpTHardpTJetFactor = %2.2f",fpTHardpTJetFactor);\r
printf(" \n");\r
\r
} \r
git://git.uio.no / u / mrichter / AliRoot.git / blobdiff