Retrofeeded recent developement in from AliHBTReaderESD (ITS track points, additional...
[u/mrichter/AliRoot.git] / HBTAN / AliHBTReaderITSv1.cxx
1 //____________________________________________________________________
2 //////////////////////////////////////////////////////////////////////
3 //                                                                  //
4 //  class AliHBTReaderITSv1                                         //
5 //                                                                  //
6 //  Reader for ITSv1 tracks. Not maintained since v1 is not         //
7 //  supposed to be used                                             //
8 //                                                                  //
9 //////////////////////////////////////////////////////////////////////
10
11 #include <Riostream.h>
12 #include <TBranch.h>
13 #include <TFile.h>
14 #include <TObjArray.h>
15 #include <TParticle.h>
16 #include <TString.h>
17 #include <TTree.h>
18
19 #include "AliHBTParticle.h"
20 #include "AliHBTParticleCut.h"
21 #include "AliHBTReaderITSv1.h"
22 #include "AliHBTRun.h"
23 #include "AliITSIOTrack.h"
24 #include "AliKalmanTrack.h"
25 #include "AliMagF.h"
26 #include "AliRun.h"
27 #include "AliStack.h"
28
29 ClassImp(AliHBTReaderITSv1)
30 /********************************************************************/
31
32 AliHBTReaderITSv1::AliHBTReaderITSv1(const Char_t* tracksfilename,const Char_t* galicefilename):
33  fITSTracksFileName(tracksfilename),
34  fGAliceFileName(galicefilename)
35 {
36  //ctor
37 }
38 /********************************************************************/
39
40 AliHBTReaderITSv1::AliHBTReaderITSv1(TObjArray* dirs, const Char_t* tracksfilename,const Char_t* galicefilename):
41  AliHBTReader(dirs),
42  fITSTracksFileName(tracksfilename),
43  fGAliceFileName(galicefilename)
44 {
45 //ctor
46 }
47 /********************************************************************/
48
49 AliHBTReaderITSv1::~AliHBTReaderITSv1()
50 {
51 //dtor
52 }
53 /********************************************************************/
54
55
56 Int_t AliHBTReaderITSv1::Read(AliHBTRun* particles, AliHBTRun *tracks)
57 {
58 //Reads data
59  Int_t nevents = 0;
60  AliITSIOTrack *iotrack=new AliITSIOTrack;
61  Int_t currentdir = 0;
62  Int_t ndirs;
63  Int_t totalnevents = 0;
64  
65  if (fDirs)
66   {
67     ndirs = fDirs->GetEntries();
68   }
69  else
70   {
71     ndirs = 0;
72   }
73  
74  do //do while is good even if 
75   {  
76    TFile* gAliceFile = OpenGAliceFile(currentdir);
77    if(gAliceFile == 0x0)
78     {
79        Error("Read","Can not open the file with gAlice");
80        delete iotrack;
81        return 1;
82     }
83    if (gAlice->TreeE())//check if tree E exists
84      {
85       nevents = (Int_t)gAlice->TreeE()->GetEntries();//if yes get number of events in gAlice
86       cout<<"________________________________________________________\n";
87       cout<<"Found "<<nevents<<" event(s) in directory "<<GetDirName(currentdir)<<endl;
88       cout<<"Setting Magnetic Field. Factor is "<<gAlice->Field()->Factor()<<endl;
89       AliKalmanTrack::SetConvConst(100/0.299792458/0.2/gAlice->Field()->Factor());
90      }
91     else
92      {//if not return an error
93        Error("Read","Can not find Header tree (TreeE) in gAlice");
94        delete iotrack;
95        return 4;
96      }
97
98    TFile *file = OpenTrackFile(currentdir);
99    if(file == 0x0)
100     {
101        Error("Read","Can not open the file with ITS tracks V1");
102        delete iotrack;
103        return 2;
104     }
105     
106    Int_t naccepted = 0;
107    char tname[30];
108    
109    for (Int_t currentEvent = 0; currentEvent < nevents; currentEvent++)
110     { 
111       cout<<"Reading Event "<<currentEvent;
112       
113       sprintf(tname,"TreeT%d",currentEvent);
114       file->cd(); 
115       TTree *tracktree=(TTree*)file->Get(tname);
116       TBranch *tbranch=tracktree->GetBranch("ITStracks");
117       tbranch->SetAddress(&iotrack);
118       
119       gAliceFile->cd();
120       gAlice->GetEvent(currentEvent);
121       gAlice->Stack()->Particles();
122
123       Int_t nentr=(Int_t)tracktree->GetEntries();
124       
125       cout<<".  Found "<<nentr<<" tracks.";
126       fflush(0);
127       
128       for (Int_t i=0; i<nentr; i++) 
129        {
130
131         tracktree->GetEvent(i);
132         if(!iotrack) continue;       
133         Int_t label = iotrack->GetLabel();
134         if (label < 0) 
135          {
136            continue;
137          }
138
139         TParticle *p = (TParticle*)gAlice->Stack()->Particle(label);
140         if(!p)
141          {
142            Warning("Read","Can not get particle with label &d",label);
143            continue;
144          }
145         if(Pass(p->GetPdgCode())) continue; //check if we are intersted with particles of this type
146                                            //if not take next partilce
147
148         AliHBTParticle* part = new AliHBTParticle(*p,i);
149         if(Pass(part)) { delete part; continue;}//check if meets all criteria of any of our cuts
150                                                 //if it does not delete it and take next good track
151         
152         Double_t px=iotrack->GetPx();
153         Double_t py=iotrack->GetPy();
154         Double_t pz=iotrack->GetPz();
155         Double_t mass = p->GetMass();
156         Double_t tEtot = TMath::Sqrt(px*px + py*py + pz*pz + mass*mass);//total energy of the track
157  
158         Double_t x= iotrack->GetX();
159         Double_t y= iotrack->GetY();
160         Double_t z= iotrack->GetZ();
161         
162         AliHBTParticle* track = new AliHBTParticle(p->GetPdgCode(), i, px, py , pz, tEtot, x, y, z, 0.);
163         if(Pass(track)) { delete  track;continue;}//check if meets all criteria of any of our cuts
164                                                   //if it does not delete it and take next good track
165
166         particles->AddParticle(totalnevents,part);//put track and particle on the run
167         tracks->AddParticle(totalnevents,track);
168         naccepted++;
169        }//end loop over tracks in the event
170
171        totalnevents++;
172        cout<<"  Accepted "<<naccepted<<" tracks"<<endl;
173      }//end of loop over events in current directory
174     
175     gAliceFile->Close();
176     delete gAliceFile;
177     gAliceFile = 0;
178     
179     file->Close(); 
180     delete file;
181     file = 0;
182     currentdir++;
183    }while(currentdir < ndirs);//end of loop over directories specified in fDirs Obj Array
184
185
186   delete iotrack;
187   fIsRead = kTRUE;
188   return 0;
189  
190  }
191 /********************************************************************/
192
193 TFile* AliHBTReaderITSv1::OpenTrackFile(Int_t ndir)
194 {
195 //opens files to be read for given directoru nomber in fDirs Array
196    const TString& dirname = GetDirName(ndir); 
197    if (dirname == "")
198     {
199       Error("OpenGAliceFile","Can not get directory name");
200       return 0x0;
201     }
202    TString filename = dirname + "/" + fITSTracksFileName;
203
204    TFile *file = TFile::Open(filename.Data());   
205    if (!file)
206     {
207       Error("Read","Can not open file %s",filename.Data());
208       return 0x0;
209     }
210    if (!file->IsOpen())
211     {
212       Error("Read","Can not open file %s",filename.Data());
213       return 0x0;
214     }
215    
216    return file;
217 }
218
219
220 /********************************************************************/
221 TFile* AliHBTReaderITSv1::OpenGAliceFile(Int_t ndir)
222 {
223 //Opens galice.root file
224   const TString& dirname = GetDirName(ndir); 
225    if (dirname == "")
226     {
227       Error("OpenGAliceFile","Can not get directory name");
228       return 0x0;
229     }
230   
231   TString filename = dirname + "/" + fGAliceFileName;
232
233   TFile* gAliceFile = TFile::Open(filename.Data());
234   if ( gAliceFile== 0x0)
235    {
236      Error("OpenFiles","Can't open file named %s",filename.Data());
237      return 0x0;
238    }
239   if (!gAliceFile->IsOpen())
240    {
241      Error("OpenFiles","Can't open file named %s",filename.Data());
242      return 0x0;
243    }
244
245   if (!(gAlice=(AliRun*)gAliceFile->Get("gAlice")))
246    {
247      Error("OpenFiles","gAlice have not been found on %s !\n",filename.Data());
248      gAliceFile->Close();
249      delete gAliceFile;
250      return 0x0;
251    }
252   
253   return gAliceFile;
254 }
255
256 /********************************************************************/
257 /********************************************************************/
258
259