1 //_________________________________________________________________________
2 // Utility Class for transverse energy studies
7 //*-- Authors: Oystein Djuvsland (Bergen)
8 //_________________________________________________________________________
9 #include "AliAnalysisEtSelector.h"
10 #include "AliAnalysisEtCuts.h"
11 #include "AliESDCaloCluster.h"
13 #include "TParticle.h"
14 #include "TParticlePDG.h"
15 #include "AliAnalysisEtCommon.h"
18 ClassImp(AliAnalysisEtSelector);
20 AliAnalysisEtSelector::AliAnalysisEtSelector(AliAnalysisEtCuts *cuts) : AliAnalysisEtCommon()
28 AliAnalysisEtSelector::AliAnalysisEtSelector() : AliAnalysisEtCommon()
35 fCuts = new AliAnalysisEtCuts();
38 AliAnalysisEtSelector::~AliAnalysisEtSelector()
45 void AliAnalysisEtSelector::SetEvent(const AliESDEvent* event)
48 if(!fInitialized) Init(event);
51 Bool_t AliAnalysisEtSelector::IsNeutralMcParticle(Int_t pIdx, AliStack& s, const TParticlePDG& pdg) const
53 return s.IsPhysicalPrimary(pIdx) &&(TMath::Abs(TMath::Abs(pdg.Charge()) - fCuts->GetMonteCarloNeutralParticle())<1e-3);
56 Bool_t AliAnalysisEtSelector::IsEmEtParticle(const Int_t pdgCode) const
58 return pdgCode == fgGammaCode || pdgCode == fgPi0Code || pdgCode == fgEtaCode || pdgCode == fgEPlusCode || pdgCode == fgEMinusCode;
62 Bool_t AliAnalysisEtSelector::PrimaryIsEmEtParticle(const Int_t pIdx, AliStack& stack) const
64 return IsEmEtParticle(stack.Particle(GetPrimary(pIdx, stack))->GetPdgCode());
66 Int_t AliAnalysisEtSelector::GetPrimary(const Int_t partIdx, AliStack& stack) const
70 Int_t mothIdx = stack.Particle(partIdx)->GetMother(0);
71 if(mothIdx < 0) return -1;
72 TParticle *mother = stack.Particle(mothIdx);
75 if(stack.IsPhysicalPrimary(mothIdx)) return mothIdx;
76 else return GetPrimary(mothIdx, stack);
80 std::cout << "WAT!" << std::endl;
86 //Bool_t AliAnalysisEtSelector::FromSecondaryInteraction(const TParticle& part, AliStack &stack) const
87 Bool_t AliAnalysisEtSelector::FromSecondaryInteraction(Int_t partID, AliStack &stack) const
89 // Bool_t partVtxSecondary = (
90 // TMath::Sqrt(part.Vx()*part.Vx() + part.Vy()*part.Vy()) > fCuts->GetPrimaryVertexCutXY()
91 // || TMath::Abs(part.Vz()) > fCuts->GetPrimaryVertexCutZ()
93 // && TMath::Sqrt(part.Vx()*part.Vx()+part.Vy()*part.Vy() + part.Vz()*part.Vz())<(fCuts->GetGeometryPhosDetectorRadius()-10);
96 //Bool_t partVtxSecondary = (TMath::Sqrt(part.Vx()*part.Vx() + part.Vy()*part.Vy()) <420);
97 //if(partVtxSecondary) return kFalse;
98 // //Let's find suspect decay (typical for secondary interaction)...
99 // if(partVtxSecondary){
100 // return SuspiciousDecayInChain(211, 111, part, stack);
101 //return stack.IsSecondaryFromMaterial(part.GetUniqueID());
102 return stack.IsSecondaryFromMaterial(partID);
113 Int_t AliAnalysisEtSelector::GetMother(Int_t partID, AliStack& stack) const {
115 TParticle *particle = stack.Particle(partID);
118 return particle->GetMother(0);
123 Bool_t AliAnalysisEtSelector::IsFromDetectorCover(Int_t partID, AliStack& stack) const{
125 TParticle *particle = stack.Particle(partID);
126 if(particle){//particle exists
127 if(stack.IsSecondaryFromMaterial(partID)){//particle is from an interaction with the material
128 //say that it's from the detector cover if its vertex is larger than 400 cm because this is where the cover is
129 return (TMath::Sqrt(particle->Vx()*particle->Vx() + particle->Vy()*particle->Vy()) >400);
137 Int_t AliAnalysisEtSelector::GetFirstMotherNotFromDetectorCover(Int_t partID, AliStack& stack) const{
138 Int_t targetParticle = partID;
139 Int_t testParticle = partID;
141 while(IsFromDetectorCover(targetParticle,stack) && testParticle>0){//if the particle is from the detector cover
142 //cout<<"Particle "<<targetParticle<<" is from detector cover. Getting mother ";
143 testParticle = GetMother(targetParticle,stack);
144 if(testParticle>0) targetParticle = testParticle;
145 //cout<<targetParticle<<" iteration "<<iteration<<endl;
148 //if(iteration>0) cout<<"iterations "<<iteration<<endl;
149 return targetParticle;
151 Bool_t AliAnalysisEtSelector::SuspiciousDecayInChain(const UInt_t suspectMotherPdg, const UInt_t suspectDaughterPdg, const TParticle &part, AliStack& stack) const
153 UInt_t partPdg = TMath::Abs(part.GetPdgCode());
154 if(part.GetFirstMother() == -1)
158 TParticle *mother = stack.Particle(part.GetFirstMother());
159 UInt_t motherPdg = TMath::Abs(mother->GetPdgCode());
160 if((suspectDaughterPdg==partPdg || 2112 == partPdg )&& suspectMotherPdg == motherPdg)
164 return SuspiciousDecayInChain(suspectMotherPdg, suspectDaughterPdg, *mother, stack);
167 Float_t AliAnalysisEtSelector::ShiftAngle(Float_t phi){//Always returns an angle in radians between -pi<phi<pi
169 while(myphi>TMath::Pi()){//angle is too high, decrease the angle
170 myphi = myphi - 2*TMath::Pi();
172 while(myphi<-TMath::Pi()){//angle is too low, increase the angle
173 myphi = myphi + 2*TMath::Pi();