1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
19 /// \file DecodeRecoCocktail.C
20 /// \brief add brief description
22 /// \author A. De Falco, H. Woehri, INFN Cagliari, July 2006
24 /// This macro reads the generation/reconstruction files in the
25 /// input directories (recodir and simdir),
26 /// builds a tree that contains, for each event, an array of muons and dimuons
27 /// (TClonesArrays of AliMUONTrackLight and AliMUONPairLight objects)
28 /// and writes the tree in an output file (outFileName)
29 /// Note that if the path for the output file is not explicitly specified,
30 /// it will be written in the directory containing the generation/reconstruction
31 /// 27-Nov-2006: modified by in order to loop on files
34 /// Updated this macro to work with new version of AliMUONRecoCheck. Also, we are
35 /// now fetching reconstructed track data from ESD and not the track tree, because
36 /// the AliMUONTrack objects for reconstructed tracks are no longer stored on disk.
37 /// - Artur Szostak <artursz@iafrica.com>
39 #if !defined(__CINT__) || defined(__MAKECINT__)
40 #include <Riostream.h>
41 #include <TClonesArray.h>
45 #include <TLorentzVector.h>
46 #include <TParticle.h>
48 #include "AliMUONRecoCheck.h"
49 #include "AliMUONTrack.h"
50 #include "AliMUONTrackLight.h"
51 #include "AliMUONPairLight.h"
52 #include "AliMUONVTrackStore.h"
53 #include "AliESDEvent.h"
54 #include "AliESDVertex.h"
55 #include "AliESDMuonTrack.h"
56 #include "AliMCEventHandler.h"
57 #include "AliMCEvent.h"
59 /*TODO: need to update this with changes made to ITS
60 #include "AliITSVertexerPPZ.h"
61 #include "AliITSLoader.h"
66 void DecodeRecoCocktail(
67 char* recodir=".", // The directory containing galice.root for reconstructed data.
68 char* simdir="generated/", // The directory containing galice.root for simulated data.
69 char* outFileName = "MuonLight.root" // The output filename containing AliMUONTrackLight and AliMUONPairLight objects.
72 /// \param recodir The directory containing galice.root for reconstructed data.
73 /// \param simdir The directory containing galice.root for simulated data.
74 /// \param outFileName The output filename containing AliMUONTrackLight and AliMUONPairLight objects.
76 char startingDir[200];
77 sprintf (startingDir,"%s",gSystem->pwd());
80 TClonesArray *muonArray = new TClonesArray("AliMUONTrackLight",10);
81 TClonesArray *dimuonArray = new TClonesArray("AliMUONPairLight",10);
82 TTree *treeOut = new TTree("tree","tree");
84 TFile *fout = new TFile(outFileName,"recreate");
87 treeOut->Branch("muons",&muonArray);
88 treeOut->Branch("dimuons",&dimuonArray);
90 AliMUONRecoCheck *rc = new AliMUONRecoCheck("AliESDs.root", simdir);
92 /*TODO: need to update this with changes made to ITS
93 AliITSLoader* ITSloader = (AliITSLoader*) runLoaderSim->GetLoader("ITSLoader");
94 AliITSVertexerPPZ *dovert = 0;
96 dovert = new AliITSVertexerPPZ("default",0,0);
97 // dovert->SetDebug(0);
98 dovert->SetDiffPhiMax(0.05);
99 dovert->SetWindow(3.);
101 AliESDVertex *vert = 0;
106 Int_t nev = rc->NumberOfEvents();
107 for(Int_t ievent = 0; ievent < nev; ievent++){ // loop over events
109 /*TODO: need to update this with changes made to ITS
110 runLoaderSim->GetHeader();
112 vert = dovert->FindVertexForCurrentEvent(ievent);
116 //printf ("Event %d of %d\n",ievent,nev);
117 muonArray->Clear(); // clean muon and dimuon arrays
118 dimuonArray->Clear();
120 AliMUONVTrackStore* recoTracks = rc->ReconstructedTracks(ievent, kFALSE); // Use tracks from actual reconstruction.
121 //AliMUONVTrackStore* trackRefs = rc->ReconstructibleTracks(ievent); // Only use reconstructible reference tracks.
122 AliMUONVTrackStore* trackRefs = rc->TrackRefs(ievent);
123 AliStack* pstack = (const_cast<AliMCEventHandler*>(rc->GetMCEventHandler()))->MCEvent()->Stack();
125 // loop over ESD tracks
126 Int_t nreftracks = 0;
127 const AliESDEvent* esd = rc->GetESDEvent();
128 Int_t nTracks = (Int_t)esd->GetNumberOfMuonTracks() ;
129 for (Int_t iTrack = 0; iTrack < nTracks; iTrack++) {
131 AliESDMuonTrack* esdTrack = esd->GetMuonTrack(iTrack);
134 if (!esdTrack->ContainTrackerData()) continue;
136 // find the corresponding MUON track
137 AliMUONTrack* trackReco = (AliMUONTrack*) recoTracks->FindObject(esdTrack->GetUniqueID());
139 // try to match the reconstructed track with a simulated one
140 Int_t nMatchClusters = 0;
141 AliMUONTrack* matchedTrackRef = rc->FindCompatibleTrack(*trackReco, *trackRefs, nMatchClusters, kFALSE, 10.);
143 if (matchedTrackRef) {
145 //store new referenced track in the muonArray
146 AliMUONTrackLight* muLight = new ((*muonArray)[nreftracks++]) AliMUONTrackLight();
148 // assign parameters concerning the reconstructed tracks
149 muLight->FillFromESD(esdTrack);
150 // muLight->FillFromAliMUONTrack(trackReco);
152 // store further information related to the simulated track
153 muLight->SetTrackPythiaLine(matchedTrackRef->GetUniqueID());
154 TParticle *part = pstack->Particle(matchedTrackRef->GetUniqueID());
155 muLight->SetTrackPDGCode(part->GetPdgCode());
156 v.SetPxPyPzE(part->Px(), part->Py(), part->Pz(), part->Energy());
158 muLight->FillMuonHistory(pstack, part);
159 // muLight.PrintInfo("A");
165 // now loop over muon pairs to build dimuons
166 Int_t nmuons = muonArray->GetEntriesFast();
168 for(Int_t itrRec1 = 0; itrRec1 < nmuons-1; itrRec1++) {
169 AliMUONTrackLight* mu1 = (AliMUONTrackLight*) muonArray->At(itrRec1);
170 for(Int_t itrRec2 = itrRec1+1; itrRec2 < nmuons; itrRec2++){
171 AliMUONTrackLight* mu2 = (AliMUONTrackLight*) muonArray->At(itrRec2);
172 AliMUONPairLight dimuLight;
173 dimuLight.SetMuons(*mu1, *mu2);
174 TClonesArray &dimuons = *dimuonArray;
175 new (dimuons[ndimuons++]) AliMUONPairLight(dimuLight);
184 gSystem->cd(startingDir);