Retrofeeded recent developement in from AliHBTReaderESD (ITS track points, additional...
[u/mrichter/AliRoot.git] / HBTAN / AliHBTReaderITSv2.cxx
CommitLineData
1b446896 1#include "AliHBTReaderITSv2.h"
2
bed069a4 3#include <Varargs.h>
4
a9bfdd7b 5#include <TString.h>
6#include <TObjString.h>
7#include <TTree.h>
a9bfdd7b 8#include <TParticle.h>
9
10#include <AliRun.h>
88cb7938 11#include <AliRunLoader.h>
12#include <AliLoader.h>
bed069a4 13#include <AliStack.h>
a9bfdd7b 14#include <AliMagF.h>
15#include <AliITStrackV2.h>
a9bfdd7b 16
17#include "AliHBTRun.h"
18#include "AliHBTEvent.h"
19#include "AliHBTParticle.h"
20#include "AliHBTParticleCut.h"
21
22
1b446896 23ClassImp(AliHBTReaderITSv2)
24
bed069a4 25AliHBTReaderITSv2::AliHBTReaderITSv2():
26 fFileName("galice.root"),
27 fRunLoader(0x0),
28 fITSLoader(0x0),
29 fMagneticField(0.0),
30 fUseMagFFromRun(kTRUE)
88cb7938 31{
32 //constructor,
33 //Defaults:
34 // galicefilename = "galice.root"
88cb7938 35}
36/********************************************************************/
37
bed069a4 38AliHBTReaderITSv2::AliHBTReaderITSv2(const Char_t* galicefilename):
39 fFileName(galicefilename),
40 fRunLoader(0x0),
41 fITSLoader(0x0),
42 fMagneticField(0.0),
43 fUseMagFFromRun(kTRUE)
a9bfdd7b 44{
45 //constructor,
46 //Defaults:
a9bfdd7b 47 // galicefilename = "galice.root"
a9bfdd7b 48}
49/********************************************************************/
50
88cb7938 51AliHBTReaderITSv2::AliHBTReaderITSv2(TObjArray* dirs, const Char_t* galicefilename):
bed069a4 52 AliHBTReader(dirs),
53 fFileName(galicefilename),
54 fRunLoader(0x0),
55 fITSLoader(0x0),
56 fMagneticField(0.0),
57 fUseMagFFromRun(kTRUE)
a9bfdd7b 58{
59 //constructor,
60 //Defaults:
a9bfdd7b 61 // galicefilename = "galice.root"
62
a9bfdd7b 63}
64/********************************************************************/
bed069a4 65
66void AliHBTReaderITSv2::Rewind()
67{
68 //rewinds reading
69 delete fRunLoader;
70 fRunLoader = 0x0;
71 fCurrentDir = 0;
72 fNEventsRead= 0;
88cb7938 73}
a9bfdd7b 74/********************************************************************/
75
bed069a4 76AliHBTReaderITSv2::~AliHBTReaderITSv2()
77{
78 //dtor
79 delete fRunLoader;
80}
a9bfdd7b 81/********************************************************************/
88cb7938 82
bed069a4 83Int_t AliHBTReaderITSv2::ReadNext()
a9bfdd7b 84{
bed069a4 85//reads data from next event
86
88cb7938 87 register Int_t i = 0; //iterator
a9bfdd7b 88
41515b05 89// AliITStrackerV2 *tracker; // ITS tracker - used for cooking labels
bed069a4 90 TTree *tracktree = 0x0; // tree for tracks
91 AliITStrackV2 *iotrack = 0x0;
a9bfdd7b 92
93 Double_t xk;
94 Double_t par[5]; //Kalman track parameters
95 Float_t phi, lam, pt;//angles and transverse momentum
96 Int_t label; //label of the current track
a9bfdd7b 97
a9bfdd7b 98
bed069a4 99 if (fParticlesEvent == 0x0) fParticlesEvent = new AliHBTEvent();
100 if (fTracksEvent == 0x0) fTracksEvent = new AliHBTEvent();
101
102 fParticlesEvent->Reset();
103 fTracksEvent->Reset();
a9bfdd7b 104 do //do while is good even if Ndirs==0 (than read from current directory)
105 {
bed069a4 106 if (fRunLoader == 0x0)
107 if (OpenNextFile()) continue;//directory counter is increased inside in case of error
108
109 if (fCurrentEvent == fRunLoader->GetNumberOfEvents())
a9bfdd7b 110 {
bed069a4 111 //read next directory
112 delete fRunLoader;//close current session
113 fRunLoader = 0x0;//assure pointer is null
114 fCurrentDir++;//go to next dir
115 continue;//directory counter is increased inside in case of error
a9bfdd7b 116 }
bed069a4 117
118 Info("ReadNext","Reading Event %d",fCurrentEvent);
119
120 fRunLoader->GetEvent(fCurrentEvent);
121
122 tracktree=fITSLoader->TreeT();
123 if (!tracktree)
88cb7938 124 {
bed069a4 125 Error("ReadNext","Can't get a tree with ITS tracks");
126 fCurrentEvent++;
88cb7938 127 continue;
128 }
bed069a4 129
130 TBranch *tbranch=tracktree->GetBranch("tracks");
131 if (!tbranch)
88cb7938 132 {
bed069a4 133 Error("ReadNext","Can't get a branch with ITS tracks");
134 fCurrentEvent++;
88cb7938 135 continue;
136 }
bed069a4 137
138 AliStack* stack = fRunLoader->Stack();
139 if (stack == 0x0)
a9bfdd7b 140 {
bed069a4 141 Error("ReadNext","Can not get stack for current event",fCurrentEvent);
142 fCurrentEvent++;
5ee8b891 143 continue;
a9bfdd7b 144 }
bed069a4 145
146 //must be here because on the beginning conv. const. is not set yet
147 if (iotrack == 0x0) iotrack = new AliITStrackV2(); //buffer track for reading data from tree
a9bfdd7b 148
bed069a4 149 Int_t ntr = (Int_t)tracktree->GetEntries();
150
151 for (i=0; i < ntr; i++) //loop over all tpc tracks
a9bfdd7b 152 {
bed069a4 153 tbranch->SetAddress(&iotrack);
154 tracktree->GetEvent(i);
a9bfdd7b 155
bed069a4 156 label=iotrack->GetLabel();
157 if (label < 0)
158 {
159 continue;
160 }
88cb7938 161
bed069a4 162 TParticle *p = stack->Particle(label);
163 if(p == 0x0) continue; //if returned pointer is NULL
164 if(p->GetPDG() == 0x0) continue; //if particle has crezy PDG code (not known to our database)
a9bfdd7b 165
bed069a4 166 if(Pass(p->GetPdgCode())) continue; //check if we are intersted with particles of this type
167 //if not take next partilce
a9bfdd7b 168
bed069a4 169 AliHBTParticle* part = new AliHBTParticle(*p,i);
170 if(Pass(part)) { delete part; continue;}//check if meets all criteria of any of our cuts
171 //if it does not delete it and take next good track
a9bfdd7b 172
bed069a4 173 iotrack->PropagateTo(3.,0.0028,65.19);
174 iotrack->PropagateToVertex();
f1e2da22 175
bed069a4 176 iotrack->GetExternalParameters(xk,par); //get properties of the track
177 phi=TMath::ASin(par[2]) + iotrack->GetAlpha();
178 if (phi<-TMath::Pi()) phi+=2*TMath::Pi();
179 if (phi>=TMath::Pi()) phi-=2*TMath::Pi();
180 lam=par[3];
181 pt=1.0/TMath::Abs(par[4]);
a9bfdd7b 182
bed069a4 183 Double_t tpx = pt * TMath::Cos(phi); //track x coordinate of momentum
184 Double_t tpy = pt * TMath::Sin(phi); //track y coordinate of momentum
185 Double_t tpz = pt * lam; //track z coordinate of momentum
186
187 Double_t mass = p->GetMass();
188 Double_t tEtot = TMath::Sqrt( tpx*tpx + tpy*tpy + tpz*tpz + mass*mass);//total energy of the track
189
190 AliHBTParticle* track = new AliHBTParticle(p->GetPdgCode(), i, tpx, tpy , tpz, tEtot, 0., 0., 0., 0.);
191 if(Pass(track))//check if meets all criteria of any of our cuts
192 //if it does not delete it and take next good track
193 {
194 delete track;
195 delete part;
196 continue;
197 }
198
199 fParticlesEvent->AddParticle(part);
200 fTracksEvent->AddParticle(track);
201 }//end of loop over tracks in the event
a9bfdd7b 202
bed069a4 203 Info("ReadNext","Read %d tracks and %d particles from event %d (event %d in dir %d).",
204 fParticlesEvent->GetNumberOfParticles(), fTracksEvent->GetNumberOfParticles(),
205 fNEventsRead,fCurrentEvent,fCurrentDir);
a9bfdd7b 206
bed069a4 207 fCurrentEvent++;
208 fNEventsRead++;
209 delete iotrack;
210 return 0;
211 }while(fCurrentDir < GetNumberOfDirs());//end of loop over directories specified in fDirs Obj Array
a9bfdd7b 212
213 delete iotrack;
bed069a4 214 return 1;
a9bfdd7b 215}
216
217/********************************************************************/
bed069a4 218Int_t AliHBTReaderITSv2::OpenNextFile()
219{
220 //opens next file
221 TString filename = GetDirName(fCurrentDir);
222 if (filename.IsNull())
223 {
224 DoOpenError("Can not get directory name");
225 return 1;
226 }
227 filename = filename +"/"+ fFileName;
e191bb57 228 fRunLoader = AliRunLoader::Open(filename,AliConfig::GetDefaultEventFolderName());
bed069a4 229 if( fRunLoader == 0x0)
230 {
231 DoOpenError("Can not open session.");
232 return 1;
233 }
234
235 if (fRunLoader->GetNumberOfEvents() <= 0)
236 {
237 DoOpenError("There is no events in this directory.");
238 return 1;
239 }
240
241 if (fRunLoader->LoadKinematics())
242 {
243 DoOpenError("Error occured while loading kinematics.");
244 return 1;
245 }
246 fITSLoader = fRunLoader->GetLoader("ITSLoader");
247 if ( fITSLoader == 0x0)
248 {
249 DoOpenError("Exiting due to problems with opening files.");
250 return 1;
251 }
252
253 Info("OpenNextSession","________________________________________________________");
254 Info("OpenNextSession","Found %d event(s) in directory %s",
255 fRunLoader->GetNumberOfEvents(),GetDirName(fCurrentDir).Data());
256 Float_t mf;
257 if (fUseMagFFromRun)
258 {
259 if (fRunLoader->LoadgAlice())
260 {
261 DoOpenError("Error occured while loading AliRun.");
262 return 1;
263 }
264 mf = fRunLoader->GetAliRun()->Field()->SolenoidField();
265 Info("OpenNextSession","Setting Magnetic Field from run: B=%fT",mf/10.);
266 fRunLoader->UnloadgAlice();
267 }
268 else
269 {
270 Info("OpenNextSession","Setting Own Magnetic Field: B=%fT",fMagneticField);
271 if (fMagneticField == 0x0)
272 {
273 Fatal("OpenNextSession","Magnetic field can not be 0.");
274 return 1;//pro forma
275 }
276 mf = fMagneticField*10.;
277 }
278 AliKalmanTrack::SetConvConst(1000/0.299792458/mf);
279
280 if (fITSLoader->LoadTracks())
281 {
282 DoOpenError("Error occured while loading TPC tracks.");
283 return 1;
284 }
285
286 fCurrentEvent = 0;
287 return 0;
288}
a9bfdd7b 289/********************************************************************/
290
bed069a4 291void AliHBTReaderITSv2::DoOpenError( const char *va_(fmt), ...)
292{
293 // Does error display and clean-up in case error caught on Open Next Session
a9bfdd7b 294
bed069a4 295 va_list ap;
296 va_start(ap,va_(fmt));
297 Error("OpenNextFile", va_(fmt), ap);
298 va_end(ap);
299
300 delete fRunLoader;
301 fRunLoader = 0x0;
302 fITSLoader = 0x0;
303 fCurrentDir++;
304}
305
306/********************************************************************/