fCurrentEvent(-1),
fDebug(-1),
fJetsList(new TList),
-fMCGenerator("PYTHIA")
+fMCGenerator(kPythia),
+fMCGeneratorString("PYTHIA"),
+fDaughMom(), fDaughMom2(),
+fMotherMom(), fGMotherMom()
{
//Ctor
}
{
// Remove all pointers.
- if (fJetsList) {
+ if (fJetsList)
+ {
fJetsList->Clear();
delete fJetsList ;
}
}
//_____________________________________________________________________________________________
-Int_t AliMCAnalysisUtils::CheckCommonAncestor(const Int_t index1, const Int_t index2,
+Int_t AliMCAnalysisUtils::CheckCommonAncestor(Int_t index1, Int_t index2,
const AliCaloTrackReader* reader,
Int_t & ancPDG, Int_t & ancStatus,
TLorentzVector & momentum, TVector3 & prodVertex)
{
- //Check the first common ancestor of 2 clusters, given the most likely labels of the primaries generating such clusters.
+ // Check the first common ancestor of 2 clusters, given the most likely labels of the primaries generating such clusters.
+
Int_t label1[100];
Int_t label2[100];
label1[0]= index1;
return ancLabel;
}
-//_____________________________________________________________________
-Int_t AliMCAnalysisUtils::CheckOrigin(const Int_t * label,
- const Int_t nlabels,
- const AliCaloTrackReader* reader)
+//________________________________________________________________________________________
+Int_t AliMCAnalysisUtils::CheckOrigin(const Int_t * label, Int_t nlabels,
+ const AliCaloTrackReader* reader, TString calorimeter)
{
- //Play with the montecarlo particles if available
+ // Play with the montecarlo particles if available.
+
Int_t tag = 0;
if(nlabels<=0) {
return kMCBadLabel;
}
+ TObjArray* arrayCluster = 0;
+ if ( calorimeter == "EMCAL" ) arrayCluster = reader->GetEMCALClusters();
+ else if ( calorimeter == "PHOS" ) arrayCluster= reader->GetPHOSClusters();
+
//Select where the information is, ESD-galice stack or AOD mcparticles branch
if(reader->ReadStack()){
- tag = CheckOriginInStack(label, nlabels, reader->GetStack());
+ tag = CheckOriginInStack(label, nlabels, reader->GetStack(), arrayCluster);
}
else if(reader->ReadAODMCParticles()){
- tag = CheckOriginInAOD(label, nlabels, reader->GetAODMCParticles());
+ tag = CheckOriginInAOD(label, nlabels, reader->GetAODMCParticles(),arrayCluster);
}
return tag ;
}
-//_____________________________________________________________________
-Int_t AliMCAnalysisUtils::CheckOrigin(const Int_t label,
- const AliCaloTrackReader* reader)
+//____________________________________________________________________________________________________
+Int_t AliMCAnalysisUtils::CheckOrigin(Int_t label, const AliCaloTrackReader* reader, TString calorimeter)
{
- //Play with the montecarlo particles if available
+ // Play with the montecarlo particles if available.
+
Int_t tag = 0;
if(label<0) {
return kMCBadLabel;
}
+ TObjArray* arrayCluster = 0;
+ if ( calorimeter == "EMCAL" ) arrayCluster = reader->GetEMCALClusters();
+ else if ( calorimeter == "PHOS" ) arrayCluster= reader->GetPHOSClusters();
+
Int_t labels[]={label};
//Select where the information is, ESD-galice stack or AOD mcparticles branch
if(reader->ReadStack()){
- tag = CheckOriginInStack(labels, 1,reader->GetStack());
+ tag = CheckOriginInStack(labels, 1,reader->GetStack(),arrayCluster);
}
else if(reader->ReadAODMCParticles()){
- tag = CheckOriginInAOD(labels, 1,reader->GetAODMCParticles());
+ tag = CheckOriginInAOD(labels, 1,reader->GetAODMCParticles(),arrayCluster);
}
return tag ;
}
-//_________________________________________________________________
-Int_t AliMCAnalysisUtils::CheckOriginInStack(const Int_t *labels,
- const Int_t nlabels,
- AliStack* stack)
+//__________________________________________________________________________________________
+Int_t AliMCAnalysisUtils::CheckOriginInStack(const Int_t *labels, Int_t nlabels,
+ AliStack* stack, const TObjArray* arrayCluster)
{
// Play with the MC stack if available. Tag particles depending on their origin.
// Do same things as in CheckOriginInAOD but different input.
Int_t mStatus = mom->GetStatusCode() ;
Int_t iParent = mom->GetFirstMother() ;
- if(fDebug > 0 && label < 8 && fMCGenerator!="") printf("AliMCAnalysisUtils::CheckOriginInStack() - Mother is parton %d\n",iParent);
+ if(fDebug > 0 && label < 8 && fMCGenerator != kBoxLike) printf("AliMCAnalysisUtils::CheckOriginInStack() - Mother is parton %d\n",iParent);
//GrandParent of the entity
TParticle * parent = NULL;
printf("AliMCAnalysisUtils::CheckOriginInStack() - First mother is directly pi0, not decayed by generator \n");
CheckOverlapped2GammaDecay(labels,nlabels, iMom, stack, tag); //set to kMCPi0 if 2 gammas in same cluster
+
}
else if(mPdg == 221)
{
if(fDebug > 2 ) printf("AliMCAnalysisUtils::CheckOriginInStack() - PYTHIA pi0 decay photon, parent pi0 with status 11 \n");
CheckOverlapped2GammaDecay(labels,nlabels, iParent, stack, tag); //set to kMCPi0 if 2 gammas in same cluster
+ // In case it did not merge, check if the decay companion is lost
+ if(!CheckTagBit(tag, kMCPi0) && !CheckTagBit(tag,kMCDecayPairInCalo))
+ CheckLostDecayPair(arrayCluster,iMom, iParent, stack, tag);
+
+ //printf("Bit set is Merged %d, Pair in calo %d, Lost %d\n",CheckTagBit(tag, kMCPi0),CheckTagBit(tag,kMCDecayPairInCalo),CheckTagBit(tag,kMCDecayPairLost));
}
else if (pPdg == 221)
{
if(fDebug > 2 ) printf("AliMCAnalysisUtils::CheckOriginInStack() - PYTHIA eta decay photon, parent pi0 with status 11 \n");
CheckOverlapped2GammaDecay(labels,nlabels, iParent, stack, tag);//set to kMCEta if 2 gammas in same cluster
+ // In case it did not merge, check if the decay companion is lost
+ if(!CheckTagBit(tag, kMCEta) && !CheckTagBit(tag,kMCDecayPairInCalo))
+ CheckLostDecayPair(arrayCluster,iMom, iParent, stack, tag);
}
else if(mStatus == 1)
{ //undecayed particle
- if(fMCGenerator == "PYTHIA")
+ if(fMCGenerator == kPythia)
{
if(iParent < 8 && iParent > 5)
{//outgoing partons
else SetTagBit(tag,kMCUnknown);
}//PYTHIA
- else if(fMCGenerator == "HERWIG")
+ else if(fMCGenerator == kHerwig)
{
if(pStatus < 197)
{//Not decay
if(fDebug > 0) printf("AliMCAnalysisUtils::CheckOriginInStack() - Checking ancestors of electrons\n");
- if (pPdg == 111) { SetTagBit(tag,kMCPi0Decay); } //Pi0 Dalitz decay
- else if (pPdg == 221) { SetTagBit(tag,kMCEtaDecay); } //Eta Dalitz decay
+ if (pPdg == 111) { SetTagBit(tag,kMCPi0Decay); SetTagBit(tag, kMCDecayDalitz) ; } //Pi0 Dalitz decay
+ else if (pPdg == 221) { SetTagBit(tag,kMCEtaDecay); SetTagBit(tag, kMCDecayDalitz) ; } //Eta Dalitz decay
else if((499 < pPdg && pPdg < 600)||(4999 < pPdg && pPdg < 6000)) { SetTagBit(tag,kMCEFromB); } //b-->e decay
else if((399 < pPdg && pPdg < 500)||(3999 < pPdg && pPdg < 5000)) { //check charm decay
if(parent)
}
-//_________________________________________________________________________
-Int_t AliMCAnalysisUtils::CheckOriginInAOD(const Int_t *labels,
- const Int_t nlabels,
- const TClonesArray *mcparticles)
+//________________________________________________________________________________________________________
+Int_t AliMCAnalysisUtils::CheckOriginInAOD(const Int_t *labels, Int_t nlabels,
+ const TClonesArray *mcparticles, const TObjArray* arrayCluster)
{
// Play with the MCParticles in AOD if available. Tag particles depending on their origin.
// Do same things as in CheckOriginInStack but different input.
+
if(!mcparticles)
{
if(fDebug >= 0)
Int_t mPdg = TMath::Abs(mPdgSign);
Int_t iParent = mom->GetMother() ;
- if(fDebug > 0 && label < 8 && fMCGenerator!="") printf("AliMCAnalysisUtils::CheckOriginInAOD() - Mother is parton %d\n",iParent);
+ if(fDebug > 0 && label < 8 && fMCGenerator != kBoxLike) printf("AliMCAnalysisUtils::CheckOriginInAOD() - Mother is parton %d\n",iParent);
//GrandParent
AliAODMCParticle * parent = NULL ;
if(fDebug > 2 ) printf("AliMCAnalysisUtils::CheckOriginInAOD() - Generator pi0 decay photon \n");
CheckOverlapped2GammaDecay(labels,nlabels, iParent, mcparticles, tag); //set to kMCPi0 if 2 gammas in same cluster
+ // In case it did not merge, check if the decay companion is lost
+ if(!CheckTagBit(tag, kMCPi0) && !CheckTagBit(tag,kMCDecayPairInCalo) && !CheckTagBit(tag,kMCDecayPairLost))
+ CheckLostDecayPair(arrayCluster,iMom, iParent, mcparticles, tag);
+
+ //printf("Bit set is Merged %d, Pair in calo %d, Lost %d\n",CheckTagBit(tag, kMCPi0),CheckTagBit(tag,kMCDecayPairInCalo),CheckTagBit(tag,kMCDecayPairLost));
}
else if (pPdg == 221)
{
if(fDebug > 2 ) printf("AliMCAnalysisUtils::CheckOriginInAOD() - Generator eta decay photon \n");
CheckOverlapped2GammaDecay(labels,nlabels, iParent, mcparticles, tag); //set to kMCEta if 2 gammas in same cluster
+ // In case it did not merge, check if the decay companion is lost
+ if(!CheckTagBit(tag, kMCEta) && !CheckTagBit(tag,kMCDecayPairInCalo))
+ CheckLostDecayPair(arrayCluster,iMom, iParent, mcparticles, tag);
}
- else if(mom->IsPhysicalPrimary() && (fMCGenerator=="PYTHIA" || fMCGenerator=="HERWIG")) //undecayed particle
+ else if( mom->IsPhysicalPrimary() && ( fMCGenerator == kPythia || fMCGenerator == kHerwig ) ) //undecayed particle
{
if(iParent < 8 && iParent > 5 )
{//outgoing partons
if(pPdg == 22) SetTagBit(tag,kMCPrompt);
else SetTagBit(tag,kMCFragmentation);
}//Outgoing partons
- else if(iParent <= 5 && (fMCGenerator=="PYTHIA" || fMCGenerator=="HERWIG"))
+ else if( iParent <= 5 && ( fMCGenerator == kPythia || fMCGenerator == kHerwig ) )
{
SetTagBit(tag, kMCISR); //Initial state radiation
}
SetTagBit(tag,kMCElectron);
- if (fDebug > 0) printf("AliMCAnalysisUtils::CheckOriginInAOD() - Checking ancestors of electrons");
- if (pPdg == 111) { SetTagBit(tag,kMCPi0Decay); } //Pi0 Dalitz decay
- else if (pPdg == 221) { SetTagBit(tag,kMCEtaDecay); } //Eta Dalitz decay
+ if(fDebug > 0) printf("AliMCAnalysisUtils::CheckOriginInAOD() - Checking ancestors of electrons");
+
+ if (pPdg == 111) { SetTagBit(tag,kMCPi0Decay); SetTagBit(tag,kMCDecayDalitz);} //Pi0 Dalitz decay
+ else if (pPdg == 221) { SetTagBit(tag,kMCEtaDecay); SetTagBit(tag,kMCDecayDalitz);} //Eta Dalitz decay
else if((499 < pPdg && pPdg < 600)||(4999 < pPdg && pPdg < 6000)) { SetTagBit(tag,kMCEFromB);} //b-hadron decay
else if((399 < pPdg && pPdg < 500)||(3999 < pPdg && pPdg < 5000))
{ //c-hadron decay check
}
-//_________________________________________________________________________
-void AliMCAnalysisUtils::CheckOverlapped2GammaDecay(const Int_t *labels,
- const Int_t nlabels,
- const Int_t mesonIndex,
- AliStack *stack,
+//_________________________________________________________________________________________
+void AliMCAnalysisUtils::CheckOverlapped2GammaDecay(const Int_t *labels, Int_t nlabels,
+ Int_t mesonIndex, AliStack *stack,
Int_t &tag)
{
- //Check if cluster is formed from the contribution of 2 decay photons from pi0 or eta. Input in stack
+ // Check if cluster is formed from the contribution of 2 decay photons from pi0 or eta. Input in stack.
if(labels[0] < 0 || labels[0] > stack->GetNtrack() || nlabels <= 1) {
if(fDebug > 2) printf("AliMCAnalysisUtils::CheckOverlapped2GammaDecay(ESD) - Exit : label[0] %d, n primaries %d, nlabels %d \n",
}
-//__________________________________________________________________________________
-void AliMCAnalysisUtils::CheckOverlapped2GammaDecay(const Int_t *labels,
- const Int_t nlabels,
- const Int_t mesonIndex,
- const TClonesArray *mcparticles,
- Int_t & tag )
+//________________________________________________________________________________________________________
+void AliMCAnalysisUtils::CheckOverlapped2GammaDecay(const Int_t *labels, Int_t nlabels, Int_t mesonIndex,
+ const TClonesArray *mcparticles, Int_t & tag )
{
- //Check if cluster is formed from the contribution of 2 decay photons from pi0 or eta. Input in AODMCParticles
+ // Check if cluster is formed from the contribution of 2 decay photons from pi0 or eta. Input in AODMCParticles.
if(labels[0] < 0 || labels[0] > mcparticles->GetEntriesFast() || nlabels <= 1)
{
}
-//_________________________________________________________________________
+//______________________________________________________________________________________________________
+void AliMCAnalysisUtils::CheckLostDecayPair(const TObjArray* arrayCluster, Int_t iMom, Int_t iParent,
+ AliStack * stack, Int_t & tag)
+{
+ // Check on ESDs if the current decay photon has the second photon companion lost.
+
+ if(!arrayCluster || iMom < 0 || iParent < 0|| !stack) return;
+
+ TParticle * parent= stack->Particle(iParent);
+
+ if(parent->GetNDaughters()!=2)
+ {
+ SetTagBit(tag, kMCDecayPairLost);
+ return ;
+ }
+
+ Int_t pairLabel = -1;
+ if ( iMom != parent->GetDaughter(0) ) pairLabel = parent->GetDaughter(0);
+ else if( iMom != parent->GetDaughter(1) ) pairLabel = parent->GetDaughter(1);
+
+ if(pairLabel<0)
+ {
+ SetTagBit(tag, kMCDecayPairLost);
+ return ;
+ }
+
+ for(Int_t iclus = 0; iclus < arrayCluster->GetEntriesFast(); iclus++)
+ {
+ AliVCluster * cluster = (AliVCluster*) arrayCluster->At(iclus);
+ for(UInt_t ilab = 0; ilab< cluster->GetNLabels(); ilab++)
+ {
+ Int_t label = cluster->GetLabels()[ilab];
+ if(label==pairLabel)
+ {
+ SetTagBit(tag, kMCDecayPairInCalo);
+ return ;
+ }
+ else if(label== iParent || label== iMom)
+ {
+ continue;
+ }
+ else // check the ancestry
+ {
+ TParticle * mother = stack->Particle(label);
+ Int_t momPDG = TMath::Abs(mother->GetPdgCode());
+ if(momPDG!=11 && momPDG!=22) continue;
+
+ //Check if "mother" of entity is converted, if not, get the first non converted mother
+ Int_t iParentClus = mother->GetFirstMother();
+ if(iParentClus < 0) continue;
+
+ TParticle * parentClus = stack->Particle(iParentClus);
+ if(!parentClus) continue;
+
+ Int_t parentClusPDG = TMath::Abs(parentClus->GetPdgCode());
+ Int_t parentClusStatus = parentClus->GetStatusCode();
+
+ if( parentClusPDG != 22 && parentClusPDG != 11 && parentClusStatus != 0) continue;
+
+ //printf("Conversion\n");
+
+ //Check if the mother is photon or electron with status not stable
+ while ((parentClusPDG == 22 || parentClusPDG == 11) && parentClusStatus != 1)
+ {
+ //New Mother
+ label = iParentClus;
+ momPDG = parentClusPDG;
+
+ iParentClus = parentClus->GetFirstMother();
+ if(iParentClus < 0) break;
+
+ parentClus = stack->Particle(iParentClus);
+ if(!parentClus) break;
+
+ parentClusPDG = TMath::Abs(parentClus->GetPdgCode());
+ parentClusStatus = parentClus->GetStatusCode() ;
+ }//while
+
+ if((momPDG == 22 || parentClusPDG ==22) && (label==pairLabel || iParentClus == pairLabel))
+ {
+ SetTagBit(tag, kMCDecayPairInCalo);
+ //printf("Conversion is paired\n");
+ return ;
+ }
+ else continue;
+
+ }
+ }
+ } // cluster loop
+
+ SetTagBit(tag, kMCDecayPairLost);
+
+}
+
+//______________________________________________________________________________________________________
+void AliMCAnalysisUtils::CheckLostDecayPair(const TObjArray * arrayCluster,Int_t iMom, Int_t iParent,
+ const TClonesArray* mcparticles, Int_t & tag)
+{
+ // Check on AODs if the current decay photon has the second photon companion lost.
+
+ if(!arrayCluster || iMom < 0 || iParent < 0|| !mcparticles) return;
+
+ AliAODMCParticle * parent = (AliAODMCParticle*) mcparticles->At(iParent);
+
+ //printf("*** Check label %d with parent %d\n",iMom, iParent);
+
+ if(parent->GetNDaughters()!=2)
+ {
+ SetTagBit(tag, kMCDecayPairLost);
+ //printf("\t ndaugh = %d\n",parent->GetNDaughters());
+ return ;
+ }
+
+ Int_t pairLabel = -1;
+ if ( iMom != parent->GetDaughter(0) ) pairLabel = parent->GetDaughter(0);
+ else if( iMom != parent->GetDaughter(1) ) pairLabel = parent->GetDaughter(1);
+
+ if(pairLabel<0)
+ {
+ //printf("\t pair Label not found = %d\n",pairLabel);
+ SetTagBit(tag, kMCDecayPairLost);
+ return ;
+ }
+
+ //printf("\t *** find pair %d\n",pairLabel);
+
+ for(Int_t iclus = 0; iclus < arrayCluster->GetEntriesFast(); iclus++)
+ {
+ AliVCluster * cluster = (AliVCluster*) arrayCluster->At(iclus);
+ //printf("\t \t ** Cluster %d, nlabels %d\n",iclus,cluster->GetNLabels());
+ for(UInt_t ilab = 0; ilab< cluster->GetNLabels(); ilab++)
+ {
+ Int_t label = cluster->GetLabels()[ilab];
+
+ //printf("\t \t label %d\n",label);
+
+ if(label==pairLabel)
+ {
+ //printf("\t \t Pair found\n");
+ SetTagBit(tag, kMCDecayPairInCalo);
+ return ;
+ }
+ else if(label== iParent || label== iMom)
+ {
+ //printf("\t \t skip\n");
+ continue;
+ }
+ else // check the ancestry
+ {
+ AliAODMCParticle * mother = (AliAODMCParticle*) mcparticles->At(label);
+ Int_t momPDG = TMath::Abs(mother->GetPdgCode());
+ if(momPDG!=11 && momPDG!=22)
+ {
+ //printf("\t \t skip, pdg %d\n",momPDG);
+ continue;
+ }
+
+ //Check if "mother" of entity is converted, if not, get the first non converted mother
+ Int_t iParentClus = mother->GetMother();
+ if(iParentClus < 0) continue;
+
+ AliAODMCParticle * parentClus = (AliAODMCParticle*) mcparticles->At(iParentClus);
+ if(!parentClus) continue;
+
+ Int_t parentClusPDG = TMath::Abs(parentClus->GetPdgCode());
+ Int_t parentClusStatus = parentClus->GetStatus();
+
+ if( parentClusPDG != 22 && parentClusPDG != 11 && parentClusStatus != 0)
+ {
+ //printf("\t \t skip, not a conversion, parent: pdg %d, status %d\n",parentClusPDG,parentClusStatus);
+ continue;
+ }
+
+ //printf("\t \t Conversion\n");
+
+ //Check if the mother is photon or electron with status not stable
+ while ((parentClusPDG == 22 || parentClusPDG == 11) && parentClusStatus != 1)
+ {
+ //New Mother
+ label = iParentClus;
+ momPDG = parentClusPDG;
+
+ iParentClus = parentClus->GetMother();
+ if(iParentClus < 0) break;
+
+ parentClus = (AliAODMCParticle*) mcparticles->At(iParentClus);
+ if(!parentClus) break;
+
+ parentClusPDG = TMath::Abs(parentClus->GetPdgCode());
+ parentClusStatus = parentClus->GetStatus() ;
+ }//while
+
+ if((momPDG == 22 || parentClusPDG ==22) && (label==pairLabel || iParentClus == pairLabel))
+ {
+ SetTagBit(tag, kMCDecayPairInCalo);
+ //printf("\t \t Conversion is paired: mom %d, parent %d\n",label,iParentClus);
+ return ;
+ }
+ else
+ {
+ //printf("\t \t Skip, finally label %d, pdg %d, parent label %d, pdg %d, status %d\n",label,momPDG,iParentClus,parentClusPDG,parentClusStatus);
+ continue;
+ }
+
+ }
+ }
+ } // cluster loop
+
+
+ SetTagBit(tag, kMCDecayPairLost);
+
+
+}
+
+//_____________________________________________________________________
TList * AliMCAnalysisUtils::GetJets(const AliCaloTrackReader * reader)
{
- //Return list of jets (TParticles) and index of most likely parton that originated it.
+ // Return list of jets (TParticles) and index of most likely parton that originated it.
+
AliStack * stack = reader->GetStack();
Int_t iEvent = reader->GetEventNumber();
AliGenEventHeader * geh = reader->GetGenEventHeader();
//Get the jet, different way for different generator
//PYTHIA
- if(fMCGenerator == "PYTHIA")
+ if(fMCGenerator == kPythia)
{
TParticle * jet = 0x0;
AliGenPythiaEventHeader* pygeh= (AliGenPythiaEventHeader*) geh;
if(fDebug > 1)
printf("AliMCAnalysisUtils::GetJets() - PythiaEventHeader: Njets: %d\n",nTriggerJets);
- Int_t iparton = -1;
- for(Int_t i = 0; i< nTriggerJets; i++){
- iparton=-1;
+ for(Int_t i = 0; i< nTriggerJets; i++)
+ {
pygeh->TriggerJet(i, tmpjet);
jet = new TParticle(94, 21, -1, -1, -1, -1, tmpjet[0],tmpjet[1],tmpjet[2],tmpjet[3], 0,0,0,0);
//Assign an outgoing parton as mother
}
}//Pythia triggered jets
//HERWIG
- else if (fMCGenerator=="HERWIG"){
+ else if (fMCGenerator == kHerwig)
+ {
Int_t pdg = -1;
//Check parton 1
TParticle * tmp = parton1;
- if(parton1->GetPdgCode()!=22){
+ if(parton1->GetPdgCode()!=22)
+ {
while(pdg != 94){
if(tmp->GetFirstDaughter()==-1) return fJetsList;
tmp = stack->Particle(tmp->GetFirstDaughter());
//Check parton 2
pdg = -1;
tmp = parton2;
- Int_t i = -1;
- if(parton2->GetPdgCode()!=22){
- while(pdg != 94){
+ if(parton2->GetPdgCode()!=22)
+ {
+ while(pdg != 94)
+ {
if(tmp->GetFirstDaughter()==-1) return fJetsList;
- i = tmp->GetFirstDaughter();
tmp = stack->Particle(tmp->GetFirstDaughter());
pdg = tmp->GetPdgCode();
}//while
+
//Add found jet to list
TParticle *jet2 = new TParticle(*tmp);
jet2->SetFirstMother(7);
}
-//________________________________________________________________________________________________________
-TLorentzVector AliMCAnalysisUtils::GetDaughter(const Int_t idaugh, const Int_t label,
+//__________________________________________________________________________________________________________
+TLorentzVector AliMCAnalysisUtils::GetDaughter(Int_t idaugh, Int_t label,
const AliCaloTrackReader* reader,
Int_t & pdg, Int_t & status, Bool_t & ok, Int_t & daughlabel)
{
- //Return the kinematics of the particle that generated the signal, its pdg and its status and its label mother
-
- TLorentzVector daughter(0,0,0,0);
+ // Return the kinematics of the particle that generated the signal, its pdg and its status and its label mother.
+ fDaughMom.SetPxPyPzE(0,0,0,0);
if(reader->ReadStack())
{
printf("AliMCAnalysisUtils::GetDaughter() - Stack is not available, check analysis settings in configuration file, STOP!!\n");
ok=kFALSE;
- return daughter;
+ return fDaughMom;
}
- if(label >= 0 && label < reader->GetStack()->GetNtrack())
+
+ Int_t nprimaries = reader->GetStack()->GetNtrack();
+ if(label >= 0 && label < nprimaries)
{
TParticle * momP = reader->GetStack()->Particle(label);
daughlabel = momP->GetDaughter(idaugh);
+ if(daughlabel < 0 || daughlabel >= nprimaries)
+ {
+ ok = kFALSE;
+ return fDaughMom;
+ }
+
TParticle * daughP = reader->GetStack()->Particle(daughlabel);
- daughP->Momentum(daughter);
+ daughP->Momentum(fDaughMom);
pdg = daughP->GetPdgCode();
status = daughP->GetStatusCode();
}
else
{
ok = kFALSE;
- return daughter;
+ return fDaughMom;
}
}
else if(reader->ReadAODMCParticles())
printf("AliMCAnalysisUtils::GetDaughter() - AODMCParticles is not available, check analysis settings in configuration file!!\n");
ok=kFALSE;
- return daughter;
+ return fDaughMom;
}
Int_t nprimaries = mcparticles->GetEntriesFast();
AliAODMCParticle * momP = (AliAODMCParticle *) mcparticles->At(label);
daughlabel = momP->GetDaughter(idaugh);
+ if(daughlabel < 0 || daughlabel >= nprimaries)
+ {
+ ok = kFALSE;
+ return fDaughMom;
+ }
+
AliAODMCParticle * daughP = (AliAODMCParticle *) mcparticles->At(daughlabel);
- daughter.SetPxPyPzE(daughP->Px(),daughP->Py(),daughP->Pz(),daughP->E());
+ fDaughMom.SetPxPyPzE(daughP->Px(),daughP->Py(),daughP->Pz(),daughP->E());
pdg = daughP->GetPdgCode();
status = daughP->GetStatus();
}
else
{
ok = kFALSE;
- return daughter;
+ return fDaughMom;
}
}
ok = kTRUE;
- return daughter;
+ return fDaughMom;
}
-//_______________________________________________________________________________________________
-TLorentzVector AliMCAnalysisUtils::GetMother(const Int_t label, const AliCaloTrackReader* reader,
- Bool_t & ok)
+//______________________________________________________________________________________________________
+TLorentzVector AliMCAnalysisUtils::GetMother(Int_t label, const AliCaloTrackReader* reader, Bool_t & ok)
{
- //Return the kinematics of the particle that generated the signal
+ // Return the kinematics of the particle that generated the signal.
Int_t pdg = -1; Int_t status = -1; Int_t momlabel = -1;
return GetMother(label,reader,pdg,status, ok,momlabel);
}
-//_______________________________________________________________________________________________
-TLorentzVector AliMCAnalysisUtils::GetMother(const Int_t label, const AliCaloTrackReader* reader,
+//_________________________________________________________________________________________
+TLorentzVector AliMCAnalysisUtils::GetMother(Int_t label, const AliCaloTrackReader* reader,
Int_t & pdg, Int_t & status, Bool_t & ok)
{
- //Return the kinematics of the particle that generated the signal
+ // Return the kinematics of the particle that generated the signal.
Int_t momlabel = -1;
return GetMother(label,reader,pdg,status, ok,momlabel);
}
-//_______________________________________________________________________________________________
-TLorentzVector AliMCAnalysisUtils::GetMother(const Int_t label, const AliCaloTrackReader* reader,
+//______________________________________________________________________________________________________
+TLorentzVector AliMCAnalysisUtils::GetMother(Int_t label, const AliCaloTrackReader* reader,
Int_t & pdg, Int_t & status, Bool_t & ok, Int_t & momlabel)
{
- //Return the kinematics of the particle that generated the signal, its pdg and its status and its label mother
+ // Return the kinematics of the particle that generated the signal, its pdg and its status and its label mother.
- TLorentzVector mom(0,0,0,0);
+ fMotherMom.SetPxPyPzE(0,0,0,0);
if(reader->ReadStack())
{
printf("AliMCAnalysisUtils::GetMother() - Stack is not available, check analysis settings in configuration file, STOP!!\n");
ok=kFALSE;
- return mom;
+ return fMotherMom;
}
if(label >= 0 && label < reader->GetStack()->GetNtrack())
{
TParticle * momP = reader->GetStack()->Particle(label);
- momP->Momentum(mom);
- pdg = momP->GetPdgCode();
- status = momP->GetStatusCode();
+ momP->Momentum(fMotherMom);
+ pdg = momP->GetPdgCode();
+ status = momP->GetStatusCode();
momlabel = momP->GetFirstMother();
}
else
{
ok = kFALSE;
- return mom;
+ return fMotherMom;
}
}
else if(reader->ReadAODMCParticles())
printf("AliMCAnalysisUtils::GetMother() - AODMCParticles is not available, check analysis settings in configuration file!!\n");
ok=kFALSE;
- return mom;
+ return fMotherMom;
}
Int_t nprimaries = mcparticles->GetEntriesFast();
if(label >= 0 && label < nprimaries)
{
AliAODMCParticle * momP = (AliAODMCParticle *) mcparticles->At(label);
- mom.SetPxPyPzE(momP->Px(),momP->Py(),momP->Pz(),momP->E());
- pdg = momP->GetPdgCode();
- status = momP->GetStatus();
+ fMotherMom.SetPxPyPzE(momP->Px(),momP->Py(),momP->Pz(),momP->E());
+ pdg = momP->GetPdgCode();
+ status = momP->GetStatus();
momlabel = momP->GetMother();
}
else
{
ok = kFALSE;
- return mom;
+ return fMotherMom;
}
}
ok = kTRUE;
- return mom;
+ return fMotherMom;
}
-//_____________________________________________________________________________________
-TLorentzVector AliMCAnalysisUtils::GetMotherWithPDG(const Int_t label, const Int_t pdg, const AliCaloTrackReader* reader,
+//___________________________________________________________________________________
+TLorentzVector AliMCAnalysisUtils::GetMotherWithPDG(Int_t label, Int_t pdg,
+ const AliCaloTrackReader* reader,
Bool_t & ok, Int_t & momlabel)
{
- //Return the kinematics of the particle that generated the signal
-
- TLorentzVector grandmom(0,0,0,0);
+ // Return the kinematics of the particle that generated the signal.
+ fGMotherMom.SetPxPyPzE(0,0,0,0);
if(reader->ReadStack())
{
printf("AliMCAnalysisUtils::GetMotherWithPDG() - Stack is not available, check analysis settings in configuration file, STOP!!\n");
ok = kFALSE;
- return grandmom;
+ return fGMotherMom;
}
if(label >= 0 && label < reader->GetStack()->GetNtrack())
if(grandmomPDG==pdg)
{
momlabel = grandmomLabel;
- grandmom.SetPxPyPzE(grandmomP->Px(),grandmomP->Py(),grandmomP->Pz(),grandmomP->Energy());
+ fGMotherMom.SetPxPyPzE(grandmomP->Px(),grandmomP->Py(),grandmomP->Pz(),grandmomP->Energy());
break;
}
printf("AliMCAnalysisUtils::GetMotherWithPDG() - AODMCParticles is not available, check analysis settings in configuration file!!\n");
ok=kFALSE;
- return grandmom;
+ return fGMotherMom;
}
Int_t nprimaries = mcparticles->GetEntriesFast();
{
//printf("AliMCAnalysisUtils::GetMotherWithPDG(AOD) - mother with PDG %d FOUND! \n",pdg);
momlabel = grandmomLabel;
- grandmom.SetPxPyPzE(grandmomP->Px(),grandmomP->Py(),grandmomP->Pz(),grandmomP->E());
+ fGMotherMom.SetPxPyPzE(grandmomP->Px(),grandmomP->Py(),grandmomP->Pz(),grandmomP->E());
break;
}
ok = kTRUE;
- return grandmom;
+ return fGMotherMom;
}
-//____________________________________________________________________________________________________
-TLorentzVector AliMCAnalysisUtils::GetGrandMother(const Int_t label, const AliCaloTrackReader* reader,
+//______________________________________________________________________________________________
+TLorentzVector AliMCAnalysisUtils::GetGrandMother(Int_t label, const AliCaloTrackReader* reader,
Int_t & pdg, Int_t & status, Bool_t & ok,
Int_t & grandMomLabel, Int_t & greatMomLabel)
{
- //Return the kinematics of the particle that generated the signal
+ // Return the kinematics of the particle that generated the signal.
- TLorentzVector grandmom(0,0,0,0);
+ fGMotherMom.SetPxPyPzE(0,0,0,0);
if(reader->ReadStack())
{
printf("AliMCAnalysisUtils::GetGrandMother() - Stack is not available, check analysis settings in configuration file, STOP!!\n");
ok = kFALSE;
- return grandmom;
+ return fGMotherMom;
}
if(label >= 0 && label < reader->GetStack()->GetNtrack())
pdg = grandmomP->GetPdgCode();
status = grandmomP->GetStatusCode();
- grandmom.SetPxPyPzE(grandmomP->Px(),grandmomP->Py(),grandmomP->Pz(),grandmomP->Energy());
+ fGMotherMom.SetPxPyPzE(grandmomP->Px(),grandmomP->Py(),grandmomP->Pz(),grandmomP->Energy());
greatMomLabel = grandmomP->GetFirstMother();
}
printf("AliMCAnalysisUtils::GetGrandMother() - AODMCParticles is not available, check analysis settings in configuration file!!\n");
ok=kFALSE;
- return grandmom;
+ return fGMotherMom;
}
Int_t nprimaries = mcparticles->GetEntriesFast();
pdg = grandmomP->GetPdgCode();
status = grandmomP->GetStatus();
- grandmom.SetPxPyPzE(grandmomP->Px(),grandmomP->Py(),grandmomP->Pz(),grandmomP->E());
+ fGMotherMom.SetPxPyPzE(grandmomP->Px(),grandmomP->Py(),grandmomP->Pz(),grandmomP->E());
greatMomLabel = grandmomP->GetMother();
}
ok = kTRUE;
- return grandmom;
+ return fGMotherMom;
}
-//___________________________________________________________________________________________________________________________
-void AliMCAnalysisUtils::GetMCDecayAsymmetryAngleForPDG(const Int_t label, const Int_t pdg, const AliCaloTrackReader* reader,
+//_______________________________________________________________________________________________________________
+void AliMCAnalysisUtils::GetMCDecayAsymmetryAngleForPDG(Int_t label, Int_t pdg, const AliCaloTrackReader* reader,
Float_t & asym, Float_t & angle, Bool_t & ok)
{
- //In case of an eta or pi0 decay into 2 photons, get the asymmetry in the energy of the photons
-
+ // In case of an eta or pi0 decay into 2 photons, get the asymmetry in the energy of the photons.
if(reader->ReadStack())
{
if(d1->GetPdgCode() == 22 && d1->GetPdgCode() == 22)
{
asym = (d1->Energy()-d2->Energy())/grandmomP->Energy();
- TLorentzVector lvd1, lvd2;
- d1->Momentum(lvd1);
- d2->Momentum(lvd2);
- angle = lvd1.Angle(lvd2.Vect());
+ d1->Momentum(fDaughMom );
+ d2->Momentum(fDaughMom2);
+ angle = fDaughMom.Angle(fDaughMom2.Vect());
}
}
else
printf("AliMCAnalysisUtils::GetMCDecayAsymmetryForPDG() - AODMCParticles is not available, check analysis settings in configuration file!!\n");
ok=kFALSE;
+ return;
}
Int_t nprimaries = mcparticles->GetEntriesFast();
if(d1->GetPdgCode() == 22 && d1->GetPdgCode() == 22)
{
asym = (d1->E()-d2->E())/grandmomP->E();
- TLorentzVector lvd1, lvd2;
- lvd1.SetPxPyPzE(d1->Px(),d1->Py(),d1->Pz(),d1->E());
- lvd2.SetPxPyPzE(d2->Px(),d2->Py(),d2->Pz(),d2->E());
- angle = lvd1.Angle(lvd2.Vect());
+ fDaughMom .SetPxPyPzE(d1->Px(),d1->Py(),d1->Pz(),d1->E());
+ fDaughMom2.SetPxPyPzE(d2->Px(),d2->Py(),d2->Pz(),d2->E());
+ angle = fDaughMom.Angle(fDaughMom2.Vect());
}
}
else
}
-
-//_______________________________________________________________________________________________________
-Int_t AliMCAnalysisUtils::GetNDaughters(const Int_t label, const AliCaloTrackReader* reader, Bool_t & ok)
+//_________________________________________________________________________________________________
+Int_t AliMCAnalysisUtils::GetNDaughters(Int_t label, const AliCaloTrackReader* reader, Bool_t & ok)
{
- // Return the the number of daughters of a given MC particle
-
+ // Return the the number of daughters of a given MC particle.
if(reader->ReadStack())
{
return -1;
}
-//_______________________________________________________________________________
-Int_t AliMCAnalysisUtils::GetNOverlaps(const Int_t * label, const UInt_t nlabels,
- const Int_t mctag, const Int_t mesonLabel,
+//_________________________________________________________________________________
+Int_t AliMCAnalysisUtils::GetNOverlaps(const Int_t * label, UInt_t nlabels,
+ Int_t mctag, Int_t mesonLabel,
AliCaloTrackReader * reader, Int_t *overpdg)
{
// Compare the primary depositing more energy with the rest,
- // if no photon/electron (conversion) or neutral meson as comon ancestor, consider it as other particle contributing
- // Give as input the meson label in case it was a pi0 or eta merged cluster
- // Init overpdg with nlabels
+ // if no photon/electron (conversion) or neutral meson as comon ancestor, consider it as other particle contributing.
+ // Give as input the meson label in case it was a pi0 or eta merged cluster.
+ // Init overpdg with nlabels.
Int_t ancPDG = 0, ancStatus = -1;
- TLorentzVector momentum; TVector3 prodVertex;
+ TVector3 prodVertex;
Int_t ancLabel = 0;
Int_t noverlaps = 0;
Bool_t ok = kFALSE;
for (UInt_t ilab = 1; ilab < nlabels; ilab++ )
{
- ancLabel = CheckCommonAncestor(label[0],label[ilab],reader,ancPDG,ancStatus,momentum,prodVertex);
+ ancLabel = CheckCommonAncestor(label[0],label[ilab],reader,ancPDG,ancStatus,fMotherMom,prodVertex);
//printf("Overlaps, i %d: Main Label %d, second label %d, ancestor: Label %d, pdg %d - tag %d \n",
// ilab,label[0],label[ilab],ancLabel,ancPDG, mctag);
Int_t mpdg = -999999, gpdg = -1;
Int_t mstatus = -1, gstatus = -1;
Int_t gLabel = -1, ggLabel = -1;
- TLorentzVector mother = GetMother (label[ilab],reader,mpdg,mstatus,mOK);
- TLorentzVector grandmother = GetGrandMother(label[ilab],reader,gpdg,gstatus,gOK, gLabel,ggLabel);
+
+ GetMother (label[ilab],reader,mpdg,mstatus,mOK);
+ fGMotherMom =
+ GetGrandMother(label[ilab],reader,gpdg,gstatus,gOK, gLabel,ggLabel);
//printf("\t Overlap!, mother pdg %d; grand mother pdg %d",mpdg,gpdg);
while( ( gpdg == 22 || TMath::Abs(gpdg==11) ) && gLabel >=0 )
{
mpdg=gpdg;
- grandmother = GetGrandMother(labeltmp,reader,gpdg,gstatus,ok, gLabel,ggLabel);
+ fGMotherMom = GetGrandMother(labeltmp,reader,gpdg,gstatus,ok, gLabel,ggLabel);
labeltmp=gLabel;
}
}
//________________________________________________________
void AliMCAnalysisUtils::Print(const Option_t * opt) const
{
- //Print some relevant parameters set for the analysis
+ // Print some relevant parameters set for the analysis.
if(! opt)
return;
printf("***** Print: %s %s ******\n", GetName(), GetTitle() ) ;
printf("Debug level = %d\n",fDebug);
- printf("MC Generator = %s\n",fMCGenerator.Data());
+ printf("MC Generator = %s\n",fMCGeneratorString.Data());
printf(" \n");
}
-//________________________________________________________
-void AliMCAnalysisUtils::PrintMCTag(const Int_t tag) const
+//__________________________________________________
+void AliMCAnalysisUtils::PrintMCTag(Int_t tag) const
{
- // print the assigned origins to this particle
+ // Print the assigned origins to this particle.
- 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, other %d, unk %d, bad %d\n",
+ 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",
tag,
CheckTagBit(tag,kMCPhoton),
CheckTagBit(tag,kMCConversion),
CheckTagBit(tag,kMCKaon),
CheckTagBit(tag,kMCAntiProton),
CheckTagBit(tag,kMCAntiNeutron),
- CheckTagBit(tag,kMCOther),
CheckTagBit(tag,kMCUnknown),
CheckTagBit(tag,kMCBadLabel)
);
-
}
+//__________________________________________________
+void AliMCAnalysisUtils::SetMCGenerator(Int_t mcgen)
+{
+ // Set the generator type.
+
+ fMCGenerator = mcgen ;
+ if (mcgen == kPythia) fMCGeneratorString = "PYTHIA";
+ else if(mcgen == kHerwig) fMCGeneratorString = "HERWIG";
+ else if(mcgen == kHijing) fMCGeneratorString = "HIJING";
+ else
+ {
+ fMCGeneratorString = "";
+ fMCGenerator = kBoxLike ;
+ }
+
+}
+
+//__________________________________________________
+void AliMCAnalysisUtils::SetMCGenerator(TString mcgen)
+{
+ // Set the generator type.
+
+ fMCGeneratorString = mcgen ;
+
+ if (mcgen == "PYTHIA") fMCGenerator = kPythia;
+ else if(mcgen == "HERWIG") fMCGenerator = kHerwig;
+ else if(mcgen == "HIJING") fMCGenerator = kHijing;
+ else
+ {
+ fMCGenerator = kBoxLike;
+ fMCGeneratorString = "" ;
+ }
+}
+
+
git://git.uio.no / u / mrichter / AliRoot.git / blobdiff