]>
Commit | Line | Data |
---|---|---|
1b446896 | 1 | #include "AliHBTReaderPPprod.h" |
2 | ||
3 | #include <iostream.h> | |
4 | #include <fstream.h> | |
5 | #include <TString.h> | |
6 | #include <TTree.h> | |
7 | #include <TFile.h> | |
8 | ||
9 | #include <AliTPCtrack.h> | |
10 | #include <AliTPCParam.h> | |
11 | #include <AliTPCtracker.h> | |
12 | ||
13 | #include "AliRun.h" | |
14 | #include "AliHBTRun.h" | |
15 | #include "AliHBTEvent.h" | |
16 | #include "AliHBTParticle.h" | |
17 | #include "AliHBTParticleCut.h" | |
18 | ||
19 | ||
20 | ClassImp(AliHBTReaderPPprod) | |
21 | ||
22 | AliHBTReaderPPprod:: | |
23 | AliHBTReaderPPprod(const Char_t* trackfilename,const Char_t* clusterfilename, | |
24 | const Char_t* goodtracksfilename,const Char_t* galicefilename): | |
25 | fTrackFileName(trackfilename),fClusterFileName(clusterfilename), | |
26 | fGAliceFileName(galicefilename), | |
27 | fGoodTPCTracksFileName(goodtracksfilename) | |
28 | { | |
29 | //constructor, only file names are set | |
30 | //Defaults: | |
31 | // trackfilename = "AliTPCtracks.root" | |
32 | // clusterfilename = "AliTPCclusters.root" | |
33 | // goodtracksfilename = "good_tracks_tpc" | |
34 | // galicefilename = "" - this means: Do not open gAlice file - | |
35 | // just leave the global pointer untached | |
36 | ||
37 | fParticles = new AliHBTRun(); | |
38 | fTracks = new AliHBTRun(); | |
39 | ||
40 | fTracksFile = 0x0; //files are opened during reading only | |
41 | fClustersFile = 0x0; | |
42 | ||
43 | fIsRead = kFALSE; | |
44 | } | |
45 | /********************************************************************/ | |
46 | ||
47 | AliHBTReaderPPprod::~AliHBTReaderPPprod() | |
48 | { | |
49 | delete fParticles; | |
50 | delete fTracks; | |
51 | } | |
52 | /********************************************************************/ | |
53 | ||
54 | AliHBTEvent* AliHBTReaderPPprod::GetParticleEvent(Int_t n) | |
55 | { | |
56 | //returns Nth event with simulated particles | |
57 | if (!fIsRead) Read(fParticles,fTracks); | |
58 | return fParticles->GetEvent(n); | |
59 | } | |
60 | /********************************************************************/ | |
61 | AliHBTEvent* AliHBTReaderPPprod::GetTrackEvent(Int_t n) | |
62 | { | |
63 | //returns Nth event with reconstructed tracks | |
64 | if (!fIsRead) Read(fParticles,fTracks); | |
65 | return fTracks->GetEvent(n); | |
66 | } | |
67 | /********************************************************************/ | |
68 | ||
69 | Int_t AliHBTReaderPPprod::GetNumberOfPartEvents() | |
70 | { | |
71 | //returns number of events of particles | |
72 | if (!fIsRead) Read(fParticles,fTracks); | |
73 | return fParticles->GetNumberOfEvents(); | |
74 | } | |
75 | ||
76 | /********************************************************************/ | |
77 | Int_t AliHBTReaderPPprod::GetNumberOfTrackEvents() | |
78 | { | |
79 | //returns number of events of tracks | |
80 | if (!fIsRead) Read(fParticles,fTracks); | |
81 | return fTracks->GetNumberOfEvents(); | |
82 | } | |
83 | /********************************************************************/ | |
84 | ||
85 | ||
86 | Int_t AliHBTReaderPPprod::Read(AliHBTRun* particles, AliHBTRun *tracks) | |
87 | { | |
88 | //reads data and puts put to the particles and tracks objects | |
89 | //reurns 0 if everything is OK | |
90 | // | |
91 | Int_t i; //iterator and some temprary values | |
92 | Int_t Nevents; | |
93 | if (!particles) //check if an object is instatiated | |
94 | { | |
95 | Error("AliHBTReaderPPprod::Read"," particles object must instatiated before passing it to the reader"); | |
96 | } | |
97 | if (!tracks) //check if an object is instatiated | |
98 | { | |
99 | Error("AliHBTReaderPPprod::Read"," tracks object must instatiated before passing it to the reader"); | |
100 | } | |
101 | particles->Reset();//clear runs == delete all old events | |
102 | tracks->Reset(); | |
103 | ||
104 | if( (i=OpenFiles()) ) | |
105 | { | |
106 | Error("AliHBTReaderPPprod::Read","Exiting due to problems with opening files. Errorcode %d",i); | |
107 | return i; | |
108 | } | |
109 | ||
110 | AliGoodTracksPP *goodTPCTracks = new AliGoodTracksPP(fGoodTPCTracksFileName); | |
111 | if (!goodTPCTracks) | |
112 | { | |
113 | Error("AliHBTReaderPPprod::Read","Exiting due to problems with opening files. Errorcode %d",i); | |
114 | return 1; | |
115 | } | |
116 | ||
117 | Nevents = 100; | |
118 | ||
119 | fClustersFile->cd(); | |
120 | AliTPCParam *TPCParam= (AliTPCParam*)fClustersFile->Get("75x40_100x60"); | |
121 | if (!TPCParam) | |
122 | { | |
123 | Error("AliHBTReaderPPprod::Read","TPC parameters have not been found !\n"); | |
124 | return 1; | |
125 | } | |
126 | ||
127 | TObjArray *tarray = new TObjArray(5000); | |
128 | tarray->SetOwner();// this causes memory leak, but in some cases deleting is infinite loop | |
129 | ||
130 | ||
131 | for(Int_t currentEvent =0; currentEvent<Nevents;currentEvent++) | |
132 | { | |
133 | cout<<"Reading Event "<<currentEvent<<endl; | |
134 | /**************************************/ | |
135 | /**************************************/ | |
136 | /**************************************/ | |
137 | fTracksFile->cd(); | |
138 | Char_t treename[100]; | |
139 | sprintf(treename,"TreeT_TPC_%d",currentEvent); | |
140 | ||
141 | TTree *tracktree=0; | |
142 | ||
143 | tracktree=(TTree*)fTracksFile->Get(treename); | |
144 | if (!tracktree) | |
145 | { | |
146 | Error("AliHBTReaderPPprod::Read","Can't get a tree with TPC tracks !\n"); | |
147 | return 1; | |
148 | } | |
149 | ||
150 | TBranch *trackbranch=tracktree->GetBranch("tracks"); | |
151 | if (!trackbranch) | |
152 | { | |
153 | Error("AliHBTReaderPPprod::Read","Can't get a branch with TPC tracks !\n"); | |
154 | return 2; | |
155 | } | |
156 | Int_t NTPCtracks=(Int_t)tracktree->GetEntries(); | |
157 | cout<<"Found "<<NTPCtracks<<" TPC tracks.\n"; | |
158 | //Copy tracks to array | |
159 | ||
160 | AliTPCtrack *iotrack=0; | |
161 | ||
162 | fClustersFile->cd(); | |
163 | AliTPCtracker *tracker = new AliTPCtracker(TPCParam,currentEvent); | |
164 | if (!tracker) | |
165 | { | |
166 | Error("AliHBTReaderPPprod::Read","Can't get a tracker !\n"); | |
167 | return 3; | |
168 | } | |
169 | tracker->LoadInnerSectors(); | |
170 | tracker->LoadOuterSectors(); | |
171 | ||
172 | for (i=0; i<NTPCtracks; i++) | |
173 | { | |
174 | iotrack=new AliTPCtrack; | |
175 | trackbranch->SetAddress(&iotrack); | |
176 | tracktree->GetEvent(i); | |
177 | tracker->CookLabel(iotrack,0.1); | |
178 | tarray->AddLast(iotrack); | |
179 | } | |
180 | ||
181 | ||
182 | fTracksFile->Delete(treename);//delete tree from memmory (and leave untached on disk)- we do not need it any more | |
183 | fTracksFile->Delete("tracks"); | |
184 | ||
185 | delete tracker; | |
186 | ||
187 | tracker = 0x0; | |
188 | trackbranch = 0x0; | |
189 | tracktree = 0x0; | |
190 | ||
191 | Int_t & ngood = goodTPCTracks->fGoodInEvent[currentEvent]; | |
192 | ||
193 | Double_t xk; | |
194 | Double_t par[5]; | |
195 | Float_t phi, lam, pt; | |
196 | Int_t label; | |
197 | Bool_t found; | |
198 | ||
199 | for (i=0; i<ngood; i++) | |
200 | { | |
201 | const struct GoodTrack & gt = goodTPCTracks->GetTrack(currentEvent,i); | |
202 | ||
203 | if(Pass(gt.code)) continue; | |
204 | ||
205 | label = gt.lab; | |
206 | found = kFALSE; //guard in case we don't find track with such a label | |
207 | for (Int_t j=0;j<NTPCtracks;j++) | |
208 | { | |
209 | iotrack = (AliTPCtrack*)tarray->At(j); | |
210 | if (iotrack->GetLabel() == label) | |
211 | { | |
212 | found = kTRUE; | |
213 | break; | |
214 | } | |
215 | } | |
216 | if(!found) | |
217 | { | |
218 | Warning("Read", | |
219 | "Sth is going wrong with tracks - there is no TPC track corresponding to goodtrack.\nGood tack label %d",label); | |
220 | continue; //put comunicate on the screen and continue loop | |
221 | } | |
222 | ||
223 | Double_t mass = TDatabasePDG::Instance()->GetParticle(gt.code)->Mass(); | |
224 | Double_t pEtot = TMath::Sqrt(gt.px*gt.px + gt.py*gt.py + gt.pz*gt.pz + mass*mass); | |
225 | ||
226 | AliHBTParticle* part = new AliHBTParticle(gt.code, gt.px, gt.py, gt.pz, pEtot, gt.x, gt.y, gt.z, 0.0); | |
227 | if(Pass(part)) continue; | |
228 | ||
229 | ||
230 | iotrack->PropagateTo(gt.x); | |
231 | iotrack->GetExternalParameters(xk,par); | |
232 | phi=TMath::ASin(par[2]) + iotrack->GetAlpha(); | |
233 | if (phi<-TMath::Pi()) phi+=2*TMath::Pi(); | |
234 | if (phi>=TMath::Pi()) phi-=2*TMath::Pi(); | |
235 | lam=par[3]; | |
236 | pt=1.0/TMath::Abs(par[4]); | |
237 | ||
238 | Double_t tpx = pt * TMath::Cos(phi); //track x coordinate of momentum | |
239 | Double_t tpy = pt * TMath::Sin(phi); //track x coordinate of momentum | |
240 | Double_t tpz = pt * lam; | |
241 | ||
242 | Double_t tEtot = TMath::Sqrt( tpx*tpx + tpy*tpy + tpz*tpz + mass*mass); | |
243 | ||
244 | AliHBTParticle* track = new AliHBTParticle(gt.code, tpx, tpy , tpz, tEtot, 0., 0., 0., 0.); | |
245 | if(Pass(track)) continue; | |
246 | ||
247 | particles->AddParticle(currentEvent,part); | |
248 | tracks->AddParticle(currentEvent,track); | |
249 | ||
250 | } | |
251 | tarray->Clear(); | |
252 | /**************************************/ | |
253 | /**************************************/ | |
254 | /**************************************/ | |
255 | } | |
256 | ||
257 | CloseFiles(); | |
258 | delete tarray; | |
259 | delete goodTPCTracks; | |
260 | fIsRead = kTRUE; | |
261 | return 0; | |
262 | } | |
263 | ||
264 | /********************************************************************/ | |
265 | Int_t AliHBTReaderPPprod::OpenFiles() | |
266 | { | |
267 | ||
268 | fTracksFile = 0; | |
269 | fTracksFile=TFile::Open("AliTPCtracks.root"); | |
270 | if (!fTracksFile->IsOpen()) | |
271 | { | |
272 | Error("AliHBTReaderPPprod::OpenFiles","Can't open AliTPCtracks.root"); | |
273 | return 1; | |
274 | } | |
275 | ||
276 | fClustersFile = 0; | |
277 | ||
278 | fClustersFile=TFile::Open("AliTPCclusters.root"); | |
279 | if (!fClustersFile->IsOpen()) | |
280 | { | |
281 | Error("AliHBTReaderPPprod::OpenFiles","Can't open AliTPCclusters.root"); | |
282 | return 2; | |
283 | } | |
284 | ||
285 | ||
286 | ||
287 | return 0; | |
288 | } | |
289 | ||
290 | ||
291 | ||
292 | /********************************************************************/ | |
293 | ||
294 | void AliHBTReaderPPprod::CloseFiles() | |
295 | { | |
296 | fTracksFile->Close(); | |
297 | fClustersFile->Close(); | |
298 | } | |
299 | ||
300 | /********************************************************************/ | |
301 | ||
302 | /********************************************************************/ | |
303 | /********************************************************************/ | |
304 | /********************************************************************/ | |
305 | ||
306 | ||
307 | AliGoodTracksPP::~AliGoodTracksPP() | |
308 | { | |
309 | delete [] fGoodInEvent; | |
310 | for (Int_t i = 0;i<fNevents;i++) | |
311 | delete [] fData[i]; | |
312 | delete [] fData; | |
313 | } | |
314 | /********************************************************************/ | |
315 | AliGoodTracksPP::AliGoodTracksPP(const TString& infilename) | |
316 | { | |
317 | ||
318 | fNevents = 100; | |
319 | ifstream in(infilename.Data()); | |
320 | ||
321 | if(!in) | |
322 | { | |
323 | cerr<<"Can not open file with Good TPC Tracks named:"<<infilename.Data()<<endl; | |
324 | delete this; | |
325 | return; | |
326 | } | |
327 | ||
328 | ||
329 | fGoodInEvent = new Int_t[fNevents]; | |
330 | fData = new struct GoodTrack* [fNevents]; | |
331 | ||
332 | Int_t i; | |
333 | for( i = 0;i<fNevents;i++) | |
334 | { | |
335 | fGoodInEvent[i] =0; | |
336 | fData[i] = new struct GoodTrack[500]; | |
337 | } | |
338 | Float_t tmp; | |
339 | Int_t evno; | |
340 | while(in>>evno) | |
341 | { | |
342 | if(fGoodInEvent[evno]>=500) | |
343 | { | |
344 | cerr<<"AliGoodTracksPP::AliGoodTracksPP() : Not enough place in the array\n"; | |
345 | continue; | |
346 | } | |
347 | in>>fData[evno][fGoodInEvent[evno]].lab; | |
348 | in>>fData[evno][fGoodInEvent[evno]].code; | |
349 | in>>fData[evno][fGoodInEvent[evno]].px; | |
350 | in>>fData[evno][fGoodInEvent[evno]].py; | |
351 | in>>fData[evno][fGoodInEvent[evno]].pz; | |
352 | in>>fData[evno][fGoodInEvent[evno]].x; | |
353 | in>>fData[evno][fGoodInEvent[evno]].y; | |
354 | in>>fData[evno][fGoodInEvent[evno]].z; | |
355 | ||
356 | in>>tmp; | |
357 | in>>tmp; | |
358 | in>>tmp; | |
359 | ||
360 | /* cout<<evno<<" "; | |
361 | cout<<fData[evno][fGoodInEvent[evno]].lab; | |
362 | cout<<" ";cout<<fData[evno][fGoodInEvent[evno]].code; | |
363 | cout<<" ";cout<<fData[evno][fGoodInEvent[evno]].px; | |
364 | cout<<" ";cout<<fData[evno][fGoodInEvent[evno]].py; | |
365 | cout<<" ";cout<<fData[evno][fGoodInEvent[evno]].pz; | |
366 | cout<<" ";cout<<fData[evno][fGoodInEvent[evno]].x; | |
367 | cout<<" ";cout<<fData[evno][fGoodInEvent[evno]].y; | |
368 | cout<<" ";cout<<fData[evno][fGoodInEvent[evno]].z; | |
369 | cout<<"\n"; | |
370 | */ | |
371 | fGoodInEvent[evno]++; | |
372 | } | |
373 | in.close(); | |
374 | cout<<"AliGoodTracksPP::AliGoodTracksPP() .... Done\n"; | |
375 | } | |
376 | ||
377 | ||
378 | ||
379 | const GoodTrack& AliGoodTracksPP::GetTrack(Int_t event, Int_t n) const | |
380 | { | |
381 | ||
382 | if( (event>fNevents) || (event<0)) | |
383 | { | |
384 | gAlice->Fatal("AliGoodTracksPP::GetTrack","No such Event %d",event); | |
385 | } | |
386 | if( (n>fGoodInEvent[event]) || (n<0)) | |
387 | { | |
388 | gAlice->Fatal("AliGoodTracksPP::GetTrack","No such Good TPC Track %d",n); | |
389 | } | |
390 | return fData[event][n]; | |
391 | ||
392 | } |