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 **************************************************************************/
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
28 #include <TClonesArray.h>
32 #include <TLorentzVector.h>
33 #include <TParticle.h>
35 #include <AliMUONRecoCheck.h>
36 #include <AliMUONTrack.h>
37 #include <AliMUONTrackParam.h>
38 #include <AliRunLoader.h>
39 #include "AliMUONTrackLight.h"
40 #include "AliMUONPairLight.h"
42 void DecodeRecoCocktail(char* dirname=".", char* outFileName = "MuonLight.root"){
43 const char *startingDir = gSystem->pwd();
45 TFile *fout = new TFile(outFileName,"recreate");
47 AliMUONRecoCheck *rc = new AliMUONRecoCheck("galice.root");
48 AliRunLoader *runLoader = rc->GetRunLoader();
50 TClonesArray *muonArray = new TClonesArray("AliMUONTrackLight",10);
51 TClonesArray *dimuonArray = new TClonesArray("AliMUONPairLight",10);
52 TTree *treeOut = new TTree("tree","tree");
55 treeOut->Branch("muons",&muonArray);
56 treeOut->Branch("dimuons",&dimuonArray);
58 runLoader->LoadKinematics("READ");;
59 Int_t nev = runLoader->GetNumberOfEvents();
62 for(Int_t ievent = 0; ievent < nev; ievent++){ // loop over events
63 printf ("Event %d of %d\n",ievent,nev);
64 muonArray->Clear(); // clean muon and dimuon arrays
66 runLoader->GetEvent(ievent);
68 rc->MakeTrackRef(); // make reconstructable tracks
70 TClonesArray * trackRecoArray = rc->GetTrackReco();
71 TClonesArray * trackRefArray = rc->GetMuonTrackRef();
72 Int_t nTrackReco = trackRecoArray->GetEntriesFast();
74 for (Int_t itrRec = 0; itrRec<nTrackReco; itrRec++) {
75 // assign parameters concerning the reconstructed tracks
76 AliMUONTrack *trackReco = (AliMUONTrack *)trackRecoArray->At(itrRec);
77 AliMUONTrackLight *muLight = new AliMUONTrackLight();
78 AliMUONTrackParam *trPar = trackReco->GetTrackParamAtVertex();
79 muLight->SetCharge(Int_t(TMath::Sign(1.,trPar->GetInverseBendingMomentum())));
80 muLight->SetPxPyPz(trPar->Px(),trPar->Py(), trPar->Pz());
81 muLight->SetTriggered(trackReco->GetMatchTrigger());
82 Double_t xyz[3] = { trPar->GetNonBendingCoor(),
83 trPar->GetBendingCoor(),
85 muLight->SetVertex(xyz);
86 // find the reference track and store further information
87 TParticle *part = muLight->FindRefTrack(trackReco,trackRefArray,runLoader);
89 v.SetPxPyPzE(part->Px(), part->Py(), part->Pz(), part->Energy());
91 muLight->FillMuonHistory(runLoader, part);
92 muLight->PrintInfo("A");
93 //store the referenced track in the muonArray:
94 TClonesArray &muons = *muonArray;
95 new (muons[nreftracks++]) AliMUONTrackLight(*muLight);
99 // now loop over muon pairs to build dimuons
100 Int_t nmuons = muonArray->GetEntriesFast();
102 for(Int_t itrRec1 = 0; itrRec1 < nmuons-1; itrRec1++) {
103 AliMUONTrackLight* mu1 = (AliMUONTrackLight*) muonArray->At(itrRec1);
104 for(Int_t itrRec2 = itrRec1+1; itrRec2 < nmuons; itrRec2++){
105 AliMUONTrackLight* mu2 = (AliMUONTrackLight*) muonArray->At(itrRec2);
106 AliMUONPairLight *dimuLight = new AliMUONPairLight();
107 dimuLight->SetMuons(*mu1, *mu2);
108 // if(dimuLight->GetCreationProcess() == 2){
109 // printf ("#dimuon = %d (%d, %d) \n", ndimuons, itrRec1, itrRec2);
110 // dimuLight->PrintInfo("A");
112 //store the referenced track in the dimuonArray:
113 TClonesArray &dimuons = *dimuonArray;
114 new (dimuons[ndimuons++]) AliMUONPairLight(*dimuLight);
117 Int_t ndimu2 = dimuonArray->GetEntriesFast();
118 // printf ("dimuonArray has %d entries\n",ndimu2);
119 // dimuonArray->Dump();
120 // printf ("filling tree\n");
123 // printf ("done\n");
128 gSystem->cd(startingDir);