Adding track references at decay points (M.Ivanov)
[u/mrichter/AliRoot.git] / STEER / AliTrackMapper.cxx
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 ////////////////////////////////////////////////////////////////////////
19 //
20 // AliTrackMapper.cxx
21 // description: 
22 //        create an AliTrackMap - a relation between track label and 
23 //        it's index in TreeH. Check all branches with hits.
24 //        Output is in the root file.
25 //
26 ////////////////////////////////////////////////////////////////////////
27
28 #include <Riostream.h>
29
30 #include "TTree.h"
31 #include "TFile.h"
32 #include "TStopwatch.h"
33 #include "TError.h"
34
35 #include "AliDetector.h"
36 #include "AliTrackMapper.h"
37 #include "AliTrackMap.h"
38 #include "AliRun.h"
39 #include "AliRunLoader.h"
40 #include "AliLoader.h"
41 #include "AliHit.h"
42
43 ClassImp(AliTrackMapper)
44
45 //_______________________________________________________________________
46 AliTrackMapper::AliTrackMapper(): 
47   fDEBUG(0)
48 {
49   //
50   // Default ctor
51   //
52 }
53
54
55 //_______________________________________________________________________
56 void AliTrackMapper::CreateMap(Int_t nEvents, Int_t firstEventNr,
57                                const char* fnMap, const char* fnHits)
58 {
59   //
60   // method to create a track map for a given number of events starting 
61   // from event firstEventNr. This method just opens and closes files and
62   // loops over events, the creation of map is delegated to 
63   // CreateMap(Int_t eventNr, TFile* fileMap) method
64   //
65   TStopwatch timer;
66   timer.Start();
67   
68   TFile *fileMap=TFile::Open(fnMap,"new");
69   if (!fileMap->IsOpen()) {cerr<<"Can't open output file "<<fnMap<<"!\n"; return;}
70
71   AliRunLoader* rl = AliRunLoader::Open(fnHits);
72   if (rl == 0x0) {cerr<<"Can't open input file "<<fnHits<<"!\n"; return;}
73   if (rl->LoadgAlice())
74     {
75       ::Error("AliTrackMapper::CreateMap","Error occured while loading AliRun");
76       return;
77     }
78     
79   if (!(gAlice=rl->GetAliRun())) {
80     cerr<<"gAlice have not been found on session !\n";
81     return;
82   }
83   
84   rl->LoadKinematics();
85   
86   for (Int_t eventNr = firstEventNr; eventNr < firstEventNr+nEvents;
87        eventNr++) {
88     CreateMap(eventNr,fileMap,rl);
89   } // end loop over events
90   
91   delete rl;
92   fileMap->Close();
93   delete fileMap;
94   timer.Stop();
95   if (fDEBUG > 0) timer.Print();
96 }
97
98 //_______________________________________________________________________
99 Int_t  AliTrackMapper::CreateMap(Int_t eventNr, TFile* fileMap,AliRunLoader* rl) 
100 {
101   //
102   // create an AliTrackMap for a given event
103   // correct gAlice must be already present in memory
104   //
105
106   rl->GetEvent(eventNr);
107
108   TTree *treeK = rl->TreeK();
109   if (!treeK) {
110     cerr<<"Error: Event "<<eventNr<<", treeK not found."<<endl;
111     return -1;
112   }
113   Int_t nAllParticles = static_cast<Int_t>(treeK->GetEntries());
114   Int_t *trackMap = new Int_t[nAllParticles];
115   for (Int_t i = 0; i<nAllParticles; i++) {trackMap[i] = kNoEntry;}
116
117
118   TObjArray *modules = gAlice->Detectors();
119   if (!modules) {
120     cerr<<"TObjArray with modules not found."<<endl;
121     return -1;
122   }
123   Int_t nModules = static_cast<Int_t>(modules->GetEntries());
124   AliHit* hit;
125   for (Int_t iModule = 0; iModule < nModules; iModule++) 
126    {
127     AliDetector * detector = dynamic_cast<AliDetector*>(modules->At(iModule));
128     if (!detector) continue;
129     AliLoader* loader = detector->GetLoader();
130     if (loader == 0x0)
131      {
132        ::Warning("AliTrackMapper::CreateMap",
133                  "Can not get loader from detector %s.",detector->GetName());
134        continue;
135      }
136     Int_t retval = loader->LoadHits();
137     if (retval) {
138       ::Error("AliTrackMapper::CreateMap",
139             "Event %d: error occured while loading hits for %s",
140             eventNr,detector->GetName());
141       return -1;
142      }
143     
144     TTree *treeH = loader->TreeH();
145     if (!treeH) {
146       ::Error("AliTrackMapper::CreateMap","Event %d: Can not get TreeH for %s",
147              eventNr,detector->GetName());
148       return -1;
149      }
150     Int_t treeHEntries = static_cast<Int_t>(treeH->GetEntries());
151     if (fDEBUG > 1) cout<<"treeHEntries "<<treeHEntries<<endl;
152      
153     detector->ResetHits();
154     
155     for (Int_t treeHIndex = 0; treeHIndex < treeHEntries; treeHIndex++)  
156      { // process only detectors with shunt = 0
157       treeH->GetEvent(treeHIndex);
158       if (detector->GetIshunt()) continue; 
159
160       hit=dynamic_cast<AliHit*>(detector->FirstHit(-1));
161       Int_t lastLabel=-1, label;
162       for( ; hit; hit=dynamic_cast<AliHit*>(detector->NextHit()) ) {
163        label=hit->Track();       
164        if (lastLabel != label) {
165          if (label < 0 || label >= nAllParticles) {
166            cerr<<"Error: label out of range. ";
167            cerr<<"Event "<<eventNr<<" treeHIndex "<<treeHIndex<<" label = "<<label<<endl;
168            return -2;
169          }
170          if (trackMap[label] >=0 && trackMap[label] != treeHIndex) {
171            cerr<<"Error: different treeHIndex for label "<<label
172               <<" indeces: "<<trackMap[label]<<" != "<<treeHIndex;
173            cerr<<" event "<<eventNr<<" detector "<<detector->GetName()<<endl;
174            return -3;
175          }
176          trackMap[label] = treeHIndex;
177          if (fDEBUG > 2) cout<<"treeHIndex, label = "<<treeHIndex<<" "<<label<<endl;
178          lastLabel = label;
179        }
180       }
181     }//loop over hits in module
182     loader->UnloadHits();
183   }//loop over modules
184
185   if (fDEBUG > 2) {
186     for (Int_t i = 0; i < nAllParticles; i++) {
187       cout<<eventNr<<"\t"<<i<<"\t"<<trackMap[i]<<endl;
188     }
189   }
190   fileMap->cd();
191   AliTrackMap* trackMapObject = new AliTrackMap(nAllParticles, trackMap);
192   trackMapObject->SetEventNr(eventNr);
193   trackMapObject->Write();
194   delete trackMapObject;
195
196   delete [] trackMap;
197   return 0;
198 }
199   
200 //_______________________________________________________________________
201 AliTrackMap* AliTrackMapper::LoadTrackMap(Int_t eventNr, const char* fnMap, TFile* &fileMap) {
202   //
203   // read an AliTrackMap object for the given event eventNr from 
204   // the file fileMap
205   //
206   fileMap=TFile::Open(fnMap);
207   if (!fileMap->IsOpen()) {cerr<<"Can't open file "<<fnMap<<" with map!\n"; return 0;}
208   char mapName[20];
209   sprintf(mapName,"AliTrackMap_%5.5d",eventNr);
210   AliTrackMap* trackMapObject = dynamic_cast<AliTrackMap*>(fileMap->Get(mapName));
211   if (!trackMapObject) {
212     cerr<<"Error: map named "<<mapName<<" not found."<<endl;
213     return 0;
214   }
215   return trackMapObject;
216 }
217
218 //_______________________________________________________________________
219 void AliTrackMapper::CheckTrackMap(Int_t eventNr, const char* fnMap) {
220   //
221   // 
222   //
223   TFile *fileMap;
224   AliTrackMap* trackMapObject = LoadTrackMap(eventNr, fnMap, fileMap);
225   if (!trackMapObject) return;
226   
227   trackMapObject->PrintValues();
228   
229   delete trackMapObject;
230   fileMap->Close();
231   delete fileMap;
232 }
233