--- /dev/null
+//
+// Select MC tracks, feed them via PicoTrack-like list
+
+#include "AliMCTruthTrackMaker.h"
+
+#include <TClonesArray.h>
+#include <TString.h>
+
+#include "AliAnalysisManager.h"
+#include "AliVEventHandler.h"
+#include "AliVEvent.h"
+#include "AliMCEvent.h"
+#include "AliVParticle.h"
+#include "AliAODMCParticle.h"
+#include "AliMCParticle.h"
+#include "AliLog.h"
+#include "AliESDMuonTrack.h"
+#include "AliAODTrack.h"
+#include "AliESDEvent.h"
+#include "AliStack.h"
+
+ClassImp(AliMCTruthTrackMaker)
+
+//________________________________________________________________________
+AliMCTruthTrackMaker::AliMCTruthTrackMaker() :
+AliAnalysisTaskSE("AliMCTruthTrackMaker"),
+fTracksOutName("MCTruthTracks"),
+fChargedMC(kFALSE),
+fFillMuMothers(kFALSE),
+fTriggerMatch(kTRUE),
+fInit(kFALSE),
+fEsdMode(kFALSE),
+fTracksIn(0x0),
+fTracksOut(0x0),
+fESD(0x0),
+fMC(0x0),
+fStack(0x0)
+{
+ // Constructor.
+}
+
+//________________________________________________________________________
+AliMCTruthTrackMaker::AliMCTruthTrackMaker(const char *name) :
+AliAnalysisTaskSE(name),
+fTracksOutName("MCTruthTracks"),
+fChargedMC(kFALSE),
+fFillMuMothers(kFALSE),
+fTriggerMatch(kTRUE),
+fInit(kFALSE),
+fEsdMode(kFALSE),
+fTracksIn(0x0),
+fTracksOut(0x0),
+fESD(0x0),
+fMC(0x0),
+fStack(0x0)
+{
+ // Constructor.
+}
+
+//________________________________________________________________________
+AliMCTruthTrackMaker::~AliMCTruthTrackMaker()
+{
+ // Destructor.
+}
+
+//________________________________________________________________________
+void AliMCTruthTrackMaker::UserCreateOutputObjects()
+{
+ // Create my user objects.
+
+ AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
+ if (!mgr) {
+ AliFatal("No analysis manager!");
+ return;
+ }
+
+ AliVEventHandler *handler = mgr->GetInputEventHandler();
+ if (!handler) {
+ AliFatal("No input handler!");
+ return;
+ }
+
+ if (handler->InheritsFrom("AliESDInputHandler"))
+ fEsdMode = kTRUE;
+ else
+ fEsdMode = kFALSE;
+
+ if (fEsdMode)
+ fTracksOut = new TClonesArray("AliMCParticle");
+ else
+ fTracksOut = new TClonesArray("AliAODMCParticle");
+ fTracksOut->SetName(fTracksOutName);
+
+}
+
+//________________________________________________________________________
+void AliMCTruthTrackMaker::UserExec(Option_t *)
+{
+ // Main loop, called for each event.
+ if (!InputEvent()) {
+ AliError("Could not retrieve event! Returning");
+ return;
+ }
+
+ if (fEsdMode && !MCEvent()) {
+ AliError("Could not retrieve ESD MC event! Returning");
+ return;
+ }
+
+ if (!fInit) {
+ // add tracks to event if not yet there
+ if (!(InputEvent()->FindListObject(fTracksOutName)))
+ InputEvent()->AddObject(fTracksOut);
+
+ if (!fEsdMode) {
+ fTracksIn = static_cast<TClonesArray*>(InputEvent()->FindListObject(AliAODMCParticle::StdBranchName()));
+ if (!fTracksIn) {
+ AliError("Could not retrieve AOD MC particles! Returning");
+ return;
+ }
+ TClass *cl = fTracksIn->GetClass();
+ if (!cl->GetBaseClass("AliAODMCParticle")) {
+ AliError(Form("%s: Collection %s does not contain AliAODMCParticle!", GetName(), AliAODMCParticle::StdBranchName()));
+ fTracksIn = 0;
+ return;
+ }
+ }
+
+ fInit = kTRUE;
+ }
+
+ // clear container (normally a null operation as the event should clean it already)
+ fTracksOut->Delete();
+
+ Int_t nacc = 0;
+
+ // first, fill MC tracks up to fEtaMax
+ const Int_t Ntracks = GetNumberOfTracks();
+ for (Int_t iTracks = 0; iTracks < Ntracks; ++iTracks) {
+ AliVParticle* track = GetTrack(iTracks);
+
+ if (!track)
+ continue;
+
+ if (TMath::Abs(track->Eta()) > fEtaMax)
+ continue;
+
+ if (fChargedMC && track->Charge() == 0) continue;
+
+ AddTrack(track, nacc);
+
+ ++nacc;
+ }
+
+ // now loop over reconstructed tracks and if accepted by muon cuts, add their primary mother
+ if (fFillMuMothers) {
+ if (!fEsdMode) {
+ AliError("Muon mother particles can only be filled when running on ESD! Returning");
+ return;
+ }
+
+ fESD = dynamic_cast<AliESDEvent*>(InputEvent());
+ if (!fESD) { AliError("ESD event not found. Nothing done!"); return; }
+ Int_t ntrks = fESD->GetNumberOfMuonTracks();
+
+ fMC = MCEvent();
+ if (!fMC) { AliError("MC event not available."); return; }
+ fStack = fMC->Stack();
+
+ for (Int_t iTrack = 0; iTrack<ntrks; iTrack++)
+ {
+ Int_t label = 0;
+ AliESDMuonTrack* muonTrack = fESD->GetMuonTrack(iTrack);
+ if(muonTrack)
+ {
+ if (IsGoodMUONtrack(*muonTrack)) {
+ label = TMath::Abs(muonTrack->GetLabel());
+ if (label>=fMC->GetNumberOfTracks()) {
+ AliError(Form("Label %d larger than number of particles on stack %d\n",label,fMC->GetNumberOfTracks()));
+ continue;
+ }
+ AliVParticle * fpmMu = (AliVParticle *) fMC->GetTrack(GetFirstPrimaryMother(label));
+ if (fpmMu) {
+ AddTrack(fpmMu, nacc);
+ ++nacc;
+ } else {
+ AliError("didn't receive a valid first primary mother here");
+ }
+ }
+ }
+ }
+
+ }
+}
+
+//________________________________________________________________________
+Int_t AliMCTruthTrackMaker::GetNumberOfTracks() const
+{
+ if (fEsdMode)
+ return MCEvent()->GetNumberOfTracks();
+ else
+ return fTracksIn->GetEntries();
+}
+
+//________________________________________________________________________
+AliVParticle* AliMCTruthTrackMaker::GetTrack(Int_t i)
+{
+ if (fEsdMode) {
+ if (!MCEvent()->IsPhysicalPrimary(i))
+ return 0;
+
+ return MCEvent()->GetTrack(i);
+ }
+ else {
+ AliAODMCParticle *part = static_cast<AliAODMCParticle*>(fTracksIn->At(i));
+ if (!part->IsPhysicalPrimary())
+ return 0;
+
+ return part;
+ }
+}
+
+//________________________________________________________________________
+void AliMCTruthTrackMaker::AddTrack(AliVParticle *track, Int_t nacc)
+{
+ if (fEsdMode) {
+ AliMCParticle *part = static_cast<AliMCParticle*>(track);
+ new ((*fTracksOut)[nacc]) AliMCParticle(part->Particle(), 0, part->Label());
+ }
+ else {
+ AliAODMCParticle *part = static_cast<AliAODMCParticle*>(track);
+ new ((*fTracksOut)[nacc]) AliAODMCParticle(*part);
+ }
+}
+
+//________________________________________________________________________
+Bool_t AliMCTruthTrackMaker::IsGoodMUONtrack(AliESDMuonTrack &track)
+{
+ // Applying track cuts for MUON tracks
+
+ if (!track.ContainTrackerData())
+ return kFALSE;
+ Double_t thetaTrackAbsEnd = TMath::ATan(track.GetRAtAbsorberEnd()/505.) * TMath::RadToDeg();
+ if ((thetaTrackAbsEnd < 2.) || (thetaTrackAbsEnd > 10.))
+ return kFALSE;
+ Double_t eta = track.Eta();
+ if ((eta < -4.) || (eta > -2.5))
+ return kFALSE;
+ if (fTriggerMatch) {
+ if (!track.ContainTriggerData())
+ return kFALSE;
+ if (track.GetMatchTrigger() < 0.5)
+ return kFALSE;
+ }
+ return kTRUE;
+}
+
+//________________________________________________________________________
+Bool_t AliMCTruthTrackMaker::IsGoodMUONtrack(AliAODTrack &track)
+{
+ // Applying track cuts for MUON tracks
+
+ if (!track.IsMuonTrack())
+ return kFALSE;
+ Double_t dThetaAbs = TMath::ATan(track.GetRAtAbsorberEnd()/505.)* TMath::RadToDeg();
+ if ((dThetaAbs<2.) || (dThetaAbs>10.))
+ return kFALSE;
+ Double_t dEta = track.Eta();
+ if ((dEta<-4.) || (dEta>-2.5))
+ return kFALSE;
+ if (fTriggerMatch) {
+ if (track.GetMatchTrigger()<0.5)
+ return kFALSE;
+ }
+ return kTRUE;
+}
+
+//________________________________________________________________________
+Int_t AliMCTruthTrackMaker::GetFirstPrimaryMother(Int_t muonlabel)
+{
+ AliMCParticle *McParticle = (AliMCParticle*)fMC->GetTrack(muonlabel);
+ if(McParticle->GetMother()<fStack->GetNprimary()) return McParticle->GetMother();
+ else
+ {
+ Int_t motherlabel = McParticle->GetMother();
+ while(motherlabel > -1)
+ {
+ AliMCParticle *MotherParticle = (AliMCParticle*)fMC->GetTrack(motherlabel);
+ if(MotherParticle->GetMother()<fStack->GetNprimary()) break;
+ else motherlabel = MotherParticle->GetMother();
+ }
+ AliMCParticle *FirstSecondaryMotherParticle = (AliMCParticle*)fMC->GetTrack(motherlabel);
+ return FirstSecondaryMotherParticle->GetMother();
+ }
+}
--- /dev/null
+#ifndef ALIMCTRUTHTRACKMAKER_H
+#define ALIMCTRUTHTRACKMAKER_H
+
+class TClonesArray;
+class AliESDEvent;
+class AliESDMuonTrack;
+class AliMCEvent;
+class AliAODTrack;
+class AliVParticle;
+class AliStack;
+
+#include "TString.h"
+#include "AliAnalysisTaskSE.h"
+
+class AliMCTruthTrackMaker : public AliAnalysisTaskSE {
+ public:
+ AliMCTruthTrackMaker();
+ AliMCTruthTrackMaker(const char *name);
+ virtual ~AliMCTruthTrackMaker();
+
+ void UserCreateOutputObjects();
+ void UserExec(Option_t *option);
+
+ void SetChargedMC(Bool_t c = kTRUE) { fChargedMC = c ; }
+ void SetEtaMax(Double_t e) { fEtaMax = e ; }
+ void SetFillMuMothers(Bool_t c = kTRUE) { fFillMuMothers = c ; }
+ void SetTriggerMatch(Bool_t b) { fTriggerMatch = b ; }
+ void SetTracksOutName(const char *name) { fTracksOutName = name ; }
+
+ protected:
+ Int_t GetNumberOfTracks() const;
+ AliVParticle *GetTrack(Int_t i);
+ void AddTrack(AliVParticle *track, Int_t nacc);
+ Bool_t IsGoodMUONtrack(AliESDMuonTrack &track);
+ Bool_t IsGoodMUONtrack(AliAODTrack &track);
+ Int_t GetFirstPrimaryMother(Int_t muonlabel);
+
+ TString fTracksOutName; // name of output track array
+ Bool_t fChargedMC; // true = only charged particles
+ Bool_t fFillMuMothers; // true = fill primary mother of reconstructed muon tracksmuon cuts
+ Bool_t fTriggerMatch; // require trigger match as muon cut?
+ Double_t fEtaMax; // maximum eta to accept tracks
+ Bool_t fInit; //!true = task initialized
+ Bool_t fEsdMode; //!switch for ESD/AOD mode
+ TClonesArray *fTracksIn; //!track array in (AOD only)
+ TClonesArray *fTracksOut; //!track array out
+ AliESDEvent *fESD; //! ESD object
+ AliMCEvent *fMC; //! MC object
+ AliStack *fStack; //! MC stack
+
+ private:
+ AliMCTruthTrackMaker(const AliMCTruthTrackMaker&); // not implemented
+ AliMCTruthTrackMaker &operator=(const AliMCTruthTrackMaker&); // not implemented
+
+ ClassDef(AliMCTruthTrackMaker, 1); // Task to select tracks in MC events
+};
+#endif
--- /dev/null
+AliMCTruthTrackMaker* AddTaskMCTruthTrackMaker(
+ const char *outname = "MCTruthTracks",
+ Bool_t mumothers = kFALSE,
+ Bool_t ch = kFALSE,
+ Double_t etamax = 1.0
+)
+{
+ // Get the pointer to the existing analysis manager via the static access method.
+ //==============================================================================
+ AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
+ if (!mgr)
+ {
+ ::Error("AddTaskMCTrackSelector", "No analysis manager to connect to.");
+ return NULL;
+ }
+
+ // Check the analysis type using the event handlers connected to the analysis manager.
+ //==============================================================================
+ if (!mgr->GetInputEventHandler())
+ {
+ ::Error("AddTaskMCTrackSelector", "This task requires an input event handler");
+ return NULL;
+ }
+
+ //-------------------------------------------------------
+ // Init the task and do settings
+ //-------------------------------------------------------
+ TString name("AliMCTruthTrackMaker");
+ name += outname;
+ AliMCTruthTrackMaker *eTask = new AliMCTruthTrackMaker(name);
+ eTask->SetFillMuMothers(mumothers);
+ eTask->SetTracksOutName(outname);
+ eTask->SetChargedMC(ch);
+ eTask->SetEtaMax(etamax);
+
+ //-------------------------------------------------------
+ // Final settings, pass to manager and set the containers
+ //-------------------------------------------------------
+ mgr->AddTask(eTask);
+
+ // Create containers for input/output
+ AliAnalysisDataContainer *cinput1 = mgr->GetCommonInputContainer();
+ mgr->ConnectInput (eTask, 0, cinput1 );
+
+ return eTask;
+}