- ESD can be used instead of AliMUONTrack objects to access the reconstructed variables.
[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
20// input directory (dirname),
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
26
27#include <iostream>
28#include <TClonesArray.h>
29#include <TFile.h>
30#include <TTree.h>
31#include <TMath.h>
32#include <TLorentzVector.h>
33#include <TParticle.h>
34#include <TSystem.h>
59830a87 35#include <TGeoManager.h>
fa6e7866 36#include <AliMUONRecoCheck.h>
37#include <AliMUONTrack.h>
38#include <AliMUONTrackParam.h>
39#include <AliRunLoader.h>
40#include "AliMUONTrackLight.h"
41#include "AliMUONPairLight.h"
22ccc301 42#include "AliMUONTrackExtrap.h"
59830a87 43#include "AliTracker.h"
44#include "AliMagFMaps.h"
fa6e7866 45
59830a87 46void DecodeRecoCocktail(char* dirname=".", char* outFileName = "MuonLight.root", char* geoFilename = "geometry.root"){
fa6e7866 47 const char *startingDir = gSystem->pwd();
48 gSystem->cd(dirname);
49 TFile *fout = new TFile(outFileName,"recreate");
50
59830a87 51 // Import TGeo geometry (needed by AliMUONTrackExtrap::ExtrapToVertex)
52 if (!gGeoManager) {
53 TGeoManager::Import(geoFilename);
54 if (!gGeoManager) {
55 Error("MUONmass_ESD", "getting geometry from file %s failed", geoFilename);
56 return;
57 }
58 }
59
60 // set mag field
61 // waiting for mag field in CDB
62 printf("Loading field map...\n");
63 AliMagFMaps* field = new AliMagFMaps("Maps","Maps", 1, 1., 10., AliMagFMaps::k5kG);
64 AliTracker::SetFieldMap(field, kFALSE);
65 // set the magnetic field for track extrapolations
66 AliMUONTrackExtrap::SetField(AliTracker::GetFieldMap());
67
fa6e7866 68 AliMUONRecoCheck *rc = new AliMUONRecoCheck("galice.root");
69 AliRunLoader *runLoader = rc->GetRunLoader();
70
71 TClonesArray *muonArray = new TClonesArray("AliMUONTrackLight",10);
72 TClonesArray *dimuonArray = new TClonesArray("AliMUONPairLight",10);
73 TTree *treeOut = new TTree("tree","tree");
74
75 fout->cd();
76 treeOut->Branch("muons",&muonArray);
77 treeOut->Branch("dimuons",&dimuonArray);
78
79 runLoader->LoadKinematics("READ");;
80 Int_t nev = runLoader->GetNumberOfEvents();
81
82 TLorentzVector v;
83 for(Int_t ievent = 0; ievent < nev; ievent++){ // loop over events
84 printf ("Event %d of %d\n",ievent,nev);
85 muonArray->Clear(); // clean muon and dimuon arrays
86 dimuonArray->Clear();
87 runLoader->GetEvent(ievent);
88 rc->ResetTracks();
89 rc->MakeTrackRef(); // make reconstructable tracks
90
91 TClonesArray * trackRecoArray = rc->GetTrackReco();
92 TClonesArray * trackRefArray = rc->GetMuonTrackRef();
93 Int_t nTrackReco = trackRecoArray->GetEntriesFast();
94 Int_t nreftracks = 0;
95 for (Int_t itrRec = 0; itrRec<nTrackReco; itrRec++) {
96 // assign parameters concerning the reconstructed tracks
97 AliMUONTrack *trackReco = (AliMUONTrack *)trackRecoArray->At(itrRec);
82724a77 98 AliMUONTrackLight muLight;
22ccc301 99 AliMUONTrackParam trPar(*((AliMUONTrackParam*) (trackReco->GetTrackParamAtHit()->First())));
100 AliMUONTrackExtrap::ExtrapToVertex(&trPar,0.,0.,0.);
101 muLight.SetCharge(Int_t(TMath::Sign(1.,trPar.GetInverseBendingMomentum())));
102 muLight.SetPxPyPz(trPar.Px(),trPar.Py(), trPar.Pz());
82724a77 103 muLight.SetTriggered(trackReco->GetMatchTrigger());
22ccc301 104 Double_t xyz[3] = { trPar.GetNonBendingCoor(),
105 trPar.GetBendingCoor(),
106 trPar.GetZ()};
82724a77 107 muLight.SetVertex(xyz);
fa6e7866 108 // find the reference track and store further information
82724a77 109 TParticle *part = muLight.FindRefTrack(trackReco,trackRefArray,runLoader);
fa6e7866 110 if (part) {
111 v.SetPxPyPzE(part->Px(), part->Py(), part->Pz(), part->Energy());
82724a77 112 muLight.SetPGen(v);
113 muLight.FillMuonHistory(runLoader, part);
114 muLight.PrintInfo("A");
fa6e7866 115 //store the referenced track in the muonArray:
116 TClonesArray &muons = *muonArray;
82724a77 117 new (muons[nreftracks++]) AliMUONTrackLight(muLight);
fa6e7866 118 }
119 }
120
121 // now loop over muon pairs to build dimuons
122 Int_t nmuons = muonArray->GetEntriesFast();
123 Int_t ndimuons = 0;
124 for(Int_t itrRec1 = 0; itrRec1 < nmuons-1; itrRec1++) {
125 AliMUONTrackLight* mu1 = (AliMUONTrackLight*) muonArray->At(itrRec1);
126 for(Int_t itrRec2 = itrRec1+1; itrRec2 < nmuons; itrRec2++){
127 AliMUONTrackLight* mu2 = (AliMUONTrackLight*) muonArray->At(itrRec2);
82724a77 128 AliMUONPairLight dimuLight;
129 dimuLight.SetMuons(*mu1, *mu2);
130 // if(dimuLight.GetCreationProcess() == 2){
fa6e7866 131 // printf ("#dimuon = %d (%d, %d) \n", ndimuons, itrRec1, itrRec2);
82724a77 132 // dimuLight.PrintInfo("A");
fa6e7866 133 // }
134 //store the referenced track in the dimuonArray:
135 TClonesArray &dimuons = *dimuonArray;
82724a77 136 new (dimuons[ndimuons++]) AliMUONPairLight(dimuLight);
fa6e7866 137 }
138 }
82724a77 139 //Int_t ndimu2 = dimuonArray->GetEntriesFast();
fa6e7866 140 // printf ("dimuonArray has %d entries\n",ndimu2);
141 // dimuonArray->Dump();
142 // printf ("filling tree\n");
143 // // fill the tree
144 treeOut->Fill();
145 // printf ("done\n");
146
147 }
148 fout->cd();
149 treeOut->Write();
150 gSystem->cd(startingDir);
82724a77 151 fout->Close();
152 delete fout;
153 delete rc;
154 delete muonArray;
155 delete dimuonArray;
fa6e7866 156}