* about the suitability of this software for any purpose. It is *
* provided "as is" without express or implied warranty. *
**************************************************************************/
-#include <TROOT.h>
//#include <TBranch.h>
//#include <TCanvas.h>
//#include <TChain.h>
//#include "AliAnalysisManager.h"
#include "AliAODEvent.h"
//#include "AliAODHandler.h"
-#include "AliAODInputHandler.h"
//#include "AliAODJet.h"
#include "AliAODMCParticle.h"
#include "AliAODTrack.h"
//#include "AliGenPythiaEventHeader.h"
#include "AliInputEventHandler.h"
//#include "AliKFVertex.h"
-#include "AliLog.h"
+//#include "AliLog.h"
#include "AliMCEvent.h"
-//#include "AliMCEventHandler.h"
-//#include "AliStack.h"
#include "AliVParticle.h"
+#include "AliAnalysisManager.h"
+#include "AliMCEventHandler.h"
+#include "AliStack.h"
+
+
////////////////////////////////////////////////
//---------------------------------------------
// Class for transverse regions analysis
fDebug(0),
fFilterBit(16),
fOnlyHadrons(kFALSE),
- fTrackEtaCut(0.8)
+ fTrackEtaCut(0.8),
+ fTrackPtMin(0),
+ fEsdTrackCuts(0x0),
+ fEsdTrackCutsSPD(0x0),
+ fEsdTrackCutsSDD(0x0)
{
// constructor
}
//____________________________________________________________________
-Bool_t AliAnalyseLeadingTrackUE::ApplyCuts(TObject* track, Int_t filterbit)
+Bool_t AliAnalyseLeadingTrackUE::ApplyCuts(TObject* track)
{
- // Reproduces the cuts of the corresponding bit in the ESD->AOD filtering
- // (see $ALICE_ROOT/ANALYSIS/macros/AddTaskESDFilter.C)
-
- AliESDtrackCuts* esdTrackCuts = new AliESDtrackCuts();
- switch (filterbit){
- /*case 16:
- // tighter cuts on primary particles for high pT tracks
- // needed as input for jetfinder
- esdTrackCuts->SetMinNClustersTPC(50);
- esdTrackCuts->SetMaxChi2PerClusterTPC(3.5);
- esdTrackCuts->SetRequireTPCRefit(kTRUE);
- esdTrackCuts->SetMaxDCAToVertexXY(2.4);
- esdTrackCuts->SetMaxDCAToVertexZ(3.2);
- esdTrackCuts->SetDCAToVertex2D(kTRUE);
- esdTrackCuts->SetRequireSigmaToVertex(kFALSE);
- esdTrackCuts->SetAcceptKinkDaughters(kFALSE);
- esdTrackCuts->SetRequireITSRefit(kTRUE); // additional cut
- break;
- */
-
- case 16:
- esdTrackCuts->GetStandardITSTPCTrackCuts2009(kTRUE);
- break;
-
- default:
- if (fDebug > 1) AliFatal("Set of cuts not defined");
- break;
- }
// select track according to set of cuts
- if (! esdTrackCuts->IsSelected(track) )return kFALSE;
-
+ if (!fEsdTrackCuts->IsSelected(track) )return kFALSE;
+ if (fEsdTrackCutsSPD && fEsdTrackCutsSDD && !fEsdTrackCutsSPD->IsSelected(track) && fEsdTrackCutsSDD->IsSelected(track)) return kFALSE;
+
return kTRUE;
}
+//____________________________________________________________________
+void AliAnalyseLeadingTrackUE::DefineESDCuts(Int_t filterbit) {
+
+ // Reproduces the cuts of the corresponding bit in the ESD->AOD filtering
+ // (see $ALICE_ROOT/ANALYSIS/macros/AddTaskESDFilter.C)
+
+ if (filterbit == -1)
+ filterbit = fFilterBit;
-
+ if (filterbit == 128)
+ {
+ fEsdTrackCuts = AliESDtrackCuts::GetStandardTPCOnlyTrackCuts();
+ fEsdTrackCuts->SetMinNClustersTPC(70);
+ }
+ else if (filterbit == 256)
+ {
+ // syst study
+ fEsdTrackCuts = AliESDtrackCuts::GetStandardTPCOnlyTrackCuts();
+ fEsdTrackCuts->SetMinNClustersTPC(80);
+ fEsdTrackCuts->SetMaxChi2PerClusterTPC(3);
+ fEsdTrackCuts->SetMaxDCAToVertexZ(2.7);
+ fEsdTrackCuts->SetMaxDCAToVertexXY(1.9);
+ }
+ else if (filterbit == 512)
+ {
+ // syst study
+ fEsdTrackCuts = AliESDtrackCuts::GetStandardTPCOnlyTrackCuts();
+ fEsdTrackCuts->SetMinNClustersTPC(60);
+ fEsdTrackCuts->SetMaxChi2PerClusterTPC(5);
+ fEsdTrackCuts->SetMaxDCAToVertexZ(3.7);
+ fEsdTrackCuts->SetMaxDCAToVertexXY(2.9);
+ }
+ else if (filterbit == 1024)
+ {
+ fEsdTrackCuts = AliESDtrackCuts::GetStandardTPCOnlyTrackCuts();
+ fEsdTrackCuts->SetMinNClustersTPC(-1);
+ fEsdTrackCuts->SetMinNCrossedRowsTPC(70);
+ fEsdTrackCuts->SetMinRatioCrossedRowsOverFindableClustersTPC(0.8);
+ }
+ else
+ {
+ fEsdTrackCuts = AliESDtrackCuts::GetStandardITSTPCTrackCuts2010();
+ fEsdTrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD, AliESDtrackCuts::kNone);
+
+ // Add SPD requirement
+ fEsdTrackCutsSPD = new AliESDtrackCuts("SPD", "Require 1 cluster in SPD");
+ fEsdTrackCutsSPD->SetClusterRequirementITS(AliESDtrackCuts::kSPD,AliESDtrackCuts::kAny);
+
+ // Add SDD requirement
+ fEsdTrackCutsSDD = new AliESDtrackCuts("SDD", "Require 1 cluster in first layer SDD");
+ fEsdTrackCutsSDD->SetClusterRequirementITS(AliESDtrackCuts::kSDD,AliESDtrackCuts::kFirst);
+ }
+}
//____________________________________________________________________
TObjArray* AliAnalyseLeadingTrackUE::FindLeadingObjects(TObject *obj)
}
+//-------------------------------------------------------------------
+TObjArray* AliAnalyseLeadingTrackUE::GetAcceptedParticles(TObject* obj, TObject* arrayMC, Bool_t onlyprimaries, Int_t particleSpecies, Bool_t useEtaPtCuts)
+{
+ // Returns an array of particles that pass the cuts, if arrayMC is given each reconstructed particle is replaced by its corresponding MC particles, depending on the parameter onlyprimaries only for primaries
+ // particleSpecies: -1 all particles are returned
+ // 0 (pions) 1 (kaons) 2 (protons) 3 (others) particles
+
+ Int_t nTracks = NParticles(obj);
+ TObjArray* tracks = new TObjArray;
+
+ // for TPC only tracks
+ Bool_t hasOwnership = kFALSE;
+ if ((fFilterBit == 128 || fFilterBit == 256 || fFilterBit == 512 || fFilterBit == 1024) && obj->InheritsFrom("AliESDEvent"))
+ hasOwnership = kTRUE;
+
+ if (!arrayMC)
+ tracks->SetOwner(hasOwnership);
+
+ // Loop over tracks or jets
+ for (Int_t ipart=0; ipart<nTracks; ++ipart) {
+ AliVParticle* part = ParticleWithCuts( obj, ipart, onlyprimaries, particleSpecies );
+ if (!part) continue;
+
+ if (useEtaPtCuts)
+ if (TMath::Abs(part->Eta()) > fTrackEtaCut || part->Pt() < fTrackPtMin)
+ {
+ if (hasOwnership)
+ delete part;
+ continue;
+ }
+
+ if (arrayMC) {
+ Int_t label = part->GetLabel();
+ if (hasOwnership)
+ delete part;
+ // re-define part as the matched MC particle
+ part = ParticleWithCuts(arrayMC, TMath::Abs(label),onlyprimaries, particleSpecies);
+ if (!part)continue;
+ }
+
+ tracks->Add(part);
+ }
+
+ return tracks;
+}
+
//-------------------------------------------------------------------
TObjArray* AliAnalyseLeadingTrackUE::GetMinMaxRegion(TList *transv1, TList *transv2)
{
Int_t nTracks;
if (obj->InheritsFrom("TClonesArray")){ // MC particles
- TClonesArray *arrayMC = dynamic_cast<TClonesArray*>(obj);
+ TClonesArray *arrayMC = static_cast<TClonesArray*>(obj);
nTracks = arrayMC->GetEntriesFast();
}else if (obj->InheritsFrom("TObjArray")){ // list of AliVParticle
- TObjArray *array = dynamic_cast<TObjArray*>(obj);
+ TObjArray *array = static_cast<TObjArray*>(obj);
nTracks = array->GetEntriesFast();
}else if (obj->InheritsFrom("AliAODEvent")){ // RECO AOD tracks
- AliAODEvent *aodEvent =
dynamic_cast<AliAODEvent*>(obj);
+ AliAODEvent *aodEvent =
static_cast<AliAODEvent*>(obj);
nTracks = aodEvent->GetNTracks();
}else if (obj->InheritsFrom("AliESDEvent")){ // RECO ESD tracks
- AliESDEvent *esdEvent =
dynamic_cast<AliESDEvent*>(obj);
+ AliESDEvent *esdEvent =
static_cast<AliESDEvent*>(obj);
nTracks = esdEvent->GetNumberOfTracks();
+ }else if (obj->InheritsFrom("AliMCEvent")){ // RECO ESD tracks
+ AliMCEvent *mcEvent = static_cast<AliMCEvent*>(obj);
+ nTracks = mcEvent->GetNumberOfTracks();
}else {
if (fDebug > 1) AliFatal(" Analysis type not defined !!! ");
return 0;
//-------------------------------------------------------------------
-AliVParticle* AliAnalyseLeadingTrackUE::ParticleWithCuts(TObject* obj, Int_t ipart, Bool_t onlyprimaries)
+AliVParticle* AliAnalyseLeadingTrackUE::ParticleWithCuts(TObject* obj, Int_t ipart, Bool_t onlyprimaries, Int_t particleSpecies)
{
// Returns track or MC particle at position "ipart" if passes selection criteria
+ // particleSpecies: -1 all particles are returned
+ // 0 (pions) 1 (kaons) 2 (protons) 3 (others) particles
AliVParticle *part=0;
if (obj->InheritsFrom("TClonesArray")){ // AOD-MC PARTICLE
- TClonesArray *arrayMC = dynamic_cast<TClonesArray*>(obj);
+ TClonesArray *arrayMC = static_cast<TClonesArray*>(obj);
part = (AliVParticle*)arrayMC->At( ipart );
if (!part)return 0;
// eventually only primaries
TMath::Abs(pdgCode)==321; // Kaon
if (!isHadron) return 0;
}
+ if (particleSpecies != -1) {
+ // find the primary mother
+ AliVParticle* mother = part;
+ while (!((AliAODMCParticle*)mother)->IsPhysicalPrimary())
+ {
+ if (((AliAODMCParticle*)mother)->GetMother() < 0)
+ {
+ mother = 0;
+ break;
+ }
+
+ mother = (AliVParticle*) arrayMC->At(((AliAODMCParticle*)mother)->GetMother());
+ if (!mother)
+ break;
+ }
+
+ if (mother)
+ {
+ Int_t pdgCode = ((AliAODMCParticle*)mother)->GetPdgCode();
+ if (particleSpecies == 0 && TMath::Abs(pdgCode)!=211)
+ return 0;
+ if (particleSpecies == 1 && TMath::Abs(pdgCode)!=321)
+ return 0;
+ if (particleSpecies == 2 && TMath::Abs(pdgCode)!=2212)
+ return 0;
+ if (particleSpecies == 3 && (TMath::Abs(pdgCode)==211 || TMath::Abs(pdgCode)==321 || TMath::Abs(pdgCode)==2212))
+ return 0;
+ }
+ else
+ {
+ // if mother not found, accept particle only in case of particleSpecies == 3. To include it in all or no sample is no solution
+ Printf("WARNING: No mother found for particle %d:", part->GetLabel());
+ part->Print();
+
+ /*
+ // this code prints the details of the mother that is missing in the AOD
+ AliMCEventHandler* fMcHandler = dynamic_cast<AliMCEventHandler*> (AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler());
+
+ AliMCEvent* fMcEvent = fMcHandler->MCEvent();
+
+ fMcEvent->Stack()->Particle(fMcEvent->Stack()->Particle(part->GetLabel())->GetMother(0))->Print();
+ fMcEvent->Stack()->Particle(fMcEvent->Stack()->Particle(fMcEvent->Stack()->Particle(part->GetLabel())->GetMother(0))->GetMother(0))->Print();
+ Printf("eta = %f", fMcEvent->Stack()->Particle(fMcEvent->Stack()->Particle(fMcEvent->Stack()->Particle(part->GetLabel())->GetMother(0))->GetMother(0))->Eta());
+ */
+
+ if (particleSpecies != 3)
+ return 0;
+ }
+ }
}else if (obj->InheritsFrom("TObjArray")){ // list of AliVParticle
- TObjArray *array = dynamic_cast<TObjArray*>(obj);
+ TObjArray *array = static_cast<TObjArray*>(obj);
part = (AliVParticle*)array->At( ipart );
if (!part)return 0;
}else if (obj->InheritsFrom("AliMCEvent")){ // MC PARTICLE
- AliMCEvent* mcEvent =
dynamic_cast<AliMCEvent*>(obj);
+ AliMCEvent* mcEvent =
static_cast<AliMCEvent*>(obj);
part = mcEvent->GetTrack( ipart );
if (!part) return 0;
// eventually only primaries
if (!isHadron) return 0;
}
*/
-
+ if (particleSpecies != -1) {
+ // find the primary mother
+ AliMCParticle* mother = (AliMCParticle*) part;
+// Printf("");
+ while (!mcEvent->IsPhysicalPrimary(mother->GetLabel()))
+ {
+// Printf("pdg = %d; mother = %d", mother->PdgCode(), mother->GetMother());
+ if (mother->GetMother() < 0)
+ {
+ mother = 0;
+ break;
+ }
+
+ mother = (AliMCParticle*) mcEvent->GetTrack(mother->GetMother());
+ if (!mother)
+ break;
+ }
+
+ if (mother)
+ {
+ Int_t pdgCode = mother->PdgCode();
+ if (particleSpecies == 0 && TMath::Abs(pdgCode)!=211)
+ return 0;
+ if (particleSpecies == 1 && TMath::Abs(pdgCode)!=321)
+ return 0;
+ if (particleSpecies == 2 && TMath::Abs(pdgCode)!=2212)
+ return 0;
+ if (particleSpecies == 3 && (TMath::Abs(pdgCode)==211 || TMath::Abs(pdgCode)==321 || TMath::Abs(pdgCode)==2212))
+ return 0;
+ }
+ else
+ {
+ // if mother not found, accept particle only in case of particleSpecies == 3. To include it in all or no sample is no solution
+ Printf("WARNING: No mother found for particle %d:", part->GetLabel());
+ //part->Dump();
+ //part->Print();
+
+ /*
+ // this code prints the details of the mother that is missing in the AOD
+ AliMCEventHandler* fMcHandler = dynamic_cast<AliMCEventHandler*> (AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler());
+
+ AliMCEvent* fMcEvent = fMcHandler->MCEvent();
+
+ fMcEvent->Stack()->Particle(fMcEvent->Stack()->Particle(part->GetLabel())->GetMother(0))->Print();
+ fMcEvent->Stack()->Particle(fMcEvent->Stack()->Particle(fMcEvent->Stack()->Particle(part->GetLabel())->GetMother(0))->GetMother(0))->Print();
+ Printf("eta = %f", fMcEvent->Stack()->Particle(fMcEvent->Stack()->Particle(fMcEvent->Stack()->Particle(part->GetLabel())->GetMother(0))->GetMother(0))->Eta());
+ */
+
+ if (particleSpecies != 3)
+ return 0;
+ }
+ }
}else if (obj->InheritsFrom("AliAODEvent")){ // RECO AOD TRACKS
- AliAODEvent *aodEvent =
dynamic_cast<AliAODEvent*>(obj);
+ AliAODEvent *aodEvent =
static_cast<AliAODEvent*>(obj);
part = aodEvent->GetTrack(ipart);
// track selection cuts
- if (!( ((AliAODTrack*)part)->TestFilterBit(fFilterBit)) )return 0;
+ if ( !(((AliAODTrack*)part)->TestFilterBit(fFilterBit)) ) return 0;
+ //if ( !(((AliAODTrack*)part)->TestFilterBit(fFilterBit)) && !(((AliAODTrack*)part)->TestFilterBit(32)) ) return 0;
// only primary candidates
//if ( ((AliAODTrack*)part)->IsPrimaryCandidate() )return 0;
// eventually only hadrons
}
}else if (obj->InheritsFrom("AliESDEvent")){ // RECO ESD TRACKS
- AliESDEvent *esdEvent =
dynamic_cast<AliESDEvent*>(obj);
+ AliESDEvent *esdEvent =
static_cast<AliESDEvent*>(obj);
part = esdEvent->GetTrack(ipart);
if (!part)return 0;
// track selection cuts
- if (!( ApplyCuts(part, fFilterBit)) )return 0;
- // only primary candidates (does not exist for ESD tracks??????)
- //if ( ((AliAODTrack*)part)->IsPrimaryCandidate() )return 0;
+ if (!( ApplyCuts(part)) )
+ return 0;
+
+ if (fFilterBit == 128 || fFilterBit == 256 || fFilterBit == 512 || fFilterBit == 1024)
+ {
+ // create TPC only tracks constrained to the SPD vertex
+
+ const AliESDVertex *vtxSPD = esdEvent->GetPrimaryVertexSPD();
+
+ AliESDtrack* track = AliESDtrackCuts::GetTPCOnlyTrack(esdEvent, ipart);
+ if(!track) return 0;
+
+ if(track->Pt()>0.){
+ // only constrain tracks above threshold
+ AliExternalTrackParam exParam;
+ // take the B-feild from the ESD, no 3D fieldMap available at this point
+ Bool_t relate = kFALSE;
+ relate = track->RelateToVertexTPC(vtxSPD,esdEvent->GetMagneticField(),kVeryBig,&exParam);
+ if(!relate)
+ {
+// Printf("relating failed");
+ delete track;
+ return 0;
+ }
+ track->Set(exParam.GetX(),exParam.GetAlpha(),exParam.GetParameter(),exParam.GetCovariance());
+ }
+
+ part = track;
+ }
// eventually only hadrons
//TO-DO
TList *transverse2 = new TList();
TObjArray *regionParticles = new TObjArray;
+ regionParticles->SetOwner(kTRUE);
+
regionParticles->AddLast(toward);
regionParticles->AddLast(away);
regionParticles->AddLast(transverse1);
if (TMath::Abs(leadVect.DeltaPhi(partVect)) > k120rad ) region = -2; //backward
// skip leading particle
- if (leading == part)
+ if (leading == part)
continue;
if (!region)continue;
//____________________________________________________________________
Bool_t AliAnalyseLeadingTrackUE::TriggerSelection(const TObject* obj)
{
+ if (!obj) // MC
+ return kFALSE;
- //Use AliPhysicsSelection to select good events
- if( !obj->InheritsFrom("AliAODInputHandler") ) { // not needed for AOD input
- if (!(((AliInputEventHandler*)obj)->IsEventSelected()))return kFALSE;
- }
+ // Use AliPhysicsSelection to select good events, works for ESD and AOD
+ if (!(((AliInputEventHandler*)obj)->IsEventSelected()&(AliVEvent::kMB|AliVEvent::kUserDefined)))
+ return kFALSE;
return kTRUE;
-
}
//____________________________________________________________________
Int_t nVertex = ((AliAODEvent*)obj)->GetNumberOfVertices();
if( nVertex > 0 ) {
AliAODVertex* vertex = (AliAODVertex*)((AliAODEvent*)obj)->GetPrimaryVertex();
- Int_t nTracksPrim = vertex->GetNContributors();
+ Int_t nTracksPrim = vertex->GetNContributors();
Double_t zVertex = vertex->GetZ();
if (fDebug > 1)AliInfo(Form(" Vertex in = %f with %d particles by %s data ...",zVertex,nTracksPrim,vertex->GetName()));
- // Select a quality vertex by number of tracks?
+ // Reject TPC only vertex
+ TString name(vertex->GetName());
+ if (name.CompareTo("PrimaryVertex") && name.CompareTo("SPDVertex"))return kFALSE;
+
+ // Select a quality vertex by number of tracks?
if( nTracksPrim < ntracks || TMath::Abs(zVertex) > zed ) {
if (fDebug > 1) AliInfo(" Primary-vertex Selection: event REJECTED ...");
return kFALSE;
return kFALSE;
}
}
- // TO-DO: ESD case for DCA studies
+
+ // ESD case for DCA studies
+ if (obj->InheritsFrom("AliESDEvent")){
+ AliESDVertex* vertex = (AliESDVertex*)((AliESDEvent*)obj)->GetPrimaryVertex();
+ if ( vertex){
+ Int_t nTracksPrim = vertex->GetNContributors();
+ Double_t zVertex = vertex->GetZ();
+ if (fDebug > 1)AliInfo(Form(" Vertex in = %f with %d particles by %s data ...",zVertex,nTracksPrim,vertex->GetName()));
+ // Reject SPD or TPC only vertex
+ TString name(vertex->GetName());
+ if (name.CompareTo("PrimaryVertex") && name.CompareTo("SPDVertex"))return kFALSE;
+
+ // Select a quality vertex by number of tracks?
+ if( nTracksPrim < ntracks || TMath::Abs(zVertex) > zed ) {
+ if (fDebug > 1) AliInfo(" Primary-vertex Selection: event REJECTED ...");
+ return kFALSE;
+ }
+ // TODO remove vertexer Z events with dispersion > 0.02: Doesn't work for AOD at present
+ //if (strcmp(vertex->GetTitle(), "AliVertexerZ") == 0 && vertex->GetDispersion() > 0.02)
+ // return kFALSE;
+ if (fDebug > 1) AliInfo(" Primary-vertex Selection: event ACCEPTED...");
+ } else {
+ if (fDebug > 1) AliInfo(" Primary-vertex Selection: event REJECTED ...");
+ return kFALSE;
+ }
+ }
return kTRUE;
}
git://git.uio.no / u / mrichter / AliRoot.git / blobdiff