SetDisplayInfo added
[u/mrichter/AliRoot.git] / HBTAN / AliHBTReaderKineTree.cxx
CommitLineData
81f7fe9f 1#include "AliHBTReaderKineTree.h"
2
d0c23b58 3#include <Riostream.h>
4//#include <Riostream.h>
81f7fe9f 5#include <TString.h>
6#include <TObjString.h>
7#include <TTree.h>
8#include <TFile.h>
9#include <TParticle.h>
10
11#include <AliRun.h>
12#include <AliStack.h>
13#include <AliHeader.h>
14
15#include "AliHBTRun.h"
16#include "AliHBTEvent.h"
17#include "AliHBTParticle.h"
18#include "AliHBTParticleCut.h"
19
20ClassImp(AliHBTReaderKineTree)
21/**********************************************************/
22
23AliHBTReaderKineTree::
24AliHBTReaderKineTree():fFileName("galice.root")
25{
26 fParticles = new AliHBTRun();
27 fIsRead = kFALSE;
28}
29
30AliHBTReaderKineTree::
31AliHBTReaderKineTree(TString& fname):fFileName(fname)
32{
33 fParticles = new AliHBTRun();
34 fIsRead = kFALSE;
35}
36
37
38/**********************************************************/
39AliHBTReaderKineTree::
40AliHBTReaderKineTree(TObjArray* dirs,const Char_t *filename):AliHBTReader(dirs),fFileName(filename)
41{
42 fParticles = new AliHBTRun();
43 fIsRead = kFALSE;
44}
45/**********************************************************/
46
47AliHBTEvent* AliHBTReaderKineTree::GetParticleEvent(Int_t n)
48 {
49 //returns Nth event with simulated particles
50 if (!fIsRead)
51 if(Read(fParticles,0x0))
52 {
53 Error("GetParticleEvent","Error in reading");
54 return 0x0;
55 }
56
57 return fParticles->GetEvent(n);
58 }
59
60/********************************************************************/
61
62Int_t AliHBTReaderKineTree::GetNumberOfPartEvents()
63 {
64 //returns number of events of particles
65 if (!fIsRead)
66 if(Read(fParticles,0x0))
67 {
68 Error("GetNumberOfPartEvents","Error in reading");
69 return 0;
70 }
71 return fParticles->GetNumberOfEvents();
72 }
73
74
75/**********************************************************/
76Int_t AliHBTReaderKineTree::
77Read(AliHBTRun* particles, AliHBTRun *tracks)
78 {
0b46db79 79 cout<<"AliHBTReaderKineTree::Read()"<<endl;
81f7fe9f 80 if (!particles) //check if an object is instatiated
81 {
82 Error("Read"," particles object must instatiated before passing it to the reader");
83 return 1;
84 }
85 particles->Reset();//clear runs == delete all old events
86
87 Int_t Ndirs;//total number of directories specified in fDirs array
88 Int_t Nevents; //number of events read in current directory
89 Int_t totalNevents = 0; //total number of read events
90 Int_t currentdir = 0; //number of current directory name is fDirs array
91
92 if (fDirs) //if array with directories is supplied by user
93 {
94 Ndirs = fDirs->GetEntries(); //get the number if directories
95 }
96 else
97 {
98 Ndirs = 0; //if the array is not supplied read only from current directory
99 }
100
101 do //do{}while; is OK even if 0 dirs specified. In that case we try to read from "./"
102 {
103 cout<<"________________________________________________________\n";
104 TFile* fileK = OpenFile(currentdir);
105 TTree* treeE = (TTree*)fileK->Get("TE");
106 if(treeE)
107 {
108 Nevents = (Int_t)treeE->GetEntries();
109 cout<<"Found "<<Nevents<<" in directory "<<GetDirName(currentdir)<<endl;
110 fileK->Delete("TE");
111 }
112 else
113 {
114 Error("Read","Cannot find Header Tree (TE)");
115 Nevents = 0;
116 }
117
118 for(Int_t currentEvent =0; currentEvent<Nevents;currentEvent++)//loop over all events
119 {
120// cout<<"processing event "<<currentEvent<<" in current dir."<<endl;
121 AliStack* stack = GetStack(currentEvent,fileK);
122 if (!stack)
123 {
124 Error("Read","Can not get stack for event %d",currentEvent);
125 continue;
126 }
127 Int_t npart = stack->GetNtrack();
128 Int_t nnn = 0;
129 for (Int_t i = 0;i<npart; i++)
130 {
131
132 TParticle * p = stack->Particle(i);
133 if (p->GetFirstMother() >= 0) continue;
134
135 if(Pass(p->GetPdgCode())) continue; //check if we are intersted with particles of this type
136 //if not take next partilce
137
138 AliHBTParticle* part = new AliHBTParticle(*p);
139 if(Pass(part)) { delete part; continue;}//check if meets all criteria of any of our cuts
140 //if it does not delete it and take next good track
141 particles->AddParticle(totalNevents,part);//put track and particle on the run
142
143 if ( (nnn%100) == 0)
144 {
145 cout<<nnn<<"\r";
146 fflush(0);
147 }
148 nnn++;
149 }
150 cout<<"Total read "<<nnn<<endl;
151 delete stack;
152 totalNevents++;
153 }
154 delete gAlice;
155 gAlice = 0;
156 if (fileK)
157 {
158 fileK->Close();
159 delete fileK;
160 fileK = 0;
161 }
162
163 currentdir++;
164 }while(currentdir < Ndirs);//end of loop over directories specified in fDirs Obj Array
165 fIsRead = kTRUE;
166 return 0;
167}
168
169/**********************************************************/
170AliStack* AliHBTReaderKineTree::GetStack(Int_t n, TFile* file)
171{
172 AliHeader *header = new AliHeader();
173 TTree* treeE = (TTree*)file->Get("TE");
174 if(treeE)
175 {
176 treeE->SetBranchAddress("Header", &header);
177 if (!treeE->GetEntry(n))
178 {
179 Error("GetEvent","Cannot find event:%dn",n);
180 delete header;
181 treeE->Delete("TE");
182 return 0x0;
183 }
184 }
185 else
186 {
187 Error("GetStack","Cannot find Header Tree (TE)n");
188 delete header;
189 return 0x0;
190 }
191
192 AliStack* stack = header->Stack();
193 if (stack)
194 {
195 if (!stack->GetEvent(n))
196 {
197 delete header;
198 treeE->Delete("TE");
199 return 0x0;
200 }
201 }
202 delete header;
203 treeE->Delete("TE");
204 return stack;
205}
206
207TFile* AliHBTReaderKineTree::OpenFile(Int_t n)
208{
209//opens file with kine tree
210
211 const TString& dirname = GetDirName(n);
212 if (dirname == "")
213 {
214 Error("OpenFiles","Can not get directory name");
215 return 0x0;
216 }
217 TString filename = dirname +"/"+ fFileName;
218 TFile *ret = TFile::Open(filename.Data());
219
220 if ( ret == 0x0)
221 {
222 Error("OpenFiles","Can't open file %s",filename.Data());
223 return 0x0;
224 }
225 if (!ret->IsOpen())
226 {
227 Error("OpenFiles","Can't open file %s",filename.Data());
228 return 0x0;
229 }
230
231 return ret;
232}