- Adding a new group for macros
[u/mrichter/AliRoot.git] / MUON / DecodeRecoCocktail.C
CommitLineData
fa6e7866 1/**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3 * *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
6 * *
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 **************************************************************************/
15
16/* $Id$ */
17
18// A. De Falco, H. Woehri, INFN Cagliari, July 2006
19// This macro reads the generation/reconstruction files in the
9bdda5f6 20// input directories (recodir and simdir),
fa6e7866 21// builds a tree that contains, for each event, an array of muons and dimuons
22// (TClonesArrays of AliMUONTrackLight and AliMUONPairLight objects)
23// and writes the tree in an output file (outFileName)
24// Note that if the path for the output file is not explicitly specified,
25// it will be written in the directory containing the generation/reconstruction
4d309390 26// 27-Nov-2006: modified by in order to loop on files
fa6e7866 27
9bdda5f6 28// 13 Nov 2007:
29// Updated this macro to work with new version of AliMUONRecoCheck. Also, we are
30// now fetching reconstructed track data from ESD and not the track tree, because
31// the AliMUONTrack objects for reconstructed tracks are no longer stored on disk.
32// - Artur Szostak <artursz@iafrica.com>
33
4d309390 34#if !defined(__CINT__) || defined(__MAKECINT__)
35#include <Riostream.h>
fa6e7866 36#include <TClonesArray.h>
37#include <TFile.h>
38#include <TTree.h>
39#include <TMath.h>
40#include <TLorentzVector.h>
41#include <TParticle.h>
42#include <TSystem.h>
59830a87 43#include <TGeoManager.h>
4d309390 44#include "AliMUONRecoCheck.h"
45#include "AliMUONTrack.h"
46#include "AliMUONTrackParam.h"
fa6e7866 47#include "AliMUONTrackLight.h"
48#include "AliMUONPairLight.h"
9bdda5f6 49#include "AliMUONVTrackStore.h"
22ccc301 50#include "AliMUONTrackExtrap.h"
9bdda5f6 51#include "AliESDEvent.h"
52#include "AliESDVertex.h"
53#include "AliMCEventHandler.h"
54#include "AliMCEvent.h"
55/*TODO: need to update this with changes made to ITS
4d309390 56#include "AliITSVertexerPPZ.h"
57#include "AliITSLoader.h"
9bdda5f6 58*/
59830a87 59#include "AliTracker.h"
60#include "AliMagFMaps.h"
4d309390 61#endif
fa6e7866 62
4d309390 63
9bdda5f6 64void DecodeRecoCocktail(
65 char* recodir=".", // The directory containing galice.root for reconstructed data.
66 char* simdir="generated/", // The directory containing galice.root for simulated data.
67 char* outFileName = "MuonLight.root", // The output filename containing AliMUONTrackLight and AliMUONPairLight objects.
61fed964 68 char* geoFilename = "geometry.root" // The filename containing the geometry.
9bdda5f6 69 )
70{
4d309390 71 char startingDir[200];
72 sprintf (startingDir,"%s",gSystem->pwd());
9bdda5f6 73 gSystem->cd(recodir);
4d309390 74
75 TClonesArray *muonArray = new TClonesArray("AliMUONTrackLight",10);
76 TClonesArray *dimuonArray = new TClonesArray("AliMUONPairLight",10);
77 TTree *treeOut = new TTree("tree","tree");
78
fa6e7866 79 TFile *fout = new TFile(outFileName,"recreate");
4d309390 80 fout->cd();
fa6e7866 81
4d309390 82 treeOut->Branch("muons",&muonArray);
83 treeOut->Branch("dimuons",&dimuonArray);
84
59830a87 85 // Import TGeo geometry (needed by AliMUONTrackExtrap::ExtrapToVertex)
86 if (!gGeoManager) {
87 TGeoManager::Import(geoFilename);
88 if (!gGeoManager) {
89 Error("MUONmass_ESD", "getting geometry from file %s failed", geoFilename);
90 return;
91 }
92 }
93
94 // set mag field
95 // waiting for mag field in CDB
96 printf("Loading field map...\n");
97 AliMagFMaps* field = new AliMagFMaps("Maps","Maps", 1, 1., 10., AliMagFMaps::k5kG);
98 AliTracker::SetFieldMap(field, kFALSE);
99 // set the magnetic field for track extrapolations
100 AliMUONTrackExtrap::SetField(AliTracker::GetFieldMap());
101
61fed964 102 AliMUONRecoCheck *rc = new AliMUONRecoCheck("AliESDs.root", simdir);
9bdda5f6 103 Int_t nev = rc->NumberOfEvents();
fa6e7866 104
9bdda5f6 105 /*TODO: need to update this with changes made to ITS
a4ee7ab9 106 AliITSLoader* ITSloader = (AliITSLoader*) runLoaderSim->GetLoader("ITSLoader");
4d309390 107 AliITSVertexerPPZ *dovert = 0;
108 if (ITSloader) {
109 dovert = new AliITSVertexerPPZ("default",0,0);
21939432 110 // dovert->SetDebug(0);
4d309390 111 dovert->SetDiffPhiMax(0.05);
112 dovert->SetWindow(3.);
113 }
114 AliESDVertex *vert = 0;
9bdda5f6 115 */
fa6e7866 116
117 TLorentzVector v;
9bdda5f6 118
fa6e7866 119 for(Int_t ievent = 0; ievent < nev; ievent++){ // loop over events
9bdda5f6 120
121 /*TODO: need to update this with changes made to ITS
a4ee7ab9 122 runLoaderSim->GetHeader();
4d309390 123 if (ITSloader) {
124 vert = dovert->FindVertexForCurrentEvent(ievent);
125 }
9bdda5f6 126 */
127
a4ee7ab9 128 //printf ("Event %d of %d\n",ievent,nev);
fa6e7866 129 muonArray->Clear(); // clean muon and dimuon arrays
130 dimuonArray->Clear();
fa6e7866 131
61fed964 132 AliMUONVTrackStore* recoTracks = rc->ReconstructedTracks(ievent); // Use tracks from actual reconstruction, but this does not work anymore.
9bdda5f6 133 //AliMUONVTrackStore* trackRefs = rc->ReconstructibleTracks(ievent); // Only use reconstructible reference tracks.
134 AliMUONVTrackStore* trackRefs = rc->TrackRefs(ievent);
61fed964 135 AliStack* pstack = (const_cast<AliMCEventHandler*>(rc->GetMCEventHandler()))->MCEvent()->Stack();
9bdda5f6 136
137 TIter next(recoTracks->CreateIterator());
138 AliMUONTrack* trackReco = NULL;
139
140 Int_t nTrackReco = recoTracks->GetSize();
61fed964 141 Int_t nTracksESD = rc->GetESDEvent()->GetNumberOfMuonTracks();
9bdda5f6 142 if (nTrackReco != nTracksESD) printf ("Tracks in recoTracks (%d) and in ESD (%d) do not match!\n", nTrackReco, nTracksESD);
143 Int_t nreftracks = 0;
144 Int_t itrRec = 0;
145 while ( (trackReco = static_cast<AliMUONTrack*>(next())) != NULL )
146 {
fa6e7866 147 // assign parameters concerning the reconstructed tracks
4d309390 148 AliMUONTrackLight muLight;
a4ee7ab9 149
61fed964 150 muLight.FillFromESD(rc->GetESDEvent()->GetMuonTrack(itrRec));
9bdda5f6 151 // muLight.FillFromAliMUONTrack(trackReco);
4d309390 152
153 // find the reference track and store further information
61fed964 154 TParticle *part = muLight.FindRefTrack(trackReco, trackRefs, pstack);
fa6e7866 155 if (part) {
156 v.SetPxPyPzE(part->Px(), part->Py(), part->Pz(), part->Energy());
82724a77 157 muLight.SetPGen(v);
9bdda5f6 158 muLight.FillMuonHistory(pstack, part);
4d309390 159 // muLight.PrintInfo("A");
fa6e7866 160 //store the referenced track in the muonArray:
161 TClonesArray &muons = *muonArray;
82724a77 162 new (muons[nreftracks++]) AliMUONTrackLight(muLight);
4d309390 163 }
9bdda5f6 164
165 itrRec++;
4d309390 166 } // end reco track
fa6e7866 167
168 // now loop over muon pairs to build dimuons
169 Int_t nmuons = muonArray->GetEntriesFast();
170 Int_t ndimuons = 0;
171 for(Int_t itrRec1 = 0; itrRec1 < nmuons-1; itrRec1++) {
172 AliMUONTrackLight* mu1 = (AliMUONTrackLight*) muonArray->At(itrRec1);
173 for(Int_t itrRec2 = itrRec1+1; itrRec2 < nmuons; itrRec2++){
174 AliMUONTrackLight* mu2 = (AliMUONTrackLight*) muonArray->At(itrRec2);
4d309390 175 AliMUONPairLight dimuLight;
82724a77 176 dimuLight.SetMuons(*mu1, *mu2);
fa6e7866 177 TClonesArray &dimuons = *dimuonArray;
82724a77 178 new (dimuons[ndimuons++]) AliMUONPairLight(dimuLight);
fa6e7866 179 }
180 }
9bdda5f6 181
fa6e7866 182 treeOut->Fill();
4d309390 183 }
184
fa6e7866 185 fout->cd();
186 treeOut->Write();
187 gSystem->cd(startingDir);
82724a77 188 fout->Close();
82724a77 189 delete muonArray;
190 delete dimuonArray;
4d309390 191 delete treeOut;
192 delete fout;
193 delete rc;
fa6e7866 194}
4d309390 195