5 #include "TClonesArray.h"
6 #include "TLorentzVector.h"
11 #include "AliAnalysisTask.h"
12 #include "AliAnalysisManager.h"
14 #include "AliESDEvent.h"
15 #include "AliESDInputHandler.h"
16 #include "AliESDMuonTrack.h"
18 #include "AliMCEventHandler.h"
19 #include "AliMCEvent.h"
21 #include "AliAnalysisTaskRecoCheck.h"
23 #include "AliMUONTrackLight.h"
24 #include "AliMUONPairLight.h"
25 #include "AliMUONVTrackStore.h"
26 #include "AliMUONTrack.h"
27 #include "AliMUONESDInterface.h"
28 #include "AliMUONRecoCheck.h"
30 // analysis task for decoding reconstructed tracks and kinematics (AliMUONRecoCheck)
31 // Authors: Bogdan Vulpescu
33 ClassImp(AliAnalysisTaskRecoCheck)
35 //________________________________________________________________________
36 AliAnalysisTaskRecoCheck::AliAnalysisTaskRecoCheck(const char *name)
37 : AliAnalysisTask(name, ""), fESDEvent(0), fTree(0), fArray1Mu(0), fArray2Mu(0), fL3Current(30000.0)
41 // Define input and output slots here
42 // Input slot #0 works with a TChain
43 DefineInput(0, TChain::Class());
44 // Output slot #0 writes into a TTree container
45 DefineOutput(0, TTree::Class());
49 //________________________________________________________________________
50 void AliAnalysisTaskRecoCheck::ConnectInputData(Option_t *)
52 // Connect ESD or AOD here
55 TTree* tree = dynamic_cast<TTree*> (GetInputData(0));
57 Printf("ERROR: Could not read chain from input slot 0");
59 tree->SetBranchStatus("*", kFALSE);
60 tree->SetBranchStatus("MuonTracks.*", kTRUE);
61 tree->SetBranchStatus("AliESDHeader.*", kTRUE);
63 AliESDInputHandler *esdH = dynamic_cast<AliESDInputHandler*> (AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
66 Printf("ERROR: Could not get ESDInputHandler");
68 fESDEvent = esdH->GetEvent();
73 //________________________________________________________________________
74 void AliAnalysisTaskRecoCheck::CreateOutputObjects()
76 // Create output tree with branches of arrays
79 fArray1Mu = new TClonesArray("AliMUONTrackLight",0);
80 fArray2Mu = new TClonesArray("AliMUONPairLight",0);
81 fTree = new TTree("tree","tree",0);
82 fTree->Branch("muons", &fArray1Mu);
83 fTree->Branch("dimuons",&fArray2Mu);
84 fTree->SetDirectory(0);
88 //________________________________________________________________________
89 void AliAnalysisTaskRecoCheck::Exec(Option_t *)
92 // Called for each event
94 // Process MC truth, therefore we receive the AliAnalysisManager and ask it for the AliMCEventHandler
95 // This handler can return the current MC event
97 Int_t nTracks = fESDEvent->GetNumberOfMuonTracks();
98 if (nTracks == 0) return;
100 AliMCEventHandler* eventHandler = dynamic_cast<AliMCEventHandler*> (AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler());
102 Printf("ERROR: Could not retrieve MC event handler");
106 AliMCEvent* mcEvent = eventHandler->MCEvent();
108 Printf("ERROR: Could not retrieve MC event");
117 AliMUONRecoCheck *rc = new AliMUONRecoCheck(fESDEvent, eventHandler);
119 AliMUONVTrackStore* trackRefs = rc->TrackRefs(-1);
120 AliStack* pstack = mcEvent->Stack();
122 // loop over ESD tracks
123 Int_t nreftracks = 0;
124 for (Int_t iTrack = 0; iTrack < nTracks; iTrack++) {
126 AliESDMuonTrack* esdTrack = fESDEvent->GetMuonTrack(iTrack);
129 if (!esdTrack->ContainTrackerData()) continue;
131 // convert ESD track to MUON track (without recomputing track parameters at each clusters)
132 AliMUONTrack muonTrack;
133 AliMUONESDInterface::ESDToMUON(*esdTrack, muonTrack, kFALSE);
135 // try to match the reconstructed track with a simulated one
136 Int_t nMatchClusters = 0;
137 AliMUONTrack* matchedTrackRef = rc->FindCompatibleTrack(muonTrack, *trackRefs, nMatchClusters, kFALSE, 10.);
139 if (matchedTrackRef) {
141 //store new referenced track in the muonArray with parameters from the reconstructed tracks
142 AliMUONTrackLight* muLight = new ((*fArray1Mu)[nreftracks++]) AliMUONTrackLight(esdTrack);
144 // store further information related to the simulated track
145 muLight->SetTrackPythiaLine(matchedTrackRef->GetUniqueID());
146 TParticle *part = pstack->Particle(matchedTrackRef->GetUniqueID());
147 muLight->SetTrackPDGCode(part->GetPdgCode());
148 v.SetPxPyPzE(part->Px(), part->Py(), part->Pz(), part->Energy());
150 muLight->FillMuonHistory(pstack, part);
156 // now loop over muon pairs to build dimuons
157 Int_t nmuons = fArray1Mu->GetEntriesFast();
159 for(Int_t itrRec1 = 0; itrRec1 < nmuons-1; itrRec1++) {
160 AliMUONTrackLight* mu1 = (AliMUONTrackLight*) fArray1Mu->At(itrRec1);
161 for(Int_t itrRec2 = itrRec1+1; itrRec2 < nmuons; itrRec2++){
162 AliMUONTrackLight* mu2 = (AliMUONTrackLight*) fArray1Mu->At(itrRec2);
163 AliMUONPairLight dimuLight;
164 dimuLight.SetMuons(*mu1, *mu2);
165 //dimuLight.PrintInfo("A");
166 TClonesArray &dimuons = *fArray2Mu;
167 new (dimuons[ndimuons++]) AliMUONPairLight(dimuLight);
174 if (nreftracks != 0) {
183 //________________________________________________________________________
184 void AliAnalysisTaskRecoCheck::Terminate(Option_t *)