1 #include "AliHBTReaderPPprod.h"
9 #include <AliTPCtrack.h>
10 #include <AliTPCParam.h>
11 #include <AliTPCtracker.h>
14 #include "AliHBTRun.h"
15 #include "AliHBTEvent.h"
16 #include "AliHBTParticle.h"
17 #include "AliHBTParticleCut.h"
20 ClassImp(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)
29 //constructor, only file names are set
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
37 fParticles = new AliHBTRun();
38 fTracks = new AliHBTRun();
40 fTracksFile = 0x0; //files are opened during reading only
45 /********************************************************************/
47 AliHBTReaderPPprod::~AliHBTReaderPPprod()
52 /********************************************************************/
54 AliHBTEvent* AliHBTReaderPPprod::GetParticleEvent(Int_t n)
56 //returns Nth event with simulated particles
57 if (!fIsRead) Read(fParticles,fTracks);
58 return fParticles->GetEvent(n);
60 /********************************************************************/
61 AliHBTEvent* AliHBTReaderPPprod::GetTrackEvent(Int_t n)
63 //returns Nth event with reconstructed tracks
64 if (!fIsRead) Read(fParticles,fTracks);
65 return fTracks->GetEvent(n);
67 /********************************************************************/
69 Int_t AliHBTReaderPPprod::GetNumberOfPartEvents()
71 //returns number of events of particles
72 if (!fIsRead) Read(fParticles,fTracks);
73 return fParticles->GetNumberOfEvents();
76 /********************************************************************/
77 Int_t AliHBTReaderPPprod::GetNumberOfTrackEvents()
79 //returns number of events of tracks
80 if (!fIsRead) Read(fParticles,fTracks);
81 return fTracks->GetNumberOfEvents();
83 /********************************************************************/
86 Int_t AliHBTReaderPPprod::Read(AliHBTRun* particles, AliHBTRun *tracks)
88 //reads data and puts put to the particles and tracks objects
89 //reurns 0 if everything is OK
91 cout<<"New I/O not implemented\n";
94 Int_t i; //iterator and some temprary values
96 if (!particles) //check if an object is instatiated
98 Error("AliHBTReaderPPprod::Read"," particles object must instatiated before passing it to the reader");
100 if (!tracks) //check if an object is instatiated
102 Error("AliHBTReaderPPprod::Read"," tracks object must instatiated before passing it to the reader");
104 particles->Reset();//clear runs == delete all old events
107 if( (i=OpenFiles()) )
109 Error("AliHBTReaderPPprod::Read","Exiting due to problems with opening files. Errorcode %d",i);
113 AliGoodTracksPP *goodTPCTracks = new AliGoodTracksPP(fGoodTPCTracksFileName);
116 Error("AliHBTReaderPPprod::Read","Exiting due to problems with opening files. Errorcode %d",i);
123 AliTPCParam *TPCParam= (AliTPCParam*)fClustersFile->Get("75x40_100x60");
126 Error("AliHBTReaderPPprod::Read","TPC parameters have not been found !\n");
130 TObjArray *tarray = new TObjArray(5000);
131 tarray->SetOwner();// this causes memory leak, but in some cases deleting is infinite loop
134 for(Int_t currentEvent =0; currentEvent<Nevents;currentEvent++)
136 cout<<"Reading Event "<<currentEvent<<endl;
137 /**************************************/
138 /**************************************/
139 /**************************************/
141 Char_t treename[100];
142 sprintf(treename,"TreeT_TPC_%d",currentEvent);
146 tracktree=(TTree*)fTracksFile->Get(treename);
149 Error("AliHBTReaderPPprod::Read","Can't get a tree with TPC tracks !\n");
153 TBranch *trackbranch=tracktree->GetBranch("tracks");
156 Error("AliHBTReaderPPprod::Read","Can't get a branch with TPC tracks !\n");
159 Int_t NTPCtracks=(Int_t)tracktree->GetEntries();
160 cout<<"Found "<<NTPCtracks<<" TPC tracks.\n";
161 //Copy tracks to array
163 AliTPCtrack *iotrack=0;
166 printf("This method is not converted to the NewIO !\n"); //I.B.
168 AliTPCtracker *tracker = new AliTPCtracker(TPCParam); //I.B.
169 //AliTPCtracker *tracker = new AliTPCtracker(TPCParam,currentEvent,"");
172 Error("AliHBTReaderPPprod::Read","Can't get a tracker !\n");
175 tracker->LoadClusters(0);//I.Belikov, "0" should be a pointer to a tree
177 for (i=0; i<NTPCtracks; i++)
179 iotrack=new AliTPCtrack;
180 trackbranch->SetAddress(&iotrack);
181 tracktree->GetEvent(i);
182 tracker->CookLabel(iotrack,0.1);
183 tarray->AddLast(iotrack);
187 fTracksFile->Delete(treename);//delete tree from memmory (and leave untached on disk)- we do not need it any more
188 fTracksFile->Delete("tracks");
196 Int_t & ngood = goodTPCTracks->fGoodInEvent[currentEvent];
200 Float_t phi, lam, pt;
204 for (i=0; i<ngood; i++)
206 const struct GoodTrack & gt = goodTPCTracks->GetTrack(currentEvent,i);
208 if(Pass(gt.code)) continue;
211 found = kFALSE; //guard in case we don't find track with such a label
212 for (Int_t j=0;j<NTPCtracks;j++)
214 iotrack = (AliTPCtrack*)tarray->At(j);
215 if (iotrack->GetLabel() == label)
224 "Sth is going wrong with tracks - there is no TPC track corresponding to goodtrack.\nGood tack label %d",label);
225 continue; //put comunicate on the screen and continue loop
228 Double_t mass = TDatabasePDG::Instance()->GetParticle(gt.code)->Mass();
229 Double_t pEtot = TMath::Sqrt(gt.px*gt.px + gt.py*gt.py + gt.pz*gt.pz + mass*mass);
231 AliHBTParticle* part = new AliHBTParticle(gt.code,i, gt.px, gt.py, gt.pz, pEtot, gt.x, gt.y, gt.z, 0.0);
232 if(Pass(part)) continue;
235 iotrack->PropagateTo(gt.x);
236 iotrack->GetExternalParameters(xk,par);
237 phi=TMath::ASin(par[2]) + iotrack->GetAlpha();
238 if (phi<-TMath::Pi()) phi+=2*TMath::Pi();
239 if (phi>=TMath::Pi()) phi-=2*TMath::Pi();
241 pt=1.0/TMath::Abs(par[4]);
243 Double_t tpx = pt * TMath::Cos(phi); //track x coordinate of momentum
244 Double_t tpy = pt * TMath::Sin(phi); //track x coordinate of momentum
245 Double_t tpz = pt * lam;
247 Double_t tEtot = TMath::Sqrt( tpx*tpx + tpy*tpy + tpz*tpz + mass*mass);
249 AliHBTParticle* track = new AliHBTParticle(gt.code,i, tpx, tpy , tpz, tEtot, 0., 0., 0., 0.);
250 if(Pass(track)) continue;
252 particles->AddParticle(currentEvent,part);
253 tracks->AddParticle(currentEvent,track);
257 /**************************************/
258 /**************************************/
259 /**************************************/
264 delete goodTPCTracks;
269 /********************************************************************/
270 Int_t AliHBTReaderPPprod::OpenFiles()
274 fTracksFile=TFile::Open("AliTPCtracks.root");
275 if (!fTracksFile->IsOpen())
277 Error("AliHBTReaderPPprod::OpenFiles","Can't open AliTPCtracks.root");
283 fClustersFile=TFile::Open("AliTPCclusters.root");
284 if (!fClustersFile->IsOpen())
286 Error("AliHBTReaderPPprod::OpenFiles","Can't open AliTPCclusters.root");
297 /********************************************************************/
299 void AliHBTReaderPPprod::CloseFiles()
301 fTracksFile->Close();
302 fClustersFile->Close();
305 /********************************************************************/
307 /********************************************************************/
308 /********************************************************************/
309 /********************************************************************/
312 AliGoodTracksPP::~AliGoodTracksPP()
314 delete [] fGoodInEvent;
315 for (Int_t i = 0;i<fNevents;i++)
319 /********************************************************************/
320 AliGoodTracksPP::AliGoodTracksPP(const TString& infilename)
324 ifstream in(infilename.Data());
328 cerr<<"Can not open file with Good TPC Tracks named:"<<infilename.Data()<<endl;
334 fGoodInEvent = new Int_t[fNevents];
335 fData = new struct GoodTrack* [fNevents];
338 for( i = 0;i<fNevents;i++)
341 fData[i] = new struct GoodTrack[500];
347 if(fGoodInEvent[evno]>=500)
349 cerr<<"AliGoodTracksPP::AliGoodTracksPP() : Not enough place in the array\n";
352 in>>fData[evno][fGoodInEvent[evno]].lab;
353 in>>fData[evno][fGoodInEvent[evno]].code;
354 in>>fData[evno][fGoodInEvent[evno]].px;
355 in>>fData[evno][fGoodInEvent[evno]].py;
356 in>>fData[evno][fGoodInEvent[evno]].pz;
357 in>>fData[evno][fGoodInEvent[evno]].x;
358 in>>fData[evno][fGoodInEvent[evno]].y;
359 in>>fData[evno][fGoodInEvent[evno]].z;
366 cout<<fData[evno][fGoodInEvent[evno]].lab;
367 cout<<" ";cout<<fData[evno][fGoodInEvent[evno]].code;
368 cout<<" ";cout<<fData[evno][fGoodInEvent[evno]].px;
369 cout<<" ";cout<<fData[evno][fGoodInEvent[evno]].py;
370 cout<<" ";cout<<fData[evno][fGoodInEvent[evno]].pz;
371 cout<<" ";cout<<fData[evno][fGoodInEvent[evno]].x;
372 cout<<" ";cout<<fData[evno][fGoodInEvent[evno]].y;
373 cout<<" ";cout<<fData[evno][fGoodInEvent[evno]].z;
376 fGoodInEvent[evno]++;
379 cout<<"AliGoodTracksPP::AliGoodTracksPP() .... Done\n";
384 const GoodTrack& AliGoodTracksPP::GetTrack(Int_t event, Int_t n) const
387 if( (event>fNevents) || (event<0))
389 gAlice->Fatal("AliGoodTracksPP::GetTrack","No such Event %d",event);
391 if( (n>fGoodInEvent[event]) || (n<0))
393 gAlice->Fatal("AliGoodTracksPP::GetTrack","No such Good TPC Track %d",n);
395 return fData[event][n];